summaryrefslogtreecommitdiff
path: root/code/linprog
diff options
context:
space:
mode:
authorRemyaDebasis2018-07-23 20:08:46 +0530
committerRemyaDebasis2018-07-23 20:08:46 +0530
commit392bc1326ebccd63e40cb55a82116208a54f2478 (patch)
treea98a596b8c4b64baa45966e3cc1ab75651def780 /code/linprog
parent69460c03b8b53068d60fd08d3180efc91e627603 (diff)
downloadFOT_Examples-392bc1326ebccd63e40cb55a82116208a54f2478.tar.gz
FOT_Examples-392bc1326ebccd63e40cb55a82116208a54f2478.tar.bz2
FOT_Examples-392bc1326ebccd63e40cb55a82116208a54f2478.zip
added code files
Diffstat (limited to 'code/linprog')
-rw-r--r--code/linprog/AirforceBombDuality.sce73
-rw-r--r--code/linprog/MoulderProblemDuality.sce87
-rw-r--r--code/linprog/MoulderProblem_RHSchange.sce78
-rw-r--r--code/linprog/MultiStage_Planning_Problem.sce165
-rw-r--r--code/linprog/ProductMix.sce87
-rw-r--r--code/linprog/Transportation_Problem.sce121
6 files changed, 611 insertions, 0 deletions
diff --git a/code/linprog/AirforceBombDuality.sce b/code/linprog/AirforceBombDuality.sce
new file mode 100644
index 0000000..3bc2ef3
--- /dev/null
+++ b/code/linprog/AirforceBombDuality.sce
@@ -0,0 +1,73 @@
+//Reference: J.K. Sharma,"OPerations Research Theory and Applications", Macmillan Publishers India Ltd., 5th edition, 2013, Chapter 4.
+
+//An Air Force is experimenting with three types of bombs P,Q and R in which three kinds of explosives, viz., A, B and C will be used. Taking the various factors into account, it has been decided to use the maximum 600 kg of explosive A, at least 480 kg of B and exactly 540 kg of explosive C. Bomb P requires 3,2,2 kg, bomb Q requires 1,4,3 kg and bomb R requires 4,2,3 kg of explosives A,B and C respectively. Bomb P is estimated to give the equivalent of a 2 ton explosion, bomb Q, a 3 ton of explosion and bomb R, a 4 ton of explosion respectively. Under what production schedule can the Air Force make the biggest bang?
+//=============================================================================
+// Copyright (C) 2018 - IIT Bombay - FOSSEE
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution. The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+// Author: Remya Kommadath
+// Organization: FOSSEE, IIT Bombay
+// Email: toolbox@scilab.in
+//============================================================================
+
+clc;
+// =============================Primal Model===================================
+nVar = 3;
+
+// Linear inequalities
+A = [3 1 4; -2 -4 -2];
+b = [600 -480]';
+nIneqCon = length(b);
+
+// Linear equality
+Aeq = [2 3 3];
+beq = 540;
+nEqCon = length(beq);
+
+lb = zeros(1,nVar);
+ub = [];
+
+Cost = -[2 3 4 ]';
+
+// =============================Dual Model===================================
+
+// Inequality constraints
+dual_A = -[A;Aeq]';
+dual_beq = Cost;
+n_dualIneq = length(dual_beq);
+Coeff_Slack = [eye(n_dualIneq,n_dualIneq);];
+dual_Aeq = [dual_A Coeff_Slack];
+dual_Cost = [b;beq; zeros(n_dualIneq,1)];
+
+dual_nVar = length(dual_Cost);
+lb = [zeros(1,nIneqCon) -%inf*ones(nEqCon) zeros(1,(dual_nVar-nIneqCon-nEqCon))];
+ub = [];
+
+[xopt,fopt,exitflag,output,lambda]=linprog(dual_Cost, [],[], dual_Aeq, dual_beq, lb,ub);
+clc
+
+select exitflag
+case 0
+ disp(xopt(1:nIneqCon+nEqCon)',"Dual values of the primal constraints")
+ disp(xopt(nIneqCon+nEqCon+1:dual_nVar)',"Dual values of the primal variables")
+ TopRow = ["Variable " "Rate of change in optima"];
+ Table = string([(1:length(xopt))' -xopt]);
+ disp("The change in optima due to one unit increase in each variables are given as below")
+ disp([TopRow;Table])
+ disp(["The optimal cost is ", string(fopt)])
+case 1
+ disp("Primal Infeasible")
+case 2
+ disp("Dual Infeasible")
+case 3
+ disp("Maximum Number of Iterations Exceeded. Output may not be optimal")
+case 4
+ disp("Solution Abandoned")
+case 5
+ disp("Primal objective limit reached")
+case 6
+ disp("Dual objective limit reached")
+end
diff --git a/code/linprog/MoulderProblemDuality.sce b/code/linprog/MoulderProblemDuality.sce
new file mode 100644
index 0000000..c7a13c3
--- /dev/null
+++ b/code/linprog/MoulderProblemDuality.sce
@@ -0,0 +1,87 @@
+// Example problem for sensitivity analysis - using duality theorem to determine the reduced cost
+
+//Reference : Bradley, Hax, and Magnanti,"Applied Mathematical Programming", Addison-Wesley, 1977, chapter 3
+
+//A custom molder has one injection-molding machine and two different dies to fit the machine. Due to differences in number of cavities and cycle times, withthe first die he can produce 100 cases of six-ounce juice glasses in six hours,while with the second die he can produce 100 cases of ten-ounce fancy cocktail glasses in five hours. The molder is approached by a new customer to produce a champagne glass. The production time for the champagne glass is 8 hours per hundred cases. He prefers to operate only on a schedule of 60 hours of production per week. He stores the week’s production in his own stockroom where he has an effective capacity of 15,000 cubic feet. A case of six-ounce juice glasses requires 10 cubic feet of storage space, while a case of ten-ounce cocktail glasses requires 20 cubic feet due to special packaging.The storage space required for the champagne glasses is 10 cubic feet per one case. The contribution of the six-ounce juice glasses is $5.00 per case; however, the only customer available will not accept more than 800 cases per week. The contribution of the ten-ounce cocktail glasses is $4.50 per case and there is no limit on the amount that can be sold where as the contribution of champagne glasses is $6.00 per case with no limit on the bound. How many cases of each type of glass should be produced each week in order to maximize the total contribution?
+//=============================================================================
+// Copyright (C) 2018 - IIT Bombay - FOSSEE
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution. The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+// Author: Remya Kommadath
+// Organization: FOSSEE, IIT Bombay
+// Email: toolbox@scilab.in
+//============================================================================
+
+clc;
+// -----------------------------Primal model-------------------------------
+nProducts = 3;
+TotalProcessingTime = 60;
+RequiredTime = [6 5 8];
+StokeRoomCapacity = 15000;
+RequiredSpace = [10 20 10]*100;
+Demand = 8;
+Revenue = [5 4.5 6];
+nVar = 3 // Each variable indicates number of hundreds of units of each product
+
+// Linear constraints
+// Processing time constraint
+A1 = RequiredTime;
+b1 = TotalProcessingTime;
+// Space constraint
+A2 = RequiredSpace;
+b2 = StokeRoomCapacity;
+// Demand constraint
+nDemand = length(Demand);
+A3 = eye(nDemand,nVar);
+b3 = Demand';
+
+A = [A1;A2;A3]; b = [b1;b2;b3];
+
+//Addition of slack variables
+nInEqConstraints = length(b);
+Aeq = eye(nInEqConstraints,nInEqConstraints);
+Aeq = [A Aeq];
+beq = b;
+
+Cost = -[Revenue zeros(1,nInEqConstraints)]';
+
+// -------------------------------Dual model------------------------------
+dualA = -A';
+dualb = Cost(1:nVar);
+
+dual_nInEqConstraints = length(dualb);
+dualAeq = eye(dual_nInEqConstraints,dual_nInEqConstraints);
+dualAeq = [dualA dualAeq];
+dualbeq = dualb;
+dualCost = beq;
+dualCost = [dualCost; zeros(dual_nInEqConstraints,1)];
+dual_nVar = length(dualCost);
+dual_lb = zeros(1,size(dualAeq,2));
+
+[dual_xopt,dual_fopt,dual_exitflag,dual_output,dual_lambda] = linprog(dualCost,[],[], dualAeq, dualbeq, dual_lb,[]);
+clc;
+select dual_exitflag
+case 0
+ disp(dual_xopt(1:nInEqConstraints)',"Dual values of the primal constraints")
+ disp(dual_xopt(nInEqConstraints+1:dual_nVar)',"Dual values of the primal variables")
+ TopRow = ["Variable " "Rate of change in optima"];
+ Table = string([(1:length(dual_xopt))' -dual_xopt]);
+ disp("The change in optima due to one unit increase in each variables are given as below")
+ disp([TopRow;Table])
+ disp(["The optimal cost is ", string(dual_fopt)])
+case 1
+ disp("Primal Infeasible")
+case 2
+ disp("Dual Infeasible")
+case 3
+ disp("Maximum Number of Iterations Exceeded. Output may not be optimal")
+case 4
+ disp("Solution Abandoned")
+case 5
+ disp("Primal objective limit reached")
+case 6
+ disp("Dual objective limit reached")
+end
diff --git a/code/linprog/MoulderProblem_RHSchange.sce b/code/linprog/MoulderProblem_RHSchange.sce
new file mode 100644
index 0000000..95721c5
--- /dev/null
+++ b/code/linprog/MoulderProblem_RHSchange.sce
@@ -0,0 +1,78 @@
+// Example problem of sensitivity analysis - changes in RHS
+
+//Reference : Bradley, Hax, and Magnanti,"Applied Mathematical Programming", Addison-Wesley, 1977, chapter 3
+
+//A custom molder has one injection-molding machine and two different dies to fit the machine. Due to differences in number of cavities and cycle times, withthe first die he can produce 100 cases of six-ounce juice glasses in six hours,while with the second die he can produce 100 cases of ten-ounce fancy cocktail glasses in five hours. The molder is approached by a new customer to produce a champagne glass. The production time for the champagne glass is 8 hours per hundred cases. He prefers to operate only on a schedule of 60 hours of production per week. He stores the week’s production in his own stockroom where he has an effective capacity of 15,000 cubic feet. A case of six-ounce juice glasses requires 10 cubic feet of storage space, while a case of ten-ounce cocktail glasses requires 20 cubic feet due to special packaging.The storage space required for the champagne glasses is 10 cubic feet per one case. The contribution of the six-ounce juice glasses is $5.00 per case; however, the only customer available will not accept more than 800 cases per week. The contribution of the ten-ounce cocktail glasses is $4.50 per case and there is no limit on the amount that can be sold where as the contribution of champagne glasses is $6.00 per case with no limit on the bound. How many cases of each type of glass should be produced each week in order to maximize the total contribution?
+//=============================================================================
+// Copyright (C) 2018 - IIT Bombay - FOSSEE
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution. The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+// Author: Remya Kommadath
+// Organization: FOSSEE, IIT Bombay
+// Email: toolbox@scilab.in
+//============================================================================
+
+clc;
+
+nProducts = 3;
+TotalProcessingTime = 60;
+RequiredTime = [6 5 8];
+StokeRoomCapacity = 15000;
+RequiredSpace = [10 20 10]*100;
+Demand = 8;
+Revenue = [5 4.5 6];
+nVar = nProducts // Each variable indicates number of hundreds of units of each product
+
+// Linear constraints
+// Processing time constraint
+A1 = RequiredTime;
+b1 = TotalProcessingTime;
+// Space constraint
+A2 = RequiredSpace;
+b2 = StokeRoomCapacity;
+// Demand constraint
+nDemand = length(Demand);
+A3 = eye(nDemand,nVar);
+b3 = Demand';
+
+A = [A1;A2;A3]; b = [b1;b2;b3];
+
+//Addition of slack variables
+nInEqConstraints = length(b);
+Aeq = eye(nInEqConstraints,nInEqConstraints);
+Aeq = [A Aeq];
+beq = b;
+
+Cost = -[Revenue zeros(1,nInEqConstraints)]';
+lb = zeros(1,size(Aeq,2));
+
+[xopt,fopt,exitflag,output,lambda]=linprog(Cost, [],[], Aeq, beq, lb,[]);
+clc;
+select exitflag
+case 0
+ mprintf('The optimum solution has been found \n')
+ disp(xopt(1:nVar)', "The optimum number of six-ounce juice glasses,ten-ounce fancy cocktail glasses and champagne glass")
+ SlackVariables = xopt(nVar+1:length(xopt));
+ disp(SlackVariables',"The surplus resource available")
+ Index = find(SlackVariables>0);
+ newbeq = beq(Index) - SlackVariables(Index);
+ changedbeq = beq;
+ changedbeq(Index) = newbeq;
+ mprintf(' The value beyond %f of the RHS value of constraint number %d will not change the optima. \n ',newbeq,Index)
+ mprintf('\n The optimal cost is %f',fopt)
+case 1
+ disp("Primal Infeasible")
+case 2
+ disp("Dual Infeasible")
+case 3
+ disp("Maximum Number of Iterations Exceeded. Output may not be optimal")
+case 4
+ disp("Solution Abandoned")
+case 5
+ disp("Primal objective limit reached")
+case 6
+ disp("Dual objective limit reached")
+end
diff --git a/code/linprog/MultiStage_Planning_Problem.sce b/code/linprog/MultiStage_Planning_Problem.sce
new file mode 100644
index 0000000..51d37dd
--- /dev/null
+++ b/code/linprog/MultiStage_Planning_Problem.sce
@@ -0,0 +1,165 @@
+//Reference : Bradley, Hax, and Magnanti,"Applied Mathematical Programming", Addison-Wesley, 1977, chapter 5
+
+// The data given in the Multi-stage planning problem are as follows
+//An automobile tire company has the ability to produce both nylon and fiberglass tires. During the next three months they have agreed to deliver tires as follows
+// =========================================
+// Date Nylon Fiberglass
+// -----------------------------------------
+// June 30 4000 1000
+// July 31 8000 5000
+// August 31 3000 5000
+// =========================================
+//The company has two presses, a Wheeling machine and a Regal machine, and appropriate molds that can be used to produce these tires, with the following production hours available in the upcoming months
+// ================================================
+// Wheeling Machine Regal Machine
+// ------------------------------------------------
+// June 700 1500
+// July 300 400
+// August 1000 300
+// ================================================
+//The production rates for each machine-and-tire combination, in terms of hours per tire, are as follows:
+//======================================================
+// Wheeling Machine Regal Machine
+//------------------------------------------------------
+//Nylon 0.15 0.16
+//Fiberglass 0.12 0.140
+//======================================================
+//The variable costs of producing tires are $5.00 per operating hour, regardlessof which machine is being used or which tire is being produced. There is also an inventory-carrying charge of $0.10 per tire per month. How should the production be scheduled in order to meet the delivery requirements at minimum costs?
+//Reported solution:
+//X = [1867 7633 3500 0 5500 2500 0 2500 2500 0 0 2667 333 5000 0];
+//===========================================================================
+// Copyright (C) 2015 - IIT Bombay - FOSSEE
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution. The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+// Author: Remya Kommadath
+// Organization: FOSSEE, IIT Bombay
+// Email: toolbox@scilab.in
+//============================================================================
+clc;
+nProducts = 2;
+nMachines = 2;
+nPeriods = 3;
+Demand = [4000 1000;8000 5000; 3000 5000];
+AvailMachineTime = [700 1500;300 400; 1000 300];
+ProdRate = [0.15 0.16;0.12 0.14];
+OperatingCost = 5;
+InventoryCost = 0.1;
+
+nVar = nProducts*nMachines*nPeriods + (nPeriods-1)*nProducts; // Dimension of the problem is determined
+
+// Linear equality constraints
+// Demand constraints
+nEqConstraints = nProducts*nPeriods
+Aeq = zeros(nEqConstraints,nVar);
+// Demand constraints for period 1
+Aeq1 = zeros(nProducts,nVar);
+for i = 1:nProducts
+ index1 = (i-1)*nProducts+1:i*nProducts;
+ Aeq1(i,index1) = 1;
+ index2 = nMachines*nProducts+i;
+ Aeq1(i,index2) = -1;
+ beq1(i,1) = Demand(1,i);
+end
+// Demand constraints for period 2 to (nPeriods-1)th period
+Aeq2 = zeros(nProducts,nVar);
+for i = 2:nPeriods-1
+ for j = 1:nProducts
+ index3 = (i-1)*(nProducts*nMachines+nProducts)+(j-1)*nMachines+1:(i-1)*(nProducts*nMachines+nProducts)+(j-1)*nMachines+nMachines;
+ Aeq2(j,index3) = 1;
+ index4 = (i-1)*(nProducts*nMachines+nProducts) - nProducts+j;
+ Aeq2(j,index4) = 1;
+ index5 = i*(nProducts*nMachines+nProducts)- nProducts+j
+ Aeq2(j,index5) = -1;
+ beq2(j,1) = Demand(i,j);
+ end
+end
+
+// Demand constraints for last period
+Aeq3 = zeros(nProducts,nVar);
+for i = 1:nProducts
+ index6 = (nProducts*nMachines+nProducts)*(nPeriods-1)+(i-1)*nProducts+1:(nProducts*nMachines+nProducts)*(nPeriods-1)+i*nProducts
+ Aeq3(i,index6) = 1;
+ index7 = (nProducts*nMachines+nProducts)*(nPeriods-1) - nProducts+i;
+ Aeq3(i,index7) = 1;
+ beq3(i,1) = Demand(nPeriods,i);
+end
+Aeq = [Aeq1;Aeq2;Aeq3];
+beq = [beq1;beq2;beq3];
+
+// Linear inequality constraints
+// Machine time constraints
+for i = 1:nPeriods
+ for j = 1:nProducts
+ Cindex = (i-1)*(nProducts*nMachines+nProducts)+j:nProducts:(i-1)*(nProducts*nMachines+nProducts)+nMachines+j;
+ Rindex = (i-1)*nProducts+j;
+ A(Rindex,Cindex) = ProdRate(:,j)';
+ b(Rindex,1) = AvailMachineTime(i,j);
+ end
+end
+
+// Objective function
+TotalProductionCost = [];
+for j = 1:nProducts
+ TotalProductionCost = [TotalProductionCost ProdRate(j,:)*OperatingCost];
+end
+
+for i = 1:nPeriods
+ index = (i-1)*(nProducts*nMachines+nProducts)+1:(i-1)*(nProducts*nMachines+nProducts)+nProducts*nMachines;
+ nindex = length(index);
+ cost(index,1) = TotalProductionCost';
+ cost(index(nindex)+1:index(nindex)+nProducts,1) = InventoryCost;
+end
+cost(nVar+1:nVar+nProducts,1) = [];
+lb = zeros(1,nVar);
+
+[xopt,fopt,exitflag,output,lambda]=linprog(cost, A, b, Aeq, beq, lb,[]);
+
+//Result representation
+
+select exitflag
+case 0
+ disp(" Optimal Solution Found")
+ M = [" "];
+ for m = 1:nMachines
+ M = [M strcat(["Machine",string(m)])];
+ end
+ P = [];
+ for p = 1:nProducts
+ P = [P;strcat(["Product ",string(p)])];
+ end
+
+ for i = 1:nPeriods
+ Sol = [];
+ for j = 1:nProducts
+ Ind1 = (i-1)*(nProducts*nMachines + nProducts)+(j-1)*nProducts+1:(i-1)*(nProducts*nMachines + nProducts)+j*nProducts;
+ Sol = [Sol;xopt(Ind1)'];
+ end
+
+ disp(strcat(["Production schedule for the Period ", string(i)]));
+ disp([M; [P string(Sol)]]);
+ end
+
+ for i = 1:nPeriods-1
+ ind = i*(nProducts*nMachines+1:nProducts*nMachines+nProducts);
+ inventory = xopt(ind);
+ disp(strcat(["Inventory at the Period ", string(i)]));
+ disp([P string(inventory)]);
+ end
+
+ disp(["The optimal cost is ", string(fopt)])
+case 1
+ disp("Primal Infeasible")
+case 2
+ disp("Dual Infeasible")
+case 3
+ disp("Maximum Number of Iterations Exceeded. Output may not be optimal")
+case 4
+ disp("Solution Abandoned")
+case 5
+ disp("Primal objective limit reached")
+case 6
+ disp("Dual objective limit reached")
+end
diff --git a/code/linprog/ProductMix.sce b/code/linprog/ProductMix.sce
new file mode 100644
index 0000000..5d66fed
--- /dev/null
+++ b/code/linprog/ProductMix.sce
@@ -0,0 +1,87 @@
+//Refernce: H.P. Williams,"Model building in mathematical programming",Fifth Edition,AJohn Wiley & Sons, Ltd., Publication,2013, chapter 1.2
+
+//An engineering factory can produce five types of product (PROD 1, PROD 2,... , PROD 5) by using two production processes: grinding and drilling. After deducting raw material costs, each unit of each product yields the following contributions to profit:
+// =====================================================
+// PROD1 PROD2 PROD3 PROD4 PROD5
+// 550 600 350 400 200
+// =====================================================
+//Each unit requires a certain time on each process. These are given below (in
+//hours). A dash indicates when a process is not needed.
+// =================================================================
+// PROD1 PROD2 PROD3 PROD4 PROD5
+// -----------------------------------------------------------------
+// Grinding 12 20 - 25 15
+// Drilling 10 8 16 - -
+// =================================================================
+//In addition, the final assembly of each unit of each product uses 20 hours of
+//an employee’s time. The factory has three grinding machines and two drilling machines and works a six-day week with two shifts of 8 hours on each day. Eight workers are employed in assembly, each working one shift a day. The problem is tofind how much of each product is to be manufactured so as to maximize the totalprofit contribution. The minimum demand for each product is added to this problem and details are as follows
+// =================================================================
+// PROD1 PROD2 PROD3 PROD4 PROD5
+// -----------------------------------------------------------------
+// Demand 5 5 2 7 2
+// =================================================================
+
+//=============================================================================
+// Copyright (C) 2015 - IIT Bombay - FOSSEE
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution. The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+// Author: Remya Kommadath
+// Organization: FOSSEE, IIT Bombay
+// Email: toolbox@scilab.in
+//============================================================================
+
+clc;
+nProducts = 5;
+nMachines = [3 2];
+nWorkers = 8;
+nDays = 6;
+nHours = 8;
+nMachineShift = 2;
+nWorkerShift = 1;
+prodProfit = [550 600 350 400 200]';
+ProcessingTime = [12 20 0 25 15; 10 8 16 0 0; 20*ones(1,nProducts)];
+
+// Linear inequalities
+A1 = ProcessingTime;
+b1 = (nDays*nHours)*[nMachines*nMachineShift nWorkers*nWorkerShift]';
+// Demand constraints
+A2 = [-1 0 0 0 0;0 -1 0 0 0;0 0 -1 0 0; 0 0 0 -1 0; 0 0 0 0 -1];
+b2 = [-5 -5 -2 -7 -2]';
+
+A = [A1;A2];
+b = [b1;b2];
+lb = zeros(1,nProducts);
+Cost = -prodProfit;
+[xopt,fopt,exitflag,output,lambda]=linprog(Cost, A, b, [], [], lb,[]);
+
+// Result representation
+select exitflag
+case 0
+ disp(" Optimal Solution Found")
+ P = [];
+ for p = 1:nProducts
+ P = [P strcat(["PROD",string(p)])];
+ end
+
+ disp([P;string(xopt')]);
+ disp(["The optimal cost is ", string(fopt)])
+case 1
+ disp("Primal Infeasible")
+ nViolatedConstraints = sum(A*xopt-b>0)
+ if nViolatedConstraints
+ disp(strcat(["No of constraints violated : ", string(nViolatedConstraints)]));
+ end
+case 2
+ disp("Dual Infeasible")
+case 3
+ disp("Maximum Number of Iterations Exceeded. Output may not be optimal")
+case 4
+ disp("Solution Abandoned")
+case 5
+ disp("Primal objective limit reached")
+case 6
+ disp("Dual objective limit reached")
+end
diff --git a/code/linprog/Transportation_Problem.sce b/code/linprog/Transportation_Problem.sce
new file mode 100644
index 0000000..9df7c49
--- /dev/null
+++ b/code/linprog/Transportation_Problem.sce
@@ -0,0 +1,121 @@
+//Reference : Dantzig G. B., Linear Programming and Extensions, Princeton University Press, Princeton, New Jersey, 1963, Chapter 3-3.
+
+//The objective of the problem is to minimize the transportation cost of goods from the canneries to the warehouses,by satisfying the supply and demand constraints. Inorder to simplify the problem, it is assumed that there are two canneries and three warehouses. The supply data of canneries are as follows
+// =============================
+// Supply
+// Cannery I 350
+// Cannery II 650
+// =============================
+//The demands at the warehouses are as given
+// =============================
+// Demand
+// Warehouse A 300
+// Warehouse B 300
+// Warehouse C 300
+// =============================
+//ShippingCost between the canneries and warehouses are
+// =========================================================
+// Warehouse A Warehouse B Warehouse C
+// =========================================================
+// Cannery I 2.5 1.7 1.8
+// Cannery II 2.5 1.8 1.4
+// =========================================================
+// It should be noted that, this case is solved without considering the inventory.
+//
+//===========================================================================
+// Copyright (C) 2015 - IIT Bombay - FOSSEE
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution. The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+// Author: Remya Kommadath
+// Organization: FOSSEE, IIT Bombay
+// Email: toolbox@scilab.in
+//============================================================================
+
+clc;
+
+nCanneries = 2;
+nWarehouses = 3;
+ShippingCost = [2.5 1.7 1.8; 2.5 1.8 1.4];
+demand = [300 300 300]';
+supply = [350 650]';
+
+nVar = nCanneries*nWarehouses;
+lb = zeros(1,nVar); // lower bounds of the amount of commodity to be shipped from each canneries to warehouses;
+
+// Linear inequality constraints
+// Total supply to the warehouses from a cannery should not exceed the corresponding supply limit of the same cannery
+for i = 1:nCanneries
+ indeX = (i-1)*nWarehouses+1:i*nWarehouses;
+ As(i,indeX) = ones(1,nWarehouses);
+end
+Bs = supply;
+As
+// Total supply to a warehouses from the canneries should satisfy the corresponding demand of the same warehouses
+//Ad = zeros()
+for i = 1:nWarehouses
+ indeX = i:nWarehouses:nVar;
+ Ad(i,indeX) = -1;
+end
+Bd = -demand;
+
+A = cat(1,As,Ad);
+b = cat(1,Bs,Bd);
+
+// Minimize the total transportation cost between all plants and markets
+Cost = zeros(nVar,1);
+for i = 1:nCanneries
+ indeX = (i-1)*nWarehouses+1:i*nWarehouses;
+ Cost(indeX) = (ShippingCost(i,:)');
+end
+
+[xopt,fopt,exitflag,output,lambda] = linprog(Cost,A,b,[],[],lb,[]);
+
+// Display of results
+select exitflag
+case 0
+ disp(" Optimal Solution Found")
+ // xopt is transformed into a nCanneries X nWarehouses sized string matrix
+ for i = 1:nCanneries
+ indeX = (i-1)*nWarehouses+1:i*nWarehouses;
+ X(i,1:nWarehouses) = xopt(indeX,1)';
+ end
+ Values = string(X);
+
+ // Labelling each markets as M_1,M_2,...,M_nWarehouses
+ Markets(1,1) = [" "];
+ for j = 2:nWarehouses+1
+ Markets(1,j) = strcat(["M_",string(j-1)]);
+ end
+
+ // Labelling each plants as P_1,P_2,...,P_nCanneries
+ for i = 1:nCanneries
+ Plants(i) = strcat(["P_",string(i)]);
+ end
+
+ table = [Markets; [Plants Values]];
+ f = gcf();
+ clf
+ as = f.axes_size; // [width height]
+ ut = uicontrol("style","table",..
+ "string",table,..
+ "position",[150 as(2)-220 300 87]) // => @top left corner of figure
+
+
+ // Modify by hand some values in the table. Then get them back from the ui:
+ matrix(ut.string,size(table))
+case 1
+ disp("Primal Infeasible")
+case 2
+ disp("Dual Infeasible")
+case 3
+ disp("Maximum Number of Iterations Exceeded. Output may not be optimal")
+case 4
+ disp("Solution Abandoned")
+case 5
+ disp("Primal objective limit reached")
+case 6
+ disp("Dual objective limit reached")
+end