Applying wind stresses at OBC boundaries

[dougiesquire]

@ashjbarnes has noticed strange behaviour occurring at the boundaries of all our NUOPC-coupled regional configurations. The video below shows the issue in a regional Tasmanian domain (left panel).

I think I've tracked the issue down to the application of surface wind stresses as a boundary condition at OBC boundary points in the MOM_vert_friction::vertvisc routine. The right panel in the video below (showing daily speed) is with zero stress boundary condition at OBC boundaries:

  if (CS%direct_stress) then
-   do j=G%jsc,G%jec ; do I=Isq,Ieq ; if (G%mask2dCu(I,j) > 0.) then
+   do j=G%jsc,G%jec ; do I=Isq,Ieq
      surface_stress(I,j) = 0.0
+     if (G%OBCmaskCu(I,j) > 0.) then
      zDS = 0.0
      stress = dt_Rho0 * forces%taux(I,j)
      do k=1,nz
        h_a = 0.5 * (h(i,j,k) + h(i+1,j,k)) + h_neglect
        hfr = 1.0 ; if ((zDS+h_a) > Hmix) hfr = (Hmix - zDS) / h_a
        u(I,j,k) = u(I,j,k) + I_Hmix * hfr * stress
        if (associated(ADp%du_dt_str)) ADp%du_dt_str(i,J,k) = (I_Hmix * hfr * stress) * Idt
        zDS = zDS + h_a ; if (zDS >= Hmix) exit
      enddo
    endif ; enddo ; enddo
  else
    do j=G%jsc,G%jec ; do I=Isq,Ieq
-     surface_stress(I,j) = dt_Rho0 * (G%mask2dCu(I,j)*forces%taux(I,j))
+     surface_stress(I,j) = dt_Rho0 * (G%OBCmaskCu(I,j)*forces%taux(I,j))
    enddo ; enddo
  endif

(and similarly for the v-component)

compare_iss.mp4

I'm not 100% confident that this is the best/right fix and I'm hoping to get conformation, or otherwise, from a MOM6 guru before implementing. E.g. possibly there are circumstances where the surface stresses should be applied at OBC boundaries? @Hallberg-NOAA maybe you can advise?

Also, there are a couple of other places in MOM6 src where the wind stresses are used:

• In MOM_vert_friction::vertFPmix when FPMIX = True
• In MOM_set_viscosity::set_viscous_ML when DYNAMIC_VISCOUS_ML = True

We don't use these in any of our configurations at the moment, but do these also require special consideration at OBC boundaries?

[dougiesquire]

Thanks to @ashjbarnes we also have an FMS-coupled version of the same configuration. Curiously this does not show the same boundary issue because the easternmost/westernmost columns of forces%taux and northernmost/southernmost rows of forces%tauy are zero:

The zeros along the eastern and northern edges are already present in the ice_ocean_boundary_type u_flux and v_flux fields passed from the coupler. The zeros along the western and southern edges are present because this symmetric edge is never filled with data when using WIND_STAGGER = "C".

I wonder if this is why this issue hasn't been noticed before (I'm assuming most coupled regional configurations out there use the FMScoupler)?

[access-hive-bot]

This issue has been mentioned on ACCESS Hive Community Forum. There might be relevant details there:

https://forum.access-hive.org.au/t/access-om3-development-highlights-updates/5902/9

[dougiesquire]

I've added this fix in #55, but I'll leave this issue open because there are still questions about whether this is the best approach and whether other changes might also be needed

[sanAkel]

The zeros along the eastern and northern edges are already present in the ice_ocean_boundary_type u_flux and v_flux fields passed from the coupler. The zeros along the western and southern edges are present because this symmetric edge is never filled with data when using WIND_STAGGER = "C".

I wonder if this is why this issue hasn't been noticed before (I'm assuming most coupled regional configurations out there use the FMScoupler)?

We do not use the FMS coupler, e.g., the regional HAFS uses CMEPS coupler and the MOM6 configuration now uses A-stagger for both winds (input) and currents (output); it used to have C-stagger, which was recently corrected from the global UFS-weather-model.

[dougiesquire]

Our OBC settings:

