Mass Flux (ND24/25/26)

Sign and direction convention

ND24: zonal flux (EW)

• Flux across the western edge per grid box
• Positive if eastward

ND25: meridional flux (NS)

• Flux across the southern edge of the box
• Positive if northward

ND26: vertical flux (UP)

• Flux across the top edge of the grid box
• Positive if downwards

ND26 Vertical Advection Bug Fix in v11-02b

To achieve mass balance in the advection diagnostic output (ND24, ND25, and ND26), Ilya Stanevic (U. Toronto) suggested that the FZ array be used to compute the ND26 vertical mass flux diagnostic for all levels. This is similar to how FX is used for ND24 and FY for ND25. Note that this fix is only for global advection since the FZ array is not defined for the nested model. With this update, the global ND26 vertical flux diagnostic is defined as the flux into the grid box above (e.g. ND26(:,:,1) is flux from level 1 to level 2 and ND26(:,:,LLPAR) is all zeros). This fix is added into v11-02b.

Consequences of the Polar Cap on Mass Fluxes

The explanation below is intended for GEOS-Chem versions 8-01-03 and newer, that is the versions with the updated advection.

In this new advection subroutine, we consider a polar cap at both poles spanning for 2 bands of latitude. A polar cap means that the concentrations in all the grid cells over this region are averaged over the region in order to form only one circular grid cell for a given altitude. Actually, the code is presently averaging the concentrations separately by band of latitude so we keep a gradient between the 2 first and 2 last bands of latitude.

Saying that each polar cap forms one circular grid cell means that the advection subroutine only calculates transport between the polar cap and the adjacent grid cells, and does no calculation within the extent of the polar cap. This has different consequences depending on the mass flux direction. To explain what happens, let's consider the South pole. It's the same for the North pole, you simply need to change the latitudes. At the South Pole, the polar cap is formed of the first and second bands of latitudes and all longitudes: J=1 and J=2, and I=1:Imax, for each altitude.

Eastward mass flux (EW) in polar cap, fx

Within the polar cap, grid cells for a given latitude and all longitudes are considered forming one grid cell. Thus there is no longitudinal transport calculated. The EW diagnostic in the polar cap regions returns values of 0 for all grid cells. In other words:

  fx(I,J,K) = 0 for I=1:Imax, J=1 or J=2, and all K (altitudes).

Vertical mass flux (UP) in polar cap, fz

This diagnostic is calculated for each grid cell separately. Each grid cell of the first band of latitude will have the same value since the concentrations for these grid cells are averaged in the advection. It will be the same for the second band of latitude, but the values of the first and second band of latitudes will be different from each other. So to get vertical mass fluxes for the whole polar cap region at the altitude K, you simply need to add up the values for all the grid cells (South pole):

  SUM( FZ( *, 1:2, K) ) in Fortran notations
  TOTAL( FZ( *, 0:1, K) ) in IDL notations

Northward mass flux (NS) in polar cap, fy

To visualize what is happening for the northward mass flux, the best is to think in polar projection. In this document, the drawing on the left represents the polar cap and the third band of latitude in polar projection. As I said, the advection code calculates the transport only between the polar cap and the adjacent cells. This means that only the fluxes going through the border between the 2d and 3d bands of latitude are calculated (blue arrows on drawing). The value of these mass fluxes will likely be different for each grid cell since the concentrations for the 3d latitude are most likely different in each grid cell.

But, the ND25 diagnostic gives some values for the 1st and 2d borders of the grid, see red arrows on the drawing on the right. In fact these fluxes are not calculated based on winds but they are geometrical averages of the fluxes calculated at the 3d border (blue arrows). So I would advise people to only use the transport between the polar cap region and the adjacent regions (blue arrows) which are properly calculated, and not to consider finer resolution in the polar cap (red arrows). If you need the total latitudinal mass flux going in/out of the polar cap (South pole), you would calculate:

  SUM( FY( *, 3, K) ) in Fortran notations
  TOTAL( FY( *, 2, K) ) in IDL notations

--Ccarouge 02:00, 26 November 2010 (EST)

Error in NS flux estimate

When trying to rigorously balance one box, I found that the North-South flux is not 100% accurate. This is simply due to the box surface area and latitude. In GEOS-Chem, the mass gain in one box is proportional to:

( south_unit_flux - north_unit_flux ) * geometry( box ), [A]

where geometry depends on box surface and latitude. But we only save south*geometry(box). So when we try to balance the saved fluxes, the mass gain becomes:

south_unit_flux * geometry( box )  -  north_unit_flux * geometry( box above ), [B]

which can significantly differ from [A]. I got almost 2% in one time step for one box. In opposite hemispheres, the difference between the two calculations will tend to be opposite. One can imagine that errors will more or less cancel each other when balancing very large region. Because box surface does not depend on longitude, there is no similar problem with east-west flux.

Related bug fix

(1) For GEOS-3, two bugs have been found in August 2007. One is the initialization of the flux that was not correctly reset to zero when dealing with a new tracer. The other is that the UP/DOWN flux between the first and second layers was not accounted for. This has been corrected in v7-04-12 or 13.

(2) For GEOS-4, GCAP, and GEOS-5, one line of code must be uncommented. It modifies flux diagnostic for boxes in the top two layers. It is line 797 of tpcore_fvdas_mod.f90. It says:

MASSFLUP(I,J,K,IQ) = MASSFLUP(I,J,K,IQ) + DTC(I,J,K,IQ)/dt

in the loop for K=1 when ND26>0. This has been corrected in v8-01-01.

--phs 17:25, 4 April 2008 (EDT)

ND27 and variable tropopause

Dave MacKenzie wrote:

In diag3.f there are some notes mentioning how this diagnostic is no good for GEOS4. Do you remember exactly why? At one point in the file it says something like "should look the same as UP-FLX-$ for O3 at 200hPa" and at another it says "should not look the same as UP-FLX-$ for O3 at 200hPa." I think the latter is correct, but regardless I'm not sure if this is outdated or what the cause of the problem might be. I'm using GEOS-Chem (GEOS4) v8-01-01 at 2x2.5. Any thoughts would be helpful.

Bob Yantosca replied:

I think a lot of this comes down to historical baggage. When we had the annual-mean-tropopause in the code, there was also a set of flags (the IFLX in diag3.f) that were generated specifically for the annual-mean-tropopause files. These flags told the directions from which stratospheric air was coming into a given grid box. (Qinbin Li originally set this up.)

However, since we have now gone to the variable tropopause (for GEOS-4 and GEOS-5), the old methodology for ND27 is no longer valid. Theoretically it is possible to compute the same IFLX flags but the problem is now they would vary in time and space so that makes the diagnostic archiving a little more trickly. And we had more pressing things to worry about at the time so we never dealt with this rigorously.

--Bob Y. 13:23, 13 July 2009 (EDT)