MathGroup Archive 2010

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Using "If" and "NDSolve" together

  • To: mathgroup at smc.vnet.net
  • Subject: [mg106839] Re: [mg106781] Using "If" and "NDSolve" together
  • From: "David Park" <djmpark at comcast.net>
  • Date: Sun, 24 Jan 2010 05:45:07 -0500 (EST)
  • References: <11235685.1264251844444.JavaMail.root@n11>

Maxim 1: Develop calculations in steps, making sure one step works before
proceeding to the next. And specifically don't nest a complex calculation as
the first argument of a Plot statement!

Maxim 2: Start your variables with small case letters to avoid conflict with
predefined Mathematica symbols.

Step 1. Define your data so that it can be easily used and changed if
necessary.

(data = {k0 -> 10^-2, k1 -> 10^-5, k2 -> 10^-6, k3 -> 10^1}) // Column  

Step 2. Write a routine to generate your differential equations and check
that it is working properly. Here I use Piecewise as a t function because I
think it will merge with differential equation solving better than an If
statement.

Clear[a, b, c]
deqns[first_, 
  last_] := {a'[t] == 
    k0 Piecewise[{{0, t <= first}, {1, first < t < last}, {0, 
         t > last}}] - k1 a[t], b'[t] == k1 a[t] - k2 b[t], 
   c'[t] == k2 b[t] - k3 c[t], a[0] == 2 10^-6, b[0] == 0, 
   c[0] == 0} /. Data  

deqns[100, 200] // Column  

Step 3. Solve the equations and define the a, b and c functions.

abcsols = First@NDSolve[deqns[100, 200], {a, b, c}, {t, 0, 1000}];
{a[t_], b[t_], c[t_]} = {a[t], b[t], c[t]} /. abcsols

Step 4. Plot the solutions. Here I scaled the functions so they would all
fit on one plot. 

Plot[{a[t], 1000 b[t], 1/2 10^10 c[t]}, {t, 0, 1000},
 Frame -> True,
 Axes -> False]  

So maybe everything is working. The plot could be improved with labeling 


David Park
djmpark at comcast.net
http://home.comcast.net/~djmpark/  


From: Benfeitas [mailto:rui.benfeitas at gmail.com] 

Hello!

I have a tricky job to do: I want to simulate something, and want to
assign a value to a variable under certain conditions. Those are, if
time is within a certain interval, I want that variable to have that
value. Otherwise, I want it to have other value.

Until now, I have been trying to use the following command:

Funcao[first_, last_] :=
 Plot[A[t] /.
   NDSolve[
          {A'[t] == k0 A0 - k1 A[t], B'[t] == k1 A[t] - k2 B[t],
          C'[t] == k2 B[t] - k3 C[t], A[0] == 2 10^-6, B[0] == 0,
          C[0] == 0} /. k0 -> 10^-2 /. k1 -> 10^-5 /. k2 -> 10^-6 /.
          k3 -> 10^1 /. If[first < t < last, A0 -> 1, 0], {A, B, C},
{t,
     0, 1000}], {t, 0, 1000}]



for some reason, it is not working, and I cannot figure out why...

With that command I will want to plot A[t], for t->{0,1000}, and
considering that A0->1 if first<t<last, and A0->0 if t<first or
t>last. That way, A[t] should be higher if first<t<last.

Can you guys please give me some tips? Thanks




  • Prev by Date: Re: Function and object naming conventions
  • Next by Date: Re: Journals dying?, apparently rather slowly (was , I->-I)
  • Previous by thread: Re: Using "If" and "NDSolve" together
  • Next by thread: Re: Using "If" and "NDSolve" together