Changeset 5160 for LMDZ6/branches/Amaury_dev/libf/phydev/physiq_mod.F90

Timestamp:

Aug 3, 2024, 2:56:58 PM (7 weeks ago)

Author:

abarral

Message:

Put .h into modules

File:

• LMDZ6/branches/Amaury_dev/libf/phydev/physiq_mod.F90 (modified) (1 diff)

LMDZ6/branches/Amaury_dev/libf/phydev/physiq_mod.F90

@@ -2 +2 @@
 MODULE physiq_mod
 
-IMPLICIT NONE
+  IMPLICIT NONE
 
 CONTAINS
 
-      SUBROUTINE physiq(nlon,nlev, &
-              debut,lafin,pdtphys, &
-              paprs,pplay,pphi,pphis,presnivs, &
-              u,v,t,qx, &
-              flxmass_w, &
-              d_u, d_v, d_t, d_qx, d_ps)
+  SUBROUTINE physiq(nlon, nlev, &
+          debut, lafin, pdtphys, &
+          paprs, pplay, pphi, pphis, presnivs, &
+          u, v, t, qx, &
+          flxmass_w, &
+          d_u, d_v, d_t, d_qx, d_ps)
 
-      USE dimphy, ONLY: klon,klev
-      USE infotrac_phy, ONLY: nqtot
-      USE lmdz_geometry, ONLY: latitude
-      USE comcstphy, ONLY: rg
-      USE iophy, ONLY: histbeg_phy,histwrite_phy
-      USE ioipsl, ONLY: getin,histvert,histdef,histend,ymds2ju
-      USE lmdz_phys_para, ONLY: jj_nb
-      USE phys_state_var_mod, ONLY: phys_state_var_init
-      USE lmdz_grid_phy, ONLY: nbp_lon,nbp_lat
+    USE dimphy, ONLY: klon, klev
+    USE infotrac_phy, ONLY: nqtot
+    USE lmdz_geometry, ONLY: latitude
+    USE comcstphy, ONLY: rg
+    USE iophy, ONLY: histbeg_phy, histwrite_phy
+    USE ioipsl, ONLY: getin, histvert, histdef, histend, ymds2ju
+    USE lmdz_phys_para, ONLY: jj_nb
+    USE phys_state_var_mod, ONLY: phys_state_var_init
+    USE lmdz_grid_phy, ONLY: nbp_lon, nbp_lat
 
-      USE lmdz_xios, ONLY: xios_update_calendar, using_xios
-      USE lmdz_wxios, ONLY: wxios_add_vaxis, wxios_set_cal, wxios_closedef
-      USE iophy, ONLY: histwrite_phy
+    USE lmdz_xios, ONLY: xios_update_calendar, using_xios
+    USE lmdz_wxios, ONLY: wxios_add_vaxis, wxios_set_cal, wxios_closedef
+    USE iophy, ONLY: histwrite_phy
 
-      IMPLICIT NONE
+    IMPLICIT NONE
 
-! Routine argument:
+    ! Routine argument:
 
-      INTEGER,INTENT(IN) :: nlon ! number of atmospheric colums
-      INTEGER,INTENT(IN) :: nlev ! number of vertical levels (should be =klev)
-      logical,INTENT(IN) :: debut ! signals first CALL to physics
-      logical,INTENT(IN) :: lafin ! signals last CALL to physics
-      REAL,INTENT(IN) :: pdtphys ! physics time step (s)
-      REAL,INTENT(IN) :: paprs(klon,klev+1) ! interlayer pressure (Pa)
-      REAL,INTENT(IN) :: pplay(klon,klev) ! mid-layer pressure (Pa)
-      REAL,INTENT(IN) :: pphi(klon,klev) ! geopotential at mid-layer
-      REAL,INTENT(IN) :: pphis(klon) ! surface geopotential
-      REAL,INTENT(IN) :: presnivs(klev) ! pseudo-pressure (Pa) of mid-layers
-      REAL,INTENT(IN) :: u(klon,klev) ! eastward zonal wind (m/s)
-      REAL,INTENT(IN) :: v(klon,klev) ! northward meridional wind (m/s)
-      REAL,INTENT(IN) :: t(klon,klev) ! temperature (K)
-      REAL,INTENT(IN) :: qx(klon,klev,nqtot) ! tracers (.../kg_air)
-      REAL,INTENT(IN) :: flxmass_w(klon,klev) ! vertical mass flux
-      REAL,INTENT(OUT) :: d_u(klon,klev) ! physics tendency on u (m/s/s)
-      REAL,INTENT(OUT) :: d_v(klon,klev) ! physics tendency on v (m/s/s)
-      REAL,INTENT(OUT) :: d_t(klon,klev) ! physics tendency on t (K/s)
-      REAL,INTENT(OUT) :: d_qx(klon,klev,nqtot) ! physics tendency on tracers
-      REAL,INTENT(OUT) :: d_ps(klon) ! physics tendency on surface pressure
+    INTEGER, INTENT(IN) :: nlon ! number of atmospheric colums
+    INTEGER, INTENT(IN) :: nlev ! number of vertical levels (should be =klev)
+    LOGICAL, INTENT(IN) :: debut ! signals first CALL to physics
+    LOGICAL, INTENT(IN) :: lafin ! signals last CALL to physics
+    REAL, INTENT(IN) :: pdtphys ! physics time step (s)
+    REAL, INTENT(IN) :: paprs(klon, klev + 1) ! interlayer pressure (Pa)
+    REAL, INTENT(IN) :: pplay(klon, klev) ! mid-layer pressure (Pa)
+    REAL, INTENT(IN) :: pphi(klon, klev) ! geopotential at mid-layer
+    REAL, INTENT(IN) :: pphis(klon) ! surface geopotential
+    REAL, INTENT(IN) :: presnivs(klev) ! pseudo-pressure (Pa) of mid-layers
+    REAL, INTENT(IN) :: u(klon, klev) ! eastward zonal wind (m/s)
+    REAL, INTENT(IN) :: v(klon, klev) ! northward meridional wind (m/s)
+    REAL, INTENT(IN) :: t(klon, klev) ! temperature (K)
+    REAL, INTENT(IN) :: qx(klon, klev, nqtot) ! tracers (.../kg_air)
+    REAL, INTENT(IN) :: flxmass_w(klon, klev) ! vertical mass flux
+    REAL, INTENT(OUT) :: d_u(klon, klev) ! physics tendency on u (m/s/s)
+    REAL, INTENT(OUT) :: d_v(klon, klev) ! physics tendency on v (m/s/s)
+    REAL, INTENT(OUT) :: d_t(klon, klev) ! physics tendency on t (K/s)
+    REAL, INTENT(OUT) :: d_qx(klon, klev, nqtot) ! physics tendency on tracers
+    REAL, INTENT(OUT) :: d_ps(klon) ! physics tendency on surface pressure
 
-INTEGER,save :: itau=0 ! counter to count number of calls to physics
-!$OMP THREADPRIVATE(itau)
-REAL :: temp_newton(klon,klev)
-INTEGER :: k
-logical, save :: first=.TRUE.
-!$OMP THREADPRIVATE(first)
+    INTEGER, save :: itau = 0 ! counter to count number of calls to physics
+    !$OMP THREADPRIVATE(itau)
+    REAL :: temp_newton(klon, klev)
+    INTEGER :: k
+    logical, save :: first = .TRUE.
+    !$OMP THREADPRIVATE(first)
 
