最优化方法

朝歌夜弦，为秦宫人

已于 2023-07-26 19:50:04 修改

阅读量1.8k

点赞数 6

分类专栏：数学建模文章标签：图论算法深度优先

于 2023-07-26 16:56:26 首次发布

本文链接：https://blog.csdn.net/weixin_44575377/article/details/131942955

版权

数学建模专栏收录该内容

3 篇文章 0 订阅

订阅专栏

一. 图论

1.最小生成树

图的生成树是它的一颗含有其所有顶点的无环连通子图,一幅加权图的最小生成树(MST)是它的一颗权值(树中的所有边的权值之和) 最小的生成树
• 适用场景：道路规划、通讯网络规划、管道铺设、电线布设等

题目数据
kruskal算法
稀疏图，按边大小排序

#include<bits/stdc++.h>
using namespace std;
const int maxn=5e3+3,maxm=2e5+5;

int n,m;
struct node{
    int fr,to,val;
    inline bool operator<(node b)const {
        return val<b.val;
    }
}e[maxm<<1];
int k,fa[maxn];
inline void add(int x,int y,int z){e[++k]=(node){x,y,z};}

inline int find(int x){
    return fa[x]==x?x:fa[x]=find(fa[x]);
}
int main()
{
    scanf("%d %d",&n,&m);
    int x,y,z;
    for(int i=1;i<=m;i++) 
    {
      	scanf("%d %d %d",&x,&y,&z);
        add(x,y,z);
        add(x,y,z);//无向图 
    }
    sort(e+1,e+k+1); //边排序
    
    for(int i=1;i<=n;i++)fa[i]=i;//并查集初始化
    int cnt=0,ans=0;
    for(int i=1;i<=k;i++)
    {
        int x=find(e[i].fr),y=find(e[i].to);
        if(x!=y)
        {
            fa[x]=y;//连边
            ans+=e[i].val;
            ++cnt;
            if(cnt==n-1)break;
        }
    }
    if(cnt==n-1)printf("%d",ans);
    else puts("orz");//图不连通
    return 0;
}

prim
稠密图

#include<bits/stdc++.h>
#define re return
#define inc(i,l,r) for(int i=l;i<=r;++i)
using namespace std;

const int maxn=5e3+3,maxm=2e5+5;
int inf=5e8;
int n,m,d[maxn][maxn],dis[maxn],vis[maxn];
int main()
{
   // freopen("a.in","r",stdin);
   	scanf("%d %d",&n,&m);
    int x,y,z;
    inc(i,1,n)inc(j,1,n)d[i][j]=inf;
    inc(i,1,m)
    {
        scanf("%d %d %d",&x,&y,&z);
        d[x][y]=d[y][x]=min(d[x][y],z);
    }
    
    vis[1]=1;//将1号点加入图中 
    inc(i,2,n)dis[i]=d[1][i];//dis初始化为1到各个点之间的最短距离 
    
    int ans=0,cnt=n-1;
    while(cnt--)//需要进行n-1次循环 
    {
        int now,minn=inf;
        inc(i,1,n)
        if(!vis[i]&&dis[i]<minn)minn=dis[i],now=i;//找当前最短距离点 now
        
     
        if(minn==inf)
        {
        	printf("orz");
        	re 0;
		}   
		vis[now]=1;
        ans+=dis[now];
        
        inc(i,1,n)dis[i]=min(dis[i],d[now][i]);
    }
    printf("%d",ans);
    re 0;
}

2.网络流

通常可以把这些边想象成道路，流量就是这条道路的车流量，容量就是道路可承受的最大的车流量 • 适用场景：企业生产运输问题、交通拥堵优化问题等

最大流
题目数据

#include<bits/stdc++.h>
#define re return
#define inc(i,l,r) for(int i=l;i<=r;++i)
using namespace std;
template<typename T>inline void rd(T&x)
{
	char c;bool f=0;
	while((c=getchar())<'0'||c>'9')if(c=='-')f=1;
	x=c^48;
	while((c=getchar())>='0'&&c<='9')x=x*10+(c^48);
	if(f)x=-x;
}
const int maxn=1e4+4,maxm=1e5+5; 
int n,m,hd[maxn],s,t;
struct node{
	int to,nt,w;
}e[maxm<<1];
int k=1;
inline void add(int x,int y,int z)
{
	e[++k]=(node){y,hd[x],z};hd[x]=k;
	e[++k]=(node){x,hd[y],0};hd[y]=k;
}

