Changeset 4033

Timestamp:

Jan 28, 2026, 3:06:52 PM (4 weeks ago)

Author:

mmaurice

Message:

Generic PCM

Introduce molecular viscosity module. Add callconduc and callmolvis
flags to call conduction and molecular viscosity, respectively. Change
phitop in conduction to phitop_conduc (there is also phitop_molvis now).
Move thermal conductivity lambda to conc_mod.F90.

MM

Location:

trunk/LMDZ.GENERIC/libf/phystd

Files:

• callkeys_mod.F90 (modified) (1 diff)
• conc_mod.F90 (modified) (3 diffs)
• conduction.F90 (modified) (4 diffs)
• inifis_mod.F90 (modified) (1 diff)
• molvis.F90 (added)
• physiq_mod.F90 (modified) (5 diffs)

trunk/LMDZ.GENERIC/libf/phystd/callkeys_mod.F90

@@ -22 +22 @@
 !$OMP THREADPRIVATE(strictboundrayleigh)
       logical,save :: callthermos
+      logical,save :: callconduc
+      logical,save :: callmolvis
       logical,save :: force_conduction
-      real,save    :: phitop
+      real,save    :: phitop_conduc
+      real,save    :: phitop_molvis
       real,save    :: zztop
       real,save    :: a_coeff
       real,save    :: s_coeff
-!$OMP THREADPRIVATE(callthermos,phitop,zztop,force_conduction,a_coeff,s_coeff)
+!$OMP THREADPRIVATE(callthermos,callconduc,callmolvis,phitop_conduc,phitop_molvis,zztop,force_conduction,a_coeff,s_coeff)
 
       logical,save :: enertest

trunk/LMDZ.GENERIC/libf/phystd/conc_mod.F90

@@ -5 +5 @@
   real,save,allocatable :: mmean(:,:)  ! mean molecular mass of the atmosphere
   real,save,allocatable :: Akknew(:,:) ! thermal conductivity cofficient
+  real,save,allocatable :: lambda(:,:) ! effective thermal conductivity
   real,save,allocatable :: cpnew(:,:)  ! specicic heat
   real,save,allocatable :: rnew(:,:)   ! specific gas constant
@@ -18 +19 @@
     allocate(mmean(ngrid,nlayer))
     allocate(Akknew(ngrid,nlayer))
+    allocate(lambda(ngrid,nlayer))
     allocate(cpnew(ngrid,nlayer))
     allocate(rnew(ngrid,nlayer))
@@ -29 +31 @@
     if (allocated(mmean)) deallocate(mmean)
     if (allocated(Akknew)) deallocate(Akknew)
+    if (allocated(lambda)) deallocate(lambda)
     if (allocated(cpnew)) deallocate(cpnew)
     if (allocated(rnew)) deallocate(rnew)

trunk/LMDZ.GENERIC/libf/phystd/conduction.F90

@@ -9 +9 @@
     
       use comcstfi_mod, only: r, cpp, mugaz
-      use callkeys_mod, only: phitop,zztop,a_coeff,s_coeff,force_conduction
+      use callkeys_mod, only: phitop_conduc,zztop,a_coeff,s_coeff,force_conduction
+      use conc_mod,     only: lambda
       use gases_h
 
@@ -48 +49 @@
       INTEGER :: i,ig,l,igas,kgas
       REAL :: alpha(ngrid,nlayer)
-      REAL :: lambda(ngrid,nlayer)
       REAL :: muvol(ngrid,nlayer)   ! kg.m-3
       REAL :: C(ngrid,nlayer)
@@ -111 +111 @@
               molar_frac(:,:,igas) = gfrac(igas)
               here(igas) = .true.
+            elseif(igas.eq.igas_CO2) then
+              !Interpolated from Hurly et al., (2007)
+              !valid between 20 and 1000 K (max 3 percent of error)
+              akk(igas) = 3.072e-4
+              skk(igas) = 0.69
+              molar_mass(igas) = 44.010e-3
+              akk_visc(igas) = 4.5054e-7
+              skk_visc(igas) = 0.6658
+              molar_frac(:,:,igas) = gfrac(igas)
+              here(igas) = .true.
               ! Add more molecules here and the reference PLEASE!
             else 
@@ -205 +215 @@
                    * (1-D(:,nlayer-1))
       C(:,nlayer)  =C(:,nlayer-1)+zt(:,nlayer-1)-zt(:,nlayer) 
-      C(:,nlayer)  =(C(:,nlayer)*lambda(:,nlayer)+phitop) &
+      C(:,nlayer)  =(C(:,nlayer)*lambda(:,nlayer)+phitop_conduc) &
                    / den(:,nlayer) 
 

trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90

@@ -292 +292 @@
      if(callthermos) then
        if (is_master) write(*,*) trim(rname)//&
+        ": call thermospheric conduction ?"
+       callconduc=.true. ! default value
+       call getin_p("conduction",callconduc)
+       if (is_master) write(*,*) trim(rname)//&
+         ": conduction = ",callconduc
+       if (is_master) write(*,*) trim(rname)//&
+         ": call thermospheric viscosity ?"
+         callmolvis=.true. ! default value
+       call getin_p("molvis",callmolvis)
+       if (is_master) write(*,*) trim(rname)//&
+          ": molvis = ",callmolvis
+       if (is_master) write(*,*) trim(rname)//&
          ": flux from thermosphere ? W/m2"
-       phitop = 0.0 ! default value
-       call getin_p("phitop",phitop)
+       phitop_conduc = 0.0 ! default value
+       call getin_p("phitop_conduc",phitop_conduc)
        if (is_master) write(*,*) trim(rname)//&
-         ": phitop = ",phitop
+         ": phitop_conduc = ",phitop_conduc
+        if (is_master) write(*,*) trim(rname)//&
+         ": momentum flux from thermosphere ?"
+       phitop_molvis = 0.0 ! default value
+       call getin_p("phitop_molvis",phitop_molvis)
+        if (is_master) write(*,*) trim(rname)//&
+         ": phitop_molvis = ",phitop_molvis
        if (is_master) write(*,*) trim(rname)//&
          ": Thickness of conduction ? m"

trunk/LMDZ.GENERIC/libf/phystd/physiq_mod.F90

@@ -71 +71 @@
                               ok_slab_ocean, photochem, rings_shadow, rmean, &
                               season, sedimentation,generic_condensation, &
-                              specOLR, callthermos, callvolcano, &
+                              specOLR, callvolcano, &
+                              callthermos, callconduc, callmolvis, &
                               startphy_file, testradtimes, tlocked, &
                               tracer, UseTurbDiff, water, watercond, &
@@ -83 +84 @@
       use callcorrk_mod, only: callcorrk
       use conduction_mod, only: conduction
+      use molvis_mod, only: molvis
       use volcano_mod, only: volcano
       use pindex_mod, only: pindex
@@ -382 +384 @@
       real zdhdif(ngrid,nlayer)                             ! Turbdiff/vdifc routines.
       real zdhadj(ngrid,nlayer)                             ! Convadj routine.
+      real zdumolvis(ngrid,nlayer)                          ! Molecular viscosity (thermosphere)
+      real zdvmolvis(ngrid,nlayer)                          ! Molecular viscosity (thermosphere)
 
       ! For Pressure and Mass :
@@ -541 +545 @@
          zdtlw(:,:) = 0.0
          zdtconduc(:,:) = 0.0
+         zdumolvis(:,:) = 0.0
+         zdvmolvis(:,:) = 0.0
 
 !        Initialize fixed variable mu
@@ -2384 +2390 @@
 ! -------------------------
       if (callthermos) then
-        call conduction(ngrid,nlayer,nq,ptimestep,pplay,pplev,  &
-                          pt,tsurf,zzlev,zzlay,muvar,pq,firstcall,zdtconduc)
-        pdt(1:ngrid,1:nlayer)=pdt(1:ngrid,1:nlayer)+zdtconduc(1:ngrid,1:nlayer)
+
+        if (callconduc) then
+          call conduction(ngrid,nlayer,nq,ptimestep,pplay,pplev,  &
+                            pt,tsurf,zzlev,zzlay,muvar,pq,firstcall,zdtconduc)
+          pdt(1:ngrid,1:nlayer)=pdt(1:ngrid,1:nlayer)+zdtconduc(1:ngrid,1:nlayer)
+        endif
+      
+        if (callmolvis) then
+          ! u
+          call molvis(ngrid,nlayer,ptimestep,pplay,pplev,pt,pdt, &
+                      pu,zzlev,zzlay,zdumolvis)
+          pdu(1:ngrid,1:nlayer) = pdu(1:ngrid,1:nlayer) + zdumolvis(1:ngrid,1:nlayer)
+          ! v
+          call molvis(ngrid,nlayer,ptimestep,pplay,pplev,pt,pdt, &
+                      pv,zzlev,zzlay,zdvmolvis)
+          pdv(1:ngrid,1:nlayer) = pdv(1:ngrid,1:nlayer) + zdvmolvis(1:ngrid,1:nlayer)
+        endif
+
       endif ! of if (callthermos)