# How to take the derivative of a NDSolve solution?

I have solved and plotted a system of differential equations, shown below, using Mathematica's NDSolve function. Now I would like to take the derivative of one the resulting interpolating functions and evaluate it at a point in the domain.

Pbub = 1.00002
kl = 2.09015;
ccl20 = 0.00206957
cc2h4eq = 0.000301892
Dc2h4 = 1.48*10^(-5);
Dcl2 = 2.90*10^(-9);
t = 183.0 + 273.15;
k1 = 11493.0*Exp[-2156.58/t];
k2 = 8.517*10^9*Exp[-7282.21/t];
\[Delta] = Dc2h4/kl;
eqn1 = Dc2h4*aa''[x] == k1*aa[x]*b[x] + k2*aa[x]*b[x]^2;
eqn2 = Dcl2*b''[x] == k1*aa[x]*b[x] + 2*k2*aa[x]*b[x]^2;
eqn3 = Dcl2*c''[x] == -k1*aa[x]*b[x];
eqn4 = Dcl2*d''[x] == -2*k2*aa[x]*b[x]^2;
eqn5 = Dcl2*ee''[x] == -2*k2*aa[x]*b[x]^2;
bc1 = aa == cc2h4eq;
bc2 = aa[\[Delta]] == 0;
bc3 = b' == 0;
bc4 = b[\[Delta]] == ccl20;
bc5 = c' == 0;
bc6 = c[\[Delta]] == 1.253/98.95;
bc7 = d == 0;
bc8 = d[\[Delta]] == 0.0001;
bc9 = ee == 0;
bc10 = ee[\[Delta]] == 0.0001;

soln = NDSolve[{eqn1, eqn2, eqn3, eqn4, eqn5, bc1, bc2, bc3, bc4, bc5,
bc6, bc7, bc8, bc9, bc10}, {aa[x], b[x], c[x], d[x], ee[x]}, {x,
0, 4*\[Delta]}]


Specifically, I would like to determine: aa'. I've tried various methods using assignment operators, evaluate, etc ... However, none of them have worked thus far. How should I proceed?

• Does the above code run without errors for you? I get error NDSolve::ndsv: Cannot find starting value for the variable aa^\[Prime]. You should make sure code posted runs without errors. You do not have Pbub defined. – Nasser Nov 16 '17 at 6:22
• I left out code that defined a few of the required parameters. The post has now been edited to include them. – Tunk Nov 16 '17 at 6:28
• You can apply the differential operator D (see documentation) also to InterpolatingFunctions. – Henrik Schumacher Nov 16 '17 at 6:31
• @HenrikSchumacher I attempted to use D, however, I could not get it to work. – Tunk Nov 16 '17 at 6:34
• Include aa' in your list of functions to be solved for. – Carl Woll Nov 16 '17 at 6:46

First, better to use aa and not aa[x]

soln=NDSolve[{eqn1,eqn2,eqn3,eqn4,eqn5,bc1,bc2,bc3,bc4,bc5,bc6,bc7,bc8,bc9,bc10},
{aa,b,c,d,ee},{x,0,4*\[Delta]}]


Now

Plot[Evaluate[aa[x]/.soln],{x,0,1}] The derivative is D[aa/.soln] And

D[aa[x]/.soln,x]/.x->0 • Thanks a ton! What is the problem with regards to using aa[x] versus aa? – Tunk Nov 16 '17 at 6:39
• @Jay there is really not a problem per say, but with NDSolve it is much easier to use aa instead of aa[x] in the call (vs. with DSolve where aa[x] would be the better choice there. It just makes post processing easier of the interpolation functions, and most of the examples in help are also this way. – Nasser Nov 16 '17 at 6:40