Fork me on GitHub
image frame

My Blog

第二章

第一题

　　只需要把$x_1x_2$换成y在利用0-1规划的性质对条件进行变形即可。

第二题

　　很明显决策变量是小区，要达到全部覆盖只需要满足每个小区能选的小学数量大于1即可。
　　这题连数据都不需要。

1
2
3
4
5
6
7
8
9
10
11
12
13
model:
sets:
col/1..6/:x;
endsets
min=@sum(col:x);
x(1)+x(2)+x(3)>1;
x(2)+x(4)>1;
x(3)+x(5)>1;
x(4)+x(6)>1;
x(5)+x(6)>1;
x(1)>1;
x(2)+x(4)+x(6)>1;
end

第三题

　　没啥好说的引入0-1的限制就行。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
model:
sets:
she/1..6/;
com/1..4/;
link(she,com):c,x;
endsets
data:
c=4 2 3 4
6 4 5 5
7 6 7 6
7 8 8 6
7 9 8 6
7 10 8 6;
enddata
max=@sum(link:c*x);
@for(she(i):@sum(com(j):x(i,j))=1);
@for(com(j):@sum(she(i):x(i,j))>1);
@for(link:@bin(x));
end

第四题

　　待更新

第五题

　　首先考虑到比较小的数据范围，可以搜索出料头小于1米的切割方案（大于1米则可以切出一根1米的棍子）。
　　对于这七种切割方案进行规划，即可得到最小使用根数（注意材料根数要是整数）。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
model:
sets:
row/1..3/;
col/1..7/:x;
link(row,col):a;
endsets
data:
a=1 2 0 0 0 0 1
0 0 3 2 1 0 1
4 1 0 2 4 6 1;
enddata
min=@sum(col:x);
@for(row(i):@sum(col(j):a(i,j)*x(j))>100);
@for(col:@gin(x));
end

第六题

　　和第三题很像

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
model:
sets:
car/1..8/;
pos/1..5/;
link(pos,car):c,x;
endsets
data:
c=30 25 18 32 27 19 22 26
29 31 19 18 21 20 30 19
28 29 30 19 19 22 23 26
29 30 19 24 25 19 18 21
21 20 18 17 16 14 16 18;
enddata
min=@sum(link:c*x);
@for(pos(i):@sum(car(j):x(i,j))=1);
@for(car(j):@sum(pos(i):x(i,j))<1);
@for(link:@bin(x));
end

第一章

第一题

　　没啥好说的，直接用linprog干就完事。只要注意把maxz取负求最小值即可。大于等于的条件也要记得加负号。

1
2
3
4
5
6
7
8
clc
clear
c=[3 -1 -1]';
a=[1 -2 1;4 -1 -2];
b=[11;-3];
aeq=[-2 0 1];beq=1;
[x,fval]=linprog(-c,a,b,aeq,beq,zeros(3,1));
fval=-fval

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
model:
sets:
col/1..3/:c,x;
row/1..2/:b;
link(row,col):a;
endsets
data:
c=3 -1 1;
a=1 -2 1 4 -1 -2;
b=11 -3;
enddata
max=@sum(col:c*x);
@for(row(i):@sum(col(j):a(i,j)*x(j))<b(i));
-2*x(1)+x(3)=1;
end

第二题

　　只要注意把$\vert x \vert$变换为$u_i+v_i$，把$x$变换为$u_i-v_i$进行建模即可。

1
2
3
4
5
6
7
8
9
10
clc
clear
c=1:4;
c=[c,c]';
aeq=[1 -1 -1 1;1 -2 1 -3;1 -1 -2 3];
beq=[0;1;-0.5];
aeq=[aeq,-aeq];
[uv,fval]=linprog(c,[],[],aeq,beq,zeros(8,1));
x=uv(1:4)-uv(5:8)
fval

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
model:
sets:
col/1..4/:c,x;
row/1..3/:b;
links(row,col):a;
endsets
data:
c=1 2 3 4;
a=1 -1 -1 1 1 -1 1 -3 1 -1 -2 3;
b=0 1 -0.5;
enddata
min=@sum(col:c*@abs(x));
@for(row(i):@sum(col(j):a(i,j)*x(j))=b(i));
@for(col:@free(x)); !x的值可正可负
end

第三题

