Solving numerically a second order ODE with Maple

Question

I am trying to solve the following ODE using the program Maple:

This ODE can easily be solved analytically on the interval [0,1/2] and the solution is given by

Now, using Maple, I got

restart;
ODE := diff(y(x),x$2)*x*abs(ln(x))-2*diff(y(x), x)=0:
initialPosition := y(0) = 0;
initialVelocity := D(y)(1/2) = 1;
solution := dsolve({ODE, initialPosition, initialVelocity})

which gives me the same solution I found earlier. However, when I try to find the solution numerically: solution_num := dsolve({ODE, initialPosition, initialVelocity}, numeric), the program displays the following error

As far as I understand, the way Maple discretizes the ODE together with the initial condition at 0, y(0) = 0, creates an error as the program cannot handle the singularity from the log(0) appearing in the equation. They suggest to use another type of discretization to circumvent this singular endpoint: the Midpoint Method. As I am just starting with Maple and I'm not comfortable with ODEs discretization, I was wondering if there was a function implementing such an algorithm to solve my problem. Any help would be appreciated.

Answer 1

restart;
ODE:=diff(y(x),x$2)*x*abs(ln(x))
     -2*diff(y(x),x)=0:
initPos:=y(0)= 0:
initVel:=D(y)(1/2)=1:

sol:=dsolve({ODE,initPos,initVel})
       assuming x>0,x<1;

    y(x) = (-x/ln(x)-Ei(1,-ln(x)))*ln(2)^2

I put this style,etc on the plot of the exact solution only so that it can be nicely overlaid with the plot of the numeric solution.

P1:=plot(eval(y(x),sol),x=0 .. 0.5,
         adaptive=false,numpoints=11,
         color=blue,style=point,
         symbol=solidcircle,symbolsize=20):

Now, the numeric solution.

numsol:=dsolve({ODE,initPos,initVel},numeric,
               method=bvp[midrich],
               maxmesh=10^4,abserr=5e-4):

P2:=plot(X->eval(y(x),numsol(X)),0 .. 0.5):

Now let's see them together.

plots:-display(P1,P2);

Another (somewhat superior) way to plot that result returned by dsolve(...,numeric).

plots:-odeplot(numsol,[x,y(x)],x=0 .. 0.5);