|
|
(11 intermediate revisions by the same user not shown) |
Line 13: |
Line 13: |
| == Relaxed Arakawa-Schubert scheme == | | == Relaxed Arakawa-Schubert scheme == |
|
| |
|
| This is the convection scheme that GEOS-Chem uses with both GEOS-3 and GEOS-5 meteorology. The source code for this scheme is located in routine <tt>NFCLDMX</tt> in <tt>convection_mod.f</tt>. | | This is the convection scheme that GEOS-Chem currently uses with [[GEOS-FP]] and [[MERRA-2]] meteorology. The source code for this scheme is located in routine <tt>DO_CONVECTION</tt> in <tt>GeosCore/convection_mod.F</tt>. |
| | |
| == Hack and Zhang-McFarlane schemes ==
| |
| | |
| This is the convection scheme that GEOS-Chem uses with both GEOS-4 and GCAP meteorology.
| |
| | |
| === GEOS-4 ===
| |
| | |
| The source code for the convection scheme required for GEOS-4 meteorology is contained within the F90 module <tt>fvdas_convect_mod.f</tt>. The main driver is called <tt>FVDAS_CONVECT</tt>. This routine calls the following routines:
| |
| | |
| #<tt>ARCONVTRAN</tt>: Prepares arrays for use by <tt>CONVTRAN</tt>
| |
| #<tt>CONVTRAN</tt>: Main driver for deep convection (i.e. Zhang/McFarlane scheme)
| |
| #<tt>HACK_CONV</tt>: Main driver for shallow convection (i.e. Hack scheme)
| |
| | |
| === GISS ===
| |
| | |
| The source code for the convection scheme required for GISS/GCAP meteorology is contained within the F90 module <tt>gcap_convect_mod.f</tt>. The main driver is called <tt>GCAP_CONVECT</tt>. This routine calls the following routines:
| |
| | |
| #<tt>ARCONVTRAN</tt>: Prepares arrays for use by <tt>CONVTRAN</tt>
| |
| #<tt>CONVTRAN</tt>: Main driver for deep convection (i.e. Zhang/McFarlane scheme)
| |
| | |
| There is no equivalent routine to <tt>HACK_CONV</tt> for the GCAP meterology. Instead, <tt>ARCONVTRAN</tt> and <tt>CONVTRAN</tt> are called with entraining mass fluxes, and then again with non-entraining mass fluxes.
| |
| | |
| --[[User:Bmy|Bob Y.]] 13:11, 19 February 2010 (EST)
| |
|
| |
|
| == Validation == | | == Validation == |
Line 56: |
Line 33: |
|
| |
|
| --[[User:Bmy|Bob Y.]] 13:11, 19 February 2010 (EST) | | --[[User:Bmy|Bob Y.]] 13:11, 19 February 2010 (EST) |
|
| |
| == Previous issues that are now resolved ==
| |
|
| |
| === Resolve very high tracer concentrations in MERRA and GEOS-FP convective scavenging ===
| |
|
| |
| <span style="color:green">'''''This update was validated with [[GEOS-Chem_v11-01_benchmark_history#v11-01d|1-month benchmark simulation v11-01d]] and [[GEOS-Chem_v11-01_benchmark_history#v11-01d-Run1|1-year benchmark simulation v11-01d-Run1]]. This version was approved on 12 Dec 2015.'''''</span>
| |
|
| |
| '''''Viral Shah wrote:'''''
| |
|
| |
| :I am running [[GEOS-Chem v11-01|v10.01]], nested-grid NA full-chemistry w/ SOA simulation with the GEOS-FP met fields, and have found instances of extremely high concentrations of certain tracers that develop all of a sudden. I traced this back to a bug (highlighted in <span style="color:red">RED</span>) in the <tt>DO_MERRA_CONVECTION</tt> (in module <tt>GeosCore/convection_mod.F</tt>) routine at the following:
| |
|
| |
| ! Check if the tracer is an aerosol or not
| |
| IF ( AER == .TRUE. ) THEN
| |
|
| |
| !---------------------------------------------------------
| |
| ! Washout of aerosol tracers
| |
| ! This is modeled as a kinetic process
| |
| !---------------------------------------------------------
| |
|
| |
| <span style="color:red">! Define ALPHA, the fraction of raindrops that
| |
| ! re-evaporate when falling from (I,J,L+1) to (I,J,L)
| |
| ALPHA = ( REEVAPCN(K) * AD(K) )
| |
| & / ( PDOWN(K+1) * AREA_M2 * 10e+0_fp )
| |
|
| |
| ! ALPHA2 is the fraction of the rained-out aerosols
| |
| ! that gets resuspended in grid box (I,J,L)
| |
| ALPHA2 = 0.5e+0_fp * ALPHA</span>
| |
|
| |
| :Here, the value of <tt>ALPHA</tt>, which should be less than or equal to 1, is sometimes much higher than 1, and for the attached example was >10^5. Since <tt>ALPHA</tt> is used to calculate resuspension of precipitating mass from above, these high values of <tt>ALPHA</tt> lead to a gain in tracer mass at the lower levels in excess of what is coming from above. I suggest we replace the <span style="color:red">RED</span> code with the following piece of code in <span style="color:green">GREEN</span> which is used elsewhere in the <tt>DO_MERRA_CONVECTION</tt> routine.
| |
|
| |
| ! %%%% CASE 1 %%%%
| |
| ! Partial re-evaporation. Less precip is leaving
| |
| ! the grid box then entered from above.
| |
| IF ( PDOWN(K+1) > PDOWN(K) .AND.
| |
| & PDOWN(K) > TINYNUM ) THEN
| |
|
| |
| <span style="color:green">! Define ALPHA, the fraction of raindrops that
| |
| ! re-evaporate when falling from grid box
| |
| ! (I,J,L+1) to (I,J,L)
| |
| ALPHA = ( REEVAPCN(K) * AD(K) )
| |
| & / ( PDOWN(K+1) * AREA_M2 * 10e+0_fp )
| |
|
| |
| ! For safety
| |
| ALPHA = MIN( ALPHA, 1e+0_fp )
| |
|
| |
| ! ALPHA2 is the fraction of the rained-out aerosols
| |
| ! that gets resuspended in grid box (I,J,L)
| |
| ALPHA2 = 0.5e+0_fp * ALPHA</span>
| |
|
| |
| ENDIF
| |
|
| |
| ! %%%% CASE 2 %%%%
| |
| ! Total re-evaporation. Precip entered from above,
| |
| ! but no precip is leaving grid box (ALPHA = 2 so
| |
| ! that ALPHA2 = 1)
| |
| IF ( PDOWN(K) < TINYNUM ) THEN
| |
| ALPHA2 = 1e+0_fp
| |
| ENDIF
| |
|
| |
| :A second bug I noticed is in the lines in <span style="color:red">RED</span> right below the first one. The current code is as follows:
| |
|
| |
| ! GAINED is the rained out aerosol coming down from
| |
| ! grid box (I,J,L+1) that will evaporate and re-enter
| |
| ! the atmosphere in the gas phase in grid box (I,J,L).
| |
| GAINED = T0_SUM * ALPHA2
| |
|
| |
| ! Amount of aerosol lost to washout in grid box
| |
| <span style="color:red">WETLOSS = Q(K,IC) * BMASS(K) / TCVV *
| |
| & WASHFRAC - GAINED</span>
| |
|
| |
| ! LOST is the rained out aerosol coming down from
| |
| ! grid box (I,J,L+1) that will remain in the liquid
| |
| ! phase in grid box (I,J,L) and will NOT re-evaporate.
| |
| LOST = T0_SUM - GAINED
| |
|
| |
| ! Update tracer concentration (V. Shah, mps, 5/20/15)
| |
| Q(K,IC) = Q(K,IC) - WETLOSS * TCVV / BMASS(K)
| |
|
| |
| :In this case, the amount of tracer <tt>GAINED</tt> is (mistakenly) not added to <tt>Q</tt> before <tt>WETLOSS</tt> is calculated. The fix for the bold faced line should be as follows:
| |
|
| |
| <span style="color:green">WETLOSS = ( Q(K,IC) * BMASS(K) / TCVV + GAINED ) *
| |
| & WASHFRAC - GAINED</span>
| |
|
| |
| :This will be consistent with the code (below) for the 'non-aerosol' case
| |
|
| |
| ! MASS_WASH is the total amount of non-aerosol tracer
| |
| ! that is available for washout in grid box (I,J,L).
| |
| ! It consists of the mass in the precipitating
| |
| ! part of box (I,J,L), plus the previously rained-out
| |
| ! tracer coming down from grid box (I,J,L+1).
| |
| ! (Eq. 15, Jacob et al, 2000).
| |
| MASS_WASH = ( F_WASHOUT * Q(K,IC) ) * BMASS(K) /
| |
| & TCVV + T0_SUM
| |
|
| |
| ! WETLOSS is the amount of tracer mass in
| |
| ! grid box (I,J,L) that is lost to washout.
| |
| ! (Eq. 16, Jacob et al, 2000)
| |
| WETLOSS = MASS_WASH * WASHFRAC - T0_SUM
| |
|
| |
| --[[User:Melissa Payer|Melissa Sulprizio]] ([[User talk:Melissa Payer|talk]]) 16:43, 14 September 2015 (UTC)
| |
|
| |
| === Fixed bug in DO_MERRA_CONVECTION affecting mass conservation ===
| |
|
| |
| <span style="color:green">'''''These updates were validated in the 1-month benchmark simulation [[GEOS-Chem v11-01 benchmark history#v11-01a|v11-01a]], which was approved on 07 Jul 2015.'''''</span>
| |
|
| |
| Mass conservation tests validating the change in tracer units within cloud convection revealed a discrepancy in mass conservation between using GEOS-5 and GEOS-FP. Comparison of subroutine <tt>DO_MERRA_CONVECTION</tt>, the convection routine used for GEOS-FP, and subroutine <tt>NFCLDMX</tt>, the convection routine used for GEOS-5, showed that a calculation in <tt>DO_MERRA_CONVECTION</tt> was incorrectly outside of the internal timestep loop.
| |
|
| |
| In <tt>DO_MERRA_CONVECTION</tt>, the quantity <tt>QB</tt>, the average mixing ratio in the column below the cloud base, is calculated in two steps:
| |
|
| |
| #the numerator of <tt>QB</tt> is calculated by summing the tracer concentrations weighted by pressure and assigned to variable <tt>QB_NUM</tt>
| |
| #<tt>QB</tt> is calculated by dividing <tt>QB_NUM</tt> by the pressure difference across all levels below the cloud base.
| |
|
| |
| While tracer concentrations were updated every internal timestep, <tt>QB_NUM</tt> was calculated prior to the timestep loop and therefore was not updated with the adjusted concentrations. As a result, mass was only conserved during the first internal timestep.
| |
|
| |
| To correct this problem, we made these modifications in subroutine <tt>DO_MERRA_CONVECTION</tt>, located in module <tt>GeosCore/convection_mod.F</tt>:
| |
|
| |
| #Change the declaration of <tt>QB_NUM</tt> from an array to a scalar.
| |
| #Move the initialization and calculation of <tt>QB_NUM</tt> to within the timestep loop just prior to the calculation of <tt>QB</tt>. Change <tt>QB_NUM</tt> from an array to a scalar as follows:
| |
| ! Calculate QB_NUM, the numerator for QB. QB is the
| |
| ! weighted average mixing ratio below the cloud base.
| |
| ! QB_NUM is equal to the grid box tracer concentrations
| |
| ! [kg/kg total air] weighted by the adjacent level pressure
| |
| ! differences and summed over all levels up to just
| |
| ! below the cloud base (ewl, 6/22/15)
| |
| QB_NUM = 0e+0_fp
| |
| DO K = 1, CLDBASE-1
| |
| QB_NUM = QB_NUM + Q(K,IC) * DELP(K)
| |
| ENDDO
| |
|
| |
| ! Compute QB, the weighted avg mixing ratio below
| |
| ! the cloud base [kg/kg total air]
| |
| QB = QB_NUM / ( PEDGE(1) - PEDGE(CLDBASE) )
| |
|
| |
| With these updates, the mass variance for GEOS-FP within convection is the same order of magnitude as for GEOS-5. For the case of CO2 with emissions turned off and initial values of 370 ppm globally, the mass fluctuation is a few hundred kg globally over a few days.
| |
|
| |
| --[[User:Lizzie Lundgren|Lizzie Lundgren]] ([[User talk:Lizzie Lundgren|talk]]) 21:50, 24 June 2015 (UTC)<br>--[[User:Bmy|Bob Y.]] ([[User talk:Bmy|talk]]) 17:55, 7 July 2015 (UTC)
| |
|
| |
| === Reduce time spent in routine DO_CONVECTION when using GEOS-FP or MERRA ===
| |
|
| |
| <span style="color:green">'''''These updates were validated in the 1-month benchmark simulation [[GEOS-Chem v10-01 benchmark history#v10-01c|v10-01c]] and approved on 29 May 2014.'''''</span>
| |
|
| |
| Using the [http://tau.uoregon.edu Tuning and Analysis Utilities (TAU)], we were able to make the following modifications to routine <tt>DO_CONVECTION</tt> (in module file <tt>GeosCore/convection_mod.F</tt>) that resulted in a speedup in the [[GEOS-FP]]/[[MERRA]] convection module. We made modifications in the following routines:
| |
|
| |
| #<tt>DO_CONVECTION</tt> (in module <tt>GeosCore/convection_mod.F</tt>)
| |
| #*We no longer use objects from <tt>GeosCore/gc_type_mod.F</tt> (i.e. <tt>OPTIONS</tt>, <tt>DIMINFO</tt>, etc.). These have been largely superseded by the <tt>Input_Opt</tt> object.
| |
| #*Remove all pointer references to fields of <tt>State_Met</tt> and <tt>State_Chm</tt> from the parallel DO loop.
| |
| #*Now use the proper # of tracers for APM aerosol microphysics.
| |
| #*Now declare the parallel DO loop with <tt>!$OMP+PRIVATE</tt>, which results in better load-balancing.
| |
| #<tt>DO_MERRA_CONVECTION</tt> (in module <tt>GeosCore/convection_mod.F</tt>):
| |
| #*Now make the following variables local pointer variables instead of arguments:
| |
| #**<tt>AD, BXHEIGHT, CMFMC, DQRCU, DTRAIN, REEVAPCN, T, Q</tt>. These point to the corresponding variables in <tt>State_Met</tt> or <tt>State_Chm</tt>.
| |
| #**<tt>H2O2s, SO2s</tt>: These point to the corresponding arrays in <tt>GeosCore/sulfate_mod.F</tt>.
| |
| #<tt>COMPUTE_F</tt> (in module <tt>GeosCore/wetscav_mod.F</tt>)
| |
| #*Now make the F argument of routine COMPUTE_F an assumed-shape array (i.e. <tt>F(:,:,:)</tt> ). This allows us to pass pointer-valued arguments to
| |
| #<tt>GeosCore/gc_type_mod.F</tt>:
| |
| #*This has now been moved to the <tt>obsolete/</tt> subdirectory.
| |
|
| |
| [[Image:GEOSFP_Conv_Speedup.png]]
| |
|
| |
| Also note that with this fix, the total CPU time of this run (1-day with GEOS-4x5 met) dropped from 16:l1 to 15:23 (min:sec). So this fix is making the overall simulation run a little faster for each simulation day.
| |
|
| |
| --[[User:Bmy|Bob Y.]] 17:09, 30 May 2014 (EDT)
| |
|
| |
| === Removed array temporaries in call to FVDAS_CONVECT ===
| |
|
| |
| <span style="color:green">'''''These updates were validated in the 1-month benchmark simulation [[GEOS-Chem v10-01 benchmark history#v10-01c|v10-01c]] and approved on 29 May 2014.'''''</span>
| |
|
| |
| The [[GEOS-Chem Unit Tester]] revealed the presence of array temporaries in the call to GEOS-4 convection routine <tt>FVDAS_CONVECT</tt>. To remove the array temporaries, we made these modifications in subroutine <tt>DO_GEOS4_CONVECT</tt>, located in module <tt>GeosCore/convection_mod.F</tt>:
| |
|
| |
| (1) Declare the <tt>F</tt> variable with the <tt>TARGET</tt> attribute, as follows:
| |
|
| |
| #if defined( APM )
| |
| INTEGER :: INDEXSOL(Input_Opt%N_TRACERS+N_APMTRA)
| |
| REAL*8, TARGET :: F (IIPAR,JJPAR,LLPAR,
| |
| & Input_Opt%N_TRACERS+N_APMTRA)
| |
| #else
| |
| INTEGER :: INDEXSOL(Input_Opt%N_TRACERS)
| |
| REAL*8, TARGET :: F (IIPAR,JJPAR,LLPAR,
| |
| & Input_Opt%N_TRACERS)
| |
| #endif
| |
|
| |
| (2) Add pointer variables <tt>p_STT</tt> and <tt>p_F</tt>:
| |
|
| |
| REAL*8, POINTER :: p_STT (:,:,:,:) ! Ptr to STT (flipped in vert)
| |
| REAL*8, POINTER :: p_F (:,:,:,:) ! Ptr to F (flipped in vert)
| |
|
| |
| (3) Point <tt>p_STT</tt> to the <tt>State_Chm%Tracers</tt> field and flip levels in the vertical. Likewise, point <tt>p_F</tt> to the <tt>F</tt> array, and also flip in the vertical:
| |
|
| |
| ! Point to various fields and flip in vertical. This avoids
| |
| ! array temporaries in the call to FVDAS_CONVECT (bmy, 4/14/14)
| |
| p_HKETA => State_Met%HKETA (:,:,LLPAR:1:-1 )
| |
| ... etc ...
| |
| p_STT => State_Chm%Tracers(:,:,LLPAR:1:-1,:)
| |
| p_F => F (:,:,LLPAR:1:-1,:)
| |
|
| |
| (4) Pass <tt>p_STT</tt> and <tt>p_F</tt> to subroutine <tt>FVDAS_CONVECT</tt>:
| |
|
| |
| ! Call the fvDAS convection routines (originally from NCAR!)
| |
| CALL FVDAS_CONVECT( TDT,
| |
| & N_TOT_TRC,
| |
| !----------------------------------------------------------------------------
| |
| ! Prior to 4/14/14:
| |
| ! Eliminate an array temporary (bmy, 4/14/14)
| |
| ! & STT (:,:,LLPAR:1:-1,:),
| |
| !----------------------------------------------------------------------------
| |
| & p_STT,
| |
| & RPDEL,
| |
| & p_HKETA,
| |
| & p_HKBETA,
| |
| & p_ZMMU,
| |
| & p_ZMMD,
| |
| & p_ZMEU,
| |
| & DP,
| |
| & NSTEP,
| |
| !----------------------------------------------------------------------------
| |
| ! Prior to 4/14/14:
| |
| ! Eliminate an array temporary (bmy, 4/14/14)
| |
| ! & F (:,:,LLPAR:1:-1,:),
| |
| !----------------------------------------------------------------------------
| |
| & p_F,
| |
| & TCVV,
| |
| & INDEXSOL,
| |
| & State_Met )
| |
|
| |
| (5) Nullify the pointers <tt>p_STT</tt> and <tt>p_F</tt> at the end of subroutine <tt>DO_GEOS4_CONVECT</tt>:
| |
|
| |
| ! Free pointers
| |
| NULLIFY( STT, p_STT, p_F )
| |
|
| |
| --[[User:Bmy|Bob Y.]] 17:08, 30 May 2014 (EDT)
| |
|
| |
| === Fixed minor issues in MERRA cloud convection routine ===
| |
|
| |
| Jenny Fisher discovered two minor issues in the [[MERRA]] convection routine <tt>DO_MERRA_CONVECTION</tt> (in module <tt>convection_mod.f</tt>):
| |
|
| |
| ==== ND38 diagnostic issue in routine DO_MERRA_CONVECTION ====
| |
|
| |
| <span style="color:green">'''''This update was tested in the 1-month benchmark simulation [[GEOS-Chem_v9-01-02_benchmark_history#v9-01-02j|v9-01-02j]] and approved on 16 Aug 2011.'''''</span>
| |
|
| |
| '''''[mailto:jafisher@fas.harvard.edu Jenny Fisher] wrote]:'''''
| |
|
| |
| :The bug is in routine <tt>DO_MERRA_CONVECTION</tt>. When <tt>DIAG38</tt> is first written (around line 2056 in my code, after the comment titled "(3.4) ND38 Diagnostic"), DIAG38 is formed by adding T0 (converted to appropriate units) to the diagnostic each level K. This is correct. However, the second time we add to <tt>DIAG38</tt> (below the cloud base, around line 2245 in my code), we add <tt>T0_SUM</tt> to each level. From what we can tell, <tt>T0_SUM</tt> represents not the washout from each level but instead the total wet scavenging loss, summed over the entire column. So we are effectively double counting in these levels and adding way too much. [Helen Amos and I] think this should instead be adding <tt>WETLOSS</tt>.
| |
|
| |
| The fix is below:
| |
|
| |
| IF ( OPTIONS%USE_DIAG38 .and. F(K,IC) > 0d0 ) THEN
| |
| DIAG38(K,IC) = DIAG38(K,IC)
| |
| <span style="color:red">!------------------------------------------------------------------------------
| |
| ! Prior to 8/16/11:
| |
| ! Now use WETLOSS instead of T0_SUM in the ND38 diagnostic below the cloud
| |
| ! base. WETLOSS is the loss in this level, but T0_SUM is the loss summed
| |
| ! over the entire column. Using T0_SUM leads us to over-count the tracer
| |
| ! scavenged out of the column. (jaf, hamos, bmy, 8/16/11)
| |
| ! & + ( T0_SUM * AREA_M2 / TCVV_DNS )
| |
| !------------------------------------------------------------------------------</span>
| |
| & + ( WETLOSS * AREA_M2 / TCVV_DNS )
| |
| ENDIF
| |
|
| |
| ==== Wrong IF tests used in routine DO_MERRA_CONVECTION ====
| |
|
| |
| <span style="color:green">'''''This update was tested in the 1-month benchmark simulation [[GEOS-Chem_v9-01-02_benchmark_history#v9-01-02j|v9-01-02j]] and approved on 16 Aug 2011.'''''</span>
| |
|
| |
| <span style="color:red">'''''NOTE: In [[GEOS-Chem v11-01]] and higher versions, tagged Hg indices are now stored in the species database.'''''</span>
| |
|
| |
| :Routine <tt>DO_MERRA_CONVECTION</tt> was using tracer flags <tt>IDTHg2</tt> and <tt>IDTHgP</tt> to determine whether routines <tt>ADD_Hg2_WD</tt>, <tt>ADD_HgP_WD</tt>, and <tt>ADD_Hg2_SNOWPACK</tt> should be called. However, <tt>IDTHg2</tt> and <tt>IDTHgP</tt> aren't ever defined and were set to zero, so these routines were never called even when the correct mercury species were used.
| |
|
| |
| The fix is to replace these two lines of code
| |
|
| |
| IF ( IS_Hg .and. IC == IDT%Hg2 ) THEN
| |
| ...
| |
| IF ( IS_Hg .and. IC == IDT%HgP ) THEN
| |
|
| |
| with these:
| |
|
| |
| IF ( IS_Hg .and. IS_Hg2( IC ) ) THEN
| |
| ...
| |
| IF ( IS_Hg .and. IS_HgP( IC ) ) THEN
| |
|
| |
| --[[User:Bmy|Bob Y.]] 16:47, 16 August 2011 (EDT)<br>--[[User:Bmy|Bob Yantosca]] ([[User talk:Bmy|talk]]) 19:12, 6 May 2016 (UTC)
| |
|
| |
| == Outstanding issues that are not yet resolved ==
| |