-! For I/Os
-INTEGER :: itau0
-REAL :: zjulian
-REAL :: dtime
-INTEGER :: nhori ! horizontal coordinate ID
-INTEGER,save :: nid_hist ! output file ID
-!$OMP THREADPRIVATE(nid_hist)
-INTEGER :: zvertid ! vertical coordinate ID
-INTEGER,save :: iwrite_phys ! output every iwrite_phys physics step
-!$OMP THREADPRIVATE(iwrite_phys)
-INTEGER,save :: iwrite_phys_omp ! intermediate variable to read iwrite_phys
-                                ! (must be shared by all threads)
-REAL :: t_ops ! frequency of the IOIPSL operations (eg average over...)
-REAL :: t_wrt ! frequency of the IOIPSL outputs
+    ! For I/Os
+    INTEGER :: itau0
+    REAL :: zjulian
+    REAL :: dtime
+    INTEGER :: nhori ! horizontal coordinate ID
+    INTEGER, save :: nid_hist ! output file ID
+    !$OMP THREADPRIVATE(nid_hist)
+    INTEGER :: zvertid ! vertical coordinate ID
+    INTEGER, save :: iwrite_phys ! output every iwrite_phys physics step
+    !$OMP THREADPRIVATE(iwrite_phys)
+    INTEGER, save :: iwrite_phys_omp ! intermediate variable to read iwrite_phys
+    ! (must be shared by all threads)
+    REAL :: t_ops ! frequency of the IOIPSL operations (eg average over...)
+    REAL :: t_wrt ! frequency of the IOIPSL outputs
 
-! initializations
-IF (debut) then ! Things to do only for the first CALL to physics
-! load initial conditions for physics (including the grid)
-  CALL phys_state_var_init() ! some initializations, required before calling phyetat0
-  CALL phyetat0("startphy.nc")
+    ! initializations
+    IF (debut) then ! Things to do only for the first CALL to physics
+      ! load initial conditions for physics (including the grid)
+      CALL phys_state_var_init() ! some initializations, required before calling phyetat0
+      CALL phyetat0("startphy.nc")
 
-! Initialize outputs:
-  itau0=0
-!$OMP MASTER
-  iwrite_phys_omp=1 !default: output every physics timestep
-  ! NB: getin() is not threadsafe; only one thread should CALL it.
-  CALL getin("iwrite_phys",iwrite_phys_omp)
-!$OMP END MASTER
-!$OMP BARRIER
-  iwrite_phys=iwrite_phys_omp
-  t_ops=pdtphys*iwrite_phys ! frequency of the IOIPSL operation
-  t_wrt=pdtphys*iwrite_phys ! frequency of the outputs in the file
-  ! compute zjulian for annee0=1979 and month=1 dayref=1 and hour=0.0
-  !CALL ymds2ju(annee0, month, dayref, hour, zjulian)
-  CALL ymds2ju(1979, 1, 1, 0.0, zjulian)
-  dtime=pdtphys
+      ! Initialize outputs:
+      itau0 = 0
+      !$OMP MASTER
+      iwrite_phys_omp = 1 !default: output every physics timestep
+      ! NB: getin() is not threadsafe; only one thread should CALL it.
+      CALL getin("iwrite_phys", iwrite_phys_omp)
+      !$OMP END MASTER
+      !$OMP BARRIER
+      iwrite_phys = iwrite_phys_omp
+      t_ops = pdtphys * iwrite_phys ! frequency of the IOIPSL operation
+      t_wrt = pdtphys * iwrite_phys ! frequency of the outputs in the file
+      ! compute zjulian for annee0=1979 and month=1 dayref=1 and hour=0.0
+      !CALL ymds2ju(annee0, month, dayref, hour, zjulian)
+      CALL ymds2ju(1979, 1, 1, 0.0, zjulian)
+      dtime = pdtphys
 #ifndef CPP_IOIPSL_NO_OUTPUT
-  ! Initialize IOIPSL output file
-  CALL histbeg_phy("histins.nc",itau0,zjulian,dtime,nhori,nid_hist)
+      ! Initialize IOIPSL output file
+      CALL histbeg_phy("histins.nc", itau0, zjulian, dtime, nhori, nid_hist)
 #endif
 
-!$OMP MASTER
+      !$OMP MASTER
 
 #ifndef CPP_IOIPSL_NO_OUTPUT 
