source: LMDZ6/trunk/libf/phylmd/phyredem.F90 @ 4352

Last change on this file since 4352 was 4298, checked in by pcadule, 21 months ago

modifications to ensure restartability of the model when CO2 tracer is passed to radiative code, and to add diagnostics variables
  • 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: 14.8 KB
Line 
1!
2! $Id: phyredem.F90 4298 2022-10-17 08:15:06Z evignon $
3!
4SUBROUTINE phyredem (fichnom)
5!
6!-------------------------------------------------------------------------------
7! Author: Z.X. Li (LMD/CNRS), 1993/08/18
8!-------------------------------------------------------------------------------
9! Purpose: Write restart state for physics.
10!-------------------------------------------------------------------------------
11  USE dimphy, ONLY: klon, klev
12  USE fonte_neige_mod,  ONLY : fonte_neige_final
13  USE pbl_surface_mod,  ONLY : pbl_surface_final
14  USE phys_state_var_mod, ONLY: radpas, zmasq, pctsrf,                       &
15                                ftsol, beta_aridity, delta_tsurf, falb_dir,  &
16                                falb_dif, qsol, fevap, radsol, solsw, sollw, &
17                                sollwdown, rain_fall, snow_fall, z0m, z0h,   &
18                                agesno, zmea, zstd, zsig, zgam, zthe, zpic,  &
19                                zval, rugoro, t_ancien, q_ancien,            &
20                                prw_ancien, prlw_ancien, prsw_ancien,        &
21                                ql_ancien, qs_ancien, rneb_ancien, u_ancien, &
22                                v_ancien, clwcon, rnebcon, ratqs, pbl_tke,   &
23                                wake_delta_pbl_tke, zmax0, f0, sig1, w01,    &
24                                wake_deltat, wake_deltaq, wake_s, wake_dens, &
25                                awake_dens, cv_gen,                          &
26                                wake_cstar,                                  &
27                                wake_pe, wake_fip, fm_therm, entr_therm,     &
28                                detr_therm, ale_bl, ale_bl_trig, alp_bl,     &
29                                ale_wake, ale_bl_stat,                       &
30                                du_gwd_rando, du_gwd_front, u10m, v10m, &
31                                treedrg, solswfdiff, delta_sal, ds_ns, dt_ns, &
32                                delta_sst, ratqs_inter
33
34  USE geometry_mod, ONLY : longitude_deg, latitude_deg
35  USE iostart, ONLY: open_restartphy, close_restartphy, enddef_restartphy, put_field, put_var
36  USE traclmdz_mod, ONLY : traclmdz_to_restart
37  USE infotrac_phy, ONLY: type_trac, types_trac, nqtot, tracers, nbtr
38  USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, co2_send, carbon_cycle_rad, RCO2_glo
39  USE indice_sol_mod, ONLY: nbsrf, is_oce, is_sic, is_ter, is_lic, epsfra
40  USE surface_data, ONLY: type_ocean, version_ocean
41  USE ocean_slab_mod, ONLY : nslay, tslab, seaice, tice, fsic
42  USE time_phylmdz_mod, ONLY: annee_ref, day_end, itau_phy, pdtphys
43  use config_ocean_skin_m, only: activate_ocean_skin 
44
45  IMPLICIT none
46
47  include "dimsoil.h"
48  include "clesphys.h"
49  include "alpale.h"
50  include "compbl.h"
51  !======================================================================
52  CHARACTER*(*) fichnom
53
54  ! les variables globales ecrites dans le fichier restart
55
56  REAL tsoil(klon, nsoilmx, nbsrf)
57  REAL qsurf(klon, nbsrf)
58  REAL snow(klon, nbsrf)
59  real fder(klon)
60  REAL run_off_lic_0(klon)
61  REAL trs(klon, nbtr)
62
63  INTEGER nid, nvarid, idim1, idim2, idim3
64  INTEGER ierr
65  INTEGER length
66  PARAMETER (length=100)
67  REAL tab_cntrl(length)
68
69  INTEGER isoil, nsrf,isw
70  CHARACTER (len=2) :: str2
71  CHARACTER (len=256) :: nam, lnam
72  INTEGER           :: it, iq, pass
73
74  !======================================================================
75
76  ! Get variables which will be written to restart file from module
77  ! pbl_surface_mod
78  CALL pbl_surface_final(fder, snow, qsurf,  tsoil)
79
80  ! Get a variable calculated in module fonte_neige_mod
81  CALL fonte_neige_final(run_off_lic_0)
82
83  !======================================================================
84
85  CALL open_restartphy(fichnom)
86
87 
88  DO ierr = 1, length
89     tab_cntrl(ierr) = 0.0
90  ENDDO
91  tab_cntrl(1) = pdtphys
92  tab_cntrl(2) = radpas
93  ! co2_ppm : current value of atmospheric CO2
94  tab_cntrl(3) = co2_ppm
95  tab_cntrl(4) = solaire
96  tab_cntrl(5) = iflag_con
97  tab_cntrl(6) = nbapp_rad
98
99  IF( iflag_cycle_diurne.GE.1 ) tab_cntrl( 7 ) = iflag_cycle_diurne
100  IF(   soil_model ) tab_cntrl( 8 ) = 1.
101  IF(     new_oliq ) tab_cntrl( 9 ) = 1.
102  IF(     ok_orodr ) tab_cntrl(10 ) = 1.
103  IF(     ok_orolf ) tab_cntrl(11 ) = 1.
104
105  tab_cntrl(13) = day_end
106  tab_cntrl(14) = annee_ref
107  tab_cntrl(15) = itau_phy
108
109  ! co2_ppm0 : initial value of atmospheric CO2
110  ! tab_cntrl(16) = co2_ppm0
111
112  !  PC -- initial value of RCO2 for the radiation scheme
113  !  tab_cntrl(17) = co2_ppm * 1.0e-06 * RMCO2 / RMD
114  IF (carbon_cycle_rad) tab_cntrl(17) = RCO2_glo
115  !PRINT*, "PC : phyredem RCO2_glo =",RCO2_glo
116
117  DO pass=1,2   ! pass=1 netcdf definition ; pass=2 netcdf write
118 
119    CALL put_var(pass, "controle", "Parametres de controle", tab_cntrl)
120
121    CALL put_field(pass,"longitude", &
122         "Longitudes de la grille physique", longitude_deg)
123
124    CALL put_field(pass,"latitude", "Latitudes de la grille physique", latitude_deg)
125
126    ! PB ajout du masque terre/mer
127
128    CALL put_field(pass,"masque", "masque terre mer", zmasq)
129
130    ! BP ajout des fraction de chaque sous-surface
131
132    ! Get last fractions from slab ocean
133    IF (type_ocean == 'slab' .AND. version_ocean == "sicINT") THEN
134        WHERE (1.-zmasq(:).GT.EPSFRA)
135            pctsrf(:,is_oce)=(1.-fsic(:))*(1.-zmasq(:))
136            pctsrf(:,is_sic)=fsic(:)*(1.-zmasq(:))
137        END WHERE
138    END IF
139
140    ! 1. fraction de terre
141
142    CALL put_field(pass,"FTER", "fraction de continent", pctsrf(:, is_ter))
143
144    ! 2. Fraction de glace de terre
145
146    CALL put_field(pass,"FLIC", "fraction glace de terre", pctsrf(:, is_lic))
147    ! 3. fraction ocean
148
149    CALL put_field(pass,"FOCE", "fraction ocean", pctsrf(:, is_oce))
150
151    ! 4. Fraction glace de mer
152
153    CALL put_field(pass,"FSIC", "fraction glace mer", pctsrf(:, is_sic))
154
155    IF(nbsrf>99) THEN
156      PRINT*, "Trop de sous-mailles";  CALL abort_physic("phyredem", "", 1)
157    END IF
158    IF(nsoilmx>99) THEN
159      PRINT*, "Trop de sous-surfaces"; CALL abort_physic("phyredem", "", 1)
160    END IF
161    IF(nsw>99) THEN
162      PRINT*, "Trop de bandes"; CALL abort_physic("phyredem", "", 1)
163    END IF
164
165!    Surface variables
166    CALL put_field_srf1(pass,"TS","Temperature",ftsol(:,:))
167
168    IF (iflag_pbl>1 .AND. iflag_wake>=1  .AND. iflag_pbl_split >=1) then
169       CALL put_field_srf1(pass, "DELTATS", &
170            "w-x surface temperature difference", delta_tsurf(:,:))
171       CALL put_field_srf1(pass,"BETAS","Aridity factor", beta_aridity(:,:))
172    end IF
173!    End surface variables
174
175! ================== Albedo =======================================
176    print*,'PHYREDEM NOUVEAU'
177    CALL put_field_srf2(pass,"A_dir_SW","Albedo direct",falb_dir(:,:,:))
178    CALL put_field_srf2(pass,"A_dif_SW","Albedo diffus",falb_dif(:,:,:))
179
180    CALL put_field_srf1(pass,"U10M", "u a 10m", u10m)
181
182    CALL put_field_srf1(pass,"V10M", "v a 10m", v10m)
183
184
185! ================== Tsoil =========================================
186    CALL put_field_srf2(pass,"Tsoil","Temperature",tsoil(:,:,:))
187!FC
188!  CALL put_field_srf2("treedrg","freinage arbres",treedrg(:,:,:))
189    CALL put_field(pass,"treedrg_ter","freinage arbres",treedrg(:,:,is_ter))
190
191
192    CALL put_field_srf1(pass,"QS"  , "Humidite",qsurf(:,:))
193
194    CALL put_field     (pass,"QSOL", "Eau dans le sol (mm)", qsol)
195
196    CALL put_field_srf1(pass,"EVAP", "Evaporation", fevap(:,:))
197
198    CALL put_field_srf1(pass,"SNOW", "Neige", snow(:,:))
199
200    CALL put_field(pass,"RADS", "Rayonnement net a la surface", radsol)
201
202    CALL put_field(pass,"solsw", "Rayonnement solaire a la surface", solsw)
203
204    CALL put_field(pass,"solswfdiff", "Fraction du rayonnement solaire a la surface qui est diffus", solswfdiff)
205
206    CALL put_field(pass,"sollw", "Rayonnement IF a la surface", sollw)
207
208    CALL put_field(pass,"sollwdown", "Rayonnement down IF a la surface", sollwdown)
209
210    CALL put_field(pass,"fder", "Derive de flux", fder)
211
212    CALL put_field(pass,"rain_f", "precipitation liquide", rain_fall)
213
214    CALL put_field(pass,"snow_f", "precipitation solide", snow_fall)
215
216    CALL put_field_srf1(pass,"Z0m", "rugosite", z0m(:,:))
217
218    CALL put_field_srf1(pass,"Z0h", "rugosite", z0h(:,:))
219
220    CALL put_field_srf1(pass,"AGESNO", "Age de la neige", agesno(:,:))
221
222    CALL put_field(pass,"ZMEA", "ZMEA", zmea)
223
224    CALL put_field(pass,"ZSTD", "ZSTD", zstd)
225
226    CALL put_field(pass,"ZSIG", "ZSIG", zsig)
227
228    CALL put_field(pass,"ZGAM", "ZGAM", zgam)
229
230    CALL put_field(pass,"ZTHE", "ZTHE", zthe)
231
232    CALL put_field(pass,"ZPIC", "ZPIC", zpic)
233
234    CALL put_field(pass,"ZVAL", "ZVAL", zval)
235
236    CALL put_field(pass,"RUGSREL", "RUGSREL", rugoro)
237
238    CALL put_field(pass,"TANCIEN", "TANCIEN", t_ancien)
239
240    CALL put_field(pass,"QANCIEN", "QANCIEN", q_ancien)
241
242    CALL put_field(pass,"QLANCIEN", "QLANCIEN", ql_ancien)
243
244    CALL put_field(pass,"QSANCIEN", "QSANCIEN", qs_ancien)
245
246    CALL put_field(pass,"RNEBANCIEN", "RNEBANCIEN", rneb_ancien)
247
248    CALL put_field(pass,"PRWANCIEN", "PRWANCIEN", prw_ancien)
249
250    CALL put_field(pass,"PRLWANCIEN", "PRLWANCIEN", prlw_ancien)
251
252    CALL put_field(pass,"PRSWANCIEN", "PRSWANCIEN", prsw_ancien)
253
254    CALL put_field(pass,"UANCIEN", "UANCIEN", u_ancien)
255
256    CALL put_field(pass,"VANCIEN", "VANCIEN", v_ancien)
257
258    CALL put_field(pass,"CLWCON", "Eau liquide convective", clwcon)
259
260    CALL put_field(pass,"RNEBCON", "Nebulosite convective", rnebcon)
261
262    CALL put_field(pass,"RATQS", "Ratqs", ratqs)
263
264    ! run_off_lic_0
265
266    CALL put_field(pass,"RUNOFFLIC0", "Runofflic0", run_off_lic_0)
267
268    ! DEB TKE PBL !
269
270    IF (iflag_pbl>1) then
271      CALL put_field_srf3(pass,"TKE", "Energ. Cineti. Turb.", &
272           pbl_tke(:,:,:))
273      CALL put_field_srf3(pass,"DELTATKE", "Del TKE wk/env.", &
274           wake_delta_pbl_tke(:,:,:))
275    END IF
276
277    ! FIN TKE PBL !
278    !IM ajout zmax0, f0, sig1, w01
279    !IM wake_deltat, wake_deltaq, wake_s, wake_cstar, wake_pe, wake_fip
280
281    CALL put_field(pass,"ZMAX0", "ZMAX0", zmax0)
282
283    CALL put_field(pass,"F0", "F0", f0)
284
285    CALL put_field(pass,"sig1", "sig1 Emanuel", sig1)
286
287    CALL put_field(pass,"w01", "w01 Emanuel", w01)
288
289    ! wake_deltat
290    CALL put_field(pass,"WAKE_DELTAT", "WAKE_DELTAT", wake_deltat)
291
292    CALL put_field(pass,"WAKE_DELTAQ", "WAKE_DELTAQ", wake_deltaq)
293
294    CALL put_field(pass,"WAKE_S", "Wake frac. area", wake_s)
295
296    CALL put_field(pass,"WAKE_DENS", "Wake num. /unit area", wake_dens)
297
298    CALL put_field(pass,"AWAKE_DENS", "Active Wake num. /unit area", awake_dens)
299
300    CALL put_field(pass,"CV_GEN", "CB birth rate", cv_gen)
301
302    CALL put_field(pass,"WAKE_CSTAR", "WAKE_CSTAR", wake_cstar)
303
304    CALL put_field(pass,"WAKE_PE", "WAKE_PE", wake_pe)
305
306    CALL put_field(pass,"WAKE_FIP", "WAKE_FIP", wake_fip)
307
308    ! thermiques
309
310    CALL put_field(pass,"FM_THERM", "FM_THERM", fm_therm)
311
312    CALL put_field(pass,"ENTR_THERM", "ENTR_THERM", entr_therm)
313
314    CALL put_field(pass,"DETR_THERM", "DETR_THERM", detr_therm)
315
316    CALL put_field(pass,"ALE_BL", "ALE_BL", ale_bl)
317
318    CALL put_field(pass,"ALE_BL_TRIG", "ALE_BL_TRIG", ale_bl_trig)
319
320    CALL put_field(pass,"ALP_BL", "ALP_BL", alp_bl)
321
322    CALL put_field(pass,"ALE_WAKE", "ALE_WAKE", ale_wake)
323
324    CALL put_field(pass,"ALE_BL_STAT", "ALE_BL_STAT", ale_bl_stat)
325
326
327    ! fisrtilp/clouds
328    CALL put_field(pass,"RATQS_INTER","Relative width of the lsc sugrid scale water",ratqs_inter)
329
330
331    IF (ANY(types_trac == 'co2i') .OR. ANY(types_trac == 'inco')) THEN
332       IF (carbon_cycle_cpl) THEN
333          IF (.NOT. ALLOCATED(co2_send)) THEN
334             ! This is the case of create_etat0_limit, ce0l
335             ALLOCATE(co2_send(klon))
336             co2_send(:) = co2_ppm0
337          END IF
338          CALL put_field(pass,"co2_send", "co2_ppm for coupling", co2_send)
339       END IF
340
341    ! trs from traclmdz_mod
342    ELSE IF (type_trac == 'lmdz') THEN
343       CALL traclmdz_to_restart(trs)
344       it = 0
345       DO iq = 1, nqtot
346          IF(.NOT.(tracers(iq)%isAdvected .AND. tracers(iq)%isInPhysics)) CYCLE
347          it = it+1
348          CALL put_field(pass,"trs_"//tracers(iq)%name, "", trs(:, it))
349       END DO
350    END IF
351
352    ! Restart variables for Slab ocean
353    IF (type_ocean == 'slab') THEN
354        IF (nslay.EQ.1) THEN
355          CALL put_field(pass,"tslab", "Slab ocean temperature", tslab)
356        ELSE
357          DO it=1,nslay
358            WRITE(str2,'(i2.2)') it
359            CALL put_field(pass,"tslab"//str2, "Slab ocean temperature", tslab(:,it))
360          END DO
361        END IF
362        IF (version_ocean == 'sicINT') THEN
363            CALL put_field(pass,"seaice", "Slab seaice (kg/m2)", seaice)
364            CALL put_field(pass,"slab_tice", "Slab sea ice temperature", tice)
365        END IF
366    END IF
367
368    if (ok_gwd_rando) call put_field(pass,"du_gwd_rando", &
369         "tendency on zonal wind due to flott gravity waves", du_gwd_rando)
370
371    IF (.not. ok_hines .and. ok_gwd_rando) call put_field(pass,"du_gwd_front", &
372         "tendency on zonal wind due to acama gravity waves", du_gwd_front)
373
374    if (activate_ocean_skin >= 1) then
375       if (activate_ocean_skin == 2 .and. type_ocean == 'couple') then
376          CALL put_field(pass, "delta_sal", &
377               "ocean-air interface salinity minus bulk salinity", delta_sal)
378          CALL put_field(pass, "delta_SST", &
379               "ocean-air interface temperature minus bulk SST", delta_sst)
380       end if
381       
382       CALL put_field(pass, "dS_ns", "delta salinity near surface", ds_ns)
383       CALL put_field(pass, "dT_ns", "delta temperature near surface", dT_ns)
384    end if
385   
386    IF (pass==1) CALL enddef_restartphy
387    IF (pass==2) CALL close_restartphy
388  ENDDO ! DO pass=1,2   ! pass=1 netcdf definition ; pass=2 netcdf write
389 
390  !$OMP BARRIER
391
392
393  CONTAINS
394
395
396SUBROUTINE put_field_srf1(pass,nam,lnam,field)
397
398  IMPLICIT NONE
399  INTEGER, INTENT(IN)           :: pass
400  CHARACTER(LEN=*),  INTENT(IN) :: nam, lnam
401  REAL,              INTENT(IN) :: field(:,:)
402  CHARACTER(LEN=256) :: nm, lm, str
403  DO nsrf = 1, SIZE(field,2)
404    WRITE(str, '(i2.2)') nsrf
405    nm=TRIM(nam)//TRIM(str)
406    lm=TRIM(lnam)//" de surface No. "//TRIM(str)
407    CALL put_field(pass,nm,lm,field(:,nsrf))
408  END DO
409
410END SUBROUTINE put_field_srf1
411
412
413SUBROUTINE put_field_srf2(pass,nam,lnam,field)
414
415  IMPLICIT NONE
416  INTEGER, INTENT(IN)            :: pass
417  CHARACTER(LEN=*),  INTENT(IN) :: nam, lnam
418  REAL,              INTENT(IN) :: field(:,:,:)
419  CHARACTER(LEN=256) :: nm, lm, str
420  DO nsrf = 1, SIZE(field,3)
421    DO isoil=1, SIZE(field,2)
422      WRITE(str, '(i2.2,"srf",i2.2)')isoil,nsrf
423!      WRITE(lunout,*)"PHYREDEM ",TRIM(nam)//TRIM(str)
424      nm=TRIM(nam)//TRIM(str)
425      lm=TRIM(lnam)//" du sol No. "//TRIM(str)
426      CALL put_field(pass,nm,lm,field(:,isoil,nsrf))
427    END DO
428  END DO
429
430END SUBROUTINE put_field_srf2
431
432
433SUBROUTINE put_field_srf3(pass,nam,lnam,field)
434
435  IMPLICIT NONE
436  INTEGER, INTENT(IN)            :: pass
437  CHARACTER(LEN=*),  INTENT(IN) :: nam, lnam
438  REAL,              INTENT(IN) :: field(:,:,:)
439  CHARACTER(LEN=256) :: nm, lm, str
440  DO nsrf = 1, SIZE(field,3)
441    WRITE(str, '(i2.2)') nsrf
442    nm=TRIM(nam)//TRIM(str)
443    lm=TRIM(lnam)//TRIM(str)
444    CALL put_field(pass,nm,lm,field(:,1:klev+1,nsrf))
445  END DO
446
447END SUBROUTINE put_field_srf3
448
449
450END SUBROUTINE phyredem
Note: See TracBrowser for help on using the repository browser.