canshushezhi.cpp

来自「此代码可实现最优控制理论中的单纯型寻优算法」· C++ 代码 · 共 455 行

// CanshuShezhi.cpp : implementation file
//

#include "stdafx.h"
#include "optimization.h"
#include "CanshuShezhi.h"
#include "math.h"
#include "optimizationDoc.h"
#include "optimizationView.h"

#ifdef _DEBUG
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif

/////////////////////////////////////////////////////////////////////////////
// CCanshuShezhi dialog


CCanshuShezhi::CCanshuShezhi(CWnd* pParent /*=NULL*/)
	: CDialog(CCanshuShezhi::IDD, pParent)
{
	//{{AFX_DATA_INIT(CCanshuShezhi)
	m_I = 0.0f;
	m_K = 0.0f;
	m_P = 0.0f;
	m_R = 0.0f;
	m_T = 0.0f;
	m_ObjFun = _T("");
	m_n = 0;
	//}}AFX_DATA_INIT
}


void CCanshuShezhi::DoDataExchange(CDataExchange* pDX)
{
	CDialog::DoDataExchange(pDX);
	//{{AFX_DATA_MAP(CCanshuShezhi)
	DDX_Text(pDX, IDC_I, m_I);
	DDX_Text(pDX, IDC_K, m_K);
	DDX_Text(pDX, IDC_P, m_P);
	DDX_Text(pDX, IDC_R, m_R);
	DDX_Text(pDX, IDC_T, m_T);
	DDX_CBString(pDX, IDC_COMBO_ObjFun, m_ObjFun);
	DDX_Text(pDX, IDC_n, m_n);
	//}}AFX_DATA_MAP
}


BEGIN_MESSAGE_MAP(CCanshuShezhi, CDialog)
	//{{AFX_MSG_MAP(CCanshuShezhi)
	ON_EN_CHANGE(IDC_R, OnChangeR)
	ON_EN_CHANGE(IDC_T, OnChangeT)
	ON_EN_CHANGE(IDC_n, OnChangen)
	ON_EN_CHANGE(IDC_K, OnChangeK)
	ON_EN_CHANGE(IDC_P, OnChangeP)
	ON_EN_CHANGE(IDC_I, OnChangeI)
	//}}AFX_MSG_MAP
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CCanshuShezhi message handlers

void CCanshuShezhi::OnChangeR() 
{
	// TODO: If this is a RICHEDIT control, the control will not
	// send this notification unless you override the CDialog::OnInitDialog()
	// function and call CRichEditCtrl().SetEventMask()
	// with the ENM_CHANGE flag ORed into the mask.
	
	// TODO: Add your control notification handler code here
	UpdateData(true);
	
}

void CCanshuShezhi::OnChangeT() 
{
	// TODO: If this is a RICHEDIT control, the control will not
	// send this notification unless you override the CDialog::OnInitDialog()
	// function and call CRichEditCtrl().SetEventMask()
	// with the ENM_CHANGE flag ORed into the mask.
	
	// TODO: Add your control notification handler code here
	UpdateData(true);
}

