单纯形算法(simplex-algorithm)实现

kollektor

已于 2022-06-16 20:22:39 修改

阅读量430

点赞数

分类专栏：笔记文章标签：算法

于 2022-05-29 13:55:31 首次发布

本文链接：https://blog.csdn.net/weixin_46127065/article/details/125030186

版权

笔记专栏收录该内容

1 篇文章 0 订阅

订阅专栏

本算法可求解一般线性规划和整数线性规划问题

算法流程图

C语言实现

matrix.h

#ifndef MATRIX_H
#define MATRIX_H

/* matrix.h单纯形矩阵有关操作头文件 */

#include<stdio.h>
#include<stdlib.h>
#include<string.h>
#include<math.h>

//一字节八位
#define u8 unsigned char    //无符号字符型 1字节
#define i8 char             //字符型 1字节
#define u16 unsigned short int      //无符号短整型 2字节
#define i16 short int       //短整型 2字节
#define u32 unsigned int        //无符号整型 4字节
#define i32 int     //整型 4字节

typedef struct MATRIX {      //矩阵数据结构
	int col, row;
	float **mem;		
} MATRIX, *Matrix;


//创建row行col列的矩阵mat
bool Create_Matrix(Matrix mat, int row, int col);
//删除矩阵mat
void Delete_Matrix(Matrix mat);
//设置矩阵mat的数据
void Set_Matrix(Matrix mat, float* val);
//原矩阵src_mat复制到目标矩阵des_mat
void Copy_Matrix(Matrix des_mat, MATRIX src_mat);
//矩阵增加exp_num行
void Matrix_Expand_Row(Matrix mat, u32 exp_num, float* val);
//矩阵增加exp_num列
void Matrix_Expand_Column(Matrix mat, u32 exp_num, float* val);
//矩阵删除下标为row_index的行 下标从0开始计
void Matrix_Remove_Row(Matrix mat, int row_index);
//矩阵删除下标为col_index的列 下标从0开始计
void Matrix_Remove_Column(Matrix mat, int col_index);
//矩阵打印输出
void print_matrix(MATRIX mat);

#endif

matrix.c


/* matrix.c矩阵有关操作实现 */

#include"matrix.h"

u8 Create_Matrix(Matrix mat, int row, int col) {
    int i;
    mat->mem = (float**)malloc(row * sizeof(float*));
    for (i = 0;i < row;i++) {       //为mat的每一行分配col个元素空间
        mat->mem[i] = (float*)malloc(col * sizeof(float));
    }
    if (mat->mem = NULL) {
        printf("malloc fail\n");
        return 1;       
    }
    mat->row = row;     //设置行数列数
    mat->col = col;
    return 0;
}

void Delete_Matrix(Matrix mat) {
    int i;
    for (i = 0;i < mat->row;i++) {
        free(mat -> mem[i]);        //释放mat每一行的col个元素空间
    }
    free(mat->mem);     //释放mat的空间
}

void Set_Matrix(Matrix mat, float* val) {
    int row, col;
    for (row = 0;row < mat->row;row++) {
        for (col = 0;col < mat->col;col++) {
            mat->mem[row][col] = val[col + row * mat->col];         //矩阵数据保存在一维数组中，当前行标+当前列标*总列数即为映射关系
        }
    }
}

void Copy_Matrix(Matrix des_mat, MATRIX src_mat) {
    int row, col;
    if (des_mat->row != src_mat.row || des_mat->col != src_mat.col) {
        printf("err, matrix attribute unmatched\n");
        return;
    }
    for (row = 0;row < des_mat->row;row++) {        //copy
        for (col = 0;col < des_mat->col;col++) {
            des_mat->mem[row][col] = src_mat.mem[row][col];
        }
    }
}

void Matrix_Expand_Row(Matrix mat,u32 exp_num,float *val) {
    MATRIX temp;    //缓存
    int row, col;
    u32 count = 0;      //计数器

    Create_Matrix(&temp, mat->row, mat->col);       //建立缓存矩阵
    Copy_Matrix(&temp, *mat);                   //将原矩阵复制到缓存矩阵中 目的传引用 源传指针
    Delete_Matrix(mat);     //删除原矩阵 重建添加exp_num行后的新矩阵
    Create_Matrix(mat, temp.row + exp_num, temp.col);

    for (row = 0;row < temp.row;row++) {        //转移数据
        for (col = 0;col < temp.col;col++) {
            mat->mem[row][col] = temp.mem[row][col];
        }
    }

    for (row = temp.row;row < mat->row;row++) {         //添加数据
        for (col = temp.col;col < mat->col;col++) {
            mat->mem[row][col] = val[count++];
        }
    }

    Delete_Matrix(&temp);       //删除缓存矩阵
}

