Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

update CICE from Consortium/main #59

Merged
merged 70 commits into from
Mar 29, 2023

Conversation

DeniseWorthen
Copy link
Collaborator

apcraig and others added 30 commits July 15, 2022 07:42
* Rename several variables to make compatible with C-grid names,
  strairxU, strairyU, strtltxU, strtltyU, strintxU, strintyU,
  taubxU, taubyU, strocnxU, strocnyU, fmU, TbU,
  waterxU, wateryU, forcexU, forceyU, aiU
Move iceumask, icenmask, iceemask from ice_flux to ice_grid
Make dyn_prep2, stepu, stepuv_CD, stepv_C, implicit_solver, and
  anderson_soler argument names a little more generic/specific
Inline boxslotcyl velocity calculation

* remove boxslotcyl_data_vel

* update documentation

* Additional updates to change to upper case for some variable names
Fix indentations as noted
* Add optargs "optional arguments" unit test.  This tests the ability
to pass optional arguments down the calling tree robustly whether they are
present or not.

* Add test to count optional arguments at 2nd level
…ortium#731)

When running a test suite that includes 'bfbcomp' tests, the test
script ('cice.test') queries the file 'suite.jobs' in the suite
directory to get the job ID of the test against which to compare for
bit-for-bit-ness.

If re-running the suite, 'suite.jobs' is not removed, so that the new
job IDs are appended to the existing file. This leads to syntax errors
when 'cice.test' looks for the job ID to compare against since the
'grep' call returns several matches.

Remove 'suite.jobs' at the start of the 'suite.submit' script generated
by 'cice.setup' to avoid that.
* Deprecate CESM ponds

* Namelist changes for deprecating cesmponds

* Update documentation
)

* initial 0-layer thermo deprecation

* capitalization matters for cpps

* set_nml.boxadv needs thermo turned on

* deprecate old ridging participation and redistribution functions

* Revert "deprecate old ridging participation and redistribution functions"

This reverts commit 95c289a.
Bring in deprecated 0-layer and cesmponds changes in Icepack.
…tium#742)

* Update cice.t-test.py to use cartopy instead of basemap.

* Bug fix to add gridlines to SH plots

* commented out contour section of plots. Default is pcolor. If contour is selected it will instead make a pcolor plot

* cice.t-test.py: addded individual colorbar to each plot. environment.yml: removed basemap, added cartopy

* Remove shrink option from difference plots

* Add blockall distribution_wght to set_nml.qc to address plotting issues in qc test

Co-authored-by: daveh150 <[email protected]>
…tium#745)

* cice.setup: allow command line to override suite options

Options chosen on the 'cice.setup' command line (using the '-s' flag) are
added to the options defined for each test in the test suite definition
file, when running a test suite.

This is implemented by appending the options from the test suite (in
variable 'sets_tmp') to the options from the command line ('sets_base')
in the variable 'sets', which is ultimately (via the variable 'setsx')
looped through to apply each option.

Since 'sets_tmp' is appended to 'sets_base', options on the command line
can't override options from the test suite file, which means one can't
e.g. run a test suite with 'kdyn=3' and expect all tests to use this
option if any option specified in the test suite set 'kdyn' to some
other value.

To allow options from the command line to override options from the test
suite, reverse the order in which 'sets_tmp' and 'sets_base' are used to
define 'sets'. This is in line with the common behaviour of the command
line taking precedence.

Adjust the documentation accordingly, fixing a typo along the way.

* cice.setup: name test suite cases less ambiguously

In the previous commit, we allowed options from the command line to
override options from the test suite definition file. However, test case
directories are still named using a sorted list of all active options,
both from command line and the suite definition file (variable
'setsarray'). This is nonoptimal since it is not clear from looking at
the test directory name which options have precedence in case of
conflict.

Change the naming logic so that options from the command line are last
in the test directory name, in a "last-one-wins" fashion. To do that,
let 'setsarray' be defined from the test suite options ('sets_tmp') and
add a second loop for the command line options ('sets_base').

Note that we do not check if the same option is mentioned both in the
test suite and the command line, in order not to complicate the code
further.

