]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Some clean ups, some structural modifications.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 10 Mar 1999 19:32:25 +0000 (19:32 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 10 Mar 1999 19:32:25 +0000 (19:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@991 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/solution_transfer.h
deal.II/deal.II/source/numerics/solution_transfer.cc

index cdc974f23e88febc43ae8532d0f0393ee32fc6de..fd6b57735f76813dbddea3bcdc8a47a0649b7bc7 100644 (file)
@@ -20,7 +20,7 @@ template <int dim> class DoFHandler;
  *
  * As the interpolation while
  * coarsening is much more complicated to organize 
- * (see further documentation below) than interpolation while pure refinement
+ * (see further documentation below) than interpolation while pure refinement,
  * #SolutionTransfer# offers two possible usages.
  * \begin{itemize}
  * \item If the grid will only be purely refined
@@ -31,6 +31,7 @@ template <int dim> class DoFHandler;
  *                                     // some refinement e.g.
  * tria->refine_and_coarsen_fixed_fraction(error_indicator, 0.3, 0);
  * tria->execute_coarsening_and_refinement();
+ * dof_handler->distribute_dofs (fe);
  * soltrans.refine_interpolate(solution);
  *                                     // if necessary interpolate some
  *                                     // more functions
@@ -47,23 +48,28 @@ template <int dim> class DoFHandler;
  * tria->prepare_coarsening_and_refinement();
  * soltrans.prepare_for_coarsening_and_refinement(solution);
  * tria->execute_coarsening_and_refinement ();
- * soltrans.refine_interpolate(solution);
+ * dof_handler->distribute_dofs (fe);
+ * soltrans.interpolate(solution);
  * \end{verbatim}
  * Multiple calling of this function #void interpolate (Vector<number> &out) const#
- * is NOT allowed. Interpolating
- * several functions can be performed
- * in one step by using #void 
- * interpolate (vector<Vector<number>
- * >&all_out) const#, see there.
+ * is NOT allowed. Interpolating several functions can be performed in one step
+ * by using #void interpolate (vector<Vector<number> >&all_out) const#, and
+ * using the respective #prepare_for_coarsening_and_refinement# function taking
+ * several vectors as input before actually refining and coarsening the
+ * triangulation (see there).
  * \end{itemize}
- * For deleting all stored data in and reinitializing the #SolutionTransfer#
- * use the #reinit()# function.
  *
- * Note that the #user_pointer# of the cell are used. Be sure that you don't need
+ * For deleting all stored data in and reinitializing the
+ * #SolutionTransfer# use the #clear()# function.
+ *
+ * Note that the #user_pointer# of some cells are used. Be sure that you don't need
  * them otherwise.
  *
- * @author Ralf Hartmann, 1999
- */
+ * The template argument #number# denotes the data type of the vectors you want
+ * to transfer.
+ *
+ * @author Ralf Hartmann, 1999 */
 template<int dim, typename number>
 class SolutionTransfer
 {
@@ -84,7 +90,7 @@ class SolutionTransfer
                                      * it has
                                      * directly after calling the Constructor
                                      */
-    void reinit();
+    void clear();
 
                                     /**
                                      * Prepares the #SolutionTransfer# for
@@ -111,17 +117,17 @@ class SolutionTransfer
     
                                     /**
                                      * Same as previous function
-                                     * but only for one discrete function
+                                     * but for only one discrete function
                                      * to interpolate.
                                      */
     void prepare_for_coarsening_and_refinement (const Vector<number> &in);
                      
                                     /**
                                      * This function
-                                     * interpolates the discrete function #in#
-                                     * on the grid before the
-                                     * refinement to the function #out#
-                                     * on the refined grid.
+                                     * interpolates the discrete function #in#,
+                                     * which is a vector on the grid before the
+                                     * refinement, to the function #out#
+                                     * which then is a vector on the refined grid.
                                      * It assumes the vectors having the
                                      * right sizes (i.e. in.size()==n_dofs_old,
                                      * out.size()==n_dofs_refined)
@@ -151,14 +157,13 @@ class SolutionTransfer
                                      * interpolates the discrete functions
                                      * that are stored in #all_out# onto
                                      * the refined and/or coarsenend grid.
-                                     * It assumes the vectors in #all_out#
-                                     * being idendically the same vectors
+                                     * It assumes the vectors in #all_in#
+                                     * denote the same vectors
                                      * as in #all_in# as parameter
                                      * of #prepare_for_refinement_and_coarsening
                                      * (vector<Vector<number> > &all_in)#.
-                                     * So they have the 
-                                     * right sizes (i.e. all_out[j].size()==
-                                     * n_dofs_old  for all j).
+                                     * However, there is no way of verifying
+                                     * this internally, so be careful here.
                                      *
                                      * Calling this function is allowed only
                                      * if first #Triangulation::prepare_coarsening_
@@ -172,14 +177,12 @@ class SolutionTransfer
                                      * several functions can be performed
                                      * in one step.
                                      */
-    void interpolate (vector<Vector<number> >&all_out) const;
+    void interpolate (const vector<Vector<number> >&all_in,
+                     vector<Vector<number> >      &all_out) const;
       
                                     /**
                                      * Same as the previous function.
                                      * It interpolates only one function.
-                                     * It assumes the vector having the
-                                     * right size (i.e. all_out[j].size()==
-                                     * n_dofs_old  for all j)
                                      * 
                                      * Multiple calling of this function is
                                      * NOT allowed. Interpolating
@@ -188,7 +191,8 @@ class SolutionTransfer
                                      * interpolate (vector<Vector<number>
                                      * >&all_out) const#
                                      */
-    void interpolate (Vector<number> &out) const;
+    void interpolate (const Vector<number> &in,
+                     Vector<number>       &out) const;
 
                                     /**
                                      * Exception
@@ -232,7 +236,15 @@ class SolutionTransfer
                                      * Exception
                                      */
     DeclException0(ExcInternalError);
-    
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException2 (ExcInvalidVectorSize,
+                   int, int,
+                   << "The data vector has the size " << arg1
+                   << ", but " << arg2 << " was expected.");
+                                 
   private:
 
                                     /**
@@ -318,11 +330,11 @@ class SolutionTransfer
                                      * After calling 
                                      * #prepare_for_refinement_and_coarsening
                                      * (vector<Vector<number> > &all_in)#
-                                     * this this pointer points to the vector
+                                     * this pointer points to the vector
                                      * #all_in# for later comparison with
                                      * the vector #all_out#
                                      */
-    const vector<Vector<number> > * vecs_ptr;
+//    const vector<Vector<number> > * vecs_ptr;
 };
 
 
index 51d66e93b8815de0a0dde88ee21df324722dbc19..e071f81e2ea38ab77aff7740135d8a1498f912a8 100644 (file)
 
 
 template<int dim, typename number>
-SolutionTransfer<dim, number>::SolutionTransfer(
-  const DoFHandler<dim> &dof):
+SolutionTransfer<dim, number>::SolutionTransfer(const DoFHandler<dim> &dof):
                dof_handler(&dof),
                n_dofs_old(0),
                prepared_for_pure_refinement(0),
-               prepared_for_coarsening_and_refinement(0),
-               vecs_ptr(0) {}
+               prepared_for_coarsening_and_refinement(0)
+{}
+
 
 
 template<int dim, typename number>
 SolutionTransfer<dim, number>::~SolutionTransfer()
 {
-  reinit();
+  clear ();
 }
 
 
+
 template<int dim, typename number>
-void SolutionTransfer<dim, number>::reinit()
+void SolutionTransfer<dim, number>::clear ()
 {
   if (indices_on_cell.size())
     indices_on_cell.erase(indices_on_cell.begin(), indices_on_cell.end());  
@@ -39,7 +40,7 @@ void SolutionTransfer<dim, number>::reinit()
     all_pointerstructs.erase(all_pointerstructs.begin(), all_pointerstructs.end());
   if (dof_values_on_cell.size())
     dof_values_on_cell.erase(dof_values_on_cell.begin(), dof_values_on_cell.end());
-  vecs_ptr=0;
+
   prepared_for_pure_refinement=false;
   prepared_for_coarsening_and_refinement=false;
 }
@@ -53,14 +54,14 @@ void SolutionTransfer<dim, number>::prepare_for_pure_refinement()
   Assert(!prepared_for_coarsening_and_refinement, 
         ExcAlreadyPrepForRefAndCoarse());
 
-  reinit();
+  clear();
 
-  unsigned int n_cells=dof_handler->get_tria().n_active_cells();
-  unsigned int total_dofs=dof_handler->get_fe().total_dofs;
+  const unsigned int n_active_cells = dof_handler->get_tria().n_active_cells();
+  const unsigned int total_dofs     = dof_handler->get_fe().total_dofs;
   n_dofs_old=dof_handler->n_dofs();
 
-  indices_on_cell=vector<vector<int> > (
-    n_cells, vector<int> (total_dofs));
+  indices_on_cell=vector<vector<int> > (n_active_cells,
+                                       vector<int> (total_dofs));
 
   DoFHandler<dim>::cell_iterator cell = dof_handler->begin(),
                                 endc = dof_handler->end();
@@ -69,8 +70,16 @@ void SolutionTransfer<dim, number>::prepare_for_pure_refinement()
     {
       if (cell->active())
        {
-         cell->set_user_pointer(&indices_on_cell[i]);
+                                          // on each cell store the indices
+                                          // of the dofs. after refining we
+                                          // get the values on the children
+                                          // by taking these indices, getting
+                                          // the respective values out of
+                                          // the data vectors and prolonging
+                                          // them to the children
          cell->get_dof_indices(indices_on_cell[i]);
+         cell->set_user_pointer(&indices_on_cell[i]);
+
          ++i;    
        }
       else
@@ -80,9 +89,11 @@ void SolutionTransfer<dim, number>::prepare_for_pure_refinement()
 }
 
 
+
 template<int dim, typename number>  
-void SolutionTransfer<dim, number>::refine_interpolate(
-  const Vector<number> &in, Vector<number> &out) const
+void
+SolutionTransfer<dim, number>::refine_interpolate(const Vector<number> &in,
+                                                 Vector<number>       &out) const
 {
   Assert(prepared_for_pure_refinement, ExcNotPrepared());
   Assert(in.size()==n_dofs_old, ExcWrongVectorSize(in.size(),n_dofs_old));
@@ -100,6 +111,12 @@ void SolutionTransfer<dim, number>::refine_interpolate(
   for (; cell!=endc; ++cell) 
     {
       if (cell->user_pointer())
+                                        // this cell was refined or not
+                                        // touched at all, so we can get
+                                        // the new values by just setting
+                                        // or interpolating to the children,
+                                        // which is both done by one
+                                        // function
        {
          indexptr=static_cast<vector<int> *>(cell->user_pointer());
          for (unsigned int i=0; i<total_dofs; ++i)
@@ -109,9 +126,10 @@ void SolutionTransfer<dim, number>::refine_interpolate(
     }
 }
 
+
+
 template<int dim, typename number>  
-void SolutionTransfer<dim, number>::refine_interpolate(
-  Vector<number> &vec) const
+void SolutionTransfer<dim, number>::refine_interpolate (Vector<number> &vec) const
 {
   Assert(vec.size()==n_dofs_old, ExcWrongVectorSize(vec.size(),n_dofs_old));
 
@@ -124,19 +142,21 @@ void SolutionTransfer<dim, number>::refine_interpolate(
 
 
 template<int dim, typename number>
-void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
-  const vector<Vector<number> > &all_in)
+void
+SolutionTransfer<dim, number>::
+prepare_for_coarsening_and_refinement(const vector<Vector<number> > &all_in)
 {
   Assert(!prepared_for_pure_refinement, ExcAlreadyPrepForRef());
   Assert(!prepared_for_coarsening_and_refinement, 
         ExcAlreadyPrepForRefAndCoarse());
-  unsigned int in_size=all_in.size();
+  
+  const unsigned int in_size=all_in.size();
   Assert(in_size!=0, ExcNoInVectorsGiven());
 
-  reinit();
+  clear();
 
-  unsigned int n_cells=dof_handler->get_tria().n_active_cells();
-  unsigned int total_dofs=dof_handler->get_fe().total_dofs;
+  const unsigned int n_active_cells = dof_handler->get_tria().n_active_cells();
+  const unsigned int total_dofs     = dof_handler->get_fe().total_dofs;
   n_dofs_old=dof_handler->n_dofs();
 
   for (unsigned int i=0; i<in_size; ++i)
@@ -144,7 +164,7 @@ void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
       Assert(all_in[i].size()==n_dofs_old,
             ExcWrongVectorSize(all_in[i].size(),n_dofs_old));
     }
-  vecs_ptr=&all_in;
+
                                   // first count the number
                                   // of cells that'll be coarsened
                                   // and that'll stay or be refined
@@ -159,25 +179,26 @@ void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
       else
        ++n_cells_to_stay_or_refine;
     }
-  Assert((n_cells_to_coarsen+n_cells_to_stay_or_refine)==n_cells,
+  Assert((n_cells_to_coarsen+n_cells_to_stay_or_refine)==n_active_cells,
         ExcInternalError());
   Assert(n_cells_to_coarsen%GeometryInfo<dim>::children_per_cell==0,
         ExcInternalError());
   Assert(n_cells_to_coarsen%GeometryInfo<dim>::children_per_cell==0,
         ExcTriaPrepCoarseningNotCalledBefore());
   
-  unsigned int n_coarsen_fathers=n_cells_to_coarsen/
-                                GeometryInfo<dim>::children_per_cell;
+  const unsigned int n_coarsen_fathers = n_cells_to_coarsen /
+                                        GeometryInfo<dim>::children_per_cell;
 
                                   // allocate the needed memory
-  indices_on_cell=vector<vector<int> > (
-    n_cells_to_stay_or_refine, vector<int> (total_dofs));
-  dof_values_on_cell=vector<vector<Vector<number> > > (
-       n_coarsen_fathers, vector<Vector<number> > (
-        in_size, Vector<number> (total_dofs)));
+  indices_on_cell    = vector<vector<int> > (n_cells_to_stay_or_refine,
+                                           vector<int> (total_dofs));
+  dof_values_on_cell
+    = vector<vector<Vector<number> > > (n_coarsen_fathers,
+                                       vector<Vector<number> > (in_size,
+                                                                Vector<number> (total_dofs)));
 
-  all_pointerstructs=vector<Pointerstruct> (
-    n_cells_to_stay_or_refine+n_coarsen_fathers);
+  all_pointerstructs = vector<Pointerstruct> (n_cells_to_stay_or_refine +
+                                             n_coarsen_fathers);
 
 
                                   // we need counters for
@@ -191,29 +212,43 @@ void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
     {
       if (cell->active() && !cell->coarsen_flag_set())
        {
+                                          // cell will not be coarsened,
+                                          // so we get away by storing the
+                                          // dof indices and later
+                                          // interpolating to the children
          cell->get_dof_indices(indices_on_cell[n_sr]);
          all_pointerstructs[n].indices_ptr=&indices_on_cell[n_sr];
          all_pointerstructs[n].dof_values_ptr=0;
          cell->set_user_pointer(&all_pointerstructs[n]);
-         ++n_sr, ++n;
+         ++n_sr;
+         ++n;
        }
       else if (cell->has_children() && cell->child(0)->coarsen_flag_set())
        {
+                                          // note that if one child has the
+                                          // coaresen flag, then all should
+                                          // have if Tria::prepare_* has
+                                          // worked correctly
          for (unsigned int i=1; i<GeometryInfo<dim>::children_per_cell; ++i)
            Assert(cell->child(i)->coarsen_flag_set(),
                   ExcTriaPrepCoarseningNotCalledBefore());
              
          for (unsigned int j=0; j<in_size; ++j)
            {
-             cell->get_interpolated_dof_values(
-               all_in[j], dof_values_on_cell[n_cf][j]);
+                                              // store the data of each of
+                                              // the input vectors
+             cell->get_interpolated_dof_values(all_in[j],
+                                               dof_values_on_cell[n_cf][j]);
            }
          all_pointerstructs[n].indices_ptr=0;
          all_pointerstructs[n].dof_values_ptr=&dof_values_on_cell[n_cf];
          cell->set_user_pointer(&all_pointerstructs[n]);    
-         ++n_cf, ++n;
+         ++n_cf;
+         ++n;
        }
       else
+                                        // some cell on the lower levels to
+                                        // which nothing will happen
        cell->clear_user_pointer();
     }
   Assert(n_sr==n_cells_to_stay_or_refine, ExcInternalError());
@@ -222,9 +257,11 @@ void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
   prepared_for_coarsening_and_refinement=true;
 }
 
+
+
 template<int dim, typename number>
-void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
-  const Vector<number> &in)
+void
+SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(const Vector<number> &in)
 {
   vector<Vector<number> > all_in=vector<Vector<number> >(1, in);
   prepare_for_coarsening_and_refinement(all_in);
@@ -233,25 +270,36 @@ void SolutionTransfer<dim, number>::prepare_for_coarsening_and_refinement(
 
 
 template<int dim, typename number>
-void SolutionTransfer<dim, number>::interpolate (
-  vector<Vector<number> > &all_out) const
+void SolutionTransfer<dim, number>::
+interpolate (const vector<Vector<number> > &all_in,
+            vector<Vector<number> >       &all_out) const
 {
   Assert(prepared_for_coarsening_and_refinement, ExcNotPrepared());
-  Assert(vecs_ptr==&all_out, ExcVectorsDifferFromInVectors());
+  for (unsigned int i=0; i<all_in.size(); ++i)
+    Assert (all_in[i].size() == n_dofs_old,
+           ExcInvalidVectorSize(all_in[i].size(), n_dofs_old));
+                             
   unsigned int out_size=all_out.size();
-                                  // local copy of the vectors
-  vector<Vector<number> > all_in(all_out);
-                                  // resize the out_all vectors
-  for (unsigned int j=0; j<out_size; ++j)
-    all_out[j].reinit(dof_handler->n_dofs());
 
-  unsigned int total_dofs=dof_handler->get_fe().total_dofs;  
+                                  // resize the output vector
+  if (all_out.size() != all_in.size())
+    all_out = vector<Vector<number> >(all_in.size(),
+                                     Vector<number>(dof_handler->n_dofs()));
+  else
+    {
+      for (unsigned int i=0; i<all_in.size(); ++i)
+       if (all_out[i].size() != all_in[i].size())
+         all_out[i].reinit (all_in[i].size());
+    };
+
+  
+  const unsigned int total_dofs=dof_handler->get_fe().total_dofs;  
   Vector<number> local_values(total_dofs);
 
-  vector<int> *indexptr;
-  Pointerstruct *structptr;
+  vector<int>             *indexptr;
+  Pointerstruct           *structptr;
   vector<Vector<number> > *valuesptr;
-  vector<int> dofs(total_dofs);
+  vector<int>              dofs(total_dofs);
 
   DoFHandler<dim>::cell_iterator cell = dof_handler->begin(),
                                 endc = dof_handler->end();
@@ -264,19 +312,33 @@ void SolutionTransfer<dim, number>::interpolate (
                                           // cell stayed or is refined
          if (indexptr=structptr->indices_ptr)
            {
+             Assert (structptr->dof_values_ptr == 0,
+                     ExcInternalError());
+
+                                              // get the values of each of
+                                              // the input data vectors on
+                                              // this cell and prolong it
+                                              // to its children
              for (unsigned int j=0; j<out_size; ++j)
                {
                  for (unsigned int i=0; i<total_dofs; ++i)
                    local_values(i)=all_in[j](indexptr->operator[](i));
-                 cell->set_dof_values_by_interpolation(
-                   local_values, all_out[j]);
+                 cell->set_dof_values_by_interpolation(local_values,
+                                                       all_out[j]);
                }
            }
                                           // children of cell were deleted
          else if (valuesptr=structptr->dof_values_ptr)
            {
-             Assert(!cell->has_children(), ExcInternalError());
+             Assert (!cell->has_children(), ExcInternalError());
+             Assert (structptr->indices_ptr == 0,
+                     ExcInternalError());
+
+                                              // get the local indices
              cell->get_dof_indices(dofs);
+
+                                              // distribute the stored data
+                                              // to the new vectors
              for (unsigned int j=0; j<out_size; ++j)
                {
                  for (unsigned int i=0; i<total_dofs; ++i)
@@ -291,15 +353,22 @@ void SolutionTransfer<dim, number>::interpolate (
 }
 
 
+
 template<int dim, typename number>
-void SolutionTransfer<dim, number>::interpolate(Vector<number> &out) const
+void SolutionTransfer<dim, number>::interpolate(const Vector<number> &in,
+                                               Vector<number>       &out) const
 {
-  vector<Vector<number> > all_in=vector<Vector<number> >(1, out);
-  interpolate(all_in);
-  out=all_in[0];
+  vector<Vector<number> > all_in=vector<Vector<number> >(1, in);
+  vector<Vector<number> > all_out=vector<Vector<number> >(1, out);
+  interpolate(all_in,
+             all_out);
+  out=all_out[0];
 }
 
 
-template class SolutionTransfer<2, double>;  
+
+
+template class SolutionTransfer<deal_II_dimension, float>;
+template class SolutionTransfer<deal_II_dimension, double>;
 
 /*----------------------------   solutiontransfer.cc     ----------------------*/

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.