Changeset 3584 for trunk/LMDZ.COMMON/libf

Timestamp:

Jan 17, 2025, 5:14:18 PM (6 days ago)

Author:

jbclement

Message:

PEM:

• Albedo is now updated only at the end of the PEM (and not at every iteration) + Correct way to set it taking into account CO2/H2O ice and frost.
• Cosmetic cleanings.

JBC

Location:

trunk/LMDZ.COMMON/libf/evolution

Files:

• changelog.txt (modified) (1 diff)
• deftank/PEMrun.job (modified) (1 diff)
• deftank/modify_startfi_orbit.sh (modified) (1 diff)
• pem.F90 (modified) (5 diffs)
• recomp_tend_co2_mod.F90 (modified) (2 diffs)
• reshape_XIOS_output.F90 (modified) (7 diffs)
• soil_settings_PEM_mod.F90 (modified) (1 diff)

trunk/LMDZ.COMMON/libf/evolution/changelog.txt

@@ -540 +540 @@
 == 15/01/2025 == JBC
 Improvement of the Bash script tools in the deftank with an automatic error detection which ends the script with a message.
+
+== 17/01/2025 == JBC
+- Albedo is now updated only at the end of the PEM (and not at every iteration) + Correct way to set it taking into account CO2/H2O ice and frost.
+- Cosmetic cleanings.

trunk/LMDZ.COMMON/libf/evolution/deftank/PEMrun.job

@@ -16 +16 @@
 # Modify here the parameters depending on your setup
 ####################################################
-# A few parameters that might need to be changed depending on your setup
 # Path to the arch.env to source:
 source ../trunk/LMDZ.COMMON/arch.env

trunk/LMDZ.COMMON/libf/evolution/deftank/modify_startfi_orbit.sh

@@ -10 +10 @@
 ###########################################
 # Name of the file
-name_file = "startfi.nc"
+name_file="startfi.nc"
 
 # New values for the orbital parameters

trunk/LMDZ.COMMON/libf/evolution/pem.F90

@@ -75 +75 @@
 #ifndef CPP_STD
     use comsoil_h,          only: tsoil, nsoilmx, ini_comsoil_h, inertiedat, mlayer, inertiesoil, flux_geo, nqsoil, qsoil
-    use surfdat_h,          only: tsurf, emis, qsurf, watercap, ini_surfdat_h, &
-                                  albedodat, zmea, zstd, zsig, zgam, zthe,     &
-                                  albedo_h2o_frost,frost_albedo_threshold,     &
-                                  emissiv, watercaptag, perennial_co2ice
+    use surfdat_h,          only: tsurf, qsurf, emis, emissiv, emisice, ini_surfdat_h,   &
+                                  albedodat, albedice, albedo_h2o_frost, albedo_h2o_cap, &
+                                  zmea, zstd, zsig, zgam, zthe, frost_albedo_threshold,  &
+                                  watercap, watercaptag, perennial_co2ice, albedo_perennialco2
     use dimradmars_mod,     only: totcloudfrac, albedo
     use dust_param_mod,     only: tauscaling
@@ -287 +287 @@
 
 ! Loop variables
-integer :: i, l, ig, nnq, t, islope, ig_loop, islope_loop, isoil
+integer :: i, l, ig, nnq, t, islope, ig_loop, islope_loop, isoil, icap
 
 ! Elapsed time with system clock
@@ -384 +384 @@
     endif
 
