]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Check for right sizes of in and output vectors.
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 24 Oct 2003 09:18:33 +0000 (09:18 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 24 Oct 2003 09:18:33 +0000 (09:18 +0000)
git-svn-id: https://svn.dealii.org/trunk@8153 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 2d43a7c1698750a5b870186b41831432cdb25a3f..a6510d86f8c66ad13351369d1c74eb9ee2c3c9d1 100644 (file)
  * (i.e. not locally coarsened) then use @p{SolutionTransfer} as follows
  * @begin{verbatim}
  * SolutionTransfer<dim, double> soltrans(*dof_handler);
+ *                                     // flag some cells for refinement, e.g.
+ * GridRefinement::refine_and_coarsen_fixed_fraction(
+ *   *tria, error_indicators, 0.3, 0);
+ *                                     // prepare the triangulation
+ *                                     // for refinement,
+ * tria->prepare_coarsening_and_refinement();
+ *                                     // tell the SolutionTransfer object
+ *                                     // that we intent to do pure refinement,
  * soltrans.prepare_for_pure_refinement();
- *                                     // some refinement e.g.
- * tria->refine_and_coarsen_fixed_fraction(error_indicator, 0.3, 0);
+ *                                     // actually execute the refinement,
  * tria->execute_coarsening_and_refinement();
+ *                                     // and redistribute dofs.
  * dof_handler->distribute_dofs (fe);
- * soltrans.refine_interpolate(solution, interpolated_solution);
- *                                     // if necessary interpolate some
- *                                     // more functions
- * soltrans.refine_interpolate(sol2, interpolated_sol2);
+ * @end{verbatim}
+ *
+ * Then there are three different possibilities of how to proceed. Either
+ * @begin{verbatim}
+ *                                     // resize and interpolate a solution
+ *                                     // vector `in-place'
+ * soltrans.refine_interpolate(solution);
+ * @end{verbatim}
+ * or, when the old solution vector is still needed,
+ * @begin{verbatim}
+ *                                     // take a copy of the solution vector
+ * Vector<double> solution_old(solution);
+ *                                     // resize solution vector to the correct
+ *                                     // size, as the @p{refine_interpolate}
+ *                                     // function requires the vectors to be
+ *                                     // of right sizes
+ * solution.reinit(dof_handler->n_dofs());
+ *                                     // and finally interpolate
+ * soltrans.refine_interpolate(solution_old, solution);
+ * @end{verbatim}
+ *
+ * Although the @p{refine_interpolate} functions are allowed to be
+ * called multiple times, e.g. for interpolating several solution
+ * vectors, there is following possibility of interpolating several
+ * functions simultaneously.
+ * @begin{verbatim}
+ * vector<Vector<double> > solutions_old(n_vectors, Vector<double> (n));
  * ...
+ * vector<Vector<double> > solutions(n_vectors, Vector<double> (n));
+ * soltrans.refine_interpolate(solutions_old, solutions);
  * @end{verbatim}
- * @item If the grid will be coarsenend and refined
+ * @item If the grid will be refined AND coarsened
  * then use @p{SolutionTransfer} as follows
  * @begin{verbatim}
  * SolutionTransfer<dim, double> soltrans(*dof_handler);
- *                                     // some refinement e.g.
- * tria->refine_and_coarsen_fixed_fraction(error_indicator, 0.3, 0.05);
- *                                     // very important:
+ *                                     // flag some cells for refinement
+ *                                     // and coarsening, e.g.
+ * GridRefinement::refine_and_coarsen_fixed_fraction(
+ *   *tria, error_indicators, 0.3, 0.05);
+ *                                     // prepare the triangulation,
  * tria->prepare_coarsening_and_refinement();
+ *                                     // prepare the SolutionTransfer object
+ *                                     // for coarsening and refinement and give
+ *                                     // the solution vector that we intent to
+ *                                     // interpolate later,
  * soltrans.prepare_for_coarsening_and_refinement(solution);
+ *                                     // actually execute the refinement,
  * tria->execute_coarsening_and_refinement ();
+ *                                     // redistribute dofs,
  * dof_handler->distribute_dofs (fe);
+ *                                     // and interpolate the solution
+ * Vector<double> interpolate_solution(dof_handler->n_dofs());
  * soltrans.interpolate(solution, interpolated_solution);
  * @end{verbatim}
  *
  * interpolated and set into the vector @p{out} that is at the end the
  * discrete function @p{in} interpolated on the refined mesh.
  *
- * The @p{refine_interpolate(in,out)} function can be called multiplely for
+ * The @p{refine_interpolate(in,out)} function can be called multiply for
  * arbitrary many discrete functions (solution vectors) on the original grid. 
  *
  * @item Solution transfer while coarsening and refinement. After 
@@ -267,11 +310,14 @@ class SolutionTransfer
                                      * several functions can be
                                      * performed in one step.
                                      *
-                                     * If the number of output
-                                     * vectors is different from the
-                                     * number of input vectors, or if
-                                     * their sizes are not correct,
-                                     * then this is corrected.
+                                     * The number of output vectors
+                                     * is assumed to be the same as
+                                     * the number of input
+                                     * vectors. Also, the sizes of
+                                     * the output vectors are assumed
+                                     * to be of the right size
+                                     * (@p{n_dofs_refined}). Otherwise
+                                     * an assertion will be thrown.
                                      */
     void interpolate (const std::vector<Vector<number> >&all_in,
                      std::vector<Vector<number> >      &all_out) const;
@@ -279,12 +325,15 @@ class SolutionTransfer
                                     /**
                                      * Same as the previous function.
                                      * It interpolates only one function.
+                                     * It assumes the vectors having the
+                                     * right sizes (i.e. @p{in.size()==n_dofs_old},
+                                     * @p{out.size()==n_dofs_refined})
                                      * 
                                      * Multiple calling of this function is
                                      * NOT allowed. Interpolating
                                      * several functions can be performed
                                      * in one step by using
-                                     * @p{interpolate (all_out)}
+                                     * @p{interpolate (all_in, all_out)}
                                      */
     void interpolate (const Vector<number> &in,
                      Vector<number>       &out) const;
@@ -332,7 +381,7 @@ class SolutionTransfer
     DeclException2(ExcWrongVectorSize,
                   int, int,
                   << "The size of the vector is " << arg1
-                  << "although it should be " << arg2 << ".");
+                  << " although it should be " << arg2 << ".");
                                  
   private:
 
