# Could Mathematica solve a differential equation asymptotically?

Is there any possibility that Mathematica could give asymptotic behavior(s) of a differential equation as it independent variable tends to a certain value?

Because I didn't know how to decompose this hard problem, I just tried the naive code to find a trial solution to the 3rd-order nonlinear ODE:

$$H^2(H_{\eta\eta\eta}-\eta^{-2}H_\eta+\eta^{-1}H_{\eta\eta})-\frac{\eta}{10}=0,$$

which is subject to the boundary conditions: $H(0)=1$ and $H_{\eta\eta}(0)=-c$ with $c$ being a constant. If one wants to seek a solution $H(\eta)$ that is an even function in $\eta$.

DSolve[{H[\[Eta]]^2*(H'''[\[Eta]] - \[Eta]^-2 H'[\[Eta]] +
1/\[Eta] H''[\[Eta]]) - \[Eta]/10 == 0, H == 1,
H'' == -c, H[\[Eta]] == H[-\[Eta]]}, H[\[Eta]], \[Eta]]


Mathematica complains that:

DSolve::litarg: To avoid possible ambiguity, the arguments of the dependent variable in {...} should literally match the independent variables.

Obviously, the error results from the even function condition, H[\[Eta]] == H[-\[Eta]]. If I comment this condition, the code just repeats itself without any solution. Another minor question is:

How to impose an even or odd function requirement in DSolve?

Motivation: The reason why I arise this problem is that the asymptotic behavior of $H(\eta)$ has been shown (I guess by perturbation method which is powerful math tool) as follows. That can be used to check the desired code/method with MMA. I really want to consult superiors of MMA here for some suggestion about how to solve an exactly unsolvable differential equation with the asymptotic method. Thank you in advance!

$H \sim 1-\frac{1}{2}c\eta^2 \quad \text{as} \quad \eta \to 0$

$H(\eta)=\tilde{h}F\left(\eta_0-\eta \right)$ near $\eta_0$, where $H(\eta)$ reaches its minmum value for a given constant $c$, denoted by $\tilde{h}=H(\eta_0)$, $F(x)\sim x^2$ as $x\ll0$ and $F(x)\sim x(\log x)^{1/3}$ for large $x\gg0$.

• If H[η] is well behaved near η == 0, then the symmetry requirement can be given by H' == 0. This being the case, it is straightforward to obtain a power series solution of the ODE near η == 0. However, if you want an approximate solution at a value η0 outside the radius of convergence of the power series at η == 0, then I see no way even in principle to apply the η == 0 boundary conditions to the approximate solution at that η0. Jan 27, 2018 at 17:29
• Incidentally, -(3/(4 5^(1/3))) η^(4/3) satisfies the ODE identically but not the boundary condition H == 1. Jan 27, 2018 at 17:31
• @bbgodfrey, Thanks. For comment1: At $\eta=\eta_0$, could we impose these boundary conditions for the ODE: $H_\eta(\eta_0)=0$, $\tilde{h}=H(\eta_0)$ and $H_{\eta\eta}(\eta_0)>0$, since it reaches its min $\tilde{h}$ at $\eta_0$? Jan 28, 2018 at 2:46
• @bbgodfrey, Thanks. For comment2: for your trial solution the boundary condition at $\eta=0$ is $H(0)=0$? Jan 28, 2018 at 2:48
• With respect to your two responses, (1) You can, in principle, impose the boundary conditions at η0, but how do you know what value η0 has? (2) Yes. Jan 28, 2018 at 5:44

It is possible to roll one's own asymptotic solver with the help of Series. As a demonstration, here we show how to obtain an asymptotic power series solution around zero.

The equations in OP, except the even function condition H[η] == H[-η] (which will turned out to be redundant for this particular problem), can be rearranged to (omit the == 0 parts):

eqs = {
H[η]^2*(H'''[η]-η^-2 H'[η]+1/η H''[η])-η/10,
H - 1,
H'' + c
};


If we assume the existence of power expansion about $x=0$:

$$H(x) = \sum_{k=0}^\infty h_k x^k$$

then a series representation of $H$ and its derivatives can be straightforwardly defined through following rules:

seriesRules = RightComposition[
ReplaceAll[{
H[x_] :>
(Inactive[Series][H[η], {η, 0, max}] // Inactive[ReplaceAll][η -> x]),
Derivative[s_Integer?Positive][H][x_] :>
(Inactive[Series][Derivative[s][H][η], {η, 0, max}] // Inactive[ReplaceAll][η -> x])
}],
(* the odd/even function constraint can be described as following rule *)
(* Inactive[ReplaceAll][{
(* odd: *)(* H :> 0,Derivative[s_Integer?EvenQ][H] :> 0 *)
(* even: *)Derivative[s_Integer?OddQ][H] :> 0
}], *)
Inactive[ReplaceAll][{H :> h, Derivative[s_Integer][H] :> h[s]}]
];


Applying it on eqs gives us its series version:

series = eqs // seriesRules;


For a given series order max, series can be Activated to become algebra equations serieseqs about $h_0,h_1,...$:

asymptoticOrder = 10;
serieseqs = series //
ReplaceAll[max -> asymptoticOrder] // Activate //
Map[
If[Head[#] === SeriesData, #[], #] &
] //
Flatten // Thread[# == 0] & // DeleteCases[True];


Luckily the equations we got here are all nice and easy to solve:

seriessol = Inactive[Solve][serieseqs,
Union[Cases[serieseqs, _h, ∞]][[;; UpTo[Length@serieseqs]]]
] // Activate;
seriessol // Apply[List, #, {2}] & // Map[Grid[#, Frame -> All] &] // Row Thus the corresponding approximate solution for $H$

H[η] // seriesRules //
ReplaceAll[max -> asymptoticOrder] // Activate //
ReplaceAll[seriessol] (But do be aware this is only a formal series solution. The convergence is yet to prove.)

• What is "UpTo" ? I get an error message: "1;;UpTo is not a valid Span specification" Mar 3, 2018 at 18:28
• @VaclavKotesovec It's a handy function introduced in version 10.3. Should be easy to replace with something involving Min. Mar 4, 2018 at 13:24
• @Silvia Should we end up your final equation with $O(\eta^{12})$? Mar 12, 2018 at 6:36
• @W.Robin The calculation rules on SeriesData are built-in. They decided the exponent of $O(\eta^{11})$. In this case, there maybe the coefficient for the $\eta^{11}$ term stored at position {3,12} of SeriesData[...], but it's the penultimate argument of it decides to which term to show in a StandardForm / TraditionalForm. Mar 12, 2018 at 17:37

Mathematica 11.3? Will have AsymptoticDSolve and support for these And WKB also It will also finally have series solution for DSolve  • thanks for the information. I only have access to v9 and v10. Is there any workaround for this AsymptoticDSolve on the lower version? Jan 27, 2018 at 16:31
It's not a answer of Yours question,only a info if Mathematica 11.3 can solve or not.
Using Mathematica 11.3,it seems can't find solution. 