一般线性规划求最小值+c语言
2015-09-16 15:18
645 查看
求最小值的线性规划,仍然使用旋转算法。
对于如下的线性函数,求最优解:
初始表数据:
算法步骤:
一个实例:
c语言程序运行结果:
c语言代码如下:
#include "stdio.h"
#include "stdlib.h"
#include "math.h"
#include "vector"
using namespace std;
#define INF 9999
#define varNum 4 //变量数
#define equNum 0 //等式约束数
#define MAXITER 1000 //最大迭代次数
typedef vector<double> dim1Vector;
typedef vector<dim1Vector> dim2Vector;
dim2Vector getFileINFO(char *File1, char *File2); //获取文本信息
void Cal_LP(dim2Vector info); //解线性方程组
dim2Vector Cal_Spin(dim1Vector C, dim2Vector A, dim1Vector DEV, int laber_i, int laber_j); //旋转计算
void main()
{
char *File1 = "info.txt";
char *File2 = "s.t.txt";
dim2Vector info;
//获取文件信息
info = getFileINFO(File1, File2);
//解线性方程组
Cal_LP(info);
}
//解线性方程组
void Cal_LP(dim2Vector info)
{
int i, j, k;
double sum;
dim1Vector dim1Temp;
dim1Vector DEV; //偏差
dim1Vector C; //目标函数的系数向量
dim1Vector X; //初始基本解
dim2Vector A; //约束条件的系数向量
dim1Vector Col;
dim1Vector Row;
//将info中的数据进行拆分
for(i=0; i<info[0].size(); i++)
{
C.push_back(info[0][i]);
Col.push_back(i);
}
for(i=1; i<info.size()-2; i++)
{
dim1Temp.clear();
for(j=0; j<info[i].size(); j++)
dim1Temp.push_back(info[i][j]);
A.push_back(dim1Temp);
Row.push_back(INF);
}
for(j=0; j<info[i].size(); j++)
X.push_back(info[i][j]);
for(j=0; j<info[info.size()-1].size(); j++)
DEV.push_back(info[info.size()-1][j]);
//初始化偏差
for(i=0; i<A.size(); i++)
{
sum = 0;
for(j=0; j<A[i].size(); j++)
sum += A[i][j]*X[j];
DEV[i] = sum-DEV[i];
}
printf("初始表:\n");
for(i=0; i<C.size(); i++)
printf("%lf ", C[i]);
printf("\n");
for(i=0; i<A.size(); i++)
{
for(j=0; j<A[i].size(); j++)
printf("%lf ", A[i][j]);
printf("%lf\n", DEV[i]);
}
printf("初始解向量:\n");
for(i=0; i<X.size(); i++)
printf("X%d=%lf ", i+1, X[i]);
printf("\n");
printf("\n");
/**********************************************************************/
//开始进入计算
/*********************************************************************/
int count;
int laber_i, laber_j; //记录选择的系数的下标
double maxNum, minNum, tempNum;
dim2Vector Spin_Result;
int flag;
//预处理等式约束
for(i=0; i<equNum; i++)
{
flag = 0;
laber_i = i;
count = 0;
maxNum = -INF;
minNum = INF;
//偏差大于0
if(DEV[i]>0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]<0 && A[i][j]/C[j]>maxNum)
{
maxNum = A[i][j]/C[j];
laber_j = j;
count++;
}
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
}
//偏差小于0
if(DEV[i]<0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]>0 && A[i][j]/C[j]<minNum)
{
minNum = A[i][j]/C[j];
laber_j = j;
count++;
}
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
}
//偏差等于0
if(DEV[i]==0)
{
if(flag == 0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]>0)
{
flag = 1;
count++;
break;
}
}
if(flag == 0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]<0)
{
flag = 2;
count++;
break;
}
}
//情况1
if(flag==1)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]>0 && A[i][j]/C[j]<minNum)
{
minNum = A[i][j]/C[j];
laber_j = j;
}
}
//情况2
if(flag==2)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]<0 && A[i][j]/C[j]>maxNum)
{
maxNum = A[i][j]/C[j];
laber_j = j;
}
}
//情况3
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
}
tempNum = Row[laber_i];
Row[laber_i] = Col[laber_j];
Col[laber_j] = tempNum;
//旋转计算
Spin_Result = Cal_Spin(C, A, DEV, laber_i, laber_j);
for(k=0; k<Spin_Result.size()-2; k++)
for(j=0; j<Spin_Result[k].size(); j++)
A[k][j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k].size(); j++)
C[j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k+1].size(); j++)
DEV[j] = Spin_Result[k+1][j];
}
/***********************************************/
//迭代运算
/***********************************************/
int iter=0;
while(iter<MAXITER)
{
count = 0;
minNum = INF;
for(i=0; i<DEV.size(); i++)
if(DEV[i]<0)
{
sum = 0;
for(j=0; j<A[i].size(); j++)
sum += A[i][j]*A[i][j];
sum = DEV[i]/sqrtf(sum);
if(sum<minNum)
laber_i = i;
count++;
}
//所以偏差均为正数,结束迭代
if(count==0)
break;
count = 0;
minNum = INF;
for(j=0; j<A[laber_i].size(); j++)
{
if(A[laber_i][j]<0)
continue;
if(C[j]/A[laber_i][j]<minNum)
{
minNum = C[j]/A[laber_i][j];
laber_j = j;
}
count++;
}
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
tempNum = Row[laber_i];
Row[laber_i] = Col[laber_j];
Col[laber_j] = tempNum;
//旋转计算
Spin_Result = Cal_Spin(C, A, DEV, laber_i, laber_j);
for(k=0; k<Spin_Result.size()-2; k++)
for(j=0; j<Spin_Result[k].size(); j++)
A[k][j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k].size(); j++)
C[j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k+1].size(); j++)
DEV[j] = Spin_Result[k+1][j];
//计步器
iter++;
}
for(i=0; i<Row.size(); i++)
if(Row[i]!=INF)
{
X[Row[i]] = X[Row[i]]+DEV[i];
}
printf("\n迭代次数:%d\n", iter);
for(i=0; i<C.size(); i++)
printf("%lf ", C[i]);
printf("\n");
for(i=0; i<A.size(); i++)
{
for(j=0; j<A[i].size(); j++)
printf("%lf ", A[i][j]);
printf("%lf\n", DEV[i]);
}
printf("最优解为:\n");
for(i=0; i<X.size(); i++)
printf("X%d=%lf ", i+1, X[i]);
printf("\n");
sum = 0;
for(i=0; i<X.size(); i++)
sum += info[0][i]*X[i];
printf("最优目标函数值:%lf\n", sum);
}
//旋转计算
dim2Vector Cal_Spin(dim1Vector C, dim2Vector A, dim1Vector DEV, int laber_i, int laber_j)
{
int i, j;
dim1Vector temp_C;
dim2Vector temp_A;
dim1Vector temp_DEV;
dim2Vector dst;
temp_A = A;
temp_C = C;
temp_DEV = DEV;
//计算旋转向量C
C[laber_j] = temp_C[laber_j]/temp_A[laber_i][laber_j];
for(j=0; j<C.size(); j++)
{
if(j==laber_j)
continue;
C[j] = temp_C[j]-(temp_C[laber_j]/temp_A[laber_i][laber_j])*temp_A[laber_i][j];
}
//计算旋转向量DEV
DEV[laber_i] = -temp_DEV[laber_i]/temp_A[laber_i][laber_j];
for(i=0; i<DEV.size(); i++)
{
if(i==laber_i)
continue;
DEV[i] = temp_DEV[i]-(temp_A[i][laber_j]/temp_A[laber_i][laber_j])*temp_DEV[laber_i];
}
//计算旋转矩阵A
A[laber_i][laber_j] = 1/temp_A[laber_i][laber_j];
for(j=0; j<A[laber_i].size(); j++)
{
if(j==laber_j)
continue;
A[laber_i][j] = -temp_A[laber_i][j]/temp_A[laber_i][laber_j];
}
for(i=0; i<A.size(); i++)
{
if(i==laber_i)
continue;
A[i][laber_j] = temp_A[i][laber_j]/temp_A[laber_i][laber_j];
}
for(i=0; i<A.size(); i++)
for(j=0; j<A[i].size(); j++)
{
if(i==laber_i || j==laber_j)
continue;
A[i][j] = temp_A[i][j]-(temp_A[i][laber_j]/temp_A[laber_i][laber_j])*temp_A[laber_i][j];
}
dst = A;
dst.push_back(C);
dst.push_back(DEV);
return dst;
}
//获取文本信息
dim2Vector getFileINFO(char *File1, char *File2)
{
int i=1;
double num;
dim1Vector temp;
dim2Vector dst;
FILE *fp1 = fopen(File1, "r");
FILE *fp2 = fopen(File2, "r");
if(fp1==NULL || fp2==NULL)
{
printf("Open File Error!!\n");
exit(0);
}
while(fscanf(fp1, "%lf", &num)!=EOF)
{
temp.push_back(num);
if(i%varNum==0)
{
dst.push_back(temp);
temp.clear();
}
i++;
}
temp.clear();
while(fscanf(fp2, "%lf", &num)!=EOF)
{
temp.push_back(num);
}
dst.push_back(temp);
return dst;
}
对于如下的线性函数,求最优解:
初始表数据:
算法步骤:
一个实例:
c语言程序运行结果:
c语言代码如下:
#include "stdio.h"
#include "stdlib.h"
#include "math.h"
#include "vector"
using namespace std;
#define INF 9999
#define varNum 4 //变量数
#define equNum 0 //等式约束数
#define MAXITER 1000 //最大迭代次数
typedef vector<double> dim1Vector;
typedef vector<dim1Vector> dim2Vector;
dim2Vector getFileINFO(char *File1, char *File2); //获取文本信息
void Cal_LP(dim2Vector info); //解线性方程组
dim2Vector Cal_Spin(dim1Vector C, dim2Vector A, dim1Vector DEV, int laber_i, int laber_j); //旋转计算
void main()
{
char *File1 = "info.txt";
char *File2 = "s.t.txt";
dim2Vector info;
//获取文件信息
info = getFileINFO(File1, File2);
//解线性方程组
Cal_LP(info);
}
//解线性方程组
void Cal_LP(dim2Vector info)
{
int i, j, k;
double sum;
dim1Vector dim1Temp;
dim1Vector DEV; //偏差
dim1Vector C; //目标函数的系数向量
dim1Vector X; //初始基本解
dim2Vector A; //约束条件的系数向量
dim1Vector Col;
dim1Vector Row;
//将info中的数据进行拆分
for(i=0; i<info[0].size(); i++)
{
C.push_back(info[0][i]);
Col.push_back(i);
}
for(i=1; i<info.size()-2; i++)
{
dim1Temp.clear();
for(j=0; j<info[i].size(); j++)
dim1Temp.push_back(info[i][j]);
A.push_back(dim1Temp);
Row.push_back(INF);
}
for(j=0; j<info[i].size(); j++)
X.push_back(info[i][j]);
for(j=0; j<info[info.size()-1].size(); j++)
DEV.push_back(info[info.size()-1][j]);
//初始化偏差
for(i=0; i<A.size(); i++)
{
sum = 0;
for(j=0; j<A[i].size(); j++)
sum += A[i][j]*X[j];
DEV[i] = sum-DEV[i];
}
printf("初始表:\n");
for(i=0; i<C.size(); i++)
printf("%lf ", C[i]);
printf("\n");
for(i=0; i<A.size(); i++)
{
for(j=0; j<A[i].size(); j++)
printf("%lf ", A[i][j]);
printf("%lf\n", DEV[i]);
}
printf("初始解向量:\n");
for(i=0; i<X.size(); i++)
printf("X%d=%lf ", i+1, X[i]);
printf("\n");
printf("\n");
/**********************************************************************/
//开始进入计算
/*********************************************************************/
int count;
int laber_i, laber_j; //记录选择的系数的下标
double maxNum, minNum, tempNum;
dim2Vector Spin_Result;
int flag;
//预处理等式约束
for(i=0; i<equNum; i++)
{
flag = 0;
laber_i = i;
count = 0;
maxNum = -INF;
minNum = INF;
//偏差大于0
if(DEV[i]>0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]<0 && A[i][j]/C[j]>maxNum)
{
maxNum = A[i][j]/C[j];
laber_j = j;
count++;
}
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
}
//偏差小于0
if(DEV[i]<0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]>0 && A[i][j]/C[j]<minNum)
{
minNum = A[i][j]/C[j];
laber_j = j;
count++;
}
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
}
//偏差等于0
if(DEV[i]==0)
{
if(flag == 0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]>0)
{
flag = 1;
count++;
break;
}
}
if(flag == 0)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]<0)
{
flag = 2;
count++;
break;
}
}
//情况1
if(flag==1)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]>0 && A[i][j]/C[j]<minNum)
{
minNum = A[i][j]/C[j];
laber_j = j;
}
}
//情况2
if(flag==2)
{
for(j=0; j<A[i].size(); j++)
if(A[i][j]<0 && A[i][j]/C[j]>maxNum)
{
maxNum = A[i][j]/C[j];
laber_j = j;
}
}
//情况3
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
}
tempNum = Row[laber_i];
Row[laber_i] = Col[laber_j];
Col[laber_j] = tempNum;
//旋转计算
Spin_Result = Cal_Spin(C, A, DEV, laber_i, laber_j);
for(k=0; k<Spin_Result.size()-2; k++)
for(j=0; j<Spin_Result[k].size(); j++)
A[k][j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k].size(); j++)
C[j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k+1].size(); j++)
DEV[j] = Spin_Result[k+1][j];
}
/***********************************************/
//迭代运算
/***********************************************/
int iter=0;
while(iter<MAXITER)
{
count = 0;
minNum = INF;
for(i=0; i<DEV.size(); i++)
if(DEV[i]<0)
{
sum = 0;
for(j=0; j<A[i].size(); j++)
sum += A[i][j]*A[i][j];
sum = DEV[i]/sqrtf(sum);
if(sum<minNum)
laber_i = i;
count++;
}
//所以偏差均为正数,结束迭代
if(count==0)
break;
count = 0;
minNum = INF;
for(j=0; j<A[laber_i].size(); j++)
{
if(A[laber_i][j]<0)
continue;
if(C[j]/A[laber_i][j]<minNum)
{
minNum = C[j]/A[laber_i][j];
laber_j = j;
}
count++;
}
if(count==0)
{
printf("线性方程组无解!!!\n");
exit(0);
}
tempNum = Row[laber_i];
Row[laber_i] = Col[laber_j];
Col[laber_j] = tempNum;
//旋转计算
Spin_Result = Cal_Spin(C, A, DEV, laber_i, laber_j);
for(k=0; k<Spin_Result.size()-2; k++)
for(j=0; j<Spin_Result[k].size(); j++)
A[k][j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k].size(); j++)
C[j] = Spin_Result[k][j];
for(j=0; j<Spin_Result[k+1].size(); j++)
DEV[j] = Spin_Result[k+1][j];
//计步器
iter++;
}
for(i=0; i<Row.size(); i++)
if(Row[i]!=INF)
{
X[Row[i]] = X[Row[i]]+DEV[i];
}
printf("\n迭代次数:%d\n", iter);
for(i=0; i<C.size(); i++)
printf("%lf ", C[i]);
printf("\n");
for(i=0; i<A.size(); i++)
{
for(j=0; j<A[i].size(); j++)
printf("%lf ", A[i][j]);
printf("%lf\n", DEV[i]);
}
printf("最优解为:\n");
for(i=0; i<X.size(); i++)
printf("X%d=%lf ", i+1, X[i]);
printf("\n");
sum = 0;
for(i=0; i<X.size(); i++)
sum += info[0][i]*X[i];
printf("最优目标函数值:%lf\n", sum);
}
//旋转计算
dim2Vector Cal_Spin(dim1Vector C, dim2Vector A, dim1Vector DEV, int laber_i, int laber_j)
{
int i, j;
dim1Vector temp_C;
dim2Vector temp_A;
dim1Vector temp_DEV;
dim2Vector dst;
temp_A = A;
temp_C = C;
temp_DEV = DEV;
//计算旋转向量C
C[laber_j] = temp_C[laber_j]/temp_A[laber_i][laber_j];
for(j=0; j<C.size(); j++)
{
if(j==laber_j)
continue;
C[j] = temp_C[j]-(temp_C[laber_j]/temp_A[laber_i][laber_j])*temp_A[laber_i][j];
}
//计算旋转向量DEV
DEV[laber_i] = -temp_DEV[laber_i]/temp_A[laber_i][laber_j];
for(i=0; i<DEV.size(); i++)
{
if(i==laber_i)
continue;
DEV[i] = temp_DEV[i]-(temp_A[i][laber_j]/temp_A[laber_i][laber_j])*temp_DEV[laber_i];
}
//计算旋转矩阵A
A[laber_i][laber_j] = 1/temp_A[laber_i][laber_j];
for(j=0; j<A[laber_i].size(); j++)
{
if(j==laber_j)
continue;
A[laber_i][j] = -temp_A[laber_i][j]/temp_A[laber_i][laber_j];
}
for(i=0; i<A.size(); i++)
{
if(i==laber_i)
continue;
A[i][laber_j] = temp_A[i][laber_j]/temp_A[laber_i][laber_j];
}
for(i=0; i<A.size(); i++)
for(j=0; j<A[i].size(); j++)
{
if(i==laber_i || j==laber_j)
continue;
A[i][j] = temp_A[i][j]-(temp_A[i][laber_j]/temp_A[laber_i][laber_j])*temp_A[laber_i][j];
}
dst = A;
dst.push_back(C);
dst.push_back(DEV);
return dst;
}
//获取文本信息
dim2Vector getFileINFO(char *File1, char *File2)
{
int i=1;
double num;
dim1Vector temp;
dim2Vector dst;
FILE *fp1 = fopen(File1, "r");
FILE *fp2 = fopen(File2, "r");
if(fp1==NULL || fp2==NULL)
{
printf("Open File Error!!\n");
exit(0);
}
while(fscanf(fp1, "%lf", &num)!=EOF)
{
temp.push_back(num);
if(i%varNum==0)
{
dst.push_back(temp);
temp.clear();
}
i++;
}
temp.clear();
while(fscanf(fp2, "%lf", &num)!=EOF)
{
temp.push_back(num);
}
dst.push_back(temp);
return dst;
}
相关文章推荐
- 详解C++编程中类的成员变量和成员函数的相关知识
- C语言基础复习
- 时间函数
- C++学习笔记之define的多行定义
- 剖析C++的面向对象编程思想
- 多个空格用一个空格替代
- C++中获取数据类型typeid
- 可重入函数
- c++ 智能指针的设计与使用
- c++ STL总结
- boost 最大流求解示例
- C++ <二分查找法>
- 链表应用
- C基础-标准C语言头文件
- C++ Primer第四版 笔记
- C语言编程优化运行速度
- C++重载操作符之二(输入输出,以及乘法)
- C语言算法小练习-1
- 【LeetCode从零单刷】Binary Tree Postorder Traversal
- C++结构体中sizeof(1)