Bonjour,

Je souhaite minimiser des données en utilisant la fonction lsqcurvefit sur la résolution d'une ODE.

Mes données sont un temps (xdata) et une mesure quelconques (ydata).
J'ai différentes séries de données correspondantes à différentes vitesses de refroidissements ("phi"). Les quatre paramètres que je souhaite déterminer sont dépendant de la température et sont déterminer en utilisant "phi" dans ODE.m (T=(Tmax+273)-(phi.*t)).

Jusqu'à présent, j'arrive a déterminer les paramètres qui permettent de résoudre l'ODE pour chaque série de données, en fournissant la valeur de "phi".

Néanmoins je souhaiterais déterminer les quatre paramètres qui résolvent l'ODE pour toute les vitesses de refroidissement, donc chaque série de données avec son "phi" correspondant. Pour moi, un bon moyen était de concaténer les données pour chaque vitesses de refroidissement et d'attribuer la valeur de "phi" correspondante à chaque série de donnée, comme on peut le faire pour une équation linéaire. Cependant pour une ODE le vecteur temps est fourni par la fonction elle même et je ne peux donc pas appliquer cette méthode puisque je ne peux pas savoir quand on passe d'une série à l'autre dans mon vecteur concaténé.

Quelqu'un aurait-il eu à faire au même problème, ou pourrait m'aider sur le sujet? Peut être par une toute autre méthode

J'espère que mon explication est relativement claire.

Merci d'avance

knyun

fit.m
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
%%Loading data
 
 P=xlsread('P05.xlsx');
 xdata=P(:,1);
 ydata=P(:,2);
 
%%fit function
 
 lb=[];
 ub=[];
 options = optimset('Largescale','off','Levenbergmarquardt','on');
 [values] = lsqcurvefit('ODEresol',a0,xdata, ydata,lb,ub,options);
ODEresol.m
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
 function y = ODEresol(a,l)
 global kK1 K01 kK2 K02
 kK1=a(1);
 K01=a(2);
 kK2=a(3);
 K02=a(4);
 init=[0.00000001];
 options = odeset('RelTol',1e-9, 'AbsTol',1e-5);
 tspan=[0 l(end)];
 [t,y] = ode45(@(t,l) ODE(t,l), tspan, init,options);
 y = interp1(t, y(:,1), l);
ODE.m
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
function dl = ODE(t,l)
 global n Tm01 Tm02 Tmax U  phi A1 A2 x0 dx kK1 K01 kK2 K02
coeff;
 
 T=(Tmax+273)-(phi.*t);
 T3=T-273;
 
 w=A2+((A1-A2)./(1+exp((T3-x0)/dx)));
 
 K=(w.*(exp(K01).*exp(-U./(R.*(T-303))).*exp(-kK1./(T.*(Tm01-T).*  (2*T./(Tm01+T)))))+(1-w).*(exp(K02).*exp(-U./(R.*(T-303))).*exp(- kK2./(T.*(Tm02-T).*(2*T./(Tm02+T))))));
 
 dl=zeros(1,1);
 dl(1)=real(n.*K.*(1-l(1)).*((-log(1-l(1))).^((n-1)/n)));