tutorial updated

extensions/tutorial.hpp

@@ -706,8 +706,172 @@ refinement.refine(...) call the mesh is globally refined up to the <c>initial le
 function and especially the interface can be represented by the mesh.
 
 */
+
 //-----------------------------------------------------------
 
 
+/*! \page boundary_conditions Boundary conditions
+
+For a second order elliptic partial differential equation one has to describe the geometry and prescribe boundary conditions to close
+the system. Three types of boundary conditions are very typical: Dirichlet-type, Neumann-type or the Robin-type boundary condition and 
+the periodic boundary conditions. AMDiS provides interfaces for all these types. In the macro-mesh one has to mark several edges as 
+boundary edges. During refinement the mark is handed down to the leaf elements. That's the crux of the matter. Boundary conditions are
+defined edge-wise, but what to do if two edges with different boundary conditions collide? In the macro-mesh one has to set for each
+separate boundary part a boundary number. Depending on this number the one or other condition is preferred for the colliding vertex.
+For problems with more than one equation on the same macro-mesh one can construct examples that are solvable with this strategy. Also
+there are problem with the periodic boundary condition and mixed finite elements. Neumann boundary conditions in complicated parts of
+a boundary need also special consideration. Here some new techniques are introduced to overcome some of the problems.
+- Dirichlet boundary conditions
+  - singular points
+  - implicit areas
+- Periodic boundary conditions
+- Implicit Neumann condition
+
+\section dirichlet_bc Dirichlet boundary conditions
+There are three methods to enforce Dirichlet boundary condition, i.e. \f$ u = g\; on\; \Gamma \f$. To use these you have to replace the 
+standard ProblemStat by ExtendedProblemStat. The methods <c>addSingularDirichletBC(...)</c> add a condition at only one single degree
+of freedom. For the second method <c>addImplicitDirichletBC(...)</c> you have to pass a function, that has negative values in the 
+vertices where to enforce the boundary condition, e.g. a signed distance function. The last method <c>addManualDirichletBC(...)</c>
+provides an interface to describe parts of a boundary easily.
+
+\subsection singular_dirichlet_bc Singular Dirichlet boundary condition
+
+For example we want to solve the Neumann problem, we have to define an additional condition, since \f$ u \f$ is only defined up to a constant.
+
+\f[ \Delta u = f\quad in\; \Omega,\quad\partial_n u = 0\quad on\; \partial\Omega \f]
+
+One way out of this dilemma is to set a Dirichlet boundary condition on just one point. There are two ways to do this: first you can give a
+dofindex for which to set the value, second you can specify a coordinate and the algorithm looks for the nearest degree of freedom to set
+the boundary condition there.
+
+\code
+template<typename ValueContainer>
+  void addSingularDirichletBC(DegreeOfFreedom idx, int row, int col,
+			      ValueContainer &values);
+template<typename ValueContainer>
+  void addSingularDirichletBC(WorldVector<double> &pos, int row, int col,
+			      ValueContainer &values);
+\endcode
+
+You have to the row and column where to set the boundary condition in the system of differential equations and as last argument you have to give
+the value that should be enforced. The <c>ValueContainer</c> can ether be scalar reference, i.e. <c>double&amp;</c>, an 
+<c>AbstractFunction&lt; double, WorldVector&lt;double&gt; &gt;</c>, or a <c>DOFVector&lt; double &gt;</c>.
+
+The next example shows how to use this methods:
+
+\code
+int main(int argc, char* argv[])
+{
+  AMDiS::init(argc, argv);
+
+  // ===== create and init the scalar problem ===== 
+  ExtendedProblemStat ellipt("ellipt");
+  ellipt.initialize(INIT_ALL);
+
+  // === create adapt info ===
+  AdaptInfo adaptInfo("ellipt->adapt", ellipt.getNumComponents());
+  AdaptStationary adapt("ellipt->adapt", ellipt, adaptInfo);
+  
+  // ===== create matrix operator =====
+  Operator matrixOperator(ellipt.getFeSpace());
+  matrixOperator.addTerm(new Simple_SOT);
+  ellipt.addMatrixOperator(matrixOperator, 0, 0);
+
+  // ===== create rhs operator =====
+  Operator rhsOperator(ellipt.getFeSpace());
+  rhsOperator.addTerm(new Simple_ZOT(1.0));
+  ellipt.addVectorOperator(rhsOperator, 0);
+
+  // ===== add singular dirichlet boundaty condition =====
+  WorldVector<double> p; p.set(0.0);
+  double zero = 0.0;
+  ellipt.addSingularDirichletBC(p, 0, 0, zero);
+
+  // ===== start simulation =====
+  adapt.adapt();
+  ellipt.writeFiles(adaptInfo, true);
+
+  AMDiS::finalize();
+}
+\endcode
+
+The solution will look like
+\image html elliptSingular.png "Singular Dirichlet condition in the bottom left corner" width=300px
+
+\subsection implicit_dirichlet_bc Implicit Dirichlet boundary condition
+
+Assume you want to solve an equation in a "complicated" domain, that can be described by a signed distance function. Then
+you can use a Diffuse-Domain approach or the implicit dirichlet boundary conditions defined in ExtendedProblemStat.h. Two
+interfaces are available. You can set a DOFVector that holds the values of the signed distance function, or you can specify
+an AbstractFunction:
+
+\code
+template<typename ValueContainer>
+  void addImplicitDirichletBC(DOFVector<double> &signedDist,
+			      int row, int col,
+			      ValueContainer &values);
+
+template<typename ValueContainer>
+  void addImplicitDirichletBC(AbstractFunction<double, WorldVector<double> > &signedDist,
+			      int row, int col,
+			      ValueContainer &values);
+\endcode
+
+Similar to the singular Dirichlet conditions the value container can be a scalar, an AbstractFunction, or a DOFVector. In the next Example
+the signed distance functions describes a disk area. At the boundary the Dirichlet condition should be enforced. The algorithm now not only
+set the Dirichlet values at the boundary, but in the whole area with distance values less than zero. So a simple algebraic equation is solved
+in the outer part of the complicated domain. Near the boundary the mesh should be refined, since one wants to resolve the boundary vertices 
+very near to the implicit interface.
+
+\code
+int main(int argc, char* argv[])
+{
+  AMDiS::init(argc, argv);
+
+  // ===== create and init the scalar problem ===== 
+  ExtendedProblemStat ellipt("ellipt");
+  ellipt.initialize(INIT_ALL);
+
+  // === create adapt info ===
+  AdaptInfo adaptInfo("ellipt->adapt", ellipt.getNumComponents());
+  AdaptStationary adapt("ellipt->adapt", ellipt, adaptInfo);
+
+  // ===== signedDist function that describes the geometry =====
+  AbstractFunction<double, WorldVector<double> > *signedDistFct = new InverseCircle(1.0);
+
+  // ===== refine mesh at interface =====
+  SignedDistRefinement refFunction(ellipt.getMesh());
+  RefinementLevelCoords2 refinement(
+    ellipt.getFeSpace(),
+    &refFunction,
+    signedDistFct);
+  refinement.refine(9);
+  
+  // ===== create matrix operator =====
+  Operator matrixOperator(ellipt.getFeSpace());
+  matrixOperator.addTerm(new Simple_SOT);
+  ellipt.addMatrixOperator(matrixOperator, 0, 0);
+
+  // ===== create rhs operator =====
+  Operator rhsOperator(ellipt.getFeSpace());
+  rhsOperator.addTerm(new Simple_ZOT(1.0));
+  ellipt.addVectorOperator(rhsOperator, 0);
+
+  // ===== add boundary conditions =====  
+  double zero = 0.0;
+  ellipt.addImplicitDirichletBC(*signedDistFct, 0, 0, zero);
+
+  // ===== start simulation =====
+  adapt.adapt();
+  ellipt.writeFiles(adaptInfo, true);
+
+  AMDiS::finalize();
+}
+\endcode
+
+The solution is shown in the next pictures. On the left one can see the mesh locally refined on the boundary of a circle:
+\image html elliptImplicit.png "Implicit Dirichlet condition at the boundary of a circle" width=300px
+
+*/
 
 #endif // AMDIS_EXTENSIONS_TUTORIAL_INCLUDE