-! IOIPSL
-  ! define vertical coordinate
-  CALL histvert(nid_hist,"presnivs","Vertical levels","Pa",klev, &
-                presnivs,zvertid,'down')
-  ! define variables which will be written in "histins.nc" file
-  CALL histdef(nid_hist,'temperature','Atmospheric temperature','K', &
-               nbp_lon,jj_nb,nhori,klev,1,klev,zvertid,32, &
-               'inst(X)',t_ops,t_wrt)
-  CALL histdef(nid_hist,'u','Eastward Zonal Wind','m/s', &
-               nbp_lon,jj_nb,nhori,klev,1,klev,zvertid,32, &
-               'inst(X)',t_ops,t_wrt)
-  CALL histdef(nid_hist,'v','Northward Meridional Wind','m/s', &
-               nbp_lon,jj_nb,nhori,klev,1,klev,zvertid,32, &
-               'inst(X)',t_ops,t_wrt)
-  CALL histdef(nid_hist,'ps','Surface Pressure','Pa', &
-               nbp_lon,jj_nb,nhori,1,1,1,zvertid,32, &
-               'inst(X)',t_ops,t_wrt)
-  ! end definition sequence
-  CALL histend(nid_hist)
+      ! IOIPSL
+      ! define vertical coordinate
+      CALL histvert(nid_hist, "presnivs", "Vertical levels", "Pa", klev, &
+              presnivs, zvertid, 'down')
+      ! define variables which will be written in "histins.nc" file
+      CALL histdef(nid_hist, 'temperature', 'Atmospheric temperature', 'K', &
+              nbp_lon, jj_nb, nhori, klev, 1, klev, zvertid, 32, &
+              'inst(X)', t_ops, t_wrt)
+      CALL histdef(nid_hist, 'u', 'Eastward Zonal Wind', 'm/s', &
+              nbp_lon, jj_nb, nhori, klev, 1, klev, zvertid, 32, &
+              'inst(X)', t_ops, t_wrt)
+      CALL histdef(nid_hist, 'v', 'Northward Meridional Wind', 'm/s', &
+              nbp_lon, jj_nb, nhori, klev, 1, klev, zvertid, 32, &
+              'inst(X)', t_ops, t_wrt)
+      CALL histdef(nid_hist, 'ps', 'Surface Pressure', 'Pa', &
+              nbp_lon, jj_nb, nhori, 1, 1, 1, zvertid, 32, &
+              'inst(X)', t_ops, t_wrt)
+      ! end definition sequence
+      CALL histend(nid_hist)
 #endif
 
-    IF (using_xios) THEN
-!XIOS
-      ! Declare available vertical axes to be used in output files:    
-      CALL wxios_add_vaxis("presnivs", klev, presnivs)
+      IF (using_xios) THEN
+        !XIOS
+        ! Declare available vertical axes to be used in output files:
+        CALL wxios_add_vaxis("presnivs", klev, presnivs)
 
-      ! Declare calendar and time step
-      CALL wxios_set_cal(dtime,"earth_360d",1,1,1,0.0,1,1,1,0.0)
-    
-      !Finalize the context:
-      CALL wxios_closedef()
-    ENDIF
-!$OMP END MASTER
-!$OMP BARRIER
-END IF ! of if (debut)
+        ! Declare calendar and time step
+        CALL wxios_set_cal(dtime, "earth_360d", 1, 1, 1, 0.0, 1, 1, 1, 0.0)
 
-! increment local time counter itau
-itau=itau+1
+        !Finalize the context:
+        CALL wxios_closedef()
+      ENDIF
+      !$OMP END MASTER
+      !$OMP BARRIER
+    END IF ! of if (debut)
 
-! set all tendencies to zero
-d_u(1:klon,1:klev)=0.
-d_v(1:klon,1:klev)=0.
-d_t(1:klon,1:klev)=0.
-d_qx(1:klon,1:klev,1:nqtot)=0.
-d_ps(1:klon)=0.
+    ! increment local time counter itau
+    itau = itau + 1
 
