[MITgcm-devel] stevens boundary conditions

[Jean-Michel Campin]

Hi Martin,

You are right, I switch Western & Eastern OB in my example !
Regarding gad_advection, we agree (and your masking of gtracer is fine
with me).
In gad_calc_rhs, a separated array for vertical flux from gmredi 
should be fine (we could think later on how to make this piece
simpler, but for now, seem good to me).
And regarding this point:
> and finally mask all terms that are mulitplied with "localT", and/or *advFac, 
> ie. the last term in the last loop. Is that correct?
I would also mask 
  (fZon(i+1,j)-fZon(i,j))
and
  (fMer(i,j+1)-fMer(i,j))
in addition to all what you propose (and then you don't have to modify
the computation of fZon & fMer).

Let me know if this makes sense.

Cheers,
Jean-Michel

On Mon, Dec 12, 2011 at 12:41:27PM +0100, Martin Losch wrote:
> Hi Jean-Michel,
> 
> On Dec 12, 2011, at 3:07 AM, Jean-Michel Campin wrote:
> 
> > Hi Martin,
> > 
> > To summarize, you want to add, for T & S, (mostly surface) forcing 
> > and vertical diffusion to the stevens OBC values (= those @ OB_Jn,OB_Js
> > and @ OB_Ie,OB_Iw, so outside the OB interior region).  
> Yes, it turned out that, as usual, I was to hasty in reading and implementing, and actually the Stevens 1990 paper suggests, that one should update the tracers with "everything", but modified advection. What seems to be particular important is appropriate surface forcing and subsequent vertical mixing (vertical diffusion). In my experiments with constant (Stevens-)T/S restoring in the Bering Strait (where the is inflow almost all the time), my temperature remain near the restoring values, all year round, which is far too cold in summer, because the surface forcing is missing.
> > 
> > I think vertical advection should be "zero out" like horizontal advection
> > which is treated separately, as you say. 
> Yes.
> > Regarding horizontal diffusion, it does not make sense to me to keep part 
> > of those fluxes for this partial update of OBC value (e.g., @ i=OB_Ie, 
> > we don't have any valid tracer value @ i=OB_Ie-1, so diffusive flux in
> > X dir @ i=OB_Ie is meaningless) and it's better to cancel all of them
> > when computing flux divergence (but we should not zero out fluxes
> > at i=OB_Ie+1, right ? )
> I think you mixed up east and west, but otherwise I am fine with removing horizontal diffusion contribution at i=OB_Ie, but keeping the flux between OB_Ie and the first interior point (at the "eastern" boundary that would be OB_Ie-1)
> > 
> > As I wrote earlier (outside the list), it's not yet clear to me which one of 
> > the 2 solutions is best.
> In principle, I do not like option 1 too much (the one that I followed), because it requires to code many things again, introducing sources of error. I like the 2nd solution better, but it may involve changes in the main code that might change the results (see my attempts that change solutions because of order of summation).
> > 
> > In the option 2 way:
> > I think in gad_advection (because of OBC multi-dimensional implementation),
> > all horizontal advection flux divergence are already canceled.
> > But I would argue that vertical advection should be canceled too
> > (this is also what you propose).
> one could either apply maskInC to fverT in a separate loop, or multiply gtracer at the end as I initially suggested.
> > And in gad_calc_rhs, only vertical diffusive flux should be kept
> > in the divergence, with the divergence of all the others (including 
> > gmredi vertical component) being zero out.
> > 
> > Is it what you tried ? or something different ?
> > Because this looks quite different:
> >> I need to separate advective flux from the rest (in fZon/fMer),
> Here it becomes clear, why I asked you for help. Obviously I did not understand exactly, what's going on in gad_calc_rhs.F. Now my impression is, that I need to mask fVerT/af near line 526 (to remove vertical advection), then somehow mask the GM contribution (probably need an extra df field), and finally mask all terms that are mulitplied with "localT", and/or *advFac, ie. the last term in the last loop. Is that correct?
> > 
> > Going to look at EmPmR forcing issue.
> > 
> > However, I am little bit concern that with vertical mixing scheme 
> > (KPP, GGL90) and may be also some surface forcing computation,
> > because of horizontal averaging, it might be difficult to get a 
> > clean solution in those columns (in fact, in seaice_obcs test.exp. 
> > which is using KPP, some of the western OB values already 
> > influence the Eastern side).
> 
> So I should try an implementation of option 2 this experiment?
> 
> Martin
> 
> > 
> > Cheers,
> > Jean-Michel
> > 
> > On Thu, Dec 08, 2011 at 10:20:18AM +0100, Martin Losch wrote:
> >> Hi Jean-Michel,
> >> 
> >> following up on my previous post: my initial tests show, that for option 2, adding a maskInC in gad_advection "works" in the sense that I does not affect the results, but masking the contribution from advection in gad_calc_rhs is problematic: I need to separate advective flux from the rest (in fZon/fMer), and this separation alone (naively, I added two new afZona,afMeri, but no further masking or anything) changes the results of exp4 by 4-6 digits for cg2d (reducing the agreement to 12, 11 and 10), so just the order of summation. Is this expected?
> >> This would make option 2 not very attractive (for non-multidimAdvection), would it?
> >> 
> >> Martin
> >> 
> >> On Dec 7, 2011, at 4:59 PM, Martin Losch wrote:
> >> 
> >>> Jean-Michel,
> >>> 
> >>> it turns out that the Stevens BC story continues. So far, I implemented updating the tracer (T/S) value on the boundary only based on advection (modified model "mean" velocity and some kind of phase velocity that is similar to orlanski), and a restoring term. 
> >>> 
> >>> Now it turns out that it is necessary to include surface forcing and vertical diffusion for the boundary values as well. Stepping back this means that basically all terms in the tracer advection/diffusion equation on the boundary should be used, except for the advection. The advection is added in a modified version, and the restoring term added. How can that be implemented technically? I see two options:
> >>> 1. my quick fix is to add the desired terms to OB?t/s as required. In the case of surface forcing and vertical diffusion, I add these contributions at the end of do_oceanic_physics, and everything else remains the same.
> >>> 2. Far more elegant: Do you think it's possible to have the model do all computations, except for advection, compute advection and restoring tendencies in obcs (_calc_stevens), and apply these as a correction in obcs_apply_ts (is useStevens)? As far as I can see, this requires a maskInC in gad_advection (at the very end) and gad_calc_rhs (advFac*maskInC) and that's it. Is my thinking too simple?
> >>> 
> >>> A related problem is that the surface forcing for salinity is incomplete on the open boundary for useRealFreshWaterFlux=.true., because EmPmR is masked by maskInC in this case. Do you see a way to add that back in?
> >>> 
> >>> Martin
> >>> 
> >>> PS. In Stevens (1990), the phase velocity to correct to model advection on the boundary is computed from time steps n-1 and n-2: 
> >>> c=dy/dt * ( t(i,j+1,k,n-1) - t(i,j+1,k,n-2) )/( t(i,j+2,k,n-2)-t(i,j+1,k,n-2) )
> >>> I use
> >>> c=dy/dt * ( t(i,j+1,k,n-1) - t(i,j+1,k,n-2) )/( t(i,j+2,k,n-1)-t(i,j+1,k,n-1) )
> >>> to save me from storing t(i,j+2,k,n-2)-t(i,j+1,k,n-2). In orlanski_south, things are done "properly". I am not sure if this is an issue or not. For it does not cause any problems. One could think of merging the orlanski and the stevens code. What do you think?
> >>> 
> >>> 
> >>> 
> >>> 
> >>> 
> >>> _______________________________________________
> >>> 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