]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
adapt to loss of FaceInfo... and a little doc xhange in step 3
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 May 2010 15:46:47 +0000 (15:46 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 May 2010 15:46:47 +0000 (15:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@21136 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-12/step-12.cc
deal.II/examples/step-3/doc/intro.dox
deal.II/examples/step-3/step-3.cc
deal.II/examples/step-39/step-39.cc

index a14cd67338dc2cc5cfc8aada54db1dba2fbe889c..42c8692ba1a20e2b0566f04b020443567a1b0ab7 100644 (file)
@@ -200,8 +200,7 @@ class Step12
                                     // below.
     typedef MeshWorker::DoFInfo<dim> DoFInfo;
     typedef typename MeshWorker::IntegrationWorker<dim>::CellInfo CellInfo;
-    typedef typename MeshWorker::IntegrationWorker<dim>::FaceInfo FaceInfo;
-
+    
                                     // The following three functions
                                     // are then the ones that get called
                                     // inside the generic loop over all
@@ -231,9 +230,9 @@ class Step12
                                     // in fact other arguments
                                     // already bound.
     static void integrate_cell_term (DoFInfo& dinfo, CellInfo& info);
-    static void integrate_boundary_term (DoFInfo& dinfo, FaceInfo& info);
+    static void integrate_boundary_term (DoFInfo& dinfo, CellInfo& info);
     static void integrate_face_term (DoFInfo& dinfo1, DoFInfo& dinfo2,
-                                    FaceInfo& info1, FaceInfo& info2);
+                                    CellInfo& info1, CellInfo& info2);
 };
 
 
@@ -461,7 +460,7 @@ void Step12<dim>::integrate_cell_term (DoFInfo& dinfo, CellInfo& info)
                                 // FESubfaceValues, in order to get access to
                                 // normal vectors.
 template <int dim>
-void Step12<dim>::integrate_boundary_term (DoFInfo& dinfo, FaceInfo& info)
+void Step12<dim>::integrate_boundary_term (DoFInfo& dinfo, CellInfo& info)
 {
   const FEValuesBase<dim>& fe_v = info.fe_values();
   FullMatrix<double>& local_matrix = dinfo.matrix(0).matrix;
@@ -508,7 +507,7 @@ void Step12<dim>::integrate_boundary_term (DoFInfo& dinfo, FaceInfo& info)
                                 // back and forth.
 template <int dim>
 void Step12<dim>::integrate_face_term (DoFInfo& dinfo1, DoFInfo& dinfo2,
-                                      FaceInfo& info1, FaceInfo& info2)
+                                      CellInfo& info1, CellInfo& info2)
 {
                                   // For quadrature points, weights,
                                   // etc., we use the
index 222a263a565497a32048a75680fecc8ae0bf3a8d..2610c50350b40d0024e0c9c9d009082e83253753 100644 (file)
@@ -4,9 +4,9 @@
 
 This is the first example where we actually use finite elements to compute
 something. We 
-will solve a simple version of Laplace's equation with zero boundary
-values, but a nonzero right hand side. This example is still quite
-simple, but it already shows the basic structure of most finite
+will solve a simple version of Poisson's equation with zero boundary
+values, but a nonzero right hand side (here equal to the constant one).
+This example shows the basic structure of most finite
 element programs, which are along the following lines:
 <ul>
   <li> Grid generation;
index f87d65b39b0f6e209bbfe1408bc2ef6dc38afb6c..3a3f92a8be69e5d5a18aed44933b5548b3c49509 100644 (file)
@@ -3,7 +3,7 @@
 
 /*    $Id$       */
 /*                                                                */
-/*    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2006, 2007, 2008 by the deal.II authors */
+/*    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2006, 2007, 2008, 2010 by the deal.II authors */
 /*                                                                */
 /*    This file is subject to QPL and may not be  distributed     */
 /*    without copyright and license information. Please refer     */
@@ -218,8 +218,7 @@ void LaplaceProblem::make_grid_and_dofs ()
   std::cout << "Number of active cells: "
            << triangulation.n_active_cells()
            << std::endl;
-                                   // Here, by active we mean the cells on the
-                                  // finest level, i.e. cells that aren't
+                                   // Here, by active we mean the cells that aren't
                                   // refined any further.  We stress the
                                   // adjective `active', since there are more
                                   // cells, namely the parent cells of the
@@ -597,8 +596,10 @@ void LaplaceProblem::assemble_system ()
                                       // quadrature points of the
                                       // value of the shape function
                                       // at that point times the
-                                      // right hand side function
-                                      // (i.e. 1) times the Jacobian
+                                      // right hand side function,
+                                      // here the constant function
+                                      // equal to one,
+                                      // times the Jacobian
                                       // determinant times the weight
                                       // of that quadrature point:
       for (unsigned int i=0; i<dofs_per_cell; ++i)
index 2e6186cea3ee6a8e01805a2ce0b7d9e389043879..a686a13596397fcf9770b84f110fe7c829e8a6bc 100644 (file)
@@ -87,9 +87,7 @@ Functions::SlitSingularityFunction<2> exact_solution;
                                 // interior faces, respectively. All
                                 // the information needed for the
                                 // local integration is provided by
-                                // MeshWorker::IntegrationWorker<dim>::CellInfo
-                                // and
-                                // MeshWorker::IntegrationWorker<dim>::FaceInfo. Note
+                                // MeshWorker::IntegrationWorker<dim>::CellInfo. Note
                                 // that this public interface cannot
                                 // be changed, because it is expected
                                 // by MeshWorker::integration_loop().
@@ -100,11 +98,11 @@ class MatrixIntegrator : public Subscriptor
     static void cell(MeshWorker::DoFInfo<dim>& dinfo,
                     typename MeshWorker::IntegrationWorker<dim>::CellInfo& info);
     static void bdry(MeshWorker::DoFInfo<dim>& dinfo,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info);
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info);
     static void face(MeshWorker::DoFInfo<dim>& dinfo1,
                     MeshWorker::DoFInfo<dim>& dinfo2,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2);
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info1,
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info2);
 };
 
 
