/[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 24 by guez, Wed Mar 3 13:23:49 2010 UTC trunk/phylmd/Radlwsw/swde.f revision 254 by guez, Mon Feb 5 10:39:38 2018 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        REAL*8 PGG(KDLON)   ! ASSYMETRY FACTOR      ! ------------------------------------------------------------------
36        REAL*8 PREF(KDLON)  ! REFLECTIVITY OF THE UNDERLYING LAYER      ! * ARGUMENTS:
37        REAL*8 PRMUZ(KDLON) ! COSINE OF SOLAR ZENITH ANGLE  
38        REAL*8 PTO1(KDLON)  ! OPTICAL THICKNESS      DOUBLE PRECISION pgg(kdlon) ! ASSYMETRY FACTOR
39        REAL*8 PW(KDLON)    ! SINGLE SCATTERING ALBEDO      DOUBLE PRECISION pref(kdlon) ! REFLECTIVITY OF THE UNDERLYING LAYER
40        REAL*8 PRE1(KDLON)  ! LAYER REFLECTIVITY (NO UNDERLYING-LAYER REFLECTION)      DOUBLE PRECISION prmuz(kdlon) ! COSINE OF SOLAR ZENITH ANGLE
41        REAL*8 PRE2(KDLON)  ! LAYER REFLECTIVITY      DOUBLE PRECISION pto1(kdlon) ! OPTICAL THICKNESS
42        REAL*8 PTR1(KDLON)  ! LAYER TRANSMISSIVITY (NO UNDERLYING-LAYER REFLECTION)      DOUBLE PRECISION pw(kdlon) ! SINGLE SCATTERING ALBEDO
43        REAL*8 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        REAL*8 ZFF, ZGP, ZTOP, ZWCP, ZDT, ZX1, ZWM      ! * LOCAL VARIABLES:
49        REAL*8 ZRM2, ZRK, ZX2, ZRP, ZALPHA, ZBETA, ZARG  
50        REAL*8 ZEXMU0, ZARG2, ZEXKP, ZEXKM, ZXP2P, ZXM2P, ZAP2B, ZAM2B      INTEGER jl
51        REAL*8 ZA11, ZA12, ZA13, ZA21, ZA22, ZA23      DOUBLE PRECISION zff, zgp, ztop, zwcp, zdt, zx1, zwm
52        REAL*8 ZDENA, ZC1A, ZC2A, ZRI0A, ZRI1A      DOUBLE PRECISION zrm2, zrk, zx2, zrp, zalpha, zbeta, zarg
53        REAL*8 ZRI0B, ZRI1B      DOUBLE PRECISION zexmu0, zarg2, zexkp, zexkm, zxp2p, zxm2p, zap2b
54        REAL*8 ZB21, ZB22, ZB23, ZDENB, ZC1B, ZC2B      DOUBLE PRECISION zam2b
55        REAL*8 ZRI0C, ZRI1C, ZRI0D, ZRI1D      DOUBLE PRECISION za11, za12, za13, za21, za22, za23
56  C     ------------------------------------------------------------------      DOUBLE PRECISION zdena, zc1a, zc2a, zri0a, zri1a
57  C      DOUBLE PRECISION zri0b, zri1b
58  C*         1.      DELTA-EDDINGTON CALCULATIONS      DOUBLE PRECISION zb21, zb22, zb23, zdenb, zc1b, zc2b
59  C      DOUBLE PRECISION zri0c, zri1c, zri0d, zri1d
60   100  CONTINUE      ! ------------------------------------------------------------------
61  C  
62        DO 131 JL   =   1, KDLON      ! *         1.      DELTA-EDDINGTON CALCULATIONS
63  C  
64  C*         1.1     SET UP THE DELTA-MODIFIED PARAMETERS  
65  C      DO jl = 1, kdlon
66   110  CONTINUE  
67  C         ! *         1.1     SET UP THE DELTA-MODIFIED PARAMETERS
68        ZFF = PGG(JL)*PGG(JL)  
69        ZGP = PGG(JL)/(1.+PGG(JL))  
70        ZTOP = (1.- PW(JL) * ZFF) * PTO1(JL)         zff = pgg(jl)*pgg(jl)
71        ZWCP = (1-ZFF)* PW(JL) /(1.- PW(JL) * ZFF)         zgp = pgg(jl)/(1.+pgg(jl))
72        ZDT = 2./3.         ztop = (1.-pw(jl)*zff)*pto1(jl)
73        ZX1 = 1.-ZWCP*ZGP         zwcp = (1-zff)*pw(jl)/(1.-pw(jl)*zff)
74        ZWM = 1.-ZWCP         zdt = 2./3.
75        ZRM2 =  PRMUZ(JL) * PRMUZ(JL)         zx1 = 1. - zwcp*zgp
76        ZRK = SQRT(3.*ZWM*ZX1)         zwm = 1. - zwcp
77        ZX2 = 4.*(1.-ZRK*ZRK*ZRM2)         zrm2 = prmuz(jl)*prmuz(jl)
78        ZRP=ZRK/ZX1         zrk = sqrt(3.*zwm*zx1)
79        ZALPHA = 3.*ZWCP*ZRM2*(1.+ZGP*ZWM)/ZX2         zx2 = 4.*(1.-zrk*zrk*zrm2)
80        ZBETA = 3.*ZWCP* PRMUZ(JL) *(1.+3.*ZGP*ZRM2*ZWM)/ZX2         zrp = zrk/zx1
81  CMAF      ZARG=MIN(ZTOP/PRMUZ(JL),200.)         zalpha = 3.*zwcp*zrm2*(1.+zgp*zwm)/zx2
82        ZARG=MIN(ZTOP/PRMUZ(JL),2.0d+2)         zbeta = 3.*zwcp*prmuz(jl)*(1.+3.*zgp*zrm2*zwm)/zx2
83        ZEXMU0=EXP(-ZARG)         ! MAF      ZARG=MIN(ZTOP/PRMUZ(JL),200.)
84  CMAF      ZARG2=MIN(ZRK*ZTOP,200.)         zarg = min(ztop/prmuz(jl), 2.0D+2)
85        ZARG2=MIN(ZRK*ZTOP,2.0d+2)         zexmu0 = exp(-zarg)
86        ZEXKP=EXP(ZARG2)         ! MAF      ZARG2=MIN(ZRK*ZTOP,200.)
87        ZEXKM = 1./ZEXKP         zarg2 = min(zrk*ztop, 2.0D+2)
88        ZXP2P = 1.+ZDT*ZRP         zexkp = exp(zarg2)
89        ZXM2P = 1.-ZDT*ZRP         zexkm = 1./zexkp
90        ZAP2B = ZALPHA+ZDT*ZBETA         zxp2p = 1. + zdt*zrp
91        ZAM2B = ZALPHA-ZDT*ZBETA         zxm2p = 1. - zdt*zrp
92  C         zap2b = zalpha + zdt*zbeta
93  C*         1.2     WITHOUT REFLECTION FROM THE UNDERLYING LAYER         zam2b = zalpha - zdt*zbeta
94  C  
95   120  CONTINUE         ! *         1.2     WITHOUT REFLECTION FROM THE UNDERLYING LAYER
96  C  
97        ZA11 = ZXP2P  
98        ZA12 = ZXM2P         za11 = zxp2p
99        ZA13 = ZAP2B         za12 = zxm2p
100        ZA22 = ZXP2P*ZEXKP         za13 = zap2b
101        ZA21 = ZXM2P*ZEXKM         za22 = zxp2p*zexkp
102        ZA23 = ZAM2B*ZEXMU0         za21 = zxm2p*zexkm
103        ZDENA = ZA11 * ZA22 - ZA21 * ZA12         za23 = zam2b*zexmu0
104        ZC1A = (ZA22*ZA13-ZA12*ZA23)/ZDENA         zdena = za11*za22 - za21*za12
105        ZC2A = (ZA11*ZA23-ZA21*ZA13)/ZDENA         zc1a = (za22*za13-za12*za23)/zdena
106        ZRI0A = ZC1A+ZC2A-ZALPHA         zc2a = (za11*za23-za21*za13)/zdena
107        ZRI1A = ZRP*(ZC1A-ZC2A)-ZBETA         zri0a = zc1a + zc2a - zalpha
108        PRE1(JL) = (ZRI0A-ZDT*ZRI1A)/ PRMUZ(JL)         zri1a = zrp*(zc1a-zc2a) - zbeta
109        ZRI0B = ZC1A*ZEXKM+ZC2A*ZEXKP-ZALPHA*ZEXMU0         pre1(jl) = (zri0a-zdt*zri1a)/prmuz(jl)
110        ZRI1B = ZRP*(ZC1A*ZEXKM-ZC2A*ZEXKP)-ZBETA*ZEXMU0         zri0b = zc1a*zexkm + zc2a*zexkp - zalpha*zexmu0
111        PTR1(JL) = ZEXMU0+(ZRI0B+ZDT*ZRI1B)/ PRMUZ(JL)         zri1b = zrp*(zc1a*zexkm-zc2a*zexkp) - zbeta*zexmu0
112  C         ptr1(jl) = zexmu0 + (zri0b+zdt*zri1b)/prmuz(jl)
113  C*         1.3     WITH REFLECTION FROM THE UNDERLYING LAYER  
114  C         ! *         1.3     WITH REFLECTION FROM THE UNDERLYING LAYER
115   130  CONTINUE  
116  C  
117        ZB21 = ZA21- PREF(JL) *ZXP2P*ZEXKM         zb21 = za21 - pref(jl)*zxp2p*zexkm
118        ZB22 = ZA22- PREF(JL) *ZXM2P*ZEXKP         zb22 = za22 - pref(jl)*zxm2p*zexkp
119        ZB23 = ZA23- PREF(JL) *ZEXMU0*(ZAP2B - PRMUZ(JL) )         zb23 = za23 - pref(jl)*zexmu0*(zap2b-prmuz(jl))
120        ZDENB = ZA11 * ZB22 - ZB21 * ZA12         zdenb = za11*zb22 - zb21*za12
121        ZC1B = (ZB22*ZA13-ZA12*ZB23)/ZDENB         zc1b = (zb22*za13-za12*zb23)/zdenb
122        ZC2B = (ZA11*ZB23-ZB21*ZA13)/ZDENB         zc2b = (za11*zb23-zb21*za13)/zdenb
123        ZRI0C = ZC1B+ZC2B-ZALPHA         zri0c = zc1b + zc2b - zalpha
124        ZRI1C = ZRP*(ZC1B-ZC2B)-ZBETA         zri1c = zrp*(zc1b-zc2b) - zbeta
125        PRE2(JL) = (ZRI0C-ZDT*ZRI1C) / PRMUZ(JL)         pre2(jl) = (zri0c-zdt*zri1c)/prmuz(jl)
126        ZRI0D = ZC1B*ZEXKM + ZC2B*ZEXKP - ZALPHA*ZEXMU0         zri0d = zc1b*zexkm + zc2b*zexkp - zalpha*zexmu0
127        ZRI1D = ZRP * (ZC1B*ZEXKM - ZC2B*ZEXKP) - ZBETA*ZEXMU0         zri1d = zrp*(zc1b*zexkm-zc2b*zexkp) - zbeta*zexmu0
128        PTR2(JL) = ZEXMU0 + (ZRI0D + ZDT*ZRI1D) / PRMUZ(JL)         ptr2(jl) = zexmu0 + (zri0d+zdt*zri1d)/prmuz(jl)
129  C  
130   131  CONTINUE      END DO
131        RETURN      RETURN
132        END    END SUBROUTINE swde
133    
134    end module swde_m
Legend:
Removed from v.24  
changed lines
  Added in v.254
  ViewVC Help
Powered by ViewVC 1.1.21