[MITgcm-devel] more seaice blues for the adjoint

[Jinlun Zhang]

Martin, I am in Japan, so I catch this email of yours quick.
Good to know B-grid also has this problem. But I thought if the motion 
gets smaller and smaller, then delta approaches zero, and zeta would hit 
maximum, and eventually, everything is constant, thus dp/dx is zero. 
This would not be related to zmin. Am I wrong here?
Jinlun

> Hi Jinlun,
>
> I check with the b-code and at first I got zero motion!!! But why? the 
> reason is not the grid, but this line:
> zmin = 4e8
> if I set this to
> zmin = 0
> as we have agree a few days/weeks ago, I do get motion.
>
> So it's not the C-grid, but this lower limit for zeta, that makes zeta 
> = constant if delta = 0, and thus P = zeta*2*delt2 = const, so no 
> dP/dx -> no motion. if zmin = 0, then zeta is not constant, neither P 
> and we have motion.
> so if zmin = 0, then I need to do something else to ensure P = 
> constant if delta = 0., eg. P = P * delta/delt2
>
> Martin
>
> On 8 Jun 2007, at 08:38, Jinlun Zhang wrote:
>
>> Hi Martin,
>> Before your fix, is it occurring with both LSR and EVP C-grid codes? 
>> Also, would it be possible for you to try the B-grid that I put in 
>> originally (if it is still one of the options)? It is OK for not 
>> trying if you don't have time. I remember I tested this thing with 
>> B-grid years ago without finding similar problem, but I am not sure 
>> now after so many years without dealing with it again. I just hope 
>> this is not one of the differences between B and C. Your explanation 
>> certainly indicates that this thing does not depend on grid. But a 
>> test on B-grid would clear any residual doubts.
>> Thanks, Jinlun
>>
>>
>> Martin Losch wrote:
>>> Chris,
>>> when there is no strain (in the absence of any forcing), the ice 
>>> should not move, right? Unless you expect it to spread like a pile 
>>> of sand with time, but I don't know what the time scale for that 
>>> would be, very slow.
>>> If delta = 0, we reset it to delta_min, because we have to devide by 
>>> it zeta = press/(2*delta), this effictively makes zeta <= zeta_max = 
>>> press/(2*delta_min), after capping zeta, we recompute press_replace 
>>> = zeta*2*delta, if delta has been replaced by delta_min, this give 
>>> non zero pressure and thus a forcing due to internal stress. That's 
>>> the numerical issue. I am not shure that I solved it very well (the 
>>> adjoint breaks), but this way the ice does not move in the absence 
>>> of forcing.
>>>
>>> Martin
>>> On 7 Jun 2007, at 18:05, chris hill wrote:
>>>
>>>> Martin,
>>>>
>>>>  Quick ice question.
>>>>
>>>> Is it just a numerical issue i.e. is it clear what the underlying 
>>>> equations are, since ice is neither a classical fluid or a 
>>>> classical solid. What time scale does it spread on?
>>>>
>>>> Chris
>>>> Martin Losch wrote:
>>>>> No, but it's a numerical issure. I did not pay any attention to 
>>>>> this before and maybe this was no problem in the original code 
>>>>> that you supplied but if you have zero strain, Delta = 0, and this 
>>>>> is replaced by some Delta_min (SEAICE_EPS in the code) in order to 
>>>>> avoid division by zero (and infinite viscosities). P however is 
>>>>> still non-zero and divergence(stress) end up being (dP/dx, dP/dy) 
>>>>> .NE. 0, so that you have a forcing of down the slope of P. This 
>>>>> does not make sense (why would ice spontaneously spread?), so that 
>>>>> the pressure is replaced by P = P*Delta/max(Delta,SEAICE_EPS), so 
>>>>> that in the (rare) case of no strain P = zeta = eta = 0 and thus 
>>>>> no forcing by stress. Makes sense to me (and is really only 
>>>>> relevant in idealized test case, just like the spurious motion of 
>>>>> sigma coordinates vs. z-coordinates in ocean models).
>>>>> Hope I did not misunderstand anything here.
>>>>> Martin
>>>>> On 7 Jun 2007, at 13:40, Jinlun Zhang wrote:
>>>>>> Hi Martin,
>>>>>> Would this be due to some kind of finite differencing problem?
>>>>>> Jinlun
>>>>>>
>>>>>> Martin Losch wrote:
>>>>>>> Hi Patrick,
>>>>>>>
>>>>>>> I have found (or rather, was pointed to) a problem with the 
>>>>>>> seaice solvers: The start to move spontaneously (in the absence 
>>>>>>> of forcing), if the sea ice distribution is NOT uniform.
>>>>>>> I have implemented a fix but this will cause problems with the 
>>>>>>> adjoint: I need terms like
>>>>>>> SQRT(deltaC), which used to be SQRT(MAX(deltaC,SEAICE_EPS_SQ)), 
>>>>>>> so that the derivate code will be involve 1/sqrt(deltaC). Should 
>>>>>>> I put this into #ifdefs?
>>>>>>>
>>>>>>> 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
>
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel