WET-DRY points initialization, how they go from DRY to WET?

[Eugen]

Hello, everyone,

As a continuation of this topic (viewtopic.php?f=19&t=4423), we decided to output all model variables at the 'problematic' point of model’s explosion (during long-time simulations), in order to investigate why.

It seems, that temperature goes randomly (?) at the next time step after it becomes masked due to WET-DRY subroutine (see the picture below). But, it also seems, that after it becomes the wet-point back, it still has this random temperature value (!). See the picture below: x-axis is the time-step axis, blue/red lines – temperature at the surface/bottom, green – mask (shifts from 1 to 0 and then back around the time step with a number 4080, blue points – sea surface elevation).




Here, on the picture below, the model explosion is caught. The temperature value could go as crazy as -40 (usually, it shouldn’t be a reason to explode, but we set an alarm on this). We assume, that if temperature could go down -40, it also could go up to NaN values that will lead to temperature vanishing in the whole domain, like in this topic (viewtopic.php?f=19&t=4423).




The question is, how the wet-dry points are initialized as wet-points after they were dry-points for a while?

Thanks in advance.

[kate]

Yes, I've battled this one in the past. The temperatures have to be valid values through wet, dry and restart. The trick is in not letting it get weird when dry, whether through vertical processes or through lateral processes. I think we've turned off the vertical heat fluxes when dry. That leaves lateral processes - the rotated mixing tensors are poorly-behaved between a dry point at 3 m above sea level and a wet point with 30 m depth. We have a coastal domain with both fjords and shallow sills - the fjords are best being deep with always-dry walls while the sills can be wet-dry. Also, the model is running without explicit lateral diffusion using the third-order upwind advection scheme - turning off the diffusion solved this particular problem for me. Still, I spent a great deal of time manually fussing with the land mask and the bathymetry before I got to that stage.

[jcwarner]

we have also spent an effort to resolve similar issues. Are you using bulk_fluxes or imposing your own stresses? In bulk_fluxes we have the wet/dry masks multiplied so that the dry cells are not heating up.

Also, this needs to be set in pre_step for the light penetration. we added:
! Add in incoming solar radiation at interior W-points using decay
! decay penetration function based on Jerlov water type.
!
IF (itrc.eq.itemp) THEN
DO k=1,N(ng)-1
DO i=Istr,Iend
FC(i,k)=FC(i,k)+dt(ng)*srflx(i,j) &
# ifdef WET_DRY
& *rmask_wet(i,j) &
# endif
& *swdk(i,j,k)
END DO
END DO
END IF

It could also be the lateral mixing, that could be causing issues. can you turn off some cppdefs to see if it goes away?
-j

[Eugen]

Yes, I use BULK_FLUXES subroutine. And in its code I saw the multiplication on rmask_wet for longwave, sensible, latent and heat fluxes, and the wind stress.

Here is a complete list of subroutines, which I use:

• #define ROMS_MODEL
  #undef NESTING
  #define MASKING
  
  #define UV_ADV
  #define UV_SADVECTION
  #define UV_COR
  #define UV_VIS2
  
  #define MY25_MIXING
  #ifdef MY25_MIXING
  # define N2S2_HORAVG
  # define KANTHA_CLAYSON
  #endif
  
  #define MIX_S_UV
  #define TS_U3HADVECTION
  #define UV_LDRAG
  
  #define DJ_GRADPS
  #define NONLIN_EOS
  
  #define TS_DIF2
  #define MIX_S_TS
  
  #define DIFF_GRID
  #define VISC_GRID
  
  #define SALINITY
  #define SOLVE3D
  #define SPLINES_VVISC
  #define SPLINES_VDIFF
  #define AVERAGES
  
  #define BULK_FLUXES
  #define SHORTWAVE
  #define LONGWAVE
  #define CLOUDS
  #define DIUNAL_SRFLUX
  
  #define ANA_SSFLUX
  #define ANA_BSFLUX
  #define ANA_BTFLUX
  
  #define SSH_TIDES
  #define UV_TIDES
  #define RAMP_TIDES
  #define ADD_FSOBC
  #define ADD_M2OBC
  #define WET_DRY
  
  #define MPI_REAL8 MPI_DOUBLE_PRECISION

Could you guess, please, which subroutine could cause those problems?
I tried to switch off DIFF_GRID and TS_DIF2, but it seems, that temperature anomalies are still at the place...
And if so, what is the best way to effectively turn off the diffusion? Should it be done in a header file, or by setting a diffusion coefficient to zero in .in file?

[jcwarner]

can you put that rmask_wet in pre_step3d.F that i showed above?
-j

[Eugen]

can you put that rmask_wet in pre_step3d.F that i showed above?
-j

I'm sorry, but isn't it already there?

[kate]

#define ANA_SSFLUX

Isn't this inconsistent with BULK_FLUXES? It may be ignored, but I'm not sure.

#define MIX_S_TS

This may be more stable than the rotated mixing tensors, but it does cause cross-isopycnal mixing.

As for your question about whether to undefine TS_DIF2 or to set the coefficient to zero, both should work and the former will do fewer computations.

[Eugen]

Thank you, Kate !

Btw, we found two possible ways to overcome this problem.

First way: the subroutine LIMIT_STFLX_COOLING, which prevents temperature of going too low. It works fine, but it is a rude way.

Second way: the added a condition into the subroutine step3D_t to keep temperature the same, if a grid cell is temporally masked:

# ifdef WET_DRY
! Apply Land/Sea mask.
DO k=1,N(ng)
DO j=JstrR,JendR
DO i=IstrR,IendR
IF ((rmask_wet(i,j).eq.0).AND.(rmask(i,j).eq.1)) THEN
t(i,j,k,nnew,itrc)=t(i,j,k,nstp,itrc)
END IF
END DO
END DO
END DO
# endif

It seems, that it works fine with short-term simulations.
For example, it was before:



And this is after we introduced the piece of code into the subroutine: