← Back to team overview

dolfin team mailing list archive

Re: Functions for applying boundary conditions to the RHS vector only now accept the solution vector x. Where Dirichlet boundary conditions are applied, b = x - bc.The is useful for applying boundary conditions within Newton iterations. [...]

 

On Fri, 2005-12-02 at 09:36 -0600, Anders Logg wrote:
> That sounds good.
> 
> One option that we should probably think about is to allow setting
> boundary conditions based on a Vector of degrees of freedom (similarly
> to what you do), so one does not need to supply an expression in terms
> of x, y, z. Maybe we can use the Function class for this, since
> it already has similar functionality?
> 

I don't see how a user could supply boundary conditions without an
expression in terms of x, y, z. 


> Then you would be able to have one applyBC() function but in your case
> you give it a special boundary condition that returns the residual.
> 
> I have not looked much yet at the new nonlinear solver so I don't
> understand all the issues. But I think it will become very useful to
> automate some of the solution process. I have some plans to let
> FFC generate a subclass of PDE (in addition to the classes it now
> generates), which would allow the following:
> 
>     MyFunction f; // Right-hand side
>     Function u;   // Solution
>     Poisson::PDE poisson(mesh, f, bc);
>     poisson.solve(u);
> 
> No need to create any vectors and attach to Functions 

This would be nice. 

> and no need to
> solve a linear (or nonlinear system), which is handled automatically
> by the PDE class. (It would choose the correct solver based on the
> given problem.)
> 

Sounds good, but I think that this is asking a lot for nonlinear
problems.


> By the way, any thoughts on which option is best, either
> 
>     u = poisson.solve();
>
> or
> 
>     poisson.solve(u);
> 

The first option could be good. It could then be 

u = poisson.solve(a, b, c);

where a, b, c are parameters for the solver (like intial time, dt, and
final time). 

> ?
> 
> Another thing I'm not completely happy with is that the subclass of
> PDE gets named Poisson::PDE (similarly to Poisson::BilinearForm etc in
> the Poisson namespace generated by FFC). Naturally, one would have the
> name PDE::Poisson. Another option would be to name it Equation, so
> that one can declare
> 
> Poisson::BilinearForm a;
> Poisson::LinearForm L;
> Poisson::Equation pde;
> 
Ok with me.

Garth

