source: LMDZ6/branches/IPSL-CM6A-MR/libf/phymar/swde.F90 @ 5447

Last change on this file since 5447 was 2089, checked in by Laurent Fairhead, 10 years ago

Inclusion de la physique de MAR


Integration of MAR physics

File size: 5.2 KB
Line 
1!OPTIONS XOPT(HSFUN)
2SUBROUTINE SWDE &
3  &( KIDIA, KFDIA, KLON &
4  &, PGG  , PREF , PRMUZ, PTO1, PW &
5  &, PRE1 , PRE2 , PTR1 , PTR2 &
6  &)
7
8!**** *SWDE* - DELTA-EDDINGTON IN A CLOUDY LAYER
9
10!     PURPOSE.
11!     --------
12!           COMPUTES THE REFLECTIVITY AND TRANSMISSIVITY OF A CLOUDY
13!     LAYER USING THE DELTA-EDDINGTON'S APPROXIMATION.
14
15!**   INTERFACE.
16!     ----------
17!          *SWDE* IS CALLED BY *SWR*, *SWNI*
18
19
20!        EXPLICIT ARGUMENTS :
21!        --------------------
22! PGG    : (KLON)             ; ASSYMETRY FACTOR
23! PREF   : (KLON)             ; REFLECTIVITY OF THE UNDERLYING LAYER
24! PRMUZ  : (KLON)             ; COSINE OF SOLAR ZENITH ANGLE
25! PTO1   : (KLON)             ; OPTICAL THICKNESS
26! PW     : (KLON)             ; SINGLE SCATTERING ALBEDO
27!     ==== OUTPUTS ===
28! PRE1   : (KLON)             ; LAYER REFLECTIVITY ASSUMING NO
29!                             ; REFLECTION FROM UNDERLYING LAYER
30! PTR1   : (KLON)             ; LAYER TRANSMISSIVITY ASSUMING NO
31!                             ; REFLECTION FROM UNDERLYING LAYER
32! PRE2   : (KLON)             ; LAYER REFLECTIVITY ASSUMING
33!                             ; REFLECTION FROM UNDERLYING LAYER
34! PTR2   : (KLON)             ; LAYER TRANSMISSIVITY ASSUMING
35!                             ; REFLECTION FROM UNDERLYING LAYER
36
37!        IMPLICIT ARGUMENTS :   NONE
38!        --------------------
39
40!     METHOD.
41!     -------
42
43!          STANDARD DELTA-EDDINGTON LAYER CALCULATIONS.
44
45!     EXTERNALS.
46!     ----------
47
48!          NONE
49
50!     REFERENCE.
51!     ----------
52
53!        SEE RADIATION'S PART OF THE MODEL'S DOCUMENTATION AND
54!        ECMWF RESEARCH DEPARTMENT DOCUMENTATION OF THE IFS
55
56!     AUTHOR.
57!     -------
58!        JEAN-JACQUES MORCRETTE  *ECMWF*
59
60!     MODIFICATIONS.
61!     --------------
62!        ORIGINAL : 88-12-15
63!                   96-05-30 Michel Deque (security in EXP())
64!                   08-03-28 Hubert Gallee(lower/upper limit on PTR1,2)
65   
66!     ------------------------------------------------------------------
67
68
69!     ------------------------------------------------------------------
70
71!*       0.1   ARGUMENTS
72!              ---------
73
74
75#include "tsmbkind.h"
76
77IMPLICIT NONE
78
79
80!     DUMMY INTEGER SCALARS
81INTEGER_M :: KFDIA
82INTEGER_M :: KIDIA
83INTEGER_M :: KLON
84
85REAL_B :: PGG(KLON),PREF(KLON),PRMUZ(KLON),PTO1(KLON),PW(KLON)
86REAL_B :: PRE1(KLON),PRE2(KLON),PTR1(KLON),PTR2(KLON)
87
88!     LOCAL INTEGER SCALARS
89INTEGER_M :: JL
90
91!     LOCAL REAL SCALARS
92REAL_B :: ZA11, ZA12, ZA13, ZA21, ZA22, ZA23, ZALPHA,&
93          &ZAM2B, ZAP2B, ZARG, ZARG2, ZB21, ZB22, ZB23, &
94          &ZBETA, ZC1A, ZC1B, ZC2A, ZC2B, ZDENA, ZDENB, &
95          &ZDT, ZEXKM, ZEXKP, ZEXMU0, ZFF, ZGP, ZRI0A, &
96          &ZRI0B, ZRI0C, ZRI0D, ZRI1A, ZRI1B, ZRI1C, &
97          &ZRI1D, ZRK, ZRM2, ZRP, ZTOP, ZWCP, ZWM, ZX1, &
98          &ZX2, ZXM2P, ZXP2P
99
100
101
102!     ------------------------------------------------------------------
103
104!*         1.      DELTA-EDDINGTON CALCULATIONS
105
106
107DO JL   =   KIDIA,KFDIA
108
109!*         1.1     SET UP THE DELTA-MODIFIED PARAMETERS
110
111
112  ZFF = PGG(JL)*PGG(JL)
113  ZGP = PGG(JL)/(_ONE_+PGG(JL))
114  ZTOP = (_ONE_- PW(JL) * ZFF) * PTO1(JL)
115  ZWCP = (1-ZFF)* PW(JL) /(_ONE_- PW(JL) * ZFF)
116  ZDT = _TWO_/3._JPRB
117  ZX1 = _ONE_-ZWCP*ZGP
118  ZWM = _ONE_-ZWCP
119  ZRM2 =  PRMUZ(JL) * PRMUZ(JL)
120  ZRK = SQRT(3._JPRB*ZWM*ZX1)
121  ZX2 = 4._JPRB*(_ONE_-ZRK*ZRK*ZRM2)
122  ZRP=ZRK/ZX1
123  ZALPHA = 3._JPRB*ZWCP*ZRM2*(_ONE_+ZGP*ZWM)/ZX2
124  ZBETA = 3._JPRB*ZWCP* PRMUZ(JL) *(_ONE_+3._JPRB*ZGP*ZRM2*ZWM)/ZX2
125!     ZARG=MIN(ZTOP/PRMUZ(JL),200.)
126  ZARG=MAX(-085._JPRB,MIN(ZTOP/PRMUZ(JL),085._JPRB))
127  ZEXMU0=EXP(-ZARG)
128  ZARG2=MIN(ZRK*ZTOP,085._JPRB)
129  ZEXKP=EXP(ZARG2)
130  ZEXKM = _ONE_/ZEXKP
131  ZXP2P = _ONE_+ZDT*ZRP
132  ZXM2P = _ONE_-ZDT*ZRP
133  ZAP2B = ZALPHA+ZDT*ZBETA
134  ZAM2B = ZALPHA-ZDT*ZBETA
135
136!*         1.2     WITHOUT REFLECTION FROM THE UNDERLYING LAYER
137
138
139  ZA11 = ZXP2P
140  ZA12 = ZXM2P
141  ZA13 = ZAP2B
142  ZA22 = ZXP2P*ZEXKP
143  ZA21 = ZXM2P*ZEXKM
144  ZA23 = ZAM2B*ZEXMU0
145  ZDENA = ZA11 * ZA22 - ZA21 * ZA12
146  ZC1A = (ZA22*ZA13-ZA12*ZA23)/ZDENA
147  ZC2A = (ZA11*ZA23-ZA21*ZA13)/ZDENA
148  ZRI0A = ZC1A+ZC2A-ZALPHA
149  ZRI1A = ZRP*(ZC1A-ZC2A)-ZBETA
150  PRE1(JL) = (ZRI0A-ZDT*ZRI1A)/ PRMUZ(JL)
151  ZRI0B = ZC1A*ZEXKM+ZC2A*ZEXKP-ZALPHA*ZEXMU0
152  ZRI1B = ZRP*(ZC1A*ZEXKM-ZC2A*ZEXKP)-ZBETA*ZEXMU0
153  PTR1(JL) = ZEXMU0+(ZRI0B+ZDT*ZRI1B)/ PRMUZ(JL)
154
155  PRE1(JL) = MAX(_ZERO_,MIN(PRE1(JL),_ONE_))   ! lower/upper limit (Hubert Gallee, LGGE, 28-03-2008)
156  PTR1(JL) = MAX(_ZERO_,MIN(PTR1(JL),_ONE_))   ! lower/upper limit (Hubert Gallee, LGGE, 28-03-2008)
157
158!*         1.3     WITH REFLECTION FROM THE UNDERLYING LAYER
159
160
161  ZB21 = ZA21- PREF(JL) *ZXP2P*ZEXKM
162  ZB22 = ZA22- PREF(JL) *ZXM2P*ZEXKP
163  ZB23 = ZA23- PREF(JL) *ZEXMU0*(ZAP2B - PRMUZ(JL) )
164  ZDENB = ZA11 * ZB22 - ZB21 * ZA12
165  ZC1B = (ZB22*ZA13-ZA12*ZB23)/ZDENB
166  ZC2B = (ZA11*ZB23-ZB21*ZA13)/ZDENB
167  ZRI0C = ZC1B+ZC2B-ZALPHA
168  ZRI1C = ZRP*(ZC1B-ZC2B)-ZBETA
169  PRE2(JL) = (ZRI0C-ZDT*ZRI1C) / PRMUZ(JL)
170  ZRI0D = ZC1B*ZEXKM + ZC2B*ZEXKP - ZALPHA*ZEXMU0
171  ZRI1D = ZRP * (ZC1B*ZEXKM - ZC2B*ZEXKP) - ZBETA*ZEXMU0
172  PTR2(JL) = ZEXMU0 + (ZRI0D + ZDT*ZRI1D) / PRMUZ(JL)
173
174  PRE2(JL) = MAX(_ZERO_,MIN(PRE2(JL),_ONE_))   ! lower/upper limit (Hubert Gallee, LGGE, 28-03-2008)
175  PTR2(JL) = MAX(_ZERO_,MIN(PTR2(JL),_ONE_))   ! lower/upper limit (Hubert Gallee, LGGE, 28-03-2008)
176
177ENDDO
178RETURN
179END SUBROUTINE SWDE
Note: See TracBrowser for help on using the repository browser.