!*********************************************************************** !* Calculate the Legendre polynomials for a complex argument * !* ------------------------------------------------------------------- * !* EXPLANATION: * !* * !* ========================================================== * !* Purpose: This program computes the Legendre polynomials * !* of the second kind Qn(z) and Qn'(z) for a complex * !* argument using subroutine CLQN. * !* Input : x --- Real part of z * !* y --- Imaginary part of z * !* n --- Degree of Qn(z), n = 0,1,... * !* Output: CQN(n) --- Qn(z) * !* CQD(n) --- Qn'(z) * !* Examples: * !* * !* z = 0.5 + 0.5 i * !* n Re[Qn(z)] Im[Qn(z)] Re[Qn'(z)] Im[Qn'(z)] * !* ----------------------------------------------------------- * !* 0 .402359D+00 .553574D+00 .800000D+00 .400000D+00 * !* 1 -.107561D+01 .477967D+00 .602359D+00 .115357D+01 * !* 2 -.136636D+01 -.725018D+00 -.242682D+01 .183390D+01 * !* 3 .182619D+00 -.206146D+01 -.622944D+01 -.247151D+01 * !* 4 .298834D+01 -.110022D+01 -.114849D+01 -.125963D+02 * !* 5 .353361D+01 .334847D+01 .206656D+02 -.123735D+02 * !* * !* z = 3.0 + 2.0 i * !* n Re[Qn(z)] Im[Qn(z)] Re[Qn'(z)] Im[Qn'(z)] * !* ----------------------------------------------------------- * !* 0 .229073D+00 -.160875D+00 -.250000D-01 .750000D-01 * !* 1 .896860D-02 -.244805D-01 .407268D-02 .141247D-01 * !* 2 -.736230D-03 -.281865D-02 .190581D-02 .155860D-02 * !* 3 -.264727D-03 -.227023D-03 .391535D-03 .314880D-04 * !* 4 -.430648D-04 -.443187D-05 .527190D-04 -.305592D-04 * !* 5 -.481362D-05 .265297D-05 .395108D-05 -.839883D-05 * !* ========================================================== * !* * !* ------------------------------------------------------------------- * !* SAMPLE RUN: * !* * !* Please enter Nmax, x and y (z=x+iy): 5 3 2 * !* x = 3.0, y = 2.0 * !* * !* n Re[Qn(z)] Im[Qn(z)] Re[Qn'(z)] Im[Qn'(z)] * !* ----------------------------------------------------------- * !* 0 0.229073D+00 -0.160875D+00 -0.250000D-01 0.750000D-01 * !* 1 0.896860D-02 -0.244805D-01 0.407268D-02 0.141247D-01 * !* 2 -0.736230D-03 -0.281865D-02 0.190581D-02 0.155860D-02 * !* 3 -0.264727D-03 -0.227023D-03 0.391535D-03 0.314880D-04 * !* 4 -0.430648D-04 -0.443187D-05 0.527190D-04 -0.305592D-04 * !* 5 -0.481362D-05 0.265297D-05 0.395108D-05 -0.839883D-05 * !* * !* ------------------------------------------------------------------- * !* REFERENCE: "Fortran Routines for Computation of Special Functions, * !* jin.ece.uiuc.edu/routines/routines.html". * !* * !* F90 Release By J-P Moreau, Paris. * !* (www.jpmoreau.fr) * !*********************************************************************** PROGRAM MCLQN IMPLICIT DOUBLE PRECISION (X,Y) IMPLICIT COMPLEX*16 (C,Z) DIMENSION CQN(0:100),CQD(0:100) WRITE(*,*) WRITE(*,10,advance='no') READ(*,*)N,X,Y WRITE(*,30)X,Y WRITE(*,*) CALL CLQN(N,X,Y,CQN,CQD) WRITE(*,*)' n Re[Qn(z)] Im[Qn(z)] Re[Qn''(z)] Im[Qn''(z)]' WRITE(*,*)' -----------------------------------------------------------' DO K=0,N WRITE(*,20)K,CQN(K),CQD(K) END DO WRITE(*,*) STOP 10 FORMAT(' Please enter Nmax, x and y (z=x+iy): ') 20 FORMAT(1X,I3,4D14.6) 30 FORMAT(3X,'x =',F5.1,', ','y =',F5.1) END SUBROUTINE CLQN(N,X,Y,CQN,CQD) ! =================================================== ! Purpose: Compute the Legendre functions Qn(z) and ! their derivatives Qn'(z) for a complex ! argument ! Input : x --- Real part of z ! y --- Imaginary part of z ! n --- Degree of Qn(z), n = 0,1,2,... ! Output: CQN(n) --- Qn(z) ! CQD(n) --- Qn'(z) ! =================================================== IMPLICIT DOUBLE PRECISION (X,Y) IMPLICIT COMPLEX*16 (C,Z) DIMENSION CQN(0:N),CQD(0:N) Z=CMPLX(X,Y) IF (Z.EQ.1.0D0) THEN DO K=0,N CQN(K)=(1.0D+300,0.0D0) CQD(K)=(1.0D+300,0.0D0) END DO RETURN ENDIF LS=1 IF (CDABS(Z).GT.1.0D0) LS=-1 CQ0=0.5D0*CDLOG(LS*(1.0D0+Z)/(1.0D0-Z)) CQ1=Z*CQ0-1.0D0 CQN(0)=CQ0 CQN(1)=CQ1 IF (CDABS(Z).LT.1.0001D0) THEN print *, CQ0 print *, CQ1 CQF0=CQ0 CQF1=CQ1 DO K=2,N CQF2=((2.0D0*K-1.0D0)*Z*CQF1-(K-1.0D0)*CQF0)/K CQN(K)=CQF2 CQF0=CQF1 CQF1=CQF2 END DO ELSE IF (CDABS(Z).GT.1.1D0) THEN KM=40+N ELSE KM=(40+N)*INT(-1.0-1.8*LOG(CDABS(Z-1.0))) ENDIF CQF2=0.0D0 CQF1=1.0D0 DO K=KM,0,-1 CQF0=((2*K+3.0D0)*Z*CQF1-(K+2.0D0)*CQF2)/(K+1.0D0) IF (K.LE.N) CQN(K)=CQF0 CQF2=CQF1 CQF1=CQF0 END DO DO K=0,N CQN(K)=CQN(K)*CQ0/CQF0 END DO ENDIF CQD(0)=(CQN(1)-Z*CQN(0))/(Z*Z-1.0D0) DO K=1,N CQD(K)=(K*Z*CQN(K)-K*CQN(K-1))/(Z*Z-1.0D0) END DO RETURN END !end of file mclqn.f90