2
$\begingroup$

I have data generated from two functions given below. How can I go about finding a fit of the form:

y'[x] == -(a + b/(c + d y[x] + e F[x] ))

where a,b,c,d,e are to be determined?

F[x_] := If[x < 37.5, 2.80 + 0.0036 x, 117.86 - 6.314 x + 0.0865 x^2] s = NDSolve[{y'[x] == -(0.0595 + 3716/(11780 + y[x])) y[x], y[0] == 650000}, y, {x, 0, 60}] data = Table[{x, Evaluate[y[x] /. s], F[x]}, {x, 0, 60, 1}]
Improve this question
$\endgroup$
6
  • $\begingroup$ @Nasser ...error corrected $\endgroup$ Commented Jan 8, 2014 at 6:35
  • $\begingroup$ What are the triplets in data ? Is it {x, y[x], y'[x]} ? $\endgroup$ Commented Jan 8, 2014 at 8:12
  • $\begingroup$ @b.gatessucks It is {x,y[x],F[x]} $\endgroup$ Commented Jan 8, 2014 at 8:26
  • 1
    $\begingroup$ @thils, I believe you have to approximate y'[x] based on the data using finite difference coefficients. Have you considered this? $\endgroup$ Commented Jan 9, 2014 at 9:44
  • 1
    $\begingroup$ @thils, I noticed that the solution from NDSolve was differentiable already; hence my answer :) $\endgroup$ Commented Jan 13, 2014 at 11:32

1 Answer 1

Reset to default
1
$\begingroup$

This seems to be a task for FindFit. First, get the solution from NDSolve as a function; let me call it solY.

F[x_] := If[x < 37.5, 2.80 + 0.0036 x, 117.86 - 6.314 x + 0.0865 x^2]; solution = NDSolve[{y'[x] == -(0.0595 + 3716/(11780 + y[x])) y[x], y[0] == 650000}, y, {x, 0, 60}]; solY = y /. First[First[solution]];

Now, generate the data in a format for FindFit. In this case, as you want y'[x] == -(a + b/(c + d y[x] + e F[x] )), my suggestion is to have two independent values, y[x] and F[x] giving the dependent value y'[x]. x is not required. Note that you can use the derivative solY'!

This data is generated by the following

data = Table[{solY[x], F[x], solY'[x]}, {x, 0, 60, 1}];

Now, feed that into FindFit, like

FindFit[data, -(a + b/(c + d v1 + e v2)), {a, b, c, d, e}, {v1, v2}]

and you get

{a -> 11425.8, b -> -2.63417*10^9, c -> 1.01206*10^7, d -> -98.9974, e -> -213696.}

The derivative solY' comes from an approximation (numerical solution to NDSolve), so I would double check results by other means. Also, your parameters are non-linear and "In the nonlinear case, it finds in general only a locally optimal fit." (see documentation)

Hope it helps.

Update

As mentioned by OP, FindFit does not find a good answer. F is itself problematic and certainly one would like to explore the data. Here some ideas on how to do this, but consider that the parameter space is 5D - so visualization is already hard.

You can use the solution to define a new function. This is a way to do it. Say you capture the solution into resp, like resp = FindFit[...];. Then you can do

g[v1_, v2_] := Evaluate[-(a + b/(c + d v1 + e v2)) /. resp];

Using ?q actually gives g with values from the solution rules.

Mathematica graphics

You can define then a new data to compare.

data2 = Map[({#1, #2, g[#1, #2]} &) @@ # &, data]; Show[ ListPlot3D[data2, ColorFunction -> (Directive[Red, Opacity[0.7]] &)], ListPlot3D[data, ColorFunction -> (Directive[Blue, Opacity[0.7]] &)], PlotRange -> All]

which gives you

Mathematica graphics

Again, hope it helps.

Improve this answer
$\endgroup$
2
  • $\begingroup$ Appear to work fine, except that I got "Failed to converge to the requested accuracy or precision within 100 iterations"...with slightly different set of a,b,c,d,e. Possibly I need to change the model slightly? $\endgroup$ Commented Jan 14, 2014 at 5:20
  • 1
    $\begingroup$ @thils, I am not sure about your problem space. I added few ideas you can use. Good luck! $\endgroup$ Commented Jan 14, 2014 at 8:26

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.