void CCanshuShezhi::OnChangen() 
{
	// TODO: If this is a RICHEDIT control, the control will not
	// send this notification unless you override the CDialog::OnInitDialog()
	// function and call CRichEditCtrl().SetEventMask()
	// with the ENM_CHANGE flag ORed into the mask.
	
	// TODO: Add your control notification handler code here
	UpdateData(true);
}
void CCanshuShezhi::OnOK() 
{
	// TODO: Add extra validation here
	
	//(1)单纯形初始化
	int i=0,j=0,k=0;//循环变量
	double Xh=0,Xg=0,Xl=0;//最坏值，次坏值，最好值
	double mid=0;//中间变量
	//各个点的p，i，X值
	double ph=0,pl=0,pg=0,ih=0,ig=0,il=0,pc=0,ic=0,pr=0,ir=0,pe=0,ie=0,ps=0,is=0,Xc=0,Xr=0,Xe=0,Xs=0;
	double E=0.01;//精度要求
	double ansp=0,ansi=0;//最终p，i值
	double Q=0;//最终目标函数值
	double pm,im,Xm,pmm,imm,Xmm;//中间量
	double lamda=0.4,miu=1.5;//压缩因子，扩张因子
	double pini,iini,Xini,hp=0.1,hi=5,iup=250,pup=20,idow=1,pdow=0.01;

	pini=m_P;
	iini=m_I;
	Xini=simulink(m_R,pini,iini,m_K,m_n,m_T);
	pm=m_P+hp;
	im=iini;
	Xm=simulink(m_R,pm,im,m_K,m_n,m_T);
	pmm=m_P;
	imm=iini+hi;
	Xmm=simulink(m_R,pmm,imm,m_K,m_n,m_T);

	Xl=Xini;pl=pini;il=iini;
	if(Xl>Xm)//最好值
	{
		Xl=Xm;pl=pm;il=im;
	}
	if(Xl>Xmm)
	{
		Xl=Xmm;pl=pmm;il=imm;	
	}

	Xh=Xini;ph=pini;ih=iini;
	if(Xh<Xm)//最坏值
	{
		Xh=Xm;ph=pm;ih=im;	
	}
	if(Xh<Xmm)
	{
		Xh=Xmm;ph=pmm;ih=imm;
	}
	
	Xg=Xini;pg=pini;ig=iini;
	if(Xg==Xh||Xg==Xl)//次坏值
	{
		Xg=Xm;pg=pm;ig=im;
	}
	if(Xg==Xh||Xg==Xl)
	{
		Xg=Xmm;pg=pmm;ig=imm;
	}	
	

	for(k=0;k<=10000;k++)
	{
		if (fabs(Xh-Xl)/Xl<=E) //符合精度要求则退出
		{
			ansp=(ph+pl)/2;
			ansi=(ih+il)/2;			
			break;
		}

		if(k>=100)
		{
			MessageBox("搜索失败","目标函数值",MB_OK);
			ansp=(ph+pl)/2;
			ansi=(ih+il)/2;			
			break;
			//return;
		}

		if(pl>pup) pl=pup;
		if(pg>pup) pg=pup;
		if(ph>pup) ph=pup;
		if(pl<pdow) pl=pdow;
		if(pg<pdow) pg=pdow;
		if(ph<pdow) ph=pdow;

		if(il>iup) il=iup;
		if(ig>iup) ig=iup;
		if(ih>iup) ih=iup;
		if(il<idow) il=idow;
		if(ig<idow) ig=idow;
		if(ih<idow) ih=idow;
	

		pr=pl+pg-ph;//反射点
		ir=il+ig-ih;
		if(pr>pup) pr=pup;
		if(pr<pdow) pr=pdow;
		if(ir>iup) ir=iup;
		if(ir<idow) ir=idow;
		Xr=simulink(m_R,pr,ir,m_K,m_n,m_T);

		if(Xr>Xh)//反射点和最坏点互换
		{
			pm=pr;im=ir;Xm=Xr;
			pr=ph;ir=ih;Xr=Xh;
			ph=pm;ih=im;Xh=Xm;
		}
		if(Xr>Xg)//压缩单纯形
		{
			ps=(pl+pg)/2+lamda*((pl+pg)/2-ph);
			is=(il+ig)/2+lamda*((il+ig)/2-ih);
			if(ps>pup) ps=pup;
			if(ps<pdow) ps=pdow;
			if(is>iup) is=iup;
			if(is<idow) is=idow;
			Xs=simulink(m_R,ps,is,m_K,m_n,m_T);

			if(Xg>=Xs)//单纯形收缩
			{
				pini=ps; iini=is; Xini=Xs;
				pm=pl;	im=il;	Xm=Xl;
				pmm=pg;	imm=ig;	Xmm=Xg;
				
				Xl=Xini;pl=pini;il=iini;
				if(Xl>Xm)	{Xl=Xm;pl=pm;il=im;}
				if(Xl>Xmm)	{Xl=Xmm;pl=pmm;il=imm;}					
				
				Xh=Xini;ph=pini;ih=iini;
				if(Xh<Xm)	{Xh=Xm;ph=pm;ih=im;	}
				if(Xh<Xmm)	{Xh=Xmm;ph=pmm;ih=imm;}
				
				Xg=Xini;pg=pini;ig=iini;
				if(Xg==Xh||Xg==Xl)	{Xg=Xm;pg=pm;ig=im;}
				if(Xg==Xh||Xg==Xl)	{Xg=Xmm;pg=pmm;ig=imm;}
				
				continue;
			}
			if(Xg<Xs)//单纯形收缩
			{
				ph=(pl+ph)/2;
				pg=(pl+pg)/2;
				ih=(il+ih)/2;
				ig=(il+ig)/2;

				Xh=simulink(m_R,ph,ih,m_K,m_n,m_T);
				Xg=simulink(m_R,pg,ig,m_K,m_n,m_T);
				
				pini=ph; iini=ih; Xini=Xh;
				pm=pl;	im=il;	Xm=Xl;
				pmm=pg;	imm=ig;	Xmm=Xg;
				
				Xl=Xini;pl=pini;il=iini;
				if(Xl>Xm)	{Xl=Xm;pl=pm;il=im;}
				if(Xl>Xmm)	{Xl=Xmm;pl=pmm;il=imm;}					
				
				Xh=Xini;ph=pini;ih=iini;
				if(Xh<Xm)	{Xh=Xm;ph=pm;ih=im;	}
				if(Xh<Xmm)	{Xh=Xmm;ph=pmm;ih=imm;}
				
				Xg=Xini;pg=pini;ig=iini;
				if(Xg==Xh||Xg==Xl)	{Xg=Xm;pg=pm;ig=im;}
				if(Xg==Xh||Xg==Xl)	{Xg=Xmm;pg=pmm;ig=imm;}
				
				continue;				
			}
		}
		

		if((1-miu)*Xh+miu*Xr>=Xl)//Xr,Xl,Xg重新组成单纯形
		{
			pini=pr; iini=ir; Xini=Xr;
			pm=pl;	im=il;	Xm=Xl;
			pmm=pg;	imm=ig;	Xmm=Xg;
			
			
			Xl=Xini;pl=pini;il=iini;
			if(Xl>Xm)	{Xl=Xm;pl=pm;il=im;}
			if(Xl>Xmm)	{Xl=Xmm;pl=pmm;il=imm;}	
			
			Xh=Xini;ph=pini;ih=iini;
			if(Xh<Xm)	{Xh=Xm;ph=pm;ih=im;	}
			if(Xh<Xmm)	{Xh=Xmm;ph=pmm;ih=imm;}
			
			Xg=Xini;pg=pini;ig=iini;
			if(Xg==Xh||Xg==Xl)	{Xg=Xm;pg=pm;ig=im;}
			if(Xg==Xh||Xg==Xl)	{Xg=Xmm;pg=pmm;ig=imm;}
			
			continue;
		}

		if(Xr<Xg)//扩张单纯形
		{
			pe=(1-miu)*(pl+pg)/2-miu*ph;
			ie=(1-miu)*(il+ig)/2-miu*ih;
			if(pe>pup) pe=pup;
			if(pe<pdow) pe=pdow;
			if(ie>iup) ie=iup;
			if(ie<idow) ie=idow;
			Xe=simulink(m_R,pe,ie,m_K,m_n,m_T);
			
			if(Xe<Xr)
			{
				pr=pe;
				ir=ie;
				Xr=Xe;
			}

			pini=pr; iini=ir; Xini=Xr;
			pm=pl;	im=il;	Xm=Xl;
			pmm=pg;	imm=ig;	Xmm=Xg;
			
			
			Xl=Xini;pl=pini;il=iini;
			if(Xl>Xm)	{Xl=Xm;pl=pm;il=im;}
			if(Xl>Xmm)	{Xl=Xmm;pl=pmm;il=imm;	}
			
			Xh=Xini;ph=pini;ih=iini;
			if(Xh<Xm)	{Xh=Xm;ph=pm;ih=im;	}
			if(Xh<Xmm)	{Xh=Xmm;ph=pmm;ih=imm;}
			
			Xg=Xini;pg=pini;ig=iini;
			if(Xg==Xh||Xg==Xl)	{Xg=Xm;pg=pm;ig=im;}
			if(Xg==Xh||Xg==Xl)	{Xg=Xmm;pg=pmm;ig=imm;}
			
			//continue;			
		}
}
    Q=simulink(m_R,ansp,ansi,m_K,m_n,m_T);
	m_P=ansp;
	m_I=ansi;
	UpdateData(false);
	CString sq;
	sq.Format("%f",Q);
	MessageBox(sq,"目标函数值",MB_OK);

	flag=1;
	COptimizationDoc *pDoc=
	(COptimizationDoc *)(((CFrameWnd *)AfxGetApp()->m_pMainWnd)->GetActiveDocument());
	POSITION m_pos = pDoc->GetFirstViewPosition();
	COptimizationView* pView = (COptimizationView*)(pDoc->GetNextView(m_pos));
	pView->Invalidate();

	CDialog::OnOK();
}



