[firedrake] Specify boundary conditions at a corner

Tue Aug 25 10:03:32 BST 2015

Hi Justin,

The imports in bcs.py are relative imports (which is slightly bad style
but...) In order to access correct module from outside you need to use:

from firedrake import utils

Regards,

David

On Tue, 25 Aug 2015 at 00:48 Justin Chang <jychang48 at gmail.com> wrote:

> David,
>
> I get the following error:
>
> Traceback (most recent call last):
>   File "Test_LSAug.py", line 41, in <module>
>     class PointDirichletBC(DirichletBC):
>   File "Test_LSAug.py", line 42, in PointDirichletBC
>     @utils.cached_property
> AttributeError: 'module' object has no attribute 'cached_property'
>
> When I comment out that line, i get a long error. When I try to import
> utils and all the other modules that bcs.py imports, i get errors
> saying those modules don't exist. Know what's going on?
>
> Thanks,
> Justin
>
> On Mon, Aug 24, 2015 at 3:14 AM, David Ham <David.Ham at imperial.ac.uk>
> wrote:
> > Hi Justin,
> >
> > Having subclassed DirichletBC, what you now have is a new type of
> Dirichlet
> > "boundary" object, so you just use it as a Dirichlet condition:
> >
> > bc3 = PointDirichletBC(W.sub(1), 0, 0)
> >
> > bc_all = [bc1, bc2, bc3]
> >
> > However Andrew is right, removing the nullspace is almost certainly a
> better
> > solution in this case.
> >
> > Cheers,
> >
> > David
> >
> > On Sun, 23 Aug 2015 at 05:23 Justin Chang <jychang48 at gmail.com> wrote:
> >>
> >> Andrew,
> >>
> >> That actually did the trick, thank you very much.
> >>
> >> But I still would like to know the answer to my original question. Or
> >> perhaps, instead of a single corner, I would like to set (x=0,y<0.1)
> >> and (x<0.1,y=0) to a specific value while everywhere else is set to a
> >> different value.
> >>
> >> Thanks,
> >> Justin
> >>
> >> On Sat, Aug 22, 2015 at 11:58 AM, Andrew McRae <
> a.mcrae12 at imperial.ac.uk>
> >> wrote:
> >> > This doesn't answer your main question, but setting an appropriate
> >> > nullspace
> >> > might be more appropriate than pinning a single value; see
> >> > http://www.firedrakeproject.org/solving-interface.html
> >> >
> >> > [If I'm wrong, I'll let the usual suspects correct my misinformation]
> >> >
> >> > On 22 August 2015 at 18:47, Justin Chang <jychang48 at gmail.com> wrote:
> >> >>
> >> >> David,
> >> >>
> >> >> How exactly do I use that class in my code? Say I have the following
> >> >> function spaces/discretization based on the least-squares finite
> >> >> element method:
> >> >>
> >> >> # Mesh
> >> >> mesh = UnitSquareMesh(seed,seed)
> >> >> V = VectorFunctionSpace(mesh,"CG",1)
> >> >> Q = FunctionSpace(mesh,"CG",1)
> >> >> W = V*Q
> >> >> v,p = TrialFunctions(W)
> >> >> w,q = TestFunctions(W)
> >> >>
> >> >> # Weak form
> >> >> g = Function(V)
> >> >>
> >> >>
> >> >>
> g.interpolate(Expression(("cos(pi*x[0])*sin(pi*x[1])+2*pi*cos(2*pi*x[0])*sin(2*pi*x[1])","-sin(pi*x[0])*cos(pi*x[1])+2*pi*sin(2*pi*x[0])*cos(2*pi*x[1])")))
> >> >> a = dot(v+grad(p),w+grad(q))*dx + div(v)*div(w)*dx
> >> >> L = dot(w+grad(q),g)*dx
> >> >>
> >> >> # Boundary conditions
> >> >> bc1 = DirichletBC(W.sub(0).sub(0),
> >> >> Expression("cos(pi*x[0])*sin(pi*x[1])"), (1,2))
> >> >> bc2 = DirichletBC(W.sub(0).sub(1),
> >> >> Expression("-sin(pi*x[0])*cos(pi*x[1])"), (3,4))
> >> >> bc_all = [bc1,bc2]
> >> >>
> >> >> # Solve
> >> >> cg_parameters = {
> >> >>     'ksp_type': 'cg',
> >> >>     'pc_type': 'bjacobi'
> >> >>     }
> >> >> solution = Function(W)
> >> >> A = assemble(a,bcs=bc_all)
> >> >> b = assemble(L,bcs=bc_all)
> >> >> solver =
> >> >> LinearSolver(A,solver_parameters=cg_parameters,options_prefix="cg_")
> >> >> solver.solve(solution,b)
> >> >>
> >> >> If I run the above code I get an error saying 'LinearSolver failed to
> >> >> converge after %d iterations with reason: %s', 196,
> >> >> 'DIVERGED_INDEFINITE_MAT'. Which I am guessing has to do with the
> lack
> >> >> of a boundary condition for the Q space, thus I want to ensure a
> >> >> unique solution by prescribing the bottom left constraint to a zero
> >> >> value.
> >> >>
> >> >> Thanks,
> >> >> Justin
> >> >>
> >> >>
> >> >> On Fri, Aug 21, 2015 at 4:52 AM, David Ham <David.Ham at imperial.ac.uk
> >
> >> >> wrote:
> >> >> > Hi Justin,
> >> >> >
> >> >> > The nice way of doing this would require better subdomain support
> >> >> > than
> >> >> > we
> >> >> > have right now. However there is a slightly hacky way of doing it
> >> >> > which
> >> >> > I
> >> >> > think will cover your case nicely.
> >> >> >
> >> >> > If you take a look at the DirichletBC class (in bcs.py), you'll
> >> >> > notice
> >> >> > that
> >> >> > the set of nodes at which the BC should be applied is calculated in
> >> >> > DirichletBC.nodes . So you could simply subclass DirichletBC and
> >> >> > replace
> >> >> > nodes with a function which returns the index of the zero node. For
> >> >> > example
> >> >> > (and I confess this is a sketch code which I haven't tried to run):
> >> >> >
> >> >> > class PointDirichletBC(DirichletBC):
> >> >> >     @utils.cached_property
> >> >> >     def nodes(self):
> >> >> >         # Find the array of coordinate values.
> >> >> >         x = self.function_space().mesh().coordinates.dat.data_ro
> >> >> >         # Find the location of the zero rows in that
> >> >> >         return np.where(~x.any(axis=1))[0]
> >> >> >
> >> >> > Does that work for you?
> >> >> >
> >> >> > Cheers,
> >> >> >
> >> >> > David
> >> >> >
> >> >> > On Fri, 21 Aug 2015 at 03:32 Justin Chang <jychang48 at gmail.com>
> >> >> > wrote:
> >> >> >>
> >> >> >> Hi all,
> >> >> >>
> >> >> >> If I create a mesh using UnitSquareMesh or UnitCubeMesh, is there
> a
> >> >> >> way to subject a single point (as opposed to an entire edge/face)
> to
> >> >> >> a
> >> >> >> DirichletBC? I want to subject the the location x=0,y=0 to some
> >> >> >> value.
> >> >> >>
> >> >> >> Thanks,
> >> >> >> Justin
> >> >> >>
> >> >> >> _______________________________________________
> >> >> >> firedrake mailing list
> >> >> >> firedrake at imperial.ac.uk
> >> >> >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake
> >> >> >
> >> >> >
> >> >> > _______________________________________________
> >> >> > firedrake mailing list
> >> >> > firedrake at imperial.ac.uk
> >> >> > https://mailman.ic.ac.uk/mailman/listinfo/firedrake
> >> >> >
> >> >>
> >> >> _______________________________________________
> >> >> firedrake mailing list
> >> >> firedrake at imperial.ac.uk
> >> >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake
> >> >
> >> >
> >> >
> >> > _______________________________________________
> >> > firedrake mailing list
> >> > firedrake at imperial.ac.uk
> >> > https://mailman.ic.ac.uk/mailman/listinfo/firedrake
> >> >
> >>
> >> _______________________________________________
> >> firedrake mailing list
> >> firedrake at imperial.ac.uk
> >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake
> >
> >
> > _______________________________________________
> > firedrake mailing list
> > firedrake at imperial.ac.uk
> > https://mailman.ic.ac.uk/mailman/listinfo/firedrake
> >
>
> _______________________________________________
> firedrake mailing list
> firedrake at imperial.ac.uk
> https://mailman.ic.ac.uk/mailman/listinfo/firedrake
>
-------------- next part --------------
HTML attachment scrubbed and removed