　　这里我是用的非线性的优化函数fmincon，但是用线性优化函数再加上整数限制算可能更好。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
clc
clear
f=@(x) -(x(1)+x(2)+1.65*x(8)+2.3*x(9)-(300/6000)*(5*x(1)+10*x(6))-(321/10000)*(7*x(2)+9*x(7)+12*x(9))-(250/4000)*(6*x(3)+8*x(8))-(783/7000)*(4*x(4)+11*x(9))-(200/4000)*x(5));
a=zeros(5,9);
a(1,1)=5;a(1,6)=10;
a(2,2)=7;a(2,7)=9;a(2,9)=12;
a(3,3)=6;a(3,8)=8;
a(4,4)=4;a(4,9)=11;
a(5,5)=7;
b=[6000;10000;4000;7000;4000];
aeq=[1 1 -1 -1 -1 0 0 0 0;0 0 0 0 0 1 1 -1 0];
beq=[0;0];
[x,fval]=fmincon(f,randn(9,1),a,b,aeq,beq,zeros(9,1));
fval=-fval

第四题

　　本题使用matlab求解很不方便，这里只提供lingo代码。
　　这道题关键在于对三个仓载货量和最大容许量成比例的建模。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
model:
sets:
wu/1..4/:a,b,c,y;
cang/1..3/:w,v;
links(wu,cang):x;
endsets
data:
c=3100 3800 3500 2850;
b=480 650 580 390;
a=18 15 23 12;
w=10 16 8;
v=6800 8700 5300;
enddata
max=@sum(wu(i):c(i)*@sum(cang(j):x(i,j)));
@for(wu(i):@sum(cang(j):x(i,j))<a(i));
@for(cang(j):@sum(wu(i):x(i,j))<w(j));
@for(cang(j):@sum(wu(i):b(i)*x(i,j))<v(j));
@for(cang(j):@sum(wu(i):x(i,j))/w(j)=(@sum(wu(i):x(i,1))/w(1)));
@for(wu(i):y(i)=@sum(cang(j):x(i,j)));
end

!a(i)为运输货物i的重量
!b(i)为单位货物所占的空间
!c(i)为单位货物的利润
!y(i)为四种物资的量
!w,v分别是空间和体积限制

第五题

　　这题和上一题一样，由于决策变量是二维变量使用lingo求解更为方便

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
model:
sets:
row/1..3/;
col/1..4/:b;
link(row,col):c,x;
endsets
data:
b=62 83 39 91;
c=132 0 100 103
  0 91 100 100
  106 89 100 98;
@text('ex.txt')=@table(x); !存到文本文件中;
enddata
min=@sum(link:c*x);
@for(col(j):@sum(row(i):x(i,j))=b(j));
end

97. 约数之和 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
typedef long long ll;
const int mod =9901; 
ll a,b,ans;
ll qpow(ll x,ll k)
{
    ll res=1;
    x%=mod;
    while(k)
    {
        if(k&1)
            res=(res*x)%mod;
        x=(x*x)%mod;
        k>>=1;
    }
    return res;
}
ll sum(ll p,ll k)
{
    if(k==0)
        return 1;
    if(k%2==0)
        return (p%mod*sum(p,k-1)+1)%mod;
    else
        return (1+qpow(p,k/2+1))*sum(p,k/2)%mod;
}
int main()
{
    cin >> a >>b;
    ans=1;
    for(int i=2;i<=a;i++)
    {
        int cnt=0;
        while(a%i==0)
        {
            cnt++;
            a/=i;
        }
        if(cnt)
            ans=ans*sum(i,cnt*b)%mod;
        //cout<<a<<endl;
    }
    if(!a)
        ans=0;
    cout << ans<<endl;
    return 0;
}