@@ -135,7 +133,7 @@ void MatrixIntegrator<dim>::cell(
 template <int dim>
 void MatrixIntegrator<dim>::bdry(
   MeshWorker::DoFInfo<dim>& dinfo,
-  typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info)
+  typename MeshWorker::IntegrationWorker<dim>::CellInfo& info)
 {
   const FEValuesBase<dim>& fe = info.fe_values();
   FullMatrix<double>& local_matrix = dinfo.matrix(0,false).matrix;
@@ -157,8 +155,8 @@ template <int dim>
 void MatrixIntegrator<dim>::face(
   MeshWorker::DoFInfo<dim>& dinfo1,
   MeshWorker::DoFInfo<dim>& dinfo2,
-  typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
-  typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2)
+  typename MeshWorker::IntegrationWorker<dim>::CellInfo& info1,
+  typename MeshWorker::IntegrationWorker<dim>::CellInfo& info2)
 {
   const FEValuesBase<dim>& fe1 = info1.fe_values();
   const FEValuesBase<dim>& fe2 = info2.fe_values();
@@ -208,11 +206,11 @@ class RHSIntegrator : public Subscriptor
 {
   public:
     static void cell(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::CellInfo& info);
-    static void bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info);
+    static void bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::CellInfo& info);
     static void face(MeshWorker::DoFInfo<dim>& dinfo1,
                     MeshWorker::DoFInfo<dim>& dinfo2,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2);
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info1,
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info2);
 };
 
 
@@ -222,7 +220,7 @@ void RHSIntegrator<dim>::cell(MeshWorker::DoFInfo<dim>&, typename MeshWorker::In
 
 
 template <int dim>
-void RHSIntegrator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info)
+void RHSIntegrator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::CellInfo& info)
 {
   const FEValuesBase<dim>& fe = info.fe_values();
   Vector<double>& local_vector = dinfo.vector(0).block(0);
@@ -244,8 +242,8 @@ void RHSIntegrator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWork
 template <int dim>
 void RHSIntegrator<dim>::face(MeshWorker::DoFInfo<dim>&,
                              MeshWorker::DoFInfo<dim>&,
-                             typename MeshWorker::IntegrationWorker<dim>::FaceInfo&,
-                             typename MeshWorker::IntegrationWorker<dim>::FaceInfo&)
+                             typename MeshWorker::IntegrationWorker<dim>::CellInfo&,
+                             typename MeshWorker::IntegrationWorker<dim>::CellInfo&)
 {}
 
 
@@ -254,11 +252,11 @@ class Estimator : public Subscriptor
 {
   public:
     static void cell(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::CellInfo& info);
-    static void bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info);
+    static void bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::CellInfo& info);
     static void face(MeshWorker::DoFInfo<dim>& dinfo1,
                     MeshWorker::DoFInfo<dim>& dinfo2,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
-                    typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2);
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info1,
+                    typename MeshWorker::IntegrationWorker<dim>::CellInfo& info2);
 };
 
 
@@ -278,7 +276,7 @@ void Estimator<dim>::cell(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::
 
 
 template <int dim>
-void Estimator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info)
+void Estimator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::CellInfo& info)
 {
   const FEValuesBase<dim>& fe = info.fe_values();
   
@@ -300,8 +298,8 @@ void Estimator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::
 template <int dim>
 void Estimator<dim>::face(MeshWorker::DoFInfo<dim>& dinfo1,
                          MeshWorker::DoFInfo<dim>& dinfo2,
-                         typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
-                         typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2)
+                         typename MeshWorker::IntegrationWorker<dim>::CellInfo& info1,
+                         typename MeshWorker::IntegrationWorker<dim>::CellInfo& info2)
 {
   const FEValuesBase<dim>& fe = info1.fe_values();
   const std::vector<double>& uh1 = info1.values[0][0];
@@ -333,7 +331,6 @@ class Step39
 {
   public:
     typedef typename MeshWorker::IntegrationWorker<dim>::CellInfo CellInfo;
-    typedef typename MeshWorker::IntegrationWorker<dim>::FaceInfo FaceInfo;
     
     Step39(const FiniteElement<dim>& fe);
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.