int deep[maxn],cur[maxn];
inline bool bfs()
{
	inc(i,1,n)deep[i]=0;
	deep[s]=1;
	queue<int>q;
	q.push(s);
	while(!q.empty())
	{
		int x=q.front();
		q.pop();
		for(int i=hd[x];i;i=e[i].nt)
		{
			int v=e[i].to;
			if(!deep[v]&&e[i].w)
			{
				deep[v]=deep[x]+1;
				if(v==t)re 1;
				q.push(v);
			}
		}
	}
	re 0;
}

inline int dfs(int x,int flow)
{
	if(x==t)re flow;
	int delta=flow;
	for(int &i=cur[x];i;i=e[i].nt)
	{
		int v=e[i].to;
		if(deep[v]==deep[x]+1&&e[i].w)
		{
			int d=dfs(v,min(e[i].w,delta));
			e[i].w-=d;e[i^1].w+=d;
			delta-=d;
			if(!delta)re flow;
		}
	}
	re flow-delta;
}
int main()
{
	rd(n),rd(m),rd(s),rd(t);
	int x,y,w;
	inc(i,1,m)
	{
		rd(x),rd(y),rd(w);
		add(x,y,w);
	}
	
	int ans=0;
	while(bfs())
	{
		inc(i,1,n)cur[i]=hd[i];
		ans+=dfs(s,2147483647);
	}
	printf("%d",ans);
	re 0;
}

最小费用最大流
最大流量的基础上要求最小的费用，有边权值
题目数据

#include<bits/stdc++.h>
#define re return
#define inc(i,l,r) for(int i=l;i<=r;++i)
using namespace std;
template<typename T>inline void rd(T&x)
{
	char c;bool f=0;
	while((c=getchar())<'0'||c>'9')if(c=='-')f=1;
	x=c^48;
	while((c=getchar())>='0'&&c<='9')x=x*10+(c^48);
	if(f)x=-x;
}
const int maxn=5e3+3;

int k=1,n,m,s,t;
int dis[maxn],ord[maxn],hd[maxn],flow[maxn],vis[maxn];
int smoney,sflow;

struct node{
	int flow,to,nt,cost;
}e[100005];
inline void add(int u,int v,int w,int f)
{
	e[++k].to=v;e[k].nt=hd[u];e[k].flow=w;e[k].cost=f;hd[u]=k;
}

inline bool spfa()
{
	queue<int>q;
	memset(dis,0x3f3f3f3f,sizeof dis);
	memset(vis,0,sizeof vis);
	q.push(s);
	dis[s]=0;
	flow[s]=0x3f3f3f3f;
	while(!q.empty())
	{
		int x=q.front();
		q.pop();
		vis[x]=0;
		for(int i=hd[x];i;i=e[i].nt)
		{
			int v=e[i].to,w=e[i].flow,f=e[i].cost;
			if(w&&dis[v]>dis[x]+f)
			{
				if(!vis[v])
				{
					vis[v]=1;
					q.push(v);
				}
				flow[v]=min(flow[x],w);
				dis[v]=dis[x]+f;
				ord[v]=i;
			}
		}
	}
	re dis[t]!=0x3f3f3f3f;
}

inline void vivi()
{
	int x=t;
	while(x!=s)
	{
		int i=ord[x];
		e[i].flow-=flow[t];
		e[i^1].flow+=flow[t];
		x=e[i^1].to;
	}
	sflow+=flow[t];
	smoney+=flow[t]*dis[t];
}
int main()
{
//	freopen("a.in","r",stdin);
	rd(n),rd(m),rd(s),rd(t);
	int u,v,w,f;
	inc(i,1,m)
	{
		rd(u),rd(v),rd(w),rd(f);
		add(u,v,w,f);
		add(v,u,0,-f);
	}
	while(spfa())
		vivi();
	printf("%d %d",sflow,smoney);
	re 0;
}

3.最短路

主要包括Dijkstra算法和Floyd算法两种，用于求解两点间的最短距离 • 适用场景：路径规划问题，如修建道路、设定救援路线等
题目数据
在这里插入图片描述

spfa

单源最短路
容易被卡

#include<bits/stdc++.h>
#define re return
#define inc(i,l,r) for(int i=l;i<=r;++i)
using namespace std;
template<typename T>inline void rd(T&x)
{
    char c;bool f=0;
    while((c=getchar())<'0'||c>'9')if(c=='-')f=1;
    x=c^48;
    while((c=getchar())>='0'&&c<='9')x=x*10+(c^48);
    if(f)x=-x;
}

