Changeset 1390


Ignore:
Timestamp:
Mar 5, 2015, 12:44:03 PM (10 years ago)
Author:
flefevre
Message:
  • changement de la composition via newstart:
    • amelioration pour eviter des valeurs negatives de CO2 dans la thermosphere
    • mise a jour des valeurs d'exemples mesurees par MSL
File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.MARS/libf/phymars/newstart.F

    r1383 r1390  
    3535      use comgeomfi_h, only: ini_fillgeom
    3636
     37
    3738      implicit none
    3839
     
    168169      real :: Mair_old,Mair_new,vmr_old,vmr_new
    169170      real,allocatable :: coefvmr(:)  ! Correction coefficient when changing composition
     171      integer :: iloc(1), iqmax
    170172
    171173c sortie visu pour les champs dynamiques
     
    10211023
    10221024              if (iq.eq.igcm_n2) then
    1023                 write(*,*) "New vmr(n2)? (MSL: 1.89E-02 at Ls~186)"
     1025                write(*,*) "New vmr(n2)? (MSL: 2.03e-02 at Ls~184)"
    10241026              endif
    10251027              if (iq.eq.igcm_ar) then
    1026                 write(*,*) "New vmr(ar)? (MSL: 1.93E-02 at Ls~186)"
     1028                write(*,*) "New vmr(ar)? (MSL: 2.07e-02 at Ls~184)"
    10271029              endif
    10281030              if (iq.eq.igcm_o2) then
    1029                 write(*,*) "New vmr(o2)? (MSL: 1.46E-03) at Ls~186"
     1031                write(*,*) "New vmr(o2)? (MSL: 1.73e-03 at Ls~184)"
    10301032              endif
    10311033              if (iq.eq.igcm_co) then
    1032                 write(*,*) "New vmr(co)? (~8.E-04)"
     1034                write(*,*) "New vmr(co)? (MSL: 7.49e-04 at Ls~184)"
    10331035              endif
    10341036 302          read(*,*,iostat=ierr) vmr_new
     
    10551057          end do
    10561058
    1057         ! Recompute mass mixing ratios everywhere, and adjust mmr(CO2) to keep sum of mmr constant.
     1059        ! Recompute mass mixing ratios everywhere, and adjust mmr of most abundant species
     1060        ! to keep sum of mmr constant.
    10581061          do l=1,llm
    10591062           do j=1,jjp1
     
    10691072                end if
    10701073              enddo
    1071               q(i,j,l,igcm_co2) = q(i,j,l,igcm_co2) + Smmr_old-Smmr_new
     1074              iloc = maxloc(q(i,j,l,:))
     1075              iqmax = iloc(1)
     1076              q(i,j,l,iqmax) = q(i,j,l,iqmax) + Smmr_old - Smmr_new
    10721077            enddo
    10731078           enddo
     
    10751080
    10761081          write(*,*)
    1077      &   "mmr(CO2) is modified everywhere to keep sum of mmr constant"
     1082     &   "The most abundant species is modified everywhere to keep "//
     1083     &   "sum of mmr constant"
    10781084          write(*,*) 'At reference site vmr(CO2)=',
    10791085     &        q(iref,jref,1,igcm_co2)*Mair_new/mmol(igcm_co2)
    1080           write(*,*) "Compared to MSL observation: vmr(CO2)= 0.9597"
    1081 
     1086          write(*,*) "Compared to MSL observation: vmr(CO2)= 0.957 "//
     1087     &   "at Ls=184"
    10821088
    10831089c      wetstart : wet atmosphere with a north to south gradient
Note: See TracChangeset for help on using the changeset viewer.