September 2015 - firedrake

Re: [firedrake] Cell center on extruded mesh
by Andrew McRae 28 Sep '15

28 Sep '15

Ah, my first instinct was that the FE way would be exactly equivalent, but of course that's not true. Anyway, lowest-order FE is similar to FV, so here's the FE approach I was thinking of: The space (mesh, "NCF", 1) has one degree of freedom per facet ( http://femtable.org/, 3rd column, bottom panel). We want to find a function in this space which is somehow the gradient of your p in DG0. It'd be nice if we could just do u = TrialFunction(W) w = TestFunction(W) a = dot(w, u)*dx L = dot(w, grad(p))*dx out = Function(W) solve(a == L, out) However, this clearly won't give the result you want, since p is in DG0, so grad(p) is 0 everywhere, so the right hand side is 0. Instead, integrate the right hand side by parts to obtain L = p_bc*dot(w, n)*ds - div(w)*p*dx, where n = FacetNormal(mesh), and p_bc is some function with the values of p that you want on the boundary (this can be just p, or you might choose to put p_bc in a continuous space, for example). Complete code for the 1D version of this is --- from firedrake import * mesh = UnitIntervalMesh(10) V = FunctionSpace(mesh, "DG", 0) f = Function(V).interpolate(Expression("2*x[0] - 1.0")) f_bc = Function(FunctionSpace(mesh, "CG", 1)).interpolate(Expression("2*x[0] - 1.0")) W = FunctionSpace(mesh, "CG", 1) u = TrialFunction(W) v = TestFunction(W) a = u*v*dx n = FacetNormal(mesh) L = f_bc*v*n[0]*ds - v.dx(0)*f*dx g = Function(W) solve(a == L, g) print g.dat.data --- Indeed, this verifies that the FE 'weak gradient' of 2x-1 is 2 everywhere. For the 3D version, W needs to be FunctionSpace(mesh, "NCF", 1), and L needs to be as given above. Andrew On 24 September 2015 at 07:15, Buesing, Henrik < HBuesing(a)eonerc.rwth-aachen.de> wrote: > > > > > -- > > Dipl.-Math. Henrik Büsing > > Institute for Applied Geophysics and Geothermal Energy > > E.ON Energy Research Center > > RWTH Aachen University > > ------------------------------------------------------ > > Mathieustr. 10 | Tel +49 (0)241 80 49907 > > 52074 Aachen, Germany | Fax +49 (0)241 80 49889 > > ------------------------------------------------------ > > http://www.eonerc.rwth-aachen.de/GGE > > hbuesing(a)eonerc.rwth-aachen.de > > ------------------------------------------------------ > > > > *Von:* firedrake-bounces(a)imperial.ac.uk [mailto: > firedrake-bounces(a)imperial.ac.uk] *Im Auftrag von *Andrew McRae > *Gesendet:* 22 September 2015 20:42 > *An:* firedrake(a)imperial.ac.uk > *Betreff:* Re: [firedrake] Cell center on extruded mesh > > > > For the distances between cell centres, yes. > > It feels like you're performing some operations that aren't particularly > finite element-y? Would you be using the approximate gradient in a > variational form immediately, or are there further calculations you want to > carry out? > > In the former case, I think there's a nice FE way to get what you want. > > > > *[Buesing, Henrik] @Andrew: So, what would be the FE way of doing this?* > > > * Henrik* > > > > > > In the latter case, we can probably hack something up if we're smart > enough... > > > > On 22 September 2015 at 18:42, Buesing, Henrik < > HBuesing(a)eonerc.rwth-aachen.de> wrote: > > > ------------------------------ > > *Von:* firedrake-bounces(a)imperial.ac.uk [firedrake-bounces(a)imperial.ac.uk]" > im Auftrag von "Andrew McRae [A.T.T.McRae(a)bath.ac.uk] > *Gesendet:* Dienstag, 22. September 2015 18:59 > *An:* firedrake(a)imperial.ac.uk > *Betreff:* Re: [firedrake] Cell center on extruded mesh > > >>Can you say a bit more about how you plan to use this? > >>I wonder to what end you would like these values? > I would like to approximate \grad(pw) on a DG0 Function space by > (p('+') - p('-'))/|x('+')-x('-')|, where x are the cell centers. On a > hexahedral mesh this should just give \Delta x, \Delta y, \Delta z, right? > > Henrik > > > > On 22 September 2015 at 17:52, Buesing, Henrik < > HBuesing(a)eonerc.rwth-aachen.de> wrote: > > Dear Firedrakers, > > on an extruded mesh: > > > > > > > > *meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = > ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", > quadrilateral, 0) vert_elt = FiniteElement("DG", interval, 0) elt = > OuterProductElement(horiz_elt, vert_elt) *How can I access the element > center and the center of the element faces? > Thank you! > > Henrik > > > > > -- > Dipl.-Math. Henrik Büsing > Applied Geophysics and Geothermal Energy > E.ON Energy Research Center > RWTH Aachen University > ----------------------------------------------- > Mathieustr. 10 | Tel +49 (0)241 80 49907 > 52074 Aachen, Germany | Fax +49 (0)241 80 49889 > ----------------------------------------------- > http://www.eonerc.rwth-aachen.de/GGE > hbuesing(a)eonerc.rwth-aachen.de > ----------------------------------------------- > > > > >

2 1

ubuntu 14.04 build problem
by osman 25 Sep '15

25 Sep '15

Hi, I just installed firedrake globally with --global --sudo options. I get this: Successfully installed firedrake Cleaning up... Creating firedrake-update script. Successfully installed Firedrake. To upgrade firedrake, run firedrake/bin/firedrake-update Now to test it: osman@osman-desktop:~/FIREDRAKE$ python helmholtz.py Traceback (most recent call last): File "helmholtz.py", line 47, in <module> from firedrake import * File "/usr/local/lib/python2.7/dist-packages/firedrake/__init__.py", line 3, in <module> from pyop2_utils import enable_mpi_prefork ImportError: cannot import name enable_mpi_prefork Have I missed anything? python is 2.7, cython0.23 Thanks in advance -osman

1 0

Re: [firedrake] Mimetic elements on deformed extruded meshes
by Colin Cotter 24 Sep '15

24 Sep '15

Hi Tuomas, Yes, you need to set u.n=0 on the bottom. In practice this is done by specifying all 3 components to zero, because only normal components are stored on the boundary. all the best --cjc On 24 September 2015 at 19:49, Tuomas Karna <tuomas.karna(a)gmail.com> wrote: > Hi Colin, > > I've been using the kinematic condition, u*nx + v*ny + w*nz = 0, but now > that you mention it this is probably wrong for Piola elements. If (u,v) is > along the bottom face by definition the w=0 at the bottom. > > - Tuomas > > > > On 09/24/2015 11:37 AM, Colin Cotter wrote: > > Hi Tuomas, > > Just checking, how are you setting this boundary condition? > > There may be pressure gradient errors over steep topography, but I would > expect the velocity not to be too bad. > > all the best > --cjc > > On 24 September 2015 at 17:39, Tuomas Karna <tuomas.karna(a)gmail.com> > wrote: > >> 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(a)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(a)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 listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake >> >> >> > > > -- > http://www.imperial.ac.uk/people/colin.cotter > > www.cambridge.org/9781107663916 > > > > > _______________________________________________ > firedrake mailing listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

2 1

Re: [firedrake] Mimetic elements on deformed extruded meshes
by Colin Cotter 24 Sep '15

24 Sep '15

Hi Tuomas, Just checking, how are you setting this boundary condition? There may be pressure gradient errors over steep topography, but I would expect the velocity not to be too bad. all the best --cjc On 24 September 2015 at 17:39, Tuomas Karna <tuomas.karna(a)gmail.com> wrote: > 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(a)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(a)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 listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

2 1

Re: [firedrake] Mimetic elements on deformed extruded meshes
by Colin Cotter 24 Sep '15

24 Sep '15

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(a)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(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

2 1

Re: [firedrake] Firedrake issue
by David Ham 24 Sep '15

24 Sep '15

Hi Tomasz, Welcome back. Let's get through a few preliminaries. First, if you haven't used firedrake for 2-3 months then it will be quite out of date. You need to update to the current version. You can do this using your existing install and update all the packages, or you can use the new installer to get a new installation which will also have an attached update script. See http://firedrakeproject.org/download.html Next, it's much better to have this sort of conversation on the Firedrake list. This keeps everyone in the loop and enables other people to intervene when I say something stupid ;). I have taken the liberty of adding you to the list and I'm replying to the list now so that this conversation ends up there. Now, we need to be talking about code, not only maths, so please push what you're working on to a git repository on github or bitbucket so we can see what's going on. Finally we can start to think about the actual problem. I think that last time we talked about this, we discussed setting up a mixed mass matrix across the two domains and checking that projection works. Does that still work (and can you post some code to a repository which does this)? Once we get past that then we need to talk about how to assemble forms on the boundary facets. This will involve doing some mucking around with maps, because at the moment the facets belong to two different meshes. We need to think a little about how to do that, but we can do so when the first steps work. Cheers, David On Tue, 22 Sep 2015 at 12:54 Tomasz Salwa [RPG] <mmtjs(a)leeds.ac.uk> wrote: > Dear David, > > > I was busy in the recent months and could not work on the fluid-structure > interaction problem. Now I returned to it and reformulated it. > > > I am getting the following error: > > > RuntimeError: Nonlinear solve failed to converge after 0 nonlinear > iterations. > Reason: > Inner linear solve failed to converge after 10000 iterations with > reason: DIVERGED_MAX_IT > > > but the problem is with the proper formulation in Firedrake. > > I include the pdf with it. Would you be able to help me with this? > > > Regards, > > Tomasz > > > > ----- > > Tomasz Salwa > > PhD Research Student > > University of Leeds > > http://www1.maths.leeds.ac.uk/~mmtjs >

2 1

Mimetic elements on deformed extruded meshes
by Tuomas Karna 24 Sep '15

24 Sep '15

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

1 0

Re: [firedrake] Cell center on extruded mesh
by Andrew McRae 24 Sep '15

24 Sep '15

For the distances between cell centres, yes. It feels like you're performing some operations that aren't particularly finite element-y? Would you be using the approximate gradient in a variational form immediately, or are there further calculations you want to carry out? In the former case, I think there's a nice FE way to get what you want. In the latter case, we can probably hack something up if we're smart enough... On 22 September 2015 at 18:42, Buesing, Henrik < HBuesing(a)eonerc.rwth-aachen.de> wrote: > > ------------------------------ > *Von:* firedrake-bounces(a)imperial.ac.uk [firedrake-bounces(a)imperial.ac.uk]" > im Auftrag von "Andrew McRae [A.T.T.McRae(a)bath.ac.uk] > *Gesendet:* Dienstag, 22. September 2015 18:59 > *An:* firedrake(a)imperial.ac.uk > *Betreff:* Re: [firedrake] Cell center on extruded mesh > > >>Can you say a bit more about how you plan to use this? > >>I wonder to what end you would like these values? > I would like to approximate \grad(pw) on a DG0 Function space by > (p('+') - p('-'))/|x('+')-x('-')|, where x are the cell centers. On a > hexahedral mesh this should just give \Delta x, \Delta y, \Delta z, right? > > Henrik > > > On 22 September 2015 at 17:52, Buesing, Henrik < > HBuesing(a)eonerc.rwth-aachen.de> wrote: > >> Dear Firedrakers, >> >> on an extruded mesh: >> >> >> >> >> >> >> >> *meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = >> ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", >> quadrilateral, 0) vert_elt = FiniteElement("DG", interval, 0) elt = >> OuterProductElement(horiz_elt, vert_elt) *How can I access the element >> center and the center of the element faces? >> Thank you! >> >> Henrik >> >> >> >> >> -- >> Dipl.-Math. Henrik Büsing >> Applied Geophysics and Geothermal Energy >> E.ON Energy Research Center >> RWTH Aachen University >> ----------------------------------------------- >> Mathieustr. 10 | Tel +49 (0)241 80 49907 >> 52074 Aachen, Germany | Fax +49 (0)241 80 49889 >> ----------------------------------------------- >> http://www.eonerc.rwth-aachen.de/GGE >> hbuesing(a)eonerc.rwth-aachen.de >> ----------------------------------------------- >> > >

2 2

Re: [firedrake] Dirichlet BC on extruded mesh with discontinuous function space
by Andrew McRae 23 Sep '15

23 Sep '15

Yes, but it should still be possible to set the value of the function in those cells. Didn't Tuomas make some progress on this functionality? I don't remember what the current status is. As a workaround, Miklos's suggestion may suffice. On 23 September 2015 at 11:09, Miklos Homolya <m.homolya14(a)imperial.ac.uk> wrote: > Dear Henrik, > > You’ve got a DG0 space, i.e. all your nodes are in the cell centres. > Neither topologically nor geometrically have you any nodes on the boundary > faces. > > Maybe you want to impose your BC weakly? > > Regards, > Miklos > > On 23 Sep 2015, at 12:05, Buesing, Henrik <HBuesing(a)eonerc.rwth-aachen.de> > wrote: > > Dear Firedrakers, > > I’m constructing a DG space on an extruded mesh. > > > > > > *meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = > ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", > quadrilateral, 0) vert_elt = FiniteElement("DG", interval, 0) elt = > OuterProductElement(horiz_elt, vert_elt)* > *DG = FunctionSpace(mesh, elt)* > *W=DG*DG* > > Now I would like to impose Dirchlet BC on one of the boundaries. > I tried > *bc0 = DirichletBC(W.sub(0), Expression("1.0e7"), 2,method="geometric")* > > and get a > > *ValueError: Geometric boundary conditions are not yet supported on > extruded meshes* > > Using the standard topological way, seems to have no effect at all. > > Is there a way around this? I would be fine with fixing the value on the > most outer cell layer, if that is possible… > > Henrik > > > > > > -- > Dipl.-Math. Henrik Büsing > Institute for Applied Geophysics and Geothermal Energy > E.ON Energy Research Center > RWTH Aachen University > ------------------------------------------------------ > Mathieustr. 10 | Tel +49 (0)241 80 49907 > 52074 Aachen, Germany | Fax +49 (0)241 80 49889 > ------------------------------------------------------ > http://www.eonerc.rwth-aachen.de/GGE > hbuesing(a)eonerc.rwth-aachen.de > ------------------------------------------------------ > > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > > >

2 1

Re: [firedrake] Cell center on extruded mesh
by Colin Cotter 23 Sep '15

23 Sep '15

Yes, it seems like you might be better off using a package that is more focussed on FV? all the best --cjc On 22 September 2015 at 19:42, Andrew McRae <A.T.T.McRae(a)bath.ac.uk> wrote: > For the distances between cell centres, yes. > > It feels like you're performing some operations that aren't particularly > finite element-y? Would you be using the approximate gradient in a > variational form immediately, or are there further calculations you want to > carry out? > > In the former case, I think there's a nice FE way to get what you want. > In the latter case, we can probably hack something up if we're smart > enough... > > On 22 September 2015 at 18:42, Buesing, Henrik < > HBuesing(a)eonerc.rwth-aachen.de> wrote: > >> >> ------------------------------ >> *Von:* firedrake-bounces(a)imperial.ac.uk [firedrake-bounces(a)imperial.ac.uk]" >> im Auftrag von "Andrew McRae [A.T.T.McRae(a)bath.ac.uk] >> *Gesendet:* Dienstag, 22. September 2015 18:59 >> *An:* firedrake(a)imperial.ac.uk >> *Betreff:* Re: [firedrake] Cell center on extruded mesh >> >> >>Can you say a bit more about how you plan to use this? >> >>I wonder to what end you would like these values? >> I would like to approximate \grad(pw) on a DG0 Function space by >> (p('+') - p('-'))/|x('+')-x('-')|, where x are the cell centers. On a >> hexahedral mesh this should just give \Delta x, \Delta y, \Delta z, right? >> >> Henrik >> >> >> On 22 September 2015 at 17:52, Buesing, Henrik < >> HBuesing(a)eonerc.rwth-aachen.de> wrote: >> >>> Dear Firedrakers, >>> >>> on an extruded mesh: >>> >>> >>> >>> >>> >>> >>> >>> *meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = >>> ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", >>> quadrilateral, 0) vert_elt = FiniteElement("DG", interval, 0) elt = >>> OuterProductElement(horiz_elt, vert_elt) *How can I access the element >>> center and the center of the element faces? >>> Thank you! >>> >>> Henrik >>> >>> >>> >>> >>> -- >>> Dipl.-Math. Henrik Büsing >>> Applied Geophysics and Geothermal Energy >>> E.ON Energy Research Center >>> RWTH Aachen University >>> ----------------------------------------------- >>> Mathieustr. 10 | Tel +49 (0)241 80 49907 >>> 52074 Aachen, Germany | Fax +49 (0)241 80 49889 >>> ----------------------------------------------- >>> http://www.eonerc.rwth-aachen.de/GGE >>> hbuesing(a)eonerc.rwth-aachen.de >>> ----------------------------------------------- >>> >> >> > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

2 1