MATLAB returning NaN value to optimization variable while using GA.

Hello,

I am carrying out simulation based optimization by connecting MATLAB and Hysys using Hysyslib tool.

The optimization variables are stored in a spreadsheet.

MATLAB is returning NaN value to the optimization variable.

What can be the possible reason for the same.

Thank you

Regards

Ajinkya Pal

[Asvin Kumar]

I suspect you will need to provide more details for anyone from the community to help you. Maybe you can start by mentioning your setup, your optimization variables and what sort of data is involved.

[Ajinkya Pal]

Dear Asvin,

Greetings.

Thank you for you kind reply.

I figured out the problem of NaN, I was using parallel computing and had global variables, so i rectified it. However, I have new problems now.

I am carrying out simulation based optimization by connecting MATLAB and HYSYS using hysylib tool. My objective and constraint functions are blackbox in the sense that they are computed in HYSYS and transported to HYSYS spreadsheet and MATLAB reads from there and supplies the values for variables for simulations to occur.

I have one objective function and 5 non linear constraints and 16 variables. The upper,lower bounds and initial values of variables are as follows:

I converted the contrained problem to un constrained by adding penalty to objective function and deploying dynamic penalty, depending on the current iteration number using particle swarm optimization.

I am deploying try and catch to eliminate the infeasible simulations.

The problem I am facing is that, my simulation crashes after few minutes of run. I want to make sure that my MATLAB code is correct and I am using right optimization algorithm.

I am not sure about the reason of crash. I have few thoughts about it:

1. I am using PSO outputfcn to get the current value of iteration and saving it in a global variable and making it available in the objective function calculation. I am not sure if I am doing it correct or not.

function stop = outfun(optimValues,state)
stop = false;
 
   switch state
       case 'init'
           
       case 'iter'
           % Concatenate current point and objective function
           % value with history. x must be a row vector.
          history.bestval = [history.bestval; optimValues.bestval];
          history.bestx = [history.bestx; optimValues.bestx];
           %saving iteration count to use as dynamic penalty
          penalty= optimValues.iteration;    
       case 'done'
         
   end
end

2. I am using try and catch to tackle infeasible solutions. When the solution is infeasible, the HYSYS Spreadsheet has <empty>. and when MATLAB tries to access it, it shows some error, but the during optimization, the simulation just freezes.

3. Should I scale and center the upper and lower bounds of my variables ? I have listed the variables and their upper,lower bounds above.

The main function called is :

function [history,xsol] = main_computeall


%declaring the penalty
global penalty;

% to store variable array and best solution after each run
history.bestx=[];
history.bestval=[];
penalty=0;

% lower bounds of 16 variables is declared

lb = [0.01 ; 0.01 ; 0.01 ; 0.01 ; 0 ; 0 ; 0 ; 24 ; -165 ; -130 ; 5; -120 ;0.1; 0.3 ; 0.01;-168]; 
   
%upper bound of 16 variables is declared 

ub = [ 0.1 ; 0.15 ; 0.15 ; 0.2 ; 0.15 ; 0.15; 0.12 ; 35 ; -150 ; -100 ; 23 ; -90 ;0.3;0.7;0.08;-158]; 

%Initial point
x0=[1.6e-2,6e-2,8e-2,0.1750,0,0,9e-2,25.6,-156,-115,22,-100,0.2,0.6,5.5e-2,-161];  

 % calling optimization algorithm 
options=optimoptions(@particleswarm,'OutputFcn',@outfun,'Display','iter','SwarmSize',200,'InitialSwarmMatrix',x0,'MaxIterations',15);

%calling function to compute objective function
fun=@computeall;


xsol=particleswarm(fun,16,lb,ub,options);

function stop = outfun(optimValues,state)
stop = false;
 
   switch state
       case 'init'
           
       case 'iter'
           % Concatenate current point and objective function
           % value with history. x must be a row vector.
          history.bestval = [history.bestval; optimValues.bestval];
          history.bestx = [history.bestx; optimValues.bestx];
           %saving iteration count to use as dynamic penalty
          penalty= optimValues.iteration;    
       case 'done'
         
   end