-    call init_testphys1d('start1D.txt','startfi.nc',therestart1D,therestartfi,ngrid,nlayer,610.,nq,q,            &
-                         time_0,ps_start_dyn(1),ucov,vcov,teta,ndt,ptif,pks,dtphys,zqsat,dq,dqdyn,day0,day,gru,grv,w,  &
+    call init_testphys1d('start1D.txt','startfi.nc',therestart1D,therestartfi,ngrid,nlayer,610.,nq,q,                 &
+                         time_0,ps_start_dyn(1),ucov,vcov,teta,ndt,ptif,pks,dtphys,zqsat,dq,dqdyn,day0,day,gru,grv,w, &
                          play,plev,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
     ps_start_dyn(2) = ps_start_dyn(1)
@@ -844 +844 @@
                         enddo
                     enddo outer2
-                endif
-                if (h2o_ice(ig,islope) > frost_albedo_threshold) then
-                    albedo(ig,1,islope) = albedo_h2o_frost
-                    albedo(ig,2,islope) = albedo_h2o_frost
-                else
-                    albedo(ig,1,islope) = albedodat(ig)
-                    albedo(ig,2,islope) = albedodat(ig)
-                    emis(ig,islope) = emissiv
                 endif
             else if (co2_ice(ig,islope) > 1.e-10 .and. d_co2ice(ig,islope) > 1.e-10) then ! Put tsurf as tcond CO2
@@ -1115 +1107 @@
 deallocate(zplev_start0,zplev_new)
 
+! III_a.6 Albedo update for start file
+do ig = 1,ngrid
+    if (latitude(ig) < 0.) then
+        icap = 2 ! Southern hemisphere
+    else
+        icap = 1 ! Northern hemisphere
+    endif
+    do islope = 1,ngrid
+        ! Bare ground
+        albedo(ig,:,islope) = albedodat(ig)
+        emis(ig,islope) = emissiv
+
+        ! CO2 ice/frost is treated after H20 ice/frost because it is considered dominant
+        ! H2O ice
+        if (h2o_ice(ig,islope) > 0.) then
+            albedo(ig,:,islope) = albedo_h2o_cap
+            emis(ig,islope) = 1.
+        endif
+        ! CO2 ice
+        if (co2_ice(ig,islope) > 0.) then
+            albedo(ig,:,islope) = albedo_perennialco2(icap)
+            emis(ig,islope) = emisice(icap)
+        endif
+        ! H2O frost
+        if (qsurf(ig,igcm_h2o_ice,islope) > 0.) then
+            albedo(ig,:,islope) = albedo_h2o_frost
+            emis(ig,islope) = 1.
+        endif
+        ! CO2 frost
+        if (qsurf(ig,igcm_co2,islope) > 0.) then
+            albedo(ig,:,islope) = albedice(icap)
+            emis(ig,islope) = emisice(icap)
+        endif
+    enddo
+enddo
+
+! III_a.7 Orbital parameters update for start file
 if (evol_orbit_pem) call recomp_orb_param(i_myear,i_myear_leg)
 

trunk/LMDZ.COMMON/libf/evolution/recomp_tend_co2_mod.F90

@@ -16 +16 @@
 !=======================================================================
 !
-!  Routine that compute the evolution of the tendencie for co2 ice
+!  To compute the evolution of the tendencie for co2 ice
 !
 !=======================================================================
@@ -67 +67 @@
 !=======================================================================
 !
-!  Routine that compute the evolution of the tendencie for h2o ice
+!  To compute the evolution of the tendencie for h2o ice
 !
 !=======================================================================

trunk/LMDZ.COMMON/libf/evolution/reshape_XIOS_output.F90

@@ -42 +42 @@
 
 do numyear = 1,2
-    write(*,*) 'numyear',numyear
-    write(str(1:1),'(i1.1)') numyear
-
-    state = nf90_open(path = "data2reshape_Y"//str//".nc", mode = nf90_nowrite, ncid = ncid1)
+    write(str,'(i1.1)') numyear
+    write(*,*) 'Reshaping of variables from "data2reshape_Y'//str//'.nc"...'
+
+    state = nf90_open(path = "data2reshape_Y"//str//".nc",mode = nf90_nowrite,ncid = ncid1)
     if (state /= nf90_noerr) call handle_err(state)
 
@@ -57 +57 @@
         endif
     endif
-    state = nf90_create(path = "data_PCM_Y"//str//".nc", cmode=or(nf90_noclobber,nf90_64bit_offset), ncid = ncid2)
+    state = nf90_create(path = "data_PCM_Y"//str//".nc",cmode = or(nf90_noclobber,nf90_64bit_offset),ncid = ncid2)
     if (state /= nf90_noerr) call handle_err(state)
 
@@ -69 +69 @@
     allocate(varids_2(nvars))
 
-    state = nf90_inq_dimids(ncid1, ndims, dimids, include_parents)
-    if (state /= nf90_noerr) call handle_err(state)
-    state = nf90_inq_varids(ncid1, nvars, varids)
+    state = nf90_inq_dimids(ncid1,ndims,dimids,include_parents)
+    if (state /= nf90_noerr) call handle_err(state)
+    state = nf90_inq_varids(ncid1,nvars,varids)
     if (state /= nf90_noerr) call handle_err(state)
 
     do i = 1,ndims
-        state = nf90_inquire_dimension(ncid1, dimids(i), name_, len_)
+        state = nf90_inquire_dimension(ncid1,dimids(i),name_,len_)
         if (state /= nf90_noerr) call handle_err(state)
         if (name_ == "lon" .or. name_ == "longitude") then
@@ -82 +82 @@
             len_ = len_ + 1
         else if (name_ == "lat".or. name_ == "latitude") then
-            dimid_lat=dimids(i)
-            len_lat=len_
+            dimid_lat = dimids(i)
+            len_lat = len_
         else if (name_ == "time_counter".or. name_ ==  "Time") then
-            dimid_time=dimids(i)
-            len_time=len_
+            dimid_time = dimids(i)
+            len_time = len_
         else if (name_ == "soil_layers".or. name_ ==  "subsurface_layers") then
-            dimid_soil=dimids(i)
+            dimid_soil = dimids(i)
             len_soil = len_
         endif
-        state = nf90_def_dim(ncid2, name_,len_,dimid_2)
+        state = nf90_def_dim(ncid2,name_,len_,dimid_2)
         if (state /= nf90_noerr) call handle_err(state)
         dimids_2(i) = dimid_2
@@ -97 +97 @@
 
     do i = 1,nvars
-        state = nf90_inquire_variable(ncid1, varids(i),name = namevar,xtype = xtype_var,ndims = numdims,natts = numatts)
-        write(*,*) "namevar00= ", namevar
+        state = nf90_inquire_variable(ncid1,varids(i),name = namevar,xtype = xtype_var,ndims = numdims,natts = numatts)
+        write(*,*) '> Treatment of '//namevar
         if (state /= nf90_noerr) call handle_err(state)
         allocate(dimid_var(numdims))
@@ -109 +109 @@
                 state = nf90_get_var(ncid1,varids(i),tempvalues_1d)
                 if (state /= nf90_noerr) call handle_err(state)
-                state = nf90_def_var(ncid2,namevar,xtype_var, dimid_var, varids_2(i))
+                state = nf90_def_var(ncid2,namevar,xtype_var,dimid_var,varids_2(i))
                 if (state /= nf90_noerr) call handle_err(state)
                 values_1d(1:len_lon) = tempvalues_1d(:)
@@ -219 +219 @@
 
     deallocate(dimids,varids,dimids_2,varids_2)
+    write(*,*) 'Done!'
 enddo
 

trunk/LMDZ.COMMON/libf/evolution/soil_settings_PEM_mod.F90

@@ -20 +20 @@
 !
 ! Modifications: Aug.2010 EM: use NetCDF90 to load variables (enables using
-!                              r4 or r8 restarts independently of having compiled
-!                              the GCM in r4 or r8)
+!                             r4 or r8 restarts independently of having compiled
+!                             the GCM in r4 or r8)
 !                June 2013 TN: Possibility to read files with a time axis
 !