How to solve this logistic differential in sage?
y=function('y')(x)
(b,r,y0)=(30,0.8,1.5)
de=diff(y,x)==r*y*((1-y)/b) # logistic equation
soln=desolve(de,y,ics=[0,y0])mathematica solution is: y(x)=(b*exp(r*x)*y0)/(exp(r*x)*y0-y0+b)
Other questions/answers found on the logistic topic were not helpful for this problem
The posted code:
( b, r, y0 ) = ( 30, 8/10, 3/2 )
var( 'x' );
y = function('y')(x)
de = diff( y, x ) == r*y*(1-y)/b
sol = desolve( de, y, ics=[0,y0] )already leads to the solution (using implicit functions of y and x):
sage: sol
-75/2*log(y(x) - 1) + 75/2*log(y(x)) == x + 75/2*log(2) + 75/2*log(3/2)In the ode-world this is already a solution. (Although i had a similar situation in an exam, and the examiner considered
i stopped too early, since the equation can be solved explicitly.) Here, the coefficients of the two log expressions involving a non-constant function of y(x) are opposite, ±75/2, so using exp on both sides the expression (y−1)/y=1−1/y leads to an algebraic equation of first order in y, that can be solved. If we want to do this in sage, there is one more step needed, before calling solve (for a combined exponential/logarithmic and algebraic equation), we substitute with an other variable, z, the expression y(x). The adjusted code for this plan and a possible answer to the question is as follows:
( b, r, y0 ) = ( 30, 8/10, 3/2 )
var( 'x,z' );
y = function('y')(x)
de = diff( y, x ) == r*y*(1-y)/b
sol = desolve( de, y, ics=[0,y0] )
sol = sol . subs( { y(x): z } )
Y = sol . solve( z, to_poly_solve=True )[0].rhs()
print "Solution: Y=%s" % Y
print "Check is", bool( diff(Y,x) == r*Y*(1-Y)/ b )
print "Y(0) = %s" % Y.subs( {x:0} )Results:
Solution: Y=3*e^(2/75*x)/(3*e^(2/75*x) - 1)
Check is True
Y(0) = 3/2Related question, answers and discussion: logistic differential equation
let me study this. thank you
re:logistic function If I understand your reply correctly, y(x)=3e^(2/75x)/(3e^(2/75x) - 1) But y(x)=3e^(2/75x)/(3e^(2/75x) - 1) is not a logistic function Plot it out and plot out y(x)=(bexp(rx)y0)/(exp(rx)*y0-y0+b) Again thank you for responding.
The above formula corresponds to...
sage: ( b, r, y0 ) = ( 30, 8/10, 3/2 )
sage: Y = b*exp(r*x)*y0/(exp(r*x)*y0-y0+b)
sage: Y
30*e^(4/5*x)/(e^(4/5*x) + 19)But please observe the difference between the two expressions:
de=diff(y,x)==r*y*((1-y)/b) as posted, and
de=diff(y,x)==r*y*(1-y/b) as leading to the above solution. Let's see if the code fits in place:
( b, r, y0 ) = ( 30, 8/10, 3/2 )
var( 'x,z' );
y = function('y')(x)
de = diff( y, x ) == r*y*(1-y/b)
sol = desolve( de, y, ics=[0,y0] )
Y = sol . subs( { y(x): z } ) . solve( z, to_poly_solve=True )[0].rhs()
print "Solution: Y=%s" % Ywhich gives indeed:
Solution: Y=30*e^(4/5*x)/(e^(4/5*x) + 19)dan_fulea: Thank you. Solved my problem/question. I am a senior (very) not a student. Math is my hobby. My mathematica version does not run on a 64 bit machine and therefore, switching to sage and learning on my own with the help of a few texts books. I have not seen/learned the code "sol.subs({y(x):z}).solve(z,to_poly_solve=True)[0].rhs() in any of my text books. After your help I found with 'solve? evaluate' but no explanation as to why. A bit exotic for my level. The code fixes my problem. Thanks again. כול הכבוד
P.S. I miss-typed my original questions. Meant to write (1-y/b); not (1-y)/b. That error caused some unecessary confusion.
Asked: 7 years ago
Seen: 660 times
Last updated: Sep 16 '17
Copyright Sage, 2010. Some rights reserved under creative commons license. Content on this site is licensed under a Creative Commons Attribution Share Alike 3.0 license.
