[MITgcm-support] Warning and time stepping issue
Anthony Coletti
ajcolett at geo.umass.edu
Thu Jan 16 12:16:31 EST 2014
Martin,
Thank you! That worked! I thought using the default deltaT variable would work just as well, but I was wrong. Setting them individually did the trick.
Thank you all.
Anthony
Anthony J. Coletti
Climate System Research Center
Department of Geosciences
Morrill Building
611 N. Pleasant Street
233 Morrill Science Center
University of Massachusetts-Amherst
paleoclimate.org
Email: ajcolett at geo.umass.edu
http://blogs.umass.edu/ajcolett/
“For me, I am driven by two main philosophies: know more today about the world than I knew yesterday and lessen the suffering of others. You'd be surprised how far that gets you.” ― Neil deGrasse Tyson
On Jan 16, 2014, at 11:19 AM, Martin Losch <Martin.Losch at awi.de> wrote:
> Hi Anthony,
>
> I suggestion to have a look at verification/global_ocean.cs32x15 again and try these parameters:
> deltaTmom =1200.,
> deltaTtracer=86400.,
> deltaTfreesurf=86400.,
> deltaTClock =86400.,
> instead of deltaT = 5000.,This will speed up the integration (again at the cost of slowed down gravity waves etc. I refer you to the literature about asynchronous time stepping or “tracer acceleration”)
>
> Martin
>
> On Jan 16, 2014, at 11:14 AM, Anthony Coletti <ajcolett at geo.umass.edu> wrote:
>
>> Hi Christoph and Martin,
>>
>> Thank you for your responses!
>>
>> @Martin: I am looking through the manual now to see if I can tweak my run to make it run faster (5000s time step is unbearably slow for such course resolution).
>>
>> @Christoph: My ice thicknesses are not out of the ordinary at all. I think the thickest amount of ice is about 4 meters thick.
>>
>> Anthony
>>
>>
>> Anthony J. Coletti
>> Climate System Research Center
>> Department of Geosciences
>> Morrill Building
>> 611 N. Pleasant Street
>> 233 Morrill Science Center
>> University of Massachusetts-Amherst
>> paleoclimate.org
>> Email: ajcolett at geo.umass.edu
>> http://blogs.umass.edu/ajcolett/
>>
>> “For me, I am driven by two main philosophies: know more today about the world than I knew yesterday and lessen the suffering of others. You'd be surprised how far that gets you.” ― Neil deGrasse Tyson
>>
>>
>>
>>
>> On Jan 16, 2014, at 10:59 AM, Christoph Voelker <christoph.voelker at awi.de> wrote:
>>
>>> Hi Anthony and Martin,
>>>
>>> I suspect thta it is not the CFL criterium, but the sea-ice in
>>> conjunction with the nonlinear free surface, under extreme cases of ice
>>> formation. I had a similar isue in some last glacial maximum model runs
>>> that I did, where the free surface code basically complained that the
>>> submerged part of the sea ice got too thick. When that happened I also
>>> got very large negative eta values. There was a not particulaly elegant
>>> trick around it that I don't remember at the moment, I'll check. But
>>> Anthony, maybe you can check your ice thickness?
>>>
>>> Cheers, Christoph
>>>
>>> On 1/16/14 4:41 PM, Martin Losch wrote:
>>>> Hi Anthony,
>>>>
>>>> I usually find the stability “analysis” in the tutorials quite illustrative. For example, for the 4deg global ocean experiment <http://mitgcm.org/public/r2_manual/latest/online_documents/node129.html>, you can find a couple of stability criteria: CFL for advection and viscosity, but also for inertial oscillations. The Coriolis is treated explicitly and I believe, that this limits your time step for the momentum equations. The only way around it is to use asynchronous time stepping (at the usual cost of distorted physics).
>>>>
>>>> Martin
>>>>
>>>> On Jan 16, 2014, at 10:13 AM, Anthony Coletti <ajcolett at geo.umass.edu> wrote:
>>>>
>>>>> I should be more specific with my problem. Apologies.
>>>>>
>>>>> See below for the original problem:
>>>>>
>>>>> I am using the sea-ice package and external forcing packages (along with other default packages) to do a control LGM simulation of the Arctic. I am currently trying to see how much sea ice I can grow in the Arctic during this period. I have prescribed most of the atmosphere, the initial SST and Salinity fields.
>>>>>
>>>>> My data file looks like this:
>>>>>
>>>>> # ====================
>>>>> # | Model parameters |
>>>>> # ====================
>>>>> #
>>>>> # Continuous equation parameters
>>>>> &PARM01
>>>>> tRef=15*20.,
>>>>> sRef=15*35.,
>>>>> viscAh =3.E5,
>>>>> viscAr =1.E-3,
>>>>> diffKhT=0.,
>>>>> diffK4T=0.,
>>>>> diffKrT=3.E-5,
>>>>> diffKhS=0.,
>>>>> diffK4S=0.,
>>>>> diffKrS=3.E-5,
>>>>> ivdc_kappa=10.,
>>>>> implicitDiffusion=.TRUE.,
>>>>> gravity=9.81,
>>>>> rhoConst=1035.,
>>>>> rhoConstFresh=1000.,
>>>>> eosType='JMD95Z',
>>>>> staggerTimeStep=.TRUE.,
>>>>> vectorInvariantMomentum=.TRUE.,
>>>>> implicitFreeSurface=.TRUE.,
>>>>> exactConserv=.TRUE.,
>>>>> select_rStar=2,
>>>>> nonlinFreeSurf=4,
>>>>> hFacInf=0.2,
>>>>> hFacSup=2.0,
>>>>> useRealFreshWaterFlux=.TRUE.,
>>>>> hFacMin=.1,
>>>>> hFacMinDr=20.,
>>>>> #- to check that it has no impact:
>>>>> # doResetHFactors=.TRUE.,
>>>>> # tempAdvScheme=7,
>>>>> # saltAdvScheme=7,
>>>>> readBinaryPrec=32,
>>>>> writeBinaryPrec=32,
>>>>> useSingleCpuIO=.TRUE.,
>>>>> &
>>>>>
>>>>> # Elliptic solver parameters
>>>>> &PARM02
>>>>> cg2dMaxIters=200,
>>>>> #cg2dTargetResidual=1.E-9,
>>>>> cg2dTargetResWunit=1.E-14,
>>>>> &
>>>>>
>>>>> # Time stepping parameters
>>>>> &PARM03
>>>>> nIter0=0.,
>>>>> endTime=31557600000.,
>>>>> deltaT= 1500.,
>>>>> abEps = 0.1,
>>>>> forcing_In_AB=.FALSE.,
>>>>> pickupStrictlyMatch=.FALSE.,
>>>>> cAdjFreq=0.,
>>>>> pChkptFreq = 31557600.,
>>>>> chkptFreq = 15778800.,
>>>>> monitorFreq = 21600.,
>>>>> &
>>>>>
>>>>> # Gridding parameters
>>>>> &PARM04
>>>>> usingCurvilinearGrid=.TRUE.,
>>>>> horizGridFile='../../../input.files/grid_inputs/grid_cs32',
>>>>> delR= 50., 70., 100., 140., 190.,
>>>>> 240., 290., 340., 390., 440.,
>>>>> 490., 540., 590., 640., 690.,
>>>>> &
>>>>>
>>>>> # Input datasets
>>>>> &PARM05
>>>>> bathyFile ='../../../input.files/LGM_cntr/lgm+ish_cntr_bathy.bin',
>>>>> hydrogThetaFile='../../../input.files/LGM_cntr/sfc_T_coldstart_32x32x15.bin',
>>>>> hydrogSaltFile ='../../../input.files/LGM_cntr/oce_salt_coldstart_32x32x15.bin',
>>>>> &
>>>>>
>>>>>
>>>>> When I lower my time step BELOW 5000. s, it does not crash and is stable the entire run. I also have bathymetry set to LGM levels (adding 120m to the bathymetry files) to make the ocean more shallow.
>>>>>
>>>>>
>>>>> Anthony
>>>>>
>>>>>
>>>>> Anthony J. Coletti
>>>>> Climate System Research Center
>>>>> Department of Geosciences
>>>>> Morrill Building
>>>>> 611 N. Pleasant Street
>>>>> 233 Morrill Science Center
>>>>> University of Massachusetts-Amherst
>>>>> paleoclimate.org
>>>>> Email: ajcolett at geo.umass.edu
>>>>> http://blogs.umass.edu/ajcolett/
>>>>>
>>>>> “For me, I am driven by two main philosophies: know more today about the world than I knew yesterday and lessen the suffering of others. You'd be surprised how far that gets you.” ― Neil deGrasse Tyson
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On Jan 16, 2014, at 9:45 AM, Anthony Coletti <ajcolett at geo.umass.edu> wrote:
>>>>>
>>>>>> Hi all,
>>>>>>
>>>>>> I am trying to run a course 32x32x15 (6 faces) run to get some preliminary data. Unfortunately, my time step will not exceed 5000.s. ~ Two minutes into the run, it crashes and I get the corresponding error::
>>>>>>
>>>>>> WARNING: r*FacC > hFacSup at 1 pts : bi,bj,Thid,Iter= 4 1 1 106
>>>>>> e.g. at i,j= 2 1 ; rStarFac,H,eta = 2.028435 2.200000E+01 2.262556E+01
>>>>>> WARNING: r*FacC > hFacSup at 1 pts : bi,bj,Thid,Iter= 4 1 1 107
>>>>>> e.g. at i,j= 2 1 ; rStarFac,H,eta = 2.333596 2.200000E+01 2.933912E+01
>>>>>> WARNING: r*FacW < hFacInf at 2 pts : bi,bj,Thid,Iter= 2 1 1 108
>>>>>> e.g. at i,j= 5 2 ; rStarFac,H,eta = 0.094709 1.368000E+03 -6.740896E+02 -1.803211E+03
>>>>>> STOP in CALC_R_STAR : too SMALL rStarFacW !
>>>>>> WARNING: r*FacS < hFacInf at 1 pts : bi,bj,Thid,Iter= 2 1 1 108
>>>>>> e.g. at i,j= 3 2 ; rStarFac,H,eta = -0.186928 2.400000E+02 5.703576E+01 -6.277609E+02
>>>>>> STOP in CALC_R_STAR : too SMALL rStarFacS !
>>>>>> WARNING: r*FacC < hFacInf at 2 pts : bi,bj,Thid,Iter= 2 1 1 108
>>>>>> e.g. at i,j= 5 2 ; rStarFac,H,eta = -0.152948 1.564000E+03 -1.803211E+03
>>>>>> STOP in CALC_R_STAR : too SMALL rStarFacC !
>>>>>>
>>>>>> If you need any other info, please let me know. Any help would be appreciated!
>>>>>>
>>>>>>
>>>>>> Thanks!
>>>>>>
>>>>>> Anthony
>>>>>>
>>>>>> Anthony J. Coletti
>>>>>> Climate System Research Center
>>>>>> Department of Geosciences
>>>>>> Morrill Building
>>>>>> 611 N. Pleasant Street
>>>>>> 233 Morrill Science Center
>>>>>> University of Massachusetts-Amherst
>>>>>> paleoclimate.org
>>>>>> Email: ajcolett at geo.umass.edu
>>>>>> http://blogs.umass.edu/ajcolett/
>>>>>>
>>>>>> “For me, I am driven by two main philosophies: know more today about the world than I knew yesterday and lessen the suffering of others. You'd be surprised how far that gets you.” ― Neil deGrasse Tyson
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>> _______________________________________________
>>>>> 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
>>>
>>>
>>> --
>>> Christoph Voelker
>>> Alfred Wegener Institute for Polar and Marine Research
>>> Am Handelshafen 12
>>> 27570 Bremerhaven, Germany
>>> e: Christoph.Voelker at awi.de
>>> t: +49 471 4831 1848
>>>
>>>
>>> _______________________________________________
>>> 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
>
>
> _______________________________________________
> 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/20140116/26e59f67/attachment.htm>
More information about the MITgcm-support
mailing list