BVP4c not giving expected behavior, is the set up right?

I am trying to solve the following nonlinear diff eq:

With boundary conditions of .

Now for the issue, if we look at , we find: (I matched BC dropping negative exponent). This above diff eq was transformed from an orignal one with a change of variables: . If we sub this into our approximate expression, we find: .

My research advisor wants a slope of 1around the r-orgin, but my solver gives a value of (I attached some code which generates plots). Presummably this is just because , so maybe there is no issue? The issue might just be my advisor doesn't trust me, and he is too sure of what he thinks the answer should be...

Basically, does this code look ok? When treating the BC, I had to set to a semi large negative number and the ∞ boundary condition was moved to 1 (or just under it). I take all and apply an arctanh to stretch the interval. This was the only way I could find a stable solution. Can this change certain behaviors and give me a different slope? I'm guessing yes, but is the approximation "ok"? As in, does the slope of 1/2 more than likely correspond to ?