6

I'm new to mathematica and was wondering if someone could help me. I have a differential equation which is:

$\frac{dy}{dx} = 1/2 y(1 - y)(a - y + a y)$

It has initial condition $y(0) = y0$ and it's assumed $0 < a < 1/2$.

So I've been trying to figure out how to use DSolve. I know the solution $y(x)$ is meant to fall in the range $0 < y(x) < 1$ but DSolve keeps giving me complex solutions. I'm assuming this is because I haven't provided it any assumptions. Here's what I've tried:

DSolve[{y'[x] == 
   1/2 (1 - y[x]) y[x] (a - y[x] + a y[x]), y[0] == y0}, y[x], x]

Which results in:

{{y -> Function[{x}, 
    InverseFunction[(
       a Log[1 - #1] + (1 - 2 a) Log[#1] + (-1 + a) Log[
          a - #1 + a #1])/(a (-1 + 2 a)) &][-(x/2) + (
      a Log[1 - y0] + Log[y0] - 2 a Log[y0] - Log[a - y0 + a y0] + 
       a Log[a - y0 + a y0])/(a (-1 + 2 a))]]}}

There are a few questions I was wondering someone could help me with:

  • Is it possible to pass assumptions to DSolve, e.g. to tell it $0 < a < 1/2$ and $0 \leq y \leq 1$? If not how can I handle this in Mathematica
  • What is the correct ways of assigning the result of DSolve to a function whereby I can pass values of $a$, $x_0$ and $x$, e.g. $y(x, x_0, x)$?
Greg Hurst
  • 35,921
  • 1
  • 90
  • 136
user3277112
  • 63
  • 1
  • 3
  • I can add this about why M can't solve the InverseFunction problem: It's similar to http://mathematica.stackexchange.com/q/56771, but the exponential form of your equation has terms of the form y^2 and y^a. For a generic a, that's hard to solve, unless you're lucky. – Michael E2 Apr 26 '15 at 19:44

3 Answers3

5

To create a function of the type you ask for, you can do the following. Note the use of Set and DSolveValue.

Clear[a, y0, x]; 
y[a_, y0_, x_] = 
 DSolveValue[{y'[x] == 1/2 (1 - y[x]) y[x] (a - y[x] + a y[x]), 
   y[0] == y0}, y[x], x];

Plot[Evaluate@Table[y[a, 0.1, x], {a, 0.1, 0.4, 0.1}], {x, 0, 50}]

Mathematica graphics

Instead of DSolveValue one can also use y[x] /. DSolve[<ode>, y[x], x].


One can pass assumptions to the functions used by DSolve via $Assumptions, which can be temporarily set with Assuming. That sometimes helps, but not in this case.

Assuming[0 < a < 1/2 && 0 < y[x] < 1 && 0 < y0 < 1, 
 DSolveValue[{y'[x] == 1/2 (1 - y[x]) y[x] (a - y[x] + a y[x]), 
   y[0] == y0}, y[x], x]
 ]

It does not help in this case because the equation that implicitly defines the solution y[x] cannot be solved.

Assuming[0 < a < 1/2 && 0 <= y[x] <= 1 && 0 <= y0 <= 1, 
 eq = DSolveValue[{y'[x] == 1/2 (1 - y[x]) y[x] (a - y[x] + a y[x]), 
      y[1] == y0}, y[x], x] /. 
    InverseFunction[f_][f0_] :> f[y] == f0 // Simplify
 ]
(*  a + 2 a^2 x + 2 a Log[1 - y] + (2 - 4 a) Log[y] + 
      2 a Log[a - y + a y] + 4 a Log[y0] + 2 Log[a - y0 + a y0] == 
    2 a^2 + a x + 2 Log[a - y + a y] + 2 a Log[1 - y0] + 2 Log[y0] + 
      2 a Log[a - y0 + a y0]  *)

Assuming[0 < a < 1/2 && 0 <= y <= 1 && 0 <= y0 <= 1,
 Solve[eq == 0, y]
 ]

Solve::nsmet: This system cannot be solved with the methods available to Solve. >>

Michael E2
  • 235,386
  • 17
  • 334
  • 747
1

Another option for you is to use ParametricNDSolve, which numerically solves differential equations with one or more parameters.

soln = ParametricNDSolveValue[
 {y'[x] == 1/2 (1 - y[x]) y[x] (a - y[x] + a y[x]), y[0] == y0}, (* your equations*)
 y, (* the expression you want to be returned*)
 {x, 0, 100}, (* the independent variable and range to solve for *)
 {{a, 0, 0.5}, {y0, 0, 1}} (* the parameters with their ranges *)
]

ParametricNDSolveValue returns a ParametricFunction. You can call this function with a set of parameters to get y as a function of x:

soln[0.2, 0.6]  (* returns InterpolatingFunction[{{0., 100.}}, <>] *)

Note that the parameters are input in the same order as they were input to ParametricNDSolve, in this case [a, y0]. The InterpolatingFunction will allow you to find the value of y at a specific value of x:

soln[0.2, 0.6][10] (* === 0.381214 *)

Although your given system does seem to have an analytic solution, this method will be useful for systems that require numerical solution.

2012rcampion
  • 7,851
  • 25
  • 44
0

One way you can let Mathematica know what you are assuming here is

Simplify[DSolve[{y'[x]==1/2 (1-y[x]) y[x] (a-y[x]+a y[x]), y[0]==y0}, y[x], x], 0<a<1/2]

One way of assigning the result and passing value of a, y0, etc is

yf = y[x]/.DSolve[{y'[x]==1/2 (1-y[x]) y[x] (a-y[x]+a y[x]), y[0]==y0}, y[x], x][[1]];
Simplify[yf/.{a->2, y0->1/2, x->4}]

That does warn you that inverse functions are being used and you should check the results carefully, but the result from DSolve is an inverse function, so I don't think that can be avoided

Bill
  • 12,001
  • 12
  • 13
  • @2012rcampion DSolve calls functions (such as Integrate and Simplify) which do use $Assumptions. So Assuming can make a difference. See http://mathematica.stackexchange.com/a/65946 or http://mathematica.stackexchange.com/a/75243, for instance. – Michael E2 Apr 26 '15 at 17:15
  • @MichaelE2 I didn't know that... that's actually quite useful. – 2012rcampion Apr 26 '15 at 19:14
  • 1
    @2012rcampion It can be useful. So can Block[{Simplify = FullSimplify}, DSolve[ode, y, x]]. But DSolve works well enough on a generically wide range of problems that these tricks do not help that often. – Michael E2 Apr 26 '15 at 19:19