void Matrix_Expand_Column(Matrix mat,u32 exp_num,float *val) {
    MATRIX temp;
    int row, col;
    u32 count = 0;
    
    Create_Matrix(&temp, mat->row, mat->col);        //建立缓存矩阵
    Copy_Matrix(&temp, *mat);       //将原矩阵复制到缓存矩阵中 目的传引用 源传指针
    Delete_Matrix(mat);         //删除原矩阵 重建添加exp_num行后的新矩阵
    Create_Matrix(mat, temp.row, temp.col + exp_num);
    for (row = 0;row < temp.row;row++) {
        for (col = 0;col < temp.col;col++) {        //转移原数据
            mat->mem[row][col] = temp.mem[row][col];
        }
        for (;col < mat->col;col++) {       //继续在当前行中向新加入列添加数据
            mat->mem[row][col] = val[count++];
        }
    }

    Delete_Matrix(&temp);
}

void Matrix_Remove_Column(Matrix mat, int col_index) {
    MATRIX temp;
    int row, col;
    u32 count;

    if (col_index > mat->col - 1) {     //试图删除不存在的列
        printf("err,matrix_remove_column\n");
        return;
    }

    Create_Matrix(&temp, mat->row, mat->col);
    Copy_Matrix(&temp, *mat);
    Delete_Matrix(mat);     //重构列数减少1的矩阵
    Create_Matrix(mat, temp.row, temp.col - 1);

    for (row = 0;row < temp.row;row++) {        //转移原矩阵数据到新矩阵 注意跳过下标为col_index的列
        count = 0;
        for (col = 0;col < temp.col;col++) {
            if (col == col_index) {         //跳过要删除的列
                continue;
            }
            mat->mem[row][count] = temp.mem[row][col];      //count控制原矩阵列标 配合跳过要删除的列
            count++;
        }
    }

    Delete_Matrix(&temp);
}

void Matrix_Remove_Row(Matrix mat, int row_index) {
    MATRIX temp;
    int row, col;
    u32 count;

    if (row_index > mat->row - 1) {
        printf("err,matrix_remove_row\n");
        return;
    }

    Create_Matrix(&temp, mat->row, mat->col);
    Copy_Matrix(&temp, *mat);
    Delete_Matrix(mat);
    Create_Matrix(mat, temp.row - 1, temp.col);

    for (col = 0;col < temp.col;col++) {
        count = 0;
        for (row = 0;row < temp.row;row++) {
            if (row == row_index) {
                continue;
            }
            mat->mem[count][col] = temp.mem[row][col];
            count++;
        }
    }

    Delete_Matrix(&temp);
}

simplex.h

#ifndef SIMPLEX_H
#define SIMPLEX_H

/* simplex.h单纯形算法数据结构与操作头文件 */

#include"matrix.h"

typedef struct PIVOT {      //基变量数据结构
    int bvar;       //base variant variable:yi
    int nbvar;      //non-base variant variable:xi
}PIVOT;

typedef enum STRATEGY {     //LP=Linear Programming,ILP=Integer Linear Programming
    LP = 0,
    ILP = 1
}STRATEGY;

u8 initSimplexModel(MATRIX mat, MATRIX g_mat, char* var_name[]);
void deleteSimplexModel(void);
u8 SimplexRun(STRATEGY strategy);
int getMaxIntValue(void);
void getIntSolution(int* solution);

#endif

simplex.c


/* simplex.c单纯形算法操作实现 */

#include"simplex.h"

#define MAX_VAL_NUM 6144        //最大变量个数
#define MAX_VAL_NAME 8      //最大变量名个数

static MATRIX temp_mat;
static MATRIX goal_mat;
static char VarName[MAX_VAL_NUM][MAX_VAL_NAME];
static u32 VarName_Map[MAX_VAL_NUM] = { 0 };
static int BVarIff;
static float* coefbuff;
static u8 expand_flag;
static u8 using_name;
static int temp_z_col;

