Message Boards Message Boards

GROUPS:

Solve two simultaneous differential equations using W|A?

Posted 5 months ago
763 Views
|
9 Replies
|
0 Total Likes
|

In need a solution of the two simultaneous differential equations:

dy/dt = 0.2 y(t) f(t) – 0.05 y(t) 
df/dt = - 0.05 y(t) f(t)                                      
with y(0) = 0.00001 and f(0) = 0.9999

I tried the Euler method provided by WolframAlpha Pro Premium: Input :

use Euler method dy/dt = 0.2 y(t) f(t) – 0.05 y(t), df/dt = - 0.05 y(t) f(t), y(0) = 0.000012, f(0) = 1.0, from 1 to 180  

output : Interpreting as:

dy/dt = 0.2 y(t) f(t) - 0.05 y(t)

The program accepts one differential equation only and a solution is provided in the form of an integral equation. The same happened with the Runge-Kutta method. Can anybody help me to find the error ? Is the syntax wrong ? Or is there a principal missunderstanding ? I am sorry, I just started to use Wolfram. Ernst

9 Replies

May be that an analytical solution exists, or perhaps not. But a numeric approach is possible

lsg = NDSolve[
   {D[y[t], t] == 0.2 f[t] y[t] - 0.05 y[t],
    D[f[t], t] == - 0.05 y[t] f[t],
    y[0] == 0.00001,
    f[0] == 0.9999},
   {y, f},
   {t, 0, 250}];
fy = y /. Flatten[lsg];
ff = f /. Flatten[lsg];
Plot[{fy[t], ff[t]}, {t, 0, 250}]

Dear Mr. Dolhaine, thank you very much for your effort to address the problem. I have tried to run your program, however, unfortunately without success. I wonder wether I am using a suitable compiler, because often the comment appears: WolframAlpha does not understand your query. During my professional life I used Fortran IV, but times have changed, and I am now an absolute beginner. I would be most grateful if you would consider the query again. Sincerely Yours Ernst Hoinkis

Quite possibly the Mathematica command is too long for Wolfram|Alpha. In which case you might want to try it in the free version of the Wolfram Cloud (a simple web search should suffice to locate this resource).

Hallo Daniel, Hans's Code runs perfectly with Mathematica. Thank you very much and best regards ! Ernst

Posted 5 months ago

Hi Daniel. Referring to your post: "In which case you might want to try it in the free version of the Wolfram Cloud", is there a free version of the Wolfram Cloud(Mathematica)? I could only get a 14-days trial.

Check here for more details: https://www.wolfram.com/cloud/ It has a button "Sign up for free access".

Alternatively, go here: https://www.wolfram.com/programming-lab/

Then use the "Start programming now (no sign-in required)" button. That should take you to: https://www.open.wolframcloud.com/env/wpl/GetStarted.nb?funnel=WPLExplorations#sidebar=explorations

Now position your curson below the "2+2", enter input, and use <shift-enter> to evaluate it.

Hello Ernst,

I was just typing an answer to you when Daniel's post appeared. I wanted to say the same.

And as far as I know you can download a (free) test-version of Mathematica which runs a month or so.

Hallo Herr Dolhaine ( oder einfach Hans), Ihr Programm läuft mit Mathematica. Herzlichen Dank ! Bleiben Sie gesund und Grüße aus Berlin ! Ernst Hoinkis

I found it interesting that with some algebraic manipulations the system of the two differential equations can be decoupled to give an expression for y as a function of f.

yy = (a (-f + f0))/b + y0 + Log[f] - Log[f0]

(Note : you should run the code given above to have all symbols available).

With appropriate modifications this is indeed the same as fy given above, or better the difference is reasonably well near to zero:

Plot[fy[0] - .2/.05 (ff[t] - ff[0]) + Log[ff[t]/ff[0]] - fy[t], {t, 0, 240},
 PlotRange -> All]

Now you can insert this y in the 2nd equation and obtain an equation for f

fp1 = - b yy f // Simplify
dgl = fF'[t] == (fp1 /. {f -> fF[t], f0 -> ff[0], y0 -> fy[0], a -> .2,  b -> .05})

Unfortunately this cannot be integrated, but NDSolve does the job

Clear[fF]
fF = fF /. Flatten[NDSolve[ {dgl , fF[0] == .9999}, fF, {t, 0, 200}]]

fF is equivalent to ff

Plot[{ff[tt], fF[tt]}, {tt, 0, 150},
 PlotStyle -> {Blue, {Thick, Dashing[{.05, .05}], Red}}]

and we can insert it into the expression for y to obtain

fY = (yy /. {a -> .2, b -> .05, y0 -> fy[0], f0 -> ff[0], f -> fF[t]})

which as well is equivalent to fy obtained earlier :

Plot[{fY, fy[t]}, {t, 0, 180},
 PlotStyle -> {Blue, {Thick, Dashing[{.05, .05}], Red}}]

If there is the need to operate with some more explicit functions one could fit the derivative of ff (or fF) to some appropriate functions. Here I chose a Gaussian and a Lorentzian

data = Table[{t, -D[ff[x], x] /. x -> t}, {t, 50, 150}];
lsg = FindFit[
  data,
  p1 Exp[-p2 (t - p3)^2] + p4/(1 + p5 (t - p6)^2),
  {{p1, .015}, {p2, .005}, {p3, 84}, {p4, .015}, {p5, .005}, {p6, 
    84.5}},
  t]

giving

empfp = p1 Exp[-p2 (t - p3)^2] + p4/(1 + p5 (t - p6)^2) /. lsg

which looks quite reasonable :

Plot[empfp, {t, 40, 120}, Epilog -> Point /@ data]

This can be integrated to give an approximation for ff :

SetOptions[Integrate,GenerateConditions->False];
empff = Integrate[-empfp, {t, 0, x}] + 0.9999
Plot[empff, {x, 0, 140}, 
Epilog -> Point /@ Table[{t, ff[t]}, {t, 20, 140, 5}]]

Insert this in the first differential equation and integrate to get fy

empfy = Exp[  Integrate[(a empff - b /. a -> .2 /. b -> .05), {x, 0, t}]]

The leading factor must be adjusted a bit to give a nice result

Plot[1.282 10^-5 empfy, {t, 0, 250},
 Epilog -> Point /@ Table[{t, fy[t]}, {t, 0, 200, 2}],
 PlotRange -> All]
Attachments:
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract