very diferent velocitities using GLS,MY25 or LMD

[flcastej]

Dear all,

I am running a nested simulation. It have 3 child grids for 3 channels with a resolution about 10 meters. The channels are quite narrow (50/70 m width) and shallow (5 meters). I have used different vertical mixing schemes: MY25, GLS and LMD. If I use GLS or LMD my model is very stable but I get low current velocities compared to the one measured by the ADCP deployed in the area but if I use MY25 my velocities are more realistic but my model get unstable. I have been doing different test but I always get the same results (MY25 more realistic but unstable). I will really appreciated any help or suggestion.

Thanks a lot.
My configuration used are :


#define STATIONS
#undef FLOATS
#undef DIAGNOSTICS_UV
#define MASKING
#define SOLVE3D
#undef ONE_WAY

#define UV_ADV
#define UV_COR
#define UV_VIS2
#define UV_QDRAG
#define DJ_GRADPS

#define CURVGRID
#define MIX_S_UV

#define TS_U3HADVECTION
#define TS_C4VADVECTION
#define SALINITY
#define NONLIN_EOS

#define MY25_MIXING
# ifdef MY25_MIXING
# undef N2S2_HORAVG
# define KANTHA_CLAYSON
# undef UV_LOGDRAG
# define K_C2ADVECTION
# undef K_C4ADVECTION
# endif

#undef GLS_MIXING
#if defined GLS_MIXING
# define KANTHA_CLAYSON
# define N2S2_HORAVG
# define CRAIG_BANNER
# define CHARNOK
# undef ZOS_HSIG
# undef TKE_WAVEDISS
#endif

#undef LMD_MIXING
#ifdef LMD_MIXING
# define LMD_CONVEC
#define LMD_DDMIX
# define LMD_SKPP
# define LMD_BKPP
# define LMD_NONLOCAL
# define LMD_RIMIX
#endif

#define ANA_BTFLUX
#define ANA_BSFLUX

#define SSH_TIDES
#define UV_TIDES
#define ADD_FSOBC
#define ADD_M2OBC

#define BULK_FLUXES
#ifdef BULK_FLUXES
# define CLOUDS
# define ALBEDO
# define DIURNAL_SRFLUX
# define LONGWAVE
# define ANA_RAIN
# undef ANA_PAIR
# define ANA_CLOUD
# define ANA_SST
# undef ANA_TAIR
# undef ANA_HUMID
# define EMINUSP
#endif

# undef ANA_DRAG
# define UV_DRAG_GRID

Initial file:

! Model iteration loops parameters.

ERstr = 1
ERend = 1
Nouter = 1
Ninner = 1
Nintervals = 1

! Number of eigenvalues (NEV) and eigenvectors (NCV) to compute for the
! Lanczos/Arnoldi problem in the Generalized Stability Theory (GST)
! analysis. NCV must be greater than NEV (see documentation below).

NEV = 2 ! Number of eigenvalues
NCV = 10 ! Number of eigenvectors



LcycleTLM == F
NTLM == 72
NDEFTLM == 0
LcycleADJ == F
NADJ == 72
NDEFADJ == 0

LrstGST = F ! GST restart switch
MaxIterGST = 500 ! maximun number of iterations
NGST = 10 ! check pointing interval


Ritz_tol = 1.0d-15



TNU2 == 2*0.0d0 2*0.0d0 2*0.0d0 2*0.0d0 ! m2/s
TNU4 == 2*0.0d0 2*0.0d0 2*0.0d0 2*0.0d0 ! m4/s


VISC2 == 4*0.001d0 ! m2/s
VISC4 == 4*0.0d0 ! m4/s


AKT_BAK == 4*1.0d-6 4*5.0d-6 ! m2/s

AKV_BAK == 4*1.0d-2 4*1.0d-2 ! m2/s


AKK_BAK == 4*5.0d-6 ! m2/s
AKP_BAK == 4*5.0d-6 ! m2/s
TKENU2 == 4*0.0d0 ! m2/s
TKENU4 == 4*0.0d0 ! m4/s


GLS_P == 4*3.0d0 ! K-epsilon
GLS_M == 4*1.5d0
GLS_N == 4*-1.0d0
GLS_Kmin == 4*7.6d-6
GLS_Pmin == 4*1.0d-12