float accuracy = 0.01f;
#define eq(a,b)	((a-b)<accuracy&&(a-b)>-accuracy)
#define lt(a,b)	((b-a)>accuracy)
#define	le(a,b)	((b-a)>accuracy||((b-a)>-accuracy&&(b-a)<accuracy))
#define isInt(a) (fabs((fabs(a)-(int)(fabs(a)+0.5)))<accuracy)
#define toInt(a) (a>0?(int)(a+0.5f):(int)(a-0.5f))

#define RESERVE_SPACE 3
#define SHOW_RESULT 1

#define RESERVE_SPACE		3
#define SHOW_RESULT			1

void showResult(MATRIX sf_mat, u8 goal) {
	int row, col;

	if (using_name) {
		for (row = 0; row < sf_mat.row; row++) {
			if (goal)
				printf("\nU = ");
			else
				printf("\n%s = ", VarName[VarName_Map[row + BVarOff]]);
			for (col = 0; col < sf_mat.col; col++) {
				printf("%.2f%s	", sf_mat.mem[row][col], VarName[VarName_Map[col]]);
			}
		}
		printf("\n");
	} else {
		for (row = 0; row < sf_mat.row; row++) {
			if (goal)
				printf("\nU = ");
			for (col = 0; col < sf_mat.col; col++) {
				printf("%.2f	", sf_mat.mem[row][col]);
			}
		}
		printf("\n");
	}
}

/************************************************************************
	forall b_i>=0, return 0,   otherwise: return 1
************************************************************************/
u8 should_expand(MATRIX mat) {
	int row;

	for (row = 0; row < mat.row; row++) {
		/* the first column is b_i */
		if (mat.mem[row][0] < 0)
			return 1;
	}

	return 0;
}

u8 addNonBasicVarName(MATRIX sf_mat, char *name) {
	if (sf_mat.col >= BVarOff) {
		printf("err, non-basic var reserve space is not enough\n");
		return 1;
	}

	if (using_name)
		strcpy(VarName[sf_mat.col], name);
	VarName_Map[sf_mat.col] = sf_mat.col;

	return 0;
}

void resetSimplexModel(MATRIX *sf_mat) {
	int i;

	//expand_flag = 0;

	if (should_expand(*sf_mat)) {
		float *val;

		if (!expand_flag) {
			//not expand before, we need allocate name space for z
			addNonBasicVarName(*sf_mat, "z");
		} else {
			//VarName_Map[temp_z_col] = sf_mat->col-1;

			for (i = temp_z_col; i < sf_mat->col; i++) {
				VarName_Map[i] = VarName_Map[i + 1];
			}
			VarName_Map[i] = sf_mat->col;
		}

		val = (float *)malloc(sf_mat->row * sizeof(float));
		for (i = 0; i < sf_mat->row; i++)
			val[i] = 1.0f;	/* each row add a fresh variable z, for z>=0 */

		matrix_expand_column(sf_mat, 1, val);
		free(val);

		expand_flag = 1;
	} else {
		expand_flag = 0;
	}

	free(coefbuff);
	coefbuff = (float *)malloc(sf_mat->row * sizeof(float));
}

/************************************************************************
Set goal form, calculate the maximum value of goal form.
The last column of mat is constant
************************************************************************/
void setGoalForm(MATRIX mat) {
	copy_matrix(&goal_mat, mat);
	if (expand_flag) {
		/* add -z for goal form */
		float val = -1.0f;
		//float val = 0.0f;
		matrix_expand_column(&goal_mat, 1, &val);
	}
}

/************************************************************************
mat: A * x_i <= b_i (x_i >=0)
the first column of matrix stores constant b_i

Slack Form: b_i - A * x_i >= 0
************************************************************************/
void createSlackForm(MATRIX *mat) {
	int row, col;

	for (row = 0; row < mat->row; row++)
		for (col = 1; col < mat->col; col++)
			mat->mem[row][col] = -mat->mem[row][col];
}

