dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Wed, Mar 21, 2007 at 06:24:37PM +0100, Garth N. Wells wrote:
> 
> 
> Anders Logg wrote:
> > On Wed, Mar 21, 2007 at 06:16:04PM +0100, Garth N. Wells wrote:
> >>
> >> Anders Logg wrote:
> >>> On Wed, Mar 21, 2007 at 04:58:19PM +0100, Garth N. Wells wrote:
> >>>> Anders Logg wrote:
> >>>>> On Wed, Mar 21, 2007 at 02:46:48PM +0100, Garth N. Wells wrote:
> >>>>>> For the new assembly, should we create a new class for applying 
> >>>>>> Dirichlet boundary conditions? Alternatively, we could add functions to 
> >>>>>> BoundaryCondition to do the job.
> >>>>> Yes, why not, but how do we want to specify boundary conditions?
> >>>>>
> >>>>> There are some options:
> >>>>>
> >>>>>    MyBoundaryCondition bc;
> >>>>>    bc.apply(A, b);
> >>>>>
> >>>>> or
> >>>>>
> >>>>>    assemble(A, b, ..., bc):
> >>>>>
> >>>> I prefer the latter as it allows for restrained degrees of freedom to be 
> >>>> eliminated from the system, rather than zeroing rows. While not a high 
> >>>> priority, I'd like to have this in the future.
> >>>>
> >>>> Should BoundaryCondition posses the members functions to apply the 
> >>>> boundary conditions? I'm reluctant to add these to Assembler in order to 
> >>>> keep Assembler clean and just for assembly.
> >>> I'm not sure I follow, didn't you say you wanted option 2 (putting it
> >>> in the assembler rather than in the BoundaryCondition class)?
> >>>
> >> The question is which class knows how to apply boundary conditions. I 
> >> was thinking something like
> >>
> >>    Assembler::assemble(A, b, .. . , *bc)
> >>    {
> >>      .
> >>      .
> >>      .
> >>      if( bc )
> >>        bc.applyBC(A);
> >>
> >>    }
> > 
> > It looks a little strange to me. To keep it simple, I'd suggest we do
> > as you say but don't involve the assembler:
> > 
> >     assemble(A, b, ...);
> >     bc.apply(A, b);
> > 
> > or bc.apply(A), bc.apply(b) etc.
> > 
> > Then the assembler only knows how to assemble and the boundary
> > conditions are handled elsewhere. (And when we wrap this into
> > LinearSolver or NonlinearSolver this is all taken care of.)
> > 
> > Note that bc.apply(A, b) can eliminate Dirichlet dofs as an
> > alternative to putting zeros into the matrix.
> > 
> 
> OK. bc.apply(A, b) can eliminate Dirichlet dofs, but it is more complex 
> to do this after assembly rather than during assembly.

Perhaps more efficient (not needing to recreate the matrix) but my
feeling is that doing it after is easier since then the assembly does not
need to care about skipping rows or columns.

/Anders


> Garth
> 
> > /Anders
> > 
> > 
> >>>>>> Also, for the new assembly, is there anything other than boundary 
> >>>>>> conditions and updating of functions that needs to be done before it's 
> >>>>>> usable?
> >>>>> I think that's it for the basic functionality of assembly. Then there
> >>>>> is the question about how we want to design the interface.
> >>>>>
> >>>>> Right now, FFC (new UFC-based) generates one .h file for each .form
> >>>>> file and each .form file can contain exactly one form, which must have
> >>>>> the name 'a'. 
> >>>>>
> >>>>> Take a look at the following lines, this is what we have now:
> >>>>>
> >>>>> 1. Poisson.form (form file)
> >>>>>
> >>>>>     a = dot(grad(v), grad(u))*dx
> >>>>>
> >>>>> 2. Poisson.h (generated with ffc -l ufc Poisson.form)
> >>>>>
> >>>>>     class Poisson: public ufc::form
> >>>>>     { 
> >>>>>       ...
> >>>>>     }
> >>>>>
> >>>>> 3. main.cpp (includes generated header file)
> >>>>>
> >>>>>     #include "Poisson.h"
> >>>>>
> >>>>>     ...
> >>>>>
> >>>>>     Poisson poisson;
> >>>>>     assemble(A, poisson, mesh);
> >>>>>
> >>>>> We have a few things to settle (which we could do pretty quickly):
> >>>>>
> >>>>> 1. Is the principle one form file - one form ok? Or do we want to be
> >>>>> able to define several forms in one file (a, L and M)?
> >>>>>
> >>>>> With one form per file, we would have PoissonBilinearForm.form and
> >>>>> PoissonLinearForm.form and then
> >>>>>
> >>>>>     #include "PoissonBilinearForm.h"
> >>>>>     #include "PoissonLinearForm.h"
> >>>>>
> >>>>>     PoissonBilinearForm a;
> >>>>>     PoissonBilinearForm L;
> >>>>>
> >>>> I prefer having multiple forms in one file.
> >>> ok, but then the number of forms will be limited to a maximum of three
> >>> forms per file:
> >>>
> >>>   bilinear form a --> FooBilinearForm
> >>>   linear form L   --> FooLinearForm
> >>>   functional M    --> FooFunctional
> >>>
> >>> Is that ok?
> >>>
> >> Yes.
> >>
> >> Garth
> >>
> >>> FFC must generate the names for the classes from the name of the .form
> >>> file so this puts a limitation on which names are available. Perhaps
> >>> one could generate Foo_0, Foo_1, Foo_2, Foo_3 if there is more than
> >>> exactly one form of each type, but then the naming scheme starts to
> >>> get complicated...
> >>>
> >>>>> 2. How do we want to attach fixed functions to forms, like the
> >>>>> right-hand side in Poisson? Before, we had
> >>>>>
> >>>>>     Function f;
> >>>>>     PoissonLinearForm L(f);
> >>>>>
> >>>>> The problem with this is that FFC now generates UFC-compliant code
> >>>>> (only including the ufc.h header). It does not include anything from
> >>>>> DOLFIN (like Function.h) which would be necessary for this to work.
> >>>>>
> >>>>> One option that I've been thinking of is to have an option to FFC
> >>>>> -l dolfin, which would generate UFC code and at the end insert code
> >>>>> that wraps the UFC form in a DOLFIN Form:
> >>>>>
> >>>>>     class PoissonLinearForm_UFC : public ufc::form
> >>>>>     {
> >>>>>       ...
> >>>>>     };
> >>>>>
> >>>>>     #include <dolfin/Function.h>
> >>>>>     #include <dolfin/Form.h>
> >>>>>     class PoissonLinearForm : public dolfin::Form
> >>>>>     {
> >>>>>     public:
> >>>>>
> >>>>>       PoissonLinearForm(Function& f) : Form(1)
> >>>>>       {
> >>>>>           init(0, f);
> >>>>>       }
> >>>>>       
> >>>>>     };
> >>>>>
> >>>>> Opinions?
> >>>>>
> >>>> Not sure what's best here. I like the above approach as nothing changes 
> >>>> on the DOLFIN side, but it does require FFC to be more than just UFC 
> >>>> compliant.
> >>> ok, I think this is the best way to do it. It hides the UFC classes
> >>> from the user which is probably good anyway.
> >>>
> >>> I think there should be two assemble() functions:
> >>>
> >>>     assemble(..., ufc::form& form, Array<Function*> coefficients, ...);
> >>>     assemble(..., dolfin::Form& form, ...)
> >>>
> >>> One that takes a ufc::form and a set of coefficients and then another
> >>> one that takes just a Form (which contains the coefficients).
> >>>
> >>> The Form class can have two member functions form() and coefficients()
> >>> that just return the ufc::form and the array of Functions so the
> >>> second assemble function can be really short.
> >>>
> >>> /Anders
> >>> _______________________________________________
> >>> DOLFIN-dev mailing list
> >>> DOLFIN-dev@xxxxxxxxxx
> >>> http://www.fenics.org/mailman/listinfo/dolfin-dev
> >>>
> >>
> >> _______________________________________________
> >> DOLFIN-dev mailing list
> >> DOLFIN-dev@xxxxxxxxxx
> >> http://www.fenics.org/mailman/listinfo/dolfin-dev
> > _______________________________________________
> > DOLFIN-dev mailing list
> > DOLFIN-dev@xxxxxxxxxx
> > http://www.fenics.org/mailman/listinfo/dolfin-dev
> > 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev