From 1f4fd3ba4b2ac6f271947b374de3f80a9c4ae77c Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 30 Sep 2009 15:01:58 +0000 Subject: [PATCH] Apply new function in ConstraintMatrix to step-37. git-svn-id: https://svn.dealii.org/trunk@19619 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-37/doc/results.dox | 44 ++++---- deal.II/examples/step-37/step-37.cc | 138 +++++++++++------------ 2 files changed, 87 insertions(+), 95 deletions(-) diff --git a/deal.II/examples/step-37/doc/results.dox b/deal.II/examples/step-37/doc/results.dox index a1ddd0c1c0..a444ddf67a 100644 --- a/deal.II/examples/step-37/doc/results.dox +++ b/deal.II/examples/step-37/doc/results.dox @@ -9,38 +9,38 @@ output: @code Cycle 0 Number of degrees of freedom: 337 -System matrix memory consumption: 0.0224 MBytes. -Multigrid objects memory consumption: 0.04636 MBytes. +System matrix memory consumption: 0.0257 MBytes. +Multigrid objects memory consumption: 0.05071 MBytes. Convergence in 10 CG iterations. Cycle 1 Number of degrees of freedom: 1313 -System matrix memory consumption: 0.08109 MBytes. -Multigrid objects memory consumption: 0.1635 MBytes. +System matrix memory consumption: 0.0925 MBytes. +Multigrid objects memory consumption: 0.1792 MBytes. Convergence in 10 CG iterations. Cycle 2 Number of degrees of freedom: 5185 -System matrix memory consumption: 0.3141 MBytes. -Multigrid objects memory consumption: 0.6208 MBytes. +System matrix memory consumption: 0.3551 MBytes. +Multigrid objects memory consumption: 0.6775 MBytes. Convergence in 10 CG iterations. Cycle 3 Number of degrees of freedom: 20609 -System matrix memory consumption: 1.243 MBytes. -Multigrid objects memory consumption: 2.434 MBytes. +System matrix memory consumption: 1.397 MBytes. +Multigrid objects memory consumption: 2.644 MBytes. Convergence in 10 CG iterations. Cycle 4 Number of degrees of freedom: 82177 -System matrix memory consumption: 4.951 MBytes. -Multigrid objects memory consumption: 9.658 MBytes. +System matrix memory consumption: 5.544 MBytes. +Multigrid objects memory consumption: 10.46 MBytes. Convergence in 10 CG iterations. Cycle 5 Number of degrees of freedom: 328193 -System matrix memory consumption: 19.77 MBytes. -Multigrid objects memory consumption: 38.5 MBytes. +System matrix memory consumption: 22.1 MBytes. +Multigrid objects memory consumption: 41.63 MBytes. Convergence in 10 CG iterations. @endcode @@ -51,32 +51,32 @@ program in three dimensions: @code Cycle 0 Number of degrees of freedom: 517 -System matrix memory consumption: 0.09657 MBytes. -Multigrid objects memory consumption: 0.1413 MBytes. +System matrix memory consumption: 0.1 MBytes. +Multigrid objects memory consumption: 0.1462 MBytes. Convergence in 9 CG iterations. Cycle 1 Number of degrees of freedom: 3817 -System matrix memory consumption: 0.6334 MBytes. -Multigrid objects memory consumption: 0.8567 MBytes. +System matrix memory consumption: 0.6612 MBytes. +Multigrid objects memory consumption: 0.8894 MBytes. Convergence in 10 CG iterations. Cycle 2 Number of degrees of freedom: 29521 -System matrix memory consumption: 4.882 MBytes. -Multigrid objects memory consumption: 6.403 MBytes. +System matrix memory consumption: 5.099 MBytes. +Multigrid objects memory consumption: 6.653 MBytes. Convergence in 10 CG iterations. Cycle 3 Number of degrees of freedom: 232609 -System matrix memory consumption: 38.68 MBytes. -Multigrid objects memory consumption: 50.28 MBytes. +System matrix memory consumption: 40.4 MBytes. +Multigrid objects memory consumption: 52.24 MBytes. Convergence in 11 CG iterations. Cycle 4 Number of degrees of freedom: 1847617 -System matrix memory consumption: 308.4 MBytes. -Multigrid objects memory consumption: 399.6 MBytes. +System matrix memory consumption: 321.9 MBytes. +Multigrid objects memory consumption: 415.1 MBytes. Convergence in 11 CG iterations. @endcode diff --git a/deal.II/examples/step-37/step-37.cc b/deal.II/examples/step-37/step-37.cc index 0388244299..1d3a3edd84 100644 --- a/deal.II/examples/step-37/step-37.cc +++ b/deal.II/examples/step-37/step-37.cc @@ -507,62 +507,69 @@ vmult_on_subrange (const unsigned int first_cell, // OK, now we are sitting in the loop that // goes over our chunks of cells. What we - // need to do is five things: First, we - // have to give the full matrices - // containing the solution at cell dofs and - // quadrature points the correct sizes. We - // use the true argument in - // order to specify that this should be - // done fast, i.e., the field will not be - // initialized since we fill them manually - // in a second anyway. Then, we copy the + // need to do is five things: First, we have + // to give the full matrices containing the + // solution at cell dofs and quadrature + // points the correct sizes. We use the + // true argument in order to + // specify that this should be done fast, + // i.e., the field will not be initialized + // since we fill them manually in the very + // next step second anyway. Then, we copy the // source values from the global vector to // the local cell range, and we perform a // matrix-matrix product to tranform the // values to the quadrature points. It is a // bit tricky to find out how the matrices - // should be multiplied with each - // other. One way to resolve this is to + // should be multiplied with each other, + // i.e., which matrix needs to be + // transposed. One way to resolve this is to // look at the matrix dimensions: // solution_cells has // current_chunk_size rows and // matrix_sizes.m columns, // whereas small_matrix has // matrix_sizes.m rows and - // matrix_sizes.n columns, - // which is also the size of columns in the - // output matrix + // matrix_sizes.n columns, which + // is also the size of columns in the output + // matrix // solution_points. Hence, the - // columns of the first matrix are as many - // as there are rows in the second, which - // means that the product is done - // non-transposed for both matrices. + // columns of the first matrix are as many as + // there are rows in the second, which means + // that the product is done non-transposed + // for both matrices. // // Once the first product is calculated, we // apply the derivative information on all - // the cells and all the quadrature points - // by calling the transform + // the cells and all the quadrature points by + // calling the transform // operation of the // Transformation class, and - // then use a second matrix-matrix product - // to get back to the solution values at - // the support points. This time, we need - // to transpose the small matrix, indicated - // by a mTmult in the - // operations. The fifth and last step is - // to add the local data into the global - // vector, which is what we did in many - // tutorial programs when assembling right - // hand sides. Just use the - // indices_local_to_global - // field to find out how local dofs and - // global dofs are related to each other. + // then use a second matrix-matrix product to + // get back to the solution values at the + // support points. This time, we need to + // transpose the small matrix, indicated by a + // mTmult in the operations. The + // fifth and last step is to add the local + // data into the global vector, which is what + // we did in many tutorial programs when + // assembling right hand sides. We use the + // indices_local_to_global field + // to find out how local dofs and global dofs + // are related to each other. Since we + // simultaneously apply the constraints, we + // hand this task off to the ConstraintMatrix + // object. Most often, itis used to work on + // one cell at a time, but since we work on a + // whole chunk of dofs, we can do that just + // as easily for all the cells at once. solution_cells.reinit (current_chunk_size,matrix_sizes.m, true); solution_points.reinit (current_chunk_size,matrix_sizes.n, true); - for (unsigned int i=0; i::Tvmult_add (Vector &dst, - // The vmult_add function that - // multiplies the matrix with vector - // src and adds the result to - // vector dst first creates a - // copy of the source vector in order to - // apply the constraints. The reason for - // doing this is that constrained dofs are - // zero when used in a solver like CG - // (since they are not real degrees of - // freedom), but the solution at the - // respective nodes might still have - // non-zero values which is necessary to - // represent the field correctly in terms - // of the FE basis functions. Then, we call + // This is the vmult_add + // function that multiplies the matrix with + // vector src and adds the + // result to vector dst. We call // a %parallel function that applies the // multiplication on a subrange of cells - // (cf. the @ref threads module), and we - // eventually condense the constraints on - // the resulting vector. + // (cf. the @ref threads module). // // TODO: Use WorkStream for parallelization // instead of apply_to_subranges, once we @@ -657,33 +653,29 @@ void MatrixFree::vmult_add (Vector &dst, const Vector &src) const { - Vector src_copy (src); - constraints.distribute(src_copy); - parallel::apply_to_subranges (0, matrix_sizes.n_cells, std_cxx1x::bind(&MatrixFree:: template vmult_on_subrange, this, _1,_2, boost::ref(dst), - boost::cref(src_copy)), + boost::cref(src)), 200); - constraints.condense (dst); // One thing to be cautious about: The - // deal.II classes expect that the - // matrix still contains a diagonal - // entry for constrained dofs - // (otherwise, the matrix would be - // singular, which is not what we + // deal.II classes expect that the matrix + // still contains a diagonal entry for + // constrained dofs (otherwise, the matrix + // would be singular, which is not what we // want). Since the - // condense command of the - // constraint matrix sets those - // constrained elements to zero, we - // have to circumvent that problem by - // setting the diagonal to some - // non-zero value. We simply set it to - // one. + // distribute_local_to_global + // command of the constraint matrix which we + // used for add the local elements into the + // global matrix does not do write anything + // with constrained elements, we have to + // circumvent that problem by setting the + // diagonal to some non-zero value. We simply + // set it to one. for (unsigned int i=0; i