Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I have a second order reaction / diffusion type ODE of the form

$\frac{D_{o}}{r^2} \frac{d}{dr}\left(r^2 \frac{dC}{dr} \right) - \frac{aC}{C+k} = 0$

where $a, k$ and $D_{o}$ are constants and $C$ is a function of $r$. I know the following;

$C(r_{o}) = \omega$ , $C(r_{n}) = 0$ and $C'(r_{n}) = 0$ , where $r_{o}$ and $r_{n}$ are positions and $\omega$ is a known constant. To model this behaviour, I am using NDSolve, which takes two conditions for these types of ODEs. For example, this gives a potential curve between the boundaries $r_{n}$ and $r_{o}$ ;

eqn = ((Do2/(r^2))*D[(r^2*(O2'[r])), r] - (a)*((O2[r])/(O2[r] + k)));

s = NDSolve[{eqn == 0, O2[rn] == 0, O2[ro] == omega}, O2, {r, rn, ro}]

This gives me one curve which I can readily plot; however is I instead tell NDSolve to use different conditions like this -

s = NDSolve[{eqn == 0, O2'[rn] == 0, O2[ro] == omega}, O2, {r, rn, ro}]

I get a different curve. However, I suspect I should be getting the same curve regardless. If the problem is narrowing it down to a unique solution from an array of possible solutions, is there a way to make Mathematica take all three BCs? I'm also open to the possibility that there is a more fundamental mathematical error or perhaps a bad assumption on my part, so please by all means point it out to me if I've missed it . Thanks in advance - grateful for any advice on this!

share|improve this question
Ah yes - a valid point. I actually use different names for my variables and the function, my mistake - was merely trying to create a M.W.E! –  DRG Sep 9 '13 at 16:24
You have effectively a second order ordinary differential equation. Thus, 2 boundary conditions should be taken, rather than 3. The third is the extra one. I can imagine 2 situations where more than 2 boundary (or initial) conditions apply. The first is the case of a degeneration. It is not the case of the equation at hand. The second is the non-linear eigenvalue problem. That is, all the three conditions work, but only for a specific value of one of the parameters. Then this value should be found in the course of solution. Is it your case? –  Alexei Boulbitch Sep 10 '13 at 8:10
It could well be, and I was wondering the same myself - I'll play around with the mathematics more and see if I can work out the root cause. Thanks, Alexei –  DRG Sep 10 '13 at 15:37
are D0 and k "really" constants? Meybe that's the key. –  Kuba Sep 10 '13 at 16:23
Perhaps - I suspect the value of rn might be key; I can get both to match completely if rn ~= 330 um. However, I have some evidence that rn = 342 um, and if I use that they don't match completely. Maybe rn affects the model greatly? –  DRG Sep 10 '13 at 16:42

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.