Ask Your Question

Revision history [back]

click to hide/show revision 1
initial version

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, $\pm 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/2

Related question, answers and discussion: logistic differential equation