# DSolve for a huge linear inhomogeneous ODE in parallel

I have a linear inhomogeneous ODE with constant coefficients, which I need to solve symbolically. The problem is that the inhomogeneity is sum of more than 100 terms (however, each term itself is pretty nice for integration). During strightforward computation with DSolve Mathematica runs out of memory.

Since the problem is linear, it's possible to solve equation for each term of inhomogeneity separately (in parallel, in perfect) and sum the results. Do you have any suggestion how can I do this? I think I need to distinguish homogeneous and inhomogeneous parts of the equation and split the last one into inhomogeneities with less terms.

Thank you,
Ivan

EDIT

Sorry for lack of example. A trivial casee is:

(* Construct inhomogeneous terms *)
inhom = Plus @@ Table[a[i] Exp[b[i] x], {i, 1, 100}];
(* Construct ODE and boundary conditions *)
eq = {u''[x] + u[x] + inhom == 0, u[0] == 0, u[1] == 0};
(* Solve ODE *)
DSolve[eq, u, x]


In reality coefficients a[i] and b[i] are more complicated constants.

EDIT

With the help of Michael's answer I managed to solve the problem:

(* Function which returns homogeneous and inhomogeneous parts of equation *)
splitEq[eq_Equal, vars_] := With[{oneSide = Subtract @@ eq},
If[Head[oneSide] =!= Plus,(* If the equation consist of the only term *)
{oneSide, 0},
With[{inhom = Select[oneSide, FreeQ[_?(MemberQ[Flatten[{vars}], #] &)]]},
{oneSide - inhom, -inhom}]
]]

(* Solve equation for each term of inhomogenity in parallel *)
ParallelDSolve[eqn_, bcs_, vars_, x_] := Block[{eq, hom, inhom, sol},
{hom, inhom} = splitEq[eqn, vars];
eq[rhs_] := {hom == rhs}~Join~bcs;
If[Head[inhom] =!= Plus, (* If inhomogenity consist of the only term *)
sol = u -> DSolveValue[eq[inhom], vars, x],
sol = With[{body = Through[ParallelMap[DSolveValue[eq[#], vars, x] &, inhom][x]]},
u -> Function[Evaluate[{x}], body]]];
sol]

(* -------------------------- *)

(* Construct inhomogeneous terms and equation*)
inhom = Sum[a[i] Exp[b[i] x], {i, 1, 10}];
eq = u''[x] + u[x] + inhom == 0

(* Solve equation *)
ParallelDSolve[eq, {u[0] == 0, u[1] == 0}, u, x]


This example takes 0.086 sec on my machine. When straightforward DSolve[{eq, u[0] == 0, u[1] == 0}, u, x] runs out of memory after a few minutes of work.

Thanks to Michael.

• Hi ! You can add a simple example, so people would not waste time searching and typing-in equations. Mar 28 '15 at 10:12

Perhaps this:

inhom = Plus @@ Table[a[i] Exp[b[i] x], {i, 1, 100}];
eq[inhom_] := {u''[x] + u[x] + inhom == 0, u[0] == 0, u[1] == 0};

sol = u -> ParallelMap[DSolveValue[eq[#], u[x], x] &, inhom]

(*
u -> -(1/(1 + b[1]^2))
a[1] (-Cos[x] + E^(x b[1]) Cos[x]^2 + Cot[1] Sin[x] -
E^b[1] Cos[1] Cot[1] Sin[x] - E^b[1] Sin[1] Sin[x] +
E^(x b[1]) Sin[x]^2) -
... -
(1/(1 + b[100]^2))
a[100] (-Cos[x] + E^(x b[100]) Cos[x]^2 + Cot[1] Sin[x] -
E^b[100] Cos[1] Cot[1] Sin[x] - E^b[100] Sin[1] Sin[x] +
E^(x b[100]) Sin[x]^2)
*)


One could also use ParallelSum or ParallelTable instead of mapping over Plus.

• Thank you! It works! One ramark: in reallity I don't have inhom and u''[x] + u[x], but they are mixed in some eq. To use your approach I need to distinguish terms, involving unknown function, and inhomogeneous term. Do you have any idea how to do it for arbitrary equation eq?
– Ivan
Mar 28 '15 at 13:22
• Does this work?: inhom = Cases[eq2 /. Equal -> Subtract, term_ /; FreeQ[term, u]], where eq2 is the equation without initial conditions, e.g. eq2 = First[eq]. Mar 28 '15 at 14:34
• That's a good idea. I have an improvement to split the equation into homogeneous and inhomogeneous parts: splitEq[Equal[lhs_, rhs_], var_] := With[{splitted = (Plus @@@ SplitBy[List @@ #, FreeQ[var]]) & /@ {lhs, rhs}}, (Subtract @@ splitted) /. {a_, b_} :> {b, -a}]. One can use it in the following way: {hom, inhom} = splitEq[eq2, u], where eq2 is the equation without initial conditions
– Ivan
Mar 28 '15 at 15:06