From e7ebaed53a81460170040065774a94e3972369c1 Mon Sep 17 00:00:00 2001 From: wolf Date: Tue, 27 Feb 2001 08:32:31 +0000 Subject: [PATCH] Take over patches dof_tools.h:1.41->1.42, dof_tools.cc:1.53->1.54. (Was: Fix quadratic memory behaviour of compute_intergrid_constraints to linear.) git-svn-id: https://svn.dealii.org/branches/Branch-3-1@4037 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 4 +- deal.II/deal.II/source/dofs/dof_tools.cc | 135 ++++++++++++++++------- 2 files changed, 99 insertions(+), 40 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 83cf803b55..7c3f50ff48 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -899,7 +899,7 @@ class DoFTools const InterGridMap &coarse_to_fine_grid_map, const std::vector > ¶meter_dofs, const std::vector &weight_mapping, - FullMatrix &weights); + std::vector > &weights); /** * This is a function that is @@ -918,7 +918,7 @@ class DoFTools const InterGridMap &coarse_to_fine_grid_map, const std::vector > ¶meter_dofs, const std::vector &weight_mapping, - FullMatrix &weights, + std::vector > &weights, const typename DoFHandler::active_cell_iterator &begin, const typename DoFHandler::active_cell_iterator &end); }; diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 7141a1e5f0..e08626f5da 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -1202,13 +1202,31 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa // global (fine grid) parameter dof // indices to the columns // - // note that the `weights' array - // can take up huge amounts of - // memory, and in particular is - // roughly quadratic in the memory - // consumption! - FullMatrix weights (n_coarse_dofs, - n_parameters_on_fine_grid); + // in the original implementation, + // the weights array was actually + // of FullMatrix type. this + // wasted huge amounts of memory, + // but was fast. nontheless, since + // the memory consumption was + // quadratic in the number of + // degrees of freedom, this was not + // very practical, so we now use a + // vector of rows of the matrix, + // and in each row a vector of + // pairs (colnum,value). this seems + // like the best tradeoff between + // memory and speed, as it is now + // linear in memory and still fast + // enough. + // + // to save some memory and since + // the weights are usually + // (negative) powers of 2, we + // choose the value type of the + // matrix to be @p{float} rather + // than @p{double}. + std::vector > weights(n_coarse_dofs); + // this is this mapping. there is one // entry for each dof on the fine grid; // if it is a parameter dof, then its @@ -1283,11 +1301,12 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa // a dof belongs to, but this at // least tests some errors #ifdef DEBUG - for (unsigned int col=0; col1) && (sum==0)), ExcInternalError()); }; @@ -1313,15 +1332,27 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa // while all others will be // constrained to this dof (and // possibly others) - std::vector representants(weights.m(), -1); - for (unsigned int parameter_dof=0; parameter_dof representants(n_coarse_dofs, -1); + for (unsigned int parameter_dof=0; parameter_dof 0, ExcInternalError()); + + // find the column where the + // representant is mentioned + map::const_iterator i = weights[parameter_dof].begin(); + for (; i!=weights[parameter_dof].end(); ++i) + if (i->second == 1) break; - Assert (column < weights.n(), ExcInternalError()); + Assert (i!=weights[parameter_dof].end(), ExcInternalError()); + const unsigned int column = i->first; // now we know in which column of // weights the representant is, but @@ -1345,8 +1376,7 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa // coarse grid, then the // respective row must be // empty! - for (unsigned int col=0; col &coa // to be unconstrained. otherwise, // all other dofs are constrained { + const unsigned int col = weight_mapping[global_dof]; + Assert (col < n_coarse_dofs, ExcInternalError()); + unsigned int first_used_row=0; - for (; first_used_row::const_iterator col_entry; + for (; first_used_rowsecond == 1) && + (representants[first_used_row] == static_cast(global_dof))) + // dof unconstrained or + // constrained to itself + // (in case this cell is + // mapped to itself, rather + // than to children of + // itself) + continue; + }; - if ((weights(first_used_row,weight_mapping[global_dof]) == 1) && - (representants[first_used_row] == static_cast(global_dof))) - // dof unconstrained or - // constrained to itself - // (in case this cell is - // mapped to itself, rather - // than to children of - // itself) - continue; // otherwise enter all constraints constraints.add_line (global_dof); constraint_line.clear (); - for (unsigned int row=first_used_row; row::const_iterator + j = weights[row].find(col); + if ((j != weights[row].end()) && (j->second != 0)) + constraint_line.push_back (std::make_pair(representants[row], + j->second)); + }; + constraints.add_entries (global_dof, constraint_line); }; }; @@ -1415,7 +1460,7 @@ DoFTools::compute_intergrid_weights (const DoFHandler &coarse_ const InterGridMap &coarse_to_fine_grid_map, const std::vector > ¶meter_dofs, const std::vector &weight_mapping, - FullMatrix &weights) + std::vector > &weights) { // simply distribute the range of // cells to different threads @@ -1448,7 +1493,7 @@ DoFTools::compute_intergrid_weights_1 (const DoFHandler &coars const InterGridMap &coarse_to_fine_grid_map, const std::vector > ¶meter_dofs, const std::vector &weight_mapping, - FullMatrix &weights, + std::vector > &weights, const typename DoFHandler::active_cell_iterator &begin, const typename DoFHandler::active_cell_iterator &end) { @@ -1589,6 +1634,18 @@ DoFTools::compute_intergrid_weights_1 (const DoFHandler &coars // among the @p{i} for each @p{j}. this will // be done by simply taking the first // @p{i} for which @p{w_{ij}==1}. + // + // guard modification of + // the weights array by a + // Mutex. since it should + // happen rather rarely + // that there are several + // threads operating on + // different intergrid + // weights, have only one + // mutex for all of them + static Threads::ThreadMutex mutex; + mutex.acquire (); for (unsigned int i=0; i &coars { const unsigned int wi = parameter_dof_indices[local_dof], wj = weight_mapping[i]; - weights(wi,wj) = global_parameter_representation(i); + weights[wi][wj] = global_parameter_representation(i); }; } else Assert (global_parameter_representation(i) == 0, ExcInternalError()); + + mutex.release (); }; }; }; -- 2.39.5