Changeset 3289 for trunk


Ignore:
Timestamp:
Apr 3, 2024, 11:49:58 AM (8 months ago)
Author:
jliu
Message:

21/03/2024 == Jliu

Fix several bugs in photochemistry related to the submit-crash-resubmit
problem: Some of the integers in the related routines are not set in
privatethreads. Consequently, the additives such as n=n+1 are accumulated with
threads, causing systermatic problems in the routine. This update fixed the
submit-crash-resubmit problem. The simulated results are same with previous
simulations as tested.

Location:
trunk/LMDZ.MARS
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.MARS/changelog.txt

    r3282 r3289  
    45834583when the first (at every time step) variable sent to wstats() was not part of
    45844584those listed in stats.def
     4585
     4586== 21/03/2024 == Jliu
     4587Fix several bugs in photochemistry related to the submit-crash-resubmit
     4588problem: Some of the integers in the related routines are not set in
     4589privatethreads. Consequently, the additives such as n=n+1 are accumulated with
     4590threads, causing systermatic problems in the routine. This update fixed the
     4591submit-crash-resubmit problem. The simulated results are same with previous
     4592simulations as tested.

  • trunk/LMDZ.MARS/libf/aeronomars/calchim_mod.F90

    r3012 r3289  
    3232      use comcstfi_h, only: pi
    3333      use chemistrydata_mod, only: read_phototable
    34       use photolysis_mod, only: init_photolysis, nphot
     34      use photolysis_mod, only: init_photolysis !, nphot
    3535      use iono_h, only: temp_elect
    3636      use wstats_mod, only: wstats
    … …  
    168168      integer :: ig_vl1
    169169
    170       integer :: nb_reaction_3_max     ! number of quadratic reactions
    171       integer :: nb_reaction_4_max     ! number of bimolecular reactions
    172       integer :: nquench               ! number of quenching + heterogeneous reactions
    173       integer :: nphotion              ! number of photoionizations
    174       integer :: nb_reaction_4_ion     ! quadratic reactions for ionosphere
    175       integer :: nb_reaction_4_deut    ! quadratic reactions for deuterium chem
    176       integer :: nb_phot_max           ! total number of photolysis+photoionizations+quenching reactions
    177 
     170      integer, save :: nb_reaction_3_max     ! number of quadratic reactions
     171      integer, save :: nb_reaction_4_max     ! number of bimolecular reactions
     172      integer, parameter :: nb_reaction_3_max0 = 6    ! number of quadratic reactions
     173      integer, parameter :: nb_reaction_4_max0 = 31   ! number of bimolecular reactions
     174      integer, parameter :: nquench =  9            ! number of quenching + heterogeneous reactions
     175      integer, save :: nphotion            ! number of photoionizations
     176      integer, parameter :: nphotion0 = 0            ! number of photoionizations
     177      integer, save :: nphot                       ! number of photolysis
     178      integer, parameter :: nphot0 = 13           ! number of photolysis
     179                                                  ! is incremented by +2 in calchim if deuterium chemisty)
     180      integer, parameter :: nb_reaction_4_ion = 64     ! quadratic reactions for ionosphere
     181      integer, parameter :: nb_reaction_4_deut = 35    ! quadratic reactions for deuterium chem
     182      integer, save :: nb_phot_max           ! total number of photolysis+photoionizations+quenching reactions
     183      integer, parameter :: nb_phot_max0 = 0 ! total number of photolysis+photoionizations+quenching reactions
     184!$OMP THREADPRIVATE(nb_reaction_3_max,nb_reaction_4_max,nphotion,nphot,nb_phot_max)
    178185
    179186      real    :: latvl1, lonvl1
    … …  
    757764         write(*,*) 'calchim: tracer indices=',niq(:)
    758765
    759 
     766         nphot = nphot0         ! number of photolysis
     767                                ! is incremented by +2 in calchim if deuterium chemisty)
    760768         if (photochem) then
    761769            if (jonline) then
    762770               print*,'calchim: Read UV absorption cross-sections'
    763771               !Add two photodissociations in deuterium chemistry included
    764                if(deutchem) nphot = nphot + 2
     772               if(deutchem) nphot = nphot0 + 2
    765773               call init_photolysis     ! on-line photolysis
    766774            else
    … …  
    843851         if (photochem) then
    844852            ! set number of reactions, depending on ion chemistry or not
    845             nb_reaction_4_ion  = 64
    846             nb_reaction_4_deut = 35
     853           ! nb_reaction_4_ion  = 64
     854           ! nb_reaction_4_deut = 35
    847855
    848856            !Default numbers if no ion and no deuterium chemistry included
    849857
    850             nb_reaction_4_max = 31     ! set number of bimolecular reactions
    851             nb_reaction_3_max = 6      ! set number of quadratic reactions
    852             nquench           = 9      ! set number of quenching + heterogeneous
    853             nphotion          = 0      ! set number of photoionizations
    854 
     858           ! nb_reaction_4_max = 31     ! set number of bimolecular reactions
     859           ! nb_reaction_3_max = 6      ! set number of quadratic reactions
     860           ! nquench           = 9      ! set number of quenching + heterogeneous
     861           ! nphotion          = 0      ! set number of photoionizations
     862            nb_reaction_4_max = nb_reaction_4_max0     ! set number of bimolecular reactions
     863            nb_reaction_3_max = nb_reaction_3_max0      ! set number of quadratic reactions
     864            nphotion          = nphotion0      ! set number of photoionizations
     865            nb_phot_max = nb_phot_max0
    855866            if (ionchem) then
    856                nb_reaction_4_max = nb_reaction_4_max + nb_reaction_4_ion 
    857                nphotion          = 18   ! set number of photoionizations
     867               nb_reaction_4_max = nb_reaction_4_max0 + nb_reaction_4_ion 
     868               nphotion          = nphotion0 + 18   ! set number of photoionizations
    858869            endif
    859870            if(deutchem) then
    … …  
    864875!        numerically as a photolysis:
    865876
    866             nb_phot_max = nphot + nphotion + nquench
    867 
     877            nb_phot_max = nb_phot_max0 + nphot + nphotion + nquench
     878!          print*, 'nb_phot_max = ', nb_phot_max
     879!         print*, 'nb_reaction_4_max=', nb_reaction_4_max
     880!          print*, 'nb_reaction_3_max=', nb_reaction_3_max
    868881!        call main photochemistry routine
    869882

  • trunk/LMDZ.MARS/libf/aeronomars/photochemistry.F90

    r3141 r3289  
    161161!integer,parameter :: i_hdo2    = 38
    162162
    163 integer :: i_last
    164 
     163integer, save :: i_last
     164!$OMP THREADPRIVATE(i_last)
    165165!Tracer indexes for photionization coeffs
    166166
    … …  
    180180
    181181integer :: ilay
    182 integer :: ind_norec
    183 integer :: ind_orec
    184 
     182integer, save :: ind_norec
     183integer, save :: ind_orec
     184!$OMP THREADPRIVATE(ind_norec, ind_orec)
    185185real :: ctimestep           ! standard timestep for the chemistry (s)
    186186real :: dt_guess            ! first-guess timestep (s)
    … …  
    503503em_no(:)=c(:,i_o)*c(:,i_n)*v_4(:,ind_norec)   !2.8e-17*(300./temp(:)))**0.5
    504504em_o2(:)=0.75*c(:,i_o)*c(:,i_o)*c(:,i_co2)*v_3(:,ind_orec)   !2.5*9.46e-34*exp(485./temp(:))*dens(:)
    505 
     505!print*, "ind_norec= ", ind_norec
     506!print*, "ind_orec= ", ind_orec
    506507!===================================================================
    507508!     save chemical species for the gcm       
    … …  
    659660
    660661use comcstfi_h
    661 use photolysis_mod, only : nphot
     662!use photolysis_mod, only : nphot
    662663
    663664implicit none
    … …  
    675676logical, intent(in)     :: deutchem
    676677integer                 :: lswitch           ! interface level between lower
     678integer,parameter :: nphot0 =13
    677679                                             ! atmosphere and thermosphere chemistries
    678680real, dimension(nlayer) :: dens              ! total number density (molecule.cm-3)
    … …  
    695697real (kind = 8), dimension(nlayer,nb_reaction_3_max) :: v_3
    696698real (kind = 8), dimension(nlayer,nb_reaction_4_max) :: v_4
    697 integer :: ind_norec
    698 integer :: ind_orec
     699integer, intent(out) :: ind_norec
     700integer, intent(out) :: ind_orec
    699701
    700702!----------------------------------------------------------------------
    … …  
    703705
    704706integer          :: ilev
    705 integer          :: nb_phot, nb_reaction_3, nb_reaction_4
     707integer, save    :: nb_phot, nb_reaction_3, nb_reaction_4
     708!$OMP THREADPRIVATE(nb_phot, nb_reaction_3, nb_reaction_4)
    706709real :: ak0, ak1, xpo, rate, rate1, rate2
    707710real :: k1a0, k1b0, k1ainf, k1a, k1b, fc, fx, x, y, gam
    … …  
    737740!     initialisation
    738741!----------------------------------------------------------------------
    739 
    740       nb_phot       = nphot + nphotion ! initialised to the number of photolysis + number of photoionization rates
     742      nb_phot       = 0
     743      nb_phot       = nphot0 + nphotion ! initialised to the number of photolysis + number of photoionization rates
    741744      nb_reaction_3 = 0
    742745      nb_reaction_4 = 0
    … …  
    776779      v_3(:,nb_reaction_3) = a002(:)
    777780      ind_orec=nb_reaction_3
    778 
     781!      print*,"ind_orec",ind_orec
    779782!---  a003: o + o3 -> o2 + o2
    780783
    … …  
    12001203      v_4(:,nb_reaction_4) = d007(:)
    12011204      ind_norec = nb_reaction_4
    1202 
     1205!      print*,"ind_norec=",ind_norec
    12031206!---  d008: n + ho2 -> no + oh
    12041207
Note: See TracChangeset for help on using the changeset viewer.