This also allows greatly simplifying the logic used to adjust 'bfbcomp'
test names to include command line options. This part of the code is
checking if the options for the test aginst which to compare contain
duplicates and 'none', which is unnecessary since these options come
directly from the test suite definition file.
* Update Icepack
Update Version

* Set visc_method='avg_strength' for gridCD to avoid some aborts
Fix a few bugs in the test suite lists
…ium#749)

At the end of subroutine ice_grid::gridbox_corners, the arrays
'lont_bounds' and 'lonu_bounds', which contain the longitude of the
corners of each grid cell on the T and U grids, are converted to to the
[0, 360] range.

In the case of rectangular grids ('grid_type = rectangular'), at the
point where 'gridbox_corners' is called in 'init_grid2', 'lont_bounds'
is not initialized, causing the code to abort if compiling with NaN
initialization. This is due to the fact that 'gridbox_verts', which
initializes 'lont_bounds' and 'latt_bounds', is not called in
'rectgrid', whereas it is called in 'popgrid[_nc]'.

Do call 'gridbox_verts' in 'rectgrid', so that 'lont_bounds' and
'latt_bounds' are correctly initalized in that case also.

Note that these calls are also missing in 'latlongrid' and 'cpomgrid',
but since these two subroutines are not used in standalone
configuration, let's not bother for now.
* Update/improve debug_blocks output, see CICE-Consortium#718.

* Add ICE_MEMUSE cice.settings flag for batch memory use
Add set_env.memsmall, memmed, memlarge options
To use, will require changes to the env machine files.  Most machines will probably not use it.
See CICE-Consortium#674.

* Add setup_machparams.csh to compute batch/launch machine parameters
Update cice.batch.csh and cice.launch.csh to use setup_machparams.csh
See CICE-Consortium#650

* Update subroutine diagnostic_abort which calls print_state
Update ice_transport_remap and ice_transport_driver to call diagnostic_abort
  during some errors.
See also CICE-Consortium#622

* Update miniconda install information
See CICE-Consortium#547

* Code cleanup based on compile with -Wall
Code cleanup based on -std f2003 and f2008 checks
Add -stand f08 to cheyenne_intel debug flags
Add -std f2008 to cheyenne_gnu debug flags
Code consistent with Fortran 2003 except for use of contiguous in
  1d evp code.

* Remove all trailing blank space with script

* Update the cheyenne env so qc testing works
Add configuration/scripts/tests/qctest.yml file
Update documentation

* Update Icepack

* Clean up some output

* fix comments

* update print_state output
* machines: eccc: unify baseline directory

* machines: eccc: fix modules initialization

Make sure to source the Csh initialization script for environment
modules ourselves, as it is not done in all environments.

While at it, for convenience add I_MPI_LIBRARY_KIND=debug to the
commented lines.
The variable ice_calendar::nstreams, which corresponds to the number of
output history streams to use for the run, is initialized in
ice_history::init_hist depending on the number of non-'x' elements in
'histfreq' in the namelist.

However, the code does use 'nstreams' before ice_history::init_hist is
called, in ice_calendar::calendar when called from
ice_calendar::init_calendar. Both 'init_calendar' and 'init_hist' are
called from CICE_InitMod::cice_init, in that order, such that the loop
that initializes 'write_history' in 'calendar' uses 'nstreams'
uninitialized.

'calendar' ends up being called at least once more during 'cice_init', from
ice_calendar::advance_timestep, at which point 'nstreams' is correctly
defined and 'write_history' is thus correctly initialized, before its
first use in 'accum_hist'.

To avoid using 'nstreams' uninitialized in the first call to 'calendar'
from 'init_calendar', initialize it to zero in 'init_calendar' before
calling 'calendar'.

This issue was discovered by compiling using the '-init=huge' flag of
the Intel compiler.
CICE-Consortium#746)

* Added dxgrow, dygrow to facilitate variable spaced grid. Modified rectgrid to generate grid from center outward using growth factor

* Adding vargrid namelist options.