GLS_CMU0 == 4*0.5477d0
GLS_C1 == 4*1.44d0
GLS_C2 == 4*1.92d0
GLS_C3M == 4*-0.4d0
GLS_C3P == 4*1.0d0
GLS_SIGK == 4*1.0d0
GLS_SIGP == 4*1.30d0
RDRG == 4*3.0d-04 ! m/s
RDRG2 == 4*0.003d0 ! nondimensional
Zob == 4*0.0d0 ! m
Zos == 4*00.0d0 ! m



BLK_ZQ == 4*10.0d0 ! air humidity
BLK_ZT == 4*10.0d0 ! air temperature
BLK_ZW == 4*5.0d0 ! winds


DCRIT == 4*0.10d0 ! m


WTYPE == 4*1
LEVSFRC == 4*15
LEVBFRC == 4*1



Vtransform == 4*2 ! transformation equation
Vstretching == 4*3 ! stretching function


THETA_S == 4*7.0d0 ! surface stretching parameter
THETA_B == 4*0.8d0 ! bottom stretching parameter
TCLINE == 4*4.0d0 ! critical depth (m)


RHO0 = 1026.0d0 ! kg/m3
BVF_BAK = 1.0d-5 ! 1/s2


DSTART = 15622 ! days
TIDE_START = 0.0d0 ! days
TIME_REF = 19680523.00 ! yyyymmdd.dd
! TIME_REF = 00010101.00 ! yyyymmdd.dd


! TNUDG == 4*0.0d0 ! days
! ZNUDG == 4*0.01d0 ! days
! M2NUDG == 4*0.01d0 ! days
! M3NUDG == 4*0.01d0 ! days


! OBCFAC == 1.0d0 ! nondimensional




LBC(isFsur) == Clo Gra RadNud RadNud \ ! free-surface, Grid 1
Nes Nes Nes Nes \ ! free-surface, Grid 2
Nes Nes Nes Nes \ ! free-surface, Grid 3
Nes Nes Nes Nes ! free-surface, Grid 4
LBC(isUbar) == Clo FLa Fla Fla \ ! 2D U-momentum, Grid 1
Nes Nes Nes Nes \ ! 2D U-momentum, Grid 2
Nes Nes Nes Nes \ ! 2D U-momentum, Grid 3
Nes Nes Nes Nes ! 2D U-momentum, Grid 4
LBC(isVbar) == Clo FLa Fla Fla \ ! 2D V-momentum, Grid 1
Nes Nes Nes Nes \ ! 2D V-momentum, Grid 2
Nes Nes Nes Nes \ ! 2D V-momentum, Grid 3
Nes Nes Nes Nes ! 2D V-momentum, Grid 4
LBC(isUvel) == Clo Rad Rad Rad \ ! 3D U-momentum, Grid 1
Nes Nes Nes Nes \ ! 3D U-momentum, Grid 2
Nes Nes Nes Nes \ ! 3D U-momentum, Grid 3
Nes Nes Nes Nes ! 3D U-momentum, Grid 4
LBC(isVvel) == Clo Rad Rad Rad \ ! 3D V-momentum, Grid 1
Nes Nes Nes Nes \ ! 3D V-momentum, Grid 2
Nes Nes Nes Nes \ ! 3D V-momentum, Grid 3
Nes Nes Nes Nes ! 3D V-momentum, Grid 4
LBC(isMtke) == Clo Rad Rad Rad \ ! mixing TKE, Grid 1
Nes Nes Nes Nes \ ! mixing TKE, Grid 2
Nes Nes Nes Nes \ ! mixing TKE, Grid 3
Nes Nes Nes Nes ! mixing TKE, Grid 4

LBC(isTvar) == Clo Gra Gra Gra \ ! temperature, Grid 1
Clo Gra Gra Gra \ ! salinity, Grid 1
Nes Nes Nes Nes \ ! temperature, Grid 2
Nes Nes Nes Nes \ ! salinity, Grid 2
Nes Nes Nes Nes \ ! temperature, Grid 3
Nes Nes Nes Nes \ ! salinity, Grid 3
Nes Nes Nes Nes \ ! temperature, Grid 4
Nes Nes Nes Nes ! salinity, Grid 4




! Linear equation of State parameters:

R0 == 4*1026.0d0 ! kg/m3
T0 == 4*17.0d0 ! Celsius
S0 == 4*37.0d0 ! nondimensional
TCOEF == 4*1.7d-4 ! 1/Celsius
SCOEF == 4*7.6d-4

