GEOS-ESM / GEOSgcm_GridComp

Repository containing the physics and IAU code for the GEOS Earth System Model
Apache License 2.0
9 stars 7 forks source link

Cannot build with FV_PRECISION=R8 #60

Closed kgerheiser closed 4 years ago

kgerheiser commented 5 years ago

I want to build gcm_GridComp with FV_PRECISION=R8 so that I can build the ctm and gcm in one build. However, a recent update to FV has introduced errors.

First, in FV_StateMod.F90:

https://github.com/GEOS-ESM/GEOSgcm_GridComp/blob/b4446a6504226c3933fd34f4a842a44ba8b7333b/GEOSagcm_GridComp/GEOSsuperdyn_GridComp/FVdycoreCubed_GridComp/FV_StateMod.F90#L3564-L3567

The values 2.e3 and 5.e3 do not agree with the function's dummy arguments since updraft_helicity in fv_diagnostics.F90 is compiled with the -r8 flag which causes default reals to be treated as r8.

This is easy enough to fix by making those values variables of FVPRC, and then passing them.

But if you do that then compilation fails in DynCore_GridCompMod.F90 because of this line:

https://github.com/GEOS-ESM/GEOSgcm_GridComp/blob/b4446a6504226c3933fd34f4a842a44ba8b7333b/GEOSagcm_GridComp/GEOSsuperdyn_GridComp/FVdycoreCubed_GridComp/DynCore_GridCompMod.F90#L4756

The output of fv_getUpdraftHelicity is defined as FVPRC, but DynCore expects it to be r4 and the export it is filling is r4.

My solution is to change the output of the subroutine to be r4, but if you do that then you need to introduce a temporary output variable of FVPRC that you pass to updraft_helicity (because it wants FVPRC, not r4) and then at the end of the subroutine assign that intermediate value to the output which is defined as r4.

This introduces an extraneous temporary variable, but I think it's the easiest way to make it work.

mathomp4 commented 5 years ago

Added @wmputman as an assignee and @tclune and @aoloso as well as they are in FV3 land a bit at present.

If you have a branch with these changes, it would be pretty simple to test for zero-diff to develop

tclune commented 5 years ago

The first change above (making the literals variables of configurable KIND) seems fine to me.

Still investigating the second issue. Description is a bit confusing, as fv_getUpdraftHelicity() is a subroutine, not a function so it does not have an output kind. So presumably you are saying that the 1st argument of the subroutine is FPRC inside but that temp2d is declared as R4 outside.

tclune commented 5 years ago

Note: Whoever submits the pull request for this should also fix the various declarations (e.g., temp2d) that are declared with KIND=4. This is not standard conforming. E.g., with the NAG Fortran compiler KIND=1 is 32 bit real with the default settings and can be something rather weird with more aggressive debug flags.

tclune commented 5 years ago

Definitely a question for @wmputman: temp2d is used far too often in the outer layer to blindly increase the KIND for FVPRC=8.

mathomp4 commented 5 years ago

Query for @tclune: when these are fixed/corrected, do they also need to be applied to the @FVdycoreCubed_GridComp repository as well?

tclune commented 5 years ago

It's simply a race condition. It is in the FV GridComp whereever it happens to live when the PR goes in. My preference would be in the new home, as I'd like to move the location asap.

kgerheiser commented 5 years ago

@tclune

What I did to get the code to build was to change the intent(out) output to be REAL4:

https://github.com/GEOS-ESM/GEOSgcm_GridComp/blob/dd376d9ee8aecc9e54b6e9ab620c400718687a2d/GEOSagcm_GridComp/GEOSsuperdyn_GridComp/FVdycoreCubed_GridComp/FV_StateMod.F90#L3557

And then introduce an intermediate value uh25_tmp that is of FVPRC so that you call updraft_helicity with uh25_tmp. And then assign uh25_tmp to the actual output of uh25

subroutine fv_getUpdraftHelicity(uh25)
   use constants_mod, only: fms_grav=>grav
   use fv_diagnostics_mod, only: get_vorticity, updraft_helicity

! made this REAL4
   real(REAL4), intent(OUT) :: uh25(FV_Atm(1)%bd%isc:FV_Atm(1)%bd%iec,FV_Atm(1)%bd%jsc:FV_Atm(1)%bd%jec)

! made an intermediate output of FVPRC
   real(FVPRC) :: uh25_tmp(FV_Atm(1)%bd%isc:FV_Atm(1)%bd%iec,FV_Atm(1)%bd%jsc:FV_Atm(1)%bd%jec)

   integer :: sphum=1
   real(FVPRC) :: vort(FV_Atm(1)%bd%isc:FV_Atm(1)%bd%iec,FV_Atm(1)%bd%jsc:FV_Atm(1)%bd%jec,FV_Atm(1)%npz)

   ! introduced these two variables for the literals
   real(FVPRC) :: z_bot, z_top

   z_bot = 2.e3
   z_top = 5.e3
   call get_vorticity(FV_Atm(1)%bd%isc, FV_Atm(1)%bd%iec, FV_Atm(1)%bd%jsc, FV_Atm(1)%bd%jec, &
                      FV_Atm(1)%bd%isd, FV_Atm(1)%bd%ied, FV_Atm(1)%bd%jsd, FV_Atm(1)%bd%jed, &
                      FV_Atm(1)%npz, FV_Atm(1)%u, FV_Atm(1)%v, vort, &
                      FV_Atm(1)%gridstruct%dx, FV_Atm(1)%gridstruct%dy, FV_Atm(1)%gridstruct%rarea)

! call this with uh25_tmp
   call updraft_helicity(FV_Atm(1)%bd%isc, FV_Atm(1)%bd%iec, FV_Atm(1)%bd%jsc, FV_Atm(1)%bd%jec, FV_Atm(1)%ng, FV_Atm(1)%npz, &
                     zvir, sphum, uh25_tmp, &
                     FV_Atm(1)%w, vort, FV_Atm(1)%delz, FV_Atm(1)%q,   &
                     FV_Atm(1)%flagstruct%hydrostatic, FV_Atm(1)%pt, FV_Atm(1)%peln, FV_Atm(1)%phis, fms_grav, z_bot, z_top)

! cast back to r4
   uh25 = uh25_tmp
end subroutine fv_getUpdraftHelicity

An alternative would be to do the same thing in the DynCore level: introduce a new variable to get the helicity, and then assign it to the export pointer.

lizziel commented 5 years ago

I am putting these updates into FVdycoreCubed_GridComp which I am now testing with GCHP. Would you like me to do a pull request to feed back to you? If yes, which branch?