3
$\begingroup$

I have 2 differential equations with 2 variables, x and y,which are a function of t and I have the parameters k1, k2 y k3.

                       dx/dt=-k1 x2+ k2 x y
                      dy/dt=k1 x2-k2 x y- k3 y

I have to adjust the equations to the following experimental data

                       xo=70.26, x(t=720)=45.78
                       xo=71.04, x(t=720)=46.32
                       xo=37.23, x(t=720)=24.67
                       xo=37.91, x(t=720)=28.78

I tried FindFit and NMinimize. The problem with FindFit is that I have multiple initial conditions. Then I used NMinimize, I tried to create a function error only with the first data (eventually, I will use the rest of the data) but NMinimize gives me the following error This is the funtion.

Remove["Global`*"] 
f[k1_?NumericQ,k2_?NumericQ,k3_?NumericQ]:=
     Module[{x,y,out},out=Abs[45.78-x[720] /.
           NDSolve[{x'[t]==-k1 x[t]^2+k2 x[t]y[t],y'[t]==k1 x[t]^2-k2 x[t]y[t]-k3 y[t],
                    x[0]==70.26,y[0]==0},{x,y},{t,0,800}]]]
res=NMinimize[{f},{k1,k2,k3}]

And this is the error

NMinimize::nnum: The function value f is not a number at {k1,k2,k3} = {8.17269,8.09533,4.9417}. >>

If anybody can help

$\endgroup$
4
  • $\begingroup$ ParametricNDSolve[] is what seems to be used now, if you have version 9. $\endgroup$ Commented May 11, 2013 at 2:08
  • $\begingroup$ Might the problem be caused by providing NMinimize with the list {f} as argument instead of just f? $\endgroup$ Commented May 11, 2013 at 6:55
  • $\begingroup$ First, the output form of NDSolve is nested list (Such as: {{x->InterpolatingFunction[{{0.,30.}},<>]}}). So you need to add such as First@, to the NDSolve to get value not {value}. Second, you need to write the function explicitly(i.e.f[k1,k2,k3] not f,) in the NMinimize. Although after these two steps, you'll still get a bunch of warnings. $\endgroup$
    – luyuwuli
    Commented May 11, 2013 at 10:07
  • $\begingroup$ I tried by using ParametricNDSolve and FindRoot on simpler situations; It works. But I can't find appropriate initial values of k1,k2,k3 (and y0), so Mathematica returns a bunch of warning too. I think a good constrains of k1,k2,k3 and y0 is needed. (BTW, the original equation is wrong with x^2 not x2) $\endgroup$
    – luyuwuli
    Commented May 11, 2013 at 10:14

1 Answer 1

2
$\begingroup$

Could be interpreted as a nonlinear least squares problem.

pts = {{70.26, 45.78}, {71.04, 46.32}, {37.23, 24.67}, {37.91, 
28.78}};

ndsoln[{k1_, k2_, k3_}, x0_] := x[720] /. 
  NDSolve[{x'[t] == -k1 x[t]^2 + k2 x[t] y[t], 
  y'[t] == k1 x[t]^2 - k2 x[t] y[t] - k3 y[t], x[0] == x0, 
  y[0] == 0}, {x, y}, {t, 0, 720}][[1]]

ssfun[{k1_?NumericQ, k2_?NumericQ, k3_?NumericQ}] := 
  Sum[(ndsoln[{k1, k2, k3}, pts[[i, 1]]] - pts[[i, 2]])^2, 
  {i, Length[pts]}]

scaledArgs = {10^-5 k1, k2, 10^2 k3};

fm = FindMinimum[ssfun[scaledArgs], {k1, 3, 4}, {k2, 20, 21}, {k3, 5, 6}]

Show[Plot[ndsoln[scaledArgs /. fm[[2]], x0], {x0, 35, 75}], 
 ListPlot[pts, PlotStyle -> Directive[ColorData[1][2], PointSize[Medium]]], 
 PlotRange -> {20, 50}, Frame -> False, Axes -> True, 
 AxesLabel -> {"x(0)", "x(720)"}, AxesOrigin -> {35, 20}]

The result

$\endgroup$
1
  • $\begingroup$ Thank you very much! It was so helpful. $\endgroup$
    – Mai
    Commented May 14, 2013 at 0:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.