dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Thu, Jul 03, 2008 at 06:07:08PM +0100, Garth N. Wells wrote:
> 
> 
> Kent-Andre Mardal wrote:
> > On to., 2008-07-03 at 08:03 +0100, Garth N. Wells wrote:
> >> DOLFIN wrote:
> >>> One or more new changesets pushed to the primary dolfin repository.
> >>> A short summary of the last three changesets is included below.
> >>>
> >>> changeset:   4374:6e1db926186ebfe9d9316b729a7d5747be3c5115
> >>> tag:         tip
> >>> parent:      4373:8f4accd9eb039d0d13412e8d871e7f1bfce896b8
> >>> parent:      4372:3abbccc877d3bb358c76ef36dbbe252155b88183
> >>> user:        Kent-Andre Mardal <kent-and@xxxxxxxxx>
> >>> date:        Wed Jul 02 16:02:17 2008 +0200
> >>> files:       dolfin/fem/Assembler.cpp dolfin/fem/DirichletBC.cpp site-packages/dolfin/assemble.py
> >>> description:
> >>> merge
> >>>
> >>>
> >>> changeset:   4373:8f4accd9eb039d0d13412e8d871e7f1bfce896b8
> >>> parent:      4369:b272b94b43ae789e2c72c90f6a027ce9fb89f6fc
> >>> user:        Kent-Andre Mardal <kent-and@xxxxxxxxx>
> >>> date:        Wed Jul 02 16:01:09 2008 +0200
> >>> files:       demo/pde/sym-dirichlet-bc/demo.py dolfin/fem/Assembler.cpp dolfin/fem/Assembler.h dolfin/fem/DirichletBC.cpp dolfin/fem/DirichletBC.h dolfin/fem/assemble.cpp dolfin/fem/assemble.h dolfin/la/AssemblyMatrix.cpp dolfin/la/AssemblyMatrix.h dolfin/la/EpetraMatrix.cpp dolfin/la/EpetraMatrix.h dolfin/la/EpetraVector.cpp dolfin/la/EpetraVector.h dolfin/la/GenericMatrix.h dolfin/la/GenericTensor.h dolfin/la/GenericVector.h dolfin/la/Matrix.h dolfin/la/PETScMatrix.cpp dolfin/la/PETScMatrix.h dolfin/la/PETScVector.cpp dolfin/la/PETScVector.h dolfin/la/Scalar.h dolfin/la/Vector.h dolfin/la/uBlasMatrix.h dolfin/la/uBlasVector.cpp dolfin/la/uBlasVector.h dolfin/swig/dolfin_fem_pre.i dolfin/swig/dolfin_headers.i site-packages/dolfin/assemble.py
> >>> description:
> >>> implemented symmetric handling of boundary conditions
> >>>
> >> Great.
> >>
> >> I looked at the demo but could see how to specify the method for 
> >> applying boundary conditions. Could you write a short description on how 
> >> to use it and it compares in terms of performance?
> >>
> >> Garth
> >>
> > 
> > Can you explain to me what the different methods do ? I assume you refer
> > to the topological, geometric, pointwise switch in computeBC (which
> > is used). 
> > 
> 
> I haven't looked at it in detail for some time, so I have to defer to 
> Anders and Kristian for a good explanation.

We have implemented 3 different strategies for locating the dofs for
which strong boundary conditions should be enforced:

1. Topological search

This will iterate over all boundary facets, and for each facet
tabulate the dofs that are on the facet using the UFC function
tabulate_facet_dofs. The BCs will be enforced for all these dofs.

This means that BCs will be enforced for all dofs belonging to a facet
which is *entirely on the boundary*. As a consequence, BCs will not be
enforced for dofs that happen to lie on the boundary but which belong
to a facet which is not on the boundary.

For example, BCs will not be enforced for some vertex dofs that are on
the boundary. This does not matter for continuous elements, since a
vertex dof will belong to at least one facet which is on the boundary
(even if it also belongs to another facet which is not on the
boundary). However, for DG methods this does matter.

2. Geometric search

This will iterate over all boundary facets and locate all nearby dofs
which are on the facet, using the UFC function tabulate_coordinates
and test the location of the dofs geometrically against the location of
the facet. This is slower than the topological search, but catches
the dofs that the topological search may miss.

3. Pointwise

This will iterate over all cells (not only the boundary facets) of the
entire mesh and call the inside() function of a SubDomain to check
which dofs should be constrained. This allows setting boundary
conditions pointwise (even inside the domain).

-- 
Anders


> > What is done is the following: 
> > After the global matrix is computed, computeBC is used to compute the
> > boundary conditions. 
> > Then the BC is enforced as follows: 
> > 
> > for all element: 
> >    get element matrix and element vector 
> >    enforce Dirichlet BC locally 
> >    symmetrice the element matrix (and update the element rhs) by 
> >      partial Gaussian elimination on the element linear system 
> >    set the element matrix and element vector into the global matrix 
> > 
> 
> So if I pass the boundary conditions to the assembler it will do this?
> 
> Is this method of applying bcs taken into account when computing the 
> sparsity pattern?
> 
> > 
> > Hence the assembly loop is as follows: 
> > 
> > 
> > for all elements: 
> >    compute cell contributions
> > 
> > for all elements: 
> >    compute interior facet contributions
> > 
> > for all elements: 
> >    compute exterior facet  contributions
> > 
> > for all elements: 
> >    enforce symmetric bc 
> > 
> 
> Isn't this different from what's described above? I thought that the 
> Dirichlet bcs were being imposed within after the computation of each 
> cell contribution?
> 
> > 
> > 
> > which can be more efficiently done as: 
> > 
> > for all elements: 
> >    compute cell contributions
> >    compute interior facet contributions
> >    compute exterior facet  contributions
> >    enforce symmetric bc 
> > 
> > 
> > Concerning efficiency the PETSc implemententation may be slow since
> > 
> >    MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
> >    MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
> > 
> > must be called between each get and set call. 
> > 
> 
> Why is a get call needed? Don't you just assemble the reduced system for 
> each cell?
> 
> Garth
> 
> > 
> > 
> > Kent
> > 
> > 
> > 
> > 
> > 
> > 
> > 
> > 
> > 
> > 
> 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev

Attachment: signature.asc
Description: Digital signature