From: Wolfgang Bangerth Date: Fri, 22 May 1998 09:31:45 +0000 (+0000) Subject: Resolve some problems with the projection of functions. X-Git-Tag: v8.0.0~22920 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=402d5838619aa8c9f6448ab9422403e1f00a9449;p=dealii.git Resolve some problems with the projection of functions. git-svn-id: https://svn.dealii.org/trunk@334 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/numerics/vectors.h b/deal.II/deal.II/include/numerics/vectors.h index 81c4812ae6..f90258bc52 100644 --- a/deal.II/deal.II/include/numerics/vectors.h +++ b/deal.II/deal.II/include/numerics/vectors.h @@ -12,6 +12,7 @@ template class DoFHandler; template class Function; template class Quadrature; +template class QGauss2; template class FiniteElement; template class Boundary; template class StraightBoundary; @@ -82,23 +83,44 @@ enum NormType { * $f_i = \int_\Omega f(x) \phi_i(x) dx$. The solution vector $v$ then is * the projection. * - * In order to get proper results, it is necessary to treat boundary - * conditions right. This is done by $L_2$-projection of the trace of the + * In order to get proper results, it may necessary to treat boundary + * conditions right. Below are listed some cases wher this may be needed. + * If needed, this is done by $L_2$-projection of the trace of the * given function onto the finite element space restricted to the boundary * of the domain, then taking this information and using it to eliminate * the boundary nodes from the mass matrix of the whole domain, using the * #MatrixTools::apply_boundary_values# function. The projection of the * trace of the function to the boundary is done with the - * #VectorTools::project_boundary_values# (see below) function. You may - * specify a flag telling the projection that the function has zero boundary - * values, in which case the $L_2$-projection onto the boundary is not - * needed. If it is needed, the #VectorTools::project_boundary_values# is - * called with a map of boundary functions of which all boundary indicators + * #VectorTools::project_boundary_values# (see below) function, which is + * called with a map of boundary functions in which all boundary indicators * from zero to 254 (255 is used for other purposes, see the #Triangulation# * class documentation) point to the function to be projected. The projection * to the boundary takes place using a second quadrature formula given to * the #project# function. * + * The projection of the boundary values first, then eliminating them from + * the global system of equations is not needed usually. It may be necessary + * if you want to enforce special restrictions on the boundary values of the + * projected function, for example in time dependant problems: you may want + * to project the initial values but need consistency with the boundary + * values for later times. Since the latter are projected onto the boundary + * in each time step, it is necessary that we also project the boundary + * values of the initial values, before projecting them to the whole domain. + * + * The selection whether the projection to the boundary first is needed is + * done with the #project_to_boundary_first# flag passed to the function. + * If #false# is given, the additional quadrature formula for faces is + * ignored. + * + * You should be aware of the fact that if no projection is to the boundary + * is requested, a function with zero boundary values may not have zero + * boundary values after projection. There is a flag for this especially + * important case, which tells the function to enforce zero boundary values + * on the respective boundary parts. Since enforced zero boundary values + * could also have been reached through projection, but are more economically + * obtain using other methods, the #project_to_boundary_first# flag is + * ignored if the #enforce_zero_boundary# flag is set. + * * The solution of the linear system is presently done using a simple CG * method without preconditioning and without multigrid. This is clearly not * too efficient, but sufficient in many cases and simple to implement. This @@ -241,6 +263,15 @@ class VectorTools { * Compute the projection of * #function# to the finite element space. * + * By default, projection to the boundary + * and enforcement of zero boundary values + * are disabled. The ordering of arguments + * to this function is such that you need + * not give a second quadrature formula if + * you don't want to project to the + * boundary first, but that you must if you + * want to do so. + * * See the general documentation of this * class for further information. */ @@ -248,11 +279,12 @@ class VectorTools { const ConstraintMatrix &constraints, const FiniteElement &fe, const Quadrature &q, - const Quadrature &q_boundary, const Boundary &boundary, const Function &function, - const bool has_zero_boundary, - dVector &vec); + dVector &vec, + const bool enforce_zero_boundary = false, + const Quadrature &q_boundary = QGauss2(), + const bool project_to_boundary_first = false); /** * Make up the list of node subject diff --git a/deal.II/deal.II/source/numerics/vectors.cc b/deal.II/deal.II/source/numerics/vectors.cc index c207759144..358af5ac37 100644 --- a/deal.II/deal.II/source/numerics/vectors.cc +++ b/deal.II/deal.II/source/numerics/vectors.cc @@ -75,11 +75,12 @@ void VectorTools<1>::project (const DoFHandler<1> &, const ConstraintMatrix &, const FiniteElement<1> &, const Quadrature<1> &, - const Quadrature<0> &, const Boundary<1> &, const Function<1> &, + dVector &, const bool , - dVector &) { + const Quadrature<0> &, + const bool ) { // this function should easily be implemented // using the template below. However some // changes have to be made since faces don't @@ -97,28 +98,19 @@ void VectorTools::project (const DoFHandler &dof, const ConstraintMatrix &constraints, const FiniteElement &fe, const Quadrature &q, - const Quadrature &q_boundary, const Boundary &boundary, const Function &function, - const bool has_zero_boundary, - dVector &vec) { + dVector &vec, + const bool enforce_zero_boundary, + const Quadrature &q_boundary, + const bool project_to_boundary_first) { // make up boundary values map boundary_values; - if (has_zero_boundary == false) - { - // set up a list of boundary functions for - // the different boundary parts. We want the - // #function# to hold on all parts of the - // boundary - FunctionMap boundary_functions; - for (unsigned char c=0; c<255; ++c) - boundary_functions[c] = &function; - project_boundary_values (dof, boundary_functions, fe, q_boundary, - boundary, boundary_values); - } - else - // no need to project boundary values + if (enforce_zero_boundary == true) + // no need to project boundary values, but + // enforce homogeneous boundary values + // anyway { DoFHandler::active_face_iterator face = dof.begin_active_face(), endf = dof.end_face(); @@ -131,7 +123,23 @@ void VectorTools::project (const DoFHandler &dof, // for all boundary nodes boundary_values[face_dof_indices[i]] = 0.; }; - }; + } + else + // no homogeneous boundary values + if (project_to_boundary_first == true) + // boundary projection required + { + // set up a list of boundary functions for + // the different boundary parts. We want the + // #function# to hold on all parts of the + // boundary + FunctionMap boundary_functions; + for (unsigned char c=0; c<255; ++c) + boundary_functions[c] = &function; + project_boundary_values (dof, boundary_functions, fe, q_boundary, + boundary, boundary_values); + }; + // set up mass matrix and right hand side