(*This program solves Δu + λf(u) = 0 for r<1, u=0 for r=1, in n dimensions. Draws the bifurcation curve (λ, u(0))*)
Clear["`*"]
f[u_]=u+ 0.5*u *Sin[u];
n=3;
u00=0;
delu=0.07;
nsteps=400;
tend= 1000;
h=$MachineEpsilon;
dirichlet={};
For[k=1,k ≤ nsteps,k++,
 u0=u00+k*delu;
sol=NDSolve[{x''[t]==-(n-1)*x'[t]/t-f[x[t]],x[h]==u0,x'[h]==0,WhenEvent[{x[t]*x'[t]==0},If[Abs[x[t]]<10^(-8),AppendTo[dirichlet,{t^2,u0}]];"StopIntegration"]},x,{t,h,tend}];
]

ListPlot[dirichlet,Joined -> True,  AxesLabel -> {"λ","u(0)"}, PlotRange -> All, PlotStyle -> Thickness[0.007],PlotLabel->"Bifurcation Diagram for Dirichlet problem"]