Advanced Search
All Posts : Product Tips & Techniques

Interesting numeric bug

Schivnorr 526
October 25 2007
0
The student is supposed to use the logistic differential equation given to them (with carrying capacity .52e12), dsolve() the differential equation with a specific initial condition, and then calculate the population at 3 different times. Notice what happens: 

[> restart;
[> my_deq := (diff(P(t), t))/P(t) = .789*(1-P(t)/(.52*10^12));

                   d
                   -- P(t)
                   dt                               -11
         my_deq := ------- = 0.789 - 0.1517307692 10    P(t)
                    P(t)

[> dsolve({my_deq,P(2290)=.83*10^10});

                                                -180681
                    16371750000000000000000 exp(-------)
                                                  100
     P(t) = -----------------------------------------------------
                            -180681                        789 t
            31484134609 exp(-------) + 1941015865391 exp(- -----)
                              100                          1000

[> a1 := evalf(16371750000000000000000*exp(-180681/100)/(31484134609*exp(-180681/100)+1941015865391*exp(-(789/1000)*10)),3);

                                           -772
                             a1 := 0.190 10

[> b1 := evalf[3](16371750000000000000000*exp(-180681/100)/(31484134609*exp(-180681/100)+1941015865391*exp(-(789/1000)*2400)));

                                            12
                              b1 := 0.520 10

[> c1 := evalf[3](16371750000000000000000*exp(-180681/100)/(31484134609*exp(-180681/100)+1941015865391*exp(-(789/1000)*2500)));

                                            12
                              c1 := 0.520 10

[> restart;
[> my_deq := (diff(P(t), t))/P(t) = .789*(1-P(t)/(.52*10^12));

                   d
                   -- P(t)
                   dt                               -11
         my_deq := ------- = 0.789 - 0.1517307692 10    P(t)
                    P(t)

[> dsolve({my_deq,P(0)=.83*10^10});

                          16371750000000000000000
           P(t) = ----------------------------------------
                                                    789 t
                  31484134609 + 1941015865391 exp(- -----)
                                                    1000

[> a2 := evalf[3](16371750000000000000000/(31484134609+1941015865391*exp(-(789/1000)*10)));

                                            12
                              a2 := 0.510 10

[> b2 := evalf[3](16371750000000000000000*exp(-180681/100)/(31484134609*exp(-180681/100)+1941015865391*exp(-(789/1000)*2400)));

                                            12
                              b2 := 0.522 10

[> c2 := evalf[3](16371750000000000000000*exp(-180681/100)/(31484134609*exp(-180681/100)+1941015865391*exp(-(789/1000)*2500)));

                                            12
                              c2 := 0.522 10

[>
Although a1 is incorrect, b1 and c1 are correct. Then, using the exact same commands to calculate b2 and c2, they are incorrect. Any idea why the my_deq := or dsolve() lines cause this numeric difference (if you take those lines out, b1 c1 b2 and c2 all take on the value .522e12)? Re-running these lines repeatedly, sometimes b1 and c1 will become .522e12, but I have not yet seen b2 or c2 become .520e12. --Schivnorr

Please Wait...