include/deal.II/base/qprojector.h

00001 //---------------------------------------------------------------------------
00002 //    @f$Id: qprojector.h 25345 2012-03-31 08:37:04Z bangerth @f$
00003 //
00004 //    Copyright (C) 2005, 2006, 2007, 2008, 2009, 2012 by the deal.II authors
00005 //
00006 //    This file is subject to QPL and may not be  distributed
00007 //    without copyright and license information. Please refer
00008 //    to the file deal.II/doc/license.html for the  text  and
00009 //    further information on this license.
00010 //
00011 //---------------------------------------------------------------------------
00012 #ifndef __deal2__qprojector_h
00013 #define __deal2__qprojector_h
00014 
00015 
00016 #include <deal.II/base/quadrature.h>
00017 #include <deal.II/base/geometry_info.h>
00018 
00019 DEAL_II_NAMESPACE_OPEN
00022 
00023 
00072 template <int dim>
00073 class QProjector
00074 {
00075   public:
00083     typedef Quadrature<dim-1> SubQuadrature;
00084 
00093     static void project_to_face (const SubQuadrature &quadrature,
00094                                  const unsigned int      face_no,
00095                                  std::vector<Point<dim> > &q_points);
00096 
00105     static Quadrature<dim>
00106     project_to_face (const SubQuadrature &quadrature,
00107                      const unsigned int      face_no);
00108 
00123     static void project_to_subface (const SubQuadrature       &quadrature,
00124                                     const unsigned int         face_no,
00125                                     const unsigned int         subface_no,
00126                                     std::vector<Point<dim> >  &q_points,
00127                                     const RefinementCase<dim-1> &ref_case=RefinementCase<dim-1>::isotropic_refinement);
00128 
00144     static Quadrature<dim>
00145     project_to_subface (const SubQuadrature       &quadrature,
00146                         const unsigned int         face_no,
00147                         const unsigned int         subface_no,
00148                         const RefinementCase<dim-1> &ref_case=RefinementCase<dim-1>::isotropic_refinement);
00149 
00180     static Quadrature<dim>
00181     project_to_all_faces (const SubQuadrature &quadrature);
00182 
00205     static Quadrature<dim>
00206     project_to_all_subfaces (const SubQuadrature &quadrature);
00207 
00227     static
00228     Quadrature<dim>
00229     project_to_child (const Quadrature<dim>  &quadrature,
00230                       const unsigned int      child_no);
00231 
00250     static
00251     Quadrature<dim>
00252     project_to_all_children (const Quadrature<dim>  &quadrature);
00253 
00261     static
00262     Quadrature<dim>
00263     project_to_line(const Quadrature<1>& quadrature,
00264                     const Point<dim>& p1,
00265                     const Point<dim>& p2);
00266 
00292     class DataSetDescriptor
00293     {
00294       public:
00304         DataSetDescriptor ();
00305 
00317         static DataSetDescriptor cell ();
00318 
00339         static
00340         DataSetDescriptor
00341         face (const unsigned int face_no,
00342               const bool         face_orientation,
00343               const bool         face_flip,
00344               const bool         face_rotation,
00345               const unsigned int n_quadrature_points);
00346 
00371         static
00372         DataSetDescriptor
00373         subface (const unsigned int face_no,
00374                  const unsigned int subface_no,
00375                  const bool         face_orientation,
00376                  const bool         face_flip,
00377                  const bool         face_rotation,
00378                  const unsigned int n_quadrature_points,
00379                  const internal::SubfaceCase<dim> ref_case=internal::SubfaceCase<dim>::case_isotropic);
00380 
00395         operator unsigned int () const;
00396 
00397       private:
00403         const unsigned int dataset_offset;
00404 
00412         DataSetDescriptor (const unsigned int dataset_offset);
00413     };
00414 
00415   private:
00429     static Quadrature<2> reflect (const Quadrature<2> &q);
00430 
00447     static Quadrature<2> rotate (const Quadrature<2> &q,
00448                                  const unsigned int n_times);
00449 };
00450 
00454 // -------------------  inline and template functions ----------------
00455 
00456 
00457 
00458 template <int dim>
00459 inline
00460 QProjector<dim>::DataSetDescriptor::
00461 DataSetDescriptor (const unsigned int dataset_offset)
00462                 :
00463                 dataset_offset (dataset_offset)
00464 {}
00465 
00466 
00467 template <int dim>
00468 inline
00469 QProjector<dim>::DataSetDescriptor::
00470 DataSetDescriptor ()
00471                 :
00472                 dataset_offset (numbers::invalid_unsigned_int)
00473 {}
00474 
00475 
00476 
00477 template <int dim>
00478 typename QProjector<dim>::DataSetDescriptor
00479 QProjector<dim>::DataSetDescriptor::cell ()
00480 {
00481   return 0;
00482 }
00483 
00484 
00485 
00486 template <int dim>
00487 inline
00488 QProjector<dim>::DataSetDescriptor::operator unsigned int () const
00489 {
00490   return dataset_offset;
00491 }
00492 
00493 
00494 /* -------------- declaration of explicit specializations ------------- */
00495 
00496 #ifndef DOXYGEN
00497 
00498 
00499 template <>
00500 void
00501 QProjector<1>::project_to_face (const Quadrature<0> &,
00502                                 const unsigned int,
00503                                 std::vector<Point<1> > &);
00504 template <>
00505 void
00506 QProjector<2>::project_to_face (const Quadrature<1>      &quadrature,
00507                                 const unsigned int        face_no,
00508                                 std::vector<Point<2> >   &q_points);
00509 template <>
00510 void
00511 QProjector<3>::project_to_face (const Quadrature<2>    &quadrature,
00512                                 const unsigned int      face_no,
00513                                 std::vector<Point<3> > &q_points);
00514 
00515 template <>
00516 Quadrature<1>
00517 QProjector<1>::project_to_all_faces (const Quadrature<0> &quadrature);
00518 
00519 
00520 template <>
00521 void
00522 QProjector<1>::project_to_subface (const Quadrature<0> &,
00523                                    const unsigned int,
00524                                    const unsigned int,
00525                                    std::vector<Point<1> > &,
00526                                    const RefinementCase<0> &);
00527 template <>
00528 void
00529 QProjector<2>::project_to_subface (const Quadrature<1>    &quadrature,
00530                                    const unsigned int      face_no,
00531                                    const unsigned int      subface_no,
00532                                    std::vector<Point<2> > &q_points,
00533                                    const RefinementCase<1> &);
00534 template <>
00535 void
00536 QProjector<3>::project_to_subface (const Quadrature<2>       &quadrature,
00537                                    const unsigned int         face_no,
00538                                    const unsigned int         subface_no,
00539                                    std::vector<Point<3> >    &q_points,
00540                                    const RefinementCase<2> &face_ref_case);
00541 
00542 template <>
00543 Quadrature<1>
00544 QProjector<1>::project_to_all_subfaces (const Quadrature<0> &quadrature);
00545 
00546 
00547 template <>
00548 bool
00549 QIterated<1>::uses_both_endpoints (const Quadrature<1> &base_quadrature);
00550 
00551 template <>
00552 QIterated<1>::QIterated (const Quadrature<1> &base_quadrature,
00553                          const unsigned int   n_copies);
00554 
00555 
00556 #endif // DOXYGEN
00557 DEAL_II_NAMESPACE_CLOSE
00558 
00559 #endif
00560 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Wed May 23 2012 06:07:31 by doxygen 1.7.3