Version 8 (modified by snguyen, 2 years ago) (diff) |
---|

# How the code for the Schwarz method implementation in NEMO is organized

## Implementation of the Schwarz method

### The time stepping loop

The essential part to understand is how the original time stepping loop of NEMO in `nemogcm.F90` is modified to allow rewinding of the state of NEMO after a coupling window has passed.

Initially the code has one loop over the command `CALL stp( istp )` for each time step.

With the Schwarz algorithm, the total amount of time steps is divided in coupling windows (called Schwarz loops in the code), each coupling window is repeated the number of schwarz iterations given as a parameter (mswr), and for each schwarz iteration a number of time steps (the length of a coupling window) is done.

For a computation of 5 days with 15 time steps per day, a coupling window of 1 day and 6 Schwarz iterations there will be 3 loops.

- One global loop over the coupling windows : 5 to iterate. Counter
`iswloop`, limits`1 to nsloops`. - One inner loop over the number of schwarz iteration for each coupling window : 6. Counter
`kswr`, limits`1 to mswr`. - One most inner loop over the time steps done in one Schwarz coupling window : 15 time steps. Counter
`istp`, limits`nit000 + (iswloop-1)*ntsinswr to nit000+iswloop*ntsinswr-1`

The counter given to subroutine `stp(istp)` is rewind to its value before entering the Schwarz loop when the Schwarz coupling window is repeated.

This is coded like this in `nemogcm.F90`:

iswloop = 1 IF (lwp) WRITE(numout,*) 'init iswloop =',iswloop DO WHILE ( iswloop <= nsloops .AND. nstop == 0 ) ! iterate on nsloops schwarz loops kswr = 1 IF (lwp) THEN WRITE(numout,*) '*** schwarz loops ***' WRITE(numout,*) 'iswloop =',iswloop WRITE(numout,*) 'kswr = ',kswr ENDIF CALL swz_store ! store Ocean state at first schwarz iteration before first time step CALL swz_store_lim3 ! store lim3 variables state CALL swz_store_limdiahsb ! store limdiahsb variables IF(lwp) WRITE(numout,*) 'store schwarz initial state ' DO WHILE ( kswr <= mswr .AND. nstop == 0 ) IF (lwp) WRITE(numout,*) ' kswr = ',kswr istp = nit000 + (iswloop-1)*ntsinswr ! set time step to first value for iswloop (current schwarz loop) IF ( kswr > 1 ) CALL swz_reinit(istp) ! reset Ocean state to first time step for new schwarz iteration IF (lwp) WRITE(numout,*) 'maxistp = ',nit000+iswloop*ntsinswr-1 DO WHILE ( istp <= nit000+iswloop*ntsinswr-1 .AND. nstop == 0 ) IF (lwp) WRITE(numout,*) ' istp = ',istp,'; istpswz = ',istpswz #if defined key_agrif CALL stp ! AGRIF: time stepping #else CALL stp( istp ) ! standard time stepping #endif istp = istp + 1 istpswz = istpswz + 1 IF( lk_mpp ) CALL mpp_max( nstop ) END DO ! istp kswr = kswr + 1 END DO ! kswr iswloop = iswloop +1 END DO ! iswloop

The variable `iswloop` is the current Schwarz loop and `nsloops` is the number of Schwarz loops for the given time of simulation and Schwarz window size. For 5 days with a window of 1 day `nsloops=5`.

Inside the Schwarz repeating loop `DO WHILE ( kswr <= mswr .AND. nstop == 0 )`, in the time stepping loop,
`istp` starts at the value `nit000 + (iswloop-1)*ntsinswr` and finishes at ` nit000+iswloop*ntsinswr-1`.

This allows `istp` to take the value it would have for the current Schwarz window if no sub-iteration was done: setting `mswr=1` gives you the original time stepping.

The Schwarz parameters are initialized (or computed) before entering the time loops, before and after the call to `nemo_init` depending on what is needed:

ntsinswr = 1 !! set to allow fldread in nemo_init before computation of schwarz indices kswr = 0 !! set here to allow day_init to set nitrst properly WRITE(numout,*) "ntsinswr = ",ntsinswr ! !-----------------------! CALL nemo_init !== Initialisations ==! ! !-----------------------! !! set integer constants for computation ! ncplfrq = 86400 !cpl_freq( 'O_QnsMix' ) read in namrun (nemoinit) ntsinswr = ncplfrq / INT(rdt) nsloops = (nitend-nit000+1)/ntsinswr !! outputs to check values IF(lwp) THEN ! control print WRITE(numout,*) WRITE(numout,*) ' computation of loop indices for schwarz computation ' WRITE(numout,*) ' nit000 = ',nit000 WRITE(numout,*) ' nitend = ',nitend WRITE(numout,*) ' rdt = ',rdt WRITE(numout,*) ' ncplfrq = ',ncplfrq WRITE(numout,*) ' ntsinswr = ',ntsinswr WRITE(numout,*) ' nsloops = ',nsloops WRITE(numout,*) ' mswr = ',mswr WRITE(numout,*) ' nitrst = ',nitrst ENDIF

