文档章节

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

贝克街的亡灵sf
 贝克街的亡灵sf
发布于 2017/04/24 15:28
字数 1315
阅读 42
收藏 0
 

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; // tableaux

private int M; // number of constraints

private int N; // number of original variables

 

private int[] basis; // basis[i] = basic variable corresponding to row i

// only needed to print out solution, not book

 

// sets up the simplex tableaux

public 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 conditions

assert check(A, b, c);//初始化判断输入的线性规划是否存在初始基本可行解、对偶最优解判断、无限解

}

 

// run simplex algorithm starting from initial BFS

private void solve() {

while (true) {

 

// find entering column q

int q = bland();

if (q == -1) break; // optimal

 

// find leaving row p

int p = minRatioRule(q);

if (p == -1) throw new RuntimeException("Linear program is unbounded");

 

// pivot

pivot(p, q);

 

// update basis

basis[p] = q;

}

}

 

// lowest index of a non-basic column with a positive cost

private 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 cost

private 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; // optimal

else 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 q

for (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 q

for (int i = 0; i <= M; i++)

if (i != p) a[i][q] = 0.0;

 

// scale row p

for (int j = 0; j <= M + N; j++)

if (j != q) a[p][j] /= a[p][q];

a[p][q] = 1.0;

}

 

// return optimal objective value

public double value() {

return -a[M][M+N];

}

 

// return primal solution vector

public 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 vector

public 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 >= 0

for (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 <= b

for (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 >= 0

for (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 >= c

for (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 = yb

private boolean isOptimal(double[] b, double[] c) {

double[] x = primal();

double[] y = dual();

double value = value();

 

// check that value = cx = yb

double 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 tableaux

public 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 = 800

public 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);

}

 

// unbounded

public 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 coefficient

public 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 client

public 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());*/

}

 

}

 

//***************************************

结果

 

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

本文转载自:http://blog.csdn.net/fjssharpsword/article/details/53198213

下一篇: 单例模式
贝克街的亡灵sf
粉丝 2
博文 42
码字总数 21597
作品 0
松江
程序员
私信 提问
【二叉搜索树】Java版本,完美版本二叉搜索树,功能齐全

正文之前 经过几天下午四点到现在的功夫,然后结合算法导论和自己的一些小尝试,总算把Java版本的二叉树给做的完美了~ emm 不信的话大家伙尽管测试。。。终于!!!树形结构!用Java实现出来...

HustWolf
2018/07/14
0
0
飞龙的程序员书单 – 数据结构、算法

入门向 啊哈!算法 这本书真心简洁易懂,dijkstra我是看课本怎么看也看不懂,最后看这本书才懂的。真心推荐。 大话数据结构 工程向 算法 Java实现 C实现 C++实现 普林斯顿的算法课程教材,Cou...

apachecn_飞龙
2016/01/15
0
0
为什么会有C和JAVA,C#抑或其他语言之争呢?

刚工作的人...真的不理解这个... 就老听人说搞C的完全看不起JAVAer,我实在不知道,为啥呢?喷子我问他原因,他就给我个结论... 啥语言,不都是翻译自己的思想而已么?还是,争论的并不是语言...

黑狗
2012/08/16
481
7
数据结构算法书籍推荐

学计算机的人是幸福的,因为在这个领域中有如此多的通俗易懂(相对来说)的经典好书,你需要做的只是坚持把它们一本一本读下去而已。在这里列出一些我看过或者准备看的算法书籍,以供参考。 ...

modernizr
2014/05/22
31.3K
13
线性规划算法原理介绍

线性规划定义: 求满足约束的最优目标,目标是变量的线性函数,约束是变量的相等或不等表达式。 单纯形算法 1 松弛变量 为将不等式转化为等式添加的非负变量 比如 将f(xi) >0 变成 xj= f(xi...

158114527090
2017/06/20
0
0

没有更多内容

加载失败,请刷新页面

加载更多

抽象同步队列AQS——AbstractQueuedSynchronizer锁详解

AQS——锁的底层支持 谈到并发,不得不谈ReentrantLock;而谈到ReentrantLock,不得不谈AbstractQueuedSynchronizer(AQS)! 类如其名,抽象的队列式的同步器,AQS定义了一套多线程访问共享资...

须臾之余
今天
2
0
springboot配置百度UEditor 富文本详解

富文本简介 UEditor是由百度web前端研发部开发所见即所得富文本web编辑器,具有轻量,可定制,注重用户体验等特点,开源基于MIT协议,允许自由使用和修改代码... 准备工作 ueditor需要单独文...

wotrd
昨天
3
0
mysql 5.7之my.cnf配置大全

[client]port = 3306socket = /tmp/mysql.sock[mysqld]###############################基础设置######################################Mysql服务的唯一编号 每个mysql服务...

Online_Reus
昨天
2
0
MAVEN打包时引入外部链接的包

1.项目引入了ORACLE的jar包,MAVEN配置如下 2.打jar包的时候需要指定下main入口函数mainClass <dependency> <groupId>com.oracle</groupId> <artifactId>ojdbc6</artifactId> ......

Cobbage
昨天
2
0

没有更多内容

加载失败,请刷新页面

加载更多

返回顶部
顶部