/************************************************************************
Find a non-basic variable, which has a positive coefficient
return: -1:doesn't have positive coefficient  otherwise: index of coefficient
************************************************************************/
int findNextNonBasicVar(void) {
	int i;

	/* the first column is constant */
	for (i = 1; i < goal_mat.col; i++) {
		/* ignore z */
		if (expand_flag && VarName_Map[i] == goal_mat.col - 1)
			continue;
		if (goal_mat.mem[0][i] > accuracy)
			return i;
	}

	return -1;
}

int findNextNonBasicVarWithZ(MATRIX sf_mat, int row) {
	int col;

	/* we calculate the min(z) here, so find negative coefficient */
	for (col = 1; col < sf_mat.col; col++) {
		//TODO: need goal_mat.mem[0][col]>0?
		//if( lt(sf_mat.mem[row][col], 0.0f) && lt(0.0f, goal_mat.mem[0][col]) ){
		if (lt(sf_mat.mem[row][col], 0.0f)) {
			return col;
		}
	}

	return -1;
}

/************************************************************************
	 find most negative b_i to be pivot with z
************************************************************************/
PIVOT findMostNegativePivot(MATRIX sf_mat) {
	int row;
	float neg_b;
	PIVOT pivot;

	pivot.bvar = -1;
	pivot.nbvar = sf_mat.col - 1;	 /* swap z and basic variable */

	neg_b = -accuracy;
	for (row = 0; row < sf_mat.row; row++) {
		if (sf_mat.mem[row][0] < neg_b) {
			neg_b = sf_mat.mem[row][0];
			pivot.bvar = row;
		}
	}

	return pivot;
}

PIVOT findPivot(MATRIX sf_mat, int nbvar) {
	int row;
	float tight_constrain = 100000000.0f;
	float temp_constrain;
	PIVOT pivot;

	pivot.bvar = -1;	/* -1 means doesn't find a valid pivot */
	pivot.nbvar = nbvar;

	for (row = 0; row < sf_mat.row; row++) {
		if (eq(0.0f, sf_mat.mem[row][nbvar]) || sf_mat.mem[row][nbvar] > 0.0f)
			continue;
		temp_constrain = sf_mat.mem[row][0] / -sf_mat.mem[row][nbvar];
		if (temp_constrain < tight_constrain) {
			tight_constrain = temp_constrain;
			pivot.bvar = row;
		}
	}

	return pivot;
}

void swapNameSpace(PIVOT pivot) {
	int map_index;

	map_index = VarName_Map[pivot.nbvar];
	VarName_Map[pivot.nbvar] = VarName_Map[pivot.bvar + BVarOff];
	VarName_Map[pivot.bvar + BVarOff] = map_index;
}

/************************************************************************
		swap basic var and non-basic var
************************************************************************/
void swapPivot(MATRIX *sf_mat, PIVOT pivot) {
	int col;
	float temp;

	temp = sf_mat->mem[pivot.bvar][pivot.nbvar];
	for (col = 0; col < sf_mat->col; col++) {
		if (col == pivot.nbvar)
			sf_mat->mem[pivot.bvar][col] = 1.0f / temp;
		else
			sf_mat->mem[pivot.bvar][col] /= -temp;
	}
}

void replaceNewBasicVar(MATRIX *sf_mat, PIVOT pivot) {
	int row, col;

	for (row = 0; row < sf_mat->row; row++) {
		if (row == pivot.bvar)
			continue;

		coefbuff[row] = sf_mat->mem[row][pivot.nbvar];

		for (col = 0; col < sf_mat->col; col++) {
			if (col == pivot.nbvar)
				sf_mat->mem[row][col] = coefbuff[row] * sf_mat->mem[pivot.bvar][col];
			else
				sf_mat->mem[row][col] += coefbuff[row] * sf_mat->mem[pivot.bvar][col];
		}
	}
}

void updateGoalForm(MATRIX *gf_mat, MATRIX *sf_mat, PIVOT pivot) {
	int col;
	float temp;

	temp = gf_mat->mem[0][pivot.nbvar];

	for (col = 0; col < gf_mat->col; col++) {
		if (col == pivot.nbvar)
			gf_mat->mem[0][col] = temp * sf_mat->mem[pivot.bvar][col];
		else
			gf_mat->mem[0][col] += temp * sf_mat->mem[pivot.bvar][col];
	}
}

