# NDSolve: Couple PDE and ODE involving integral

I am trying to numerically solve an equation with NDSolve, where there is a ODE coupled to a PDE, like the following:

    NDSolve[{
Derivative[1, 0][u][t,x] ==
-0.8*u[t, x] - (5.*n[t]*Derivative[0, 1][u][t, x])/(5. + n[t]) + (50.*n[t]^2*Derivative[0, 2][u][t, x])/(5. + n[t])^2,
Derivative[n][t] ==
100. - 0.8*n[t] - (0.5*Integrate[u[t, x], {x, 0, 2*Pi}]*n[t])/(5. + n[t]),
u[0, x] == 200.,
n == 50.,
u[t, 0] == u[t, 2*Pi]},
{u, n},
{t, 0., 2.},
{x, 0., 6.28}]


Unfortunately Mathematica tells me "NDSolve::ndode: Input is not an ordinary differential equation".

When i change n[t] to n[t,x] it calculates something, but n[t,x] doesn't stay uniform over time, which it should because n wasn't a function of x in the first place.

Would anyone know a way around this?

• Can't you just integrate your n function by hand and explicitly add it to the equation?
– gpap
Apr 15, 2014 at 9:55
• The actual problem is a little bit more complicated, so there is no formal solution for n[t]
– Lufu
Apr 15, 2014 at 9:57
• Can you post a minimal example capturing the intricacies of the original problem then?
– gpap
Apr 15, 2014 at 9:58
• Sorry I can't make sense of this as there are a bunch of brackets missing, all the variables have been switched to [Theta], [Omega] in place of \[Theta], \[Omega] and a constant isn't defined. You could edit the question with the full problem - that would help people get interested in this.
– gpap
Apr 15, 2014 at 10:17
• NDSolve[{Derivative[1, 0][u][t,x] == -0.8*u[t, x] - (5.*n[t]*Derivative[0, 1][u][t, x])/ (5. + n[t]) + (50.*n[t]^2*Derivative[0, 2][u][t,x])/(5. + n[t])^2, Derivative[n][t] == 100. - 0.8*n[t] - (0.5*Integrate[u[t,x], {x, 0, 2*Pi}]*n[t])/(5. + n[t]), u[0,x] == 200., n == 50., u[t, 0] == u[t, 2*Pi]}, {u, n}, {t, 0., 2.}, {x, 0., 6.28}] This should be correct now.
– Lufu
Apr 15, 2014 at 11:39

As mentioned in the comment above, NDSolve doesn't know how to disretize your system, so we do it ourselves. I'll use pdetoode for the disretization of the derivative term and Gaussian quadrature formula for the discretization of integral:

With[{u = u[t, x], n = n[t]},
eq@1 = D[u, t] == -0.8 u - (5. n D[u, x])/(5. + n) + (50. n^2 D[u, x, x])/(5. + n)^2;
eq@2 = D[n, t] == 100. - 0.8 n - (0.5 (Integrate[u, {x, 0, 2 Pi}]) n)/(5. + n);
{ic@1, ic@2} = {u == 200, n == 5} /. t -> 0];

points = 25; difforder = 4; domain = {0, 2 Pi};

{nodes, weights} = Most[NIntegrateGaussRuleData[points, MachinePrecision]];
midgrid = Rescale[nodes, {0, 1}, domain];

intrule = HoldPattern@
Integrate[__] -> -Subtract @@ domain weights.Map[u[#][t] &, midgrid];

grid = Flatten[{domain[], midgrid, domain[[-1]]}];

(*Definition of pdetoode isn't included in this post,
ptoofunc = pdetoode[u[t, x], t, grid, difforder, True];
ode@1 = ptoofunc@eq@1;
ode@2 = eq@2 /. intrule;
{odeic@1, odeic@2} = {ptoofunc@ic@1, ic@2};

tend = 2;
{usollst, nsol} =
NDSolveValue[{ode /@ {1, 2}, odeic /@ {1, 2}}, {u /@ grid, n}, {t, 0, tend}];
usol = rebuild[usollst, grid];


Notice the periodic b.c. is set with the 5th argument of pdetoode.

Plot3D[usol[t, x], {t, 0, tend}, {x, ##}] & @@ domain Plot[nsol[t], {t, 0, tend}, PlotRange -> All]
` 