]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Also pass the locations of evaluation points through DataPostprocessor.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Oct 2010 12:41:43 +0000 (12:41 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Oct 2010 12:41:43 +0000 (12:41 +0000)
git-svn-id: https://svn.dealii.org/trunk@22292 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/data_out.h
deal.II/deal.II/include/numerics/data_out_faces.h
deal.II/deal.II/include/numerics/data_out_rotation.h
deal.II/deal.II/include/numerics/data_postprocessor.h
deal.II/deal.II/source/numerics/data_out.cc
deal.II/deal.II/source/numerics/data_out_faces.cc
deal.II/deal.II/source/numerics/data_out_rotation.cc
deal.II/deal.II/source/numerics/data_postprocessor.cc
deal.II/doc/news/changes.h
deal.II/examples/step-29/step-29.cc
deal.II/examples/step-33/step-33.cc

index 2fe312fdbe98337f549d7417ee8779e5859af74f..dffd98ea68c877e1d35d27b508640bf3f63a7008 100644 (file)
@@ -122,6 +122,8 @@ namespace internal
        const dealii::hp::QCollection<dim> q_collection;
        const dealii::hp::MappingCollection<dim,spacedim> mapping_collection;
        dealii::hp::FEValues<dim,spacedim> x_fe_values;
+       
+       std::vector<Point<spacedim> > patch_evaluation_points;
 
        const std::vector<std::vector<unsigned int> > *cell_to_patch_index_map;
     };
index a4cdc00d91ed57dffee3c305128faa790bbbef4d..c26f4da8ec6b05b4fa4c6b0656b81480e0898612 100644 (file)
@@ -54,6 +54,7 @@ namespace internal
        dealii::hp::FEFaceValues<dim> x_fe_values;
 
        std::vector<Point<dim> > patch_normals;
+       std::vector<Point<spacedim> > patch_evaluation_points;
     };
   }
 }
index af8d501a9e573c9687ee2cc700c79a788c4905a9..f293f9f73a9583677f62b1eef001b7518751b9f0 100644 (file)
@@ -53,6 +53,8 @@ namespace internal
        
        const dealii::hp::QCollection<dim> q_collection;
        dealii::hp::FEValues<dim,spacedim> x_fe_values;
+       
+       std::vector<Point<spacedim> > patch_evaluation_points;
     };
   }
 }
index 8fd4567a9c572cf9dd4f772f489850ff83d7fa86..910eb4204eaa05602da6a18a6fbe2444f582286e 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2007, 2008 by the deal.II authors
+//    Copyright (C) 2007, 2008, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -80,6 +80,28 @@ class DataPostprocessor: public Subscriptor
                                      */
     virtual ~DataPostprocessor ();
 
+                                    /**
+                                     * @deprecated
+                                     *
+                                     * This function only exists for backward
+                                     * compatibility as this is the interface
+                                     * provided by previous versions of the
+                                     * library. The default implementation of
+                                     * the other function of same name below
+                                     * calls this function by simply dropping
+                                     * the argument that denotes the
+                                     * evaluation points. Since this function
+                                     * might at one point go away, you should
+                                     * overload the other function instead.
+                                     */
+    virtual
+    void
+    compute_derived_quantities_scalar (const std::vector<double>         &uh,
+                                      const std::vector<Tensor<1,dim> > &duh,
+                                      const std::vector<Tensor<2,dim> > &dduh,
+                                      const std::vector<Point<dim> >    &normals,
+                                      std::vector<Vector<double> >      &computed_quantities) const;
+                                      
                                     /**
                                      * This is the main function which actually
                                      * performs the postprocessing. The last
@@ -113,16 +135,22 @@ class DataPostprocessor: public Subscriptor
                                       const std::vector<Tensor<1,dim> > &duh,
                                       const std::vector<Tensor<2,dim> > &dduh,
                                       const std::vector<Point<dim> >    &normals,
-                                      std::vector<Vector<double> >      &computed_quantities) const;
+                                      const std::vector<Point<dim> >    &evaluation_points,
+                                      std::vector<Vector<double> >      &computed_quantities) const;                          
     
                                     /**
-                                     * Same as above function, but
-                                     * this function is called when
-                                     * the original data vector
-                                     * represents vector data,
-                                     * i.e. the finite element in use
-                                     * has multiple vector
-                                     * components.
+                                     * @deprecated
+                                     *
+                                     * This function only exists for backward
+                                     * compatibility as this is the interface
+                                     * provided by previous versions of the
+                                     * library. The default implementation of
+                                     * the other function of same name below
+                                     * calls this function by simply dropping
+                                     * the argument that denotes the
+                                     * evaluation points. Since this function
+                                     * might at one point go away, you should
+                                     * overload the other function instead.
                                      */
     virtual
     void
