Solving a second order coupled BVP with boundaries at infinity using DifferentialEquations.jl

New to Julia
,
1

I’m trying to solve the following coupled differential equations in Julia by using DifferentialEquations.jl and BoundaryValueDiffEq.jl

\begin{align} \frac{d^2f}{du^2}+\frac{1}{u}\frac{df}{du}-\frac{n^2}{u^2}(1-a)^2f-\frac{\lambda}{e^2}(f-1)(f-\tilde{V}_1)(f-\tilde{V}_m)&=0\\ \frac{d^2a}{du^2}-\frac{1}{u}\frac{da}{du}+2f^2(1-a)&=0 \end{align}

With boundary conditions a[0]=f[0]=0, a[\infty]=f[\infty]=1.

Indeed, I already got some solutions using the GeneralMIRK4() solver.The problem arise when I increase the interval of integration and the compilation takes for ever to end. This might be happening because I’m not implementing the best method to solve the equations. Any advice in how to speed up my code or how to implement another method to solve this system of differential equations is appreciated.

function vortexgg!(dx,x,p,t)
    dx[1] = x[2]
    dx[2] = -x[2]/t+x[1]*(1.0-x[3])^2*(1/t^2)+
    (14.0)*(x[1]-1)*(x[1]-(0.3))*(x[1]-(0.6))

    dx[3] = x[4]
    dx[4] = x[4]/t-(2.0*x[1]^2)*(1.0-x[3])
 end


 function bc1!(residual, u, p, t)
    residual[1] = u[1][1] - 0.0
    residual[2] = u[end][1] - 1.0
    residual[3] = u[1][3] - 0.0
    residual[4] = u[end][3] - 1.0

 end
 bvp1 = BVProblem(vortexgg!, bc1!, [0.00,0.0,0.0,0.00], tspan)
  st = solve(bvp1, GeneralMIRK4(),dt=0.1,alg_hints = [:stiff])
2

Try doing a change of variables like u = t/(1-t) that maps [0,\infty) to [0,1).

3 Likes
3

Thanks, I already tried doing that but somehow it makes the system more instable, and the functions don’t converge to the appropiate solution. I’ll give it another try, I’m guessing that this could speed up the code.

4

I am not sure what BoundaryValueDiffEq.jl uses, but I would consider simply using collocation on either [0,1] or [0, \infty] (which is a bit more tricky), see Boyd, J. P. (2001). Chebyshev and fourier spectral methods.

1 Like
5

I am not sure this is a good choice on large intervals.

Do you know where (f,a) converge when u goes to infinity?

In any cases, I would give a go to Shooting (see docs) with appropriate ode stepper.

Finally, the factor 1/u is not very nice from an ODE point of vue. Keep in mind that the solvers implement existence results from theorem. I see that your vector field is singular, a well chosen change of time variable could make it more well posed?

1 Like
6

So, this topic kinda piqued my interest and kept simmering in the back of my head. :grin: I finally sat down to work on an asymptotic expansion of the problem at the boundaries.
The behaviour of the functions, especially as u\rightarrow \infty in the present problem, can either guide you to an appropriate coordinate transformation (e.g. to make the functions decay nicely in the new coordinates), or you can use them to shift the right boundary to a location before \infty if you get a singularity there. I used both of these techniques to solve a BVP during my PhD work.

It’s been quite a while since I last did that, so I am pretty rusty. I started at u\rightarrow 0 first, and found

f_L(u\rightarrow 0) \sim u\cdot f_1 - u^2 \cdot f_2 + \ldots\\ a_L(u\rightarrow0) \sim u^2 \cdot a_1 - u^4 \cdot a_2 + \ldots

Is this useful to you? Values for the coefficients f_1 and a_1 are unknown, but that’s not the primary result of the analysis anyway. Interestingly, I have expressions for f_2 and a_2 (with some confidence).

As expected, the right boundary expansion turned out to be trickier, so I don’t have a result to show, yet – the effort passed my threshold for “nice evening fun”, the first part was ~8 pages of paper already.