/************************************************************************
swap pivot and reformulate slack form
************************************************************************/
void reformulation(MATRIX *sf_mat, PIVOT pivot) {
	/* step1: swap basic var and non-basic var */
	swapNameSpace(pivot);
	swapPivot(sf_mat, pivot);
	/* step2: replace new basic var with non-basic var */
	replaceNewBasicVar(sf_mat, pivot);
	/* step3: update goal form */
	updateGoalForm(&goal_mat, sf_mat, pivot);
}

void dispSolution(MATRIX sf_mat, u8 show_res) {
	int row, col;

	if (show_res) {
		showResult(goal_mat, 1);
		showResult(temp_mat, 0);
		printf("\nSolution:\n");
	}
	for (col = 1; col < sf_mat.col; col++) {
		float solution = 0.0f;

		for (row = 0; row < sf_mat.row; row++) {
			if (VarName_Map[row + BVarOff] == col) {
				solution = sf_mat.mem[row][0];
			}
		}
		if (show_res) {
			if (using_name)
				printf("%s=%.2f, ", VarName[col], solution);
			else
				printf("%.2f, ", solution);
		}
	}
	if (show_res)
		printf("U=%.2f\n", goal_mat.mem[0][0]);
}

void dispIntSolution(MATRIX sf_mat, u8 show_res) {
	int row, col;

	if (show_res) {
		showResult(goal_mat, 1);
		showResult(temp_mat, 0);
		printf("\nInterger Solution:\n");
	}

	for (col = 1; col < temp_mat.col - expand_flag; col++) {
		int sol = 0;

		for (row = 0; row < temp_mat.row; row++) {
			if (VarName_Map[row + BVarOff] == col) {
				sol = toInt(temp_mat.mem[row][0]);
			}
		}
		if (show_res) {
			if (using_name)
				printf("%s=%d, ", VarName[col], sol);
			else
				printf("%d, ", sol);
		}
	}

	if (show_res)
		printf("U=%.2f\n", goal_mat.mem[0][0]);
}

int getNonIntIndex(MATRIX sf_mat, int row_index) {
	int row, col;

	if (row_index >= sf_mat.row) {
		return -1;
	}

	for (row = (int)row_index; row < sf_mat.row; row++) {
		if (!isInt(sf_mat.mem[row][0])) {
			for (col = 1; col < sf_mat.col; col++) {
				if (!isInt(sf_mat.mem[row][col]))
					return row;	/* find a basic var which has non-interger bi and xi */
			}
		}
	}

	return -100;
}

/* add gomory constraint to calculate interger solution */
void addGomoryConstraint(MATRIX *s_mat, int row_index) {
	static int cons_cnt = 1;
	float *gomory_cons;
	int col;

	gomory_cons = (float *)malloc(s_mat->col * sizeof(float));

	/* extract basic var fractional part, with b>=0 */
	gomory_cons[0] = -(s_mat->mem[row_index][0] - (int)s_mat->mem[row_index][0]);

	/* extract non-basic var fractional part */
	for (col = 1; col < s_mat->col; col++) {
		if (eq(s_mat->mem[row_index][col], 0.0f) || isInt(s_mat->mem[row_index][col])) {
			/* no fraction part */
			gomory_cons[col] = 0.0f;
			continue;
		}
		gomory_cons[col] = -(s_mat->mem[row_index][col] - (float)ceil(s_mat->mem[row_index][col]));
	}

	/* add gomory constraint */
	matrix_expand_row(s_mat, 1, gomory_cons);

	if (using_name)
		sprintf(VarName[s_mat->row - 1 + BVarOff], "g%d", cons_cnt++);
	VarName_Map[s_mat->row - 1 + BVarOff] = s_mat->row - 1 + BVarOff;

	free(gomory_cons);
}

/* remove z column */
int removeZ(MATRIX *s_mat) {
	int z_col = -1;
	if (!expand_flag)
		return -1;

	for (z_col = 0; z_col < s_mat->col; z_col++) {
		if (VarName_Map[z_col] == s_mat->col - 1) {
			break;
		}
	}

	if (z_col >= 0) {
		matrix_remove_column(s_mat, (int)z_col);
		matrix_remove_column(&goal_mat, (int)z_col);
	}

	return z_col;
}