@@ -131,6 +159,24 @@ class DataPostprocessor: public Subscriptor
                                       const std::vector<std::vector<Tensor<2,dim> > > &dduh,
                                       const std::vector<Point<dim> >                  &normals,
                                       std::vector<Vector<double> >                    &computed_quantities) const;
+                                      
+                                    /**
+                                     * Same as the
+                                     * compute_derived_quantities_scalar()
+                                     * function, but this function is called
+                                     * when the original data vector
+                                     * represents vector data, i.e. the
+                                     * finite element in use has multiple
+                                     * vector components.
+                                     */
+    virtual
+    void
+    compute_derived_quantities_vector (const std::vector<Vector<double> >              &uh,
+                                      const std::vector<std::vector<Tensor<1,dim> > > &duh,
+                                      const std::vector<std::vector<Tensor<2,dim> > > &dduh,
+                                      const std::vector<Point<dim> >                  &normals,
+                                      const std::vector<Point<dim> >                  &evaluation_points,
+                                      std::vector<Vector<double> >                    &computed_quantities) const;                            
 
                                     /**
                                      * Return the vector of strings describing
index 845d23e2e321044806ffbdd57e8124ad4d0a6fe2..0926a72d9971eb804048187f35b4ae1272ea022e 100644 (file)
@@ -768,12 +768,18 @@ build_one_patch (const std::pair<cell_iterator, unsigned int> *cell_and_index,
                  if (update_flags & update_hessians)
                    this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                    data.patch_hessians);
+                                                                   
+                 if (update_flags & update_quadrature_points)
+                   data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+                   
+                   
                  std::vector<Point<DH::space_dimension> > dummy_normals;
                  postprocessor->
                    compute_derived_quantities_scalar(data.patch_values,
                                                      data.patch_gradients,
                                                      data.patch_hessians,
                                                      dummy_normals,
+                                                     data.patch_evaluation_points,
                                                      data.postprocessed_values[dataset]);
                }
              else
@@ -791,12 +797,18 @@ build_one_patch (const std::pair<cell_iterator, unsigned int> *cell_and_index,
                  if (update_flags & update_hessians)
                    this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                    data.patch_hessians_system);
+                 
+                 if (update_flags & update_quadrature_points)
+                   data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+                 
                  std::vector<Point<DH::space_dimension> > dummy_normals;
+                 
                  postprocessor->
                    compute_derived_quantities_vector(data.patch_values_system,
                                                      data.patch_gradients_system,
                                                      data.patch_hessians_system,
                                                      dummy_normals,
+                                                     data.patch_evaluation_points,
                                                      data.postprocessed_values[dataset]);
                }
 
index 22b0716203d4b8347f0d95b02ffdc3013b73f3cc..3b4b62e340ba6f94b7b1692d658afef2eda0b43f 100644 (file)
@@ -170,11 +170,16 @@ build_one_patch (const FaceDescriptor *cell_and_face,
                  if (update_flags & update_hessians)
                    this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                    data.patch_hessians);
+                 
+                 if (update_flags & update_quadrature_points)
+                   data.patch_evaluation_points = fe_patch_values.get_quadrature_points();                                                
+                                                                   
                  postprocessor->
                    compute_derived_quantities_scalar(data.patch_values,
                                                      data.patch_gradients,
                                                      data.patch_hessians,
                                                      data.patch_normals,
+                                                     data.patch_evaluation_points,
                                                      data.postprocessed_values[dataset]);
                }
              else
@@ -192,11 +197,16 @@ build_one_patch (const FaceDescriptor *cell_and_face,
                  if (update_flags & update_hessians)
                    this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                    data.patch_hessians_system);
+                                                                   
+                 if (update_flags & update_quadrature_points)
+                   data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+                   
                  postprocessor->
                    compute_derived_quantities_vector(data.patch_values_system,
                                                      data.patch_gradients_system,
                                                      data.patch_hessians_system,
                                                      data.patch_normals,
+                                                     data.patch_evaluation_points,
                                                      data.postprocessed_values[dataset]);
                }
 
index d456b165d01665e479b3f86ce1320c39ea647c1c..148c5b86276780c52ad63fc67685f4771a9b63f6 100644 (file)
@@ -220,6 +220,9 @@ build_one_patch (const cell_iterator *cell,
                      if (update_flags & update_hessians)
                        this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                        data.patch_hessians);
+                                                                       
+                       if (update_flags & update_quadrature_points)
+                     data.patch_evaluation_points = fe_patch_values.get_quadrature_points();                                           
 
                      std::vector<Point<DH::space_dimension> > dummy_normals;
                      postprocessor->
@@ -227,6 +230,7 @@ build_one_patch (const cell_iterator *cell,
                                                          data.patch_gradients,
                                                          data.patch_hessians,
                                                          dummy_normals,
+                                                         data.patch_evaluation_points,
                                                          data.postprocessed_values[dataset]);
                    }
                  else
@@ -244,12 +248,17 @@ build_one_patch (const cell_iterator *cell,
                      if (update_flags & update_hessians)
                        this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                        data.patch_hessians_system);
+                                                                       
+                       if (update_flags & update_quadrature_points)
+                     data.patch_evaluation_points = fe_patch_values.get_quadrature_points();                                           
+                                                                       
                      std::vector<Point<DH::space_dimension> > dummy_normals;
                      postprocessor->
                        compute_derived_quantities_vector(data.patch_values_system,
                                                          data.patch_gradients_system,
                                                          data.patch_hessians_system,
                                                          dummy_normals,
+                                                         data.patch_evaluation_points,
                                                          data.postprocessed_values[dataset]);
                    }
                      
index 217b201479431e4c752a439351afced92f1392e9..502b24452ceb2c99759ba988cfe69f64a888ba74 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2007 by the deal.II authors
+//    Copyright (C) 2007, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -35,11 +35,25 @@ compute_derived_quantities_scalar (const std::vector<double>         &/*uh*/,
 }
 
 
