[MITgcm-support] periodic but "closed" channel flow
Patrick Rosendahl
Patrick.Rosendahl at zmaw.de
Tue Mar 11 08:25:08 EDT 2008
Hi,
I intended to "close" the channel, so that a return flow develops. Since
this is a periodic = infinite channel, this is a problem.
Yes, the OBCS package provides the Orlansky-type of BC together with the
balance-switch. But because it has to provide the BC at futureTime,
there is always a little discrepancy between the real flow at futureTime
and the BC. This results in a tilted surface, which is (if I look at my
infinte channel) a saw-tooth shaped surface.
Hence I copy/pasted the OBCSbalance routine and applied it to the entire
flowfield, while keeping the true periodic calculation of the model.
If I look at the flow in a moving frame, the model will not incorporate
the friction due to the return flow. Also, the velocity at the bottom
will not be zero anymore. Again, my apology, I forgot to mention my
no-slip condition.
Not dumb, just not appropriate for my specific problem setting (e.g.
no-slip condition at the bottom).
Patrick
David Ferreira wrote:
> Hi all,
> Not sure I fully understand but isn't your problem
> equivalent to look at the flow field in a frame moving
> at the mean zonal speed?
> At equilibirum, because of conservation, the mean zonal flow
> (your Tr_T correction) should be constant longitudinally, it
> looks like a change of frame.
> What I am suggesting is that there is no correction to
> apply and no code to modify (or obsc package to use), but
> just to diagnose Tr_T a posteriori and remove it from your
> zonal speed outputs.
> Dumb ?
> david
>
>
> Quoting Patrick Rosendahl <Patrick.Rosendahl at zmaw.de>:
>
>> Martin,
>>
>> I forgot to mention that I closed the N+S boundary. Yes, I am
>> expecting a static stationary solution.
>> The E+W boundaries are open, and windstresss will generate some
>> profile into the wind direction, causing a netflow in that direction.
>> This patch will act as if there was a constant pressure pushing as
>> much water back until the netflow over the boundary is zero (in this
>> case 5 pixels left of the E boundary), e.g. the channel is "closed".
>>
>> I dont expect this to work with time-varying solutions. In this case
>> you should average over 1 waveperiod, if that length can be determined
>> at all.
>>
>> Thank you for the answer, and yes, it helps. So I will apply the
>> correction in the momentum-correction-step routine only (also the
>> obcs_apply_uv is called from there in line 104).
>>
>> Patrick
>>
>> Martin Losch wrote:
>>
>>> Patrick,
>>> what are you trying to do? a zonally periodic channel that still has
>>> open boundaries in the zonal direction? Why?
>>> On average, a zonal channel, unless you have a localized unbalanced
>>> mass sources somewhere (e.g. freshwater flux at the surface), should
>>> have no zonal netflow anywhere. However, if you have a time varying
>>> solution, there will be fluctuations in the mass flux (with zero
>>> mean). Removing these mass fluxes as you seem to have done will
>>> introduce shocks to the system at every timestep. Also I don't know
>>> if obcs is the best place to do that. If you have to change the
>>> entire flow field, I would really try to do it in the correction
>>> step (momentum_correction_step.F), all you need to do is include
>>> GRID.h into that file and then you can do what you want to do.
>>>
>>> Does that help?
>>>
>>> Martin
>>>
>>> On 5 Mar 2008, at 14:53, Patrick Rosendahl wrote:
>>>
>>>> Hello all,
>>>>
>>>> this is code I use for simulating a infinite channel (periodic
>>>> domain in EW-direction) with return flow (netflow in any section is
>>>> 0).
>>>>
>>>> It is currently hooked up in obcs_apply_uv.F via
>>>>
>>>> CALL OBCS_MODIFY( bi, bj, uFld, vFld, myThid )
>>>>
>>>> and implemented like this in obcs_calc.F inside the IF-
>>>> useOBCSbalance-THEN switch
>>>>
>>>> C
>>>> C calculate the mass transfer across some section
>>>> C
>>>> Tr_T = 0. _d 0
>>>> Ar_T = 0. _d 0
>>>> DO K=1,Nr
>>>> DO J=1-Oly,sNy+Oly
>>>> I_obc = sNx-5
>>>> IF (I_obc.ne.0) THEN
>>>> Ar = drF(k)*hFacC(I_obc,j,k,bi,bj)*dyG(I_obc,j,bi,bj)
>>>> * _maskW(I_obc,J,K,bi,bj)
>>>> Ar_T = Ar_T + Ar
>>>> Tr_T = Tr_T + Ar * uVel(I_obc,J,K,bi,bj)
>>>> ENDIF
>>>> ENDDO
>>>> ENDDO
>>>> _GLOBAL_SUM_R8( Ar_T , myThid )
>>>> _GLOBAL_SUM_R8( Tr_T , myThid )
>>>> Tr_T = (0. - Tr_T)/Ar_T
>>>> C
>>>> C correct the entire u-field
>>>> C
>>>> DO K=1,Nr
>>>> DO J=1-Oly,sNy+Oly
>>>> DO I=1-Olx,sNx+Olx
>>>> uVel(I,J,K,bi,bj) = uVel(I,J,K,bi,bj) + Tr_T * _maskW
>>>> (I_obc,J,K,bi,bj)
>>>> ENDDO
>>>> ENDDO
>>>> ENDDO
>>>>
>>>>
>>>>
>>>>
>>>> I first tried the Orlanski-type BC, but it does not work very well.
>>>> Probably there is a proper way to hook the code.
>>>>
>>>> Patrick
>>>> _______________________________________________
>>>> MITgcm-support mailing list
>>>> MITgcm-support at mitgcm.org
>>>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>>>
>>>
>>>
>>> _______________________________________________
>>> MITgcm-support mailing list
>>> MITgcm-support at mitgcm.org
>>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>>
>>
>> _______________________________________________
>> MITgcm-support mailing list
>> MITgcm-support at mitgcm.org
>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>>
>
>
> _______________________________________________
> MITgcm-support mailing list
> MITgcm-support at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-support
More information about the MITgcm-support
mailing list