! Slipperiness parameter: 1.0 (free slip) or -1.0 (no slip)

GAMMA2 == 4*1.0d0

! Starting (DstrS) and ending (DendS) day for adjoint sensitivity forcing.
! DstrS must be less or equal to DendS. If both values are zero, their
! values are reset internally to the full range of the adjoint integration.

DstrS == 0.0d0 ! starting day
DendS == 0.0d0

[kate]

When it goes unstable, have you tried changing the timestep to something shorter? What is the nature of the instability?

[flcastej]

Dear Kate,

Thanks a lot for your answer and your help. The model get unstable after 4 hours, the wind blowing at this time is not very strong (5 m/s) and the sea level variations quite small. I have tried to use smaller time steps, different values for RDRG2 and several sigma layer distributions but I always get the same results. If I use other vertical schemes (GLS or LMD) my model is stable but the effect of the wind is very small and I don’t get realistic values for the current. I am only able to get them if I use the MY25 scheme. Below you will find two plots: the bathymetry and the u-component where the instability can be seen clearly.

[lanerolle]

Your model bathymetry in general is quite shallow and it blows up in a region where the water depth < 2m. Does the model have tides? If so what is the tidal range? Have you tried switching on the wetting-drying in ROMS?

[flcastej]

Dear Lyon Lanerolle,

Thanks a lot for your answer I am working in a quite shallow area. I have tides but the tidal range is very small (max. 20 cm) so the area never get dry.

Regards,

-Francisco

[lanerolle]

So the domain of computation is deeper than 20 cm (0.2 m) everywhere and anything shallower is under the permanent dry mask (mask_rho, mask_u, mask_v)? If possible, could you please show us a plot of your permanent dry mask (say mask_rho)?

[flcastej]

The white area in the pictures showed in my post correspond to the mask_rho.

Yes, the minimum depth in my domain is 0.5 meters to ensure that I am not going to have any dry area.

My problem come using the MY25 vertical scheme, if I use other scheme my model is stable but I get low velocities.

Thanks a lot for you help.

[lanerolle]

I apply my GLS and MY2.5 schemes in ROMS as follows and I find that the GLS & MY2.5 give similar results with the exception that the former provides better stratification and is usually ~40% more expensive computationally.

#if defined MY25_MIXING || defined GLS_MIXING
# define N2S2_HORAVG
# define KANTHA_CLAYSON
#endif

As you can see, I always use the N2S2_HORAVG option but I do not specify an advection scheme choice for the mixing scheme - I let ROMS pick the default.

[flcastej]

Thanks a lot for your suggestion, I will try to do so.

[kate]

You do appear to need some more lateral smoothing. Another option would be to increase:

Code: Select all

VISC2 == 4*0.001d0 ! m2/s

at least for the coarser grids. You don't show us where and how fine each of your four grids are.

[flcastej]

Dear Kate,

Thanks a lot for your help. Attached you will find a plot of the coarser grid and the three nested grids. Two of the nested grids are refined to a scale of 3 and another one to 4. In the plot is showed this values. I have changed the values of VISC2 increasing the value for the coarser grid and re-launch the model.



Dear Lyon,

I have launched the model with your configuration but I am getting the same problem.

[lanerolle]

If MY2.5 still blows up and GLS does not and gives diminished currents then the only thing I can conclude is that the bottom stress is coming into the two turbulent viscosity formulations differently - see my_corstep.F for example where (bustr, bvstr) and applied at the bottom most layer as a boundary condition. So it appears that the bottom stress which comes into the MY2.5 is less and causes the model to become unstable whereas what comes into the GLS scheme is at a lower level and allows the model to be numerically stable.

Have you tried both the quadratic (e.g. Cd=0.005) and logarithmic (e.g. Zob=0.005) bottom stress formulations? Do they both cause MY2.5 to blow-up?

[flcastej]

Dear Lyon,

Thanks a lot for your advice. I have launched the model following Kate suggestion and changing the value of VISC2 and it seem that it´s working. Once it finished I will follow your suggestion and test to run the model with different bottom stress formulation.

[flcastej]

Dear Kate and Lyon,

Thanks a lot for your help. Finally it seem that I have been able to get rid off the instability trying different values for VISC2 and Zob. Now I have another problem but I think that I am going to be able to handle it by my own.

Thanks !!