Description

MATLAB code

Results

Table 1. Solution report

 

 

Description

Case study Stochastic Utility Problem (see Formal Problem Statement) in MATLAB Environment is solved with mpsg_solver PSG function.

 

Main MATLAB code is in file CS_Stochastic_Utility_Problem_Scen_Lin_Blocks_Toolbox.m.

Data are saved in file CS_Utility_data.zip:

 

MATLAB code

Let us describe the main operations. To run case study you need to do the following main steps:

 

In file CS_Stochastic_Utility_Problem_Scen_Lin_Blocks_Toolbox.m:

 

Reading and checking problem data:

 

fid = fopen('l500.ins','r');

if fid<0, disp(' '); disp(' DATA_file does not exist'); break;

end  

 

Read n:

 

rline=fgetl(fid);

if rline == -1, disp(' '); disp(['Error in ', file_name, ': Keyword "n=" is not found']); terminate %End  %quit  %  Stop

else; rline=strtrim(lower(rline)); rl=length(rline);

if(rl>0); nw = sscanf(rline, '%i'); n=nw(1);

end  

end  

 

Read V_k, S_k:

 

k=0;

while true, rline=fgetl(fid);

if rline == -1, break;

else; rl=length(rline);

  if(rl>2), [v(k+1),s(k+1)] = strread(rline); k=k+1; else, break;

  end  

end  

end

fclose(fid);

 

 

Read scenarios for scenarios problem statement:

 

ps = fopen('Utility_senario_forL500_Short.txt', 'r');

 

J=0; tic;

while true, rline=fgetl(ps);

if rline == -1, break; else; J=J+1; end  

end

frewind(ps); disp(['counting scenarios ',num2str(toc)]);

 

matrix_scen(J,n)=0; l=0; tic;

while true, rline=fgetl(ps);

if rline == -1, break;

else; rl=length(rline); l=l+1;

   ip = findstr(rline(8:14),':');

   matrix_scen(l,1:n)=strread(rline(8+ip:rl), '%f' )';

end

end

disp(['reading scenarios ',num2str(toc)]);

 

Make pmatrix_scenario:

 

pmatrix_scen((n+1)*J,3)=0.; id=0; tic;

for j=1:J

for i=1:n; id=id+1;

   pmatrix_scen(id,1) = j;

   pmatrix_scen(id,2) = i;

   pmatrix_scen(id,3) = matrix_scen(j,i);

end

 id=id+1;

 pmatrix_scen(id,1) = j;

 pmatrix_scen(id,2) = n+j;                   % variable T_j

 pmatrix_scen(id,3) = -1;                      

end

Spmatrix_scen = sparse(pmatrix_scen(:,1),pmatrix_scen(:,2),pmatrix_scen(:,3));

clear matrix_scen pmatrix_scen;

disp(['filling pmatrix of scenarios ',num2str(toc)]);

 

Pack matrices and point:

 

variables_X='x1'; variables_T='T_1';

for j=2:n; variables_X =[variables_X, char(9),'x',int2str(j)]; end

for j=2:J; variables_T =[variables_T, char(9),'T_',int2str(j)]; end

count=1;

toolboxstruc_arr(count) = tbpsg_pmatrix_pack('pmatrix_scen', Spmatrix_scen, [variables_X char(9) variables_T]);

 

 matrix_block_k(1,2)=0.;

for l=1:k

   count=count+1;

   toolboxstruc_arr(count) = tbpsg_matrix_pack(['matrix_',int2str(l)], -s(l),['T_1'],v(l));

end

 avg_matrices=['block_',int2str(J),'(matrix_1)'];

for l=2:k; avg_matrices=[avg_matrices ', block_',int2str(J),'(matrix_',int2str(l),')'];

end

 

count=count+1;

toolboxstruc_arr(count) = tbpsg_matrix_pack('matrix_budget', ones(1,n),variables_X);

 

count=count+1;

toolboxstruc_arr(count) = tbpsg_point_pack('point_lowerbounds',zeros(1,n),variables_X);

 

Form problem statement:

 

problem_statement=sprintf('%s\n',...

  'problem: problem_Utility_Scen_Block, type = minimize',...

  'objective: objective_avg_multi, linearize = 0',...  

   ['  avg_max_risk_(', avg_matrices ,')'],...

  'constraint: constraint_linear_mult, lower_bound = 0, upper_bound = 0',...

  '  linearmulti_1(pmatrix_scen)',...

  'constraint: constraint_budget, lower_bound = 1, upper_bound = 1',...

  '  linear_budget(matrix_budget)',...

  'box_of_variables: lowerbounds = 0',...

  'solver: car, precision = 7',...

  '');

 

tic

 

%Uncomment the following line to open the problem in Toolbox Window:

%tbpsg_toolbox(problem_statement,toolboxstruc_arr);

 

Optimize problem:

 

[solution_str, outargstruc_arr] = tbpsg_run(problem_statement, toolboxstruc_arr);

 

Display solution and optimal point:

 

disp(solution_str);

optimal_point_data = tbpsg_optimal_point_data(solution_str, outargstruc_arr);

optimal_point_vars = tbpsg_optimal_point_vars(solution_str, outargstruc_arr);

 

Results

Table 1. Solution report

 

Problem: solution_status = optimal

Timing: Data_loading_time = 0.02, Preprocessing_time = 0.02, Solving_time = 0.02

Variables: optimal_point = point_problem_utility_scen_block

Objective: objective_avg_multi = -8.707287523413

Constraint: constraint_linear_mult = vector_constraint_linear_mult, [6.661338147751e-016]

Constraint: constraint_budget = 1.000000000000e+000, [-2.220446049250e-016]

Function: avg_max_risk_(block_13(matrix_1), block_13(matrix_2), block_13(matrix_3), block_13(matrix_4), block_13(matrix_5), block_13(matrix_6), block_13(matrix_7), block_13(matrix_8), block_13(matrix_9), block_13(matrix_10), block_13(matrix_11)) = -8.707287523413e+000

Function: linearmulti_1(pmatrix_scen) = vector_constraint_linear_mult

Function: linear_budget(matrix_budget) = 1.000000000000e+000