From: Timo Heister Date: Tue, 27 Nov 2012 15:51:16 +0000 (+0000) Subject: merge from mainline X-Git-Tag: v8.0.0~120^2~107 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=02bce20a003bfddff26b8a963fd08ca7d6d345b2;p=dealii.git merge from mainline git-svn-id: https://svn.dealii.org/branches/branch_unify_linear_algebra@27697 0785d39b-7218-0410-832d-ea1e28bc413d --- 02bce20a003bfddff26b8a963fd08ca7d6d345b2 diff --cc deal.II/examples/step-13/step-13.cc index 63e62cc64b,b85cb3d50a..b180dc8325 --- a/deal.II/examples/step-13/step-13.cc +++ b/deal.II/examples/step-13/step-13.cc @@@ -292,32 -203,19 +203,19 @@@ namespace Step1 void PointValueEvaluation:: operator () (const DoFHandler &dof_handler, - const Vector &solution) const + const Vector &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::active_cell_iterator cell = dof_handler.begin_active(), diff --cc deal.II/examples/step-14/step-14.cc index c1d4613872,6242b228e3..3dc49df692 --- a/deal.II/examples/step-14/step-14.cc +++ b/deal.II/examples/step-14/step-14.cc @@@ -216,19 -198,15 +198,15 @@@ namespace Step1 void PointXDerivativeEvaluation:: operator () (const DoFHandler &dof_handler, - const Vector &solution) const + const Vector &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 vertex_quadrature; FEValues fe_values (dof_handler.get_fe(), vertex_quadrature, diff --cc deal.II/examples/step-20/step-20.cc index fd816b8be1,49bb88bd09..c7c9fca010 --- a/deal.II/examples/step-20/step-20.cc +++ b/deal.II/examples/step-20/step-20.cc @@@ -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 - double RightHandSide::value (const Point & /*p*/, + double RightHandSide::value (const Point &/*p*/, const unsigned int /*component*/) const { return 0; diff --cc deal.II/examples/step-28/step-28.cc index 0363c2dad6,dc4a31b387..049d63568e --- a/deal.II/examples/step-28/step-28.cc +++ b/deal.II/examples/step-28/step-28.cc @@@ -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 ¶meters; + const Parameters ¶meters; const MaterialData material_data; FE_Q fe; diff --cc deal.II/examples/step-29/step-29.cc index 4e4bff6dd9,253a0c92f6..b97a0558b7 --- a/deal.II/examples/step-29/step-29.cc +++ b/deal.II/examples/step-29/step-29.cc @@@ -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 - UltrasoundProblem::UltrasoundProblem (ParameterHandler ¶m) + UltrasoundProblem::UltrasoundProblem (ParameterHandler ¶m) : prm(param), dof_handler(triangulation), diff --cc deal.II/examples/step-33/step-33.cc index 1721b8b098,77199f884b..6b4b5cd418 --- a/deal.II/examples/step-33/step-33.cc +++ b/deal.II/examples/step-33/step-33.cc @@@ -765,26 -587,19 +587,19 @@@ namespace Step3 EulerEquations::Postprocessor:: compute_derived_quantities_vector (const std::vector > &uh, const std::vector > > &duh, - const std::vector > > & /*dduh*/, - const std::vector > & /*normals*/, - const std::vector > & /*evaluation_points*/, + const std::vector > > &/*dduh*/, + const std::vector > &/*normals*/, + const std::vector > &/*evaluation_points*/, std::vector > &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 duh vector only - // contains data if we really need it (the - // system knows about this because we say - // so in the - // get_needed_update_flags() - // 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 duh vector only + // contains data if we really need it (the system knows about this because + // we say so in the get_needed_update_flags() 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) diff --cc deal.II/examples/step-40/step-40.cc index 0dc140fb0e,9b613a9505..ffc9cd8723 --- a/deal.II/examples/step-40/step-40.cc +++ b/deal.II/examples/step-40/step-40.cc @@@ -13,28 -13,11 +13,26 @@@ // @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 #include + +#include + +#define USE_PETSC_LA + +namespace LA +{ +#ifdef USE_PETSC_LA + using namespace dealii::LinearAlgebraPETSc; +#else + using namespace dealii::LinearAlgebraTrilinos; +#endif +} + + #include #include #include @@@ -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); @@@ -621,14 -436,11 +459,13 @@@ 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)); diff --cc deal.II/examples/step-6/step-6.cc index 0947c81dfa,b62734b212..3708b7c482 --- a/deal.II/examples/step-6/step-6.cc +++ b/deal.II/examples/step-6/step-6.cc @@@ -622,25 -426,19 +440,25 @@@ void Step6::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); } + +#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 - // solve function. + // 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 solve function. } @@@ -1076,97 -709,55 +729,58 @@@ void Step6::run ( // @sect3{The main 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 - // Assert macro which we have - // used to test against programming - // errors. If uncaught, these - // exceptions propagate the call tree - // up to the main 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 try 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 Assert macro which we have used to test + // against programming errors. If uncaught, these exceptions propagate the + // call tree up to the main 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 try 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 exception, then we can - // use the what 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 + // exception, then we can use the what 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 exc.what() - // function will return - // approximately the same string as - // would be generated if the - // exception was thrown using the - // Assert 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 exc.what() function will return + // approximately the same string as would be generated if the exception was + // thrown using the Assert 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 - // return 1; does): + // Apart from this, there isn't much that we can do except exiting the + // program with an error code (this is what the return 1; + // does): catch (std::exception &exc) { std::cerr << std::endl << std::endl diff --cc deal.II/examples/step-9/step-9.cc index 3ae105a770,34c1acbb41..25b77c3c2b --- a/deal.II/examples/step-9/step-9.cc +++ b/deal.II/examples/step-9/step-9.cc @@@ -1459,23 -976,15 +976,15 @@@ namespace Step template void GradientEstimation::estimate (const DoFHandler &dof_handler, - const Vector &solution, + const Vector &solution, Vector &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(), diff --cc deal.II/include/deal.II/base/exceptions.h index ed9b1615e1,822a0f7fe7..d74bbb592b --- a/deal.II/include/deal.II/base/exceptions.h +++ b/deal.II/include/deal.II/base/exceptions.h @@@ -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;