u8 solveSimplex(void) {
	u8 res = 0;
	int next_nbvar;
	PIVOT pivot;

	/* we should first calculate the min(z), if min(z)=0, then change z to non-basic var */
	if (expand_flag) {
		int row_z;

		pivot = findMostNegativePivot(temp_mat);
		row_z = (int)pivot.bvar;

		if (pivot.bvar == -1) {
			//printf("NA, findMostNegativePivot\n");
			return 0;
		}
		reformulation(&temp_mat, pivot);

		while (1) {
			if (eq(temp_mat.mem[row_z][0], 0.0f)) {
				/* min(z)=0 */
				int col;
				u8 res = 0;
				for (col = 1; col < temp_mat.col; col++) {
					if (!eq(temp_mat.mem[row_z][col], 0.0f)) {
						pivot.bvar = row_z;
						pivot.nbvar = col;
						reformulation(&temp_mat, pivot);
						res = 1;
						break;
					}
				}
				if (res)
					break;
				else {
					printf("some err here\n");
					while (1);
				}

			}

			next_nbvar = findNextNonBasicVarWithZ(temp_mat, row_z);
			if (next_nbvar == -1) {
				//printf("NA, findNextNonBasicVarForZ\n");
				//showResult(goal_mat, 1);
				//showResult(temp_mat, 0);
				return 0;
			}
			pivot = findPivot(temp_mat, (int)next_nbvar);
			reformulation(&temp_mat, pivot);
			if (pivot.bvar == row_z) {
				//TODO: is it correct?
				/* z becomes non-basic variable now, we can continue */
				break;
			}
		}
	}

	while (1) {
		next_nbvar = findNextNonBasicVar();
		if (next_nbvar < 0) {
			/* non positive coefficient, we get the maximum value */
			res = 1;
			break;
		}

		pivot = findPivot(temp_mat, (int)next_nbvar);
		if (pivot.bvar == -1) {
			/* goal formula can be infinite, so doesn't have maximum value */
			//printf("NA, inf result\n");
			res = 0;
			break;
		} else {
			next_nbvar = 1;
		}
		reformulation(&temp_mat, pivot);
	}

	return res;
}

u8 incise(void) {
	int nonIntIndex = -1;
	u8 res;

	while (1) {
		//nonIntIndex = getNonIntIndex(temp_mat, nonIntIndex+1);
		nonIntIndex = getNonIntIndex(temp_mat, 0);
		/* Handle ILP */
		if (nonIntIndex >= 0) {
			temp_z_col = removeZ(&temp_mat);
			addGomoryConstraint(&temp_mat, (int)nonIntIndex);
			resetSimplexModel(&temp_mat);
			setGoalForm(goal_mat);

			res = solveSimplex();
			if (!res)	//NA
				return res;
		} else if (nonIntIndex == -1) {
			//TODO: is it possible for this situation?
			/* used up all gomory constraints and can't find interger solution */
			//printf("NA, can't find interger solution\n");
			return 0;
		} else {
			/* Get interger soluiton */
			return 1;
		}
	}
}

int getMaxIntValue(void) {
	return toInt(goal_mat.mem[0][0]);
}

void getIntSolution(int *solution) {
	int row, col;

	for (col = 1; col < temp_mat.col - expand_flag; col++) {
		float sol = 0.0f;

		for (row = 0; row < temp_mat.row; row++) {
			if (VarName_Map[row + BVarOff] == col) {
				sol = temp_mat.mem[row][0];
			}
		}
		solution[col - 1] = toInt(sol);
	}
}