> Any thoughts?
> 
> /Anders
> 
> 
> On Fri, Dec 02, 2005 at 10:53:29AM +0100, Garth N. Wells wrote:
> > Now that I'm thinking properly, a function like FEM::assembleBCresidual should
> > return the residual for both Neumann and Dirichlet BCs (since non-zero Neumann
> > BCs are not yet implemented, I'm starting to forget that they even exist!)
> > 
> > I'll revert to the old applyBC to keep things consistent and add
> > FEM::assembleBCresidual which we can discuss if needed.
> > 
> > Garth  
> > 
> > 
> > 
> > Quoting Johan Hoffman <jhoffman@xxxxxxxxxxx>:
> > 
> > > I think that would be a good idea. And the name is an accurate description
> > > of what is going on.
> > > 
> > > /Johan
> > > 
> > > 
> > > > A possibility would be a function "FEM::assembleBCresidual".
> > > >
> > > > Garth
> > > >
> > > > On Fri, 2005-12-02 at 09:11 +0100, Garth N. Wells wrote:
> > > >> Quoting Anders Logg <logg@xxxxxxxxx>:
> > > >>
> > > >> > This looks strange to me, but maybe I don't understand it.
> > > >> >
> > > >> > What if you don't even have a vector x? What vector should x should
> > > >> > you then give as an argument?
> > > >>
> > > >> If you don't have a vector x, you don't want to use this function :-).
> > > >> It's
> > > >> basically returning the Dirichlet boundary condition "residual".
> > > >>
> > > >> >
> > > >> > Maybe it is better to just call the standard version of applyBC()
> > > >> > without x and then do
> > > >> >
> > > >> > b *= -1.0;
> > > >> > b += x;
> > > >>
> > > >> This would add x to all elements of b, not just to elements where
> > > >> Dirichlet
> > > >> boundary conditions are applied. I originally used the old applyBC for a
> > > >> vector
> > > >> to use an analogous procedure to the above, but this required two
> > > >> temporary
> > > >> vectors to be created, initialised and have boundary conditions applied
> > > >> at each
> > > >> Newton iteration (in addtion to the usual assemble of b) - not very
> > > >> nice.
> > > >>
> > > >> >
> > > >> > This would be two lines of code extra (one if you put them on the same
> > > >> > line... :-). It's probably better to put this in the Newton solver
> > > >> > since the Newton solver knows why this operation should be performed.
> > > >>
> > > >> The Newton solver knows nothing about the boundary coniditons. The BC's
> > > >> need to
> > > >> be applied by the user in the resdiudal vector b. This is done in the
> > > >> class
> > > >> NonlinearFunction.
> > > >>
> > > >> A possibility is a function that inserts components of x into b (rather
> > > >> than the
> > > >> BC value), where Dirichlet BCs are applied. PETSc has a function to the
> > > >> set the
> > > >> vector R, where F(u) = R is being solved, so R would contain the
> > > >> Dirichlet BCs.
> > > >>
> > > >>
> > > >> Garth
> > > >>
> > > >> >
> > > >> > /Anders
> > > >> >
> > > >> >
> > > >> > On Thu, Dec 01, 2005 at 02:16:12PM +0100, dolfin@xxxxxxxxxx wrote:
> > > >> > > Commit from garth (2005-12-01 14:16 CET)
> > > >> > > -----------------
> > > >> > >
> > > >> > > Functions for applying boundary conditions to the RHS vector only
> > > >> > now accept the solution vector x. Where Dirichlet boundary conditions
> > > >> > are applied, b = x - bc.The is useful for applying boundary conditions
> > > >> > within Newton iterations.
> > > >> > >
> > > >> > > This is not consistent with the other applyBC functions which return
> > > >> b = bc
> > > >> > (no minus sign). This is due to the definition of the nonlinear
> > > >> function
> > > >> > F(u)=0.
> > > >> > > For Poisson's equation, I've been using the format F(u) = (grad
> > > >> v,grad u) -
> > > >> > (v,f), hence the need for the minus. Is there a preference as to how
> > > >> we
> > > >> > set-up nonlinear problems?
> > > >> > >
> > > >> > > Would it be better to use a different name for funtions that return
> > > >> the
> > > >> > difference between the current approximate solution and the Dirichlet
> > > >> BC?
> > > >> > >
> > > >> > >   dolfin  src/kernel/fem/FEM.cpp       1.43
> > > >> > >   dolfin  src/kernel/fem/dolfin/FEM.h  1.23
> > > >> > >
> > > >> > > _______________________________________________
> > > >> > > DOLFIN-dev mailing list
> > > >> > > DOLFIN-dev@xxxxxxxxxx
> > > >> > > http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> > > >> > >
> > > >> >
> > > >> >
> > > >> > _______________________________________________
> > > >> > DOLFIN-dev mailing list
> > > >> > DOLFIN-dev@xxxxxxxxxx
> > > >> > http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> > > >> >
> > > >>
> > > >> _______________________________________________
> > > >> DOLFIN-dev mailing list
> > > >> DOLFIN-dev@xxxxxxxxxx
> > > >> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> > > >
> > > > tel.     +31 15 278 7922
> > > > fax.     +31 15 278 6383
> > > > e-mail   g.n.wells@xxxxxxxxxxxxxxx
> > > > url      http://www.mechanics.citg.tudelft.nl/~garth
> > > >
> > > >
> > > >
> > > > _______________________________________________
> > > > DOLFIN-dev mailing list
> > > > DOLFIN-dev@xxxxxxxxxx
> > > > http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> > > >
> > > 
> > > 
> > > 
> > > _______________________________________________
> > > DOLFIN-dev mailing list
> > > DOLFIN-dev@xxxxxxxxxx
> > > http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> > > 
> > 
> > _______________________________________________
> > DOLFIN-dev mailing list
> > DOLFIN-dev@xxxxxxxxxx
> > http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> > 
> 





Follow ups

References