[MITgcm-devel] bloc and bloc in shelfice_thermodynamics

Martin Losch Martin.Losch at awi.de
Thu Oct 1 03:15:01 EDT 2015


Hi Dimitris,

the only reason for this is to save one floating point operation. I agree that is this can be confusing and I myself stumble last week over something like this in the seaice_evp code that I wrote myself.

I guess performance isn’t an issue here, so go ahead and change it for clarity

M.

> On 01 Oct 2015, at 01:08, Menemenlis, Dimitris (329D) <Dimitris.Menemenlis at jpl.nasa.gov> wrote:
> 
> Hi Martin, Yoshi and I wonder why you don’t multiply uLoc and vLoc by 0.5 in lines below:
> 
>          uLoc(I,J) = recip_hFacC(I,J,K,bi,bj) *
>     &         ( uVel(I,  J,K,bi,bj) * _hFacW(I,  J,K,bi,bj)
>     &         + uVel(I+1,J,K,bi,bj) * _hFacW(I+1,J,K,bi,bj) )
>          vLoc(I,J) = recip_hFacC(I,J,K,bi,bj) *
>     &         ( vVel(I,J,  K,bi,bj) * _hFacS(I,J,  K,bi,bj)
>     &         + vVel(I,J+1,K,bi,bj) * _hFacS(I,J+1,K,bi,bj) )
> 
>            uLoc(I,J) = ( uLoc(I,J) * drF(K)*_hFacC(I,J,K,bi,bj)
>     &           + drKp1 * recip_hFacC(I,J,Kp1,bi,bj) *
>     &           ( uVel(I,  J,Kp1,bi,bj) * _hFacW(I,  J,Kp1,bi,bj)
>     &           + uVel(I+1,J,Kp1,bi,bj) * _hFacW(I+1,J,Kp1,bi,bj) )
>     &           ) * recip_drLoc
>            vLoc(I,J) = ( vLoc(I,J) * drF(K)*_hFacC(I,J,K,bi,bj)
>     &           + drKp1 * recip_hFacC(I,J,Kp1,bi,bj) *
>     &           ( vVel(I,J,  Kp1,bi,bj) * _hFacS(I,J,  Kp1,bi,bj)
>     &           + vVel(I,J+1,Kp1,bi,bj) * _hFacS(I,J+1,Kp1,bi,bj) )
>     &           ) * recip_drLoc
> 
> You do it later on in line:
> 
>            ustarSq = shiCdrag * MAX( 1.D-6,
>     &           0.25 _d 0 *(uLoc(I,J)*uLoc(I,J)+vLoc(I,J)*vLoc(I,J)) )
> 
> The problem with this formulation is that if one is trying to modify
> drag formulation, it’s very easy to forget that uLoc = 2 * u
> (Yoshi spent days tracking this missing factor :-(
> 
> Shall we change for clarity?
> 
> Cheers, Dimitris
> 
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel




More information about the MITgcm-devel mailing list