source: dynamico_lmdz/simple_physics/phyparam/physics/radiative_sw.F90 @ 4241

Last change on this file since 4241 was 4233, checked in by dubos, 6 years ago

simple_physics : enforce F2003 strictly
File size: 7.2 KB
Line 
1MODULE radiative_sw
2
3#include "use_logging.h"
4
5  IMPLICIT NONE
6  SAVE
7
8  PRIVATE
9
10  PUBLIC :: sw
11
12CONTAINS
13
14  PURE SUBROUTINE monGATHER(n,a,b,index)
15    INTEGER, INTENT(IN) ::  n,index(n)
16    REAL, INTENT(IN) :: b(n)
17    REAL, INTENT(OUT) :: a(n)
18    INTEGER :: i
19
20    DO i=1,n
21       a(i)=b(index(i))
22    END DO
23  END SUBROUTINE monGATHER
24
25  PURE subroutine monscatter(n,a,index,b)
26    INTEGER, INTENT(IN) :: n,index(n)
27    REAL, INTENT(IN) :: b(n)
28    REAL, INTENT(OUT) :: a(n)
29    INTEGER :: i
30    DO i=1,n
31       a(index(i))=b(i)
32    END DO
33  end subroutine monscatter
34
35  SUBROUTINE sw(ngrid,nlayer,ldiurn, &
36       coefvis,albedo, &
37       plevel,ps_rad,pmu,pfract,psolarf0, &
38       fsrfvis,dtsw, &
39       lwrite)
40    USE phys_const, ONLY : cpp, g
41    !=======================================================================
42    !
43    !   Rayonnement solaire en atmosphere non diffusante avec un
44    !   coefficient d absoprption gris.
45    !
46    !=======================================================================
47    !
48    !   declarations:
49    !   -------------
50    !
51    !
52    !   arguments:
53    !   ----------
54    !
55    INTEGER, INTENT(IN) :: ngrid,nlayer
56    LOGICAL, INTENT(IN) :: lwrite,ldiurn
57    REAL, INTENT(IN)    :: albedo(ngrid), coefvis
58    REAL, INTENT(IN)    :: pmu(ngrid), pfract(ngrid)
59    REAL, INTENT(IN)    :: plevel(ngrid,nlayer+1), ps_rad
60    REAL, INTENT(IN)    :: psolarf0
61    REAL, INTENT(OUT)   :: fsrfvis(ngrid),dtsw(ngrid,nlayer)
62
63    REAL zalb(ngrid),zmu(ngrid),zfract(ngrid)
64    REAL zplev(ngrid,nlayer+1)
65    REAL zflux(ngrid),zdtsw(ngrid,nlayer)
66
67    INTEGER ig,l,nlevel,index(ngrid),ncount,igout
68    REAL ztrdir(ngrid,nlayer+1),ztrref(ngrid,nlayer+1)
69    REAL zfsrfref(ngrid)
70    REAL z1(ngrid)
71    REAL zu(ngrid,nlayer+1)
72    REAL tau0
73
74    !-----------------------------------------------------------------------
75    !   1. initialisations:
76    !   -------------------
77
78    nlevel=nlayer+1
79
80    !-----------------------------------------------------------------------
81    !   Definitions des tableaux locaux pour les points ensoleilles:
82    !   ------------------------------------------------------------
83
84    IF (ldiurn) THEN
85       ncount=0
86       DO ig=1,ngrid
87          index(ig)=0
88       ENDDO
89       DO ig=1,ngrid
90          IF(pfract(ig).GT.1.e-6) THEN
91             ncount=ncount+1
92             index(ncount)=ig
93          ENDIF
94       ENDDO
95       CALL monGATHER(ncount,zfract,pfract,index)
96       CALL monGATHER(ncount,zmu,pmu,index)
97       CALL monGATHER(ncount,zalb,albedo,index)
98       DO l=1,nlevel
99          CALL monGATHER(ncount,zplev(1,l),plevel(1,l),index)
100       ENDDO
101    ELSE
102       ncount=ngrid
103       zfract(:)=pfract(:)
104       zmu(:)=pmu(:)
105       zalb(:)=albedo(:)
106       zplev(:,:)=plevel(:,:)
107    ENDIF
108
109    !-----------------------------------------------------------------------
110    !   calcul des profondeurs optiques integres depuis p=0:
111    !   ----------------------------------------------------
112
113    tau0=-.5*log(coefvis)
114
115    ! calcul de la partie homogene de l opacite
116    tau0=tau0/ps_rad
117    DO l=1,nlayer+1
118       DO ig=1,ncount
119          zu(ig,l)=tau0*zplev(ig,l)
120       ENDDO
121    ENDDO
122
123    !-----------------------------------------------------------------------
124    !   2. calcul de la transmission depuis le sommet de l atmosphere:
125    !   -----------------------------------------------------------
126
127    DO l=1,nlevel
128       DO ig=1,ncount
129          ztrdir(ig,l)=exp(-zu(ig,l)/zmu(ig))
130       ENDDO
131    ENDDO
132
133    IF (lwrite) THEN
134       igout=ncount/2+1
135       WRITELOG(*,*)
136       WRITELOG(*,*) 'Diagnostique des transmission dans le spectre solaire'
137       WRITELOG(*,*) 'zfract, zmu, zalb'
138       WRITELOG(*,*) zfract(igout), zmu(igout), zalb(igout)
139       WRITELOG(*,*) 'Pression, quantite d abs, transmission'
140       DO l=1,nlayer+1
141          WRITELOG(*,*) zplev(igout,l),zu(igout,l),ztrdir(igout,l)
142       ENDDO
143    ENDIF
144
145    !-----------------------------------------------------------------------
146    !   3. taux de chauffage, ray. solaire direct:
147    !   ------------------------------------------
148
149    DO l=1,nlayer
150       DO ig=1,ncount
151          zdtsw(ig,l)=g*psolarf0*zfract(ig)*zmu(ig)*      &
152               (ztrdir(ig,l+1)-ztrdir(ig,l))/    &
153               (cpp*(zplev(ig,l)-zplev(ig,l+1)))
154       ENDDO
155    ENDDO
156    IF (lwrite) THEN
157       WRITELOG(*,*)
158       WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires:'
159       WRITELOG(*,*) ' 1 taux de chauffage lie au ray. solaire  direct'
160       DO l=1,nlayer
161          WRITELOG(*,*) zdtsw(igout,l)
162       ENDDO
163    ENDIF
164
165
166    !-----------------------------------------------------------------------
167    !   4. calcul du flux solaire arrivant sur le sol:
168    !   ----------------------------------------------
169
170    DO ig=1,ncount
171       z1(ig)=zfract(ig)*zmu(ig)*psolarf0*ztrdir(ig,1)
172       zflux(ig)=(1.-zalb(ig))*z1(ig)
173       zfsrfref(ig)=    zalb(ig)*z1(ig)
174    ENDDO
175    IF (lwrite) THEN
176       WRITELOG(*,*)
177       WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires:'
178       WRITELOG(*,*) ' 2 flux solaire net incident sur le sol'
179       WRITELOG(*,*) zflux(igout)
180    ENDIF
181
182
183    !-----------------------------------------------------------------------
184    !   5.calcul des traansmissions depuis le sol, cas diffus:
185    !   ------------------------------------------------------
186
187    DO l=1,nlevel
188       DO ig=1,ncount
189          ztrref(ig,l)=exp(-(zu(ig,1)-zu(ig,l))*1.66)
190       ENDDO
191    ENDDO
192
193    IF (lwrite) THEN
194       WRITELOG(*,*)
195       WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires'
196       WRITELOG(*,*) ' 3 transmission avec les sol'
197       WRITELOG(*,*) 'niveau     transmission'
198       DO l=1,nlevel
199          WRITELOG(*,*) l,ztrref(igout,l)
200       ENDDO
201    ENDIF
202
203    !-----------------------------------------------------------------------
204    !   6.ajout a l echauffement de la contribution du ray. sol. reflechit:
205    !   -------------------------------------------------------------------
206
207    DO l=1,nlayer
208       DO ig=1,ncount
209          zdtsw(ig,l)=zdtsw(ig,l)+ &
210               g*zfsrfref(ig)*(ztrref(ig,l+1)-ztrref(ig,l))/ &
211               (cpp*(zplev(ig,l+1)-zplev(ig,l)))
212       ENDDO
213    ENDDO
214
215    !-----------------------------------------------------------------------
216    !   10. sorites eventuelles:
217    !   ------------------------
218
219    IF (lwrite) THEN
220       WRITELOG(*,*)
221       WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires:'
222       WRITELOG(*,*) ' 3 taux de chauffage total'
223       DO l=1,nlayer
224          WRITELOG(*,*) zdtsw(igout,l)
225       ENDDO
226    ENDIF
227
228    IF (ldiurn) THEN
229       fsrfvis(:)=0.
230       CALL monscatter(ncount,fsrfvis,index,zflux)
231       dtsw(:,:)=0.
232       DO l=1,nlayer
233          CALL monscatter(ncount,dtsw(1,l),index,zdtsw(1,l))
234       ENDDO
235    ELSE
236       WRITELOG(*,*) 'NOT DIURNE'
237       fsrfvis(:)=zflux(:)
238       dtsw(:,:)=zdtsw(:,:)
239    ENDIF
240    !        call dump2d(iim,jjm-1,zflux(2),'ZFLUX      ')
241    !        call dump2d(iim,jjm-1,fsrfvis(2),'FSRVIS     ')
242    !        call dump2d(iim,jjm-1,ztrdir(2,1),'ztrdir     ')
243    !        call dump2d(iim,jjm-1,pmu(2),'pmu        ')
244    !        call dump2d(iim,jjm-1,pfract(2),'pfract     ')
245    !        call dump2d(iim,jjm-1,albedo(2),'albedo     ')
246    !        call dump2d(iim,jjm-1,ztrdir(2,1),'ztrdir     ')
247    LOG_INFO('rad_sw')
248
249  END SUBROUTINE sw
250
251END MODULE radiative_sw
Note: See TracBrowser for help on using the repository browser.