Evaluate time-dependent tendencies at correct time for each RK4 stage

[sbrus89]

This PR corrects a longstanding issue with RK4 and tendency terms which have an explicit time dependence (see https://github.com/E3SM-Project/E3SM/issues/5364 and https://github.com/MPAS-Dev/MPAS-Model/issues/375). Big thanks to @gcapodag and @jeremy-lilly for providing the fix in their LTS development branch linked in https://github.com/E3SM-Project/E3SM/issues/5364.

Here, I have migrated these changes from their MPAS-Dev branch to E3SM and added the modifications to correct the convergence issues for the manufactured solution test case: https://github.com/E3SM-Project/polaris/pull/72

[BFB] (for all standard E3SM tests)

[sbrus89]

Testing:

With this PR, the manufactured solution test case in Polaris converges at 2nd order when refining in space/time:

When refining in time only, 4th order convergence is achieved:

However, I have not been able to reproduce the 4th order convergence for the hurricane test case from Lilly et al. 2023 (Figure 4). The best I have been able to achieve is first order: For this result, I ran from 2012-10-13_00:00:00 to 2012-10-13_01:00:00 starting from a 2012-10-13_00:00:00 restart file. I used dt=0.1 as a reference solution.

[xylar]

Relevant discussion of this PR can be seen in https://github.com/E3SM-Ocean-Discussion/E3SM/pull/55

[mark-petersen]

Tested with gnu on chicoma for both debug and optimized, compared to the master branch point, for MPAS-Ocean standalone and the compass nightly suite. I am able to get a bfb match between this branch and master. In particular, the RK4 tests are BFB:

00:02 PASS ocean_global_ocean_Icos240_WOA23_RK4_performance_test
00:07 PASS ocean_global_ocean_Icos240_WOA23_RK4_restart_test
00:05 PASS ocean_global_ocean_Icos240_WOA23_RK4_decomp_test
00:05 PASS ocean_global_ocean_Icos240_WOA23_RK4_threads_test

That makes sense, because config_use_time_varying_atmospheric_forcing and config_use_time_varying_land_ice_forcing are both false for these tests, so the PR should change nothing.

One small note: The step ocean_global_ocean_Icos240_WOA23_init didn't match between this PR and master, which seems unrelated to this PR. But then I reran it and it matched. That is a mystery, but maybe I just made a mistake in my run sequence.

I also tested E3SM on perlmutter. This does not exercise the RK4 code, but it at least compiles it. E3SM passes:

./create_test SMS_Ld3.T62_oQU240.GMPAS-NYF.pm-cpu_intel -p m4572 -q debug --walltime 00:30:00
./create_test SMS_Ld3.T62_oQU240.GMPAS-NYF.pm-cpu_gnu  -p m4572 -q debug --walltime 00:30:00

[xylar]

@mark-petersen, how are you able to do standalone testing on Chicoma? There shouldn't be compatible spack packages until https://github.com/MPAS-Dev/compass/pull/865 goes in.

[mark-petersen]

Hmm... I'm on the head of compass/main and created the load_dev* this morning, then compiled and it ran fine. Here are my steps:

git log -n 1
commit cab53e2c8dc9be0efbd4a7a6c2150b0c599c1870 (HEAD -> main, origin/main, origin/HEAD)
Merge: 0e79f61e3 8b35df7f3
Author: Xylar Asay-Davis <[email protected]>
Date:   Sat Oct 19 11:04:22 2024 -0600
    Merge pull request #836 from xylar/add-normalized-dib-dismf

./conda/configure_compass_env.py --conda ${HOMEDIR}/miniforge3 --compiler gnu
source load_dev_compass_1.5.0-alpha.1_chicoma-cpu_gnu_mpich.sh
compass suite -s -c ocean -t nightly -p $e/master/components/mpas-ocean -w $n/$DIR

etc. And I had compiled stand-alone in separate directories, not in the compass submodule. It all worked fine, so I don't know.

[xylar]

I see. I'm a little surprised that works. If you were to try to recreate the load script, I don't think it would work -- at least that's what others have reported.

[mark-petersen]

I think the first two plots above show that this PR is working as intended. Those manufactured solutions provide the forcing at a requested moment in time. The hurricane test case (third plot) shows 1st rather than expected 4th order convergence. That was tricky to obtain in Lilly et al. 2023, and there must be some small difference in the way that the forcing data is interpolated in time. I'm happy for this PR to be merged despite that difference. The manufactured solution shows that this code in this PR is correct. In practice the wind data is updated at 1 or 6 hour intervals, so the model does not need interpolation that satisfies fourth-order convergence at very small time intervals.

[rljacob]

@cbegeman @xylar any updates?

[jonbob]

Passes:

• SMS_D_Ld1.ne30pg2_r05_IcoswISC30E3r5.WCYCL1850.chrysalis_intel.allactive-wcprod

merged to next

[jonbob]

merged to master