+template <int dim>
+void
+DataPostprocessor<dim>::
+compute_derived_quantities_scalar (const std::vector<double>         &uh,
+                                  const std::vector<Tensor<1,dim> > &duh,
+                                  const std::vector<Tensor<2,dim> > &dduh,
+                                  const std::vector<Point<dim> >    &normals,
+                                  const std::vector<Point<dim> >    &/*evaluation_points*/,
+                                  std::vector<Vector<double> >      &computed_quantities) const
+{
+  compute_derived_quantities_scalar(uh, duh, dduh, normals, computed_quantities);
+}
+
+
 
 template <int dim>
 void
 DataPostprocessor<dim>::
-compute_derived_quantities_vector (const std::vector<Vector<double> >              &/*uh*/,
+compute_derived_quantities_vector (const std::vector<Vector<double> > &/*uh*/,
                                   const std::vector<std::vector<Tensor<1,dim> > > &/*duh*/,
                                   const std::vector<std::vector<Tensor<2,dim> > > &/*dduh*/,
                                   const std::vector<Point<dim> >                  &/*normals*/,
@@ -51,6 +65,21 @@ compute_derived_quantities_vector (const std::vector<Vector<double> >
 
 
 
+template <int dim>
+void
+DataPostprocessor<dim>::
+compute_derived_quantities_vector (const std::vector<Vector<double> > &uh,
+                                  const std::vector<std::vector<Tensor<1,dim> > > &duh,
+                                  const std::vector<std::vector<Tensor<2,dim> > > &dduh,
+                                  const std::vector<Point<dim> >                  &normals,
+                                  const std::vector<Point<dim> >                  &/*evaluation_points*/,
+                                  std::vector<Vector<double> >                    &computed_quantities) const
+{
+  compute_derived_quantities_vector(uh, duh, dduh, normals, computed_quantities);
+}
+
+
+
 template <int dim>
 std::vector<DataComponentInterpretation::DataComponentInterpretation>
 DataPostprocessor<dim>::get_data_component_interpretation () const
index fdfca2971f1769d31b9f9d11fd155fe38b228d3c..7f332bef92ccfb1429a15eb9c360e7e2d4858609 100644 (file)
@@ -280,6 +280,15 @@ through DoFHandler::get_tria() and DoFHandler::get_fe(), respectively.
 <h3>deal.II</h3>
 
 <ol>
+  <li><p>Changed: The DataPostprocessor functions now take an additional
+  argument that indicates the location of an evaluation point. For backward
+  compatibility, the old function signature still exists so that applications
+  that overload one of the existing functions continue to work.
+  The old signature has been deprecated, however, and will be removed in a
+  future version.
+  <br>
+  (Scott Miller 2010/10/08)
+  </p></li>
 
   <li><p>Changed: FETools is now a namespace rather than a class with only
   static member functions.
index 80e692366faca6cfa7e1b6cc59481074e633c5fb..2ca1113d73b4ad6f64b6b42503f6357bdb3f87b8 100644 (file)
@@ -2,7 +2,7 @@
 /* Author: Moritz Allmaras, Texas A&M University, 2007 */
 
 /*                                                                */
-/*    Copyright (C) 2007, 2008 by the deal.II authors and M. Allmaras   */
+/*    Copyright (C) 2007, 2008, 2010 by the deal.II authors and M. Allmaras   */
 /*                                                                */
 /*    This file is subject to QPL and may not be  distributed     */
 /*    without copyright and license information. Please refer     */
@@ -414,6 +414,7 @@ class ComputeIntensity : public DataPostprocessor<dim>
       const std::vector< std::vector< Tensor< 1, dim > > > &, 
       const std::vector< std::vector< Tensor< 2, dim > > > &, 
       const std::vector< Point< dim > > &,
+      const std::vector<Point<dim> > &,
       std::vector< Vector< double > > &
     ) const;
 
@@ -518,6 +519,7 @@ ComputeIntensity<dim>::compute_derived_quantities_vector (
   const std::vector< std::vector< Tensor< 1, dim > > >  & /*duh*/,
   const std::vector< std::vector< Tensor< 2, dim > > >  & /*dduh*/,
   const std::vector< Point< dim > >                     & /*normals*/,
+  const std::vector<Point<dim> >                        & /*evaluation_points*/,
   std::vector< Vector< double > >                       & computed_quantities
 ) const
 {
index bd72e790f5904159f593ad91ba94bf52ba3ed66f..1d1d4cc3de70fa7c53a8d47257a36162c1c6e1c8 100644 (file)
@@ -711,6 +711,7 @@ struct EulerEquations
                                           const std::vector<std::vector<Tensor<1,dim> > > &duh,
                                           const std::vector<std::vector<Tensor<2,dim> > > &dduh,
                                           const std::vector<Point<dim> >                  &normals,
+                                          const std::vector<Point<dim> >                  &evaluation_points,
                                           std::vector<Vector<double> >                    &computed_quantities) const;
 
        virtual std::vector<std::string> get_names () const;
@@ -765,6 +766,7 @@ compute_derived_quantities_vector (const std::vector<Vector<double> >
                                   const std::vector<std::vector<Tensor<1,dim> > > &duh,
                                   const std::vector<std::vector<Tensor<2,dim> > > &/*dduh*/,
                                   const std::vector<Point<dim> >                  &/*normals*/,
+                                  const std::vector<Point<dim> >                  &/*evaluation_points*/,
                                   std::vector<Vector<double> >                    &computed_quantities) const
 {
                                   // At the beginning of the function, let us

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.