Hybrid vertical coordinate for ARW

[davegill]

TYPE: new feature

KEYWORDS: hybrid eta isobaric coordinate

SOURCE: internal

DESCRIPTION OF CHANGES:

1. The build of WRF now requires a modification to the CPP since the variable argument lists macros is utilized. We can no longer use "-traditional-cpp", which does not permit the variable argument macro. Getting rid of -traditional-cpp means that the apostrophe is assumed to be part of a string. The build now gets rid of apostrophe situations within a Fortran comment. The "%" character (the separator between a derived data object and a field (such as grid%u_2) cannot easily be checked-on or parsed within CPP. Part of the build turns the grid%mu strings into a parseable single string, the CPP is run, then the build turns the parseable single string back into grid%mu. We now must include -C so that the Fortran concatenation operator "//" is interpreted correctly.
2. There are two types of major changes to the source code to accommodate the hybrid coordinate. The first involves mu (still defined as δpd/δη) and the second involves pressure. The 3d dry pressure is now defined as:
   pd = B(η) ( pds – pt ) + ( η – B(η) ) ( p0 – pt ) + pt
   Note that when B(η) = η then the vertical coordinate simplifies to the previous terrain following coordinate. When B(η) = 0, then the vertical coordinate is isobaric. The B(η) array is computed internally based on a couple of nameless options.
   C3 = B(η)
   C4 = ( η – B(η) ) ( p0 – pt )

The mu field, is no longer defined as the column pressure:
μ = δpd/δη = (δB(η)/δη) pc + ( 1 – δB(η)/δη ) ( p0 – pt )

We set
C1 = δB(η)/δη
and
C2 = ( 1 – δB(η)/δη ) ( p0 – pt )

This the functional form of total (or reference) mu as:
μ = C1( η ) pc + C2 ( η )
and the perturbation mu as:
μ = C1( η ) pc

All of the C1, C2, C3, C4 arrays are defined on both half and full levels.

1. Following is a before and after of the CPP styled modifications introduced into some of the files:

BEFORE:

define mu(...) (c1(k)*XXPCTXX(VA_ARGS)+c2(k))

define XXPCTXX(...) mu(VA_ARGS)

define muavg(...) (c1(k)*XXPCAVGXX(VA_ARGS)+c2(k))

define XXPCAVGXX(...) muavg(VA_ARGS)

define mub(...) (c1(k)*XXPCBXX(VA_ARGS)+c2(k))

define XXPCBXX(...) mub(VA_ARGS)

REAL, dimension(ims:ime,jms:jme) , intent(in) :: mu, mub

CALL abcd ( ru_tendf, mu, mu_1, mub, muavg, config_flags

DO j = jts, jte
DO i = its,ite

phb(i,k,j) = (r_d/p1000)_(t_init(i,k,j)+t0)_mub(i,j)

tendency(i,k,j)=tendency(i,k,j)- &
(mrdx_0.5_((mu(i+1,j)+mu(i,j))_H1(i+1,k,j)- &
(mu(i-1,j)+mu(i,j))H1(i ,k,j))+ &
mrdy_0.5((mu(i,j+1)+mu(i,j))_H2(i,k,j+1)- &
(mu(i,j-1)+mu(i,j))_H2(i,k,j ))- &
msfty(i,j)mu(i,j)(H1avg(i,k+1,j)-H1avg(i,k,j)+ &
H2avg(i,k+1,j)-H2avg(i,k,j) &
)*rdzw(i,k,j) &

DO j = j_start, j_end
DO i = i_start, i_end
  MUAVG(i,j) = 0.25 * ( MU(i-1,j  ) + MU(i,j  ) +  &
                        MU(i-1,j-1) + MU(i,j-1) )

END DO
END DO

DO j = j_start, j_end
DO k = kts, ktf
DO i = i_start, i_end

titau(i,k,j) = muavg(i,j) * mtau(i,k,j)

END DO
END DO
END DO

AFTER:
REAL, dimension(ims:ime,jms:jme) , intent(in) :: mu, mub

CALL abcd ( ru_tendf, mu, mu_1, mub, muavg, config_flags

DO j = jts, jte
DO i = its,ite

phb(i,k,j) = (r_d/p1000)(t_init(i,k,j)+t0)(c1(k)*mub(i,j)+c2(k))

tendency(i,k,j)=tendency(i,k,j)- &
(mrdx_0.5_(((c1(k)_mu(i+1,j)+c2(k))+(c1(k)_mu(i,j)+c2(k)))_H1(i+1,k,j)- &
((c1(k)_mu(i-1,j)+c2(k))+(c1(k)_mu(i,j)+c2(k)))H1(i ,k,j))+ &
mrdy_0.5(((c1(k)_mu(i,j+1)+c2(k))+(c1(k)_mu(i,j)+c2(k)))_H2(i,k,j+1)- &
((c1(k)_mu(i,j-1)+c2(k))+(c1(k)_mu(i,j)+c2(k)))H2(i,k,j ))- &
msfty(i,j)(c1(k)mu(i,j)+c2(k))(H1avg(i,k+1,j)-H1avg(i,k,j)+ &
H2avg(i,k+1,j)-H2avg(i,k,j) &
)_rdzw(i,k,j) &

DO j = j_start, j_end
DO i = i_start, i_end
  MUAVG(i,j) = 0.25 * ( MU(i-1,j ) + MU(i,j ) + &
                        MU(i-1,j-1) + MU(i,j-1) )

END DO
END DO

DO j = j_start, j_end
DO k = kts, ktf
DO i = i_start, i_end

titau(i,k,j) = (c1(k)*muavg(i,j)+c2(k)) * mtau(i,k,j)

END DO
END DO
END DO

Note that the variable argument macro ONLY modified variables that have a set of following parentheses, and that the arguments within the parens are maintained.

1. For testing purposes with the hill2d_x case, additional namelists and input soundings are provided.
2. The hybrid option is activated via the namelist.input file:
   &dynamics
   hybrid_opt = 2
   etac = 0.3
   /
   Users may turn off the hybrid coordinate completely so that the coordinate is terrain following only. When hybrid_opt=0, the model is in terrain following mode. When the hybrid option is activated, the etac value (the critical eta value) indicates at which level the model becomes isobaric.

LIST OF MODIFIED FILES (annotated if not obvious, not required to be on a single line):
M Registry/Registry.EM
D Registry/registry.hyb_coord
M arch/configure_new.defaults
M arch/noopt_exceptions
M arch/noopt_exceptions_f
M arch/postamble_new
M dyn_em/couple_or_uncouple_em.F
M dyn_em/module_advect_em.F
M dyn_em/module_bc_em.F
M dyn_em/module_big_step_utilities_em.F
M dyn_em/module_diffusion_em.F
M dyn_em/module_em.F
M dyn_em/module_first_rk_step_part1.F
M dyn_em/module_first_rk_step_part2.F
M dyn_em/module_initialize_hill2d_x.F
M dyn_em/module_initialize_real.F
M dyn_em/module_polarfft.F
M dyn_em/module_small_step_em.F
M dyn_em/nest_init_utils.F
M dyn_em/solve_em.F
M dyn_em/start_em.F
M external/RSL_LITE/module_dm.F
M external/io_esmf/makefile
M main/Makefile
M main/ndown_em.F
M main/real_em.F
M phys/module_microphysics_driver.F
M phys/module_mp_morr_two_moment.F
M share/dfi.F
M share/mediation_integrate.F
M share/module_bc.F
D test/em_hill2d_x/input_sounding-U=10,N=0.01
D test/em_hill2d_x/input_sounding-U=15,N=0.01
D test/em_hill2d_x/input_sounding-layers-20mps
M test/em_hill2d_x/namelist.input
D test/em_hill2d_x/namelist.input-100m_hill-10mps-N^2=0.0001-30km_deep-20km_damping-dampcoef=0.1-etac=0.2
D test/em_hill2d_x/namelist.input-700m_hill-15mps-n^2=0.0001-25km_deep-15km_damping-dampcoef=0.08-etac=0.2
D test/em_hill2d_x/namelist.input-HILL
D test/em_hill2d_x/namelist.input-HILL-51
D test/em_hill2d_x/namelist.input-HILL-schar
D test/em_hill2d_x/test.csh
M test/em_real/namelist.input

TESTS CONDUCTED (explicitly state mandatory, voluntary, and assigned tests, not required to be on a single line):

1. Several hill2d cases indicate reasonable functioning of the hybrid coordinate
2. Several real data cases point to reasonable behavior (larger cases from 13 km RAP)
3. With a couple of 1-line mods to flip multiplications and divides for mu, when the option is disabled, bit-for-bit results are available for before vs after for regression test
4. Bit for bit results for DFI options 0, 1, 2, 3 (small test case, single processor)

Known limitations

1. Still working on parallel bit for bit with GSD domain
2. Nesting (standard, moving, vertical), FDDA are not enabled
3. Little to no safety net is in place for users. For example, setting etac=0.3 over the US with a 10 mb stop would be OK, while the same settings over the Himalayas would likely be troublesome.
4. All idealized cases need weighting terms c1, c2, c3, c4 introduced in module_initialize_<some_ideal_case>.F
5. Chem and DA need hybrid registry included

[llpcarson]

Dave - Have you done a full reg-test yet? I thought that parts of HWRF needed the -traditional flag for CPP, but that was some years ago, so perhaps it works OK without it now.

[davegill]

Laurie,
The NMM code came out unscathed in the regression test - honestly, about
the only thing that did! It is likely that the sed commands to remove the
apostrophes are what is allowing us to no longer use -traditional-cpp.
Dave

On Fri, Sep 16, 2016 at 12:05 PM, Laurie Carson notifications@github.com
wrote:

Dave - Have you done a full reg-test yet? I thought that parts of HWRF
needed the -traditional flag for CPP, but that was some years ago, so
perhaps it works OK without it now.

—
You are receiving this because you authored the thread.
Reply to this email directly, view it on GitHub
#5 (comment), or mute
the thread
https://github.com/notifications/unsubscribe-auth/AMFFegJH71C16qcKm9leGUe-uVMxrF8Jks5qqtp1gaJpZM4J_MYm
.

[jamiebresch]

We need hybrid_opt in global attributes of the wrfinput and wrfout to know what coordinate is used.

[davegill]

Jamie,
Yes, the global attributes are required.
I just started a list of things to do, this is #2.
Dave

On Fri, Sep 16, 2016 at 12:27 PM, Jamie Bresch notifications@github.com
wrote:

We need hybrid_opt in global attributes of the wrfinput and wrfout to know
what coordinate is used.

—
You are receiving this because you authored the thread.
Reply to this email directly, view it on GitHub
#5 (comment), or mute
the thread
https://github.com/notifications/unsubscribe-auth/AMFFeovrRAQ3QF1oFwpOfFlPEeGk1iemks5qqt95gaJpZM4J_MYm
.

[mkavulich]

@davegill At the beginning of dyn_em/module_bc_em.F you've added a whole bunch of empty lines; can you remove these?

[davegill]

Mike,
I initially went through and added quite a few blank lines in the new code,
and then did some "grep"ing on the various MU strings, with the associated
new line number (including black lines). Then I went back and started
adding in the required CPP directives for each of the necessary MU
variables. This purposefully left quite a few blank lines at the top of
every Fortran file so that my grep line numbers would not change. I'll go
back and review all of the modified files for these blank lines near the
top.
Dave

On Mon, Sep 19, 2016 at 6:05 PM, Michael Kavulich notifications@github.com
wrote:

@mkavulich requested changes on this pull request.

@davegill https://github.com/davegill At the beginning of
dyn_em/module_bc_em.F you've added a whole bunch of empty lines; can you
remove these?

—
You are receiving this because you authored the thread.
Reply to this email directly, view it on GitHub
#5 (review), or mute
the thread
https://github.com/notifications/unsubscribe-auth/AMFFesID_G_6DR_dMMILabqoqFX0cXD0ks5qryNkgaJpZM4J_MYm
.

[jamiebresch]

@davegill The CPP stuff is basically reverting what was done in the 41f94bb (Nov 1, 2013) commit.

[davegill]

Jamie,

Good catch for those Linux users who were originally hit with this problem.

For the hybrid coordinate cpp-automatic substitution of mu(i,j) into
(c1(k)*mu(i,j)+c2(k), we need to use the variable argument option. That
means that we cannot use -traditional-cpp (because the variable argument
option is not available with -traditional-cpp). When we get rid of
-traditional-cpp, then we starting interpreting C++ styled comments (the
"//" string) as C++ styled comments, and by default cpp throws comments
away. Unfortunately, those "//" strings are the Fortran concatenation
operator, which we do not want thrown away inside of the cpp text
processing. To remove the problem of the misinterpreted Fortran
concatenation operator, we need to stick the "-C" flag back onto the cpp
command (which basically says DO NOT GET RID OF COMMENTS).

Let me see if I can turn the "//" string inside of Fortran files into
something unprocessed bycpp, and then swap it back in after cpp has done
its processing. Here's what I have in mind, swapping in the string
"fort_concat" for the "//" string (and then back again).

/Users/gill>cat foo1

abd // abd

/Users/gill>sed -e "s////fort_concat/g" foo1

abd fort_concat abd

/Users/gill>sed -e "s////fort_concat/g" foo1 > foo2

/Users/gill>cat foo2

abd fort_concat abd

/Users/gill>sed -e "s/fort_concat////g" foo2

abd // abd

If I can get a full build to work this way, I can remove the -C option from
the cpp.

Dave

On Tue, Sep 20, 2016 at 1:55 PM, Jamie Bresch notifications@github.com
wrote:

@davegill https://github.com/davegill The CPP stuff is basically
reverting what was done in the 41f94bb
41f94bb
(Nov 1, 2013) commit.

—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
#5 (comment), or mute
the thread
https://github.com/notifications/unsubscribe-auth/AMFFeotOdE_OxqAmNBoC-86XJihBF3zSks5qsDpLgaJpZM4J_MYm
.

[mkavulich]

Closed in favor of #104