Hi Colin,

Thanks for your reply. On a couple of test cases I see grid-scale noise in the velocity field, which may be significant over steep bathymetry gradients. But thus far I have been imposing the bottom impermeability boundary condition in the conventional way, which may be the origin of the noise. I'll keep testing this.

- Tuomas

On 09/24/2015 01:18 AM, Colin Cotter wrote:

Hi Tuomas,

  So the way I have been advocating mimetic/compatible spaces for atmosphere/ocean with terrain-following meshes is to use the Piola transform, this does indeed make the horizontal space no longer horizontal as you say. This allows an easy implementation of the u.n=0 boundary condition on sloping surfaces. For our tests for flows over mountains, this doesn't cause any problems as far as I can see. This is slightly different from the classical approach which is to keep storing horizontal velocities, but this leads to a more complicated coupled boundary condition on sloping surfaces.

My suggestion is to just try it and see if there are any problems.

all the best

--cjc

On 24 September 2015 at 07:46, Tuomas Karna <tuomas.karna@gmail.com> wrote:

Hi all,

I have a question related to mimetic element spaces on extruded meshes.

Hydrostatic ocean models solve a momentum equation for the horizontal
velocity components uv = (u, v, 0), while w is computed from the
continuity equation diagnostically.

I've been using the following spaces for uv and w:

# for horizontal velocity components
Uh_elt = FiniteElement('RT', triangle, self.order+1)
Uv_elt = FiniteElement('DG', interval, self.order)
U_elt = HDiv(OuterProductElement(Uh_elt, Uv_elt))
# for vertical velocity component
Wh_elt = FiniteElement('DG', triangle, self.order)
Wv_elt = FiniteElement('CG', interval, self.order+1)
W_elt = HDiv(OuterProductElement(Wh_elt, Wv_elt))
U = FunctionSpace(self.mesh, U_elt)  # uv
W = FunctionSpace(self.mesh, W_elt)  # w

This works fine for orthogonal extrusion, where all horizontal faces are
strictly horizontal. However, when you extrude the mesh over variable
topography this is no longer true. Then the U space is no longer
strictly horizontal which breaks the model. For example projecting a
vector (1, 0, 0) on U space on highly deformed mesh doesn't give you (1,
0, 0).

Earlier Andrew suggested enriching the U space with the W elements to
better represent a true 3d vector:

UW_elt = EnrichedElement(U_elt, W_elt)
U = FunctionSpace(self.mesh, UW_elt)

This has couple of issues though:
- the U space no longer has an equivalent 2d vector space, so 2d-3d
coupling breaks
- this is fairly complex space, and hence slow.

Do you see any way around this? Would it make more sense to just use DG
space for U and W instead?

Cheers,

Tuomas


_______________________________________________
firedrake mailing list
firedrake@imperial.ac.uk
https://mailman.ic.ac.uk/mailman/listinfo/firedrake

--

http://www.imperial.ac.uk/people/colin.cotter

www.cambridge.org/9781107663916

_______________________________________________
firedrake mailing list
firedrake@imperial.ac.uk
https://mailman.ic.ac.uk/mailman/listinfo/firedrake