Changeset 3838


Ignore:
Timestamp:
Jul 7, 2025, 3:11:29 PM (3 days ago)
Author:
jbclement
Message:

Mars PCM 1D:
Fixing formula for prescribing atmospheric water profile in the 1D model (already under discretization) + Addition of more explicit comments.
JBC

Location:
trunk/LMDZ.MARS
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.MARS/changelog.txt

    r3831 r3838  
    48944894
    48954895== 04/07/2025 == JBC
    4896 Imposing the prescribed atmospheric water profile though 'dq' instead of 'q' directly (more correct).
     4896Imposing the prescribed atmospheric water profile though 'dq' instead of 'q' directly for the 1D model (more correct).
     4897
     4898== 04/07/2025 == JBC
     4899More explicit code and comments about the prescribtion of the atmospheric water profile in the 1D model.
     4900
     4901== 07/07/2025 == JBC
     4902Fixing formula for prescribing atmospheric water profile in the 1D model (already under discretization) + Addition of more explicit comments.
  • trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F90

    r3836 r3838  
    195195    ! Wind increment: specific for 1D
    196196    ! -------------------------------
    197     ! The physics compute the tendencies on u and v,
    198     ! here we just add Coriolos effect
     197    ! The physics computes the tendencies on u and v, here we just add Coriolos effect
    199198    !do ilayer = 1,nlayer
    200199    !    du(ilayer) = du(ilayer) + ptif*(v(ilayer) - grv)
     
    207206    !endif
    208207
     208    ! Water tracer increment: specific for 1D
     209    ! ---------------------------------------
     210    ! The physics computes the tendency on q(1,:,igcm_h2o_vap), here we can mimic an effect of the "3D dynamics"
     211    ! either by forcing a profile or by relaxing towards a prescribed profile
     212    if (water .and. prescribed_h2ovap) then ! If atmospheric water profile is prescribed
     213        if (h2ovap_relax_time < 0.) then ! Forcing
     214            ! For some tests: unless it reaches saturation (maximal value)
     215            !call watersat(nlayer,temp,play,zqsat)
     216            !dq(1,:,igcm_h2o_vap) = (min(zqsat(:),q_prescribed_h2o_vap(:)) - q(1,:,igcm_h2o_vap))/dttestphys
     217            dq(1,:,igcm_h2o_vap) = (q_prescribed_h2o_vap(:) - q(1,:,igcm_h2o_vap))/dttestphys
     218        else ! Relaxation
     219            dq(1,:,igcm_h2o_vap) = dq(1,:,igcm_h2o_vap) - (q(1,:,igcm_h2o_vap) - q_prescribed_h2o_vap(:))/h2ovap_relax_time
     220        endif
     221     endif
     222
    209223    ! Compute time for next time step
    210224    ! -------------------------------
     
    231245    play = aps + psurf*bps
    232246
    233     ! Prescription of atmospheric water if asked
    234     ! ------------------------------------------
    235     if (water .and. prescribed_h2ovap) then ! If atmospheric water profile is prescribed
    236         if (h2ovap_relax_time < 0.) then ! Forcing
    237             !call watersat(nlayer,temp,play,zqsat)
    238             dq(1,:,igcm_h2o_vap) = (q_prescribed_h2o_vap(:) - q(1,:,igcm_h2o_vap))/dttestphys
    239         else ! Relaxation
    240             dq(1,:,igcm_h2o_vap) = (q_prescribed_h2o_vap(:) + (q(1,:,igcm_h2o_vap) - q_prescribed_h2o_vap(:))*dexp(-dttestphys/h2ovap_relax_time) - q(1,:,igcm_h2o_vap))/dttestphys
    241         endif
    242      endif
    243 
    244247    ! Compute tracers for next time step
    245248    ! ----------------------------------
Note: See TracChangeset for help on using the changeset viewer.