const int maxn=1e4+5,maxm=5e5+5;
int n,m,hd[maxn];
struct node{
    int to,nt,val;
}e[maxm<<1];
int k,s;
inline void add(int x,int y,int z){
    e[++k]=(node){y,hd[x],z};hd[x]=k;
}

#define ll long long
ll inf=2147483647,dis[maxn],vis[maxn];
inline void spfa()
{
    inc(i,1,n)dis[i]=inf;
    dis[s]=0;
    queue<int>q;
    q.push(s);
    while(!q.empty())
    {
        int x=q.front();
        q.pop();
        vis[x]=0;
        for(int i=hd[x];i;i=e[i].nt)
        {
            int v=e[i].to;
            if(dis[v]>dis[x]+e[i].val)
            {
                dis[v]=dis[x]+e[i].val;
                if(!vis[v])
                {
                    vis[v]=1;
                        q.push(v);    
                }
            }
        }
    }
}

int main()
{
   // freopen("a.in","r",stdin);
    rd(n),rd(m),rd(s);
    int x,y,z;
    inc(i,1,m)
    {
        rd(x),rd(y),rd(z);
        add(x,y,z);
    }
    
    spfa();
    
    inc(i,1,n)printf("%lld ",dis[i]);
    re 0;
}

dijsktra

单源最短路
无负边

#include<bits/stdc++.h>
#define re return
#define inc(i,l,r) for(int i=l;i<=r;++i)
using namespace std;
template<typename T>inline void rd(T&x)
{
    char c;bool f=0;
    while((c=getchar())<'0'||c>'9')if(c=='-')f=1;
    x=c^48;
    while((c=getchar())>='0'&&c<='9')x=x*10+(c^48);
    if(f)x=-x;
}

const int maxn=1e4+5,maxm=5e5+5;
int n,m,hd[maxn];
struct node{
    int to,nt,val;
}e[maxm<<1];
int k,s;
inline void add(int x,int y,int z){
    e[++k]=(node){y,hd[x],z};hd[x]=k;
}

#define ll long long
ll inf=2147483647,dis[maxn];
struct KKK{
int x,val;
inline bool operator<(KKK u)const 
{
    re val>u.val;
}
};
inline void dij()
{
    inc(i,1,n)dis[i]=inf;
    dis[s]=0;
    priority_queue<KKK>q;
    q.push((KKK){s,0});
    while(!q.empty())
    {
        KKK u=q.top();
        q.pop();
        int x=u.x;
        if(dis[x]!=u.val)continue;
        for(int i=hd[x];i;i=e[i].nt)
        {
            int v=e[i].to;
            if(dis[v]>dis[x]+e[i].val)
            {
                dis[v]=dis[x]+e[i].val;
                q.push((KKK){v,dis[v]});
            }
        }
    }
}

int main()
{
    //freopen("a.in","r",stdin);
    rd(n),rd(m),rd(s);
    int x,y,z;
    inc(i,1,m)
    {
        rd(x),rd(y),rd(z);
        add(x,y,z);
    }
    
    dij();
    
    inc(i,1,n)printf("%lld ",dis[i]);
    re 0;
}

floyd

O(n^3)
多源最短路

#include<bits/stdc++.h>
using namespace std;
#define inf 1234567890
#define maxn 10005
inline int read()
{
    int x=0,k=1; char c=getchar();
    while(c<'0'||c>'9'){if(c=='-')k=-1;c=getchar();}
    while(c>='0'&&c<='9')x=(x<<3)+(x<<1)+(c^48),c=getchar();
    return x*k;
}//快读
int a[maxn][maxn],n,m,s;
inline void floyd()
{
    for(int k=1;k<=n;k++)
    //这里要先枚举k（可以理解为中转点）
	{
        for(int i=1;i<=n;i++)
		{
            if(i==k||a[i][k]==inf)
            {
                continue;
			}
			for(int j=1;j<=n;j++)
			{
                a[i][j]=min(a[i][j],a[i][k]+a[k][j]);
                //松弛操作，即更新每两个点之间的距离
                //松弛操作有三角形的三边关系推出
                //即两边之和大于第三边
            }
        }
    }
}
int main()
{
    n=read(),m=read(),s=read();
    for(int i=1;i<=n;i++)
	{
        for(int j=1;j<=n;j++)
		{
            a[i][j]=inf;
        }
    }
	//初始化，相当于memset(a,inf,sizeof(a))
    for(int i=1,u,v,w;i<=m;i++)
	{
        u=read(),v=read(),w=read();
        a[u][v]=min(a[u][v],w);
        //取min可以对付重边
    }
    floyd();
    a[s][s]=0;
    for(int i=1;i<=n;i++)
    {
        printf("%d ",a[s][i]);
    }
    return 0;
}

