Unstable with MY25_MIXING, but GLS works

General scientific issues regarding ROMS

Moderators: arango, robertson

Post Reply
Message
Author
chughes
Posts: 2
Joined: Mon May 29, 2017 4:22 pm
Location: Dalhousie University

Unstable with MY25_MIXING, but GLS works

#1 Post by chughes » Sun Jun 16, 2019 2:36 pm

Dear all,

I am looking to eventually setup a coupled ROMS-SWAN model within COAWST for an idealized hurricane in deep water with no land and a uniform bathymetry of 5000 m. However, before I can go on to the coupling, I have to get ROMS to work stand-alone. When I run the simulation using GLS_MIXING (either k-epsilon or kkl) it is stable and the results make sense. However, any time that I try to run a simulation with this particular setup using MY25_MIXING, the currents become unstable and this starts being visible after a simulation time of 12 hours. The simulation eventually blows up after 26.5 hours.

Oddly enough, when I try run ROMS with the same .h file and .in file for a different grid that also has no land and a uniform bathymetry of 5000 m, and a different wind setup, MY25_MIXING is fine.

The differences between the stable and unstable MY25_MIXING cases are the following:
-different wind (unstable uses an asymmetric wind with an inflow angle included, stable is axisymmetric tangential wind with no inflow angle)
-different initial S profiles (the unstable S profile is based on a WOA 2009 September climatology profile located near the idealized grid center (latitude = 20 degrees N, longitude = =-55 degrees) (the stable initial S is 35 psu everywhere)
-different initial T profiles (both from the same climatology), the unstable one is at the same location as the S profile, the stable T profile is a mid-latitude location with a much shallower MLD)
-unstable MY25 grid dimensions: 10 to 30 degrees latitude, -70 to -40 degrees longitude, 1/12 degree resolution, N=64
-stable MY25 grid dimensions: 6 to 34 degrees latitude, -59 to -41 degrees longitude, 1/10 degree resolution, N=32
-unstable hurricane moves westward at a constant latitude (20 degrees) with a translation speed of 4.5 m/s
-stable hurricane moves northward at a constant longitude (-51 degrees) with a translation speed of 10 m/s

Things that I have tried to fix the MY25_MIXING instability, but that have not changed it:
-reducing the 3D timestep from 60s to 30s and the barotropic time step from 3s (NDTFAST=20) to 1s (NDTFAST=30)
-reducing the number of vertical levels from N=64 to N=32
-setting the heat fluxes to zero (I added my own switch ZERO_HEATFLUX to bulk_flux.F)
-setting the wind stress to zero (by setting the drag in bulk_flux.F to zero via my own switch CD_ZERO)
-setting both the heat fluxes and wind stress to zero produced the null outcome (zero currents instead of being unstable)
-running the simulation on both the 2015 and 2019 versions of COAWST (ROMS stand-alone in both)

If the problem was the drag coefficient change that I made to bulk_flux.F (following the methodology of Oey et al., 2007), then GLS should be unstable and so should the results from the other grid. If the wind input is the problem, then GLS should be unstable as well and the heat fluxes only simulation should be stable.

Thank you for your time and consideration,
Colin

Here is my .h file (that I called roms.h):
/*
** svn $Id: inlet_test.h 838 2008-11-17 04:22:18Z jcwarner $
*******************************************************************************
** Copyright (c) 2002-2008 The ROMS/TOMS Group **
** Licensed under a MIT/X style license **
** See License_ROMS.txt **
*******************************************************************************
**
** Options for Inlet Test Case, waves-ocean (SWAN/ROMS) two-way coupling.
**
** Application flag: INLET_TEST
** Input script: ocean_inlet_test.in
** coupling_inlet_test.in
** sediment_inlet_test.in
*/

#define ROMS_MODEL

#define UV_VIS2
#define MIX_S_UV
#define MASKING
#define UV_ADV
#define UV_COR
#define TS_MPDATA
#define DJ_GRADPS
#define SPLINES_VDIFF
#define SPLINES_VVISC
#define RI_SPLINES

#define SOLVE3D
#undef ANA_INITIAL /*initial conditions*/
#undef ANA_SMFLUX /* analytical wind stress*/
#define ANA_FSOBC
#define ANA_M2OBC

/* define only one of the following 5 */
#define UV_LOGDRAG
#undef UV_QDRAG
#undef MB_BBL
#undef SG_BBL
#undef SSW_BBL
#ifdef SSW_BBL
# define SSW_CALC_ZNOT
#endif

#ifdef SOLVE3D
# define MY25_MIXING
# ifdef MY25_MIXING
# define N2S2_HORAVG
# define KANTHA_CLAYSON
# endif
# ifdef GLS_MIXING
# define KANTHA_CLAYSON
# define N2S2_HORAVG
# endif
#endif

# undef SEDIMENT
# ifdef SEDIMENT
# define SUSPLOAD
# undef BEDLOAD_SOULSBY
# undef BEDLOAD_MPM
# define SED_MORPH
# endif
# if defined SEDIMENT || defined SG_BBL || defined MB_BBL || defined SSW_BBL
# undef ANA_SEDIMENT
# endif

