GPmain1

 Calls all functions to define the MIGP problem  
--------------------------------------------------------------------------
 Date: July 30, 2010, ver01
 By: Eugene Ben-Awuah
--------------------------------------------------------------------------

0001 % Calls all functions to define the MIGP problem
0002 %--------------------------------------------------------------------------
0003 % Date: July 30, 2010, ver01
0004 % By: Eugene Ben-Awuah
0005 %--------------------------------------------------------------------------
0006 
0007 function GPmain1
0008 
0009 %Put together the objective function matrix, 'c' and save for optimization. The
0010 %inputToMILPBlocks file contains details of the block model and mining
0011 %details
0012 clear
0013 load inputToGPBlocks;
0014 
0015 %the model assumes that there exists a pre-emptive penalty cost for deviating from
0016 %the goals. the most costly goal deviation is mining capacity, C1; then processing
0017 %capacity, C2; and then dyke material requirement, C3. This penalty cost in
0018 %deviating from the target is to convert the deviational tonnes to dollar
0019 %values. NB: cost should be POSTIVE to have a negative impact on the NPV; refer your notes.
0020 %C1=4; C2=3; C3=2; %Order of deviational cost using abitrary values.
0021 
0022 %the model assumes that there exists a pre-emptive priority structure among
0023 %the goals. the most important goal is mining capacity, P1; then processing
0024 %capacity, P2; and then dyke material requirement, P3; refer your notes.
0025 %P1=1; P2=1; P3=1; %Order of goal importance using abitrary values.
0026 Blocks =Blocks120; %numOfPeriods =5; interestRate = 0.1;
0027 [c normVPQ] = objective_function(Blocks, numOfPeriods, interestRate,...
0028                                  C1, C2, C3, P1, P2, P3);
0029 normC = normVPQ;
0030 save C c normC;
0031 clear c normVPQ;
0032 
0033 %Put together the constraint matrix, 'A', by combining the goals matrix and
0034 %the constraint matrix as well as the upper bound matrix. The miningSlack
0035 %is used to introduce a lower bound to find a feasible optimal solution.
0036 %The goals first.
0037 %miningTarget =1950000; miningSlack =0;
0038 load inputToGPBlocks;
0039 [A_miningG miningGoal miningGoal_LB]=mining_capacity_goal(Blocks, numOfPeriods,...
0040                        miningTarget, miningSlack);
0041 
0042 %processingTarget =1350000; processingSlack =0;
0043 [A_processingG processingGoal processingGoal_LB]=...
0044          processing_capacity_goal(Blocks, numOfPeriods, processingTarget,... 
0045                                                          processingSlack);
0046 
0047 %dykemTarget =520000; dykemSlack =0;
0048 [A_dykemG dykemGoal dykemGoal_LB]=...
0049          dykem_capacity_goal(Blocks, numOfPeriods, dykemTarget, dykemSlack);
0050 
0051 %Put the goals matrix together and its upper bound
0052 A_goals = [A_miningG; A_processingG; A_dykemG];
0053 b_goals = [miningGoal; processingGoal; dykemGoal];
0054 b_L_goals = [miningGoal_LB; processingGoal_LB; dykemGoal_LB];
0055 clear A_miningG A_processingG A_dykemG miningGoal processingGoal dykemGoal
0056 
0057 %Constraints matrix start from here
0058 %Maxgrade =0.15; Mingrade =0.06;
0059 [A_grade b_Ugrade] = grade_blending_constraint1(Blocks, numOfPeriods,...
0060                      Maxgrade, Mingrade);
0061 
0062 %MaxfinesOre =0.30; MinfinesOre =0;
0063 [A_finesOregrade b_finesOreUgrade] = fines_ore_blending_constraint1...
0064                                      (Blocks, numOfPeriods,...
0065                                      MaxfinesOre, MinfinesOre);
0066 
0067 %MaxfinesDyke =0.5; MinfinesDyke=0;
0068 [A_finesDykegrade b_finesDykeUgrade] = fines_dyke_blending_constraint1...
0069                                        (Blocks, numOfPeriods,...
0070                                         MaxfinesDyke, MinfinesDyke);
0071                                      
0072 [A_ore_mining b_ore_mining] = ore_mining_constraint(Blocks, numOfPeriods);  
0073 
0074 [A_mining_d b_mining_d] = deviation_positivity_constraint(Blocks, numOfPeriods); 
0075      
0076 [A_slope b_Uslope] = slope_constraints(blocksAbove,...
0077                                        numBlocks,...
0078                                        numOfPeriods); 
0079 
0080 [A_slope2 b_Uslope2] = slope2_constraints(numBlocks,...
0081                                           numOfPeriods); 
0082 
0083 %Put the constraint matrix together and its upper bound
0084 A_constraints = [A_grade; A_finesOregrade; A_finesDykegrade; A_ore_mining;...
0085                 A_mining_d; A_slope; A_slope2];
0086 b_constraints = [b_Ugrade; b_finesOreUgrade; b_finesDykeUgrade;...
0087                 b_ore_mining; b_mining_d; b_Uslope; b_Uslope2];
0088             
0089 clear A_grade A_finesOregrade A_finesDykegrade A_ore_mining...
0090                 A_mining_d A_slope A_slope2;
0091 clear b_Ugrade b_finesOreUgrade b_finesDykeUgrade...
0092                 b_ore_mining b_mining_d b_Uslope b_Uslope2;
0093             
0094 %Put the goal matrix and the constraint matrix together to get the
0095 %'A MATRIX' and its upper bound and save
0096 A = [A_goals; A_constraints];
0097 b_U = [b_goals; b_constraints];
0098 save A A 
0099 save b_U b_U
0100 clear b_U
0101 clear A_goals A_constraints b_constraints
0102 
0103 %Put together the lower bound. The equality constraint first and then the
0104 %inequality constraint. The inequality constraint was set to negative
0105 %infinity and the equality constraint is set as same as the upper bound
0106 [m,n] = size(A);
0107 [k1,k2] =size(b_goals);
0108 
0109 %For optimization solvency, a slack of 1% in the equality constraint lower
0110 %bound was used to enable a feasible solution.
0111 %eq_slack = 0.027;
0112 %b_L_goals = b_goals-eq_slack;
0113 b_L     = -inf*ones(m-k1,1);    %lower bound; k1 is the size of the equality
0114 b_L     = [b_L_goals;b_L];      %constraint which must be set as the same
0115 save b_L b_L                    %as the upper bound to enforce equality.
0116 clear
0117 
0118 % call tomlab/cplex sovler on the problem
0119 load inputToGPBlocks;
0120 [Result]= goal_program();
0121 
0122 %plots
0123 %numOfPeriods =5; miningTarget =1950000; processingTarget =1350000; dykemTarget =520000;
0124 %Plots for the goals that was set prior to optimization
0125 miningGoal = ([miningTarget]')/1000000;
0126 processingGoal = ([processingTarget]')/1000000;
0127 dykemGoal = ([dykemTarget]')/1000000;
0128 
0129 %Plots for acceptable deviations that was set prior to optimization
0130 %miningd1 = (repmat(miningd1,numOfPeriods,1))/1000000;
0131 miningd1 = ([miningd1]')/1000000;
0132 %processingd2 = (repmat(processingd2,numOfPeriods,1))/1000000;
0133 processingd2 = ([processingd2]')/1000000;
0134 %dykemd3 = (repmat(dykemd3,numOfPeriods,1))/1000000;
0135 dykemd3 = ([dykemd3]')/1000000;
0136 %dykemd4 = (repmat(dykemd4,numOfPeriods,1))/1000000;
0137 dykemd4 = ([dykemd4]')/1000000;
0138 
0139 load Solution;
0140 load C normC;
0141 load Discounted_Values;
0142 output_plot_blocks1(Blocks120, Result, numOfPeriods, numBlocks,...
0143                           normC, xdim, ydim, zdim, v, p, q,...
0144                           miningGoal, processingGoal, dykemGoal,...
0145                           miningd1, processingd2, dykemd3, dykemd4) 
0146 
0147 %close all
0148 %close figure 1 figure 2 figure 3 figure 4 figure 5 figure 6 figure 7...
0149       %figure 8 figure 9 figure 10 figure 11 figure 12 figure 13 figure 14...
0150       %figure 15 figure 16 figure 17 figure 18 figure 19 figure 20;
0151 clear;                      
0152      load A; load b_U; load C; load Discounted_Values
0153      save problem;
0154      
0155      clear;
0156      load Solution; load Period; load scheduledBlocks; load graphs;
0157      save Solution;
0158      clear;
0159      
0160      delete A.mat; delete A_slope.mat; delete b_U.mat; delete Period.mat;
0161      delete C.mat; delete b_L.mat; delete xA.mat; delete xT.mat;
0162      delete yA.mat; delete scheduledBlocks.mat;
0163      delete graphs.mat; delete Discounted_Values.mat
0164 
0165 end