<div dir="ltr"><div>Hi MITgcmers,</div><div><br></div><div>I am trying to implement a simple configuration to model some nonlinear properties of internal gravity waves in the ocean later on. To do so, I am generating internal waves in a confined box (inviscid ocean with linear EOS in stratified temperature, no stratified salt) using specific boundary conditions. However, it seems that the total energy exhibits some spurious periodical time dependencies that are not supposed to appear according to the linear theory of internal gravity waves. My guess would be that I might get a numerical instability somewhere. Here is a detail of my issue.<br></div><div><br></div><div>My reference configuration consists of a 2D x-z stratified non-rotating domain with nx = 800; ny = 1, nz = 350; dx = dy = 30; dz = 10 (see below this mail for more details about the configuration extracted from my data file), with periodic boundary conditions along the x direction and w=0 at the top and bottom along the z direction, and with the following simulations parameters:</div><div><div>Nz = 1e-3; %buoyancy frequency<br>g = 9.81;<br>alpha = 2e-4; %thermal expansion coefficient<br></div><div>Tref = Nz^2 *(z+Hz/2) /(alpha*g); %stratified temperature<br></div><div>My simulations are run during a five day period (nTimeSteps =  28800, deltaT = 15.)<br></div> </div><div><br></div><div>The initial setup for the generated wave calls the following definitions (indicated in the gendata.m file):</div><div>nkx0 = 10; %integer for horizontal BC<br>nkz0 = 20; %integer for vertical BC<br>kx0 = 2*pi*nkx0/Lx; %periodic BC<br>kz0 = pi*nkz0/Hz; %imposed by w=0 at top and bottom<br>k0 = sqrt(kx0^2 + kz0^2); %norm of the wavevector of the generated internal gravity wave<br>omega0 = Nz*kx0/k0; %wave frequency</div><div><br></div><div>psi0 = sqrt(4*energyfactor)/k0 %wave amplitude with energyfactor = 10^(-4);<br></div><div><div>U0 = kz0*psi0; %x velocity amplitude<br></div><div>theta0 = - 1/(alpha*g) * (kx0/omega0) * Nz^2 * psi0; %potential temperature amplitude<br></div></div><div><br></div><div>In this box, the generated wave has initial condition (indicated in the gendata.m file):</div><div>U = U0.*cos(kz0*z).*cos(kx0*x)<br>V = zeros(nz,nx);<br>T = Tref + theta0.*sin(kz0*z).*cos(kx0*x) <br></div><div><br></div><div>The linear theory predicts a resulting stationary wave with time evolution: <br></div><div>U = U0.*cos(kz0*z).*cos(kx0*x-omega0*t)<br>V = zeros(nz,nx);</div><div>W = U0.*kx0./kz0.*sin(kz0*z).*sin(kx0*x-omega0*t)<br>T = Tref + theta0.*sin(kz0*z).*cos(kx0*x-omega0*t)</div><div>W being obtained using the continuity equation.<br></div><div><br></div><div>The mechanical energy of the wave is then given by:</div><div>E = KE + PE;</div><div>KE = 0.5*(U.^2+V.^2+W.^2) %kinetic energy;</div><div>PE = 0.5*(b/N)^2 %potential energy with b the buoyancy linked with the potential temperature by b= alpha.*g.*theta = alpha.*g.*theta0.*sin(kz0*z).*cos(kx0*x-omega0*t) (linear EOS);<br></div><div>So all calculations done, the total energy is equal to:</div><div>E = 0.5.*psi0.^2.*[kx0.^2.*sin(kz0*z).^2+kz0.^2.*cos(kx0*x-omega0*t)).^2]</div><div>Then averaging over the 2D space x-z:</div><div>e = <E> = (1/4).*k0.^2.*psi0.^2</div><div>Hence the total energy of the wave averaged over space e = <E> should not depend on time. The same commentary applies for the space mean kinetic and potential energy ke = <KE> and pe = <PE> and we have equipartition of kinetic and potential energy as ke = pe = (1/8).*k0.^2.*psi0.^2 (non rotating domain). <br></div><div><br></div><div>However, my simulations suggest a time dependency of (e,ke,pe), confirmed when looking at the Fourier transform of the mean space energy e (see screen capture <u>EnergyTimeSeries.png</u> and <u>EnergyFFT.png</u>). The equipartition is also affected even though it is respected in space+time average (but we are interested in the space average only). Can someone tell me what could be wrong apparently in my setup? Any ideas on how to correct that?<br></div><div><br></div><div>You will also find attached the time series of U, V, T and W (<u>timeSeries.png</u>) with their Fourier transform (<u>FFTtime</u>), and a video showing the time evolution of these quantities and of their space Fourier transform in a (kx,kz) domain (<u>FFTkxkzrun.avi</u>)<br></div><div><br></div><div><div>Thank you in advance,<br></div><div><br></div><div>Kevin<br></div></div><div><br></div><div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1UNsA7ut1bbj2-ZA-5rnXzpx7fByfZ3s0/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">EnergyFFT.png</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="opacity: 0.55; cursor: pointer; float: right; position: relative; top: -1px; display: none;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1r5ZvJwehTolIpQkdvWVw9-jJPormdMmJ/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">EnergyTimeSeries.png</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="opacity: 0.55; cursor: pointer; float: right; position: relative; top: -1px; display: none;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/1vi9cVrFJWvrp3PibZlpP61AiKxKTqbe3/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">FFTkxkzrun.avi</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="opacity: 0.55; cursor: pointer; float: right; position: relative; top: -1px; display: none;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:rgb(245,245,245);padding:5px;color:rgb(34,34,34);font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid rgb(221,221,221);line-height:1"><a href="https://drive.google.com/file/d/13dDtlq0GigcYrqOA4I8HwNaK-l-gtThQ/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:rgb(17,85,204);text-decoration:none;vertical-align:bottom">FFTtime.png</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="opacity: 0.55; cursor: pointer; float: right; position: relative; top: -1px; display: none;"></div><div class="gmail_chip gmail_drive_chip" style="width:396px;height:18px;max-height:18px;background-color:#f5f5f5;padding:5px;color:#222;font-family:arial;font-style:normal;font-weight:bold;font-size:13px;border:1px solid #ddd;line-height:1"><a href="https://drive.google.com/file/d/1-AuODYPDz24O13AkI6u9rigz8FRb4qoI/view?usp=drive_web" target="_blank" style="display:inline-block;max-width:366px;overflow:hidden;text-overflow:ellipsis;white-space:nowrap;text-decoration:none;padding:1px 0;border:none"><img style="vertical-align: bottom; border: none;" src="https://ssl.gstatic.com/docs/doclist/images/icon_10_generic_list.png"> <span dir="ltr" style="color:#15c;text-decoration:none;vertical-align:bottom">timeSeries.png</span></a><img src="//ssl.gstatic.com/ui/v1/icons/common/x_8px.png" style="display:none; opacity: 0.55; cursor: pointer; float: right; position: relative; top: -1px; "></div></div><div><br></div><div>#########<br></div><div><br></div><div>Config options taken from the data file:</div><div><br></div><div> no_slip_sides  = .TRUE.,<br> no_slip_bottom = .TRUE.,<br><br> viscAh  = 0.E-3,<br> viscAz  = 0.E-3,<br><br> diffKhT = 0.E-4,<br> diffKzT = 0.E-4,<br><br> f0   = 0.E-4,<br> beta = 0.E-11,<br><br> tAlpha = 2.E-4,<br> sBeta  = 0.E-4,<br><br> gravity = 9.81,<br> gBaro   = 9.81,<br><br> rigidLid = .FALSE.,<br> implicitFreeSurface=.TRUE.,<br> exactConserv = .TRUE.<br><br> nonHydrostatic = .TRUE.,<br> useSingleCpuIO = .TRUE.,<br><br> readBinaryPrec  = 64,<br> writeBinaryPrec = 64,<br> writeStatePrec  = 64,<br><br> staggerTimeStep=.TRUE.<br><br> saltStepping  = .FALSE.,<br><br> tempAdvScheme = 7, !7th Order One Step method with Monotonicity Preserving Limiter ???? GAD.h<br><br> hFacMin=0.1,<br><br> &<br><br># Elliptic solver parameters<br> &PARM02<br> cg2dMaxIters       =  10000,<br> cg2dTargetResidual = 1.E-14,<br><br> cg3dMaxIters       =    400,<br> cg3dTargetResidual = 1.E-14,<br> &<br><br># Time stepping parameters<br> &PARM03<br><br> niter0 = 0,<br><br> nTimeSteps  =  28800, <br> deltaT      =    15.,<br> dumpFreq    =  900.,<br> monitorFreq = 43200.,<br><br> abEps       =    0.1,<br> pChkptFreq  = 86400.,<br> chkptFreq   =    0.0,<br> &<br><br># Gridding parameters<br> &PARM04<br> usingCartesianGrid=.TRUE.,<br> usingSphericalPolarGrid=.FALSE.,</div></div>