Java实现算法导论中线性规划单纯形算法

需在理解算法导论中线性规划单纯性算法基础上理解Java实现的代码，结合http://blog.csdn.net/fjssharpsword/article/details/53195556理解。

具体代码如下：

package sk.mlib;import java.util.Random;/**************************************************************************  Compilation:  javac Simplex.java*  Execution:    java Simplex**  Given an M-by-N matrix A, an M-length vector b, and an*  N-length vector c, solve the  LP { max cx : Ax <= b, x >= 0 }.*  Assumes that b >= 0 so that x = 0 is a basic feasible solution.**  Creates an (M+1)-by-(N+M+1) simplex tableaux with the *  RHS in column M+N, the objective function in row M, and*  slack variables in columns M through M+N-1.**************************************************************************/public class Simplex {private static final double EPSILON = 1.0E-10;private double[][] a;   // tableauxprivate int M;          // number of constraintsprivate int N;          // number of original variablesprivate int[] basis;    // basis[i] = basic variable corresponding to row i// only needed to print out solution, not book// sets up the simplex tableauxpublic Simplex(double[][] A, double[] b, double[] c) {M = b.length;N = c.length;a = new double[M+1][N+M+1];for (int i = 0; i < M; i++)for (int j = 0; j < N; j++)a[i][j] = A[i][j];for (int i = 0; i < M; i++) a[i][N+i] = 1.0;for (int j = 0; j < N; j++) a[M][j]   = c[j];for (int i = 0; i < M; i++) a[i][M+N] = b[i];basis = new int[M];for (int i = 0; i < M; i++) basis[i] = N + i;solve();// check optimality conditionsassert check(A, b, c);//初始化判断输入的线性规划是否存在初始基本可行解、对偶最优解判断、无限解}// run simplex algorithm starting from initial BFSprivate void solve() {while (true) {// find entering column qint q = bland();if (q == -1) break;  // optimal// find leaving row pint p = minRatioRule(q);if (p == -1) throw new RuntimeException("Linear program is unbounded");// pivotpivot(p, q);// update basisbasis[p] = q;}}// lowest index of a non-basic column with a positive costprivate int bland() {for (int j = 0; j < M + N; j++)if (a[M][j] > 0) return j;return -1;  // optimal}// index of a non-basic column with most positive costprivate int dantzig() {int q = 0;for (int j = 1; j < M + N; j++)if (a[M][j] > a[M][q]) q = j;if (a[M][q] <= 0) return -1;  // optimalelse return q;}// find row p using min ratio rule (-1 if no such row)private int minRatioRule(int q) {int p = -1;for (int i = 0; i < M; i++) {if (a[i][q] <= 0) continue;else if (p == -1) p = i;else if ((a[i][M+N] / a[i][q]) < (a[p][M+N] / a[p][q])) p = i;}return p;}// pivot on entry (p, q) using Gauss-Jordan elimination，主元操作，交换变量private void pivot(int p, int q) {// everything but row p and column qfor (int i = 0; i <= M; i++)for (int j = 0; j <= M + N; j++)if (i != p && j != q) a[i][j] -= a[p][j] * a[i][q] / a[p][q];// zero out column qfor (int i = 0; i <= M; i++)if (i != p) a[i][q] = 0.0;// scale row pfor (int j = 0; j <= M + N; j++)if (j != q) a[p][j] /= a[p][q];a[p][q] = 1.0;}// return optimal objective valuepublic double value() {return -a[M][M+N];}// return primal solution vectorpublic double[] primal() {double[] x = new double[N];for (int i = 0; i < M; i++)if (basis[i] < N) x[basis[i]] = a[i][M+N];return x;}// return dual solution vectorpublic double[] dual() {double[] y = new double[M];for (int i = 0; i < M; i++)y[i] = -a[M][N+i];return y;}// is the solution primal feasible?构造辅助线性规划private boolean isPrimalFeasible(double[][] A, double[] b) {double[] x = primal();// check that x >= 0for (int j = 0; j < x.length; j++) {if (x[j] < 0.0) {System.out.println("x[" + j + "] = " + x[j] + " is negative");return false;}}// check that Ax <= bfor (int i = 0; i < M; i++) {double sum = 0.0;for (int j = 0; j < N; j++) {sum += A[i][j] * x[j];}if (sum > b[i] + EPSILON) {System.out.println("not primal feasible");System.out.println("b[" + i + "] = " + b[i] + ", sum = " + sum);return false;}}return true;}// is the solution dual feasible?构造对偶线性规划private boolean isDualFeasible(double[][] A, double[] c) {double[] y = dual();// check that y >= 0for (int i = 0; i < y.length; i++) {if (y[i] < 0.0) {System.out.println("y[" + i + "] = " + y[i] + " is negative");return false;}}// check that yA >= cfor (int j = 0; j < N; j++) {double sum = 0.0;for (int i = 0; i < M; i++) {sum += A[i][j] * y[i];}if (sum < c[j] - EPSILON) {System.out.println("not dual feasible");System.out.println("c[" + j + "] = " + c[j] + ", sum = " + sum);return false;}}return true;}// check that optimal value = cx = ybprivate boolean isOptimal(double[] b, double[] c) {double[] x = primal();double[] y = dual();double value = value();// check that value = cx = ybdouble value1 = 0.0;for (int j = 0; j < x.length; j++)value1 += c[j] * x[j];double value2 = 0.0;for (int i = 0; i < y.length; i++)value2 += y[i] * b[i];if (Math.abs(value - value1) > EPSILON || Math.abs(value - value2) > EPSILON) {System.out.println("value = " + value + ", cx = " + value1 + ", yb = " + value2);return false;}return true;}private boolean check(double[][]A, double[] b, double[] c) {return isPrimalFeasible(A, b) && isDualFeasible(A, c) && isOptimal(b, c);}// print tableauxpublic void show() {System.out.println("M = " + M);System.out.println("N = " + N);for (int i = 0; i <= M; i++) {for (int j = 0; j <= M + N; j++) {System.out.printf("%7.2f ", a[i][j]);}System.out.println();}System.out.println("value = " + value());for (int i = 0; i < M; i++)if (basis[i] < N) System.out.println("x_" + basis[i] + " = " + a[i][M+N]);System.out.println();}public static void test(double[][] A, double[] b, double[] c) {Simplex lp = new Simplex(A, b, c);System.out.println("value = " + lp.value());double[] x = lp.primal();for (int i = 0; i < x.length; i++)System.out.println("x[" + i + "] = " + x[i]);double[] y = lp.dual();for (int j = 0; j < y.length; j++)System.out.println("y[" + j + "] = " + y[j]);}public static void test1() {double[][] A = {{ -1,  1,  0 },{  1,  4,  0 },{  2,  1,  0 },{  3, -4,  0 },{  0,  0,  1 },};double[] c = { 1, 1, 1 };double[] b = { 5, 45, 27, 24, 4 };test(A, b, c);}// x0 = 12, x1 = 28, opt = 800public static void test2() {double[] c = {  13.0,  23.0 };double[] b = { 480.0, 160.0, 1190.0 };double[][] A = {{  5.0, 15.0 },{  4.0,  4.0 },{ 35.0, 20.0 },};test(A, b, c);}// unboundedpublic static void test3() {double[] c = { 2.0, 3.0, -1.0, -12.0 };double[] b = {  3.0,   2.0 };double[][] A = {{ -2.0, -9.0,  1.0,  9.0 },{  1.0,  1.0, -1.0, -2.0 },};test(A, b, c);}// degenerate - cycles if you choose most positive objective function coefficientpublic static void test4() {double[] c = { 10.0, -57.0, -9.0, -24.0 };double[] b = {  0.0,   0.0,  1.0 };double[][] A = {{ 0.5, -5.5, -2.5, 9.0 },{ 0.5, -1.5, -0.5, 1.0 },{ 1.0,  0.0,  0.0, 0.0 },};test(A, b, c);}// test clientpublic static void main(String[] args) {try                 { test1();             }catch (Exception e) { e.printStackTrace(); }System.out.println("--------------------------------");/*try                 { test2();             }catch (Exception e) { e.printStackTrace(); }System.out.println("--------------------------------");try                 { test3();             }catch (Exception e) { e.printStackTrace(); }System.out.println("--------------------------------");try                 { test4();             }catch (Exception e) { e.printStackTrace(); }System.out.println("--------------------------------");int M = Integer.parseInt(args[0]);int N = Integer.parseInt(args[1]);double[] c = new double[N];double[] b = new double[M];double[][] A = new double[M][N];for (int j = 0; j < N; j++)c[j] = new Random().nextInt(1000);for (int i = 0; i < M; i++)b[i] = new Random().nextInt(1000);for (int i = 0; i < M; i++)for (int j = 0; j < N; j++)A[i][j] = new Random().nextInt(100);Simplex lp = new Simplex(A, b, c);System.out.println(lp.value());*/}}

test1测试数据执行结果如下：

value = 22.0
x[0] = 9.000000000000002
x[1] = 9.0
x[2] = 4.0
y[0] = -0.0
y[1] = 0.1428571428571428
y[2] = 0.4285714285714286
y[3] = -0.0
y[4] = 1.0
--------------------------------

Java实现算法导论中线性规划单纯形算法相关推荐

机器学习的算法和普通《算法导论》里的算法有什么本质上的异同
机器学习的算法和普通<算法导论>里的算法有什么本质上的异同? 本人非计算机专业出身,对这些方向感兴趣,所以有此一问.曾经问过一些人,说是机器学习全是数学,是用数学的方式试图去描述和理解我们 ...
算法导论NO.1：算法基础
目录算法导论NO.1:算法基础前言一.算法是什么? 二.算法的特征 1.有穷性 2.确切性 3.输入 4.输出 5.可行性三.算法的评价总结前言新坑今天来填第一天,最近返校后培训数学建模 ...
【算法导论】内部排序算法总结
排序名称时间复杂度空间复杂度稳定性直接插入排序 O(n^2) O(1) 稳定折半插入排序 O(n^2) O(1) 稳定希尔排序 O(n^2) O(1) 不稳定冒泡排序 O(n^2) O( ...
路径规划;a*算法 demo_路径规划A*算法
(做封面图片,可真是费了很长时间) 前言路径规划涉及到很多内容. 但路径规划的起点,应该都是A*这个来自斯坦福大学算法. 说起A*,网上有很多相关的内容.各种博客,各种解释. 但是我愣是没看明白. ...
【算法导论】地图染色算法
地图染色问题可以根据四色定理来解决.所谓四色定理,就是指可以用不多于四种的颜色对地图着色,使相邻的行政区域不重色,因此我们可以用四色定理的结论,用回溯算法对一幅给定的地图染色. 算法的 ...
算法导论 — 8.1 排序算法的下界
笔记我们所熟知的插入排序.归并排序.快速排序等排序算法,它们的排序过程都依赖于比较元素间的大小,我们称这些算法为比较排序.本节讨论比较排序算法的运行时间的下界. 给定一个输入序列 < a 1 ...
算法导论06--红黑树构建算法
一.目的 1．熟悉算法设计的基本思想 2．掌握构建红黑树的方法二.内容与设计思想编写随机整数生成算法,生成S到T范围内的N个随机整数并输出: 编写红黑树构建算法,中序遍历各节点,输出颜色和值: 随 ...
【算法导论】python实现算法导论中的算法(传送门)
1.python实现插入排序: 点击跳转. 2.python实现归并排序: 点击跳转. 3.python实现堆排序: 点击跳转. 4.python实现快速排序:点击跳转. 5.python实现计数排序 ...
《算法导论》红黑树详解（一）：概念
在学习红黑树之前,读者应先掌握二叉查找树的相关知识.学习红黑树或者二叉查找树,推荐大家看<算法导论>.<算法导论>原书第3版高清PDF 带详细书签目录下载密码:acis & ...

Java实现算法导论中线性规划单纯形算法

Java实现算法导论中线性规划单纯形算法相关推荐

最新文章

热门文章