Changeset 3060

Timestamp:

Sep 28, 2023, 3:15:23 PM (14 months ago)

Author:

jbclement

Message:

Mars PCM 1D:
Related to commit r3056, correction of bugs and some adaptations of the subroutine 'init_testphys1d'.
JBC

Location:

trunk/LMDZ.MARS

Files:

• changelog.txt (modified) (1 diff)
• libf/phymars/dyn1d/init_testphys1d_mod.F90 (modified) (14 diffs)
• libf/phymars/dyn1d/testphys1d.F90 (modified) (6 diffs)

trunk/LMDZ.MARS/changelog.txt

@@ -4232 +4232 @@
 - update callphys.def.MCD6 to disable Frost metamorphism option
 
-
+== 28/09/2023 == JBC
+Related to commit r3056, correction of bugs and some adaptations of the subroutine 'init_testphys1d'.

trunk/LMDZ.MARS/libf/phymars/dyn1d/init_testphys1d_mod.F90

@@ -5 +5 @@
 contains
 
-SUBROUTINE init_testphys1d(ngrid,nq,nlayer,odpref,ndt,ptif,pks,dttestphys,startfiles_1D,q,zqsat,qsurf,dq,dqdyn, &
-                           day0,day,time,psurf,tsurf,gru,grv,u,v,w,q2,play,plev,tsoil,temp,albedo,emis,         &
-                           latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
+SUBROUTINE init_testphys1d(startfiles_1D,therestart1D,therestartfi,ngrid,nlayer,odpref,nq,ndt,ptif,pks,dttestphys, &
+                           q,zqsat,qsurf,dq,dqdyn,day0,day,time,psurf,tsurf,gru,grv,u,v,w,q2,play,plev,tsoil,temp, &
+                           albedo,emis,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
 
 use ioipsl_getincom,          only: getin ! To use 'getin'
@@ -47 +47 @@
 ! Arguments
 !=======================================================================
+logical, intent(in) :: startfiles_1D, therestart1D, therestartfi ! Use of "start1D.txt" and "startfi.nc" files
 integer, intent(in) :: ngrid, nlayer
 real,    intent(in) :: odpref ! DOD reference pressure (Pa)
@@ -55 +56 @@
 real,                              intent(out) :: ptif, pks
 real,                              intent(out) :: dttestphys    ! testphys1d timestep
-logical,                           intent(out) :: startfiles_1D ! Use of "start1D.txt" and "startfi.nc" files
 real, dimension(:,:), allocatable, intent(out) :: q             ! tracer mixing ratio (e.g. kg/kg)
 real, dimension(:),   allocatable, intent(out) :: zqsat         ! useful for (atm_wat_profile=2)
@@ -62 +62 @@
 integer,                           intent(out) :: day0          ! initial (sol ; =0 at Ls=0) and final date
 real,                              intent(out) :: day           ! date during the run
-real,                              intent(out) :: time          ! time (0<time<1 ; time=0.5 a midi)
+real,                              intent(out) :: time          ! time (0<time<1; time=0.5 at noon)
 real,                              intent(out) :: psurf         ! Surface pressure
 real, dimension(1),                intent(out) :: tsurf         ! Surface temperature
@@ -70 +70 @@
 real, dimension(nlayer),           intent(out) :: play          ! Pressure at the middle of the layers (Pa)
 real, dimension(nlayer + 1),       intent(out) :: plev          ! intermediate pressure levels (pa)
-real, dimension(nsoilmx),          intent(out) :: tsoil         ! subsurface soik temperature (K)
+real, dimension(nsoilmx),          intent(out) :: tsoil         ! subsurface soil temperature (K)
 real, dimension(nlayer),           intent(out) :: temp          ! temperature at the middle of the layers
 real, dimension(1,1),              intent(out) :: albedo        ! surface albedo
@@ -90 +90 @@
 
 ! RV & JBC: Use of "start1D.txt" and "startfi.nc" files
-logical                :: found, therestart1D, therestartfi
+logical                :: found
 character(len = 30)    :: header
 real, dimension(100)   :: tab_cntrl
@@ -118 +118 @@
 
 !------------------------------------------------------
-! Loading run parameters from "run.def" file
-!------------------------------------------------------
-! check if 'run.def' file is around (otherwise reading parameters
-! from callphys.def via getin() routine won't work.
-inquire(file = 'run.def',exist = there)
-if (.not. there) then
-    write(*,*) 'Cannot find required file "run.def"'
-    write(*,*) '  (which should contain some input parameters along with the following line: INCLUDEDEF=callphys.def)'
-    write(*,*) ' ... might as well stop here ...'
-    stop
-endif
-
-write(*,*)'Do you want to use "start1D.txt" and "startfi.nc" files?'
-startfiles_1D = .false.
-call getin("startfiles_1D",startfiles_1D)
-write(*,*) " startfiles_1D = ", startfiles_1D
-
-if (startfiles_1D) then
-    inquire(file = 'start1D.txt',exist = therestart1D)
-    if (.not. therestart1D) then
-        write(*,*) 'There is no "start1D.txt" file!'
-        write(*,*) 'Initialization is done with default values.'
-    endif
-    inquire(file = 'startfi.nc',exist = therestartfi)
-    if (.not. therestartfi) then
-        write(*,*) 'There is no "startfi.nc" file!'
-        write(*,*) 'Initialization is done with default values.'
-    endif
-endif
-
-!------------------------------------------------------
 ! Prescribed constants to be set here
 !------------------------------------------------------
@@ -155 +124 @@
 ! Mars planetary constants
 ! ------------------------
-rad = 3397200.                   ! mars radius (m) ~3397200 m
-omeg = 4.*asin(1.)/(daysec)      ! rotation rate (rad.s-1)
-g = 3.72                         ! gravity (m.s-2) ~3.72
-mugaz = 43.49                    ! atmosphere mola mass (g.mol-1) ~43.49
-rcp = .256793                    ! = r/cp ~0.256793
+rad = 3397200.            ! mars radius (m) ~3397200 m
+daysec = 88775.           ! length of a sol (s) ~88775 s
+omeg = 4.*asin(1.)/daysec ! rotation rate (rad.s-1)
+g = 3.72                  ! gravity (m.s-2) ~3.72
+mugaz = 43.49             ! atmosphere mola mass (g.mol-1) ~43.49
+rcp = .256793             ! = r/cp ~0.256793
 r = 8.314511*1000./mugaz
 cpp = r/rcp
-daysec = 88775.                  ! length of a sol (s) ~88775 s
-year_day = 669                   ! length of year (sols) ~668.6
-periheli = 206.66                ! minimum sun-mars distance (Mkm) ~206.66
-aphelie = 249.22                 ! maximum sun-mars distance (Mkm) ~249.22
-halfaxe = 227.94                 ! demi-grand axe de l'ellipse
-peri_day = 485.                  ! perihelion date (sols since N. Spring)
-obliquit = 25.2                  ! Obliquity (deg) ~25.2
-eccentric = 0.0934               ! Eccentricity (0.0934)
+year_day = 669            ! length of year (sols) ~668.6
+periheli = 206.66         ! minimum sun-mars distance (Mkm) ~206.66
+aphelie = 249.22          ! maximum sun-mars distance (Mkm) ~249.22
+halfaxe = 227.94          ! demi-grand axe de l'ellipse
+peri_day = 485.           ! perihelion date (sols since N. Spring)
+obliquit = 25.2           ! Obliquity (deg) ~25.2
+eccentric = 0.0934        ! Eccentricity (0.0934)
 
 ! Planetary Boundary Layer and Turbulence parameters
@@ -202 +171 @@
     stop
 else
-    write(*,*) "testphys1d: Reading file traceur.def"
+    write(*,*) "init_testphys1d: Reading file traceur.def"
     ! read number of tracers:
     read(90,*,iostat = ierr) nq
     nqtot = nq ! set value of nqtot (in infotrac module) as nq
     if (ierr /= 0) then
-        write(*,*) "testphys1d: error reading number of tracers"
+        write(*,*) "init_testphys1d: error reading number of tracers"
         write(*,*) "   (first line of traceur.def) "
         stop
     endif
     if (nq < 1) then
-        write(*,*) "testphys1d: error number of tracers"
+        write(*,*) "init_testphys1d: error number of tracers"
         write(*,*) "is nq=",nq," but must be >=1!"
         stop
@@ -225 +194 @@
     read(90,'(80a)',iostat = ierr) line ! store the line from traceur.def
     if (ierr /= 0) then
-        write(*,*) 'testphys1d: error reading tracer names...'
+        write(*,*) 'init_testphys1d: error reading tracer names...'
         stop
     endif
@@ -271 +240 @@
 
 ! Initialize tracers here:
-write(*,*) "testphys1d: initializing tracers"
+write(*,*) "init_testphys1d: initializing tracers"
 if (.not. therestart1D) then
     call read_profile(nq,nlayer,qsurf,q)
@@ -434 +403 @@
 ! ovverride iphysiq value that has been set by conf_phys
 if (iphysiq /= 1) then
-    write(*,*) "testphys1d: setting iphysiq=1"
+    write(*,*) "init_testphys1d: setting iphysiq=1"
     iphysiq = 1
 endif
@@ -538 +507 @@
 enddo
 if (igcm_co2 == 0) then
-    write(*,*) "testphys1d error, missing co2 tracer!"
+    write(*,*) "init_testphys1d error, missing co2 tracer!"
     stop
 endif
@@ -725 +694 @@
 
 END MODULE init_testphys1d_mod
-

trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F90

@@ -1 +1 @@
 PROGRAM testphys1d
 
+use ioipsl_getincom,         only: getin ! To use 'getin'
 use comsoil_h,               only: inertiedat, inertiesoil, nsoilmx
 use surfdat_h,               only: albedodat, perenial_co2ice, watercap
@@ -94 +95 @@
 integer                 :: nq = 1 ! number of tracers
 real, dimension(1)      :: latitude, longitude, cell_area
-
-character(len = 2)  :: str2
-character(len = 7)  :: str7
-character(len = 44) :: txt
+logical                 :: there
+character(len = 2)      :: str2
+character(len = 7)      :: str7
+character(len = 44)     :: txt
 
 ! RV & JBC: Use of "start1D.txt" and "startfi.nc" files
-logical :: startfiles_1D
+logical :: startfiles_1D, therestart1D, therestartfi
 
 ! JN & JBC: Force atmospheric water profiles
@@ -120 +121 @@
 !call initcomgeomphy
 
-call init_testphys1d(ngrid,nq,nlayer,odpref,ndt,ptif,pks,dttestphys,startfiles_1D,q,zqsat,qsurf,dq,dqdyn, &
-                     day0,day,time,psurf,tsurf,gru,grv,u,v,w,q2,play,plev,tsoil,temp,albedo,emis,         &
-                     latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
+!------------------------------------------------------
+! Loading run parameters from "run.def" file
+!------------------------------------------------------
+! check if 'run.def' file is around. Otherwise reading parameters
+! from callphys.def via getin() routine won't work.
+inquire(file = 'run.def',exist = there)
+if (.not. there) then
+    write(*,*) 'Cannot find required file "run.def"'
+    write(*,*) '  (which should contain some input parameters along with the following line: INCLUDEDEF=callphys.def)'
+    write(*,*) ' ... might as well stop here ...'
+    stop
+endif
+
+write(*,*)'Do you want to use "start1D.txt" and "startfi.nc" files?'
+startfiles_1D = .false.
+therestart1D = .false.
+therestartfi = .false.
+call getin("startfiles_1D",startfiles_1D)
+write(*,*) " startfiles_1D = ", startfiles_1D
+
+if (startfiles_1D) then
+    inquire(file = 'start1D.txt',exist = therestart1D)
+    if (.not. therestart1D) then
+        write(*,*) 'There is no "start1D.txt" file!'
+        write(*,*) 'Initialization is done with default values.'
+    endif
+    inquire(file = 'startfi.nc',exist = therestartfi)
+    if (.not. therestartfi) then
+        write(*,*) 'There is no "startfi.nc" file!'
+        write(*,*) 'Initialization is done with default values.'
+    endif
+endif
+
+call init_testphys1d(startfiles_1D,therestart1D,therestartfi,ngrid,nlayer,odpref,nq,ndt,ptif,pks,dttestphys, &
+                     q,zqsat,qsurf,dq,dqdyn,day0,day,time,psurf,tsurf,gru,grv,u,v,w,q2,play,plev,tsoil,temp, &
+                     albedo,emis,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
 
 ! Write a "startfi" file
@@ -128 +162 @@
 ! This file will be read during the first call to "physiq".
 ! It is needed to transfert physics variables to "physiq"...
-if (.not. startfiles_1D) then
+if (.not. therestartfi) then
     call physdem0("startfi.nc",longitude,latitude,nsoilmx,ngrid, &
                   llm,nq,dttestphys,float(day0),0.,cell_area,    &
@@ -135 +169 @@
                   tsurf,tsoil,inertiesoil,albedo,emis,q2,qsurf,tauscaling, &
                   totcloudfrac,wstar,watercap,perenial_co2ice)
-endif !(.not. startfiles_1D )
+endif !(.not. therestartfi)
 
 !=======================================================================
@@ -260 +294 @@
 !***********************************************************************
 !***********************************************************************
-