source: LMDZ6/branches/contrails/libf/phylmd/phystokenc_mod.f90 @ 5452

Last change on this file since 5452 was 5268, checked in by abarral, 2 months ago

.f90 <-> .F90 depending on cpp key use

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.8 KB
Line 
1!
2! $Id: phystokenc_mod.f90 5268 2024-10-23 17:02:39Z aborella $
3!
4MODULE phystokenc_mod
5
6  IMPLICIT NONE
7
8  LOGICAL,SAVE :: offline
9!$OMP THREADPRIVATE(offline)
10  INTEGER,SAVE :: istphy
11!$OMP THREADPRIVATE(istphy)
12
13
14CONTAINS
15
16  SUBROUTINE init_phystokenc(offline_dyn,istphy_dyn)
17    IMPLICIT NONE
18    LOGICAL,INTENT(IN) :: offline_dyn
19    INTEGER,INTENT(IN) :: istphy_dyn
20
21    offline=offline_dyn
22    istphy=istphy_dyn
23
24  END SUBROUTINE init_phystokenc
25
26SUBROUTINE phystokenc (nlon,nlev,pdtphys,rlon,rlat, &
27     pt,pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
28     pfm_therm,pentr_therm, &
29     cdragh, pcoefh,yu1,yv1,ftsol,pctsrf, &
30     frac_impa,frac_nucl, &
31     pphis,paire,dtime,itap, &
32     psh, pda, pphi, pmp, pupwd, pdnwd)
33 
34  USE ioipsl
35  USE dimphy
36  USE infotrac_phy, ONLY : nqtot
37  USE iophy
38  USE indice_sol_mod
39  USE print_control_mod, ONLY: lunout
40  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
41 
42  IMPLICIT NONE
43 
44!======================================================================
45! Auteur(s) FH
46! Objet: Ecriture des variables pour transport offline
47!
48!  Note (A Cozic - July 2024): when coupled with inca, offline fields are no
49!  longer calculated in this routine but directly in the physics code.
50!======================================================================
51
52! Arguments:
53!
54  REAL,DIMENSION(klon,klev), INTENT(IN)     :: psh   ! humidite specifique
55  REAL,DIMENSION(klon,klev), INTENT(IN)     :: pda
56  REAL,DIMENSION(klon,klev,klev), INTENT(IN):: pphi
57  REAL,DIMENSION(klon,klev), INTENT(IN)     :: pmp
58  REAL,DIMENSION(klon,klev), INTENT(IN)     :: pupwd ! saturated updraft mass flux
59  REAL,DIMENSION(klon,klev), INTENT(IN)     :: pdnwd ! saturated downdraft mass flux
60
61!   EN ENTREE:
62!   ==========
63!
64!   divers:
65!   -------
66!
67  INTEGER nlon ! nombre de points horizontaux
68  INTEGER nlev ! nombre de couches verticales
69  REAL pdtphys ! pas d'integration pour la physique (seconde)
70  INTEGER itap
71  INTEGER, SAVE :: physid
72!$OMP THREADPRIVATE(physid)
73
74!   convection:
75!   -----------
76!
77  REAL pmfu(klon,klev)  ! flux de masse dans le panache montant
78  REAL pmfd(klon,klev)  ! flux de masse dans le panache descendant
79  REAL pen_u(klon,klev) ! flux entraine dans le panache montant
80  REAL pde_u(klon,klev) ! flux detraine dans le panache montant
81  REAL pen_d(klon,klev) ! flux entraine dans le panache descendant
82  REAL pde_d(klon,klev) ! flux detraine dans le panache descendant
83  REAL pt(klon,klev)
84  REAL,ALLOCATABLE,SAVE :: t(:,:)
85!$OMP THREADPRIVATE(t)
86!
87  REAL rlon(klon), rlat(klon), dtime
88  REAL zx_tmp_3d(nbp_lon,nbp_lat,klev),zx_tmp_2d(nbp_lon,nbp_lat)
89
90!   Couche limite:
91!   --------------
92!
93  REAL cdragh(klon)          ! cdrag
94  REAL pcoefh(klon,klev)     ! coeff melange CL
95  REAL pcoefh_buf(klon,klev) ! coeff melange CL + cdrag
96  REAL yv1(klon)
97  REAL yu1(klon),pphis(klon),paire(klon)
98
99!   Les Thermiques : (Abderr 25 11 02)
100!   ---------------
101  REAL, INTENT(IN) ::  pfm_therm(klon,klev+1)
102  REAL pentr_therm(klon,klev)
103 
104  REAL,ALLOCATABLE,SAVE :: entr_therm(:,:)
105  REAL,ALLOCATABLE,SAVE :: fm_therm(:,:)
106!$OMP THREADPRIVATE(entr_therm)
107!$OMP THREADPRIVATE(fm_therm)
108!
109!   Lessivage:
110!   ----------
111!
112  REAL frac_impa(klon,klev)
113  REAL frac_nucl(klon,klev)
114!
115! Arguments necessaires pour les sources et puits de traceur
116!
117  REAL ftsol(klon,nbsrf)  ! Temperature du sol (surf)(Kelvin)
118  REAL pctsrf(klon,nbsrf) ! Pourcentage de sol f(nature du sol)
119!======================================================================
120!
121  INTEGER i, k, kk
122  REAL,ALLOCATABLE,SAVE :: mfu(:,:)  ! flux de masse dans le panache montant
123  REAL,ALLOCATABLE,SAVE :: mfd(:,:)  ! flux de masse dans le panache descendant
124  REAL,ALLOCATABLE,SAVE :: en_u(:,:) ! flux entraine dans le panache montant
125  REAL,ALLOCATABLE,SAVE :: de_u(:,:) ! flux detraine dans le panache montant
126  REAL,ALLOCATABLE,SAVE :: en_d(:,:) ! flux entraine dans le panache descendant
127  REAL,ALLOCATABLE,SAVE :: de_d(:,:) ! flux detraine dans le panache descendant
128  REAL,ALLOCATABLE,SAVE :: coefh(:,:) ! flux detraine dans le panache descendant
129 
130  REAL,ALLOCATABLE,SAVE :: pyu1(:)
131  REAL,ALLOCATABLE,SAVE :: pyv1(:)
132  REAL,ALLOCATABLE,SAVE :: pftsol(:,:)
133  REAL,ALLOCATABLE,SAVE :: ppsrf(:,:)
134!$OMP THREADPRIVATE(mfu,mfd,en_u,de_u,en_d,de_d,coefh)
135!$OMP THREADPRIVATE(pyu1,pyv1,pftsol,ppsrf)
136
137
138  REAL,DIMENSION(:,:), ALLOCATABLE,SAVE     :: sh 
139  REAL,DIMENSION(:,:), ALLOCATABLE,SAVE     :: da
140  REAL,DIMENSION(:,:,:), ALLOCATABLE,SAVE   :: phi
141  REAL,DIMENSION(:,:), ALLOCATABLE,SAVE     :: mp
142  REAL,DIMENSION(:,:), ALLOCATABLE,SAVE     :: upwd
143  REAL,DIMENSION(:,:), ALLOCATABLE,SAVE     :: dnwd
144 
145  REAL, SAVE :: dtcum
146  INTEGER, SAVE:: iadvtr=0
147!$OMP THREADPRIVATE(dtcum,iadvtr)
148  REAL zmin,zmax
149  LOGICAL ok_sync
150  CHARACTER(len=12) :: nvar
151  logical, parameter :: lstokenc=.FALSE.
152!
153!======================================================================
154
155  iadvtr=iadvtr+1
156
157! Dans le meme vecteur on recombine le drag et les coeff d'echange
158  pcoefh_buf(:,1)      = cdragh(:)
159  pcoefh_buf(:,2:klev) = pcoefh(:,2:klev)
160 
161  ok_sync = .TRUE.
162
163! Initialization done only once
164!======================================================================
165  IF (iadvtr==1) THEN
166     ALLOCATE( t(klon,klev))
167     ALLOCATE( mfu(klon,klev)) 
168     ALLOCATE( mfd(klon,klev)) 
169     ALLOCATE( en_u(klon,klev))
170     ALLOCATE( de_u(klon,klev))
171     ALLOCATE( en_d(klon,klev))
172     ALLOCATE( de_d(klon,klev))
173     ALLOCATE( coefh(klon,klev))
174     ALLOCATE( entr_therm(klon,klev))
175     ALLOCATE( fm_therm(klon,klev))
176     ALLOCATE( pyu1(klon))
177     ALLOCATE( pyv1(klon))
178     ALLOCATE( pftsol(klon,nbsrf))
179     ALLOCATE( ppsrf(klon,nbsrf))
180     
181     ALLOCATE(sh(klon,klev))
182     ALLOCATE(da(klon,klev))
183     ALLOCATE(phi(klon,klev,klev))
184     ALLOCATE(mp(klon,klev))
185     ALLOCATE(upwd(klon,klev))
186     ALLOCATE(dnwd(klon,klev))
187
188     CALL initphysto('phystoke', dtime, dtime*istphy,dtime*istphy,physid)
189     
190     ! Write field phis and aire only once
191     CALL histwrite_phy(physid,lstokenc,"phis",itap,pphis)
192     CALL histwrite_phy(physid,lstokenc,"aire",itap,paire)
193     CALL histwrite_phy(physid,lstokenc,"longitudes",itap,rlon)
194     CALL histwrite_phy(physid,lstokenc,"latitudes",itap,rlat)
195
196  END IF
197 
198 
199! Set to zero cumulating fields
200!======================================================================
201  IF (MOD(iadvtr,istphy)==1.OR.istphy==1) THEN
202     WRITE(lunout,*)'reinitialisation des champs cumules a iadvtr=',iadvtr
203     mfu(:,:)=0.
204     mfd(:,:)=0.
205     en_u(:,:)=0.
206     de_u(:,:)=0.
207     en_d(:,:)=0.
208     de_d(:,:)=0.
209     coefh(:,:)=0.
210     t(:,:)=0.
211     fm_therm(:,:)=0.
212     entr_therm(:,:)=0.
213     pyv1(:)=0.
214     pyu1(:)=0.
215     pftsol(:,:)=0.
216     ppsrf(:,:)=0.
217     sh(:,:)=0.
218     da(:,:)=0.
219     phi(:,:,:)=0.
220     mp(:,:)=0.
221     upwd(:,:)=0.
222     dnwd(:,:)=0.
223     dtcum=0.
224  ENDIF
225 
226
227! Cumulate fields at each time step
228!======================================================================
229  DO k=1,klev
230     DO i=1,klon
231        mfu(i,k)=mfu(i,k)+pmfu(i,k)*pdtphys
232        mfd(i,k)=mfd(i,k)+pmfd(i,k)*pdtphys
233        en_u(i,k)=en_u(i,k)+pen_u(i,k)*pdtphys
234        de_u(i,k)=de_u(i,k)+pde_u(i,k)*pdtphys
235        en_d(i,k)=en_d(i,k)+pen_d(i,k)*pdtphys
236        de_d(i,k)=de_d(i,k)+pde_d(i,k)*pdtphys
237        coefh(i,k)=coefh(i,k)+pcoefh_buf(i,k)*pdtphys
238        t(i,k)=t(i,k)+pt(i,k)*pdtphys
239        fm_therm(i,k)=fm_therm(i,k)+pfm_therm(i,k)*pdtphys
240        entr_therm(i,k)=entr_therm(i,k)+pentr_therm(i,k)*pdtphys
241        sh(i,k) = sh(i,k) + psh(i,k)*pdtphys
242        da(i,k) = da(i,k) + pda(i,k)*pdtphys
243        mp(i,k) = mp(i,k) + pmp(i,k)*pdtphys
244        upwd(i,k) = upwd(i,k) + pupwd(i,k)*pdtphys
245        dnwd(i,k) = dnwd(i,k) + pdnwd(i,k)*pdtphys
246     ENDDO
247  ENDDO
248
249  DO kk=1,klev
250     DO k=1,klev
251        DO i=1,klon
252           phi(i,k,kk) = phi(i,k,kk) + pphi(i,k,kk)*pdtphys
253        END DO
254     END DO
255  END DO
256
257  DO i=1,klon
258     pyv1(i)=pyv1(i)+yv1(i)*pdtphys
259     pyu1(i)=pyu1(i)+yu1(i)*pdtphys
260  END DO
261  DO k=1,nbsrf
262     DO i=1,klon
263        pftsol(i,k)=pftsol(i,k)+ftsol(i,k)*pdtphys
264        ppsrf(i,k)=ppsrf(i,k)+pctsrf(i,k)*pdtphys
265     ENDDO
266  ENDDO
267 
268! Add time step to cumulated time
269  dtcum=dtcum+pdtphys
270 
271
272! Write fields to file, if it is time to do so
273!======================================================================
274  IF(MOD(iadvtr,istphy)==0) THEN
275
276     ! normalize with time period
277     DO k=1,klev
278        DO i=1,klon
279           mfu(i,k)=mfu(i,k)/dtcum
280           mfd(i,k)=mfd(i,k)/dtcum
281           en_u(i,k)=en_u(i,k)/dtcum
282           de_u(i,k)=de_u(i,k)/dtcum
283           en_d(i,k)=en_d(i,k)/dtcum
284           de_d(i,k)=de_d(i,k)/dtcum
285           coefh(i,k)=coefh(i,k)/dtcum
286           t(i,k)=t(i,k)/dtcum 
287           fm_therm(i,k)=fm_therm(i,k)/dtcum
288           entr_therm(i,k)=entr_therm(i,k)/dtcum
289           sh(i,k)=sh(i,k)/dtcum
290           da(i,k)=da(i,k)/dtcum
291           mp(i,k)=mp(i,k)/dtcum
292           upwd(i,k)=upwd(i,k)/dtcum
293           dnwd(i,k)=dnwd(i,k)/dtcum
294        ENDDO
295     ENDDO
296     DO kk=1,klev
297        DO k=1,klev
298           DO i=1,klon
299              phi(i,k,kk) = phi(i,k,kk)/dtcum
300           END DO
301        END DO
302     END DO
303     DO i=1,klon
304        pyv1(i)=pyv1(i)/dtcum
305        pyu1(i)=pyu1(i)/dtcum
306     END DO
307     DO k=1,nbsrf
308        DO i=1,klon
309           pftsol(i,k)=pftsol(i,k)/dtcum
310           ppsrf(i,k)=ppsrf(i,k)/dtcum
311        ENDDO
312     ENDDO
313
314     ! write fields
315     CALL histwrite_phy(physid,lstokenc,"t",itap,t)
316     CALL histwrite_phy(physid,lstokenc,"mfu",itap,mfu)
317     CALL histwrite_phy(physid,lstokenc,"mfd",itap,mfd)
318     CALL histwrite_phy(physid,lstokenc,"en_u",itap,en_u)
319     CALL histwrite_phy(physid,lstokenc,"de_u",itap,de_u)
320     CALL histwrite_phy(physid,lstokenc,"en_d",itap,en_d)
321     CALL histwrite_phy(physid,lstokenc,"de_d",itap,de_d)
322     CALL histwrite_phy(physid,lstokenc,"coefh",itap,coefh)     
323     CALL histwrite_phy(physid,lstokenc,"fm_th",itap,fm_therm)
324     CALL histwrite_phy(physid,lstokenc,"en_th",itap,entr_therm)
325     CALL histwrite_phy(physid,lstokenc,"frac_impa",itap,frac_impa)
326     CALL histwrite_phy(physid,lstokenc,"frac_nucl",itap,frac_nucl)
327     CALL histwrite_phy(physid,lstokenc,"pyu1",itap,pyu1)
328     CALL histwrite_phy(physid,lstokenc,"pyv1",itap,pyv1)
329     CALL histwrite_phy(physid,lstokenc,"ftsol1",itap,pftsol(:,1))
330     CALL histwrite_phy(physid,lstokenc,"ftsol2",itap,pftsol(:,2))
331     CALL histwrite_phy(physid,lstokenc,"ftsol3",itap,pftsol(:,3))
332     CALL histwrite_phy(physid,lstokenc,"ftsol4",itap,pftsol(:,4))
333     CALL histwrite_phy(physid,lstokenc,"psrf1",itap,ppsrf(:,1))
334     CALL histwrite_phy(physid,lstokenc,"psrf2",itap,ppsrf(:,2))
335     CALL histwrite_phy(physid,lstokenc,"psrf3",itap,ppsrf(:,3))
336     CALL histwrite_phy(physid,lstokenc,"psrf4",itap,ppsrf(:,4))
337     CALL histwrite_phy(physid,lstokenc,"sh",itap,sh)
338     CALL histwrite_phy(physid,lstokenc,"da",itap,da)
339     CALL histwrite_phy(physid,lstokenc,"mp",itap,mp)
340     CALL histwrite_phy(physid,lstokenc,"upwd",itap,upwd)
341     CALL histwrite_phy(physid,lstokenc,"dnwd",itap,dnwd)
342
343
344! phi
345     DO k=1,klev
346        IF (k<10) THEN
347           WRITE(nvar,'(i1)') k
348        ELSE IF (k<100) THEN
349           WRITE(nvar,'(i2)') k
350        ELSE
351           WRITE(nvar,'(i3)') k
352        END IF
353        nvar='phi_lev'//trim(nvar)
354       
355        CALL histwrite_phy(physid,lstokenc,nvar,itap,phi(:,:,k))
356     END DO
357     
358     ! Syncronize file
359!$OMP MASTER
360     IF (ok_sync) CALL histsync(physid)
361!$OMP END MASTER
362     
363     
364     ! Calculate min and max values for some fields (coefficients de lessivage)
365     zmin=1e33
366     zmax=-1e33
367     DO k=1,klev
368        DO i=1,klon
369           zmax=MAX(zmax,frac_nucl(i,k))
370           zmin=MIN(zmin,frac_nucl(i,k))
371        ENDDO
372     ENDDO
373     WRITE(lunout,*)'------ coefs de lessivage (min et max) --------'
374     WRITE(lunout,*)'facteur de nucleation ',zmin,zmax
375     zmin=1e33
376     zmax=-1e33
377     DO k=1,klev
378        DO i=1,klon
379           zmax=MAX(zmax,frac_impa(i,k))
380           zmin=MIN(zmin,frac_impa(i,k))
381        ENDDO
382     ENDDO
383     WRITE(lunout,*)'facteur d impaction ',zmin,zmax
384     
385  ENDIF ! IF(MOD(iadvtr,istphy)==0)
386
387END SUBROUTINE phystokenc
388
389END MODULE phystokenc_mod
Note: See TracBrowser for help on using the repository browser.