BOOL CCanshuShezhi::OnInitDialog() 
{
	CDialog::OnInitDialog();
	
	// TODO: Add extra initialization here
	m_I =20;
	m_K = 0.53;
	m_n = 3;
	m_ObjFun = -1;
	m_P = 2;
	m_R = 1;
	m_T = 12.5;
	m_ObjFun = "绝对误差的矩的积分";
	UpdateData(false);	
	return TRUE;  // return TRUE unless you set the focus to a control
	              // EXCEPTION: OCX Property Pages should return FALSE
}



double CCanshuShezhi::simulink(double R,double P,double I,double K,int n,double T)
{
	//数据传递给文档类，用于保存，画图，in,out
	//文档类指针
COptimizationDoc *pDoc=
	(COptimizationDoc *)(((CFrameWnd *)AfxGetApp()->m_pMainWnd)->GetActiveDocument());
	//对象离散化，欧拉法仿真，计算目标函数，对象结构：（K）/[(aS+1)]n
	int i=0,j=0;
	double dt=1.0;//步长
	int st=1000;//仿真时间
	int lp;//仿真步数
	lp=int(st/dt);
	double Q=0;
	double x1=0,x2=0;
	double c,d;	
	c=exp(-dt/T);
	d=1-c;
	if (n>1)
	{		
		double *x=new double [n-1];
		for(i=0;i<n-1;i++)
		{
			x[i]=0.0;
		}
		for (i=0;i<lp;i++)//仿真
		{
			x1 = x1+dt*(R-x[n-2])/(P*I);
			x2=c*x2+K*d*(x1+R/P-x[n-2]/P);
			
			x[0]=x[0]*c+d*x2;
			for(j=1;j<n-1;j++)
			{
				x[j]=x[j]*c+d*x[j-1];
			}
			
			Q=Q+dt*dt*fabs(R-x[n-2])*i;
			pDoc->input[i]=R;
			pDoc->output[i]=x[n-2];
		}
		delete [n-1] x;
	}
	if (n==1)
	{		
		for (i=0;i<lp;i++)//仿真
		{
			x1 = x1+dt*(R-x2)/(P*I);
			x2=c*x2+K*d*(x1+R/P-x2/P);
			
			Q=Q+dt*dt*fabs(R-x2)*i;
			pDoc->input[i]=R;
			pDoc->output[i]=x2;
		
		
		}	
	}
	return Q;
}



void CCanshuShezhi::OnChangeK() 
{
	// TODO: If this is a RICHEDIT control, the control will not
	// send this notification unless you override the CDialog::OnInitDialog()
	// function and call CRichEditCtrl().SetEventMask()
	// with the ENM_CHANGE flag ORed into the mask.
	
	// TODO: Add your control notification handler code here
	UpdateData(true);
}

void CCanshuShezhi::OnChangeP() 
{
	// TODO: If this is a RICHEDIT control, the control will not
	// send this notification unless you override the CDialog::OnInitDialog()
	// function and call CRichEditCtrl().SetEventMask()
	// with the ENM_CHANGE flag ORed into the mask.
	
	// TODO: Add your control notification handler code here
	UpdateData(true);
}

void CCanshuShezhi::OnChangeI() 
{
	// TODO: If this is a RICHEDIT control, the control will not
	// send this notification unless you override the CDialog::OnInitDialog()
	// function and call CRichEditCtrl().SetEventMask()
	// with the ENM_CHANGE flag ORed into the mask.
	
	// TODO: Add your control notification handler code here
	UpdateData(true);
}