Having a problem with soling a second order ODE equation using Matlab

[kubekas]

I am having a problem coding a Matlab code that solves a second ODE equation which I give below:

x^3*(1-2*x*M)d^2J(x)/dx^2+2*(2*x^2+i*nu*x-7*x^3*M)*dJ(x)/dx
-2*(2*x+8*M*x^2+i*nu)*J(x)=0.

where

M = 1 (Mass of a black hole),
nu = 0.74734+0.17792*i,
J is a function of x,
i represents a complex number.

This equation is very complex to solve. Fistly I tried the Matlab shooting
method and it did not work because this equation is singular at the boundries i.e
0 and 1. Now I am trying to solve it using Matlab Finite Difference Method. I must point out that this problem has no analytic solution and we hope that it can be solved numerically.

Can anyone out there help me with this problem.
Thanks
Amos

 

[lalbatros]

kubekas,

You did not specify the boundary conditions. Boundary conditions, or initial values ?
I assume J is a complex function.
I am surprised by the fixed complex value for nu in your formulation. Where does that come from (I am not a computer!).
Could you state the problem a little bit more clearly and eventually give the physical meaning.

Before trying rigorous handling of this problem, have you tried some intuitive approach, like by suppressing numerically the singularities?

I don't see why x=1 is singular, but have you noticed that x=1/2 is a singularity?

Have you tried to turn this 2nd order equation in a 1st order system?

Have you tried to write analytical solutions for constant coefficients and to patch these solutions together? By this, I mean developing a simple personnal solver instead of relying on a general purpose solver: this could provide you a better understanding of the question and possibly a specialised solver that could be simpler and more efficient.

Have you read something about singularities in differential equations, have you seen other examples of that in other field of physics. I could suggest you wave propagation in inhomogeneuous plasmas, for example, where such singularites occurs when the refractive index changes sign. Physically there is often an improved model where the singularity disappears, and the solution can be physically very different. (for example: wave attenuation effects in plasmas or magnetisation/finite larmor radius/ effects that remove the singularity)

Michel

 

[kubekas]

Hi, Michel

Thanks for the replay. We derived the above equation by solving linearized Einstein equations for the Bondi-Sachs metric when the background is Schwarzschild. To do that we solved hypersurface, evolution, and contraints equations and in the end we ended up with the above equation when {\it l}(angular momentum)=2. This equation is similar to both the Regge-Wheeler and Zerilli equations in that they all discribe the same linear perturbations of the Schwarzschild black. By solving the above equation, we hope to produce quasi-normal modes of a Schwarzcshild black hole hence to see that the above equation is solverable, we are using the already calculated \nu for {\it l}=1 from the Zerilli equation.

Yes that is true at x=1/2 we got a singulary, typing error. Like I outlined, I tried to surprise these singularity when I was trieying the shooting method and I was advised by my promoter that unfortunately this method will not for this problem. In that processe we god four 1st order systems of equation and that is because J is complex as you have pointed out, that means firslty we separated the above equation into two 2nd order equations, i.e. one real and the other imaginary and then were able to reduce the oder for these two equations. Right now we think that the Finit Difference method will work since we will not have to separate this equation and that we will not have deal we singularities.

The aboundry condition are: x(0)=0 and at x(1/2)=1.


Amos

 

[lalbatros]

Another suggestion.

It might not be too complicated to find some analytical solutions around the singularities. These solutions will depends only on a few integration constants. These special solution in these finite domains around the singularities can be integrated in the full Finite Difference system of equations and solve the problem as usual by solving the FD linear system.

I think it is rather logical to study the behaviour of the solutions in the vicinity of the singularities. An analytical approach can be used numerically and can also give further insight in the physics. A fully numerical approach means a numerical analysis of the singularity and implies a capable software.

You could also see what Mathematica can offer you. I can't tell you a lot about it because I own a very old version.

Michel

 

[lalbatros]

One more suggestion,

You could try a Fourier analysis.
You would have to pay for a bunch of FFT's and solving then a large system.

I also wondered if you were not dealing with an eignevalue problem.
Would the "nu" variable not be some evaluation of this eigenvalue?
If I am not mistaken, then you would need to solve the eigenvalue problem with some ad-hoc method (one eigenvalue or the full spectrum). And you could compare with your reference model. Solving the system with an approximate eigenvalue might not lead to the correct result.

Michel