[MITgcm-devel] seaice adjoint and EVP

Jinlun Zhang zhang at apl.washington.edu
Thu May 24 11:58:34 EDT 2007 

It is not just over open water, but also in the central arctic. However, 
the noise is suppressed with 1s timestep over both open water and pack 
ice. So I start to think perhaps nothing is wrong, just needing a small 
timestep.
Jinlun

Martin Losch wrote:

> I have a new feeble theory for the noise in the evp solver over open  
> ocean:
>
> heff = 0 over open ocean, therefore seaiceMassU/V = 0.
> momentum equation in seaice_evp is discretized (in time) as
> m*duice/dt = -m dphi/dx + tau_air + cd*(uice-uocean) + m*f*vice +  
> \nabla\sigma
> m*(uice(n+1)-uice(n))/dt = -m dphi(n)/dx + tau_air(n) - cd*(uice(n+1)- 
> uocean(n)) + m*f*vice(n) + \nabla\sigma(n),
> so coriolis is explicit, ice-ocean stress is implicit. if the mass m  
> is zero (and zetaMin=0, so that zeta=eta=press = 0 over open ocen)  
> this reduces to
> cd*uice(n+1) = tau_air(n) + cd*uocean(N)
> so that uice ist a purely diagnostic quantity and not time stepped.  
> cd is a function of uice-uocean at the nth time step, averaged to  
> center points and the averaged back to velocity points.
>
> Dimitris, could that be the problem, somehow I don't think so, but  
> you can try by putting a minimum seaiceMassU/V in seaice_dynsolver.F,  
> say seaiceMassU = max(seaiceMassU,SEAICE_rhoIce*0.05)
>
> Martin
>
> On 22 May 2007, at 18:59, Jinlun Zhang wrote:
>
>> Hi Martin,
>> Yeah, we are sort of stuck, but hey it is very interesting and  
>> revealing.
>> I would vote against masking ice velocities over open water  because, 
>> as mentioned earlier, the ice velocities would be wrong at  ice edge 
>> and the ice velocity discontinuity at ice edge will get  into ocean. 
>> (o:. We don't do the masking with LSR solver, perhaps  we can avoid 
>> doing that with EVP.
>> Jinlun
>>
>> Martin Losch wrote:
>>
>>> Hi Jinlun,
>>> the evp-solver is only in place for the C-grid. I don't have the  
>>> time  to code the solver for the b-grid now. The b-grid code (for  
>>> LSR) is  still working, but I have not kept it up to date, so  there 
>>> may be a  few thing different other than the different grids.
>>>
>>> In general I though that the c-grid is perfect for evp as all the   
>>> discretizations fall in place naturally. Only for this \delta  term  
>>> one needs to average from center to corner points and vice  versa  
>>> (have a look at seaice_calc_strainrates and seaice_evp).  However,  
>>> there may be issues with the coriolis terms (commonly a  problem 
>>> with  the c-grid).
>>>
>>> Actually, Elizabeth told us that she masks ice velocities over  
>>> open  water in CICE.
>>>
>>> Now we are a little stuck, aren't we?
>>>
>>> Martin
>>>
>>> PS. I need to be able to reproduce these results myself (I  haven't  
>>> been able to, yet), maybe I can debug the stuff this way.  Via 
>>> email  etc. it's quite demanding (o:
>>>
>>> On 21 May 2007, at 19:15, Jinlun Zhang wrote:
>>>
>>>> I wouldn't think C-grid is problematic with EVP as we have seen.   
>>>> But just to make sure, is it possible to use the original B-grid   
>>>> EVP to see if the same things occur? There was a B-grid ice  model  
>>>> setup in place that may be used for doing B-grid.
>>>> Better not zap out things over open ocean. Otherwise,  
>>>> discontinuity  may occur and ocean may be screwed up.
>>>> Jinlun
>>>>
>>>>>> Hi Martin,
>>>>>>
>>>>>>> 1. Are these figures all with with zMin = 0?
>>>>>>
>>>>>>
>>>>>>
>>>>> In this case it may be worth turning of individual terms in the   
>>>>> rhs  of the momentum equations
>>>>> 1. dphiSurf/dx and dphiSurf/dy (in seaice_dynsolver)
>>>>> 2. surface wind stress (taux/y=0 in seaice_get_dynforcing)
>>>>> 3. ice-ocean stress (DWATN in seaice_evp)
>>>>> 4. Coriolis
>>>>> 5. stressDivergence
>>>>> 4 and 5 should be zero over open ocean anyway so I do not see  
>>>>> how   these terms can lead to the stripes.
>>>>> We should get to the bottom of what is causing these stripes.   
>>>>> that  way we can probably understand the noise in the ice  
>>>>> fields,  too.
>>>>>
>>>>>>
>>>>>> Yes, all the figures and results under
>>>>>> http://ecco2.jpl.nasa.gov/data1/arctic/output/tests/
>>>>>> (except for the oldtest subdirectory) are with zMin=0.
>>>>>>
>>>>>>> 2. Do you have an EVP run that does not blow up at all   
>>>>>>> (regardless  of noise)?
>>>>>>
>>>>>>
>>>>>>
>>>>>> I have not run any of the zMin=0/SEAICEuseFlooding=.true.  tests  
>>>>>> out  for very
>>>>>> long, but I am almost certain that none of these new   
>>>>>> integrations  will crash,
>>>>>> including the SEAICE_deltaTevp=60.
>>>>>> The crashes had to do with snow accumulation and could happen  
>>>>>> to   both LSR or to
>>>>>> EVP solutions.
>>>>>
>>>>>
>>>>>
>>>>> That's good news. It mean that we can (in principle) maskRHS  
>>>>> flag  and  not worry about the stripes.
>>>>>
>>>>>>
>>>>>>> 3. What's the convergence criterion for LSR, and how many    
>>>>>>> interations do you allow/do? In other words how close is the   
>>>>>>> LSR  solution to VP?
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>       LSR_ERROR          = 2e-4,
>>>>>>       SOLV_MAX_ITERS=1500
>>>>>
>>>>>
>>>>>
>>>>> That's not very much, is it? For an accurate VP solution I  would  
>>>>> put  LSR_ERROR = 1e-7 to 1e-13, right?
>>>>>
>>>>>>
>>>>>>> c. the same is true for the wind-ice/ocean-ice stress terms   
>>>>>>> which  in involve
>>>>>>>  averaging perpendicular to the stripes (unless the turning   
>>>>>>> angle  is not
>>>>>>> equal to zero, in which case there is also averaging in the   
>>>>>>> other  directions,
>>>>>>> but you don't do that, do you?).
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> No I use SEAICE_airTurnAngle=SEAICE_waterTurnAngle=0.
>>>>>
>>>>>
>>>>>
>>>>> Good.
>>>>>
>>>>>>
>>>>>>> About question 3 (is it really a VP solution): Could you   
>>>>>>> diagnose  SIsigI and SIsigII (snapshots!!!! I guess one is   
>>>>>>> enough) for all  (or some) solutions and
>>>>>>>  plot them (plot(SIsigII(:),SIsigI(:),'x')? These should be  
>>>>>>> the   principle components of sigma normalized by the  strength/ 
>>>>>>> pressure P.
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> With SEAICE_dumpFreq, SIGMA1, SIGMA2, and SIGMA12 are  diagnosed  
>>>>>> by  default for
>>>>>> the EVP solutions but not for LSR.  Are these the same as  
>>>>>> SIsigI   and SIsigII?
>>>>>> Figure for SIGMA1, SIGMA2 for EVP solution is here:
>>>>>> http://ecco2.jpl.nasa.gov/data1/arctic/output/tests/figs/  
>>>>>> SIGMA2232.ps
>>>>>> Does it look as expected?
>>>>>
>>>>>
>>>>>
>>>>> sigma1/2/12 are not the principle stress components. I have  
>>>>> added   diagnostics that are called SIsigI and SIsigII, which is  
>>>>> what you   want. In principle you could computed them yourself  
>>>>> (from snapshots):
>>>>> SIsigI = 0.5*(sigma1 + sqrt(sigma2^2 + 4*sigma12^2)/Press
>>>>> SIsigII = 0.5*(sigma1 - sqrt(sigma2^2 + 4*sigma12^2)/Press
>>>>>
>>>>> Press = max(1.e-13,Pstar * HEFF *exp( -20*(1-AREA)));
>>>>>
>>>>> see seaice_do_diags.F (and seaice_dynsolver.F)
>>>>>
>>>>>>
>>>>>>> I am also a little concerned that the LSR and EVP solutions   
>>>>>>> look  so different
>>>>>>> in the ice-covered area, can that be attributed to that   
>>>>>>> different  boundary
>>>>>>> conditons? Can you try a run with no slip for the evp solver?
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> Is LSR no slip by default?  How do you specify no slip for evp   
>>>>>> solver?
>>>>>
>>>>>
>>>>>
>>>>> LSR is half slip and that's hardwired. I didn't want to bother   
>>>>> this  the boundary conditions if EVP works, because it's so  much  
>>>>> simpler to  do that in EVP. But now I may have to  reconsider 
>>>>> this  decision.
>>>>> EVP is free slip by default. SEAICE_no_slip = .true. makes it  no  
>>>>> slip.
>>>>>
>>>>> Martin
>>>>>
>>>>> _______________________________________________
>>>>

More information about the MITgcm-devel mailing list