> restart;

We express (Laplacian)^2 of psi in polar coordinates.

 > Lpsi := diff(psi(r,theta),r\$2)+1/r*diff(psi(r,theta),r) + 1/r^2*diff(psi(r,theta),theta\$2); L2psi := diff(Lpsi,r\$2)+1/r*diff(Lpsi,r) + 1/r^2*diff(Lpsi,theta\$2): collect(L2psi,r);

We compute (Laplacian)^2 of the nth term of the Fourier series (in theta) of psi.

 > nterm := psi[n](r)*exp(I*n*theta); Lnterm := diff(nterm,r\$2)+1/r*diff(nterm,r) + 1/r^2*diff(nterm,theta\$2); L2nterm := collect(diff(Lnterm,r\$2)+1/r*diff(Lnterm,r) + 1/r^2*diff(Lnterm,theta\$2),r): collect(simplify(L2nterm*exp(-I*n*theta)),r)*exp(I*n*theta);

We solve the ODE that results from setting (Laplacian)^2 of the nth term of the Fourier series equal to zero.

 > ode1 := diff(psi[n](r),r\$4) + 2/r*diff(psi[n](r),r\$3) - (2*n^2+1)/r^2*diff(psi[n](r),r\$2)          + (2*n^2+1)/r^3*diff(psi[n](r),r) + (n^4 - 4*n^2)/r^4*psi[n](r) = 0; dsolve(ode1);

But the solution has a different form when n = 0, 1 and -1.

 > for j from 0 to 1 do sprintf("for n = %d:",j); ode1 := diff(psi[j](r),r\$4) + 2/r*diff(psi[j](r),r\$3) - (2*j^2+1)/r^2*diff(psi[j](r),r\$2)          + (2*j^2+1)/r^3*diff(psi[j](r),r) + (j^4 - 4*j^2)/r^4*psi[j](r) = 0; dsolve(ode1); end do;

So these are the general solutions:

 > psi[n] := r -> _C1*r^(-n+2)+_C2*r^(n+2)+_C3*r^(-n)+_C4*r^n; psi[0] := r -> _C1*r^2+_C2*r^2*ln(r)+_C3+_C4*ln(r); psi[1] := r -> _C1*r^3+_C2/r+_C3*r+_C4*r*ln(r);

 >

Maple TM is a registered trademark of Waterloo Maple Inc.
Math rendered by WebEQ