[MITgcm-devel] lab_sea.hb87 restart problem
Jean-Michel Campin
jmc at ocean.mit.edu
Mon Nov 12 19:10:59 EST 2007
Hi Martin,
I think I know (well, part of it) why the restart does not work
in this lab_sea test:
I printed seaice_sigma12(i,j,bi,bj) for i=3,
and clearly the value for j=1 and j=sNy+1 tend to diverge
(but both are stepped forward in seaice_evp.F, but not exchanged).
And the fact they diverge might be particular to the no-slip BC
(why ? not sure yet).
When we restart, we only write/read [1:sNy] and do an exch (1 -> j+1),
and this explains why we get differences.
2 little things:
1) I was thinking of adding iEVPstep,myTimer & myIter as argument
to S/R seaice_calc_strainrates.F . It really helps to debug.
I've done it in my own version, but could be worth to check-in.
2) The call to seaice_calc_strainrates in seaice_ocean_stress.F
is wrong. Did you already wrote about this in a previous e-mail ?
Cheers,
Jean-Michel
On Mon, Nov 12, 2007 at 11:46:10AM -0500, Jean-Michel Campin wrote:
> Hi Martin,
>
> I started to investigate into this direction (overlap & loop range)
> but did not go very far.
> Trying to do some debub print & plots, but it's a pain that
> myTime,myIter are not passed to all those S/R.
>
> If we add some EXCH, it will still be good to know why we
> need them. And also to call the right EXCH (e.g., e12 is
> at a corner point, so an EXCH_Z would be better).
>
> Cheers,
> Jean-Michel
>
> On Mon, Nov 12, 2007 at 10:20:05AM +0100, Martin Losch wrote:
> > Jean-Michel,
> > here's the solution:
> > I need to exchange the strain rates e11/e22/e12 (all three of them)
> > after they have been comuted in seaice_calc_strainrates.F
> > _EXCH_XY_R8( e11, myThid )
> > _EXCH_XY_R8( e22, myThid )
> > _EXCH_XY_R8( e12, myThid )
> > then the pickup is OK.
> > This changes the testreport results of lab_sea.hb87 (but only that
> > one and only the exchanges of e11 and e12). I do not see any reason,
> > why these extra exchanges should change lab_sea.hb87 but not any of
> > the other experiments. The computational stencil in the IF
> > (SEAICE_no_slip) THEN/ENDIF block is absolutely the same as in the
> > previous lines of code.
> >
> > Absolutely clueless,
> > Martin
> >
> > On 9 Nov 2007, at 20:44, Martin Losch wrote:
> >
> > >Hi Jean-Michel,
> > >
> > >I am stuck. With the following parameters in data.seaice
> > >>&SEAICE_PARM01
> > >> SEAICEwriteState = .TRUE.,
> > >> SEAICE_initialHEFF = 1.0,
> > >> SEAICE_deltaTtherm = 3600.,
> > >> SEAICE_deltaTdyn = 3600.,
> > >># useHB87stressCoupling = .true.,
> > >> SEAICE_no_slip = .TRUE.,
> > >> SEAICE_deltaTevp = 60.,
> > >I can get a perfect pickup, if I set hFacU to zero in the if
> > >(SEAICE_no_slip) block in seaice_calc_strainrates.F like this:
> > >> hFacU = _maskW(i,j,k,bi,bj) - _maskW(i,j-1,k,bi,bj)
> > >> hFacV = _maskS(i,j,k,bi,bj) - _maskS(i-1,j,k,bi,bj)
> > >> hFacU = 0. _d 0
> > >>CML hFacV = 0. _d 0
> > >>
> > >> e12(I,J,bi,bj) = e12(I,J,bi,bj)
> > >> & + recip_rAz(i,j,bi,bj) *
> > >> & ( hFacU * ( _dxC(i,j-1,bi,bj)*uFld(i,j ,bi,bj)
> > >> & + _dxC(i,j, bi,bj)*uFld(i,j-1,bi,bj) )
> > >> & + hFacV * ( _dyC(i-1,j,bi,bj)*vFld(i ,j,bi,bj)
> > >> & + _dyC(i, j,bi,bj)*vFld(i-1,j,bi,bj) ) )
> > >> & - hFacU
> > >> & * 0.25 _d 0 * (uFld(I,J,bi,bj)+uFld(I ,J-1,bi,bj))
> > >> & * ( _tanPhiAtU(I,J,bi,bj) + _tanPhiAtU(I,J-1,bi,bj) )
> > >> & *recip_rSphere
> > >This means, that at least the no-slip problem is related to uIce.
> > >But it cannot be directly related to uice, because the pickup (and
> > >this piece of code) is ok for many other cases (e.g. LSR).
> > >Somewhere the true problem is masked. I don't have a clue. I don't
> > >even know what to test next ...
> > >
> > >Martin
> > >
> > >
> > >On 9 Nov 2007, at 17:40, Jean-Michel Campin wrote:
> > >
> > >>Hi Martin,
> > >>
> > >>I don't know much about the cause of the problem,
> > >>but you narrow down quiet well where it sould be.
> > >>
> > >>I agree that it should not be too difficult to
> > >>make this tst1+1 script more tolerent regarding blanks
> > >>upper-case/lower case.
> > >>
> > >>I have a general question:
> > >>What about adding few "data.tst" in some verification/*/input
> > >>experiments
> > >>to facilitate those kind of test ? I am not sure it's a good
> > >>idea, but would interesting to know what other people think.
> > >>
> > >>Cheers,
> > >>Jean-Michel
> > >>
> > >>On Fri, Nov 09, 2007 at 05:21:37PM +0100, Martin Losch wrote:
> > >>>OK, there error was somewhere between keyboard and chair: your
> > >>>script
> > >>>is case sensitive and I was sloppy in editing data.tst.
> > >>>
> > >>>Now I have repeated the tests with your script and this is what I
> > >>>get:
> > >>>original lab_sea.hb87
> > >>>>csysm3::tr_run.hb87> ./tst1+1 2
> > >>>>10 12 14
> > >>>>0000000010 0000000012 0000000014
> > >>>>gcmExc=mitgcmuv
> > >>>>-- compare cg2d_init_res :
> > >>>>run 1iA:
> > >>>> 6.08551932893086E-02
> > >>>> 6.55586521348751E-02
> > >>>>run 1iB:
> > >>>> 7.02066252452855E-02
> > >>>> 7.29377130933022E-02
> > >>>>run 2it:
> > >>>> 6.08551932893086E-02
> > >>>> 6.55586521348751E-02
> > >>>> 7.01960854205512E-02
> > >>>> 7.29481919712837E-02
> > >>>lab_sea.hb87 with useHB87stressCoupling = .false.,
> > >>>>csysm3::tr_run.hb87> ./tst1+1 2
> > >>>>10 12 14
> > >>>>0000000010 0000000012 0000000014
> > >>>>gcmExc=mitgcmuv
> > >>>>-- compare cg2d_init_res :
> > >>>>run 1iA:
> > >>>> 6.08586878524103E-02
> > >>>> 6.55715983076715E-02
> > >>>>run 1iB:
> > >>>> 7.02079681169673E-02
> > >>>> 7.29487338830253E-02
> > >>>>run 2it:
> > >>>> 6.08586878524103E-02
> > >>>> 6.55715983076715E-02
> > >>>> 7.02078133431693E-02
> > >>>> 7.29489594851236E-02
> > >>>with useHB87stressCoupling = .false., and evp turned off (#
> > >>>SEAICE_deltaTevp = 60.,)
> > >>>>10 12 14
> > >>>>0000000010 0000000012 0000000014
> > >>>>gcmExc=mitgcmuv
> > >>>>-- compare cg2d_init_res :
> > >>>>run 1iA:
> > >>>> 6.08355827055190E-02
> > >>>> 1.16426574689837E-01
> > >>>>run 1iB:
> > >>>> 7.80912267704927E-02
> > >>>> 7.36159469063700E-02
> > >>>>run 2it:
> > >>>> 6.08355827055190E-02
> > >>>> 1.16426574689837E-01
> > >>>> 7.80912267704927E-02
> > >>>> 7.36159469063700E-02
> > >>>with useHB87stressCoupling = .true., and evp turned off (#
> > >>>SEAICE_deltaTevp = 60.,)
> > >>>and with useSeaice=.false.,
> > >>>>10 12 14
> > >>>>0000000010 0000000012 0000000014
> > >>>>gcmExc=mitgcmuv
> > >>>>-- compare cg2d_init_res :
> > >>>>run 1iA:
> > >>>> 6.13113158130375E-02
> > >>>> 6.57463167891239E-02
> > >>>>run 1iB:
> > >>>> 7.04816593980862E-02
> > >>>> 7.33615824658603E-02
> > >>>>run 2it:
> > >>>> 6.13113158130375E-02
> > >>>> 6.57463167891239E-02
> > >>>> 7.04816593980862E-02
> > >>>> 7.33615824658603E-02
> > >>>
> > >>>My conclusions from these numbers:
> > >>>1. Your script is great!!!
> > >>>2. The restart problem is related to EVP.
> > >>>However, I can reproduce your result in lab_sea/run, that the
> > >>>restart
> > >>>does work with EVP. Puzzling. Therefore I reran the test in lab_sea/
> > >>>run with SEAICE_no_slip = .true. (which is also used in
> > >>>tr_run.hb87),
> > >>>and then:
> > >>>>-- compare cg2d_init_res :
> > >>>>run 1iA:
> > >>>> 6.08981197561294E-02
> > >>>> 6.55972956406261E-02
> > >>>>run 1iB:
> > >>>> 7.02007081962388E-02
> > >>>> 7.29011820001111E-02
> > >>>>run 2it:
> > >>>> 6.08981197561294E-02
> > >>>> 6.55972956406261E-02
> > >>>> 7.02007451116285E-02
> > >>>> 7.29010160470697E-02
> > >>>For tr_run.lsr everything is OK even with SEAICE_no_slip = .true.
> > >>>
> > >>>So the problem is related to a combination of EVP and no-slip
> > >>>boundary conditions (and not the hb87 stress formulation). Now, the
> > >>>flag SEAICE_no_slip is not used anywhere in the EVP code, except in
> > >>>seaice_calc_strainrates.F, which is also used by the LSR code. I
> > >>>don't have no clue where this problem comes from. Why should break
> > >>>this flag the pickup for the evp code only, if it is not used in the
> > >>>code that is exclusive to EVP?
> > >>>
> > >>>Any idea?
> > >>>
> > >>>Martin
> > >>>
> > >>>
> > >>>
> > >>>
> > >>>
> > >>>_______________________________________________
> > >>>MITgcm-devel mailing list
> > >>>MITgcm-devel at mitgcm.org
> > >>>http://mitgcm.org/mailman/listinfo/mitgcm-devel
> > >>_______________________________________________
> > >>MITgcm-devel mailing list
> > >>MITgcm-devel at mitgcm.org
> > >>http://mitgcm.org/mailman/listinfo/mitgcm-devel
> > >
> > >_______________________________________________
> > >MITgcm-devel mailing list
> > >MITgcm-devel at mitgcm.org
> > >http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >
> > _______________________________________________
> > MITgcm-devel mailing list
> > MITgcm-devel at mitgcm.org
> > http://mitgcm.org/mailman/listinfo/mitgcm-devel
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel
More information about the MITgcm-devel
mailing list