-! compute tendencies to return to the dynamics:
-! "friction" on the first layer
-d_u(1:klon,1)=-u(1:klon,1)/86400.
-d_v(1:klon,1)=-v(1:klon,1)/86400.
-! newtonian relaxation towards temp_newton()
-DO k=1,klev
-  temp_newton(1:klon,k)=280.+cos(latitude(1:klon))*40.-pphi(1:klon,k)/rg*6.e-3
-  d_t(1:klon,k)=(temp_newton(1:klon,k)-t(1:klon,k))/1.e5
-END DO
+    ! set all tendencies to zero
+    d_u(1:klon, 1:klev) = 0.
+    d_v(1:klon, 1:klev) = 0.
+    d_t(1:klon, 1:klev) = 0.
+    d_qx(1:klon, 1:klev, 1:nqtot) = 0.
+    d_ps(1:klon) = 0.
 
+    ! compute tendencies to return to the dynamics:
+    ! "friction" on the first layer
+    d_u(1:klon, 1) = -u(1:klon, 1) / 86400.
+    d_v(1:klon, 1) = -v(1:klon, 1) / 86400.
+    ! newtonian relaxation towards temp_newton()
+    DO k = 1, klev
+      temp_newton(1:klon, k) = 280. + cos(latitude(1:klon)) * 40. - pphi(1:klon, k) / rg * 6.e-3
+      d_t(1:klon, k) = (temp_newton(1:klon, k) - t(1:klon, k)) / 1.e5
+    END DO
 
-PRINT*,'PHYDEV: itau=',itau
+    PRINT*, 'PHYDEV: itau=', itau
 
-! write some outputs:
-! IOIPSL
+    ! write some outputs:
+    ! IOIPSL
 #ifndef CPP_IOIPSL_NO_OUTPUT 
-IF (modulo(itau,iwrite_phys)==0) THEN
-  CALL histwrite_phy(nid_hist,.FALSE.,"temperature",itau,t)
-  CALL histwrite_phy(nid_hist,.FALSE.,"u",itau,u)
-  CALL histwrite_phy(nid_hist,.FALSE.,"v",itau,v)
-  CALL histwrite_phy(nid_hist,.FALSE.,"ps",itau,paprs(:,1))
-END IF
+    IF (modulo(itau, iwrite_phys)==0) THEN
+      CALL histwrite_phy(nid_hist, .FALSE., "temperature", itau, t)
+      CALL histwrite_phy(nid_hist, .FALSE., "u", itau, u)
+      CALL histwrite_phy(nid_hist, .FALSE., "v", itau, v)
+      CALL histwrite_phy(nid_hist, .FALSE., "ps", itau, paprs(:, 1))
+    END IF
 #endif
 
-!XIOS
-IF (using_xios) THEN
-  !$OMP MASTER
-    !Increment XIOS time
-    CALL xios_update_calendar(itau)
-  !$OMP END MASTER
-  !$OMP BARRIER
+    !XIOS
+    IF (using_xios) THEN
+      !$OMP MASTER
+      !Increment XIOS time
+      CALL xios_update_calendar(itau)
+      !$OMP END MASTER
+      !$OMP BARRIER
 
-    !Send fields to XIOS: (NB these fields must also be defined as
-    ! <field id="..." /> in iodef.xml to be correctly used
-    CALL histwrite_phy("temperature",t)
-    CALL histwrite_phy("temp_newton",temp_newton)
-    CALL histwrite_phy("u",u)
-    CALL histwrite_phy("v",v)
-    CALL histwrite_phy("ps",paprs(:,1))
-ENDIF
+      !Send fields to XIOS: (NB these fields must also be defined as
+      ! <field id="..." /> in iodef.xml to be correctly used
+      CALL histwrite_phy("temperature", t)
+      CALL histwrite_phy("temp_newton", temp_newton)
+      CALL histwrite_phy("u", u)
+      CALL histwrite_phy("v", v)
+      CALL histwrite_phy("ps", paprs(:, 1))
+    ENDIF
 
-! if lastcall, then it is time to write "restartphy.nc" file
-IF (lafin) THEN
-  CALL phyredem("restartphy.nc")
-END IF
+    ! if lastcall, then it is time to write "restartphy.nc" file
+    IF (lafin) THEN
+      CALL phyredem("restartphy.nc")
+    END IF
 
-END SUBROUTINE  physiq
+  END SUBROUTINE  physiq
 
 END MODULE physiq_mod