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

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

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

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