243. 一个简单的整数问题2 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
const int maxn =1e5+1;
typedef long long ll;
ll d[4*maxn],b[4*maxn],a[maxn],n,m,l,r,k;
string op;
void build(ll s,ll t,ll p)
{
    ll mid=(t+s)>>1;
    if(s==t)
    {
        d[p]=a[s];
        return ;
    }
    build(s,mid,p*2);
    build(mid+1,t,p*2+1);
    d[p]=d[p*2]+d[p*2+1];
}
void pushdown(ll s,ll t,ll p)
{
    ll m=(s+t)>>1;
    if(b[p])
    {
        d[p*2]+=(m-s+1)*b[p];
        d[p*2+1]+=(t-m)*b[p];
        b[p*2+1]+=b[p];b[p*2]+=b[p];
        b[p]=0;
    }
}
void add(ll l,ll r,ll c,ll s,ll t,ll p)
{
    ll mid=(s+t)>>1;
    if(l<=s&&t<=r)
    {
        d[p]+=(t-s+1)*c;
        b[p]+=c;
        return ;
    }
    pushdown(s,t,p);
    if(l<=mid)
        add(l,r,c,s,mid,p*2);
    if(r>mid)
        add(l,r,c,mid+1,t,p*2+1);
    d[p]=d[p*2]+d[p*2+1];
}
ll getsum(ll l,ll r,ll s,ll t,ll p)
{
    ll res=0;
    if(l<=s&&t<=r)
    {
        return d[p];
    }
    ll mid=(s+t)>>1;
    pushdown(s,t,p);
    if(l<=mid)
        res+=getsum(l,r,s,mid,p*2);
    if(r>mid)
        res+=getsum(l,r,mid+1,t,p*2+1);
    return res;
}
int main()
{
    cin>>n>>m;
    for(int i=1;i<=n;i++)
        scanf("%lld",&a[i]);
    build(1,n,1);
    while (m -- )
    {
        cin>>op;
        if(op=="Q")
        {
            scanf("%lld%lld",&l,&r);
            printf("%lld\n",getsum(l,r,1,n,1));
        }
        else
        {
            scanf("%lld%lld%lld",&l,&r,&k);
            add(l,r,k,1,n,1);
        }
    }
    return 0;
}

242. 一个简单的整数问题 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
typedef long long ll;
const int maxn=1e5+1; 
ll a[maxn],tree[maxn],n,m,l,r,k;
char op[2];
ll lowbit(ll x)
{
    return x&-x;
}
void add(ll x,ll k)
{
    while(x<=n)
    {
        tree[x]+=k;
        x+=lowbit(x);
    }
}
ll getsum(ll x)
{
    ll res=0;
    while(x)
    {
        res+=tree[x];
        x-=lowbit(x);
    }
    return res;
}
int main()
{
    cin>>n>>m;
    for (int i = 1; i <= n; i ++ )
        scanf("%lld",a+i);
    while(m--)
    {
        scanf("%s",op);
        if(strcmp(op,"Q")==0)
        {
            scanf("%lld",&l);
            printf("%lld\n",a[l]+getsum(l));
        }
        else if(strcmp(op,"C")==0)
        {
            scanf("%lld%lld%lld",&l,&r,&k);
            add(l,k);add(r+1,-k);
        }
    }
    return 0;
}

Matlab中的随机数生成函数

前言

　　Matlab 统计工具箱中有 27 种概率分布，工具箱对每一种分布都提供 5 类函数，其命令的字符是：
pdf 概率密度； cdf 分布函数； inv 分布函数的反函数；stat 均值与方差； rnd 随机数生成
　　这里主要介绍常用的rnd类函数

rand 函数

1
2
rand(n) %返回一个 n×n 在区间 (0,1) 内均匀分布的随机数矩阵
rand(m,n) %返回一个 m×n 在区间 (0,1) 内均匀分布的随机数矩阵

unifrnd 函数

1
unifrnd(a,b,m,n) %返回一个 m×n 在区间 (a,b) 内均匀分布的随机数矩阵

randn 函数

1
2
randn(n) %返回一个 n×n 标准正态分布的随机数矩阵
randn(m,n) %返回一个 m×n 标准正态分布的随机数矩阵

exprnd 函数

1
exprnd(lambda,m,n) %返回一个均值为lambda的 m×n 的指数分布随机数矩阵

normrnd 函数

1
normrnd(mu,sigma,m,n) %返回一个满足N(mu,sigma)的 m×n 的正态分布随机数矩阵

poissrnd 函数

1
poissrnd(lambda,m,n) %返回一个参数为lambda的 m×n 的泊松分布随机数矩阵

binornd 函数

1
binornd(n,p,m,n) %返回一个满足B(n,p)的 m×n 的二项分布随机数矩阵

randperm 函数

