/[lmdze]/trunk/phylmd/Radlwsw/swde.f
ViewVC logotype

Diff of /trunk/phylmd/Radlwsw/swde.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/phylmd/Radlwsw/swde.f revision 71 by guez, Mon Jul 8 18:12:18 2013 UTC trunk/Sources/phylmd/Radlwsw/swde.f revision 207 by guez, Thu Sep 1 10:30:53 2016 UTC
# Line 1  Line 1 
1        SUBROUTINE SWDE (PGG,PREF,PRMUZ,PTO1,PW,  module swde_m
2       S                 PRE1,PRE2,PTR1,PTR2)  
3        use dimens_m    IMPLICIT NONE
4        use dimphy  
5        use raddim  contains
6        IMPLICIT none  
7  C    SUBROUTINE swde(pgg, pref, prmuz, pto1, pw, pre1, pre2, ptr1, ptr2)
8  C     ------------------------------------------------------------------      USE dimens_m
9  C     PURPOSE.      USE dimphy
10  C     --------      USE raddim
11  C           COMPUTES THE REFLECTIVITY AND TRANSMISSIVITY OF A CLOUDY      ! ------------------------------------------------------------------
12  C     LAYER USING THE DELTA-EDDINGTON'S APPROXIMATION.      ! PURPOSE.
13  C      ! --------
14  C     METHOD.      ! COMPUTES THE REFLECTIVITY AND TRANSMISSIVITY OF A CLOUDY
15  C     -------      ! LAYER USING THE DELTA-EDDINGTON'S APPROXIMATION.
16  C  
17  C          STANDARD DELTA-EDDINGTON LAYER CALCULATIONS.      ! METHOD.
18  C      ! -------
19  C     REFERENCE.  
20  C     ----------      ! STANDARD DELTA-EDDINGTON LAYER CALCULATIONS.
21  C  
22  C        SEE RADIATION'S PART OF THE MODEL'S DOCUMENTATION AND      ! REFERENCE.
23  C        ECMWF RESEARCH DEPARTMENT DOCUMENTATION OF THE IFS      ! ----------
24  C  
25  C     AUTHOR.      ! SEE RADIATION'S PART OF THE MODEL'S DOCUMENTATION AND
26  C     -------      ! ECMWF RESEARCH DEPARTMENT DOCUMENTATION OF THE IFS
27  C        JEAN-JACQUES MORCRETTE  *ECMWF*  
28  C      ! AUTHOR.
29  C     MODIFICATIONS.      ! -------
30  C     --------------      ! JEAN-JACQUES MORCRETTE  *ECMWF*
31  C        ORIGINAL : 88-12-15  
32  C     ------------------------------------------------------------------      ! MODIFICATIONS.
33  C* ARGUMENTS:      ! --------------
34  C      ! ORIGINAL : 88-12-15
35        DOUBLE PRECISION PGG(KDLON)   ! ASSYMETRY FACTOR      ! ------------------------------------------------------------------
36        DOUBLE PRECISION PREF(KDLON)  ! REFLECTIVITY OF THE UNDERLYING LAYER      ! * ARGUMENTS:
37        DOUBLE PRECISION PRMUZ(KDLON) ! COSINE OF SOLAR ZENITH ANGLE  
38        DOUBLE PRECISION PTO1(KDLON)  ! OPTICAL THICKNESS      DOUBLE PRECISION pgg(kdlon) ! ASSYMETRY FACTOR
39        DOUBLE PRECISION PW(KDLON)    ! SINGLE SCATTERING ALBEDO      DOUBLE PRECISION pref(kdlon) ! REFLECTIVITY OF THE UNDERLYING LAYER
40        DOUBLE PRECISION PRE1(KDLON)  ! LAYER REFLECTIVITY (NO UNDERLYING-LAYER REFLECTION)      DOUBLE PRECISION prmuz(kdlon) ! COSINE OF SOLAR ZENITH ANGLE
41        DOUBLE PRECISION PRE2(KDLON)  ! LAYER REFLECTIVITY      DOUBLE PRECISION pto1(kdlon) ! OPTICAL THICKNESS
42        DOUBLE PRECISION PTR1(KDLON)  ! LAYER TRANSMISSIVITY (NO UNDERLYING-LAYER REFLECTION)      DOUBLE PRECISION pw(kdlon) ! SINGLE SCATTERING ALBEDO
43        DOUBLE PRECISION PTR2(KDLON)  ! LAYER TRANSMISSIVITY      DOUBLE PRECISION pre1(kdlon) ! LAYER REFLECTIVITY (NO UNDERLYING-LAYER REFLECTION)
44  C      DOUBLE PRECISION pre2(kdlon) ! LAYER REFLECTIVITY
45  C* LOCAL VARIABLES:      DOUBLE PRECISION ptr1(kdlon) ! LAYER TRANSMISSIVITY (NO UNDERLYING-LAYER REFLECTION)
46  C      DOUBLE PRECISION ptr2(kdlon) ! LAYER TRANSMISSIVITY
47        INTEGER jl  
48        DOUBLE PRECISION ZFF, ZGP, ZTOP, ZWCP, ZDT, ZX1, ZWM      ! * LOCAL VARIABLES:
49        DOUBLE PRECISION ZRM2, ZRK, ZX2, ZRP, ZALPHA, ZBETA, ZARG  
50        DOUBLE PRECISION ZEXMU0, ZARG2, ZEXKP, ZEXKM, ZXP2P, ZXM2P, ZAP2B      INTEGER jl
51        DOUBLE PRECISION ZAM2B      DOUBLE PRECISION zff, zgp, ztop, zwcp, zdt, zx1, zwm
52        DOUBLE PRECISION ZA11, ZA12, ZA13, ZA21, ZA22, ZA23      DOUBLE PRECISION zrm2, zrk, zx2, zrp, zalpha, zbeta, zarg
53        DOUBLE PRECISION ZDENA, ZC1A, ZC2A, ZRI0A, ZRI1A      DOUBLE PRECISION zexmu0, zarg2, zexkp, zexkm, zxp2p, zxm2p, zap2b
54        DOUBLE PRECISION ZRI0B, ZRI1B      DOUBLE PRECISION zam2b
55        DOUBLE PRECISION ZB21, ZB22, ZB23, ZDENB, ZC1B, ZC2B      DOUBLE PRECISION za11, za12, za13, za21, za22, za23
56        DOUBLE PRECISION ZRI0C, ZRI1C, ZRI0D, ZRI1D      DOUBLE PRECISION zdena, zc1a, zc2a, zri0a, zri1a
57  C     ------------------------------------------------------------------      DOUBLE PRECISION zri0b, zri1b
58  C      DOUBLE PRECISION zb21, zb22, zb23, zdenb, zc1b, zc2b
59  C*         1.      DELTA-EDDINGTON CALCULATIONS      DOUBLE PRECISION zri0c, zri1c, zri0d, zri1d
60  C      ! ------------------------------------------------------------------
61   100  CONTINUE  
62  C      ! *         1.      DELTA-EDDINGTON CALCULATIONS
63        DO 131 JL   =   1, KDLON  
64  C  
65  C*         1.1     SET UP THE DELTA-MODIFIED PARAMETERS      DO jl = 1, kdlon
66  C  
67   110  CONTINUE         ! *         1.1     SET UP THE DELTA-MODIFIED PARAMETERS
68  C  
69        ZFF = PGG(JL)*PGG(JL)  
70        ZGP = PGG(JL)/(1.+PGG(JL))         zff = pgg(jl)*pgg(jl)
71        ZTOP = (1.- PW(JL) * ZFF) * PTO1(JL)         zgp = pgg(jl)/(1.+pgg(jl))
72        ZWCP = (1-ZFF)* PW(JL) /(1.- PW(JL) * ZFF)         ztop = (1.-pw(jl)*zff)*pto1(jl)
73        ZDT = 2./3.         zwcp = (1-zff)*pw(jl)/(1.-pw(jl)*zff)
74        ZX1 = 1.-ZWCP*ZGP         zdt = 2./3.
75        ZWM = 1.-ZWCP         zx1 = 1. - zwcp*zgp
76        ZRM2 =  PRMUZ(JL) * PRMUZ(JL)         zwm = 1. - zwcp
77        ZRK = SQRT(3.*ZWM*ZX1)         zrm2 = prmuz(jl)*prmuz(jl)
78        ZX2 = 4.*(1.-ZRK*ZRK*ZRM2)         zrk = sqrt(3.*zwm*zx1)
79        ZRP=ZRK/ZX1         zx2 = 4.*(1.-zrk*zrk*zrm2)
80        ZALPHA = 3.*ZWCP*ZRM2*(1.+ZGP*ZWM)/ZX2         zrp = zrk/zx1
81        ZBETA = 3.*ZWCP* PRMUZ(JL) *(1.+3.*ZGP*ZRM2*ZWM)/ZX2         zalpha = 3.*zwcp*zrm2*(1.+zgp*zwm)/zx2
82  CMAF      ZARG=MIN(ZTOP/PRMUZ(JL),200.)         zbeta = 3.*zwcp*prmuz(jl)*(1.+3.*zgp*zrm2*zwm)/zx2
83        ZARG=MIN(ZTOP/PRMUZ(JL),2.0d+2)         ! MAF      ZARG=MIN(ZTOP/PRMUZ(JL),200.)
84        ZEXMU0=EXP(-ZARG)         zarg = min(ztop/prmuz(jl), 2.0D+2)
85  CMAF      ZARG2=MIN(ZRK*ZTOP,200.)         zexmu0 = exp(-zarg)
86        ZARG2=MIN(ZRK*ZTOP,2.0d+2)         ! MAF      ZARG2=MIN(ZRK*ZTOP,200.)
87        ZEXKP=EXP(ZARG2)         zarg2 = min(zrk*ztop, 2.0D+2)
88        ZEXKM = 1./ZEXKP         zexkp = exp(zarg2)
89        ZXP2P = 1.+ZDT*ZRP         zexkm = 1./zexkp
90        ZXM2P = 1.-ZDT*ZRP         zxp2p = 1. + zdt*zrp
91        ZAP2B = ZALPHA+ZDT*ZBETA         zxm2p = 1. - zdt*zrp
92        ZAM2B = ZALPHA-ZDT*ZBETA         zap2b = zalpha + zdt*zbeta
93  C         zam2b = zalpha - zdt*zbeta
94  C*         1.2     WITHOUT REFLECTION FROM THE UNDERLYING LAYER  
95  C         ! *         1.2     WITHOUT REFLECTION FROM THE UNDERLYING LAYER
96   120  CONTINUE  
97  C  
98        ZA11 = ZXP2P         za11 = zxp2p
99        ZA12 = ZXM2P         za12 = zxm2p
100        ZA13 = ZAP2B         za13 = zap2b
101        ZA22 = ZXP2P*ZEXKP         za22 = zxp2p*zexkp
102        ZA21 = ZXM2P*ZEXKM         za21 = zxm2p*zexkm
103        ZA23 = ZAM2B*ZEXMU0         za23 = zam2b*zexmu0
104        ZDENA = ZA11 * ZA22 - ZA21 * ZA12         zdena = za11*za22 - za21*za12
105        ZC1A = (ZA22*ZA13-ZA12*ZA23)/ZDENA         zc1a = (za22*za13-za12*za23)/zdena
106        ZC2A = (ZA11*ZA23-ZA21*ZA13)/ZDENA         zc2a = (za11*za23-za21*za13)/zdena
107        ZRI0A = ZC1A+ZC2A-ZALPHA         zri0a = zc1a + zc2a - zalpha
108        ZRI1A = ZRP*(ZC1A-ZC2A)-ZBETA         zri1a = zrp*(zc1a-zc2a) - zbeta
109        PRE1(JL) = (ZRI0A-ZDT*ZRI1A)/ PRMUZ(JL)         pre1(jl) = (zri0a-zdt*zri1a)/prmuz(jl)
110        ZRI0B = ZC1A*ZEXKM+ZC2A*ZEXKP-ZALPHA*ZEXMU0         zri0b = zc1a*zexkm + zc2a*zexkp - zalpha*zexmu0
111        ZRI1B = ZRP*(ZC1A*ZEXKM-ZC2A*ZEXKP)-ZBETA*ZEXMU0         zri1b = zrp*(zc1a*zexkm-zc2a*zexkp) - zbeta*zexmu0
112        PTR1(JL) = ZEXMU0+(ZRI0B+ZDT*ZRI1B)/ PRMUZ(JL)         ptr1(jl) = zexmu0 + (zri0b+zdt*zri1b)/prmuz(jl)
113  C  
114  C*         1.3     WITH REFLECTION FROM THE UNDERLYING LAYER         ! *         1.3     WITH REFLECTION FROM THE UNDERLYING LAYER
115  C  
116   130  CONTINUE  
117  C         zb21 = za21 - pref(jl)*zxp2p*zexkm
118        ZB21 = ZA21- PREF(JL) *ZXP2P*ZEXKM         zb22 = za22 - pref(jl)*zxm2p*zexkp
119        ZB22 = ZA22- PREF(JL) *ZXM2P*ZEXKP         zb23 = za23 - pref(jl)*zexmu0*(zap2b-prmuz(jl))
120        ZB23 = ZA23- PREF(JL) *ZEXMU0*(ZAP2B - PRMUZ(JL) )         zdenb = za11*zb22 - zb21*za12
121        ZDENB = ZA11 * ZB22 - ZB21 * ZA12         zc1b = (zb22*za13-za12*zb23)/zdenb
122        ZC1B = (ZB22*ZA13-ZA12*ZB23)/ZDENB         zc2b = (za11*zb23-zb21*za13)/zdenb
123        ZC2B = (ZA11*ZB23-ZB21*ZA13)/ZDENB         zri0c = zc1b + zc2b - zalpha
124        ZRI0C = ZC1B+ZC2B-ZALPHA         zri1c = zrp*(zc1b-zc2b) - zbeta
125        ZRI1C = ZRP*(ZC1B-ZC2B)-ZBETA         pre2(jl) = (zri0c-zdt*zri1c)/prmuz(jl)
126        PRE2(JL) = (ZRI0C-ZDT*ZRI1C) / PRMUZ(JL)         zri0d = zc1b*zexkm + zc2b*zexkp - zalpha*zexmu0
127        ZRI0D = ZC1B*ZEXKM + ZC2B*ZEXKP - ZALPHA*ZEXMU0         zri1d = zrp*(zc1b*zexkm-zc2b*zexkp) - zbeta*zexmu0
128        ZRI1D = ZRP * (ZC1B*ZEXKM - ZC2B*ZEXKP) - ZBETA*ZEXMU0         ptr2(jl) = zexmu0 + (zri0d+zdt*zri1d)/prmuz(jl)
129        PTR2(JL) = ZEXMU0 + (ZRI0D + ZDT*ZRI1D) / PRMUZ(JL)  
130  C      END DO
131   131  CONTINUE      RETURN
132        RETURN    END SUBROUTINE swde
133        END  
134    end module swde_m
Legend:
Removed from v.71  
changed lines
  Added in v.207
  ViewVC Help
Powered by ViewVC 1.1.21