index 99b15d2f9aa073b0fb27a19fbda1da6eaf9a9105..e5bb7b47140fcba48787b0592e35847329057971 100644 (file)
@@ -282,24 +282,16 @@ interpolate (const std::vector<Vector<number> > &all_in,
             std::vector<Vector<number> >       &all_out) const
 {
   Assert(prepared_for==coarsening_and_refinement, ExcNotPrepared());
-  for (unsigned int i=0; i<all_in.size(); ++i)
+  const unsigned int size=all_in.size();
+  Assert(all_out.size()==size, ExcDimensionMismatch(all_out.size(), size));
+  for (unsigned int i=0; i<size; ++i)
     Assert (all_in[i].size() == n_dofs_old,
            ExcWrongVectorSize(all_in[i].size(), n_dofs_old));
-
-  
-  const unsigned int in_size = all_in.size();
-
-                                  // resize the output vector if
-                                  // necessary
-  if (all_out.size() != in_size)
-    all_out.resize (in_size, Vector<number>(dof_handler->n_dofs()));
-  else
-    for (unsigned int i=0; i<in_size; ++i)
-      if (all_out[i].size() != dof_handler->n_dofs())
-       all_out[i].reinit (dof_handler->n_dofs());
-
-  for (unsigned int i=0; i<in_size; ++i)
-    for (unsigned int j=0; j<in_size; ++j)
+  for (unsigned int i=0; i<all_out.size(); ++i)
+    Assert (all_out[i].size() == dof_handler->n_dofs(),
+           ExcWrongVectorSize(all_out[i].size(), dof_handler->n_dofs()));
+  for (unsigned int i=0; i<size; ++i)
+    for (unsigned int j=0; j<size; ++j)
       Assert(&all_in[i] != &all_out[j],
              ExcMessage ("Vectors cannot be used as input and output"
                          " at the same time!"));
@@ -335,7 +327,7 @@ interpolate (const std::vector<Vector<number> > &all_in,
                                               // data vectors on this
                                               // cell and prolong it
                                               // to its children
-             for (unsigned int j=0; j<in_size; ++j)
+             for (unsigned int j=0; j<size; ++j)
                {
                  for (unsigned int i=0; i<dofs_per_cell; ++i)
                    local_values(i)=all_in[j](indexptr->operator[](i));
@@ -358,7 +350,7 @@ interpolate (const std::vector<Vector<number> > &all_in,
                                               // distribute the
                                               // stored data to the
                                               // new vectors
-             for (unsigned int j=0; j<in_size; ++j)
+             for (unsigned int j=0; j<size; ++j)
                for (unsigned int i=0; i<dofs_per_cell; ++i)
                  all_out[j](dofs[i])=valuesptr->operator[](j)(i);
            }
@@ -375,6 +367,11 @@ template<int dim, typename number>
 void SolutionTransfer<dim, number>::interpolate(const Vector<number> &in,
                                                Vector<number>       &out) const
 {
+  Assert (in.size()==n_dofs_old,
+         ExcWrongVectorSize(in.size(), n_dofs_old));
+  Assert (out.size()==dof_handler->n_dofs(),
+         ExcWrongVectorSize(out.size(), dof_handler->n_dofs()));
+
   std::vector<Vector<number> > all_in(1);
   all_in[0] = in;
   std::vector<Vector<number> > all_out(1);

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.