source: dynamico_lmdz/simple_physics/phyparam/physics/radiative_lw.F90 @ 4244

Last change on this file since 4244 was 4244, checked in by dubos, 4 years ago

simple_physics : output LW fluxes
File size: 7.0 KB
Line 
1MODULE radiative_lw
2
3#include "use_logging.h"
4
5  IMPLICIT NONE
6  SAVE
7
8  PRIVATE
9
10  PUBLIC :: lw
11
12  LOGICAL, PARAMETER :: lstrong=.TRUE.
13  REAL, PARAMETER :: stephan=5.67e-08
14
15CONTAINS
16
17  SUBROUTINE lw(ngrid,nlayer,coefir,emissiv, &
18       pp,ps_rad,ptsurf,pt,              &
19       pfluxir,pdtlw,                    &
20       lverbose, lwrite)
21    USE phys_const, ONLY : cpp, g
22    USE writefield_mod, ONLY : writefield
23    !=======================================================================
24    !
25    !   calcul de l evolution de la temperature sous l effet du rayonnement
26    !   infra-rouge.
27    !   Pour simplifier, les transmissions sont precalculees et ne
28    !   dependent que de l altitude.
29    !
30    !   arguments:
31    !   ----------
32    !
33    !   entree:
34    !   -------
35    !      ngrid             nombres de points de la grille horizontale
36    !      nlayer              nombre de couches
37    !      ptsurf(ngrid)     temperature de la surface
38    !      pt(ngrid,nlayer)    temperature des couches
39    !      pp(ngrid,nlayer+1)  pression entre les couches
40    !      lwrite            variable logique pour sorties
41    !
42    !   sortie:
43    !   -------
44    !      pdtlw(ngrid,nlayer) taux de refroidissement
45    !      pfluxir(ngrid)    flux infrarouge sur le sol
46    !
47    !=======================================================================
48
49    !   declarations:
50    !   -------------
51
52    !   arguments:
53    !   ----------
54
55    INTEGER, INTENT(IN)  ::  ngrid,nlayer
56    REAL,    INTENT(IN)  :: coefir,emissiv(ngrid),ps_rad
57    REAL,    INTENT(IN)  ::  ptsurf(ngrid),pt(ngrid,nlayer),pp(ngrid,nlayer+1)
58    REAL,    INTENT(OUT) ::  pdtlw(ngrid,nlayer),pfluxir(ngrid)
59    LOGICAL, INTENT(IN)  :: lwrite, lverbose
60
61    !   variables locales:
62    !   ------------------
63
64    INTEGER nlevel,ilev,ig,i,il
65    REAL zplanck(ngrid,nlayer+1),zcoef
66    REAL zfluxup(ngrid,nlayer+1),zfluxdn(ngrid,nlayer+1)
67    REAL zflux(ngrid,nlayer+1)
68    REAL zlwtr1(ngrid),zlwtr2(ngrid)
69    REAL zup(ngrid,nlayer+1),zdup(ngrid)
70
71    CHARACTER(6), PARAMETER :: tag='rad/lw'
72    !-----------------------------------------------------------------------
73    !   initialisations:
74    !   ----------------
75
76    nlevel=nlayer+1
77
78    !-----------------------------------------------------------------------
79    !   2. calcul des quantites d absorbants:
80    !   -------------------------------------
81
82    !   absorption forte
83    IF(lstrong) THEN
84       DO ilev=1,nlevel
85          DO ig=1,ngrid
86             zup(ig,ilev)=pp(ig,ilev)*pp(ig,ilev)/(2.*g)
87          ENDDO
88       ENDDO
89       IF(lverbose) THEN
90          DO ilev=1,nlayer
91             WRITELOG(*,*) ' up(',ilev,')  =  ',zup(ngrid/2+1,ilev)
92          ENDDO
93          LOG_DBG(tag)
94       ENDIF
95       zcoef=-log(coefir)/sqrt(ps_rad*ps_rad/(2.*g))
96
97       !   absorption faible
98    ELSE
99       DO ilev=1,nlevel
100          DO ig=1,ngrid
101             zup(ig,ilev)=pp(ig,ilev)
102          ENDDO
103       ENDDO
104       zcoef=-log(coefir)/ps_rad
105    ENDIF
106
107
108    !-----------------------------------------------------------------------
109    !   2. calcul de la fonction de corps noir:
110    !   ---------------------------------------
111
112    DO ilev=1,nlayer
113       DO ig=1,ngrid
114          zplanck(ig,ilev)=pt(ig,ilev)*pt(ig,ilev)
115          zplanck(ig,ilev)=stephan* &
116               zplanck(ig,ilev)*zplanck(ig,ilev)
117       ENDDO
118    ENDDO
119
120    !-----------------------------------------------------------------------
121    !   4. flux descendants:
122    !   --------------------
123
124    DO ilev=1,nlayer
125       DO ig=1,ngrid
126          zfluxdn(ig,ilev)=0.
127       ENDDO
128       DO ig=1,ngrid
129          zdup(ig)=zup(ig,ilev)-zup(ig,nlevel)
130       ENDDO
131       CALL lwtr(ngrid,zcoef,lstrong,zdup,zlwtr1)
132
133       DO il=nlayer,ilev,-1
134          zlwtr2(:)=zlwtr1(:)
135          DO ig=1,ngrid
136             zdup(ig)=zup(ig,ilev)-zup(ig,il)
137          ENDDO
138          CALL lwtr(ngrid,zcoef,lstrong,zdup,zlwtr1)
139          DO ig=1,ngrid
140             zfluxdn(ig,ilev)=zfluxdn(ig,ilev)+ &
141                  zplanck(ig,il)*(zlwtr1(ig)-zlwtr2(ig))
142          ENDDO
143       ENDDO
144    ENDDO
145
146    DO ig=1,ngrid
147       zfluxdn(ig,nlevel)=0.
148       pfluxir(ig)=emissiv(ig)*zfluxdn(ig,1)
149    ENDDO
150
151    DO ig=1,ngrid
152       zfluxup(ig,1)=ptsurf(ig)*ptsurf(ig)
153       zfluxup(ig,1)=emissiv(ig)*stephan*zfluxup(ig,1)*zfluxup(ig,1) &
154            +(1.-emissiv(ig))*zfluxdn(ig,1)
155    ENDDO
156
157    !-----------------------------------------------------------------------
158    !   3. flux montants:
159    !   ------------------
160
161    DO ilev=1,nlayer
162       DO ig=1,ngrid
163          zdup(ig)=zup(ig,1)-zup(ig,ilev+1)
164       ENDDO
165       CALL lwtr(ngrid,zcoef,lstrong,zdup,zlwtr1)
166       DO ig=1,ngrid
167          zfluxup(ig,ilev+1)=zfluxup(ig,1)*zlwtr1(ig)
168       ENDDO
169       DO il=1,ilev
170          zlwtr2(:)=zlwtr1(:)
171          DO ig=1,ngrid
172             zdup(ig)=zup(ig,il+1)-zup(ig,ilev+1)
173          ENDDO
174          CALL lwtr(ngrid,zcoef,lstrong,zdup,zlwtr1)
175          DO ig=1,ngrid
176             zfluxup(ig,ilev+1)=zfluxup(ig,ilev+1)+ &
177                  zplanck(ig,il)*(zlwtr1(ig)-zlwtr2(ig))
178          ENDDO
179       ENDDO
180
181    ENDDO
182
183    !-----------------------------------------------------------------------
184    !   5. calcul des flux nets:
185    !   ------------------------
186
187    DO ilev=1,nlevel
188       DO ig=1,ngrid
189          zflux(ig,ilev)=zfluxup(ig,ilev)-zfluxdn(ig,ilev)
190       ENDDO
191    ENDDO
192
193    !-----------------------------------------------------------------------
194    !   6. Calcul des taux de refroidissement:
195    !   --------------------------------------
196
197    DO ilev=1,nlayer
198       DO ig=1,ngrid
199          pdtlw(ig,ilev)=(zflux(ig,ilev+1)-zflux(ig,ilev))* &
200               g/(cpp*(pp(ig,ilev+1)-pp(ig,ilev)))
201       ENDDO
202    ENDDO
203
204    !-----------------------------------------------------------------------
205    !   10. sorties eventuelles:
206    !   ------------------------
207
208    IF (lverbose) THEN
209       WRITELOG(*,*) 'Diagnostique rayonnement thermique'
210       WRITELOG(*,*) 'temperature     ', &
211            'flux montant    flux desc.     taux de refroid.'
212       i=ngrid/2+1
213       WRITELOG(6,'(4e18.4)') ptsurf(i)
214       DO ilev=1,nlayer
215          WRITELOG(6,'(i4,4e18.4)') ilev,pt(i,ilev), &
216               zfluxup(i,ilev),zfluxdn(i,ilev),pdtlw(i,ilev)
217       ENDDO
218       WRITELOG(6,'(4e18.4)') zfluxup(i,nlevel),zfluxdn(i,nlevel)
219       LOG_DBG(tag)
220    ENDIF
221
222    IF(lwrite) THEN
223       CALL writefield('lwflux_up', 'Upward LW flux', 'W/m2', zfluxup)
224       CALL writefield('lwflux_down', 'Downward LW flux', 'W/m2', zfluxdn)
225    END IF
226    !-----------------------------------------------------------------------
227
228  END SUBROUTINE lw
229
230  PURE SUBROUTINE lwtr(ngrid,coef,lstrong,dup,transm)
231    INTEGER, INTENT(IN) :: ngrid
232    REAL,    INTENT(IN) :: coef
233    LOGICAL, INTENT(IN) :: lstrong
234    REAL,    INTENT(IN) :: dup(ngrid)
235    REAL,    INTENT(OUT) :: transm(ngrid)
236    INTEGER ig
237    IF(lstrong) THEN
238       DO ig=1,ngrid
239          transm(ig)=exp(-coef*sqrt(dup(ig)))
240       ENDDO
241    ELSE
242       DO ig=1,ngrid
243          transm(ig)=exp(-coef*dup(ig))
244       ENDDO
245    ENDIF
246
247  END SUBROUTINE lwtr
248
249END MODULE radiative_lw
Note: See TracBrowser for help on using the repository browser.