There is also an "absolute" time step counter : `istpswz` which is incremented at each time step including those repeated for the Schwarz algorithm. It is needed for the coupling with OASIS.

### The Schwarz algorithm: storing and restoring

The Schwarz algorithm is implemented by storing and restoring the state of the ocean and ice at the proper position inside the three loops structure.

The routines used to store and restore a state are all in `schwarz.F90` which is a new module.

The storing is done just before the schwarz looping is done: inside the `iswloop` loop, before the `kswr` loop.

CALL swz_store ! store Ocean state at first schwarz iteration before first time step CALL swz_store_lim3 ! store lim3 variables state CALL swz_store_limdiahsb ! store limdiahsb variables IF(lwp) WRITE(numout,*) 'store schwarz initial state ' DO WHILE ( kswr <= mswr .AND. nstop == 0 )

Like this we are certain to store the initial condition for a Schwarz loop.

The restoring is done at the beginning of the `kswr` loop before the `istp` loop.
If we are not in the first iteration the variables state is restored:

IF ( kswr > 1 ) CALL swz_reinit(istp) ! reset Ocean state to first time step for new schwarz iteration

The routine `swz_reinit` is reproducing the initialisation structure from `nemo_init` in `nemogcm.F90` by selecting only the relevant subroutines and coding the corresponding restoring of variables.

(In case you wish to use PISCES you will have to add the relevant routines in `swz_reinit` and module `schwarz.F90`)

`Unfortunately this is not the only place variables are stored and restored!`

Because NEMO is modular some implementation details are repeated in the modules some modules have their own storing and restoring details.

The files concerned are:

sbcmod.F90 sbcrnf.F90 traqsr.F90 trasbc.F90

## Coding details: modification of other files

The routines:

iom.F90 step.F90

have been modified to output only one time serie for a set of Schwarz loops. The parameter `ksout` selects which Schwarz loop is given for outputs.

Some routines have a particular behavior when `istp=nit000` or/and `istp=nitend`. They where modified to behave correctly when doing Schwarz iterations:

closea.F90 diaar5.F90 diafwb.F90 diahth.F90 dynspg_ts.F90 dynvor.F90 fldread.F90 limdyn.F90 limitd_me.F90 limrst.F90 limthd.F90 limtrp.F90 limupdate1.F90 limupdate2.F90 restart.F90 sbcblk_core.F90 sbcice_lim.F90 sbcssm.F90 stpctl.F90 traadv_eiv.F90 traadv_tvd.F90 trazdf_imp.F90 zdftke.F90 zdftmx.F90

Some other routines which are not compiled have also been modified for the sake of consistency. These are the routines relative to PISCES present in `ORCA1_LIM3_PISCES/MY_SRC` which get copied in `ORCA2_LIM3_PISCES`.

They are:

p4zflx.F90 p4zsed.F90 p4zsms.F90 trcrst.F90 trcsms_age.F90 trcsms_cfc.F90 trcstp.F90 trcwri.F90

You will need these updates if you want to add PISCES support to the Schwarz algorithm.

Lastly, but much more important, some files have been modified to declare the variables needed to store the parameters of the Schwarz algorithm and the storing/restoring variables.

The parameters are in:

in_out_manager.F90

and read in

domain.F90

The "Schwarz fields" named with `_s` appended at the end of the original variable name are defined in these subroutines:

dom_oce.F90 dynspg_oce.F90 ice.F90 limdiahsb.F90 oce.F90 sbc_oce.F90 sbcrnf.F90 zdf_oce.F90 zdftke.F90

## a few words on coupling with OASIS

The coupling with OASIS is straightforward.

The calls to `sbc_cpl_snd` in `step.F90` and `sbc_cpl_rcv` in `sbcmod.F90` are modified to use the absolute counter `istpswz` instead of `istp`.
In `step.F90`:

!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> ! Coupled mode !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< IF( lk_oasis ) CALL sbc_cpl_snd( istpswz ) ! coupled mode : field exchanges

In `sbcmod.F90`:

CASE( jp_purecpl ) ; CALL sbc_cpl_rcv ( istpswz, nn_fsbc, nn_ice ) ! pure coupled formulation

The only case treated is the pure coupling mode.

## a few last words on testing and validation

The NEMO part of the Schwarz algorithm was validated by running a forced configuration with given atmospheric boundary conditions. The code reproduces mswr times the same day over 5 days with a coupling and Schwarz window of one day. The result was compared to the original code with no Schwarz iterations in the same configuration. The outputs are identical for the same day and any Schwarz iteration this day. Hence the storing and restoring is reliable.