* Refactored rectgrid to compute dx/dy first. Then ULON/ULAT. Added scale_dxdy flag to check of want grid spacing scaled. Renamed dxgrow/dygrow to dxscale/dyscale.

* Added method to check for odd nx_global/ny_global when applying grid spacing scale factors

* Update comments before computing dx/dy

* Update comments when checking for even/odd

* made grid_lonref, grid_latref namelist varaibles. Removed vargrid_suite.ts. Updated the box nml to specify the default grid_lonref and grid_latref for future reference.

* Change grid_lonref/grid_latref to lonrefrect,latrefrect. Reduce default vargrid tests to 3 per B,C,CD grid.

* Make new subroutine rectgrid_scale_dxdy to implemnet grid scaling. Remove explicit latrefrec/lonrefrect from set_nml. Make dxscale,dyscale,latrefrec,lonrefrec double precition in ice_in

* Add set_nml.scale1 to test vargrid with dxscale/dyscale = 1.d0

* Remove lonrefrec/latrefrec from boxnodyn

* Add lonrefrect, latrafrect to documentation

* Inserted new scaled grid varibles in alphabetical order in documentation
…CE-Consortium#754)

* Update Pull Request Template, add question about updating Icepack
* Refactored evp sub cycling loop

* corrected indent and case for dyn_haloUptdate
…um#758)

In 3fedc78 (Allow for read of tlat, tlon, anglet with popgrid (CICE-Consortium#463),
2020-06-24), ice_grid::init_grid2 was changed so that ice_grid::Tlatlon,
which computes the TLAT and TLON arrays from ULAT and ULON, is only
called if the private module variable 'l_readCenter' is false.

The idea is that if the grid file contains a variable 'anglet', then it
is assumed that it also contains variables 'tlon' and 'tlat', and so
these fields are read directly instead of being computed.

This logic, however, was only implemented in ice_grid::popgrid_nc, which
sets 'l_readCenter' depending on the presence or absence of 'anglet' in
the grid file. This means that if 'popgrid_nc' is not called (for
example with "grid_format='bin'", in which case init_grid2 calls
'popgrid' and not 'popgrid_nc'), then 'l_readCenter' is used
uninitialized in init_grid2, and so it's possible that 'Tlatlon' is not
called, if 'l_readCenter' happens to be initialized to true.

This in turns leads to 'TLAT' and 'TLON' being uninitialized, and the
code failing when accessing these arrays if compiling with NaN
initialization.

Fix this by initializing 'l_readCenter' at the beginning of init_grid2,
such that it is initialized for all choices of 'grid_format' and
'grid_type'. Remove the initialization in 'popgrid_nc'.
* move deformation out of loop for B grid only

* Moved C and CD grid deformations

* correct location of bgrid deformation
* Update dEdd implementation

- Update Icepack with several fixes (changes answers for tr_aero)
- Update bgc aerosol table to higher precision
- Add modal aerosol tests
- Update test suites

* Update Icepack including dEdd fixes
)

* Update grid averaging for tmass, aice, uvelT, vvelT

- Update tmass and aice T2U mapping, switch from "F" to "S", F was backwards compatible but not correct (changes answers)
- Update ocean forcing T2U averaging in ocn_data_ncar_init, change "F" to "A".
- Update uvelT, vvelT averaging in step_therm1, change from 4 point average to U2TA (changes answers for highfreq=.true.)
- Remove history grids not needed (i.e. ustr3Dz)

* Refactor uvelT, vvelT implementation
Mention that CICE must be cloned with '--recurse-submodules' for Icepack
to also be cloned, formatting the 'git clone' command as a code block,
and fix the link to the Git Workflow Guide.
* Refactor strocnxT, strocnyT implementation

- add aiU to ice_state
- migrate computation of strocnxT and strocnyT to ice_step, needed for thermodynamics,
  better code reuse.
- add strocnxT_sf, strocnyT_sf as coupling field, could be computed differently than
  the thermodynanics version.  The _sf field computation should be in scale fluxes, but
  because scale_fluxes is called on a block and the _sf fields require a halo update
  among other things, the computation can't be done in scale_fluxes.