#override OBC_NUMBER_OF_SEGMENTS = 4 !None
#override OBC_FREESLIP_VORTICITY = False !None
#override OBC_FREESLIP_STRAIN = False !None
#override OBC_COMPUTED_VORTICITY = True !None
#override OBC_COMPUTED_STRAIN = True !None
#override OBC_ZERO_BIHARMONIC = True !None
#override OBC_TRACER_RESERVOIR_LENGTH_SCALE_OUT = 3.0E+04 !None
#override OBC_TRACER_RESERVOIR_LENGTH_SCALE_IN = 3000.0 !None
#override BRUSHCUTTER_MODE = True !None
#override OBC_SEGMENT_001 = "J=N,I=N:0,FLATHER,ORLANSKI,NUDGED,ORLANSKI_TAN,NUDGED_TAN" !None
#override OBC_SEGMENT_001_VELOCITY_NUDGING_TIMESCALES = 0.3, 360.0 !None
#override OBC_SEGMENT_001_DATA = "U=file:forcing_obc_segment_001.nc(u),V=file:forcing_obc_segment_001.nc(v),SSH=file:forcing_obc_segment_001.nc(eta),TEMP=file:forcing_obc_segment_001.nc(temp),SALT=file:forcing_obc_segment_001.nc(salt)" !None
#override OBC_SEGMENT_002 = "J=0,I=0:N,FLATHER,ORLANSKI,NUDGED,ORLANSKI_TAN,NUDGED_TAN" !None
#override OBC_SEGMENT_002_VELOCITY_NUDGING_TIMESCALES = 0.3, 360.0 !None
#override OBC_SEGMENT_002_DATA = "U=file:forcing_obc_segment_002.nc(u),V=file:forcing_obc_segment_002.nc(v),SSH=file:forcing_obc_segment_002.nc(eta),TEMP=file:forcing_obc_segment_002.nc(temp),SALT=file:forcing_obc_segment_002.nc(salt)" !None
#override OBC_SEGMENT_003 = "I=N,J=0:N,FLATHER,ORLANSKI,NUDGED,ORLANSKI_TAN,NUDGED_TAN" !None
#override OBC_SEGMENT_003_VELOCITY_NUDGING_TIMESCALES = 0.3, 360.0 !None
#override OBC_SEGMENT_003_DATA = "U=file:forcing_obc_segment_003.nc(u),V=file:forcing_obc_segment_003.nc(v),SSH=file:forcing_obc_segment_003.nc(eta),TEMP=file:forcing_obc_segment_003.nc(temp),SALT=file:forcing_obc_segment_003.nc(salt)" !None
#override OBC_SEGMENT_004 = "I=0,J=N:0,FLATHER,ORLANSKI,NUDGED,ORLANSKI_TAN,NUDGED_TAN" !None
#override OBC_SEGMENT_004_VELOCITY_NUDGING_TIMESCALES = 0.3, 360.0 !None
#override OBC_SEGMENT_004_DATA = "U=file:forcing_obc_segment_004.nc(u),V=file:forcing_obc_segment_004.nc(v),SSH=file:forcing_obc_segment_004.nc(eta),TEMP=file:forcing_obc_segment_004.nc(temp),SALT=file:forcing_obc_segment_004.nc(salt)" !None

[dougiesquire]

We do not use the FMS coupler, e.g., the regional HAFS uses CMEPS coupler and the MOM6 configuration now uses A-stagger for both winds (input) and currents (output); it used to have C-stagger, which was recently corrected from the global UFS-weather-model.

Thanks @sanAkel. I don't think this is due to an incorrect stagger setting in this case. In both our NUOPC- and FMS-coupled cases, the stagger is hardcoded to the correct value in the driver (here for NUOPC, here for FMS)

[sanAkel]

We do not use the FMS coupler, e.g., the regional HAFS uses CMEPS coupler and the MOM6 configuration now uses A-stagger for both winds (input) and currents (output); it used to have C-stagger, which was recently corrected from the global UFS-weather-model.

Thanks @sanAkel. I don't think this is due to an incorrect stagger setting in this case. In both our NUOPC- and FMS-coupled cases, the stagger is hardcoded to the correct value in the driver (here for NUOPC, here for FMS)

👍 got it @dougiesquire. That line sets wind stagger to A. But what about the staggering of the surface currents? From what I understood, NUOPC/ CMEPS works (correct!) only if the currents are also A.

[dougiesquire]

👍 got it @dougiesquire. That line sets wind stagger to A. But what about the staggering of the surface currents? From what I understood, NUOPC/ CMEPS works (correct!) only if the currents are also A.

Yup, we use OCEAN_SURFACE_STAGGER = "A" in the NUOPC-coupled case. We have OCEAN_SURFACE_STAGGER = "C" in the FMS-coupled configuration, which I think is correct, though the FMS coupler is a little opaque to me.

[sanAkel]

👍 got it @dougiesquire. That line sets wind stagger to A. But what about the staggering of the surface currents? From what I understood, NUOPC/ CMEPS works (correct!) only if the currents are also A.

Yup, we use OCEAN_SURFACE_STAGGER = "A" in the NUOPC-coupled case. We have OCEAN_SURFACE_STAGGER = "C" in the FMS-coupled configuration, which I think is correct, though the FMS coupler is a little opaque to me.

Good to confirm, sounds like stagger is not the case of the problem!

[manishvenu]

I thought I'd add that this PR indirectly fixed an exact OBC restart issue we were having. I don't understand it that well, but the direct stress lines of code here are reading from the halo points, which aren't saved during restarts.