]> https://gitweb.dealii.org/ - dealii.git/commitdiff
merge from mainline
authorTimo Heister <timo.heister@gmail.com>
Tue, 27 Nov 2012 15:51:16 +0000 (15:51 +0000)
committerTimo Heister <timo.heister@gmail.com>
Tue, 27 Nov 2012 15:51:16 +0000 (15:51 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_unify_linear_algebra@27697 0785d39b-7218-0410-832d-ea1e28bc413d

30 files changed:
1  2 
deal.II/examples/step-13/step-13.cc
deal.II/examples/step-14/step-14.cc
deal.II/examples/step-18/step-18.cc
deal.II/examples/step-20/step-20.cc
deal.II/examples/step-21/step-21.cc
deal.II/examples/step-22/step-22.cc
deal.II/examples/step-23/step-23.cc
deal.II/examples/step-28/step-28.cc
deal.II/examples/step-29/step-29.cc
deal.II/examples/step-31/step-31.cc
deal.II/examples/step-32/step-32.cc
deal.II/examples/step-33/step-33.cc
deal.II/examples/step-35/step-35.cc
deal.II/examples/step-37/step-37.cc
deal.II/examples/step-40/step-40.cc
deal.II/examples/step-42/step-42.cc
deal.II/examples/step-43/step-43.cc
deal.II/examples/step-46/step-46.cc
deal.II/examples/step-47/step-47.cc
deal.II/examples/step-6/step-6.cc
deal.II/examples/step-9/step-9.cc
deal.II/include/deal.II/base/exceptions.h
deal.II/include/deal.II/dofs/dof_tools.h
deal.II/include/deal.II/lac/petsc_solver.h
deal.II/include/deal.II/lac/sparse_matrix.h
deal.II/include/deal.II/matrix_free/fe_evaluation.h
deal.II/include/deal.II/numerics/vector_tools.templates.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/fe/fe_q.cc
deal.II/source/lac/petsc_solver.cc

index 63e62cc64bae6f2e2f4cf446b644e2708646b422,b85cb3d50a0a129499e1ee0b7ecda1e40aa29ea1..b180dc832547002dd249773488b8cad759b185d4
@@@ -292,32 -203,19 +203,19 @@@ namespace Step1
      void
      PointValueEvaluation<dim>::
      operator () (const DoFHandler<dim> &dof_handler,
 -                 const Vector<double> &solution) const
 +                 const Vector<double>  &solution) const
      {
-       // First allocate a variable that
-       // will hold the point
-       // value. Initialize it with a
-       // value that is clearly bogus,
-       // so that if we fail to set it
-       // to a reasonable value, we will
-       // note at once. This may not be
-       // necessary in a function as
-       // small as this one, since we
-       // can easily see all possible
-       // paths of execution here, but
-       // it proved to be helpful for
-       // more complex cases, and so we
-       // employ this strategy here as
-       // well.
+       // First allocate a variable that will hold the point value. Initialize
+       // it with a value that is clearly bogus, so that if we fail to set it
+       // to a reasonable value, we will note at once. This may not be
+       // necessary in a function as small as this one, since we can easily see
+       // all possible paths of execution here, but it proved to be helpful for
+       // more complex cases, and so we employ this strategy here as well.
        double point_value = 1e20;
  
-       // Then loop over all cells and
-       // all their vertices, and check
-       // whether a vertex matches the
-       // evaluation point. If this is
-       // the case, then extract the
-       // point value, set a flag that
-       // we have found the point of
+       // Then loop over all cells and all their vertices, and check whether a
+       // vertex matches the evaluation point. If this is the case, then
+       // extract the point value, set a flag that we have found the point of
        // interest, and exit the loop.
        typename DoFHandler<dim>::active_cell_iterator
        cell = dof_handler.begin_active(),
index c1d4613872548835754f131d47447347976ff445,6242b228e3d37e5ac894894ec8f24c0c8a53be29..3dc49df6928786b9d082f1d50fac157ffff7fe76
@@@ -216,19 -198,15 +198,15 @@@ namespace Step1
      void
      PointXDerivativeEvaluation<dim>::
      operator () (const DoFHandler<dim> &dof_handler,
 -                 const Vector<double> &solution) const
 +                 const Vector<double>  &solution) const
      {
-       // This time initialize the
-       // return value with something
-       // useful, since we will have to
-       // add up a number of
-       // contributions and take the
-       // mean value afterwards...
+       // This time initialize the return value with something useful, since we
+       // will have to add up a number of contributions and take the mean value
+       // afterwards...
        double point_derivative = 0;
  
-       // ...then have some objects of
-       // which the meaning wil become
-       // clear below...
+       // ...then have some objects of which the meaning wil become clear
+       // below...
        QTrapez<dim>  vertex_quadrature;
        FEValues<dim> fe_values (dof_handler.get_fe(),
                                 vertex_quadrature,
Simple merge
index fd816b8be1a573ee5691e16e898869f6d68dcd48,49bb88bd09c133a432b1072ce7ba4e9f63bc010f..c7c9fca01042a6cd323455b85c437aef351d2cd7
@@@ -191,15 -150,11 +150,11 @@@ namespace Step2
    };
  
  
-   // And then we also have to define
-   // these respective functions, of
-   // course. Given our discussion in
-   // the introduction of how the
-   // solution should look like, the
-   // following computations should be
-   // straightforward:
+   // And then we also have to define these respective functions, of
+   // course. Given our discussion in the introduction of how the solution
+   // should look like, the following computations should be straightforward:
    template <int dim>
 -  double RightHandSide<dim>::value (const Point<dim>  & /*p*/,
 +  double RightHandSide<dim>::value (const Point<dim>  &/*p*/,
                                      const unsigned int /*component*/) const
    {
      return 0;
Simple merge
Simple merge
Simple merge
index 0363c2dad69a11b94291487252ab9d42528b81b7,dc4a31b38764ad9ca20814b667f87b37ce17fdbb..049d63568e5c2789781e05f4a73425110413febe
@@@ -1638,21 -1219,13 +1219,13 @@@ namespace Step2
  
      // @sect5{Private member variables}
  
-     // Next, we have a few member
-     // variables. In particular,
-     // these are (i) a reference to
-     // the parameter object (owned by
-     // the main function of this
-     // program, and passed to the
-     // constructor of this class),
-     // (ii) an object describing the
-     // material parameters for the
-     // number of energy groups
-     // requested in the input file,
-     // and (iii) the finite element
-     // to be used by all energy
-     // groups:
+     // Next, we have a few member variables. In particular, these are (i) a
+     // reference to the parameter object (owned by the main function of this
+     // program, and passed to the constructor of this class), (ii) an object
+     // describing the material parameters for the number of energy groups
+     // requested in the input file, and (iii) the finite element to be used by
+     // all energy groups:
 -    const Parameters &parameters;
 +    const Parameters  &parameters;
      const MaterialData material_data;
      FE_Q<dim>          fe;
  
index 4e4bff6dd9bec4f435bf263a287a7aeb097e1062,253a0c92f6e8ad87dc65ab43a59d5c9d076a6995..b97a0558b7bcbf4125278640f1394f0294421846
@@@ -572,16 -400,12 +400,12 @@@ namespace Step2
  
  
  
-   // The constructor takes the
-   // ParameterHandler object and stores
-   // it in a reference. It also
-   // initializes the DoF-Handler and
-   // the finite element system, which
-   // consists of two copies of the
-   // scalar Q1 field, one for $v$ and
-   // one for $w$:
+   // The constructor takes the ParameterHandler object and stores it in a
+   // reference. It also initializes the DoF-Handler and the finite element
+   // system, which consists of two copies of the scalar Q1 field, one for $v$
+   // and one for $w$:
    template <int dim>
 -  UltrasoundProblem<dim>::UltrasoundProblem (ParameterHandler &param)
 +  UltrasoundProblem<dim>::UltrasoundProblem (ParameterHandler  &param)
      :
      prm(param),
      dof_handler(triangulation),
Simple merge
Simple merge
index 1721b8b098f8f5ecd5ed72bf6386dd0166e84096,77199f884b12b95241ce7375272948d24bea0546..6b4b5cd4182c1111ff4aa84c2069ce95e16d8592
@@@ -765,26 -587,19 +587,19 @@@ namespace Step3
    EulerEquations<dim>::Postprocessor::
    compute_derived_quantities_vector (const std::vector<Vector<double> >              &uh,
                                       const std::vector<std::vector<Tensor<1,dim> > > &duh,
 -                                     const std::vector<std::vector<Tensor<2,dim> > > & /*dduh*/,
 -                                     const std::vector<Point<dim> >                  & /*normals*/,
 -                                     const std::vector<Point<dim> >                  & /*evaluation_points*/,
 +                                     const std::vector<std::vector<Tensor<2,dim> > > &/*dduh*/,
 +                                     const std::vector<Point<dim> >                  &/*normals*/,
 +                                     const std::vector<Point<dim> >                  &/*evaluation_points*/,
                                       std::vector<Vector<double> >                    &computed_quantities) const
    {
-     // At the beginning of the function, let us
-     // make sure that all variables have the
-     // correct sizes, so that we can access
-     // individual vector elements without
-     // having to wonder whether we might read
-     // or write invalid elements; we also check
-     // that the <code>duh</code> vector only
-     // contains data if we really need it (the
-     // system knows about this because we say
-     // so in the
-     // <code>get_needed_update_flags()</code>
-     // function below). For the inner vectors,
-     // we check that at least the first element
-     // of the outer vector has the correct
-     // inner size:
+     // At the beginning of the function, let us make sure that all variables
+     // have the correct sizes, so that we can access individual vector
+     // elements without having to wonder whether we might read or write
+     // invalid elements; we also check that the <code>duh</code> vector only
+     // contains data if we really need it (the system knows about this because
+     // we say so in the <code>get_needed_update_flags()</code> function
+     // below). For the inner vectors, we check that at least the first element
+     // of the outer vector has the correct inner size:
      const unsigned int n_quadrature_points = uh.size();
  
      if (do_schlieren_plot == true)
Simple merge
Simple merge
index 0dc140fb0e90777bc7e9fbe48dfc7c77a1b2b0ed,9b613a950578a05e7221f77a166d71d14b621335..ffc9cd87234fe28c9f5956efbb4cbe5e7bcdeae3
  
  // @sect3{Include files}
  //
- // Most of the include files we need for this
- // program have already been discussed in
- // previous programs. In particular, all of
- // the following should already be familiar
- // friends:
+ // Most of the include files we need for this program have already been
+ // discussed in previous programs. In particular, all of the following should
+ // already be familiar friends:
  #include <deal.II/base/quadrature_lib.h>
  #include <deal.II/base/function.h>
 +
 +#include <deal.II/lac/abstract_linear_algebra.h>
 +
 +#define USE_PETSC_LA
 +
 +namespace LA
 +{
 +#ifdef USE_PETSC_LA
 +  using namespace dealii::LinearAlgebraPETSc;
 +#else
 +  using namespace dealii::LinearAlgebraTrilinos;
 +#endif
 +}
 +
 +
  #include <deal.II/lac/vector.h>
  #include <deal.II/lac/full_matrix.h>
  #include <deal.II/lac/solver_cg.h>
@@@ -339,39 -235,28 +250,28 @@@ namespace Step4
                                        locally_relevant_dofs);
      locally_relevant_solution = 0;
      system_rhs.reinit (mpi_communicator,
 -                       dof_handler.n_dofs(),
 -                       dof_handler.n_locally_owned_dofs());
 +                       locally_owned_dofs);
 +
      system_rhs = 0;
  
-     // The next step is to compute hanging node
-     // and boundary value constraints, which we
-     // combine into a single object storing all
+     // The next step is to compute hanging node and boundary value
+     // constraints, which we combine into a single object storing all
      // constraints.
      //
-     // As with all other things in %parallel,
-     // the mantra must be that no processor can
-     // store all information about the entire
-     // universe. As a consequence, we need to
-     // tell the constraints object for which
-     // degrees of freedom it can store
-     // constraints and for which it may not
-     // expect any information to store. In our
-     // case, as explained in the @ref
-     // distributed module, the degrees of
-     // freedom we need to care about on each
-     // processor are the locally relevant ones,
-     // so we pass this to the
-     // ConstraintMatrix::reinit function. As a
-     // side note, if you forget to pass this
-     // argument, the ConstraintMatrix class
-     // will allocate an array with length equal
-     // to the largest DoF index it has seen so
-     // far. For processors with high MPI
-     // process number, this may be very large
-     // -- maybe on the order of billions. The
-     // program would then allocate more memory
-     // than for likely all other operations
-     // combined for this single array.
+     // As with all other things in %parallel, the mantra must be that no
+     // processor can store all information about the entire universe. As a
+     // consequence, we need to tell the constraints object for which degrees
+     // of freedom it can store constraints and for which it may not expect any
+     // information to store. In our case, as explained in the @ref distributed
+     // module, the degrees of freedom we need to care about on each processor
+     // are the locally relevant ones, so we pass this to the
+     // ConstraintMatrix::reinit function. As a side note, if you forget to
+     // pass this argument, the ConstraintMatrix class will allocate an array
+     // with length equal to the largest DoF index it has seen so far. For
+     // processors with high MPI process number, this may be very large --
+     // maybe on the order of billions. The program would then allocate more
+     // memory than for likely all other operations combined for this single
+     // array.
      constraints.clear ();
      constraints.reinit (locally_relevant_dofs);
      DoFTools::make_hanging_node_constraints (dof_handler, constraints);
  
      SolverControl solver_control (dof_handler.n_dofs(), 1e-12);
  
 +#ifdef USE_PETSC_LA
 +
      PETScWrappers::SolverCG solver(solver_control, mpi_communicator);
  
-     // Ask for a symmetric preconditioner by
-     // setting the first parameter in
+     // Ask for a symmetric preconditioner by setting the first parameter in
      // AdditionalData to true.
 -    PETScWrappers::PreconditionBoomerAMG
 +    LA::MPI::PreconditionAMG
      preconditioner(system_matrix,
                     PETScWrappers::PreconditionBoomerAMG::AdditionalData(true));
  
Simple merge
Simple merge
Simple merge
Simple merge
index 0947c81dfa7f7be1238863041a366c41d0f50514,b62734b2123562e50f8085507466a1ee37c3b755..3708b7c48248c1dc18e45cc9979755622d18655f
@@@ -622,25 -426,19 +440,25 @@@ void Step6<dim>::assemble_system (
                              fe_values.JxW(q_point));
            }
  
+       // Finally, transfer the contributions from @p cell_matrix and
+       // @p cell_rhs into the global objects.
        cell->get_dof_indices (local_dof_indices);
-       // transfer the contributions from @p cell_matrix and @cell_rhs into the global objects.
-       constraints.distribute_local_to_global(cell_matrix, cell_rhs, local_dof_indices, system_matrix, system_rhs);
+       constraints.distribute_local_to_global(cell_matrix,
+                                            cell_rhs,
+                                            local_dof_indices,
+                                            system_matrix, system_rhs);
      }
