[MITgcm-support] Help me get rid of my noise!
Ryan Abernathey
ryan.abernathey at gmail.com
Wed May 7 14:43:09 EDT 2014
Thanks for all the helpful suggestions everyone!
I am now convinced that this isn't as "unphysical" as I thought.
In response to Christopher's suggestions, does anyone know whether I should
have convective adjustment active together with GGL90? I know for KPP you
are supposed to turn off all convective adjustment. Is that also the case
for GGL?
Thanks again,
Ryan
On Wed, May 7, 2014 at 2:10 PM, Christopher L. P. Wolfe <
christopher.wolfe at stonybrook.edu> wrote:
>
> My guess that you’re seeing hydrostatic convection. You don’t seem to have
> set either cAdjFreq or ivdc_kappa and my understanding is that they both
> default to zero. So, your model’s trying to do convection through the
> momentum equations and the convective scale is smaller than your grid
> spacing (hydrostatic convection runs to the grid scale).
>
> Christopher
>
>
> On May 7, 2014, at 1:20 PM, Ryan Abernathey <ryan.abernathey at gmail.com>
> wrote:
>
> I guess my image didn't make it through on that last email.
>
>
>
> On Wed, May 7, 2014 at 1:11 PM, Ryan Abernathey <ryan.abernathey at gmail.com
> > wrote:
>
>> Hello MITgcm Support Group,
>>
>> I am running a 2D model (depth and latitude) in a Southern Ocean
>> configuration. I am trying to eliminate what appears to be unphysical
>> grid-scale noise in my W field--columns of alternately very high and low W.
>> This is a snapshot of W after 5 years of spinup. This pattern is steady in
>> time.
>> As you can see, the noise is concentrated mainly (but not exclusively)
>> over the shelf break. A possibly related issue is that my cg2d solver is
>> relatively far from convergence:
>> cg2d: Sum(rhs),rhsMax = -5.46186194609888E-05 6.64557830841732E-02
>>
>> My model is 0.5 degree resolution in latitude, meaning that it is
>> relatively high resolution compared to the tutorials (e.g. ideal_2D_oce),
>> but of course it has no eddies or any temporal variability at all, since it
>> is 2D. This puts it in an unknown [at least to me] parameter regime.
>>
>> I have tried cranking the lateral viscosity (both laplacian and
>> biharmonic) up as high as they will go. I have enabled both linear and
>> quadratic bottom drag. I am NOT using the CD scheme; the noise was much
>> worse with it enabled. I am NOT using KPP but rather GGL90. (I tried with
>> KPP and it was not very different.)
>>
>> I would appreciate any suggestions you have about how to get rid of this
>> noise and achieve smooth fields.
>>
>> Thanks,
>> Ryan
>>
>> ####### data file ##########
>> &PARM01
>> useAreaViscLength=.FALSE.,
>> viscAr=1.0E-03,
>> viscA4=1e12,
>> viscAh=1e5,
>> tempAdvScheme=77,
>> # KrT fields are read in from a file
>> diffKrT=0.0E-5,
>> diffKhT=0.0,
>> diffK4T=0.0,
>> staggerTimeStep=.TRUE.,
>> # initial vertical profiles of T and S
>> tRef=44*0.0,
>> sRef=44*35.0,
>> # equation of state
>> eosType='LINEAR',
>> rhonil=1035.,
>> eosType='LINEAR',
>> tAlpha=2.000000E-04,
>> sBeta=0.000000E+00,
>> saltStepping=.FALSE.,
>> # boundary conditions
>> no_slip_sides=.FALSE.,
>> no_slip_bottom=.TRUE.,
>> # needed for points above the U sponge layer (shelf)
>> bottomDragLinear=1.E-3,
>> bottomDragQuadratic = 0.0012
>> # additional drag is done with RBCS
>> # physical parameters
>> gravity=9.810000E+00,
>> # implicit diffusion and convective adjustment
>> implicitDiffusion=.TRUE.,
>> implicitViscosity=.TRUE.,
>> implicitFreeSurface=.TRUE.,
>> # exact volume conservation
>> # exactConserv=.TRUE.,
>> # C-V scheme for Coriolis term
>> useCDscheme=.FALSE.,
>> # partial cells for smooth topography
>> hFacMin=5.000000E-02,
>> # file IO stuff
>> readBinaryPrec=64,
>> useSingleCpuIO=.TRUE.,
>> debugLevel=1,
>> &
>> # elliptic solver parameters
>> &PARM02
>> cg2dMaxIters=500,
>> cg2dTargetResidual=1E-09,
>> &
>> # timestepping parameters
>> &PARM03
>> nIter0=0000172800,
>> # nTimeSteps=1000,
>> # nTimeSteps=12961,
>> nTimeSteps=345600,
>> # deltaT=3600.0,
>> # deltaT=2400.,
>> # deltaT=1800.0,
>> deltaT=900.0,
>> abEps=0.1,
>> # needed for sea ice
>> forcing_In_AB = .FALSE.,
>> ChkptFreq=31104000.
>> pChkptFreq=311040000.,
>> taveFreq=311040000.0,
>> dumpFreq=31104000.,
>> monitorFreq=2592000.0,
>> tauThetaClimRelax=2592000.0,
>> dumpInitAndLast=.TRUE.,
>> pickupStrictlyMatch=.FALSE.,
>> &
>> # gridding parameters
>> &PARM04
>> usingSphericalPolarGrid=.TRUE.,
>> delX=0.5,
>> delY=100*0.5,
>> xgOrigin=0.,
>> ygOrigin=-80.,
>> delR=10., 10., 10., 10., 10., 10., 10., 12., 14.,
>> 16., 19., 22., 26., 30., 36., 42., 50., 60.,
>> 72., 85., 100., 120., 140., 166., 200., 200., 200.,
>> 200., 200., 200., 200., 200., 200., 200., 200., 200.,
>> 200., 200., 200., 200., 200., 200., 200., 200.
>> &
>> # Input datasets
>> &PARM05
>> bathyFile='bathyFile_shelf_2Dchan.bin',
>> zonalWindFile='zonalWindFile.bin',
>> thetaClimFile='thetaClimFile_WOA98_2Dchan.bin',
>> diffKrFile='diffKrFile_nosponge_2Dchan.bin',
>> hydrogThetaFile='hydrogThetaFile_WOA98_2Dchan.bin',
>> hydrogSaltFile='hydrogSaltFile_WOA98_2Dchan.bin',
>> &
>>
>> ########## data.pkg ##############
>> &PACKAGES
>> useMNC = .TRUE.,
>> useGMREDI = .TRUE.,
>> useGGL90 = .TRUE.,
>> useRBCS = .TRUE.,
>> &
>>
>> ######## data.gmredi ##############
>> &GM_PARM01
>> GM_MNC = .FALSE,
>> GM_AdvForm = .TRUE.,
>> GM_UseBVP = .TRUE.,
>> GM_background_K = 2000,
>> GM_isopycK = 2000,
>> GM_Kmin_horiz = 50,
>> GM_taper_scheme = 'dm95',
>> GM_BVP_ModeNumber = 1,
>> GM_BVP_cMin = .1,
>> &end
>>
>> ########### data.ggl90 ############
>> &GGL90_PARM01
>> GGL90writeState=.FALSE.,
>> GGL90TKEmin = 1.e-7,
>> mxlMaxFlag=2,
>> GGL90mixingLengthMin=3.,
>> &
>>
>>
>> ############ SIZE.h ##############
>> ...
>> PARAMETER (
>> & sNx = 1,
>> & sNy = 100,
>> & OLx = 4,
>> & OLy = 4,
>> & nSx = 1,
>> & nSy = 1,
>> & nPx = 1,
>> & nPy = 1,
>> & Nx = sNx*nSx*nPx,
>> & Ny = sNy*nSy*nPy,
>> & Nr = 44)
>> ...
>>
>>
> <W.png>_______________________________________________
>
> MITgcm-support mailing list
> MITgcm-support at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-support
>
>
>
> _______________________________________________
> MITgcm-support mailing list
> MITgcm-support at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-support
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140507/860def40/attachment.htm>
More information about the MITgcm-support
mailing list