1 |
|
|
MODULE SUPOL_MOD |
2 |
|
|
CONTAINS |
3 |
|
|
SUBROUTINE SUPOL(KNSMAX,DDMU,DDPOL,DDA,DDB,DDC,DDD,DDE,DDF,DDG,DDH,DDI) |
4 |
|
|
|
5 |
|
|
!**** *SUPOL * - Routine to compute the Legendre polynomials |
6 |
|
|
|
7 |
|
|
! Purpose. |
8 |
|
|
! -------- |
9 |
|
|
! For a given value of mu, computes the Legendre |
10 |
|
|
! polynomials. |
11 |
|
|
|
12 |
|
|
!** Interface. |
13 |
|
|
! ---------- |
14 |
|
|
! *CALL* *SUPOL(KNSMAX,DDMU,DDPOL,DDA,DDB,DDC,DDD,DDE |
15 |
|
|
! ,DDF,DDG,DDH,DDI) |
16 |
|
|
|
17 |
|
|
! Explicit arguments : |
18 |
|
|
! -------------------- |
19 |
|
|
! KNSMAX : Truncation (triangular) |
20 |
|
|
! DDMU : Abscissa at which the polynomials are computed (mu) |
21 |
|
|
! DDPOL : Polynomials (the first index is m and the second n) |
22 |
|
|
|
23 |
|
|
|
24 |
|
|
! Implicit arguments : None |
25 |
|
|
! -------------------- |
26 |
|
|
|
27 |
|
|
! Method. |
28 |
|
|
! ------- |
29 |
|
|
! See documentation |
30 |
|
|
|
31 |
|
|
! Externals. |
32 |
|
|
! ---------- |
33 |
|
|
|
34 |
|
|
! Reference. |
35 |
|
|
! ---------- |
36 |
|
|
! ECMWF Research Department documentation of the IFS |
37 |
|
|
|
38 |
|
|
! Author. |
39 |
|
|
! ------- |
40 |
|
|
! Mats Hamrud and Philippe Courtier *ECMWF* |
41 |
|
|
|
42 |
|
|
! Modifications. |
43 |
|
|
! -------------- |
44 |
|
|
! Original : 87-10-15 |
45 |
|
|
! K. YESSAD (MAY 1998): modification to avoid underflow. |
46 |
|
|
! ------------------------------------------------------------------ |
47 |
|
|
|
48 |
|
|
USE PARKIND1 ,ONLY : JPIM ,JPRB |
49 |
|
|
USE PARKIND2 ,ONLY : JPRH |
50 |
|
|
|
51 |
|
|
IMPLICIT NONE |
52 |
|
|
|
53 |
|
|
INTEGER(KIND=JPIM),INTENT(IN) :: KNSMAX |
54 |
|
|
REAL(KIND=JPRH) ,INTENT(IN) :: DDMU |
55 |
|
|
REAL(KIND=JPRH) ,INTENT(IN) :: DDC(0:KNSMAX,0:KNSMAX) |
56 |
|
|
REAL(KIND=JPRH) ,INTENT(IN) :: DDD(0:KNSMAX,0:KNSMAX) |
57 |
|
|
REAL(KIND=JPRH) ,INTENT(IN) :: DDE(0:KNSMAX,0:KNSMAX) |
58 |
|
|
REAL(KIND=JPRH) ,INTENT(IN) :: DDA(0:KNSMAX),DDB(0:KNSMAX),DDF(0:KNSMAX) |
59 |
|
|
REAL(KIND=JPRH) ,INTENT(IN) :: DDG(0:KNSMAX),DDH(0:KNSMAX),DDI(0:KNSMAX) |
60 |
|
|
REAL(KIND=JPRH) ,INTENT(OUT) :: DDPOL(0:KNSMAX,0:KNSMAX) |
61 |
|
|
|
62 |
|
|
REAL(KIND=JPRH) :: DLX,DLSITA,DL1SITA,DLKM2,DLKM1,DLK,DL1,DLS |
63 |
|
|
|
64 |
|
|
INTEGER(KIND=JPIM) :: JM, JN |
65 |
|
|
REAL(KIND=JPRB) :: Z |
66 |
|
|
|
67 |
|
|
! ------------------------------------------------------------------ |
68 |
|
|
|
69 |
|
|
!* 1. First two columns. |
70 |
|
|
! ------------------ |
71 |
|
|
|
72 |
|
|
DLX=DDMU |
73 |
|
|
DLSITA=SQRT(1.0_JPRB-DLX*DLX) |
74 |
|
|
|
75 |
|
|
! IF WE ARE LESS THAN 1Meter FROM THE POLE, |
76 |
|
|
IF(ABS(REAL(DLSITA,KIND(Z))) <= SQRT(EPSILON(Z)))THEN |
77 |
|
|
DLX=1._JPRB |
78 |
|
|
DLSITA=0._JPRB |
79 |
|
|
DL1SITA=0._JPRB |
80 |
|
|
ELSE |
81 |
|
|
DL1SITA=1.0_JPRB/DLSITA |
82 |
|
|
ENDIF |
83 |
|
|
DLKM2=1._JPRB |
84 |
|
|
DLKM1=DLX |
85 |
|
|
DDPOL(0,0)=DLKM2 |
86 |
|
|
DDPOL(0,1)=DLKM1*DDA(1) |
87 |
|
|
DDPOL(1,1)=DLSITA*DDB(1) |
88 |
|
|
DO JN=2,KNSMAX |
89 |
|
|
DLK=DDF(JN)*DLX*DLKM1-DDG(JN)*DLKM2 |
90 |
|
|
DL1=DDI(JN)*(DLKM1-DLX*DLK)*DL1SITA |
91 |
|
|
DDPOL(0,JN)=DLK*DDA(JN) |
92 |
|
|
DDPOL(1,JN)=DL1*DDB(JN) |
93 |
|
|
DLKM2=DLKM1 |
94 |
|
|
DLKM1=DLK |
95 |
|
|
ENDDO |
96 |
|
|
|
97 |
|
|
! ------------------------------------------------------------------ |
98 |
|
|
|
99 |
|
|
!* 2. Diagonal (the terms 0,0 and 1,1 have already been computed) |
100 |
|
|
! ----------------------------------------------------------- |
101 |
|
|
|
102 |
|
|
DLS=DL1SITA*TINY(DLS) |
103 |
|
|
|
104 |
|
|
!OCL SCALAR |
105 |
|
|
DO JN=2,KNSMAX |
106 |
|
|
DDPOL(JN,JN)=DDPOL(JN-1,JN-1)*DLSITA*DDH(JN) |
107 |
|
|
IF ( ABS(DDPOL(JN,JN)) < DLS ) DDPOL(JN,JN)=0.0_JPRB |
108 |
|
|
ENDDO |
109 |
|
|
|
110 |
|
|
! ------------------------------------------------------------------ |
111 |
|
|
|
112 |
|
|
!* 3. General recurrence. |
113 |
|
|
! ------------------- |
114 |
|
|
|
115 |
|
|
DO JN=3,KNSMAX |
116 |
|
|
!DIR$ IVDEP |
117 |
|
|
!OCL NOVREC |
118 |
|
|
DO JM=2,JN-1 |
119 |
|
|
DDPOL(JM,JN)=DDC(JM,JN)*DDPOL(JM-2,JN-2)& |
120 |
|
|
&-DDD(JM,JN)*DDPOL(JM-2,JN-1)*DLX & |
121 |
|
|
&+DDE(JM,JN)*DDPOL(JM ,JN-1)*DLX |
122 |
|
|
ENDDO |
123 |
|
|
ENDDO |
124 |
|
|
|
125 |
|
|
! ------------------------------------------------------------------ |
126 |
|
|
|
127 |
|
|
END SUBROUTINE SUPOL |
128 |
|
|
END MODULE SUPOL_MOD |
129 |
|
|
|
130 |
|
|
|