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