Changeset 3036 for trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F

Timestamp:

Sep 8, 2023, 6:47:24 PM (16 months ago)

Author:

jbclement

Message:

Mars PCM 1D:
In testphys1D.F, "dtphys" was renamed into a local variable "dttestphys" because "dtphys" is also a variable of time_phylmdz_mod.F90. Thus, it could have been modified by the physics through the module (for ex. by tabfi.F when reading startfi.nc) and affect the dynamics which is prohibited.
JBC

File:

• trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F (modified) (12 diffs)

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

@@ -23 +23 @@
      &                     obliquit, emin_turb, lmixmin
       use comcstfi_h, only: pi, rad, omeg, g, mugaz, rcp, r, cpp
-      use time_phylmdz_mod, only: daysec, dtphys, day_step,
-     &                            ecritphy, iphysiq
+      use time_phylmdz_mod, only: daysec, day_step, ecritphy, iphysiq
       use dimradmars_mod, only: tauvis,totcloudfrac
       use dust_param_mod, only: tauscaling
@@ -95 +94 @@
       real,parameter :: odpref=610. ! DOD reference pressure (Pa)
 c
-      INTEGER day0,dayn          ! date initial (sol ; =0 a Ls=0) and final
-      REAL day           ! date durant le run
-      REAL time             ! time (0<time<1 ; time=0.5 a midi)
+      INTEGER day0,dayn   ! initial (sol ; =0 at Ls=0) and final date
+      REAL day            ! date during the run
+      REAL time           ! time (0<time<1 ; time=0.5 a midi)
+      REAL dttestphys         ! testphys1d timestep
       REAL play(nlayer)   ! Pressure at the middle of the layers (Pa)
       REAL plev(nlayer+1) ! intermediate pressure levels (pa)
       REAL psurf,tsurf(1)      
-      REAL u(nlayer),v(nlayer)  ! zonal, meridional wind
-      REAL gru,grv   ! prescribed "geostrophic" background wind
-      REAL temp(nlayer)   ! temperature at the middle of the layers
-      REAL,ALLOCATABLE :: q(:,:) ! tracer mixing ratio (e.g. kg/kg)
+      REAL u(nlayer),v(nlayer) ! zonal, meridional wind
+      REAL gru,grv             ! prescribed "geostrophic" background wind
+      REAL temp(nlayer)        ! temperature at the middle of the layers
+      REAL,ALLOCATABLE :: q(:,:)   ! tracer mixing ratio (e.g. kg/kg)
       REAL,ALLOCATABLE :: qsurf(:) ! tracer surface budget (e.g. kg.m-2)
       REAL tsoil(nsoilmx)   ! subsurface soik temperature (K)
       REAL emis(1)          ! surface layer
       REAL albedo(1,1)      ! surface albedo
-      REAL :: wstar(1)=0.    ! Thermals vertical velocity
-      REAL q2(nlayer+1)   ! Turbulent Kinetic Energy
-      REAL zlay(nlayer)   ! altitude estimee dans les couches (km)
+      REAL :: wstar(1)=0.   ! Thermals vertical velocity
+      REAL q2(nlayer+1)     ! Turbulent Kinetic Energy
+      REAL zlay(nlayer)     ! altitude estimee dans les couches (km)
 
 c    Physical and dynamical tandencies (e.g.  m.s-2, K/s, Pa/s)
@@ -416 +416 @@
       dayn=day0+ndt
       ndt=ndt*day_step     
-      dtphys=daysec/day_step  
+      dttestphys=daysec/day_step
 
 c Imposed surface pressure
@@ -498 +498 @@
 !  some initializations (some of which have already been
 !  done above!) and loads parameters set in callphys.def
-!  and allocates some arrays 
-!Mars possible matter with dtphys in input and include!!!
+!  and allocates some arrays
+! Mars possible matter with dttestphys in input and include!!!
 ! Initializations below should mimick what is done in iniphysiq for 3D GCM
       call init_interface_dyn_phys
