simplex.cpp
来自「改进的单纯形算法」· C++ 代码 · 共 973 行 · 第 1/2 页
CPP
973 行
m_prolongationvalue += pow(ppobjvalue[i],2.0);
}
if(m_prolongationvalue<m_lowvalue)
{
/* 继续延伸 */
while((m_prolongationvalue<=m_lowvalue)&&(flag<4))
{
/* 去掉原来的最大值点,取而代之以原来的最小点,
并用延伸点替换单纯形中目标值最小的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_variable[m_lownumber][i];
p_variable[m_lownumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = p_fitness[m_lownumber];
m_highvalue = m_lowvalue;
p_fitness[m_lownumber] = m_prolongationvalue;
m_lowvalue = m_prolongationvalue;
for(i=0;i<m_numpoint;i++)
{
p_bestobjvalue[i] = ppobjvalue[i];
}
for(i=0; i<m_numvar; i++)
{
if(p_variable[m_highnumber][i]>p_variable[m_lownumber][i])
p_prolongation[i] = p_variable[m_lownumber][i]-m_gamma*p_center[i];
else
p_prolongation[i] = p_variable[m_lownumber][i]+m_gamma*p_center[i];
}
/* 将参数存储到文件 */
WriteToFem(p_prolongation);
if(m_blErrorExit)return;
/* 有限元计算 */
ObjFunction();
if(m_blErrorExit) return;
/* 读入测点的计算值 */
ReadFromFem(p_calculation);
if(m_blErrorExit)return;
m_prolongationvalue = 0.0;
/* 计算测点的观测值与计算值之差 */
for(i=0; i<m_numpoint; i++)
{
ppobjvalue[i] = p_observation[i]-p_calculation[i];
m_prolongationvalue += pow(ppobjvalue[i],2.0);
}
}
if(m_prolongationvalue<m_lowvalue)
{
/* 去掉原来的最大值点,取而代之以原来的最小点,
并用延伸点替换单纯形中目标值最小的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_variable[m_lownumber][i];
p_variable[m_lownumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = p_fitness[m_lownumber];
m_highvalue = m_lowvalue;
p_fitness[m_lownumber] = m_prolongationvalue;
m_lowvalue = m_prolongationvalue;
for(i=0;i<m_numpoint;i++)
{
p_bestobjvalue[i] = ppobjvalue[i];
}
}
else
{
/* 用延伸点替换单纯形中目标值最大的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = m_prolongationvalue;
m_highvalue = m_prolongationvalue;
/* 收缩 */
Constriction();
if(m_blErrorExit)return;
}
}
else
{
if(m_prolongationvalue<m_highvalue)
{
/* 用延伸点替换单纯形中目标值最大的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = m_prolongationvalue;
m_highvalue = m_prolongationvalue;
/* 收缩 */
Constriction();
if(m_blErrorExit)return;
}
else
{
/* 延伸失败,调整延伸系数,重新延伸 */
flag = 0;
while((m_prolongationvalue>m_highvalue)&&(flag<4))
{
/* 调整延伸系数 */
m_gamma /= 2.0;
for(i=0; i<m_numvar; i++)
{
if(p_variable[m_highnumber][i]>p_variable[m_lownumber][i])
p_prolongation[i] = p_variable[m_lownumber][i]-m_gamma*p_center[i];
else
p_prolongation[i] = p_variable[m_lownumber][i]+m_gamma*p_center[i];
}
/* 将参数存储到文件 */
WriteToFem(p_prolongation);
if(m_blErrorExit)return;
/* 有限元计算 */
ObjFunction();
if(m_blErrorExit) return;
/* 读入测点的计算值 */
ReadFromFem(p_calculation);
if(m_blErrorExit)return;
m_prolongationvalue = 0.0;
/* 计算测点的观测值与计算值之差 */
for(i=0; i<m_numpoint; i++)
{
ppobjvalue[i] = p_observation[i]-p_calculation[i];
m_prolongationvalue += pow(ppobjvalue[i],2.0);
}
}
m_gamma = 1.0;
if(m_prolongationvalue<m_lowvalue)
{
/* 去掉原来的最大值点,取而代之以原来的最小点,
并用延伸点替换单纯形中目标值最小的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_variable[m_lownumber][i];
p_variable[m_lownumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = p_fitness[m_lownumber];
m_highvalue = m_lowvalue;
p_fitness[m_lownumber] = m_prolongationvalue;
m_lowvalue = m_prolongationvalue;
for(i=0;i<m_numpoint;i++)
{
p_bestobjvalue[i] = ppobjvalue[i];
}
}
else
{
/* 用延伸点替换单纯形中目标值最大的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = m_prolongationvalue;
m_highvalue = m_prolongationvalue;
/* 收缩 */
Constriction();
if(m_blErrorExit)return;
}
}
}
delete [] ppobjvalue;
}
/* 收缩 */
void Simplex::Constriction()
{
int i,j;
double * pcobjvalue;
pcobjvalue = new double[m_numpoint];
/****************** 求收缩点及其目标值 ******************/
for(i=0; i<m_numvar; i++)
{
if(p_variable[m_highnumber][i]>p_variable[m_lownumber][i])
m_constriction[i] = p_variable[m_highnumber][i]-m_beta*p_center[i];
else
m_constriction[i] = p_variable[m_highnumber][i]+m_beta*p_center[i];
}
/* 将参数存储到文件 */
WriteToFem(m_constriction);
if(m_blErrorExit)return;
/* 有限元计算 */
ObjFunction();
if(m_blErrorExit) return;
/* 读入测点的计算值 */
ReadFromFem(p_calculation);
if(m_blErrorExit)return;
m_constrictionvalue = 0.0;
/* 计算测点的观测值与计算值之差 */
for(i=0; i<m_numpoint; i++)
{
pcobjvalue[i] = p_observation[i]-p_calculation[i];
m_constrictionvalue += pow(pcobjvalue[i],2.0);
}
/****************** 求收缩点及其目标值 ******************/
if(m_constrictionvalue<m_lowvalue)
{
/* 去掉原来的最大值点,取而代之以原来的最小点,
并用收缩点替换单纯形中目标值最小的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_variable[m_lownumber][i];
p_variable[m_lownumber][i] = m_constriction[i];
}
p_fitness[m_highnumber] = p_fitness[m_lownumber];
m_highvalue = m_lowvalue;
p_fitness[m_lownumber] = m_constrictionvalue;
m_lowvalue = m_constrictionvalue;
for(i=0;i<m_numpoint;i++)
{
p_bestobjvalue[i] = pcobjvalue[i];
}
}
else
{
if(m_constrictionvalue<m_highvalue)
{
/* 用收缩点替换单纯形中目标值最大的点 */
for(i=0; i<m_numvar; i++)
{
p_variable[m_highnumber][i] = p_prolongation[i];
}
p_fitness[m_highnumber] = m_prolongationvalue;
m_highvalue = m_prolongationvalue;
/* 缩小边长 */
for(i=0; i<=m_numvar; i++)
for(j=0; j<m_numvar; j++)
p_variable[i][j] = 0.5*(p_variable[m_lownumber][j]+p_variable[i][j]);
/******** 重新计算新的单纯形的目标值 ********/
for(i=0; i<=m_numvar; i++)
{
p_fitness[i] = 0.0;
/* 将参数存储到文件 */
WriteToFem(p_variable[i]);
if(m_blErrorExit)return;
/* 有限元计算 */
ObjFunction();
if(m_blErrorExit) return;
/* 读入测点的计算值 */
ReadFromFem(p_calculation);
if(m_blErrorExit)return;
/* 计算测点的观测值与计算值之差的平方 */
for(j=0; j<m_numpoint; j++)
{
pcobjvalue[j] = p_observation[j]-p_calculation[j];
p_fitness[i] += pow(pcobjvalue[j],2.0);
}
}
}
}
}
/* 有限元计算 */
void Simplex::ObjFunction()
{
if(EXCUTE_NUMERICAL_ANALYSIS()==false)
m_blErrorExit=true;
}
/* 写参数值 */
void Simplex::WriteToFem(double *pvariable)
{
CString str;
CString csBackOutFileName = mGcsFileName+_T(".btof");//back analysis to fem
ofstream cBackOutFile;
cBackOutFile.open((LPCSTR)csBackOutFileName);
if (cBackOutFile.fail())
{
AfxMessageBox("Cannot open input file of back analysis!!!",MB_OK);
m_blErrorExit = true;
return;
}
//write
int index=1;
cBackOutFile << m_numvar << endl;
for(int i=0; i<m_numvar; i++)
{
str.Format("%d %d %d %d %.6f\n",index++, p_BAKindflag[i],
p_iparameter[i], p_jparameter[i], pvariable[i]);
cBackOutFile << (LPCSTR)str;
}
//边界荷载输出
if (m_numbndload > 0)
{
index=1;
cBackOutFile << m_numbndload << endl;
for(int j=0; j<m_numbndload; j++)
{
str.Format("%d %d %d %d %d %d %.4f %.4f %.4f %.4f %.4f %.4f\n",
p_Bndload[j].index++, p_Bndload[j].wAddStage, p_Bndload[j].wAddStep,
p_Bndload[j].iNode1, p_Bndload[j].iNode2, p_Bndload[j].iCoordFlag,
p_Bndload[j].BeginLoad1, p_Bndload[j].BeginLoad2, p_Bndload[j].BeginLoad3,
p_Bndload[j].EndLoad1, p_Bndload[j].EndLoad2, p_Bndload[j].EndLoad3);
cBackOutFile << (LPCSTR)str;
}
}
cBackOutFile.close();
}
/* 读取计算值 */
void Simplex::ReadFromFem(double *calculatevalue)
{
CString csFromFemFileName = mGcsFileName+_T(".ftob");//file from fem
ifstream cFromFemFile;
cFromFemFile.open((LPCSTR)csFromFemFileName);
if (cFromFemFile.fail())
{
AfxMessageBox("Cannot open input file of fem!!!",MB_OK);
m_blErrorExit = true;
return;
}
int index;
for(int i=0;i<m_numpoint;i++)
{
cFromFemFile >> index;
cFromFemFile >> calculatevalue[i];
}
cFromFemFile.close();
}
/* 从BAA文件中读入待反演参数的性质标志以及初始值 */
void Simplex::ReadBAA()
{
char tempchar[256];
int index;
CString csBaaFileName = mGcsFileName+_T(".BAA");
ifstream cBaaFile;
cBaaFile.open((LPCSTR)csBaaFileName);
if (cBaaFile.fail())
{
AfxMessageBox("Cannot open input BAA file!!!",MB_OK);
m_blErrorExit = true;
return;
}
cBaaFile.getline(tempchar,256);
//读入各个参数的上下限,并计算出各个参数的精度
for(int i=0; i<m_numvar; i++)
{
cBaaFile >> index;
cBaaFile >> p_BAKindflag[i];
cBaaFile >> p_iparameter[i];
cBaaFile >> p_jparameter[i];
cBaaFile >> p_variable[0][i];
cBaaFile.getline(tempchar,256);
}
cBaaFile >> m_numbndload;
if(m_numbndload>0)
{
if(p_Bndload)
{
delete[]p_Bndload ;
p_Bndload = NULL ;
}
p_Bndload = new bndload[m_numbndload];
for(int j=0; j<m_numbndload; j++)
{
cBaaFile >> p_Bndload[j].index >> p_Bndload[j].wAddStage
>> p_Bndload[j].wAddStep >> p_Bndload[j].iNode1
>> p_Bndload[j].iNode2 >> p_Bndload[j].iCoordFlag
>> p_Bndload[j].BeginLoad1 >> p_Bndload[j].BeginLoad2
>> p_Bndload[j].BeginLoad3 >> p_Bndload[j].EndLoad1
>> p_Bndload[j].EndLoad2 >> p_Bndload[j].EndLoad3;
cBaaFile.getline(tempchar,256);
}
}
cBaaFile.close();
}
//从MEA文件中读取监测数据
void Simplex::ReadMEA()
{
CString csFromMeaFileName = mGcsFileName+_T(".MEA");//file from fem
ifstream cFromMeaFile;
cFromMeaFile.open((LPCSTR)csFromMeaFileName);
if (cFromMeaFile.fail())
{
AfxMessageBox("Cannot open input file of mea!!!",MB_OK);
m_blErrorExit = true;
return;
}
int iTemp;
double dTemp;
char tempch[256];
cFromMeaFile.getline(tempch,256);
for(int i=0;i<m_numpoint;i++)
{
cFromMeaFile >> iTemp >> dTemp >> dTemp >> dTemp >> iTemp >> iTemp
>> iTemp >> iTemp >> iTemp >> iTemp >> iTemp;
cFromMeaFile >> p_observation[i];
cFromMeaFile.getline(tempch,256);
}
cFromMeaFile.close();
}
/* 输出每一代演化的最佳结果及其适应值 */
void Simplex::WriteResult()
{
CString csToUserFileName = mGcsFileName+_T(".result");//file to user
ofstream cResultFile;
cResultFile.open((LPCSTR)csToUserFileName,ios::app);
if (cResultFile.fail())
{
AfxMessageBox("Cannot open output result file!!!",MB_OK);
m_blErrorExit = true;
return;
}
cResultFile << "第" << m_numlap << "步优化结果:" << endl;
cResultFile << "适应值:" << m_lowvalue << endl;
cResultFile << "反演参数值依次为:";
for(int i=0;i<m_numvar;i++)
cResultFile << p_variable[m_lownumber][i] <<" ";
cResultFile << endl;
cResultFile << "各量测点量测值与计算值之差分别为:" ;
for(i=0; i<m_numpoint; i++)
cResultFile << p_bestobjvalue[i] <<" ";
cResultFile << endl << endl;
cResultFile.close();
}
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?