source: trunk/LMDZ.COMMON/libf/evolution/planet.F90

Last change on this file was 4180, checked in by jbclement, 16 hours ago

PEM:

  • Rework layering-related logic, especially clarify interactions between surface and subsurface water tendencies and disable CO2 lag resistance (inconsistent without updating pressure and mass balance + PCM).
  • Prevent simultaneous activation of layering and ice flows.
  • Add warning when flux_geo /= 0 while soil is disabled.
  • Add new utility function "is_lvl_enabled" for displaying.
  • Replace deprecated 'minieps' with 'eps'/'tol'.

JBC
File size: 14.4 KB
Line 
1MODULE planet
2!-----------------------------------------------------------------------
3! NAME
4!     planet
5!
6! DESCRIPTION
7!     Centralized climate-state storage and lifecycle helpers for PEM.
8!
9! AUTHORS & DATE
10!     JB Clement, 03/2026
11!
12! NOTES
13!     Ownership criterion for declarations:
14!       - Declare in module "planet" variables that are persistent climate
15!         state or persistent references used across multiple time steps.
16!       - Declare in program "pem.F90" transient workflow/control varaibles,
17!         one-shot initialization buffers and per-iteration working buffers.
18!-----------------------------------------------------------------------
19
20! DEPENDENCIES
21! ------------
22use numerics,         only: dp, qp, k4
23use layered_deposits, only: layering, del_layering
24
25! DECLARATION
26! -----------
27implicit none
28
29! PARAMETERS
30! ----------
31! Pressure-related:
32real(dp), dimension(:),   allocatable :: ps_avg          ! Average surface pressure [Pa]
33real(dp), dimension(:,:), allocatable :: ps_ts           ! Surface pressure timeseries [Pa]
34real(dp), dimension(:),   allocatable :: ps_dev          ! Deviation of surface pressure [Pa]
35real(dp)                              :: ps_avg_glob_ini ! Global average pressure at initialization [Pa]
36real(dp)                              :: ps_avg_glob_old ! Global average pressure of previous time step [Pa]
37real(dp)                              :: ps_avg_glob     ! Global average pressure of current time step [Pa]
38
39! Ice-related:
40real(dp),    dimension(:,:), allocatable :: h2o_ice                    ! H2O ice [kg/m2]
41real(dp),    dimension(:,:), allocatable :: co2_ice                    ! CO2 ice [kg/m2]
42real(dp)                                 :: h2oice_sublim_coverage_ini ! Initial surface area of sublimating H2O ice [m2]
43real(dp)                                 :: co2ice_sublim_coverage_ini ! Initial surface area of sublimating CO2 ice [m2]
44logical(k4), dimension(:,:), allocatable :: is_h2oice_ini              ! Initial location of H2O ice
45logical(k4), dimension(:,:), allocatable :: is_co2ice_ini              ! Initial location of CO2 ice
46logical(k4), dimension(:,:), allocatable :: is_co2ice_disappeared      ! Flag to check if CO2 ice disappeared at the previous timestep
47
48! Surface-related:
49real(dp), dimension(:,:), allocatable :: tsurf_avg           ! Average surface temperature [K]
50real(dp), dimension(:,:), allocatable :: tsurf_dev           ! Deviation of surface temperature [K]
51real(dp), dimension(:,:), allocatable :: h2o_surfdensity_avg ! Average water surface density [kg/m^3]
52
53! Soil-related:
54real(dp), dimension(:,:,:),   allocatable :: tsoil_avg    ! Average soil temperature [K]
55real(dp), dimension(:,:,:),   allocatable :: tsoil_dev    ! Deviation of soil temperature [K]
56real(dp), dimension(:,:,:,:), allocatable :: tsoil_ts     ! Soil temperature timeseries [K]
57real(dp), dimension(:,:,:,:), allocatable :: tsoil_ts_old ! Soil temperature timeseries at the previous time step [K]
58
59! Layering-related:
60type(layering), dimension(:,:), allocatable :: layerings_map ! Layering for each grid point and slope
61
62! Sorption-related:
63real(dp), dimension(:,:,:), allocatable :: h2o_soildensity_avg ! Average of soil water soil density [kg/m^3]
64real(dp), dimension(:),     allocatable :: delta_co2_ads       ! Quantity of CO2 exchanged due to adsorption/desorption [kg/m^2]
65real(dp), dimension(:),     allocatable :: delta_h2o_ads       ! Quantity of H2O exchanged due to adsorption/desorption [kg/m^2]
66real(dp), dimension(:,:,:), allocatable :: h2o_ads_reg         ! H2O adsorbed in the regolith [kg/m^2]
67real(dp), dimension(:,:,:), allocatable :: co2_ads_reg         ! CO2 adsorbed in the regolith [kg/m^2]
68
69! Ice table-related:
70real(dp), dimension(:,:),   allocatable :: icetable_depth     ! Depth of the ice table [m]
71real(dp), dimension(:,:),   allocatable :: icetable_thickness ! Thickness of the ice table [m]
72real(dp), dimension(:,:,:), allocatable :: ice_porefilling    ! Amount of porefilling [m^3/m^3]
73real(dp), dimension(:,:),   allocatable :: icetable_depth_old ! Old depth of the ice table [m]
74real(dp), dimension(:),     allocatable :: delta_icetable     ! Total mass of the H2O exchanged with the ice table [kg]
75real(dp), dimension(:,:),   allocatable :: flux_ssice_avg     ! Average of total flux exchanged with subsurface ice [kg/m2/y]
76
77! Tracer-related:
78real(dp), dimension(:,:), allocatable :: q_co2_ts     ! CO2 mass mixing ratio in the first layer [kg/kg]
79real(dp), dimension(:,:), allocatable :: q_co2_ts_ini ! Initial CO2 mass mixing ratio in the first layer [kg/kg]
80real(dp), dimension(:,:), allocatable :: q_h2o_ts     ! H2O mass mixing ratio in the first layer [kg/kg]
81
82! Tendency-related:
83real(dp), dimension(:,:), allocatable :: d_co2ice     ! Tendency of perennial CO2 ice [kg/m2/y]
84real(dp), dimension(:,:), allocatable :: d_co2ice_ini ! Tendency of perennial CO2 ice at the beginning [kg/m2/y]
85real(dp), dimension(:,:), allocatable :: d_h2oice     ! Tendency of perennial H2O ice [kg/m2/y]
86
87contains
88!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
89
90!=======================================================================
91SUBROUTINE ini_planet()
92!-----------------------------------------------------------------------
93! NAME
94!     ini_planet
95!
96! DESCRIPTION
97!     Initialize scalar state values of module planet.
98!
99! AUTHORS & DATE
100!     JB Clement, 03/2026
101!
102! NOTES
103!
104!-----------------------------------------------------------------------
105
106! DECLARATION
107! -----------
108implicit none
109
110! CODE
111! ----
112ps_avg_glob_ini = 0._dp
113ps_avg_glob_old = 0._dp
114ps_avg_glob = 0._dp
115h2oice_sublim_coverage_ini = 0._dp
116co2ice_sublim_coverage_ini = 0._dp
117
118END SUBROUTINE ini_planet
119!=======================================================================
120
121!=======================================================================
122SUBROUTINE allocate_xios_state()
123!-----------------------------------------------------------------------
124! NAME
125!     allocate_xios_state
126!
127! DESCRIPTION
128!     Allocate arrays needed before loading PCM/XIOS data.
129!
130! AUTHORS & DATE
131!     JB Clement, 03/2026
132!
133! NOTES
134!
135!-----------------------------------------------------------------------
136
137! DEPENDENCIES
138! ------------
139use geometry, only: ngrid, nslope, nday, nsoil
140
141! DECLARATION
142! -----------
143implicit none
144
145! CODE
146! ----
147allocate(ps_avg(ngrid))
148allocate(ps_ts(ngrid,nday))
149allocate(tsurf_avg(ngrid,nslope))
150allocate(h2o_surfdensity_avg(ngrid,nslope))
151allocate(tsoil_avg(ngrid,nsoil,nslope))
152allocate(tsoil_ts(ngrid,nsoil,nslope,nday))
153allocate(h2o_soildensity_avg(ngrid,nsoil,nslope))
154allocate(q_h2o_ts(ngrid,nday))
155allocate(q_co2_ts(ngrid,nday))
156allocate(flux_ssice_avg(ngrid,nslope))
157
158ps_avg(:) = 0._dp
159ps_ts(:,:) = 0._dp
160tsurf_avg(:,:) = 0._dp
161h2o_surfdensity_avg(:,:) = 0._dp
162tsoil_avg(:,:,:) = 0._dp
163tsoil_ts(:,:,:,:) = 0._dp
164h2o_soildensity_avg(:,:,:) = 0._dp
165q_h2o_ts(:,:) = 0._dp
166q_co2_ts(:,:) = 0._dp
167flux_ssice_avg(:,:) = 0._dp
168
169END SUBROUTINE allocate_xios_state
170!=======================================================================
171
172!=======================================================================
173SUBROUTINE allocate_deviation_state()
174!-----------------------------------------------------------------------
175! NAME
176!     allocate_deviation_state
177!
178! DESCRIPTION
179!     Allocate persistent deviation fields derived from PCM averages.
180!
181! AUTHORS & DATE
182!     JB Clement, 03/2026
183!
184! NOTES
185!
186!-----------------------------------------------------------------------
187
188! DEPENDENCIES
189! ------------
190use geometry, only: ngrid, nslope, nsoil_PCM
191
192! DECLARATION
193! -----------
194implicit none
195
196! CODE
197! ----
198allocate(ps_dev(ngrid))
199allocate(tsurf_dev(ngrid,nslope))
200allocate(tsoil_dev(ngrid,nsoil_PCM,nslope))
201
202ps_dev(:) = 0._dp
203tsurf_dev(:,:) = 0._dp
204tsoil_dev(:,:,:) = 0._dp
205
206END SUBROUTINE allocate_deviation_state
207!=======================================================================
208
209!=======================================================================
210SUBROUTINE allocate_startevo_state()
211!-----------------------------------------------------------------------
212! NAME
213!     allocate_startevo_state
214!
215! DESCRIPTION
216!     Allocate arrays loaded from startevo or evolving afterward.
217!
218! AUTHORS & DATE
219!     JB Clement, 03/2026
220!
221! NOTES
222!
223!-----------------------------------------------------------------------
224
225! DEPENDENCIES
226! ------------
227use geometry, only: ngrid, nslope, nsoil
228
229! DECLARATION
230! -----------
231implicit none
232
233! LOCAL VARIABLES
234! ---------------
235integer :: i, islope
236
237! CODE
238! ----
239allocate(h2o_ice(ngrid,nslope))
240allocate(co2_ice(ngrid,nslope))
241allocate(icetable_depth(ngrid,nslope))
242allocate(icetable_thickness(ngrid,nslope))
243allocate(ice_porefilling(ngrid,nsoil,nslope))
244allocate(h2o_ads_reg(ngrid,nsoil,nslope))
245allocate(co2_ads_reg(ngrid,nsoil,nslope))
246allocate(delta_h2o_ads(ngrid))
247allocate(delta_co2_ads(ngrid))
248allocate(layerings_map(ngrid,nslope))
249
250h2o_ice(:,:) = 0._dp
251co2_ice(:,:) = 0._dp
252icetable_depth(:,:) = 0._dp
253icetable_thickness(:,:) = 0._dp
254ice_porefilling(:,:,:) = 0._dp
255h2o_ads_reg(:,:,:) = 0._dp
256co2_ads_reg(:,:,:) = 0._dp
257delta_h2o_ads(:) = 0._dp
258delta_co2_ads(:) = 0._dp
259do islope = 1,nslope
260    do i = 1,ngrid
261        layerings_map(i,islope)%nb_str = 0
262        nullify(layerings_map(i,islope)%bottom,layerings_map(i,islope)%top)
263    end do
264end do
265
266END SUBROUTINE allocate_startevo_state
267!=======================================================================
268
269!=======================================================================
270SUBROUTINE allocate_tendencies()
271!-----------------------------------------------------------------------
272! NAME
273!     allocate_tendencies
274!
275! DESCRIPTION
276!     Allocate perennial ice tendency arrays.
277!
278! AUTHORS & DATE
279!     JB Clement, 03/2026
280!
281! NOTES
282!
283!-----------------------------------------------------------------------
284
285! DEPENDENCIES
286! ------------
287use geometry, only: ngrid, nslope
288
289! DECLARATION
290! -----------
291implicit none
292
293! CODE
294! ----
295allocate(d_h2oice(ngrid,nslope))
296allocate(d_co2ice(ngrid,nslope))
297
298d_h2oice(:,:) = 0._dp
299d_co2ice(:,:) = 0._dp
300
301END SUBROUTINE allocate_tendencies
302!=======================================================================
303
304!=======================================================================
305SUBROUTINE allocate_initial_snapshots()
306!-----------------------------------------------------------------------
307! NAME
308!     allocate_initial_snapshots
309!
310! DESCRIPTION
311!     Allocate arrays storing initial-state snapshots.
312!
313! AUTHORS & DATE
314!     JB Clement, 03/2026
315!
316! NOTES
317!
318!-----------------------------------------------------------------------
319
320! DEPENDENCIES
321! ------------
322use geometry, only: ngrid, nslope, nday
323
324! DECLARATION
325! -----------
326implicit none
327
328! CODE
329! ----
330allocate(d_co2ice_ini(ngrid,nslope))
331allocate(q_co2_ts_ini(ngrid,nday))
332allocate(is_h2oice_ini(ngrid,nslope))
333allocate(is_co2ice_ini(ngrid,nslope))
334
335d_co2ice_ini(:,:) = 0._dp
336q_co2_ts_ini(:,:) = 0._dp
337is_h2oice_ini(:,:) = .false.
338is_co2ice_ini(:,:) = .false.
339
340END SUBROUTINE allocate_initial_snapshots
341!=======================================================================
342
343!=======================================================================
344SUBROUTINE allocate_loop_state()
345!-----------------------------------------------------------------------
346! NAME
347!     allocate_loop_state
348!
349! DESCRIPTION
350!     Allocate arrays that start being used during the main PEM loop.
351!
352! AUTHORS & DATE
353!     JB Clement, 03/2026
354!
355! NOTES
356!
357!-----------------------------------------------------------------------
358
359! DEPENDENCIES
360! ------------
361use geometry, only: ngrid, nslope, nsoil, nday
362
363! DECLARATION
364! -----------
365implicit none
366
367! CODE
368! ----
369allocate(delta_icetable(ngrid))
370allocate(icetable_depth_old(ngrid,nslope))
371allocate(is_co2ice_disappeared(ngrid,nslope))
372allocate(tsoil_ts_old(ngrid,nsoil,nslope,nday))
373
374delta_icetable(:) = 0._dp
375icetable_depth_old(:,:) = 0._dp
376is_co2ice_disappeared(:,:) = .false.
377tsoil_ts_old(:,:,:,:) = 0._dp
378
379END SUBROUTINE allocate_loop_state
380!=======================================================================
381
382!=======================================================================
383SUBROUTINE end_planet()
384!-----------------------------------------------------------------------
385! NAME
386!     end_planet
387!
388! DESCRIPTION
389!     Finalize module planet.
390!
391! AUTHORS & DATE
392!     JB Clement, 03/2026
393!
394! NOTES
395!     Deallocate all allocatable state arrays of module planet.
396!-----------------------------------------------------------------------
397
398! DEPENDENCIES
399! ------------
400use geometry, only: ngrid, nslope
401
402! DECLARATION
403! -----------
404implicit none
405
406! LOCAL VARIABLES
407! ---------------
408integer :: i, islope
409
410! CODE
411! ----
412if (allocated(ps_avg)) deallocate(ps_avg)
413if (allocated(ps_ts)) deallocate(ps_ts)
414if (allocated(ps_dev)) deallocate(ps_dev)
415if (allocated(h2o_ice)) deallocate(h2o_ice)
416if (allocated(co2_ice)) deallocate(co2_ice)
417if (allocated(is_h2oice_ini)) deallocate(is_h2oice_ini)
418if (allocated(is_co2ice_ini)) deallocate(is_co2ice_ini)
419if (allocated(is_co2ice_disappeared)) deallocate(is_co2ice_disappeared)
420if (allocated(tsurf_avg)) deallocate(tsurf_avg)
421if (allocated(tsurf_dev)) deallocate(tsurf_dev)
422if (allocated(h2o_surfdensity_avg)) deallocate(h2o_surfdensity_avg)
423if (allocated(tsoil_avg)) deallocate(tsoil_avg)
424if (allocated(tsoil_dev)) deallocate(tsoil_dev)
425if (allocated(tsoil_ts)) deallocate(tsoil_ts)
426if (allocated(tsoil_ts_old)) deallocate(tsoil_ts_old)
427if (allocated(layerings_map)) then
428    do islope = 1,nslope
429        do i = 1,ngrid
430            call del_layering(layerings_map(i,islope))
431        end do
432    end do
433    deallocate(layerings_map)
434end if
435if (allocated(h2o_soildensity_avg)) deallocate(h2o_soildensity_avg)
436if (allocated(delta_co2_ads)) deallocate(delta_co2_ads)
437if (allocated(delta_h2o_ads)) deallocate(delta_h2o_ads)
438if (allocated(h2o_ads_reg)) deallocate(h2o_ads_reg)
439if (allocated(co2_ads_reg)) deallocate(co2_ads_reg)
440if (allocated(icetable_depth)) deallocate(icetable_depth)
441if (allocated(icetable_thickness)) deallocate(icetable_thickness)
442if (allocated(ice_porefilling)) deallocate(ice_porefilling)
443if (allocated(icetable_depth_old)) deallocate(icetable_depth_old)
444if (allocated(delta_icetable)) deallocate(delta_icetable)
445if (allocated(q_co2_ts)) deallocate(q_co2_ts)
446if (allocated(q_co2_ts_ini)) deallocate(q_co2_ts_ini)
447if (allocated(q_h2o_ts)) deallocate(q_h2o_ts)
448if (allocated(d_co2ice)) deallocate(d_co2ice)
449if (allocated(d_co2ice_ini)) deallocate(d_co2ice_ini)
450if (allocated(d_h2oice)) deallocate(d_h2oice)
451if (allocated(flux_ssice_avg)) deallocate(flux_ssice_avg)
452
453END SUBROUTINE end_planet
454!=======================================================================
455
456END MODULE planet
Note: See TracBrowser for help on using the repository browser.