Problem with nudging time scales at open boundaries

[nencio]

Hi all,

I am running a simulation with 4 open boundaries, using radiation conditions with nudging for tracers and 3D velocity (cpp options ***_TRADIATION, ***_TNUDGING, ***_M3RADIATION, and ***_MNUDGING defined).
I also wanted the nudging coefficients to vary linearly from the value at the boundary to 0, within 14 points from the boundary. To do so I activated the cpp options TCLIM_NUDGING and M3CLIM_NUDGING, and modified the routine ana_nudgcoef.h (I attached it to this post).

With my configuration the values of the nudging coefficients at the 4 boundaries are set in ana_nudgcoef.h through (e.g. nudging coefficients for 3D velocity at the Western boundary):

Code: Select all

      IF (SOUTH_WEST_CORNER) THEN
        M3obc_out(ng,iwest)=0.5_r8*                                     &
     &                      (CLIMA(ng)%M3nudgcof(0,1)+                  &
     &                       CLIMA(ng)%M3nudgcof(1,1))
        M3obc_in (ng,iwest)=obcfac(ng)*M3obc_out(ng,iwest)
      END IF

If I run the simulation using one cpu the coefficients at the 4 boundaries look ok (they are output in the his netcdf file). However, if I run the same exact simulation with multiple processors (I tried 2x1 and 2x2 tiling), then their value at the Eastern and Northern open boundaries is set to 0.
Has anybody else encountered a similar problem?

I played around with the code for a bit, and to me it seems it could be a problem of comunication between the different cpus. The values of the coefficients are computed correctly at all four boundaries, but then coming out of the ana_nudgcoef.h routine the values of M3obc_out,M3obc_in,Tobc_out and Tobc_in are 0 at ieast and inorth.

If I don't use the TCLIM_NUDGING and M3CLIM_NUDGING (and so, no nudging coefficients in the interior), then the coefficients are set (still in ana_nudgcoef.h) through:

Code: Select all

      IF (SOUTH_WEST_TEST) THEN
        M3obc_out(ng,iwest)=M3nudg(ng)
        M3obc_in (ng,iwest)=obcfac(ng)*M3nudg(ng)
      END IF

and in this case the values in the history file are ok, either I run single or multiprocessor.

Thanks a lot,
francesco

[nencio]

Hi all,

After few days spent on the code, here's my understanding on what's happening, and the way I fixed it.
Any comment, correction or further suggestion is very welcomed.

Let's consider only the open boundary nudging coefficients for free surface(FSobc_in, FSobc_out), in a simulation in which the domain is tiled 2x2.
FSobc_in, FSobc_out are declared in Modules/mod_scalars.F:

Code: Select all

        real(r8), dimension(Ngrids,4) :: FSobc_in
        real(r8), dimension(Ngrids,4) :: FSobc_out

If Ngrids =1 (no nested grids) then they are 1x4 arrays, each element corresponding to a boundary: 1 = west; 2 = south; 3 = east; 4 = north.

The values of FSobc_in, FSobc_out are assigned in Functional/ana_nudgcoef.h, which is called for each tile of the domain by Nonlinear/initial.F:

Code: Select all

#ifdef NUDGING_COFF
!
!-----------------------------------------------------------------------
!  If appropriate, set nudging coefficiests time scales.
!-----------------------------------------------------------------------
!
!$OMP PARALLEL DO PRIVATE(thread,subs,tile) SHARED(ng,numthreads)
        DO thread=0,numthreads-1
          subs=NtileX(ng)*NtileE(ng)/numthreads
          DO tile=subs*thread,subs*(thread+1)-1
            CALL ana_nudgcoef (ng, TILE, iNLM)
          END DO
        END DO
!$OMP END PARALLEL DO
#endif

This means that each tile define its own version of both M3obc_in, M3obc_out.

If the ssh field is nudged to climatological values (cpp option ZCLIM_NUDGING defined) then the 4 values are assigned in Functional/ana_nudgcoef.h as:

Code: Select all

# ifdef WEST_FSNUDGING
#  ifdef ZCLM_NUDGING
      IF (SOUTH_WEST_CORNER) THEN
        FSobc_out(ng,iwest)=CLIMA(ng)%Znudgcof(0,1)
        FSobc_in (ng,iwest)=obcfac(ng)*FSobc_out(ng,iwest)
      END IF
