Changeset 3974


Ignore:
Timestamp:
Nov 21, 2025, 4:28:06 PM (8 weeks ago)
Author:
jbclement
Message:

Mars PCM:
Simplification of the orbital parameters initialization (resolve ticket #109): removing redundant actions in subroutine 'iniorbit' in regard of what is done in "tabfi.F" and moving the subroutine 'lsp2solp' from "tabfi.F" to "planete_h.F90" + some cleanings.
JBC

Location:
trunk
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.COMMON/libf/evolution/pem.F90

    r3967 r3974  
    8686use tracer_mod,                 only: noms, igcm_h2o_ice, igcm_co2, mmol, igcm_h2o_vap ! Tracer names and molar masses
    8787use mod_phys_lmdz_para,         only: is_parallel, is_sequential, is_mpi_root, is_omp_root, is_master
    88 use planete_h,                  only: aphelie, periheli, year_day, peri_day, obliquit, iniorbit
    8988use surfini_mod,                only: surfini
    9089use comcstfi_h,                 only: mugaz
     
    597596!    I_i Compute orbit criterion
    598597!------------------------
    599 call iniorbit(aphelie,periheli,year_day,peri_day,obliquit)
    600 
    601598n_myear_leg = Max_iter_pem
    602599if (evol_orbit_pem) call orbit_param_criterion(i_myear,n_myear_leg)
  • trunk/LMDZ.MARS/changelog.txt

    r3969 r3974  
    40924092
    40934093== 07/07/2023 == JBC
    4094 Minor fix for initialization of tracers indices and rework of "atm_wat_profile" in testphys1d.F.
     4094Minor fix for initialization of tracers indices and rework of 'atm_wat_profile' in "testphys1d.F".
    40954095
    40964096== 10/07/2023 == EM
     
    41544154
    41554155== 31/07/2023 == JBC
    4156 Addition in 1D of relaxation towards a profile for atmospheric water. It needs the flag "atm_wat_tau" (real):
     4156Addition in 1D of relaxation towards a profile for atmospheric water. It needs the flag 'atm_wat_tau' (real):
    41574157 - atm_wat_tau < 0. -> no relaxation (default);
    4158  - atm_wat_tau >= 0. -> relaxation towards the value "atm_wat_profile" with relaxation time "atm_wat_tau".
     4158 - atm_wat_tau >= 0. -> relaxation towards the value "atm_wat_profile" with relaxation time 'atm_wat_tau'.
    41594159
    41604160== 07/08/2023 == EM
     
    41634163
    41644164== 09/08/2023 == JBC
    4165 Small fixes to be able to run the Mars PCM 1D without "water".
     4165Small fixes to be able to run the Mars PCM 1D without 'water'.
    41664166
    41674167== 14/08/2023 == JBC
    4168 Related to commit r3026: improvement of error message in initracer.F (now it gives correctly the only identified tracers) + one small correction to run PCM 1D without water.
     4168Related to commit r3026: improvement of error message in "initracer.F" (now it gives correctly the only identified tracers) + one small correction to run PCM 1D without 'water'.
    41694169
    41704170== 8/09/2023 == JBC
    4171 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.
     4171In 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.
    41724172Some code cleaning in regards of tests in 1D.
    41734173
     
    42654265    - Correction of a bug in "writediagfi.F": the case of using the 1D model with parallelization was not anticipated so that the "diagfi.nc" file was filled with NaNf;
    42664266    - Addition of the file "start1D.txt" as an example in the directory deftank/;
    4267     - Some "cosmetic" modifications in "improvedclouds_mod.F", "write_output_mod.F90" and "testphys1d.F90".
     4267    - Some cosmetic modifications in "improvedclouds_mod.F", "write_output_mod.F90" and "testphys1d.F90".
    42684268
    42694269== 22/10/2023 == EM
     
    48634863
    48644864== 12/06/2025 == JBC
    4865 Move paleoclimate variables from time-independent to time-dependent writing of "stratfi.nc".
     4865Move paleoclimate variables from time-independent to time-dependent writing of "startfi.nc".
    48664866
    48674867== 16/06/2025 == JM
     
    50035003Bug fix: threadsave statement missing on variable save "isample". Could fail
    50045004when writting diagfi.
     5005
     5006== 21/11/2025 == JBC
     5007Simplification of the orbital parameters initialization (resolve ticket #109): removing redundant actions in subroutine 'iniorbit' in regard of what is done in "tabfi.F" and moving the subroutine 'lsp2solp' from "tabfi.F" to "planete_h.F90" + some cleanings.
  • trunk/LMDZ.MARS/deftank/run.def.1d

    r3943 r3974  
    5454# Prescription of atmospheric water vapour profile defined in file "profile_ref_h2o_vap" (kg/kg)
    5555ctrl_h2ovap=.false.
    56 # Relaxation time (s). If <0 then no time relaxation i.e. forcing
     56# Relaxation time (s) towards atmospheric water vapour profile. If <0 then no time relaxation i.e. forcing
    5757relaxtime_h2ovap=-1.
    5858# Prescription of atmospheric water ice profile defined in file "profile_ref_h2o_ice" (kg/kg)
    5959ctrl_h2oice=.false.
    60 # Relaxation time (s). If <0 then no time relaxation i.e. forcing
     60# Relaxation time (s) towards atmospheric water ice profile. If <0 then no time relaxation i.e. forcing
    6161relaxtime_h2oice=-1.
    6262# hybrid vertical coordinate ? (.true. for hybrid and .false. for sigma levels)
  • trunk/LMDZ.MARS/libf/phymars/phyetat0_mod.F90

    r3964 r3974  
    3535use comcstfi_h,          only: pi
    3636use geometry_mod,        only: latitude
    37 use soil_settings_mod, only: soil_settings
    38 use tabfi_mod, only: tabfi
    39 use callkeys_mod, only: startphy_file, rdstorm, hdo
    40 use callkeys_mod, only: CLFvarying, CLFfixval
     37use soil_settings_mod,   only: soil_settings
     38use tabfi_mod,           only: tabfi
     39use callkeys_mod,        only: startphy_file, rdstorm, hdo, CLFvarying, CLFfixval
    4140
    4241implicit none
  • trunk/LMDZ.MARS/libf/phymars/physiq_mod.F

    r3918 r3974  
    8585      use turb_mod, only: q2, wstar, ustar, sensibFlux,
    8686     &                    zmax_th, hfmax_th, turb_resolved
    87       use planete_h, only: aphelie, periheli, year_day, peri_day,
    88      &                     obliquit
    89       use planete_h, only: iniorbit
     87      use planete_h, only: obliquit
    9088      use orbite_mod, only: orbite
    9189      USE comcstfi_h, only: r, cpp, mugaz, g, rcp, pi, rad
     
    718716c        ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    719717         CALL surfini(ngrid,nslope,qsurf)
    720          CALL iniorbit(aphelie,periheli,year_day,peri_day,obliquit)
     718
    721719c        initialize soil
    722720c        ~~~~~~~~~~~~~~~
  • trunk/LMDZ.MARS/libf/phymars/planete_h.F90

    r3095 r3974  
    11MODULE planete_h
    2 IMPLICIT NONE
     2
     3implicit none
    34
    45REAL,SAVE :: aphelie   ! Aphelion, in Mkm
     
    2324
    2425CONTAINS
     26!-----------------------------------------------------------------------
     27! Gives sol at perihelion for ls at perihelion (for precession cycles)
     28  SUBROUTINE lsp2solp(lsp,solp)
     29
     30  use comcstfi_h, only: pi
     31
     32  implicit none
     33
     34  ! Arguments:
     35  real, intent(in) :: lsp ! ls at perihelion
     36  real, intent(out) :: solp ! sol at perihelion
     37
     38  ! Locals:
     39  real :: zx0 ! eccentric anomaly at Ls=0
     40  real, parameter :: degrad = 180.d0/pi
     41
     42  ! Compute orbit geometrical parameters
     43  e_elips = (aphelie - periheli)/(aphelie + periheli)
     44  zx0 = -2.0*datan(dtan(0.5*lsp/degrad)*dsqrt((1. - e_elips)/(1. + e_elips)))
     45  if (zx0 <= 0.) zx0 = zx0 + 2.*pi
     46
     47  ! Compute sol at perihelion
     48  solp = year_day*(1. - (zx0 - e_elips*dsin(zx0))/(2.*pi))
     49
     50  END SUBROUTINE lsp2solp
     51
     52!-----------------------------------------------------------------------
     53! Initialize orbital parameters
     54  SUBROUTINE iniorbit()
     55
     56  use comcstfi_h, only: pi
     57
     58  implicit none
     59
     60  ! Locals:
     61  character(8), parameter :: myname = "iniorbit"
     62  real    :: zxref, zanom, zz, zx0, zdx
     63  integer :: iter
     64
     65  ! Info about orbital values
     66  write(*,*)myname,': Aphelion in Mkm                 =',aphelie
     67  write(*,*)myname,': Perihelion in Mkm               =',periheli
     68  write(*,*)myname,': Number of days in the year      =',year_day
     69  write(*,*)myname,': Date of perihelion in days      =',peri_day
     70  write(*,*)myname,': Obliquity in degrees            =',obliquit
     71
     72  ! Compute orbit geometrical parameters
     73  e_elips = (aphelie - periheli)/(aphelie + periheli)
     74  p_elips = 0.5*(periheli + aphelie)*(1. - e_elips*e_elips)/unitastr
     75  write(*,*)myname,': Orbit eccentricity              =',e_elips
     76  write(*,*)myname,': Ellipse semi-latus rectum       =',p_elips
    2577     
    26   SUBROUTINE iniorbit(paphelie,pperiheli,pyear_day,pperi_day,pobliq)
    27   use comcstfi_h, only: pi
    28   IMPLICIT NONE
    29   ! initialize orbital parameters
    30  
    31   REAL,INTENT(IN) :: paphelie,pperiheli,pyear_day,pperi_day,pobliq
     78  ! Compute mean anomaly zanom
     79  zz = (year_day - peri_day)/year_day
     80  zanom = 2.*pi*(zz - nint(zz))
     81  zxref = abs(zanom)
     82  write(*,*)myname,': zanom                           =',zanom
    3283
    33   CHARACTER(LEN=8),PARAMETER :: myname="iniorbit"
    34   REAL :: zxref,zanom,zz,zx0,zdx
    35   INTEGER :: iter
    36 
    37   ! copy over input values
    38   aphelie =paphelie
    39   periheli=pperiheli
    40   year_day=pyear_day
    41   obliquit=pobliq
    42   peri_day=pperi_day
    43  
    44   write(*,*)myname,': Perihelion in Mkm  ',periheli
    45   write(*,*)myname,': Aphelion  in Mkm  ',aphelie
    46   write(*,*)myname,': Obliquity in degrees  :',obliquit
    47      
    48   ! compute orbit geometrical parameters
    49   e_elips=(aphelie-periheli)/(periheli+aphelie)
    50   p_elips=0.5*(periheli+aphelie)*(1-e_elips*e_elips)/unitastr
    51      
    52   write(*,*)myname,': e_elips',e_elips
    53   write(*,*)myname,': p_elips',p_elips
    54      
    55   !  compute mean anomaly zanom
    56 
    57   zz=(year_day-pperi_day)/year_day
    58   zanom=2.*pi*(zz-nint(zz))
    59   zxref=abs(zanom)
    60   write(*,*)myname,': zanom  ',zanom
    61 
    62   ! solve equation  zx0 - e * sin (zx0) = zxref for eccentric anomaly zx0
     84  ! Solve equation  zx0 - e * sin (zx0) = zxref for eccentric anomaly zx0
    6385  ! using Newton method
    6486
    65   zx0=zxref+e_elips*sin(zxref)
    66   DO iter=1,100
    67     zdx=-(zx0-e_elips*sin(zx0)-zxref)/(1.-e_elips*cos(zx0))
    68     if(abs(zdx).le.(1.e-12)) exit
    69     zx0=zx0+zdx
     87  zx0 = zxref + e_elips*sin(zxref)
     88  DO iter = 1,100
     89    zdx = -(zx0 - e_elips*sin(zx0) - zxref)/(1. - e_elips*cos(zx0))
     90    if (abs(zdx).le.(1.e-12)) exit
     91    zx0 = zx0 + zdx
    7092  ENDDO
    7193
    72   zx0=zx0+zdx
    73   if(zanom.lt.0.) zx0=-zx0
    74   write(*,*)myname,': zx0   ',zx0
     94  zx0 = zx0 + zdx
     95  if (zanom < 0.) zx0 = -zx0
     96  write(*,*)myname,': Eccentric anomaly zx0           =',zx0
    7597
    76   lsperi=-2.*atan(sqrt((1.+e_elips)/(1.-e_elips))*tan(zx0/2.))
     98  lsperi = -2.*atan(sqrt((1. + e_elips)/(1. - e_elips))*tan(zx0/2.))
    7799  lsperi = modulo(lsperi,2.*pi)
    78   write(*,*)myname,': Perihelion solar long. Ls (deg)=', &
    79                           lsperi*180./pi
     100  write(*,*)myname,': Perihelion solar long. Ls (deg) =', lsperi*180./pi
    80101
    81102  END SUBROUTINE iniorbit
  • trunk/LMDZ.MARS/libf/phymars/tabfi.F

    r3902 r3974  
    6161      use time_phylmdz_mod, only: daysec, dtphys
    6262      use planete_h, only: aphelie, emin_turb, lmixmin, obliquit,
    63      &                     peri_day, periheli, year_day
     63     &                     peri_day, periheli, year_day,
     64     &                     lsp2solp, iniorbit
     65
    6466      implicit none
    6567 
     
    460462          write(*,*) 'perihelion =',periheli
    461463          write(*,*) 'and',year_day,'sols/year'
    462           call lsp2solp(peri_ls,peri_day,aphelie,periheli,year_day)
     464          call lsp2solp(peri_ls,peri_day)
    463465          write(*,*) 'peri_day (new value):',peri_day
    464466
     
    543545
    544546c-----------------------------------------------------------------------
     547c       Initialization of orbital parameters
     548c-----------------------------------------------------------------------
     549      call iniorbit()
     550
     551c-----------------------------------------------------------------------
    545552c       Case when using a start file from before March 1996 (without iceradius...
    546553c-----------------------------------------------------------------------
     
    555562       end if
    556563
    557 c-----------------------------------------------------------------------
    558564      END SUBROUTINE tabfi
    559565
    560 
    561 
    562      
    563 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    564 ! gives sol at perihelion for ls at perihelion (for precession cycles)
    565       subroutine lsp2solp(lsp,solp,aphelie,periheli,year_day)
    566 
    567       implicit none
    568 !  Arguments:
    569       real lsp     ! Input: ls at perihelion
    570       real solp    ! Output: sol at perihelion
    571       real aphelie,periheli,year_day ! Input: parameters
    572  
    573 !  Local:
    574       double precision zx0 ! eccentric anomaly at Ls=0
    575       double precision e_elips
    576       double precision pi,degrad
    577      
    578       parameter (pi=4.d0*atan(1.d0))
    579       parameter (degrad=180.d0/pi)
    580 
    581       e_elips=(aphelie-periheli)/(aphelie+periheli)     
    582       zx0 = -2.0*datan(dtan(0.5*lsp/degrad)
    583      .          *dsqrt((1.-e_elips)/(1.+e_elips)))
    584       if (zx0 .le. 0.) zx0 = zx0 + 2.*pi
    585      
    586       solp  = year_day*(1.-(zx0-e_elips*dsin(zx0))/(2.*pi))
    587 
    588 
    589       end subroutine lsp2solp
    590 
    591 
    592566      END MODULE tabfi_mod
    593 
Note: See TracChangeset for help on using the changeset viewer.