二、动态规划

dp数学建模的应用不多，感觉主要还是各种环境下的一个背包（多维dp），对于状压，数位都不怎么涉及
总的思想来说偏向背包较为单一，找递归方程时的问题主要在于结合其他信息，可能涉及到概率（马尔科夫链），图（偏向树形dp之类的）

背包
采药
一维

#include<bits/stdc++.h>
using namespace std;
int a[100001];
int o;
int main()
{
	int t,m;
	cin>>t>>m;
int w[m+1],v[m+1];
for(int i=1;i<=m;i++)
cin>>w[i]>>v[i];
	
	for(int j=1;j<=t;j++)	
	for(int k=1;k<=m;k++)
		{
		    a[j]=o;
		    if(j>=w[k])
		    {a[j]=max(a[j],a[j-w[k]]+v[k]); 
		    o=max(o,a[j]);}
		}
			
	cout<<a[t];
	return 0;
}

二维

#include<bits/stdc++.h>
using namespace std;
int a[101][1001];
int main()
{
int n,m;
cin>>n>>m;
int w[m+1];
int v[m+1];
for(int i=1;i<=m;i++)
cin>>w[i]>>v[i];

for(int i=1;i<=m;i++)
for(int j=1;j<=n;j++)
{ a[i][j]=a[i-1][j];
if(j>=w[i])a[i][j]=max(a[i-1][j-w[i]]+v[i],a[i][j]); 

}
cout<<a[m][n];
return 0;
}

2.树形dp
在连通图上dp
没有上司的舞会

#include<bits/stdc++.h>
using namespace std;
int f[6005][3],head[6005],cc[6005],in[6005];
int n,a,b,c;

void dfs(int x)
{   
	for(int i=head[x];i!=0;i=cc[i])
 	{	
 	  dfs(i);
	 f[x][1]=f[x][1]+f[i][0];
	 f[x][0]+=max(f[i][0],f[i][1]);
	}
}

int main()
{  

scanf("%d",&n);
for(int i=1;i<=n;i++)
scanf("%d",&f[i][1]);

for(int i=1;i<n;i++)
{scanf("%d%d",&a,&b); 
cc[a]=head[b];
head[b]=a;in[a]=1;}

for(int i=1;i<=n;i++)
if(in[i]==0){c=i;break;}
dfs(c);
printf("%d",max(f[c][1],f[c][0]));

return 0;
}

三、启发式算法

模拟退火
模拟退火算法（Simulate Anneal，SA）是一种通用概率演算法，用来在一个大的搜寻空间内找寻命题的最优解。
https://www.cnblogs.com/flashhu/p/8884132.html

#include<bits/stdc++.h>
#define re return
#define st static
#define mem(a,b) memset((a),(b),sizeof(a))
#define inc(i,l,r) for(register int i=l;i<=r;++i)
#define dec(i,l,r) for(register int i=l;i>=r;--i)

using namespace std;
int n,x[1005],y[1005],w[1005];
double EPS=10e-15,D=0.975;

double calc(double x0,double y0)
{
	double res=0;
	inc(i,1,n)res+=sqrt((x[i]-x0)*(x[i]-x0)+(y[i]-y0)*(y[i]-y0))*w[i];
	re res;
}

int main()
{
//	freopen("in.txt","r",stdin);
	double x1=0,y1=0,x0,y0,bx,by,res,ans,best,time=2;
	scanf("%d",&n);
	inc(i,1,n){
		scanf("%d%d%d",&x[i],&y[i],&w[i]);
		bx+=x[i];by+=y[i];
	}
	
	best=calc(bx/=n,by/=n);//用平均值作为最初解
	srand(20130317);
	while(time--)
	{
		x1=bx,y1=by,ans=best;
		for(double T=100000;T>EPS;T*=D)
		{
			x0=x1+T*(2*rand()-RAND_MAX);y0=y1+T*(2*rand()-RAND_MAX);
			res=calc(x0,y0);
			if(res<best){best=res;bx=x0;by=y0;}
			if(res<ans||exp((ans-res)/T)>(double)(rand())/RAND_MAX)
			{ans=res;x1=x0;y1=y0;}
		}
	}
	printf("%.3lf %.3lf",bx,by);
	re 0;
}