- // Now we are done assembling the linear
-   // system.  The constrained nodes are still
-   // in the linear system (there is a one on
-   // the diagonal of the matrix and all other
-   // entries for this line are set to zero)
-   // but the computed values are invalid. We
-   // compute the correct values for these
-   // nodes at the end of the
-   // <code>solve</code> function.
 +
 +#ifdef USE_PETSC_LA
 +  system_matrix.compress();
 +  system_rhs.compress();
 +#endif
 +
+   // Now we are done assembling the linear system.  The constrained nodes are
+   // still in the linear system (there is a one on the diagonal of the matrix
+   // and all other entries for this line are set to zero) but the computed
+   // values are invalid. We compute the correct values for these nodes at the
+   // end of the <code>solve</code> function.
  }
  
  
@@@ -1076,97 -709,55 +729,58 @@@ void Step6<dim>::run (
  
  // @sect3{The <code>main</code> function}
  
- // The main function is unaltered in
- // its functionality from the
- // previous example, but we have
- // taken a step of additional
- // caution. Sometimes, something goes
- // wrong (such as insufficient disk
- // space upon writing an output file,
- // not enough memory when trying to
- // allocate a vector or a matrix, or
- // if we can't read from or write to
- // a file for whatever reason), and
- // in these cases the library will
- // throw exceptions. Since these are
- // run-time problems, not programming
- // errors that can be fixed once and
- // for all, this kind of exceptions
- // is not switched off in optimized
- // mode, in contrast to the
- // <code>Assert</code> macro which we have
- // used to test against programming
- // errors. If uncaught, these
- // exceptions propagate the call tree
- // up to the <code>main</code> function, and
- // if they are not caught there
- // either, the program is aborted. In
- // many cases, like if there is not
- // enough memory or disk space, we
- // can't do anything but we can at
- // least print some text trying to
- // explain the reason why the program
- // failed. A way to do so is shown in
- // the following. It is certainly
- // useful to write any larger program
- // in this way, and you can do so by
- // more or less copying this function
- // except for the <code>try</code> block that
- // actually encodes the functionality
- // particular to the present
- // application.
+ // The main function is unaltered in its functionality from the previous
+ // example, but we have taken a step of additional caution. Sometimes,
+ // something goes wrong (such as insufficient disk space upon writing an
+ // output file, not enough memory when trying to allocate a vector or a
+ // matrix, or if we can't read from or write to a file for whatever reason),
+ // and in these cases the library will throw exceptions. Since these are
+ // run-time problems, not programming errors that can be fixed once and for
+ // all, this kind of exceptions is not switched off in optimized mode, in
+ // contrast to the <code>Assert</code> macro which we have used to test
+ // against programming errors. If uncaught, these exceptions propagate the
+ // call tree up to the <code>main</code> function, and if they are not caught
+ // there either, the program is aborted. In many cases, like if there is not
+ // enough memory or disk space, we can't do anything but we can at least print
+ // some text trying to explain the reason why the program failed. A way to do
+ // so is shown in the following. It is certainly useful to write any larger
+ // program in this way, and you can do so by more or less copying this
+ // function except for the <code>try</code> block that actually encodes the
+ // functionality particular to the present application.
 -int main ()
 +int main (int argc, char *argv[])
  {
  
-   // The general idea behind the
-   // layout of this function is as
-   // follows: let's try to run the
-   // program as we did before...
+   // The general idea behind the layout of this function is as follows: let's
+   // try to run the program as we did before...
    try
      {
        deallog.depth_console (0);
 +      Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv);
  
 -      Step6<2> laplace_problem_2d;
 -      laplace_problem_2d.run ();
 +      {
 +        Step6<2> laplace_problem_2d;
 +        laplace_problem_2d.run ();
 +      }
      }
-   // ...and if this should fail, try
-   // to gather as much information as
-   // possible. Specifically, if the
-   // exception that was thrown is an
-   // object of a class that is
-   // derived from the C++ standard
-   // class <code>exception</code>, then we can
-   // use the <code>what</code> member function
-   // to get a string which describes
-   // the reason why the exception was
+   // ...and if this should fail, try to gather as much information as
+   // possible. Specifically, if the exception that was thrown is an object of
+   // a class that is derived from the C++ standard class
+   // <code>exception</code>, then we can use the <code>what</code> member
+   // function to get a string which describes the reason why the exception was
    // thrown.
    //
-   // The deal.II exception classes
-   // are all derived from the
-   // standard class, and in
-   // particular, the <code>exc.what()</code>
-   // function will return
-   // approximately the same string as
-   // would be generated if the
-   // exception was thrown using the
-   // <code>Assert</code> macro. You have seen
-   // the output of such an exception
-   // in the previous example, and you
-   // then know that it contains the
-   // file and line number of where
-   // the exception occured, and some
-   // other information. This is also
-   // what the following statements
-   // would print.
+   // The deal.II exception classes are all derived from the standard class,
+   // and in particular, the <code>exc.what()</code> function will return
+   // approximately the same string as would be generated if the exception was
+   // thrown using the <code>Assert</code> macro. You have seen the output of
+   // such an exception in the previous example, and you then know that it
+   // contains the file and line number of where the exception occured, and
+   // some other information. This is also what the following statements would
+   // print.
    //
-   // Apart from this, there isn't
-   // much that we can do except
-   // exiting the program with an
-   // error code (this is what the
-   // <code>return 1;</code> does):
+   // Apart from this, there isn't much that we can do except exiting the
+   // program with an error code (this is what the <code>return 1;</code>
+   // does):
    catch (std::exception &exc)
      {
        std::cerr << std::endl << std::endl
index 3ae105a770b68b72ee3efb9b99f3125c382840d8,34c1acbb41aeef387753bb3180df14ba4011ceda..25b77c3c2b1b7d05c15205df46e4f830dcd65080
@@@ -1459,23 -976,15 +976,15 @@@ namespace Step
    template <int dim>
    void
    GradientEstimation::estimate (const DoFHandler<dim> &dof_handler,
 -                                const Vector<double> &solution,
 +                                const Vector<double>  &solution,
                                  Vector<float>         &error_per_cell)
    {
-     // Before starting with the work,
-     // we check that the vector into
-     // which the results are written,
-     // has the right size. It is a
-     // common error that such
-     // parameters have the wrong size,
-     // but the resulting damage by not
-     // catching these errors are very
-     // subtle as they are usually
-     // corruption of data somewhere in
-     // memory. Often, the problems
-     // emerging from this are not
-     // reproducible, and we found that
-     // it is well worth the effort to
+     // Before starting with the work, we check that the vector into which the
+     // results are written, has the right size. It is a common error that such
+     // parameters have the wrong size, but the resulting damage by not
+     // catching these errors are very subtle as they are usually corruption of
+     // data somewhere in memory. Often, the problems emerging from this are
+     // not reproducible, and we found that it is well worth the effort to
      // check for such things.
      Assert (error_per_cell.size() == dof_handler.get_tria().n_active_cells(),
              ExcInvalidVectorLength (error_per_cell.size(),
index ed9b1615e1f708b6f28ff769ec19fb4c4843bc24,822a0f7fe757eb5134fd9d516acec85f53675905..d74bbb592bfd9e99628227acfefbf441d0ed33bf
@@@ -155,12 -134,11 +134,11 @@@ protected
    /**
     * Name of the exception and call sequence.
     */
 -  const char *exc;
 +  const char  *exc;
  
    /**
-    * A backtrace to the position
-    * where the problem happened, if
-    * the system supports this.
+    * A backtrace to the position where the problem happened, if the
+    * system supports this.
     */
    char **stacktrace;
  
Simple merge
Simple merge
Simple merge

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.