1
2
randperm(n) %返回行向量，其中包含从 1 到 n 没有重复元素的整数随机排列
randperm(n,k) %返回行向量，其中包含在 1 到 n 之间随机选择的 k 个唯一整数(1-n随机选择k个互不相同的整数)

perms 函数

1
2
perms(n) %返回的矩阵包含了向量 v 中元素按字典序反序的所有排列。  
% P 的每一行包含 v 中 n 个元素的一个不同排列。矩阵 P 具有与 v 相同的数据类型，包含 n! 行和 n 列。

附录

　　以下介绍几种绘制数据分布图的方法，建议配合以上的随机数生成函数学习：

hist 函数

1
2
3
4
[N,X] = hist(Y,M); 
%得到数组（行、列均可）Y 的频数表。
%它将区间[min(Y),max(Y)]等分为 M 份（缺省时M 设定为 10）。
%N 返回 M 个小区间的频数，X 返回 M 个小区间的中点。

cdfplot 函数

1
2
3
cdfplot(X) %作样本X（向量）的经验累积分布(也就是F(X))函数图形
%如果表头(title)出现矩形则是中文显示异常，可以通过设置字体解决
%set(0,'defaultAxesFontName', 'Microsft YaHei UI');实测有效，就是每次启动都要设置

参考文献以及资料

1、司守奎数学建模算法课件（2版含源程序）

2、https://max.book118.com/html/2017/0702/119736142.shtm

3、南信大的数学建模课件（Monte-Carlo 方法介绍及其建模应用）

96. 奇怪的汉诺塔 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
int d[15],f[15];
int main()
{
    d[1]=1;
    for(int i=2;i<=12;i++)
        d[i]=d[i-1]*2+1;
    memset(f,0x3f,sizeof(f));
    f[1]=1;
    for(int i=2;i<=12;i++)
    {
        for(int j=1;j<i;j++)
        {
            f[i]=min(f[i],f[j]*2+d[i-j]);
        }
    }
    for (int i = 1; i <= 12; i ++ )
        cout << f[i]<<endl;
    return 0;
}

95. 费解的开关 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;
char g[10][10];
int n,t;
int dx[]={0,0,0,1,-1};
int dy[]={0,1,-1,0,0};
void turn(int x,int y)
{
    for(int i=0;i<5;i++)
    {
        int xx=x+dx[i];
        int yy=y+dy[i];
        if(xx>=0&&xx<5&&yy>=0&&yy<5)
            g[xx][yy]^=1;
    }
}
void print()
{
    for(int i=0;i<5;i++)
        cout << g[i]<<endl;
}
int work()
{
    int ans=0x3f3f3f3f;
    for(int k=0;k<(1<<5);k++)
    {
        char buf[10][10];
        memcpy(buf,g,sizeof(g));
        int res=0;
        for(int i=0;i<5;i++)
        {
            if((k>>i)&1)
            {
                res++;
                turn(0,i);
                // print();
                // cout<<endl;
            }
        }
        for(int i=0;i<4;i++)
        {
            for(int j=0;j<5;j++)
            {
                if(g[i][j]=='0')
                {
                    res++;
                    turn(i+1,j);
                    // print();
                    // cout<<endl;
                }
            }
        }
        bool ok=1;
        for(int i=0;i<5;i++)
        {
            if(g[4][i]=='0')
            {
                ok=0;
                break;
            }
        }
        //cout<<res<<" 01 "<<ok<<endl;
        //print();
        if(ok)
            ans=min(ans,res);
        memcpy(g,buf,sizeof(g));
    }
    //cout << ans<<" "<<res<<endl;
    if(ans>6)
        return -1;
    else
        return ans;
}
int main()
{
    cin >> t;
    while(t--)
    {
        for(int i=0;i<5;i++)
            cin>>g[i];
        // print();
        // cout<<endl;
        cout << work()<<endl;
    }
    return 0;
}

94. 递归实现排列型枚举 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#include <bits/stdc++.h>
using namespace std;
int n,a[10];
int main()
{
    cin >> n;
    for(int i=1;i<=n;i++)
        a[i]=i;
    do
    {
        for(int i=1;i<=n;i++)
            printf("%d ",a[i]);
        printf("\n");
    }while(next_permutation(a+1,a+1+n));
    return 0;
}
  • Copyrights © 2022-2026 CPY
  • Visitors: | Views:

请我喝杯咖啡吧~

支付宝
微信