u8 initSimplexModel(MATRIX mat, MATRIX g_mat, char *var_name[]) {
	u8 res = 0;
	int i;

	if (var_name == NULL)
		using_name = 0;
	else
		using_name = 1;
	expand_flag = 0;

	res |= create_matrix(&goal_mat, 1, mat.col);
	/* we will do calculation on temp_mat later */
	res |= create_matrix(&temp_mat, mat.row, mat.col);
	copy_matrix(&temp_mat, mat);

	if (should_expand(temp_mat)) {
		float *val;
		int i;
		expand_flag = 1;

		val = (float *)malloc(temp_mat.row * sizeof(float));
		for (i = 0; i < temp_mat.row; i++)
			val[i] = -1.0f;	/* each row add a fresh variable z, for z>=0 */

		matrix_expand_column(&temp_mat, 1, val);
		free(val);
	}

	coefbuff = (float *)malloc(mat.row * sizeof(float));

	for (i = 0; i < mat.col; i++) {
		if (using_name)
			strcpy(VarName[i], var_name[i]);
		VarName_Map[i] = i;
	}
	if (expand_flag) {
		if (using_name)
			strcpy(VarName[mat.col], "z");
		VarName_Map[mat.col] = mat.col;
		BVarOff = mat.col + RESERVE_SPACE;
	} else {
		BVarOff = mat.col + RESERVE_SPACE;	/* basic value offset of VarName array */
	}
	for (i = 0; i < mat.row; i++) {
		if (using_name)
			sprintf(VarName[i + BVarOff], "y%d", i + 1);
		VarName_Map[i + BVarOff] = i + BVarOff;
	}

	setGoalForm(g_mat);
	createSlackForm(&temp_mat);

	return res;
}

void deleteSimplexModel(void) {
	delete_matrix(&goal_mat);
	delete_matrix(&temp_mat);
	free(coefbuff);
}

/***************************************************************************
Input:
	c_mat: constraint matrix
	g_mat: goal matrix
	var_name: variable name
	strategy:
				LP:	Linear Program
			   ILP:	Interger Linear Program
Return:
	0: unsatisfy	otherwise: satisfy
***************************************************************************/
u8 SimplexRun(STRATEGY strategy) {
	u8 res = 0;

	res = solveSimplex();

	if (res) {
		if (strategy == LP) {
			dispSolution(temp_mat, SHOW_RESULT);
		} else {
			res = incise();
			dispIntSolution(temp_mat, SHOW_RESULT);
		}
	}

	return res;
}

main.c

/* main.c */

#include "matrix.h"
#include "simplex.h"

MATRIX mat;
MATRIX goal_mat;

void Test1(STRATEGY strategy) {
	char* var_name[] = { "","x1","x2", "x3", "x4" };
	//	 x1 +     + 2x3       <= 18
	//      + 2x2       - 7x4 <= 0
	//   x1 +  x2 +  x3 +  x4 <= 9			[ x1+x2=b <=> (x1+x2>=b && x1+x2<=b) ]
	// - x1 -  x2 -  x3 -  x4 <= -9
	//        -x2 +  x3 - 2x4 <= -1
	float val[] =
	{
		18,  1,  0,  2,  0,
		 0,  0,  2,  0, -7,
		 9,  1,  1,  1,  1,
		-9, -1, -1, -1, -1,
		-1,  0,-1,  1, -2,
	};
	//	maximize U = x1 + x2 + 3x3 - x4
	float goal[] = { 0, 1, 1, 3, -1 };
	int max;
	int solution[4];
	int i;

	create_matrix(&mat, 5, 5);
	set_matrix(&mat, val);

	create_matrix(&goal_mat, 1, 5);
	set_matrix(&goal_mat, goal);

	initSimplexModel(mat, goal_mat, var_name);
	SimplexRun(strategy);
	max = getMaxIntValue();
	getIntSolution(solution);
	deleteSimplexModel();

	delete_matrix(&mat);
	delete_matrix(&goal_mat);
}

void Test2(STRATEGY strategy) {
	u8 i;
	char* var_name[] = { "", "x1", "x2" };
	//  9x1  +  7x2  <=   56
	// 7x1  +  20x2  <=  70
	float val[] =
	{
		56,   9,  7,
	   70,  7, 20,
	};
	//	maximize U = 40x1 + 90x2
	float goal[] = { 0, 40, 90 };
	int max;
	int solution[2];

	create_matrix(&mat, 2, 3);
	set_matrix(&mat, val);

	create_matrix(&goal_mat, 1, 3);
	set_matrix(&goal_mat, goal);

	initSimplexModel(mat, goal_mat, var_name);
	SimplexRun(strategy);
	max = getMaxIntValue();
	getIntSolution(solution);
	deleteSimplexModel();

	delete_matrix(&mat);
	delete_matrix(&goal_mat);
}

void initSolver(void);
void solverTest(void);

int main() {
	Test1(LP);
	Test2(ILP);

	printf("\nPress any key to exit");
	getchar();
	return 0;
}