end

end

The function to calculate the objective function is:

function J = computeall(x)

 global penalty;
 
% this function passes the value of variable vector to hsysys and the
% objective and contraints is calculated.
% try and catch is used 
try
    
 pen=penalty;
 %connecting MATLAB and HYSYS
 h=hyconnect;
%hycase=h.Activedocument;
CONS=hyspread(h,'CONSTRAINTS');
OBJ=hyspread(h,'OBJ');
   
    %putting solver on hold
    hysolvertoggle(h);
    
    %Assigning values to the variables in hysys to run
   OBJ.Cell('I1').CellValue=x(1);
   OBJ.Cell('I2').CellValue=x(2);
   OBJ.Cell('I3').CellValue=x(3);
   OBJ.Cell('I4').CellValue=x(4);
   OBJ.Cell('I5').CellValue=x(5);
   OBJ.Cell('I6').CellValue=x(6);
   OBJ.Cell('I7').CellValue=x(7);
   OBJ.Cell('I8').CellValue=x(8);
   OBJ.Cell('I9').CellValue=x(9);
   OBJ.Cell('I10').CellValue=x(10);
   OBJ.Cell('I11').CellValue=x(11);
   
   OBJ.Cell('K1').CellValue=x(12);
   OBJ.Cell('K5').CellValue=x(13);
   OBJ.Cell('K2').CellValue=x(14);
   OBJ.Cell('K3').CellValue=x(15);
   OBJ.Cell('K4').CellValue=x(16);
    
    %putting solver active
   hysolvertoggle(h);


    
  
    %pausing matlab till the simulation runs, initially when the was not
    %used the matlab was going to subsequent lines of code before the
    %objective and constraint was calculated and getting empty values.
    %tw1=hyissolving(h);
    %try
    
    while hyissolving(h)==-1
    pause (1.2); 
    end
    
    %releasing HYSYS object to free memory 
 hyrelease(h);
 
 if pen==0
     k=1;
 else
     k = pen+1;
 end
 

hk= sqrt(k);

    
    c=zeros(5,1);
    theta=zeros(5,1);
    gama=zeros(5,1);
    capH = zeros(5,1);
   
    %Assigning parameter
    MTA=2;
    
    total_fuel_need = OBJ.Cell('e18').CellValue;
   total_turb_pow = OBJ.Cell('h26').CellValue;
  reqgen =  CONS.Cell('e8').CellValue;
    wi_fg = CONS.Cell('B10').CellValue;
    CALC_MTA=CONS.Cell('B11').CellValue;
   hhv_prod = CONS.Cell('C9').CellValue;
    %passing value of objective from computed hsyys 
  
   Jorg=CONS.Cell('e6').CellValue;
   
   
  c(1)=max(0, (MTA - CALC_MTA));
  c(2)=max(0,(total_fuel_need - total_turb_pow));
  c(3)=max(0,(0.99- reqgen));
  c(4)= max(0,(1060 - hhv_prod)) ;
  c(5)= max(0,(28 - wi_fg));
  
for i=1:5
    if c(i)<= 0.001
        theta(i)=10;
    else
        if c(i)<=0.1
            theta(i)=20;
        else
            if c(i)<=1
                theta(i)=100;
            else
                theta(i)=300;
            end
        end
    end
end

for i=1:5
    if c(i)<1
        gama(i)=1;
    else
        gama(i)=2;
    end
end

for i=1:5
    capH(i)= theta(i) * ((c(i))^gama(i));
end

conH = capH(1)+ capH(2)+ capH(3)+ capH(4)+ capH(5);

J = Jorg + hk*conH;
     
catch
    J=100000;
    
    
end

Thank You

Regards

Ajinkya Pal

[Fang Yang]

Dear Ajinkya,

I am working with Small scale lng plant optimization with N2 expander which combine HYSYS with Matlab GAslover, can you please share the corrected codings with me? I am a beginner, I just want to get an idea,

Thanks.

Hope you have a lovely day.