Re: Struggling with FindFit for a model which takes results from
- To: mathgroup at smc.vnet.net
- Subject: [mg108288] Re: Struggling with FindFit for a model which takes results from
- From: Valeri Astanoff <astanoff at gmail.com>
- Date: Fri, 12 Mar 2010 07:13:43 -0500 (EST)
- References: <hn5b2i$6q2$1@smc.vnet.net>
Good day,
I tried and fix a few things :
(I truncated data to 10 minutes)
In[1]:= Lp := 2.64*^-27;(*mol^2/(min atm um^5)*)
Rad = 63;(*um*)A = 4 \[Pi] Rad^2;(*um^2*)
R = 0.08206;(*L atm/(mol K)*)
T = 286;
minutes = 10;
Vb = 0.19;
Vbv = 4/3 \[Pi] Rad^3*Vb (*um^3*);
volconv = 1*^15 (*converts um^3 to L;1x10^15 um^3 per L*);
molconv = 1*^6 (*converts moles to micromoles;1x10^6 umol per mol*);
PMVWater = 1.8*^7;(*um^3/umol*)PMVDMSO = 7.13*^7;
PCPA1 := 4.38*^-31;
BSalt = 2.749 (*salt*);
BCPA1 = 2.423 (*DMSO*);
insecondvirialcoefficientlist = {BSalt, BCPA1};
exsecondvirialcoefficientlist = {BSalt, BCPA1};
NinWater0 = 4/3 \[Pi] Rad^3 (1 - Vb)/PMVWater;
NinSalt =
1.663*2.062*^-4 (*Salts;assumes 0.3 osmolal=8A 0.176 molal*);
NinCPA10 = 0.1*^-100(*DMSO*);
In[19]:= totmolin[t_] := (NinWater[t] + NinSalt + NinCPA1[t]);
XinSalt[t_] := NinSalt/totmolin[t];
XinCPA1[t_] := NinCPA1[t]/totmolin[t];
insolutemolefxlist[t_] := {XinSalt[t], XinCPA1[t]};
waterin[t_] :=
Total[insolutemolefxlist[t]] (1 +
insecondvirialcoefficientlist.insolutemolefxlist[t]);
NexWater = 0.056*^6 (*water*);
NexSalt = 1.663*245;(*assumes mole fraction of 0.003*)
NexCPA1 = 1694 (*DMSO*);
totmolex = NexWater + NexSalt + NexCPA1;
XexWater = NexWater/totmolex;
XexSalt = NexSalt/totmolex;
XexCPA1 = NexCPA1/totmolex;
exsolutemolefxlist = {XexSalt, XexCPA1} ;
waterex =
Total[exsolutemolefxlist]*(1 +
exsecondvirialcoefficientlist.exsolutemolefxlist);
V0 = NinWater0*PMVWater + Vbv + NinCPA10*PMVDMSO ;
tref = T + 273.16;
volconv = 1*^15;
In[35]:= (*The 2 parameters I wish to find are Lp and PCPA1 here in \
the 2 ODEs*)
solution[Lp_?NumericQ, PCPA1_?NumericQ] :=
NDSolve[{NinWater'[t] ==
8 A - Lp*volconv*molconv*A*R*T (waterex - waterin[t]),
NinCPA1'[t] == 8 A PCPA1*volconv*molconv*A*R*T (Log[XexCPA1] +
(0.5 - BCPA1) XexWater (1 - XexCPA1) -
(0.5 - BSalt) XexWater XexSalt -
Log[(NinCPA1[t]/(NinWater[t] +
NinSalt + NinCPA1[t]))] - (0.5 - BCPA1) (NinWater[t]/
(NinWater[t] + NinSalt + NinCPA1[t])) (1 - (NinCPA1[t]/
(NinWater[t] + NinSalt + NinCPA1[t]))) + (0.5 -
BSalt) (NinWater[t]/
(NinWater[t] + NinSalt +
NinCPA1[t])) (NinSalt/(NinWater[t] + NinSalt +
NinCPA1[t]))),
NinWater[0] == 8 A NinWater0, NinCPA1[0] == 8 A NinCPA10},
{NinWater, NinCPA1}, {t, 0, minutes}]
(*End of ODEs*)
In[36]:= fNinWater[t_?NumericQ, Lp_?NumericQ, PCPA1_?NumericQ] :=
(NinWater /. solution[Lp, PCPA1][[1]])[t];
fNinCPA1[t_?NumericQ, Lp_?NumericQ, PCPA1_?NumericQ] :=
(NinCPA1 /. solution[Lp, PCPA1][[1]])[t];
Check :
In[38]:= fNinWater[1, 1, 1]
Out[38]= 0.0471314
In[39]:= (*Start of Curve Fitting*)
Vtotal[t_?NumericQ, Lp_?NumericQ, PCPA1_?NumericQ] :=
(Vbv + fNinWater[t, Lp, PCPA1]*PMVWater +
fNinCPA1[t, Lp, PCPA1]*PMVDMSO);
(*So I am making this calculation using output from the ODEs
to determine the total volume,which includes a fixed volume (Vbv),
the result from NinWater'[t]*molar volume of water,and likewise for \
the DMSO.
The sum of all 3 components is the total volume,which is what I am \
trying to fit to the data (see below).*)
Check :
In[40]:= Vtotal[1, 1, 1]
Out[40]= 1.14902*10^6
In[41]:= data = {1.047397 10^6, 1.028567106, 1.010057 10^6, 991863,
974012, 956510, 939369, 922600, 906215, 890226};
In[42]:= ClearAll[Lp, PCPA1];
FindFit[data,
Vtotal[t, Lp, PCPA1], {{Lp, 2.64*^-27}, {PCPA1, 4.38*^-31}}, t]
During evaluation of In[42]:= NDSolve::nderr: Error test failure at t
\
== 4.159268807727241`*^-6; unable to continue. >>
Out[43]= {Lp -> 8.55856*10^-28, PCPA1 -> -3.67101*10^-17}
hth
--
V.Astanoff