/* use if bulk fluxes computation */
# define BULK_FLUXES /* use if bulk fluxes computation */
# define ZERO_HEATFLUX /* use if setting heat flux to zero */
# undef BULK_FLUXES2D /* use if defining fluxes in 2D*/
# undef CCSM_FLUXES /* use if you want 2D Cd (Large and Yeager,2009)*/
# if defined BULK_FLUXES
# undef CD_Z12 /* use the Zijlema et al. 2012 Cd */
# define CD_Oey07 /* use the Oey et al., 2007 Cd*/
# undef CD_ZERO /* sets the drag to zeros */
# endif
# undef LONGWAVE_OUT /* use if computing outgoing longwave radiation */
# define ANA_BTFLUX
# define ANA_BSFLUX
/* surface fluxes expect the wind*/
# define ANA_HUMIDITY
# define ANA_PAIR
# define ANA_RAIN
# define ANA_SRFLUX
# define ANA_TAIR

# define WIND_MINUS_CURRENT /* add relative wind effect */

# define ATM_PRESS
# undef EMINUSP
# undef SOLAR_SOURCE

#define AVERAGES
#undef STATIONS
#define DIAGNOSTICS_UV
#define DIAGNOSTICS_TS

Note, the white dot is the storm center location and the dashed magenta line is the idealized storm track. (MWS = 50 m/s, RMW = 55 km)

Surface currents at simulation time 12 hours (zoomed in to see the subtle instabilities):
surface_currents_MY25_coawst2019_simtime_12hours.png
surface_currents_MY25_coawst2019_simtime_12hours.png (35.13 KiB) Viewed 411 times
Surface currents at simulation time 18 hours (entire domain):
surface_currents_MY25_coawst2019_simtime_18hours_fullgrid.png
surface_currents_MY25_coawst2019_simtime_18hours_fullgrid.png (38.52 KiB) Viewed 411 times
Surface currents at simulation time 12 hours (GLS kkl for comparison):
surface_currents_GLS_kkl_coawst2019_simtime_12hours.png
surface_currents_GLS_kkl_coawst2019_simtime_12hours.png (36.64 KiB) Viewed 411 times
Surface currents at simulation time 18 hours (GLS kkl for comparison, entire domain):
surface_currents_GLS_kkl_coawst2019_simtime_18hours_entire_domain.png
surface_currents_GLS_kkl_coawst2019_simtime_18hours_entire_domain.png (30.15 KiB) Viewed 411 times
.in file:
ocean_holland.in
(142.95 KiB) Downloaded 16 times
log file (a redirected output file):
sim.log
(359.65 KiB) Downloaded 19 times
bulk_flux.F:
bulk_flux.F
(72.98 KiB) Downloaded 17 times

jcwarner
Posts: 798
Joined: Wed Dec 31, 2003 6:16 pm
Location: USGS, USA

Re: Unstable with MY25_MIXING, but GLS works

#2 Post by jcwarner » Mon Jun 17, 2019 1:01 pm

it has been a while since I was in this part of the code, but one quick thought is that we have surface flux bc's for the gls term that is not included in the my25 method.

# if defined CHARNOK
Zos_eff(i)=MAX(charnok_alpha(ng)/g*0.5_r8* &
& SQRT((sustr(i,j)+sustr(i+1,j))**2+ &
& (svstr(i,j)+svstr(i,j+1))**2), &
& Zos_min)
# elif defined ZOS_HSIG
Zos_eff(i)=MAX(zos_hsig_alpha(ng)*Hwave(i,j), Zos_min)
# else
Zos_eff(i)=Zos_min
# endif

so the Zos_eff[ective] is set by the zo_surface in the ocean.in. If you look int Kumar's 2012 OM vortex force paper, we describe the several options for the surface fluxes of tke (and gls). Eqnts 43-47.

Also look in Sandro's paper:
Carniel, S., Warner, J.C., Chiggiato, J., Sclavo, M., 2009. Investigating the impact of
surface wave breaking on modelling the trajectories of drifters in the Northern Adriatic Sea during a wind-storm event. Ocean Model. 30, 225–239.
doi:10.1016/j.ocemod.2009.07.001.
Eqtns 8-13 for the gls bc.

-john

chughes
Posts: 2
Joined: Mon May 29, 2017 4:22 pm
Location: Dalhousie University

Re: Unstable with MY25_MIXING, but GLS works

#3 Post by chughes » Tue Jun 18, 2019 8:26 pm

Sorry for the trouble, but the problem was the initial temperature profile. I used a monthly climatology T profile and then added the bottom part from an annual climatology T profile. The result was an unstable in-situ temperature profile instead of an idealized potential temperature profile. GLS_MIXING did not have any trouble with this, but it turns out that MY25_MIXING is very sensitive to the initial temperature profile.

Here is a comparison of the initial T profile that led to MY25_MIXING being unstable vs my new corrected initial T profile.
Initial_T_comparison_entire_water_column_Colin_Grid_N64.png
Initial_T_comparison_entire_water_column_Colin_Grid_N64.png (21.22 KiB) Viewed 211 times
The new surface currents look much better and resemble GLS kkl. There would likely be more agreement if I ran another simulation of GLS kkl with the updated T profile.
surface_currents_MY25_coawst2019_simtime_12hours_fixed_initial_T_profile.png
surface_currents_MY25_coawst2019_simtime_12hours_fixed_initial_T_profile.png (62.42 KiB) Viewed 211 times

Post Reply