]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
remove FEVALUESBASE template argument in IntegrationInfo
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 12 Mar 2010 20:46:18 +0000 (20:46 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 12 Mar 2010 20:46:18 +0000 (20:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@20815 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/mesh_worker.h
deal.II/deal.II/include/numerics/mesh_worker_info.h
deal.II/deal.II/include/numerics/mesh_worker_info.templates.h
deal.II/deal.II/include/numerics/mesh_worker_loop.h
deal.II/deal.II/source/numerics/mesh_worker_info.cc
deal.II/deal.II/source/numerics/mesh_worker_info.inst.in
deal.II/examples/step-38/step-38.cc
deal.II/examples/step-39/step-39.cc

index 4dc8718f46caf7e87c15baa749a46bb84a85abbd..3602b4331f6199be1178bbae8b6b3b42eb62db7a 100644 (file)
@@ -169,13 +169,13 @@ namespace MeshWorker
                                        * The info type expected by a
                                        * cell integrator.
                                        */
-      typedef IntegrationInfo<dim, FEValuesBase<dim, dim>, dim> CellInfo;
+      typedef IntegrationInfo<dim, dim> CellInfo;
 
                                       /**
                                        * The info type expected by a
                                        * face integrator.
                                        */
-      typedef IntegrationInfo<dim, FEFaceValuesBase<dim, dim>, dim> FaceInfo;
+      typedef IntegrationInfo<dim, dim> FaceInfo;
 
                                       /**
                                        * Initialize default values.
index f48fe4b9e19b825048f6a6cd3297a2ada013d422..767a5c4142718d85fc027e6891409b308abea43e 100644 (file)
@@ -474,15 +474,15 @@ namespace MeshWorker
  * @ingroup MeshWorker
  * @author Guido Kanschat, 2009
  */
-  template<int dim, class FEVALUESBASE, int spacedim = dim>
+  template<int dim, int spacedim = dim>
   class IntegrationInfo
   {
     private:
                                       /// vector of FEValues objects
-      std::vector<boost::shared_ptr<FEVALUESBASE> > fevalv;
+      std::vector<boost::shared_ptr<FEValuesBase<dim, spacedim> > > fevalv;
     public:
-      static const unsigned int dimension = FEVALUESBASE::dimension;
-      static const unsigned int space_dimension = FEVALUESBASE::space_dimension;
+      static const unsigned int dimension = dim;
+      static const unsigned int space_dimension = spacedim;
       
                                       /**
                                        * Constructor.
@@ -494,7 +494,7 @@ namespace MeshWorker
                                        * clone to be used by
                                        * WorksTream::run().
                                        */
-      IntegrationInfo(const IntegrationInfo<dim, FEVALUESBASE, spacedim>& other);
+      IntegrationInfo(const IntegrationInfo<dim, spacedim>& other);
       
                                       /**
                                        * Build all internal
@@ -560,7 +560,7 @@ namespace MeshWorker
                                        * exception, if applied to a
                                        * vector of elements.
                                        */
-      const FEVALUESBASE& fe_values () const;
+      const FEValuesBase<dim, spacedim>& fe_values () const;
 
                                       /// Access to finite elements
                                       /**
@@ -571,7 +571,7 @@ namespace MeshWorker
                                        *
                                        * @see DGBlockSplitApplication
                                        */
-      const FEVALUESBASE& fe_values (const unsigned int i) const;
+      const FEValuesBase<dim, spacedim>& fe_values (const unsigned int i) const;
 
                                       /**
                                        * The vector containing the
@@ -686,10 +686,10 @@ namespace MeshWorker
     public:
 
 /// The type of the info object for cells
-      typedef IntegrationInfo<dim, FEValuesBase<dim, spacedim>, spacedim> CellInfo;
+      typedef IntegrationInfo<dim, spacedim> CellInfo;
 
 /// The type of the info objects for faces.
-      typedef IntegrationInfo<dim, FEFaceValuesBase<dim, spacedim>, spacedim> FaceInfo;
+      typedef IntegrationInfo<dim, spacedim> FaceInfo;
       
       template <class WORKER>
       void initialize(const WORKER&,
@@ -1047,31 +1047,31 @@ namespace MeshWorker
 
 //----------------------------------------------------------------------//
 
-  template <int dim, class FVB, int spacedim>
-  inline const FVB&
-  IntegrationInfo<dim,FVB,spacedim>::fe_values() const
+  template <int dim, int spacedim>
+  inline const FEValuesBase<dim, spacedim>&
+  IntegrationInfo<dim,spacedim>::fe_values() const
   {
     AssertDimension(fevalv.size(), 1);
     return *fevalv[0];
   }
 
 
-  template <int dim, class FVB, int spacedim>
-  inline const FVB&
-  IntegrationInfo<dim,FVB,spacedim>::fe_values(unsigned int i) const
+  template <int dim, int spacedim>
+  inline const FEValuesBase<dim, spacedim>&
+  IntegrationInfo<dim,spacedim>::fe_values(unsigned int i) const
   {
     Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
     return *fevalv[i];
   }
 
 
-  template <int dim, class FVB, int spacedim>
+  template <int dim, int spacedim>
   inline void
-  IntegrationInfo<dim,FVB,spacedim>::reinit(const DoFInfo<dim, spacedim>& info)
+  IntegrationInfo<dim,spacedim>::reinit(const DoFInfo<dim, spacedim>& info)
   {
     for (unsigned int i=0;i<fevalv.size();++i)
       {
-       FVB& febase = *fevalv[i];
+       FEValuesBase<dim, spacedim>& febase = *fevalv[i];
        if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
          {
                                             // This is a subface
index 55c3bff7e17ea6587602c2b19c54a33293a5fd89..28c1afe4f976e9c5146473c62eb5910361ce4842 100644 (file)
@@ -72,8 +72,8 @@ namespace MeshWorker
 
 //----------------------------------------------------------------------//
 
-  template<int dim, class FVB, int sdim>
-  IntegrationInfo<dim,FVB,sdim>::IntegrationInfo()
+  template<int dim, int sdim>
+  IntegrationInfo<dim,sdim>::IntegrationInfo()
                  :
                  fevalv(0),
                  multigrid(false),
@@ -81,8 +81,8 @@ namespace MeshWorker
   {}
 
   
-  template<int dim, class FVB, int sdim>
-  IntegrationInfo<dim,FVB,sdim>::IntegrationInfo(const IntegrationInfo<dim,FVB,sdim>& other)
+  template<int dim, int sdim>
+  IntegrationInfo<dim,sdim>::IntegrationInfo(const IntegrationInfo<dim,sdim>& other)
                  :
                  multigrid(other.multigrid),
                  values(other.values),
@@ -94,20 +94,21 @@ namespace MeshWorker
     fevalv.resize(other.fevalv.size());
     for (unsigned int i=0;i<other.fevalv.size();++i)
       {
-       const FVB& p = *other.fevalv[i];
+       const FEValuesBase<dim,sdim>& p = *other.fevalv[i];
        const FEValues<dim,sdim>* pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
        const FEFaceValues<dim,sdim>* pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
        const FESubfaceValues<dim,sdim>* ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
        
        if (pc != 0)
-         fevalv[i] = typename boost::shared_ptr<FVB> (
-           reinterpret_cast<FVB*>(new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
-                                                          pc->get_quadrature(), pc->get_update_flags())));
+         fevalv[i] = typename boost::shared_ptr<FEValuesBase<dim,sdim> > (
+           reinterpret_cast<FEFaceValuesBase<dim,sdim>*>(
+             new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
+                                     pc->get_quadrature(), pc->get_update_flags())));
        else if (pf != 0)
-         fevalv[i] = typename boost::shared_ptr<FVB> (
+         fevalv[i] = typename boost::shared_ptr<FEValuesBase<dim,sdim> > (
            new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
        else if (ps != 0)
-         fevalv[i] = typename boost::shared_ptr<FVB> (
+         fevalv[i] = typename boost::shared_ptr<FEValuesBase<dim,sdim> > (
            new FESubfaceValues<dim,sdim> (ps->get_mapping(), ps->get_fe(), ps->get_quadrature(), ps->get_update_flags()));
        else
          Assert(false, ExcInternalError());
@@ -115,10 +116,10 @@ namespace MeshWorker
   }
   
 
-  template<int dim, class FVB, int sdim>
+  template<int dim, int sdim>
   template <class FEVALUES>
   void
-  IntegrationInfo<dim,FVB,sdim>::initialize(
+  IntegrationInfo<dim,sdim>::initialize(
     const FiniteElement<dim,sdim>& el,
     const Mapping<dim,sdim>& mapping,
     const Quadrature<FEVALUES::integral_dimension>& quadrature,
@@ -128,7 +129,7 @@ namespace MeshWorker
     if (block_info == 0 || block_info->local().size() == 0)
       {
        fevalv.resize(1);             
-       fevalv[0] = boost::shared_ptr<FVB> (
+       fevalv[0] = boost::shared_ptr<FEValuesBase<dim,sdim> > (
          new FEVALUES (mapping, el, quadrature, flags));
       }
     else
@@ -136,7 +137,7 @@ namespace MeshWorker
        fevalv.resize(el.n_base_elements());
        for (unsigned int i=0;i<fevalv.size();++i)
          {
-           fevalv[i] = boost::shared_ptr<FVB> (
+           fevalv[i] = boost::shared_ptr<FEValuesBase<dim,sdim> > (
              new FEVALUES (mapping, el.base_element(i), quadrature, flags));
          }
       }
@@ -144,9 +145,9 @@ namespace MeshWorker
   }
   
 
-  template<int dim, class FVB, int sdim>
+  template<int dim, int sdim>
   void
-  IntegrationInfo<dim,FVB,sdim>::initialize_data(
+  IntegrationInfo<dim,sdim>::initialize_data(
     const boost::shared_ptr<VectorDataBase<dim,sdim> > data)
   {
     global_data = data;
@@ -198,18 +199,18 @@ namespace MeshWorker
   }
 
 
-  template<int dim, class FVB, int sdim>
+  template<int dim, int sdim>
   void
-  IntegrationInfo<dim,FVB,sdim>::clear()
+  IntegrationInfo<dim,sdim>::clear()
   {
     fevalv.resize(0);
   }
 
 
 
-  template<int dim, class FVB, int sdim>
+  template<int dim, int sdim>
   void
-  IntegrationInfo<dim,FVB,sdim>::fill_local_data(const DoFInfo<dim, sdim>& info, bool split_fevalues)
+  IntegrationInfo<dim,sdim>::fill_local_data(const DoFInfo<dim, sdim>& info, bool split_fevalues)
   {
      if (split_fevalues)
        {
index 0e63154685bc3f1ba62bcfc7c205d61d20ebf606..a423fc2af627ee6f7af142b5b34ceb48063461cd 100644 (file)
@@ -253,17 +253,17 @@ namespace MeshWorker
       }
     
                                     // Loop over all cells
-    WorkStream::run(begin, end,
-                   std_cxx1x::bind(cell_action<INFOBOX, dim, spacedim, ITERATOR>, _1, _3, _2,
-                                   cell_worker, boundary_worker, face_worker, cells_first, true),
-                   std_cxx1x::bind(internal::assemble<dim,spacedim,ASSEMBLER>, _1, assembler),
-                   info, dof_info);
+//     WorkStream::run(begin, end,
+//                 std_cxx1x::bind(cell_action<INFOBOX, dim, spacedim, ITERATOR>, _1, _3, _2,
+//                                 cell_worker, boundary_worker, face_worker, cells_first, true),
+//                 std_cxx1x::bind(internal::assemble<dim,spacedim,ASSEMBLER>, _1, assembler),
+//                 info, dof_info);
     
-//     for (ITERATOR cell = begin; cell != end; ++cell)
-//       {
-//     cell_action(cell, dof_info, info, cell_worker, boundary_worker, face_worker, cells_first);
-//     dof_info.assemble(assembler);
-//       }
+    for (ITERATOR cell = begin; cell != end; ++cell)
+      {
+       cell_action(cell, dof_info, info, cell_worker, boundary_worker, face_worker, cells_first, true);
+       dof_info.assemble(assembler);
+      }
   }
 
 /**
index 9d418fd23bc3cea5c7bf2f4888ea8f2ebd068b22..663b29e61b39cede7924584b790efcd779d8d4bd 100644 (file)
@@ -28,18 +28,17 @@ namespace MeshWorker
 {
   template class LocalResults<double>;
   template class DoFInfo<deal_II_dimension,deal_II_dimension>;
+  template class IntegrationInfo<deal_II_dimension>;
   
-#include "mesh_worker_info.inst"
-  
-  template void IntegrationInfo<deal_II_dimension, FEValuesBase<deal_II_dimension> >
+  template void IntegrationInfo<deal_II_dimension>
   ::initialize<FEValues<deal_II_dimension> >(
     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
     const Quadrature<FEValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
-  template void IntegrationInfo<deal_II_dimension, FEFaceValuesBase<deal_II_dimension> >
+  template void IntegrationInfo<deal_II_dimension>
   ::initialize<FEFaceValues<deal_II_dimension> >(
     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
     const Quadrature<FEFaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
-  template void IntegrationInfo<deal_II_dimension, FEFaceValuesBase<deal_II_dimension> >
+  template void IntegrationInfo<deal_II_dimension>
   ::initialize<FESubfaceValues<deal_II_dimension> >(
     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
     const Quadrature<FESubfaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
index 8e123d5444717a093e6f10fad43bdff6428c5120..d79ce7b9bee2c48cd02b639bff236402ab57453c 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id: fe_field_function.inst.in 19046 2009-07-08 19:30:23Z bangerth $
 //
-//    Copyright (C) 2009 by the deal.II authors
+//    Copyright (C) 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -13,7 +13,6 @@
 
 for (FEV : FEVALUES_BASES)
 {
-  template class IntegrationInfo<deal_II_dimension, FEV>;
 }
 
 for (FEV : FEVALUES_BASES)
index 8eff89bddf586ce3be512f2ebccee5d2a13fef5c..e672ba2a29915ca6130b5af022e64fd6f97b12b9 100644 (file)
@@ -463,7 +463,7 @@ void Step12<dim>::integrate_cell_term (DoFInfo& dinfo, CellInfo& info)
 template <int dim>
 void Step12<dim>::integrate_boundary_term (DoFInfo& dinfo, FaceInfo& info)
 {
-  const FEFaceValuesBase<dim>& fe_v = info.fe_values();
+  const FEValuesBase<dim>& fe_v = info.fe_values();
   FullMatrix<double>& local_matrix = dinfo.matrix(0).matrix;
   Vector<double>& local_vector = dinfo.vector(0).block(0);
 
@@ -514,12 +514,12 @@ void Step12<dim>::integrate_face_term (DoFInfo& dinfo1, DoFInfo& dinfo2,
                                   // etc., we use the
                                   // FEFaceValuesBase object of the
                                   // first argument.
-  const FEFaceValuesBase<dim>& fe_v = info1.fe_values();
+  const FEValuesBase<dim>& fe_v = info1.fe_values();
 
                                   // For additional shape functions,
                                   // we have to ask the neighbors
                                   // FEFaceValuesBase.
-  const FEFaceValuesBase<dim>& fe_v_neighbor = info2.fe_values();
+  const FEValuesBase<dim>& fe_v_neighbor = info2.fe_values();
 
                                   // Then we get references to the
                                   // four local matrices. The letters
index c2edd8872066bed42dec81baf0a5cb8a17306b94..558e32645e149b043d4e441ea7d499278563f704 100644 (file)
@@ -137,7 +137,7 @@ void MatrixIntegrator<dim>::bdry(
   MeshWorker::DoFInfo<dim>& dinfo,
   typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info)
 {
-  const FEFaceValuesBase<dim>& fe = info.fe_values();
+  const FEValuesBase<dim>& fe = info.fe_values();
   FullMatrix<double>& local_matrix = dinfo.matrix(0,false).matrix;
   
   const unsigned int deg = fe.get_fe().tensor_degree();
@@ -160,8 +160,8 @@ void MatrixIntegrator<dim>::face(
   typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
   typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2)
 {
-  const FEFaceValuesBase<dim>& fe1 = info1.fe_values();
-  const FEFaceValuesBase<dim>& fe2 = info2.fe_values();
+  const FEValuesBase<dim>& fe1 = info1.fe_values();
+  const FEValuesBase<dim>& fe2 = info2.fe_values();
   FullMatrix<double>& matrix_v1u1 = dinfo1.matrix(0,false).matrix;
   FullMatrix<double>& matrix_v1u2 = dinfo1.matrix(0,true).matrix;
   FullMatrix<double>& matrix_v2u1 = dinfo2.matrix(0,true).matrix;
@@ -224,7 +224,7 @@ void RHSIntegrator<dim>::cell(MeshWorker::DoFInfo<dim>&, typename MeshWorker::In
 template <int dim>
 void RHSIntegrator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info)
 {
-  const FEFaceValuesBase<dim>& fe = info.fe_values();
+  const FEValuesBase<dim>& fe = info.fe_values();
   Vector<double>& local_vector = dinfo.vector(0).block(0);
   
   std::vector<double> boundary_values(fe.n_quadrature_points);
@@ -279,7 +279,7 @@ void Estimator<dim>::cell(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::
 template <int dim>
 void Estimator<dim>::bdry(MeshWorker::DoFInfo<dim>& dinfo, typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info)
 {
-  const FEFaceValuesBase<dim>& fe = info.fe_values();
+  const FEValuesBase<dim>& fe = info.fe_values();
   
   std::vector<double> boundary_values(fe.n_quadrature_points);
   exact_solution.value_list(fe.get_quadrature_points(), boundary_values);
@@ -301,7 +301,7 @@ void Estimator<dim>::face(MeshWorker::DoFInfo<dim>& dinfo1,
                          typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info1,
                          typename MeshWorker::IntegrationWorker<dim>::FaceInfo& info2)
 {
-  const FEFaceValuesBase<dim>& fe = info1.fe_values();
+  const FEValuesBase<dim>& fe = info1.fe_values();
   const std::vector<double>& uh1 = info1.values[0][0];
   const std::vector<double>& uh2 = info2.values[0][0];
   const std::vector<Tensor<1,dim> >& Duh1 = info1.gradients[0][0];
@@ -708,7 +708,8 @@ Step39<dim>::run(unsigned int n_steps)
 
 int main()
 {
+  deallog.log_execution_time(true);
   FE_DGQ<2> fe1(3);
   Step39<2> test1(fe1);
-  test1.run(13);
+  test1.run(20);
 }

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.