link ODE45 solver with the optimization toolbox

Question

0 votes

How can I simultaneously solve the ode45 and pass the output in the Fitness function of the ga toolbox? I am having a problem to link ode45 solver with the optimization toolbox.

I have a set of coupled nonlinear ODEs; I need to solve these ODE's, then optimize two parameters using genetic algorithm for maximizing one output of ode45.

   function dy = systemEQ(t,y, Trf, Vscrap) 
   global c1 c2 Lc K2 PI a viscsw ex alfa D Tsat beta eps V S 
   global condsw rosw Cpsw  Nb
    %Trf = -7; %-6/-7/-8/-9/-10/
   %Vscrap = 78.54 ;%26.179/57.6/78.54
  dy(1) =  ((alfa*D*PI*(Tsat-Trf)^2)/V) + (eps*Vscrap*y(2));
  dy(2) =  (( beta*(Tsat-y(5)))*y(1))+ (((alfa*D*PI*(Tsat-Trf)^2)/V)*Lc) +...
            (c1*eps*Vscrap*y(3)); 
  dy(3) =  (2*( beta*(Tsat-y(5)))*y(2))+ (((alfa*D*PI*(Tsat-Trf)^2)/V)*Lc*Lc)...
           + (c2*eps*Vscrap*y(4));
  dy(4) =  (PI/6 )*((3*( beta*(Tsat-y(5)))*y(3))+ ...
           (((alfa*D*PI*(Tsat-Trf)^2)/V)*Lc*Lc*Lc));  
  dy(5) =  (a*S*((condsw*rosw*Cpsw*Vscrap*Nb/15*PI)^(0.5))*(Trf-y(5)))+...
           (a*(viscsw*(1+(2.5*y(4))+(10.05*y(4)^2)+(2.73*0.001*exp(16.6*y(4))))*...
           (2*PI*ex*Vscrap)^2))+(K2*(3*( beta*(Tsat-y(5)))*y(3)+ ...
           ((alfa*D*PI*(Tsat-Trf)^2)/V)*Lc*Lc*Lc));
  dy= dy(:)
  Solution with ode45
global dsw V Nb eps rosw H D viscsw velosw  
global deltaH Prsw C0 condsw Cps Cpw S roi ex Lc beta alfa Tsat Cpsw a a2 K2 c1 c2 PI
% Données geométriques
            V= 0.000387;%m3
            S= 135.5; %m-1
            H= 0.4;%m
            D= 0.05;%m
            Nb= 2;
% Données caractéristiques de l'eau de mer
            rosw=1028.11; %kg/m3
            viscsw=0.000001826;
            velosw=1449;%m/s
            deltaH= 334700; %j/kg
            roi=917;%kg/m3
            Prsw=13.4;
            C0=0.035;
            condsw=0.563;%j/m.k.s
            Cps=880;%j/kg.k
            Cpw=4182;%j/kg.k
            Tsat= -1.9;%k
% Données cinétiques        
           ex = 17;
           Lc = 0.000005;%m
           beta = 0.000003; % m s-1 K-1
           alfa = 10^8; % s-1m-2K-
           eps = 5; %m-1
%Point de fonctionnement
         %Vscrap = 47.12;%rad/s
         dsw = 0.014;%kg/s  25/50/75
         %Trf= -6; %k
% Constantes        
        PI= 3.14;
        c1= 0.5874;
        c2=0.2599;
% Calcul des parametres et constantes
            Cpsw =(C0*Cps)+((1-C0)*Cpw); %j/kg.k
              %he=(condsw*rosw*Cpsw*Vscrap*Nb/15*PI)^(0.5);%j/m2.k.s
               a=1/(rosw*Cpsw);% k.m3/j
              a2=a*deltaH*roi; %k
              %K1= a*S*he;
              K2= a2*PI/6;%k
  % Conditions initiales
           y01 = 0;       
           y02 = 0;
           y03 = 0;
           y04 = 0;
           y05 = Tsat;
[t,y]= ode45(@systemEQ,0:0.1:30,[y01 y02 y03 y04 y05]);
y1 = y(:,1);
y2 = y(:,2);
y3 = y(:,3);
y4 = y(:,4);
y5 = y(:,5);
%plot (t,y(:,4));
plot (t,y(:,5));

I need to optimize Trf in [-10, -6] and Vscrap in [26, 74] for maximizing Y4 with the constraint Y5<-6. i have solved the ODE's and I need to pass the output in the Fitness function of ga toolbox, so I ask how to do?