...
# ifdef EAST_FSNUDGING
#  ifdef ZCLM_NUDGING
      IF (NORTH_EAST_CORNER) THEN
        FSobc_out(ng,ieast)=CLIMA(ng)%Znudgcof(Lm(ng)+1,Mm(ng))
        FSobc_in (ng,ieast)=obcfac(ng)*FSobc_out(ng,ieast)
      END IF
...
# ifdef SOUTH_FSNUDGING
#  ifdef ZCLM_NUDGING
      IF (SOUTH_WEST_CORNER) THEN
        FSobc_out(ng,isouth)=CLIMA(ng)%Znudgcof(1,0)
        FSobc_in (ng,isouth)=obcfac(ng)*FSobc_out(ng,isouth)
      END IF
...
# ifdef NORTH_FSNUDGING
#  ifdef ZCLM_NUDGING
      IF (NORTH_EAST_CORNER) THEN
        FSobc_out(ng,inorth)=CLIMA(ng)%Znudgcof(Lm(ng),Mm(ng)+1)
        FSobc_in (ng,inorth)=obcfac(ng)*FSobc_out(ng,inorth)
      END IF

According to this, in the case of a 2x2 tiling, the Bottom Left tile updates the values of its FSobc_in, FSobc_out at iwest and isouth; The Upper Right tile updates the values of its owns at ieast and inorth; wheras the FSobc_in, FSobc_out of the Upper Left and Bottom Right tile remain unchanged.

To send the updated value to the other tiles, the function mp_bcastf is called at every boundary. Considering the western boundary as an example:

Code: Select all

# ifdef WEST_FSNUDGING
#  ifdef ZCLM_NUDGING
      IF (SOUTH_WEST_CORNER) THEN
        FSobc_out(ng,iwest)=CLIMA(ng)%Znudgcof(0,1)
        FSobc_in (ng,iwest)=obcfac(ng)*FSobc_out(ng,iwest)
      END IF
...
#   ifdef DISTRIBUTE
      IF (ng.eq.Ngrids) THEN
        CALL mp_bcastf (ng, model, FSobc_out(:,iwest), Ngrids)
        CALL mp_bcastf (ng, model, FSobc_in (:,iwest), Ngrids)
      ENDIF
#   endif
#  else

I don't know too much about MPI, but I think that in principle this function is supposed to exchanged the updated values from one tile to the others, in order to have updated FSobc_in, FSobc_out for all 4 tiles. However it doesn't work, since the resulting nudging coefficients are ok for the western and southern boundary, but are 0 for northern and eastern boundary (it seems that instead of exchanging the updated values, it assings the values from the Bottom Left tile to the others).

To fix this problem I changed the code above into:

Code: Select all

# ifdef WEST_FSNUDGING
#  ifdef ZCLM_NUDGING
        FSobc_out(ng,iwest)=FSnudg(ng)
        FSobc_in (ng,iwest)=obcfac(ng)*FSobc_out(ng,iwest)
...
#  else

which is exactly the same way open boundary nudging coefficients are defined in the case of undefined ZCLIM_NUDGING. Values of FSnudg and obcfac are defined from the External/*.in file.

It's definitely not the most elegant of the solutions, but it works in my case.
Again, any further comment/correction/suggestion is welcomed.

francesco

[arango]

Yes, good catch This is a parallel (distributed-memory) bug. Your solution is also incorrect. The problem here is that the MPI broadcasting function mp_bcastf (master node broadcast to all children) is used We need to do reduction operation instead. The CLIMA(ng)%Znudgcof is dimensioned different in all the tile partitions. It is dimensioned to the local size of the tile but with global indices. For example, in a 2x2 partition, CLIMA(ng)%Znudgcof(0,1) is only bounded in tile 0 (lower left corner tile). If you compile with arrays bounds you will get an error.

I will fix this problem

[nencio]

...I wasn't sure which one of the mp_* functions in Utility/distribute.F to use to fix the problem...

Another thing I noticed (that might give some problems using mp_reduce) looking at Modules/mod_scalars.F is that only Tobc_out and Tobc_in are initialized to IniVal. Open boundary nudging coefficients for FS, M2 and M3 are not.

francesco

[arango]

Yes, we need to initialize those values in mod_scalars.F and use mp_collect. Notice that, currently, mp_collect uses mpi_allreduce.

Please update your code. This bug was corrected in the following trac ticket.