From: Wolfgang Bangerth Date: Sat, 25 Jan 2014 13:20:09 +0000 (+0000) Subject: More cleanups in SolutionTransfer. X-Git-Tag: v8.2.0-rc1~957 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7219f8d62c06951c31bb21a466f2a5162ecb8fe2;p=dealii.git More cleanups in SolutionTransfer. git-svn-id: https://svn.dealii.org/trunk@32312 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/source/dofs/dof_accessor_set.cc b/deal.II/source/dofs/dof_accessor_set.cc index 77fdbe8faf..b4fef1b51f 100644 --- a/deal.II/source/dofs/dof_accessor_set.cc +++ b/deal.II/source/dofs/dof_accessor_set.cc @@ -70,10 +70,14 @@ set_dof_values_by_interpolation (const Vector &local_values, this->get_fe().get_interpolation_matrix (this->dof_handler->get_fe()[fe_index], interpolation); + // do the interpolation to the target space. for historical + // reasons, matrices are set to size 0x0 internally even + // we reinit as 4x0, so we have to treat this case specially Vector tmp (this->get_fe().dofs_per_cell); - interpolation.vmult (tmp, local_values); + if ((tmp.size() > 0) && (local_values.size() > 0)) + interpolation.vmult (tmp, local_values); - //now set the dof values in the global vector + // now set the dof values in the global vector this->set_dof_values (tmp, values); } } diff --git a/deal.II/source/numerics/solution_transfer.cc b/deal.II/source/numerics/solution_transfer.cc index d15d6768dd..2a6ba4b849 100644 --- a/deal.II/source/numerics/solution_transfer.cc +++ b/deal.II/source/numerics/solution_transfer.cc @@ -436,13 +436,15 @@ interpolate (const std::vector &all_in, const std::vector > *const valuesptr =pointerstruct->second.dof_values_ptr; - // cell stayed or is - // refined + // cell stayed as it was or was refined if (indexptr) { Assert (valuesptr == 0, ExcInternalError()); + const unsigned int old_fe_index = + pointerstruct->second.active_fe_index; + // get the values of // each of the input // data vectors on this @@ -451,135 +453,17 @@ interpolate (const std::vector &all_in, unsigned int in_size = indexptr->size(); for (unsigned int j=0; jactive_fe_index(); - bool different_elements = false; - if (cell->has_children()) - { - new_fe_index = cell->child(0)->active_fe_index(); - for (unsigned int child=1; childn_children(); ++child) - if (cell->child(child)->active_fe_index() != new_fe_index) - { - different_elements = true; - break; - } - } - if (different_elements == true || - new_fe_index != pointerstruct->second.active_fe_index) - { - const unsigned int old_index = - pointerstruct->second.active_fe_index; - tmp.reinit (in_size, true); - for (unsigned int i=0; ihas_children() ? - cell->child(0)->get_fe().dofs_per_cell - : cell->get_fe().dofs_per_cell, true); - // do the interpolation. we get into trouble if the - // interpolation_hp(new,old) matrix hasn't been computed. - // this can happen if the respective elements don't support - // the corresponding interpolation; if that's the case, then - // the computation of the matrix simply sets the matrix - // back to size zero. so if we get here and that is - // the wrong size, then this may be because the elements - // haven't implemented the correct function yet - // - // there is one wrinkle. we would like to only error out if - // the size of the matrix is 0 times 0 but at least one - // of the elements has more than one dof per cell. the - // problem is that if you reinit a matrix to 4x0, it automatically - // sets its size to 0x0. so we can only execute the following - // test if *both* elements have dofs_per_cell>0, not if *at - // least one* have. - Assert (! ((interpolation_hp(new_fe_index,old_index).m() == 0) - && - (interpolation_hp(new_fe_index,old_index).n() == 0) - && - ((dof_handler->get_fe()[new_fe_index].dofs_per_cell > 0) - && - (dof_handler->get_fe()[old_index].dofs_per_cell > 0))), - ExcMessage ("The interpolation between two different " - "elements you are trying to use here has " - "not been implemented for this pair of " - "elements!")); - - // simple case where all children have the - // same FE index: just interpolate to their FE - // first and then use the standard routines - if (tmp.size() > 0) - interpolation_hp(new_fe_index,old_index).vmult (local_values, tmp); - else - local_values = 0; - } - - if (cell->has_children() == false) - cell->set_dof_values (local_values, all_out[j]); - else - for (unsigned int child=0; childn_children(); ++child) - { - if (different_elements == true) - { - const unsigned int c_index = - cell->child(child)->active_fe_index(); - if (c_index != old_index) - { - AssertDimension (tmp.size(), - interpolation_hp(c_index,old_index).n()); - local_values.reinit(cell->child(child)->get_fe().dofs_per_cell, true); - - // do the interpolation. same problem as above - Assert (! ((interpolation_hp(c_index,old_index).m() == 0) - && - (interpolation_hp(c_index,old_index).n() == 0) - && - ((dof_handler->get_fe()[c_index].dofs_per_cell > 0) - && - (dof_handler->get_fe()[old_index].dofs_per_cell > 0))), - ExcMessage ("The interpolation between two different " - "elements you are trying to use here has " - "not been implemented for this pair of " - "elements!")); - - if (tmp.size() > 0) - interpolation_hp(c_index,old_index).vmult (local_values, tmp); - else - local_values = 0; - } - else - local_values = tmp; - } - tmp2.reinit (cell->child(child)->get_fe().dofs_per_cell, true); - cell->child(child)->get_fe().get_prolongation_matrix(child, cell->refinement_case()) - .vmult (tmp2, local_values); - cell->child(child)->set_dof_values(tmp2,all_out[j]); - } - } - else - { - const unsigned int dofs_per_cell - = cell->get_fe().dofs_per_cell; - - AssertDimension (dofs_per_cell, indexptr->size()); - local_values.reinit (dofs_per_cell); - for (unsigned int i=0; iset_dof_values_by_interpolation(local_values, - all_out[j], - cell->active_fe_index()); - } + tmp.reinit (in_size, true); + for (unsigned int i=0; iset_dof_values_by_interpolation (tmp, all_out[j], + old_fe_index); } } - // children of cell were - // deleted else if (valuesptr) + // the children of this cell were + // deleted { Assert (!cell->has_children(), ExcInternalError()); Assert (indexptr == 0,