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

Last change on this file since 4239 was 4233, checked in by dubos, 5 years ago

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