@@ -513 +513 @@
       call phys_state_var_init(1,llm,nq,tname,
      .          day0,dayn,time,
-     .          daysec,dtphys,
+     .          daysec,dttestphys,
      .          rad,g,r,cpp,
      .          nqperes,nqfils)! MVals: variables isotopes
@@ -838 +838 @@
       atm_wat_tau = -1. ! Default: no time relaxation
       if (water) then
-      print*, 'Relax atmospheric water vapor profile?'
+      write(*,*) 'Relax atmospheric water vapor profile?'
       call getin("atm_wat_tau",atm_wat_tau)
       write(*,*) "atm_wat_tau = ", atm_wat_tau
@@ -864 +864 @@
 
       call physdem0("startfi.nc",longitude,latitude,nsoilmx,ngrid,
-     &              llm,nq,dtphys,float(day0),0.,cell_area,
+     &              llm,nq,dttestphys,float(day0),0.,cell_area,
      &              albedodat,inertiedat,def_slope,subslope_dist)
       call physdem1("startfi.nc",nsoilmx,ngrid,llm,nq,
-     &              dtphys,time,
+     &              dttestphys,time,
      &              tsurf,tsoil,inertiesoil,albedo,emis,
      &              q2,qsurf,tauscaling,
@@ -922 +922 @@
           else ! Relaxation towards the value atm_wat_profile with relaxation time atm_wat_tau
         q(:,igcm_h2o_vap) = atm_wat_profile*g/psurf + (q(:,igcm_h2o_vap)
-     &              - atm_wat_profile*g/psurf)*dexp(-dtphys/atm_wat_tau)
+     &          - atm_wat_profile*g/psurf)*dexp(-dttestphys/atm_wat_tau)
             q(:,igcm_h2o_vap) = min(zqsat(:),q(:,igcm_h2o_vap))
             q(:,igcm_h2o_ice) = 0. ! reset h2o ice
@@ -932 +932 @@
 c       call physics
 c       --------------------
-      CALL physiq (1,llm,nq,
-     ,     firstcall,lastcall,
-     ,     day,time,dtphys,
-     ,     plev,play,phi,
-     ,     u, v,temp, q,  
-     ,     w,
+      CALL physiq (1,llm,nq,firstcall,lastcall,
+     ,     day,time,dttestphys,plev,play,phi,
+     ,     u,v,temp,q,w,
 C - outputs
      s     du, dv, dtemp, dq,dpsurf)
@@ -963 +960 @@
 c       ---------------------------------------
         firstcall=.false.
-        time=time+dtphys/daysec
+        time=time+dttestphys/daysec
         IF (time.gt.1.E+0) then
             time=time-1.E+0
@@ -973 +970 @@
 
         DO ilayer=1,nlayer
-           u(ilayer)=u(ilayer)+dtphys*du(ilayer)
-           v(ilayer)=v(ilayer)+dtphys*dv(ilayer)
-           temp(ilayer)=temp(ilayer)+dtphys*dtemp(ilayer)
+           u(ilayer)=u(ilayer)+dttestphys*du(ilayer)
+           v(ilayer)=v(ilayer)+dttestphys*dv(ilayer)
+           temp(ilayer)=temp(ilayer)+dttestphys*dtemp(ilayer)
         ENDDO
 
 c       compute pressure for next time step 
 c       ----------------------------------------------------------
-           psurf=psurf+dtphys*dpsurf(1)   ! surface pressure change
+           psurf=psurf+dttestphys*dpsurf(1)   ! surface pressure change
            DO ilevel=1,nlevel
              plev(ilevel)=ap(ilevel)+psurf*bp(ilevel)
@@ -991 +988 @@
         DO iq = 1, nq
           DO ilayer=1,nlayer
-             q(ilayer,iq)=q(ilayer,iq)+dtphys*dq(ilayer,iq)
+             q(ilayer,iq)=q(ilayer,iq)+dttestphys*dq(ilayer,iq)
           ENDDO
         ENDDO