Ignore:
Timestamp:
Oct 22, 2015, 4:31:35 PM (9 years ago)
Author:
mturbet
Message:

Cleaning of condense_co2

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.GENERIC/libf/phystd/condense_co2.F90

    r1484 r1485  
    11      subroutine condense_co2(ngrid,nlayer,nq,ptimestep, &
    22          pcapcal,pplay,pplev,ptsrf,pt,                  &
    3           pphi,pdt,pdu,pdv,pdtsrf,pu,pv,pq,pdq,          &
    4           piceco2,pdqsurfc,albedo,pemisurf,              &
     3          pdt,pdtsrf,pq,pdq,                             &
     4          pqsurf,pdqsurfc,albedo,pemisurf,               &
    55          albedo_bareground,albedo_co2_ice_SPECTV,       &
    6           pdtc,pdtsrfc,pdpsrf,pduc,pdvc,                 &
    7           pdqc)
     6          pdtc,pdtsrfc,pdpsrfc,pdqc)
    87
    98      use radinc_h, only : L_NSPECTV, naerkind
    … …  
    2120!     Purpose
    2221!     -------
    23 !     Condense and/or sublime CO2 ice on the ground and in the
    24 !     atmosphere, and sediment the ice.
     22!     Condense and/or sublime CO2 ice on the ground and in the atmosphere, and sediment the ice.
    2523! 
    2624!     Inputs
    2725!     ------
    28 !     ngrid                 Number of vertical columns
    29 !     nlayer                Number of layers
    30 !     pplay(ngrid,nlayer)   Pressure layers
    31 !     pplev(ngrid,nlayer+1) Pressure levels
    32 !     pt(ngrid,nlayer)      Temperature (in K)
    33 !     ptsrf(ngrid)          Surface temperature
     26!     ngrid                 Number of vertical columns.
     27!     nlayer                Number of vertical layers.
     28!     nq                    Number of tracers.
     29!     ptimestep             Duration of the physical timestep (s).
     30!     pplay(ngrid,nlayer)   Pressure layers (Pa).
     31!     pplev(ngrid,nlayer+1) Pressure levels (Pa).
     32!     pt(ngrid,nlayer)      Atmospheric Temperatures (K).
     33!     ptsrf(ngrid)          Surface temperatures (K).
     34!     pq(ngrid,nlayer,nq)   Atmospheric tracers mixing ratios (kg/kg of air).
     35!     pqsurf(ngrid,nq)      Surface tracers (kg/m2).
    3436!     
    35 !     pdt(ngrid,nlayer)     Time derivative before condensation/sublimation of pt
    36 !     pdtsrf(ngrid)         Time derivative before condensation/sublimation of ptsrf
    37 !     pqsurf(ngrid,nq)      Sedimentation flux at the surface (kg.m-2.s-1)
     37!     pdt(ngrid,nlayer)     Time derivative before condensation/sublimation of pt.
     38!     pdtsrf(ngrid)         Time derivative before condensation/sublimation of ptsrf.
     39!     pdq(ngrid,nlayer,nq)  Time derivative before condensation/sublimation of
     40!
     41!     albedo_bareground(ngrid)           Albedo of the bare ground.
     42!     albedo_co2_ice_SPECTV(L_NSPECTV)   Spectral albedo of CO2 ice.
    3843!     
    3944!     Outputs
    4045!     -------
    41 !     pdpsrf(ngrid)       \  Contribution of condensation/sublimation
    42 !     pdtc(ngrid,nlayer)  /  to the time derivatives of Ps, pt, and ptsrf
    43 !     pdtsrfc(ngrid)     /
     46!     pdpsrfc(ngrid)          \       Contribution of condensation/sublimation
     47!     pdtc(ngrid,nlayer)       \            to the time derivatives of
     48!     pdtsrfc(ngrid)           /     Surface Pressure, Atmospheric Temperatures,
     49!     pdqsurfc(ngrid)         /         Surface Temperatures, Surface Tracers,
     50!     pdqc(ngrid,nlayer,nq)  /                and Atmospheric Tracers.*
     51!
     52!     pemisurf(ngrid)              Emissivity of the surface.
    4453!     
    4554!     Both
    4655!     ----
    47 !     piceco2(ngrid)               CO2 ice at the surface (kg/m2)
    48 !     albedo(ngrid,L_NSPECTV)      Spectral albedo at the surface
    49 !     pemisurf(ngrid)              Emissivity of the surface
     56!     albedo(ngrid,L_NSPECTV)      Spectral albedo of the surface.
    5057!
    5158!     Authors
    … …  
    5764!==================================================================
    5865
    59 !-----------------------------------------------------------------------
    60 !     Arguments
     66!--------------------------
     67!        Arguments
     68!--------------------------
     69
    6170
    6271      INTEGER,INTENT(IN) :: ngrid
    … …  
    6978      REAL,INTENT(IN) :: ptsrf(ngrid)
    7079      REAL,INTENT(IN) :: pt(ngrid,nlayer)
    71       REAL,INTENT(IN) :: pphi(ngrid,nlayer)
    7280      REAL,INTENT(IN) :: pdt(ngrid,nlayer)
    73       REAL,INTENT(IN) :: pdu(ngrid,nlayer)
    74       REAL,INTENT(IN) :: pdv(ngrid,nlayer)
    7581      REAL,INTENT(IN) :: pdtsrf(ngrid)
    76       REAL,INTENT(IN) :: pu(ngrid,nlayer)
    77       REAL,INTENT(IN) :: pv(ngrid,nlayer)
    7882      REAL,INTENT(IN) :: pq(ngrid,nlayer,nq)
     83      REAL,INTENT(IN) :: pqsurf(ngrid,nq)
    7984      REAL,INTENT(IN) :: pdq(ngrid,nlayer,nq)
    8085      REAL,INTENT(IN) :: albedo_bareground(ngrid)
    8186      REAL,INTENT(IN) :: albedo_co2_ice_SPECTV(L_NSPECTV)
    82       REAL,INTENT(INOUT) :: piceco2(ngrid)
    83       REAL,INTENT(OUT) :: albedo(ngrid,L_NSPECTV)
     87      REAL,INTENT(INOUT) :: albedo(ngrid,L_NSPECTV)
    8488      REAL,INTENT(OUT) :: pemisurf(ngrid)
    8589      REAL,INTENT(OUT) :: pdtc(ngrid,nlayer)
    8690      REAL,INTENT(OUT) :: pdtsrfc(ngrid)
    87       REAL,INTENT(OUT) :: pdpsrf(ngrid)
    88       REAL,INTENT(OUT) :: pduc(ngrid,nlayer)
    89       REAL,INTENT(OUT) :: pdvc(ngrid,nlayer)
     91      REAL,INTENT(OUT) :: pdpsrfc(ngrid)
    9092      REAL,INTENT(OUT) :: pdqc(ngrid,nlayer,nq)
    9193      REAL,INTENT(OUT) :: pdqsurfc(ngrid)
    9294
    93 !-----------------------------------------------------------------------
    94 !     Local variables
    95 
    96       INTEGER l,ig,icap,ilay,it,iq,nw
    97 
    98       REAL reffrad(ngrid,nlayer) ! radius (m) of the co2 ice particles
    99       REAL*8 zt(ngrid,nlayer)
    100       REAL zq(ngrid,nlayer,nq)
    101       REAL zcpi
    102       REAL ztcond (ngrid,nlayer)
    103       REAL ztnuc (ngrid,nlayer)
    104       REAL ztcondsol(ngrid)
    105       REAL zcondicea(ngrid,nlayer), zcondices(ngrid)
    106       REAL zfallice(ngrid), Mfallice(ngrid)
    107       REAL zmflux(nlayer+1)
    108       REAL zu(nlayer),zv(nlayer)
    109       REAL ztsrf(ngrid)
    110       REAL ztc(nlayer), ztm(nlayer+1)
    111       REAL zum(nlayer+1) , zvm(nlayer+1)
    112       LOGICAL condsub(ngrid)
    113       REAL subptimestep
    114       Integer Ntime
    115       real masse (ngrid,nlayer), w(ngrid,nlayer,nq)
    116       real wq(ngrid,nlayer+1)
    117       real vstokes,reff
    118 
    119 !     Special diagnostic variables
    120       real tconda1(ngrid,nlayer)
    121       real tconda2(ngrid,nlayer)
    122       real zdtsig (ngrid,nlayer)
    123       real zdt (ngrid,nlayer)
    124 
    125 !-----------------------------------------------------------------------
    126 !     Saved local variables
     95!------------------------------
     96!       Local variables
     97!------------------------------
     98
     99      INTEGER l,ig,icap,ilay,iq,nw,igas,it
     100
     101      REAL reffrad(ngrid,nlayer) ! Radius (m) of the CO2 ice particles.
     102      REAL*8 zt(ngrid,nlayer)    ! Updated Atmospheric Temperatures (K).
     103      REAL ztsrf(ngrid)          ! Updated Surface Temperatures (K).
     104      REAL zq(ngrid,nlayer,nq)   ! Updated Atmospheric tracers mixing ratios (kg/kg of air).
     105      REAL piceco2(ngrid)        ! Updated Surface Tracer (kg/m2).
     106      REAL ztcond (ngrid,nlayer) ! Atmospheric Temperatures of condensation of CO2.
     107      REAL ztnuc (ngrid,nlayer)  ! Atmospheric Nucleation Temperatures.
     108      REAL ztcondsol(ngrid)      ! Temperatures of condensation of CO2 at the surface.
     109      REAL zcondices(ngrid)      ! Condensation rate on the ground (kg/m2/s).
     110      REAL zfallice(ngrid)       ! Flux of ice falling on the surface (kg/m2/s).
     111      REAL Mfallice(ngrid)       ! Total amount of ice fallen to the ground during the timestep (kg/m2).
     112      REAL wq(ngrid,nlayer+1)    ! Total amount of ice fallen to the ground during the timestep (kg/m2).
     113      REAL subptimestep          ! Duration of the subtimestep (s) for the sedimentation.
     114      Integer Ntime              ! Number of subtimesteps.
     115      REAL masse (ngrid,nlayer)  ! Mass of atmospheric layers (kg/m2)
     116      REAL w(ngrid,nlayer,nq)    !
     117      REAL vstokes,reff          !
     118      REAL ppco2                 !
     119
     120
     121!------------------------------------------
     122!         Saved local variables
     123!------------------------------------------
     124
    127125
    128126      REAL,SAVE :: latcond=5.9e5
    129127      REAL,SAVE :: ccond
    130       REAL,SAVE :: cpice=1000.
    131128      REAL,SAVE,ALLOCATABLE,DIMENSION(:) :: emisref
    132 !$OMP THREADPRIVATE(latcond,ccond,cpice,emisref)
    133 
     129!$OMP THREADPRIVATE(latcond,ccond,emisref)
    134130      LOGICAL,SAVE :: firstcall=.true.
    135131!$OMP THREADPRIVATE(firstcall)
    136       REAL,EXTERNAL :: SSUM
    137 
    138       REAL,EXTERNAL :: CBRT
    139 
    140132      INTEGER,SAVE :: i_co2ice=0      ! co2 ice
    141133!$OMP THREADPRIVATE(i_co2ice)
    142134      CHARACTER(LEN=20) :: tracername ! to temporarily store text
    143135
    144       integer igas
    145 
    146       real ppco2
    147 
    148 !-----------------------------------------------------------------------
    149 !     Initializations
    150 
    151       pdqc(1:ngrid,1:nlayer,1:nq)=0
    152       pdtc(1:ngrid,1:nlayer)=0
    153       zq(1:ngrid,1:nlayer,1:nq)=0
    154       zt(1:ngrid,1:nlayer)=0
    155 
    156 !     Initialisation
     136
     137!------------------------------------------------
     138!       Initialization at the first call
     139!------------------------------------------------
     140
     141
    157142      IF (firstcall) THEN
    158143
    159          ALLOCATE(emisref(ngrid)) !! this should be deallocated in lastcall...
    160 
    161          ! find CO2 ice tracer
     144         ALLOCATE(emisref(ngrid))
     145         ! Find CO2 ice tracer.
    162146         do iq=1,nq
    163147            tracername=noms(iq)
    … …  
    167151         enddo
    168152         
    169         write(*,*) "condense_co2: i_co2ice=",i_co2ice       
    170 
    171         if((i_co2ice.lt.1))then
    172            print*,'In condens_cloud but no CO2 ice tracer, exiting.'
    173            print*,'Still need generalisation to arbitrary species!'
    174            stop
    175         endif
     153         write(*,*) "condense_co2: i_co2ice=",i_co2ice       
     154
     155         if((i_co2ice.lt.1))then
     156            print*,'In condens_cloud but no CO2 ice tracer, exiting.'
     157            print*,'Still need generalisation to arbitrary species!'
     158            stop
     159         endif
    176160
    177161         ccond=cpp/(g*latcond)
    178162         print*,'In condens_cloud: ccond=',ccond,' latcond=',latcond
    179163
    180 !          Prepare special treatment if gas is not pure CO2
    181              !if (addn2) then
    182              !   m_co2   = 44.01E-3 ! CO2 molecular mass (kg/mol)   
    183              !   m_noco2 = 28.02E-3 ! N2 molecular mass (kg/mol) 
    184 !               Compute A and B coefficient use to compute
    185 !               mean molecular mass Mair defined by
    186 !               1/Mair = q(ico2)/m_co2 + (1-q(ico2))/m_noco2
    187 !               1/Mair = A*q(ico2) + B
    188              !   A = (1/m_co2 - 1/m_noco2)
    189              !   B = 1/m_noco2
    190              !endif
    191 
    192 !          Minimum CO2 mixing ratio below which mixing occurs with layer above:
    193            !qco2min =0.75 
     164         ! Prepare special treatment if gas is not pure CO2
     165         ! if (addn2) then
     166         !    m_co2   = 44.01E-3 ! CO2 molecular mass (kg/mol)   
     167         !    m_noco2 = 28.02E-3 ! N2 molecular mass (kg/mol) 
     168         ! Compute A and B coefficient use to compute
     169         ! mean molecular mass Mair defined by
     170         ! 1/Mair = q(ico2)/m_co2 + (1-q(ico2))/m_noco2
     171         ! 1/Mair = A*q(ico2) + B
     172         ! A = (1/m_co2 - 1/m_noco2)
     173         ! B = 1/m_noco2
     174         ! endif
     175
     176         ! Minimum CO2 mixing ratio below which mixing occurs with layer above : qco2min =0.75 
    194177
    195178           firstcall=.false.
    196179      ENDIF
    197       zcpi=1./cpp
    198 
    199 !-----------------------------------------------------------------------
    200 !     Calculation of CO2 condensation / sublimation
    201 !
    202 !     Variables used:
    203 !     piceco2(ngrid)       amount of co2 ice on the ground  (kg/m2)
    204 !     zcondicea(ngrid,l)   condensation rate in layer l     (kg/m2/s)
    205 !     zcondices(ngrid)     condensation rate on the ground  (kg/m2/s)
    206 !     zfallice(ngrid)      flux of ice falling on surface   (kg/m2/s)
    207 !     pdtc(ngrid,nlayer) dT/dt due to phase changes       (K/s)
    208      
    209 
    210 !     Tendencies initially set to 0 (except pdtc)
    211       DO l=1,nlayer
    212          DO ig=1,ngrid
    213             zcondicea(ig,l) = 0.
    214             pduc(ig,l)  = 0
    215             pdvc(ig,l)  = 0
    216             pdqc(ig,l,i_co2ice)  = 0
    217          END DO
    218       END DO
    219      
    220       DO ig=1,ngrid
    221          Mfallice(ig) = 0.
    222          zfallice(ig) = 0.
    223          zcondices(ig) = 0.
    224          pdtsrfc(ig) = 0.
    225          pdpsrf(ig) = 0.
    226          condsub(ig) = .false.
    227          pdqsurfc(ig) = 0.
    228       ENDDO
    229 
    230 !-----------------------------------------------------------------------
     180
     181
     182!------------------------------------------------
     183!        Tendencies initially set to 0
     184!------------------------------------------------
     185
     186
     187      pdqc(1:ngrid,1:nlayer,1:nq)     = 0.
     188      pdtc(1:ngrid,1:nlayer)          = 0.
     189      zq(1:ngrid,1:nlayer,1:nq)       = 0.
     190      zt(1:ngrid,1:nlayer)            = 0.
     191      Mfallice(1:ngrid)               = 0.     
     192      zfallice(1:ngrid)               = 0.
     193      zcondices(1:ngrid)              = 0.
     194      pdtsrfc(1:ngrid)                = 0.
     195      pdpsrfc(1:ngrid)                = 0.
     196      pdqsurfc(1:ngrid)               = 0.
     197
     198
     199!----------------------------------
    231200!     Atmospheric condensation
     201!----------------------------------
    232202
    233203
    … …  
    238208!            DO ig=1,ngrid
    239209!              qco2=pq(ig,l,ico2)+pdq(ig,l,ico2)*ptimestep
    240 !!             Mean air molecular mass = 1/(q(ico2)/m_co2 + (1-q(ico2))/m_noco2)
     210!              Mean air molecular mass = 1/(q(ico2)/m_co2 + (1-q(ico2))/m_noco2)
    241211!              mmean=1/(A*qco2 +B)
    242212!              vmr_co2(ig,l) = qco2*mmean/m_co2
    … …  
    251221!      end if
    252222
    253 !     Forecast the atmospheric frost temperature 'ztcond' and nucleation temperature 'ztnuc'
     223
     224      ! Forecast the atmospheric frost temperature 'ztcond' and nucleation temperature 'ztnuc'.
    254225      DO l=1,nlayer
    255226         DO ig=1,ngrid
    … …  
    260231      ENDDO
    261232     
    262 !     Initialize zq and zt at the beginning of the sub-timestep loop
    263       DO l=1,nlayer
    264          DO ig=1,ngrid
     233      ! Initialize zq and zt at the beginning of the sub-timestep loop and qsurf.
     234      DO ig=1,ngrid
     235         piceco2(ig)=pqsurf(ig,i_co2ice)
     236         DO l=1,nlayer
    265237            zt(ig,l)=pt(ig,l)
    266238            zq(ig,l,i_co2ice)=pq(ig,l,i_co2ice)
    … …  
    274246      ENDDO
    275247
    276 !    Calculate the mass of each atmospheric layer (kg.m-2)
     248      ! Calculate the mass of each atmospheric layer (kg.m-2)
    277249      do  ilay=1,nlayer
    278250         DO ig=1,ngrid
    … …  
    281253      end do
    282254
    283 !     -----------------------------------------------
     255
     256!-----------------------------------------------------------
    284257!     START CONDENSATION/SEDIMENTATION SUB-TIME LOOP
    285 !     -----------------------------------------------
    286       Ntime =  20               ! number of sub-timestep
     258!-----------------------------------------------------------
     259
     260
     261      Ntime =  20  ! number of sub-timestep
    287262      subptimestep = ptimestep/float(Ntime)           
    288263
     264      ! Add the tendencies from other physical processes at each subtimstep.
    289265      DO it=1,Ntime
    290 
    291 !     Add the tendencies from other physical processes at each subtimstep
    292266         DO l=1,nlayer
    293267            DO ig=1,ngrid
    … …  
    297271         END DO
    298272
    299 
    300 !     Gravitational sedimentation
     273         ! Gravitational sedimentation starts.
    301274           
    302 !     sedimentation computed from radius computed from q in module radii_mod
     275         ! Sedimentation computed from radius computed from q in module radii_mod.
    303276         call co2_reffrad(ngrid,nlayer,nq,zq,reffrad)
    304277         
    305          do  ilay=1,nlayer
     278         DO  ilay=1,nlayer
    306279            DO ig=1,ngrid
    307280
    … …  
    316289                   pplev(ig,ilay)/(r*pt(ig,ilay))
    317290
    318             end do
    319          end do
    320 
    321 !     Computing q after sedimentation
    322 
     291            END DO
     292         END DO
     293
     294         ! Computing q after sedimentation
    323295         call vlz_fi(ngrid,nlayer,zq(1,1,i_co2ice),2.,masse,w(1,1,i_co2ice),wq)
    324296
    325297
    326 !     Progressively accumulating the flux to the ground
    327 !     Mfallice is the total amount of ice fallen to the ground
     298         ! Progressively accumulating the flux to the ground.
     299         ! Mfallice is the total amount of ice fallen to the ground.
    328300         DO ig=1,ngrid
    329301            Mfallice(ig) =  Mfallice(ig) + wq(ig,i_co2ice)
    330          end do
    331 
    332 
    333 !     Condensation / sublimation in the atmosphere
    334 !     --------------------------------------------
    335 !     (calculation of zcondicea, zfallice and pdtc)
    336 !     (MODIFICATIONS FOR EARLY MARS: falling heat neglected, condensation
    337 !     of CO2 into tracer i_co2ice)
     302         END DO
     303
     304!----------------------------------------------------------
     305!       Condensation / sublimation in the atmosphere
     306!----------------------------------------------------------
     307!     (MODIFICATIONS FOR EARLY MARS: falling heat neglected, condensation of CO2 into tracer i_co2ice)
     308
    338309
    339310         DO l=nlayer , 1, -1
    … …  
    341312               pdtc(ig,l)=0.
    342313
    343 
    344    !     ztcond-> ztnuc in test beneath to nucleate only when super saturation occurs(JL 2011)
     314               ! ztcond-> ztnuc in test beneath to nucleate only when super saturation occurs(JL 2011)
    345315               IF ((zt(ig,l).LT.ztnuc(ig,l)).or.(zq(ig,l,i_co2ice).gt.1.E-10)) THEN
    346                   condsub(ig)=.true.
    347316                  pdtc(ig,l)   = (ztcond(ig,l) - zt(ig,l))/subptimestep
    348317                  pdqc(ig,l,i_co2ice) = pdtc(ig,l)*ccond*g
    349318
    350 !    Case when the ice from above sublimes entirely
     319                  ! Case when the ice from above sublimes entirely
    351320                  IF ((zq(ig,l,i_co2ice).lt.-pdqc(ig,l,i_co2ice)*subptimestep) &
    352321                      .AND. (zq(ig,l,i_co2ice).gt.0)) THEN
    … …  
    357326                  END IF
    358327
    359 !    Temperature and q after condensation
     328                  ! Temperature and q after condensation
    360329                  zt(ig,l)   = zt(ig,l)   + pdtc(ig,l)   * subptimestep
    361330                  zq(ig,l,i_co2ice) = zq(ig,l,i_co2ice) + pdqc(ig,l,i_co2ice) * subptimestep
    … …  
    364333            ENDDO
    365334         ENDDO
    366       ENDDO                     ! end of subtimestep loop
    367 
    368 !     Computing global tendencies after the subtimestep
     335         
     336      ENDDO! end of subtimestep loop.
     337
     338      ! Computing global tendencies after the subtimestep.
    369339      DO l=1,nlayer
    370340         DO ig=1,ngrid
    … …  
    381351
    382352!-----------------------------------------------------------------------
    383 !     Condensation/sublimation on the ground
    384 !     (calculation of zcondices and pdtsrfc)
    385 
    386 !     forecast of ground temperature ztsrf and frost temperature ztcondsol
     353!              Condensation/sublimation on the ground
     354!-----------------------------------------------------------------------
     355
     356
     357      ! Forecast of ground temperature ztsrf and frost temperature ztcondsol.
    387358      DO ig=1,ngrid
    388359         ppco2=gfrac(igas_CO2)*pplay(ig,1)
    … …  
    406377     
    407378      DO ig=1,ngrid
     379     
    408380         IF(ig.GT.ngrid/2+1) THEN
    409381            icap=2
    … …  
    412384         ENDIF
    413385
    414 !    Loop over where we have condensation / sublimation
     386         ! Loop over where we have condensation / sublimation
    415387         IF ((ztsrf(ig) .LT. ztcondsol(ig)) .OR.           &        ! ground condensation
    416388                    (zfallice(ig).NE.0.) .OR.              &        ! falling snow
    417389                    ((ztsrf(ig) .GT. ztcondsol(ig)) .AND.  &        ! ground sublimation
    418390                    ((piceco2(ig)+zfallice(ig)*ptimestep) .NE. 0.))) THEN
    419             condsub(ig) = .true.
    420 
    421 !    Condensation or partial sublimation of CO2 ice
     391
     392
     393            ! Condensation or partial sublimation of CO2 ice
    422394            zcondices(ig)=pcapcal(ig)*(ztcondsol(ig)-ztsrf(ig)) &
    423395                /(latcond*ptimestep)         
    424396            pdtsrfc(ig) = (ztcondsol(ig) - ztsrf(ig))/ptimestep
    425397
    426 !    If the entire CO_2 ice layer sublimes
    427 !    (including what has just condensed in the atmosphere)
     398            ! If the entire CO_2 ice layer sublimes
     399            ! (including what has just condensed in the atmosphere)
    428400            IF((piceco2(ig)/ptimestep+zfallice(ig)).LE. &
    429401                -zcondices(ig))THEN
    … …  
    433405            END IF
    434406
    435 !     Changing CO2 ice amount and pressure
    436 
    437             pdqsurfc(ig) = zcondices(ig) + zfallice(ig)
    438             piceco2(ig)  = piceco2(ig) + pdqsurfc(ig)*ptimestep
    439             pdpsrf(ig)   = -pdqsurfc(ig)*g
    440 
    441             IF(ABS(pdpsrf(ig)*ptimestep).GT.pplev(ig,1)) THEN
     407            ! Changing CO2 ice amount and pressure
     408            piceco2(ig)  =  piceco2(ig)   + pdqsurfc(ig)*ptimestep
     409            pdqsurfc(ig) =  zcondices(ig) + zfallice(ig)
     410            pdpsrfc(ig)  = -pdqsurfc(ig)*g
     411
     412            IF(ABS(pdpsrfc(ig)*ptimestep).GT.pplev(ig,1)) THEN
    442413               PRINT*,'STOP in condens'
    443414               PRINT*,'condensing more than total mass'
    444415               PRINT*,'Grid point ',ig
    445416               PRINT*,'Ps = ',pplev(ig,1)
    446                PRINT*,'d Ps = ',pdpsrf(ig)
     417               PRINT*,'d Ps = ',pdpsrfc(ig)
    447418               STOP
    448419            ENDIF
    449420         END IF
    450       ENDDO
    451 
     421         
     422      ENDDO ! end of ngrid loop.
     423
     424
     425!---------------------------------------------------------------------------------------------
    452426!     Surface albedo and emissivity of the ground below the snow (emisref)
    453 !     --------------------------------------------------------------------
     427!---------------------------------------------------------------------------------------------
     428
     429
    454430      DO ig=1,ngrid
     431     
    455432         IF(lati(ig).LT.0.) THEN
    456433            icap=2 ! Southern Hemisphere
    … …  
    477454         end if
    478455         
    479          piceco2(ig)  = piceco2(ig) - pdqsurfc(ig)*ptimestep ! This line was added so that tendencies are added outside the routine. MT2015.
    480          if(.not.piceco2(ig).ge.0.) THEN
    481             if(piceco2(ig).le.-1.e-8) print*,   &
    482               'WARNING 2 : in condense_co2cloud: piceco2(',ig,')=', piceco2(ig)
    483             piceco2(ig)=0.
    484          endif
    485          
    486       end do
     456      END DO
    487457
    488458      return
    489     end subroutine condense_co2
    490 
    491 !-------------------------------------------------------------------------
    492     subroutine get_tcond_co2(p,tcond)
    493 !   Calculates the condensation temperature for CO2
     459     
     460      end subroutine condense_co2
     461
     462
     463
     464!-------------------------------------------------------------------------
     465!-------------------------------------------------------------------------
     466!-------------------------------------------------------------------------
     467!-------------------------------------------------------------------------
     468!-------------------------------------------------------------------------
     469!-------------------------------------------------------------------------
     470
     471
     472
     473      subroutine get_tcond_co2(p,tcond)  ! Calculates the condensation temperature for CO2
     474     
    494475
    495476      implicit none
    … …  
    500481      peff=p
    501482
    502       if(peff.lt.ptriple)then
    503          tcond = (-3167.8)/(log(.01*peff)-23.23) ! Fanale's formula
     483      if(peff.lt.ptriple) then
     484         tcond = (-3167.8)/(log(.01*peff)-23.23) ! Fanale's formula.
    504485      else
    505          tcond = 684.2-92.3*log(peff)+4.32*log(peff)**2
    506          ! liquid-vapour transition (based on CRC handbook 2003 data)
     486         tcond = 684.2-92.3*log(peff)+4.32*log(peff)**2 ! liquid-vapour transition (based on CRC handbook 2003 data)       
    507487      endif
    508488      return
    509489
    510     end subroutine get_tcond_co2
    511    
    512    
    513    
    514  
    515 !-------------------------------------------------------------------------
    516     subroutine get_tnuc_co2(p,tnuc)
    517 !   Calculates the nucleation temperature for CO2, based on a simple super saturation criterion
    518 !     (JL 2011)
     490      end subroutine get_tcond_co2
     491
     492
     493
     494!-------------------------------------------------------------------------
     495!-------------------------------------------------------------------------
     496!-------------------------------------------------------------------------
     497!-------------------------------------------------------------------------
     498!-------------------------------------------------------------------------
     499!-------------------------------------------------------------------------
     500
     501
     502
     503      subroutine get_tnuc_co2(p,tnuc)
     504      ! Calculates the nucleation temperature for CO2, based on a simple super saturation criterion. JL 2011.
    519505
    520506      use callkeys_mod, only: co2supsat
    … …  
    527513      peff=p/co2supsat
    528514
    529       if(peff.lt.ptriple)then
     515      if(peff.lt.ptriple) then
    530516         tnuc = (-3167.8)/(log(.01*peff)-23.23) ! Fanale's formula
    531517      else
    … …  
    533519         ! liquid-vapour transition (based on CRC handbook 2003 data)
    534520      endif
     521     
    535522      return
    536523
    537     end subroutine get_tnuc_co2
     524      end subroutine get_tnuc_co2
Note: See TracChangeset for help on using the changeset viewer.