ROMS 3.0 possible bug with volume conservation

[lanerolle]

For a particular application that I am running, I need to use the volume conservation CPP options in ROMS (I use EAST_VOLCONS, WEST_VOLCONS and SOUTH_VOLCONS).

I test to see whether the volume is indeed being conserved by integrating the (normal) velocities along each of the open boundary faces of the model domain and then summing the contributions up; I do this inside the ROMS itself (in ROMS/Nonlinear/diag.F).

I find that in a SERIAL model configuration, I get near perfect volume conservation (the residual volume ~ 1.0 x 10^-17) but with MPI the volume conservation is destroyed and I get a large net residual volume. The ROMS however does not blow-up but continues to run in a stable fashion.

I have looked at the obc_volcons.F routine (in ROMS/Nonlinear/) but I could not figure out how to fix this MPI-related bug. Could someone please help me out?

Thanks a lot,

Lyon.

[wilkin]

Lyon,

It sounds to me like a bug in your MPI modifications to the diagnostics, but you'd need to post this to the forum for us to be sure. To perform a summation over multiple tiles requires some care in the implementation in the coding since you can't simply accumulate a sum over i,j into a scalar variable subtotal because you must control what each processor knows about the variable.

That said, some caution is required if you use only some of NORTH_, SOUTH_, WEST_, EAST_ VOLCONS. The logic on the OBC_VOLCONS options is only correct if all open boundaries are active. If you look at obc_volcons.F it's pretty obvious what is going on. A cumulative summation of west,east,south,north volumes fluxes and AREAS is made. Any net volume flux is accommodated by adjusting the velocity on ALL boundaries by dividing the imbalance by the total area. If a user selects, say, only #define WEST_VOLCONS when both west and south are open, what does this mean? If there is a net inflow across the south that *mostly* goes out the west, this will be turned into a large net inflow when WEST_VOLCONS effectively cancels out the outflow.

Logically, I believe the calculation should be done such that the volume flux should always be over ALL boundaries, but the area summation and velocity adjustment be made over those #define xxxx_VOLCONS. The NORTH_, SOUTH_, WEST_, EAST_ VOLCONS options then indicate the user's choice of which boundaries to ADJUST to balance volume flux calculated for the complete domain.

One final point that some people miss so I'll mention it here: it makes no sense to apply the VOLCONS options in a regional model with TIDES because as the tide phase changes so to with the volume. While hope the sub-tidally frequency volume would hold constant, this is not the intent of the VOLCONS option.

John.

[Barbara]

Now I define the OBC of an idealized upwelling case in a rectangular grid （water depth is parallel to the E-W direction) as

EW_PERIODIC,
SOUTH_FSGRADIENT
SOUTH_M2RADIATION
SOUTH_M3RADIATION
SOUTH_TRADIATION

North is masked land.

Then I activate SOUTH_VOLCONS, but the sum of the zeta in the whole domain is NOT ZERO, except at t = 0.
(If the mass is conserved, the sum of zeta should be ZERO at any time, am I right?)

Thanks for your comments!

[wilkin]

See my earlier post to this thread.

There can be a problem when only some boundaries have the VOLCONS flag active. However, I would not have thought this is one of them. What goes in east comes out west, and north is closed, so only south contributes to the volume imbalance and distributing it over south only ought to work.

This is one of those cases where scrutinizing what is actually active in the f90 code might betray the problem.

You need to look at how the net volume transport is calculated, and then how it gets distributed across the open boundaries.

[Barbara]

Thanks for your reply!

In my idealized case, before I add SOUTH_VOLCONS, the summation of zeta in the whole domain is decreasing with time (mass lose). When I activate South_VOLCONS only or West, East, South _VOLCONS simultaneously, the summation of zeta is increasing with time (mass gain).

I have checked the obc_volcons.f90, and it looks fine for the 3 main steps,

Code: Select all

      IF (Jstr.eq.1) THEN
        DO i=Istr,Iend
          cff=0.5_r8*(zeta(i,Jstr-1,kinp)+h(i,Jstr-1)+                  &
     &                zeta(i,Jstr  ,kinp)+h(i,Jstr  ))*om_v(i,Jstr)
          cff=cff*vmask(i,Jstr)
          my_area=my_area+cff
          my_flux=my_flux+cff*vbar(i,JstrV-1,kinp)
        END DO
      END IF

However, I don't understand why define JstrV=Jstr+1 in the south boundary and use vbar(i,JstrV-1,kinp) here, instead of vbar(i,Jstr,kinp)?

Code: Select all

      IF (Jstr.eq.1) THEN
        DO i=-2+IstrU,MIN(Iend+1,Lm(ng))+1
          Dvom(i,Jstr)=0.5_r8*(Drhs(i,Jstr)+Drhs(i,Jstr-1))*            &
     &                 (vbar(i,Jstr,kinp)-ubar_xs)*                     &
     &                 om_v(i,Jstr)
          Dvom(i,Jstr)=Dvom(i,Jstr)*vmask(i,Jstr)
        END DO
      END IF

Code: Select all

      IF (Jstr.eq.1) THEN
        DO i=Istr,Iend
          vbar(i,Jstr,kinp)=(vbar(i,Jstr,kinp)-ubar_xs)
          vbar(i,Jstr,kinp)=vbar(i,Jstr,kinp)*vmask(i,Jstr)
        END DO
      END IF

Any suggestions on which f90 code to check, or if I should try other kinds of OBC?
Thanks a lot!

[LauraB]

Hi Barbara and hi all,

I'm having the same problem you described in October 15, 2007: An idealized configuration with periodic alongshore conditions and only one open boundary. Even when the VOLCONS option is on at the open boundary, total volume is not conserved. I also tried having a wall instead of an open boundary, and although in this case the NET_VOLUME in the output looks constant, the mean zeta at each time is not zero (though it is small, order 1e-4 meters).

I was wondering if you found the problem in the simulation you described last year. I'm especially concerned about this problem because my model does not conserve heat either, although there are no surface nor bottom fluxes. I thought it was due to my OBCs, but since it loses heat even when I have a wall instead of an open boundary, I think it might be related to the mass unbalance....

Thank you very much in advance,

Laura

[wilkin]

If you look in Nonlinear/obc_volcons.F you will see that, indeed, the logic therein is not adapted to the case of mixing periodic boundary conditions and open boundaries in an application. So, yes, your experience is consistent with what I see in the code, and the code needs modification to handle your situation.

The periodic part of the boundary flow needs no adjustment for volume conservation because it balances itself exactly. So what is required is:

(a) exclusion of any periodic open boundaries from the my_flux and my_area cumulative summations (probably simply by bracketing those sections with #ifdef WEST_VOLCONS and #ifdef EAST_VOLCONS inside a #if !defined EW_PERIODIC block (and similarly for north/south), and

(b) similarly restricting the velocity adjustment to the non-periodic boundaries in conserve_mass_tile.

In fact, the logic is not correct for the case of a user choosing to apply volume conservation adjustment on only some of the non-periodic open boundaries. (I confess this is something I changed myself long ago in my copy of version 2.3 and never communicated back to the forum.)

Consider the case that, e.g., a domain has 3 open boundaries west, south, east, but the user wants only to adjust south and east to achieve volume conservation. The flux through west must be included in the my_flux sum (as it is now) but the west area must not be added to my_area. This is because all nonperiodic boundaries must be integrated to determine the mass flux imbalance, but only the area of the boundaries to adjust should go into the calculation of the adjustment average velocity. Also, this adjustment must not be applied to the west boundary in this example.

John Wilkin.