- Update the coupling layers to use the _sf version of the fields.
- CICE-Consortium#761 suggests the values of strocnxT,
  strocnyT should not be scaled for use in thermodynamics.  This commit does not make
  that change yet, but allows for that change to be made easily.
- These changes are bit-for-bit for a full test suite on cheyenne with 3 compilers.

* Update computation of strocnxT, strocnyT passed into icepack_step_therm1

- No longer divided by aice
- strocnxT_sf, strocnyT_sf are still computed in the same way as before

* Rename strocn[x,y]T_sf to strocn[x,y]T_iavg

Revert strocn[x,y]T passed into thermodynamics to be the version divided
by aice, specifically strocn[x,y]T_iavg.  This is identical to earlier
implementations.
When the 'default_season' namelist setting was added in 01494c7 (Nml
settings (CICE-Consortium#208), 2018-10-19) to replace 'l_winter' and 'l_spring', a
call to 'broadcast_scalar' was forgotten, such that the 'default_season'
value from the namelist is only used on the first MPI process; all other
processes get the hardcoded default value 'winter' defined in
'ice_init::input_data', resulting in different initialization across the
grid for several variables if anything other than 'winter' is used in
the namelist.

Fix that by broadcasting 'default_season' to all MPI procs.
- refactor options boxsym and boxislands to boxopen, boxclosed, and boxforce
- update box test names as needed
- update calculation of test suite failures in cice.results.csh
- add documentation to user guide for rectangular grids
- Compute and used cdn_ocn on U, E, and N cell locations as needed for dynamics.
- Add halo updates in dynamics before calling grid_average.  This doesn't
  change answers, but is the safe thing to do in general.  Performance does
  not seem to be affected.
…sk (CICE-Consortium#773)

* Change icetmask to logical consistent with iceumask, icenmask, iceemask

- Add icetmask as logical array to ice_grid.F90, was integer array
- Update use of icetmask in code for consistency with new type
- Add ice_HaloUpdate2DL1 to support halo updates for logical fields in both mpi and serial ice_boundary.F90
- Modify some capital T,U,N,E in ice_dyn_shared.F90 to t,u,n,e for better consistency in code

* Update cicecore/cicedynB/infrastructure/comm/mpi/ice_boundary.F90

* Update comment in code

* Revert changes to T,U,N,E in ice_dyn_shared.F90, working toward additional changes

* Move ice[T,U,N,E}mask from ice_grid to ice_dyn_shared

* rename icell[t,u,n,e] to icell[T,U,N,E], rename indx[t,u,n,e] to indx[T,U,N,E]

* remove ice[t,u,n,e]grid from ice_grid
…CE-Consortium#774)

* doc: fix typo in index (bfbflag)

* doc: correct default value of 'maxits_nonlin'

The "Table of namelist options" in the user guide lists 'maxits_nonlin'
as having a default value of 1000, whereas its actual default is 4, both
in the namelist and in 'ice_init.F90'. This has been the case since the
original implementation of the implicit solver in f7fd063 (dynamics: add
implicit VP solver (CICE-Consortium#491), 2020-09-22).

Fix the documentation.

* doc: VP solver is validated with OpenMP

When the implicit VP solver was added in f7fd063 (dynamics: add implicit
VP solver (CICE-Consortium#491), 2020-09-22), it had not yet been tested with OpenMP
enabled.

The OpenMP implementation was carefully reviewed and then fixed in
d1e972a (Update OMP (CICE-Consortium#680), 2022-02-18), which lead to all runs of the
'decomp' suite completing and all restart tests passing. The 'bfbcomp'
tests are still failing, but this is due to the code not using the CICE
global sum implementation correctly, which will be fixed in the next
commits.

Update the documentation accordingly.

* ice_dyn_vp: activate OpenMP in 'dyn_prep2' loop

When the OpenMP implementation was reviewed and fixed in d1e972a (Update
OMP (CICE-Consortium#680), 2022-02-18), the 'PRIVATE' clause of the OpenMP directive
for the loop where 'dyn_prep2' is called in 'implicit_solver' was
corrected in line with what was done in 'ice_dyn_evp', but OpenMP was
left unactivated for this loop (the 'TCXOMP' was not changed to a real
'OMP' directive).

Activate OpenMP for this loop. All runs and restart tests of the
'decomp_suite' still pass with this change.

* machines: eccc : add ICE_MACHINE_MAXRUNLENGTH to ppp[56]

* machines: eccc: use PBS-enabled OpenMPI for 'ppp6_gnu'

The system installation of OpenMPI at /usr/mpi/gcc/openmpi-4.1.2a1/ is
not compiled with support for PBS. This leads to failures as the MPI
runtime does not have the same view of the number of available processors
as the job scheduler.

Use our own build of OpenMPI, compiled with PBS support, for the
'ppp6_gnu'  environment, which uses OpenMPI.

* machines: eccc: set I_MPI_FABRICS=ofi

Intel MPI 2021.5.1, which comes with oneAPI 2022.1.2, seems to have an
intermittent bug where a call to 'MPI_Waitall' fails with:

    Abort(17) on node 0 (rank 0 in comm 0): Fatal error in PMPI_Waitall: See the MPI_ERROR field in MPI_Status for the error code

and no core dump is produced. This affects at least these cases of the
'decomp' suite:

- *_*_restart_gx3_16x2x1x1x800_droundrobin
- *_*_restart_gx3_16x2x2x2x200_droundrobin

This was reported to Intel and they suggested setting the variable
'I_MPI_FABRICS' to 'ofi' (the default being 'shm:ofi' [1]). This
disables shared memory transport and indeeds fixes the failures.

Set this variable for all ECCC machine files using Intel MPI.

[1] https://www.intel.com/content/www/us/en/develop/documentation/mpi-developer-reference-linux/top/environment-variable-reference/environment-variables-for-fabrics-control/communication-fabrics-control.html

* machines: eccc: set I_MPI_CBWR for BASEGEN/BASECOM runs

Intel MPI, in contrast to OpenMPI (as far as I was able to test, and see
[1], [2]), does not (by default) guarantee that repeated runs of the same
code on the same machine with the same number of MPI ranks yield the
same results when collective operations (e.g. 'MPI_ALLREDUCE') are used.

Since the VP solver uses MPI_ALLREDUCE in its algorithm, this leads to
repeated runs of the code giving different answers, and baseline
comparing runs with code built from the same commit failing.

When generating a baseline or comparing against an existing baseline,
set the environment variable 'I_MPI_CBWR' to 1 for ECCC machine files
using Intel MPI [3], so that (processor) topology-aware collective
algorithms are not used and results are reproducible.

Note that we do not need to set this variable on robert or underhill, on
which jobs have exclusive node access and thus job placement (on
processors) is guaranteed to be reproducible.

[1] https://stackoverflow.com/a/45916859/
[2] https://scicomp.stackexchange.com/a/2386/
[3] https://www.intel.com/content/www/us/en/develop/documentation/mpi-developer-reference-linux/top/environment-variable-reference/i-mpi-adjust-family-environment-variables.html#i-mpi-adjust-family-environment-variables_GUID-A5119508-5588-4CF5-9979-8D60831D1411

* ice_dyn_vp: fgmres: exit early if right-hand-side vector is zero

If starting a run with with "ice_ic='none'" (no ice), the linearized
problem for the ice velocity A x = b will have b = 0, since all terms in
the right hand side vector will be zero:

- strint[xy] is zero because the velocity is zero
- tau[xy] is zero because the ocean velocity is also zero
- [uv]vel_init is zero
- strair[xy] is zero because the concentration is zero
- strtlt[xy] is zero because the ocean velocity is zero

We thus have a linear system A x = b with b=0, so we
must have x=0.

In the FGMRES linear solver, this special case is not taken into
account, and so we end up with an all-zero initial residual since
workspace_[xy] is also zero because of the all-zero initial guess
'sol[xy]', which corresponds to the initial ice velocity. This then
leads to a division by zero when normalizing the first Arnoldi vector.

Fix this special case by computing the norm of the right-hand-side
vector before starting the iterations, and exiting early if it is zero.
This is in line with the GMRES implementation in SciPy [1].

[1] https://github.com/scipy/scipy/blob/651a9b717deb68adde9416072c1e1d5aa14a58a1/scipy/sparse/linalg/_isolve/iterative.py#L620-L628

Close: phil-blain#42

* ice_dyn_vp: add global_norm, global_dot_product functions

The VP solver uses a linear solver, FGMRES, as part of the non-linear
iteration. The FGMRES algorithm involves computing the norm of a
distributed vector field, thus performing global sums.

These norms are computed by first summing the squared X and Y components
of a vector field in subroutine 'calc_L2norm_squared', summing these
over the local blocks, and then doing a global (MPI) sum using
'global_sum'.

This approach does not lead to reproducible results when the MPI
distribution, or the number of local blocks, is changed, for reasons
explained in the "Reproducible sums" section of the Developer Guide
(mostly, floating point addition is not associative). This was partly
pointed out in [1] but I failed to realize it at the time.

Introduce a new function, 'global_dot_product', to encapsulate the
computation of the dot product of two grid vectors, each split into two
arrays (for the X and Y components).

Compute the reduction locally as is done in 'calc_L2norm_squared', but
throw away the result and use the existing 'global_sum' function when
'bfbflag' is active, passing it the temporary array used to compute the
element-by-element product.

This approach avoids a performance regression from the added work done
in 'global_sum', such that non-bfbflag runs are as fast as before.

Note that since 'global_sum' loops on the whole array (and not just ice
points as 'global_dot_product'), make sure to zero-initialize the 'prod'
local array.

Also add a 'global_norm' function implemented using
'global_dot_product'. Both functions will be used in subsequent commits
to ensure bit-for-bit reproducibility.

* ice_dyn_vp: use global_{norm,dot_product} for bit-for-bit output reproducibility

Make the results of the VP solver reproducible if desired by refactoring
the code to use the subroutines 'global_norm' and 'global_dot_product'
added in the previous commit.

The same pattern appears in the FGMRES solver (subroutine 'fgmres'), the
preconditioner 'pgmres' which uses the same algorithm, and the
Classical and Modified Gram-Schmidt algorithms in 'orthogonalize'.

These modifications do not change the number of global sums in the
fgmres, pgmres and the MGS algorithm. For the CGS algorithm, there is
(in theory) a slight performance impact as 'global_dot_product' is
called inside the loop, whereas previously we called
'global_allreduce_sum' after the loop to compute all 'initer' sums at
the same time.

To keep that optimization, we would have to implement a new interface
'global_allreduce_sum' which would take an array of shape
(nx_block,ny_block,max_blocks,k) and sum over their first three
dimensions before performing the global reduction over the k dimension.

We choose to not go that route for now mostly because anyway the CGS
algorithm is (by default) only used for the PGMRES preconditioner, and
so the cost should be relatively low as 'initer' corresponds to
'dim_pgmres' in the namelist, which should be kept low for efficiency
(default 5).

These changes lead to bit-for-bit reproducibility (the decomp_suite
passes) when using 'precond=ident' and 'precond=diag' along with
'bfbflag=reprosum'. 'precond=pgmres' is still not bit-for-bit because
some halo updates are skipped for efficiency. This will be addressed in
a following commit.

[1] CICE-Consortium#491 (comment)

* ice_dyn_vp: do not skip halo updates in 'pgmres' under 'bfbflag'

The 'pgmres' subroutine implements a separate GMRES solver and is used
as a preconditioner for the FGMRES linear solver. Since it is only a
preconditioner, it was decided to skip the halo updates after computing
the matrix-vector product (in 'matvec'), for efficiency.

This leads to non-reproducibility since the content of the non-updated
halos depend on the block / MPI distribution.

Add the required halo updates, but only perform them when we are
explicitely asking for bit-for-bit global sums, i.e. when 'bfbflag' is
set to something else than 'not'.

Adjust the interfaces of 'pgmres' and 'precondition' (from which
'pgmres' is called) to accept 'halo_info_mask', since it is needed for
masked updates.

Closes CICE-Consortium#518

* ice_dyn_vp: use global_{norm,dot_product} for bit-for-bit log reproducibility

In the previous commits we ensured bit-for-bit reproducibility of the
outputs when using the VP solver.

Some global norms computed during the nonlinear iteration still use the
same non-reproducible pattern of summing over blocks locally before
performing the reduction. However, these norms are used only to monitor
the convergence in the log file, as well as to exit the iteration when
the required convergence level is reached ('nlres_norm'). Only
'nlres_norm' could (in theory) influence the output, but it is unlikely
that a difference due to floating point errors would influence the 'if
(nlres_norm < tol_nl)' condition used to exist the nonlinear iteration.

Change these remaining cases to also use 'global_norm', leading to
bit-for-bit log reproducibility.

* ice_dyn_vp: remove unused subroutine and cleanup interfaces

The previous commit removed the last caller of 'calc_L2norm_squared'.
Remove the subroutine.

Also, do not compute 'sum_squared' in 'residual_vec', since the variable
'L2norm' which receives this value is also unused in 'anderson_solver'
since the previous commit. Remove that variable, and adjust the
interface of 'residual_vec' accordingly.

* ice_global_reductions: remove 'global_allreduce_sum'

In a previous commit, we removed the sole caller of
'global_allreduce_sum' (in ice_dyn_vp::orthogonalize). We do not
anticipate that function to be ued elsewhere in the code, so remove it
from ice_global_reductions. Update the 'sumchk' unit test accordingly.

* doc: mention VP solver is only reproducible using 'bfbflag'

The previous commits made sure that the model outputs as well as the log
file output are bit-for-bit reproducible when using the VP solver by
refactoring the code to use the existing 'global_sum' subroutine.

Add a note in the documentation mentioning that 'bfbflag' is required to
get bit-for-bit reproducible results under different decompositions /
MPI counts when using the VP solver.

Also, adjust the doc about 'bfbflag=lsum8' being the same as
'bfbflag=off' since this is not the case for the VP solver: in the first
case we use the scalar version of 'global_sum', in the second case we
use the array version.

* ice_dyn_vp: improve default parameters for VP solver

During QC testing of the previous commit, the 5 years QC test with the
updated VP solver failed twice with "bad departure points" after a few
years of simulation. Simply bumping the number of nonlinear iterations
(maxits_nonlin) from 4 to 5 makes these failures disappear and allow the
simulations to run to completion, suggesting the solution is not
converged enough with 4 iterations.

We also noticed that in these failing cases, the relative tolerance for
the linear solver (reltol_fmgres = 1E-2) is too small to be reached in
less than 50 iterations (maxits_fgmres), and that's the case at each
nonlinear iteration. Other papers mention a relative tolerance of 1E-1
for the linear solver, and using this value also allows both cases to
run to completion (even without changing maxits_nonlin).

Let's set the default tolerance for the linear solver to 1E-1, and let's
be conservative and bump the number of nonlinear iterations to 10. This
should give us a more converged solution and add robustness to the
default settings.
apcraig and others added 19 commits December 7, 2022 15:00
Update Icepack to latest release version
Remove trailing blank space
CICE-Consortium#801)

* Adding method to check namelist in any order. Use subroutine in ice_namelist_mod.F90 to search for namelist in ice_in.

* Moved goto_nml subroutine to ice_fileunits.F90. Removed ice_namelist_mod.F90

* Cleanup indentations with tmpstr2 use

* Cleanup spacing and intentation

* For namelist check, remove extra continuation after making ice_abort string.

Co-authored-by: Tony Craig <[email protected]>
* Fix OMP setup

* Update meshgrid
allow nudiag_set to be available outside of cesm; may prefer
to fix in coupling interface
…rtium#810)

* Modified doc for Dupont et al ref

* Minor modif to v_i calculation in seabed_stress_factor_prob to prevent rare instability
* Fix OMP setup

* Update meshgrid

* Small updates for CESM

* Add change for UFS
* Add time_period_freq to history file metadata
… variables (CICE-Consortium#817)

* doc/source/conf.py: adjust for sphinxcontrib.bibtex 2.0

The sphinxcontrib.bibtex Sphinx extension used for the bibliography now
wants the bibliography file to be configured in the Sphinx configuration
file (conf.py) instead of in the source file where the bibliography is
included. This is new in sphinxcontrib.bibtex 2.0 [1], so let's do that.

Keeping the filename also in zreferences.rst does not hurt and lets us
stay compatible with earlier versions of sphinxcontrib.bibtex, so let's
keep it there also.

[1] https://sphinxcontrib-bibtex.readthedocs.io/en/latest/changes.html#id5

* ice_history: correct units for 'sigP'

The intenal ice pressure 'sigP' is is units of N/m, as can be seen in
ice_dyn_shared::principal_stress. However, the corresponding history
variable is wrongly defined in ice_history::init_hist with unit '1'
(dimensionless). This means the wrong unit is written to the NetCDF
history output. This dates back to the introduction of that variable in
6ed2359 (Added pressure, modified norm of principal stresses and made
small modifs to basal stress following Till's comments, 2018-03-02).

Fix the unit. While at it, add an entry for 'sigP' in the index, from
which this variable is missing.

Reported-by: Frederic Dupont <[email protected]>
Reported-by: Jean-Francois Lemieux <[email protected]>

* doc: clarify stress variables

Try to make the doc a little less confusing by cross-referencing the
code variables used for stress computations with the corresponding
variables in the science guide a little bit more, and vice-versa:

- mention the doc variables sigma_1, sigma_2 in the index entries for
stressp, stressm
- mention the code variables stressp, stressm when the doc variables
sigma_1, sigma_2 are introduced
- introduce new doc variables sigma_n,1 and sigma_n,2 to denote the
normalized principal stresses, and add the equation for those. This
allows mentioning that they are normalized by the ice strength, which
was not mentioned elsewhere.
- mention these new doc variables in the index entry for sig1, sig2
- refer to the normal stress sigma_11, sigma_22 by their variable names
when mentioning them in the sentence that introduces the ice pressure
- mention the code variables sig1, sig2 in the "Implementation" part of
the user guide when mentioning the 'principal_stresses' subroutine.

Helped-by: Jean-Francois Lemieux <[email protected]>
* ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).

* drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.
* Initial halochk unit test implementation

- Add halochk unit test
- Add "unknown" and "noupdate" checks to ice_boundary
- Remove field_loc_Wface, not used anywhere, not supported
- Update cice_decomp.csh script

To Do: validate tripole and tripoleT, add unit test to test suite

* - Fix bug in serial/ice_boundary.F90 tripoleT halo update
- Reduce redundant tripole buffer copies in serial/ice_boundary.F90
- Generalize iSrc wraparound calculation in ice_boundary.F90
- Add open, cyclic, tripole, and tripoleT set_nml files
- Update unittest suite

* - Add haloUpdate_stress test to halo unit test
- Add tripoleT support to haloUpdate_stress
- Add abort check that nx_global is even for tripole grids
- Update documentation

* update documentation

* update documentation

* update documentation

* update documentation

* update documentation

* update documentation

* Update halochk test to make haloupdate_stress test more robust, less chance for false positive

* update documentation

* update documentation

* update documentation
to fix the latex/pdf errors recently trapped by readthedocs.
* send wlat to and fro for FSD

* Update Icepack to include FSDmods

---------

Co-authored-by: cmbitz <[email protected]>
@zach1221
Copy link

This PR has been reviewed, and all tests are done on PR #1672. @NeilBarton-NOAA are you able to merge this PR?

@NeilBarton-NOAA
Copy link
Collaborator

I'm not authorized to merge this. maybe @DeniseWorthen?

@DeniseWorthen
Copy link
Collaborator Author

Thanks @NeilBarton-NOAA.

@zach1221 I'm the code manager for CICE and generally do the merges. I assume we're ready now, right?

@zach1221
Copy link

@DeniseWorthen, understood. Yes, please proceed to merge this PR.

@DeniseWorthen DeniseWorthen merged commit 5840cd1 into NOAA-EMC:emc/develop Mar 29, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

updates from Consortium/main