00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__mapping_q1_h
00013 #define __deal2__mapping_q1_h
00014
00015
00016 #include <deal.II/base/derivative_form.h>
00017 #include <deal.II/base/config.h>
00018 #include <deal.II/base/table.h>
00019 #include <deal.II/base/qprojector.h>
00020 #include <deal.II/grid/tria_iterator.h>
00021 #include <deal.II/dofs/dof_accessor.h>
00022 #include <deal.II/fe/mapping.h>
00023
00024 #include <cmath>
00025
00026 DEAL_II_NAMESPACE_OPEN
00027
00030
00031
00051 template <int dim, int spacedim=dim>
00052 class MappingQ1 : public Mapping<dim,spacedim>
00053 {
00054 public:
00058 MappingQ1 ();
00059
00060 virtual Point<spacedim>
00061 transform_unit_to_real_cell (
00062 const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00063 const Point<dim> &p) const;
00064
00082 virtual Point<dim>
00083 transform_real_to_unit_cell (
00084 const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00085 const Point<spacedim> &p) const;
00086
00087 virtual void
00088 transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
00089 VectorSlice<std::vector<Tensor<1,spacedim> > > output,
00090 const typename Mapping<dim,spacedim>::InternalDataBase &internal,
00091 const MappingType type) const;
00092
00093 virtual void
00094 transform (const VectorSlice<const std::vector<DerivativeForm<1, dim,spacedim> > > input,
00095 VectorSlice<std::vector<Tensor<2,spacedim> > > output,
00096 const typename Mapping<dim,spacedim>::InternalDataBase &internal,
00097 const MappingType type) const;
00098
00099 virtual
00100 void
00101 transform (const VectorSlice<const std::vector<Tensor<2, dim> > > input,
00102 VectorSlice<std::vector<Tensor<2,spacedim> > > output,
00103 const typename Mapping<dim,spacedim>::InternalDataBase &internal,
00104 const MappingType type) const;
00105
00106
00107 protected:
00113 template < int rank >
00114 void
00115 transform_fields(const VectorSlice<const std::vector<Tensor<rank,dim> > > input,
00116 VectorSlice< std::vector<Tensor<rank,spacedim> > > output,
00117 const typename Mapping<dim,spacedim>::InternalDataBase &internal,
00118 const MappingType type) const;
00122 template < int rank >
00123 void
00124 transform_gradients(const VectorSlice<const std::vector<Tensor<rank,dim> > > input,
00125 VectorSlice< std::vector<Tensor<rank,spacedim> > > output,
00126 const typename Mapping<dim,spacedim>::InternalDataBase &internal,
00127 const MappingType type) const;
00131 template < int rank >
00132 void
00133 transform_differential_forms(
00134 const VectorSlice<const std::vector<DerivativeForm<rank, dim, spacedim> > > input,
00135 VectorSlice<std::vector<DerivativeForm<rank, spacedim, spacedim> > > output,
00136 const typename Mapping<dim,spacedim>::InternalDataBase &internal,
00137 const MappingType type) const;
00138
00139 public:
00140
00141
00142
00143
00149 virtual
00150 Mapping<dim,spacedim> * clone () const;
00151
00156 class InternalData : public Mapping<dim,spacedim>::InternalDataBase
00157 {
00158 public:
00163 InternalData(const unsigned int n_shape_functions);
00164
00172 double shape (const unsigned int qpoint,
00173 const unsigned int shape_nr) const;
00174
00179 double &shape (const unsigned int qpoint,
00180 const unsigned int shape_nr);
00181
00187 Tensor<1,dim> derivative (const unsigned int qpoint,
00188 const unsigned int shape_nr) const;
00189
00195 Tensor<1,dim> &derivative (const unsigned int qpoint,
00196 const unsigned int shape_nr);
00197
00203 Tensor<2,dim> second_derivative (const unsigned int qpoint,
00204 const unsigned int shape_nr) const;
00205
00211 Tensor<2,dim> &second_derivative (const unsigned int qpoint,
00212 const unsigned int shape_nr);
00213
00220 virtual std::size_t memory_consumption () const;
00221
00229 std::vector<double> shape_values;
00230
00238 std::vector<Tensor<1,dim> > shape_derivatives;
00239
00248 std::vector<Tensor<2,dim> > shape_second_derivatives;
00249
00269 std::vector<DerivativeForm<1,dim, spacedim > > covariant;
00270
00282 std::vector< DerivativeForm<1,dim,spacedim> > contravariant;
00283
00306 std::vector<std::vector<Tensor<1,dim> > > unit_tangentials;
00307
00311 std::vector<std::vector<Tensor<1,spacedim> > > aux;
00312
00318 std::vector<Point<spacedim> > mapping_support_points;
00319
00325 typename Triangulation<dim,spacedim>::cell_iterator cell_of_current_support_points;
00326
00334 bool is_mapping_q1_data;
00335
00349 unsigned int n_shape_functions;
00350 };
00351
00359 typedef
00360 typename QProjector<dim>::DataSetDescriptor
00361 DataSetDescriptor;
00362
00367 virtual void
00368 fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00369 const Quadrature<dim> &quadrature,
00370 typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
00371 typename std::vector<Point<spacedim> > &quadrature_points,
00372 std::vector<double> &JxW_values,
00373 std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
00374 std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
00375 std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
00376 std::vector<Point<spacedim> > &cell_normal_vectors,
00377 CellSimilarity::Similarity &cell_similarity) const;
00378
00383 virtual void
00384 fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00385 const unsigned int face_no,
00386 const Quadrature<dim-1> &quadrature,
00387 typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
00388 typename std::vector<Point<spacedim> > &quadrature_points,
00389 std::vector<double> &JxW_values,
00390 typename std::vector<Tensor<1,spacedim> > &boundary_form,
00391 typename std::vector<Point<spacedim> > &normal_vectors) const ;
00392
00397 virtual void
00398 fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00399 const unsigned int face_no,
00400 const unsigned int sub_no,
00401 const Quadrature<dim-1>& quadrature,
00402 typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
00403 typename std::vector<Point<spacedim> > &quadrature_points,
00404 std::vector<double> &JxW_values,
00405 typename std::vector<Tensor<1,spacedim> > &boundary_form,
00406 typename std::vector<Point<spacedim> > &normal_vectors) const ;
00407
00420 void compute_shapes (const std::vector<Point<dim> > &unit_points,
00421 InternalData &data) const;
00422
00431 void compute_data (const UpdateFlags flags,
00432 const Quadrature<dim> &quadrature,
00433 const unsigned int n_orig_q_points,
00434 InternalData &data) const;
00435
00448 void compute_face_data (const UpdateFlags flags,
00449 const Quadrature<dim> &quadrature,
00450 const unsigned int n_orig_q_points,
00451 InternalData &data) const;
00452
00457 void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00458 const unsigned int npts,
00459 const DataSetDescriptor data_set,
00460 const CellSimilarity::Similarity cell_similarity,
00461 InternalData &data,
00462 std::vector<Point<spacedim> > &quadrature_points) const;
00463
00468 void compute_fill_face (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00469 const unsigned int face_no,
00470 const unsigned int subface_no,
00471 const unsigned int npts,
00472 const DataSetDescriptor data_set,
00473 const std::vector<double> &weights,
00474 InternalData &mapping_data,
00475 std::vector<Point<spacedim> > &quadrature_points,
00476 std::vector<double> &JxW_values,
00477 std::vector<Tensor<1,spacedim> > &boundary_form,
00478 std::vector<Point<spacedim> > &normal_vectors) const;
00479
00484 virtual void compute_shapes_virtual (const std::vector<Point<dim> > &unit_points,
00485 InternalData &data) const;
00486
00513 Point<spacedim> transform_unit_to_real_cell_internal (const InternalData &mdata) const;
00514
00549 void transform_real_to_unit_cell_internal (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00550 const Point<spacedim> &p,
00551 InternalData &mdata,
00552 Point<dim> &p_unit) const;
00553
00554
00555
00556
00557
00558
00559
00560
00565 virtual
00566 bool preserves_vertex_locations () const;
00567
00578 DeclException0(ExcTransformationFailed);
00579
00580
00581 protected:
00582
00583 template<int dim_>
00584 void transform_real_to_unit_cell_internal_codim1
00585 (const typename Triangulation<dim_,dim_+1>::cell_iterator &cell,
00586 const Point<dim_+1> &p,
00587 InternalData &mdata,
00588 Point<dim_> &p_unit) const;
00589
00629 Point<dim>
00630 transform_real_to_unit_cell_initial_guess (const std::vector<Point<spacedim> > &vertex,
00631 const Point<spacedim> &p) const;
00632
00633
00634 private:
00663 virtual UpdateFlags update_once (const UpdateFlags flags) const;
00664
00709 virtual UpdateFlags update_each (const UpdateFlags flags) const;
00710
00711 virtual
00712 typename Mapping<dim,spacedim>::InternalDataBase *
00713 get_data (const UpdateFlags,
00714 const Quadrature<dim>& quadrature) const;
00715
00716 virtual
00717 typename Mapping<dim,spacedim>::InternalDataBase *
00718 get_face_data (const UpdateFlags flags,
00719 const Quadrature<dim-1>& quadrature) const;
00720
00721 virtual
00722 typename Mapping<dim,spacedim>::InternalDataBase *
00723 get_subface_data (const UpdateFlags flags,
00724 const Quadrature<dim-1>& quadrature) const;
00725
00742 virtual void compute_mapping_support_points(
00743 const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00744 std::vector<Point<spacedim> > &a) const;
00745
00751 static const unsigned int n_shape_functions = GeometryInfo<dim>::vertices_per_cell;
00752 };
00753
00754
00755 template<>
00756 void
00757 MappingQ1<2,3>::
00758 transform_real_to_unit_cell_internal
00759 (const Triangulation<2,3>::cell_iterator &cell,
00760 const Point<3> &p,
00761 InternalData &mdata,
00762 Point<2> &p_unit) const;
00763
00764
00765 template<>
00766 void
00767 MappingQ1<1,2>::
00768 transform_real_to_unit_cell_internal (const Triangulation<1,2>::cell_iterator &cell,
00769 const Point<2> &p,
00770 InternalData &mdata,
00771 Point<1> &p_unit) const;
00772
00773
00774 template<>
00775 void
00776 MappingQ1<1,3>::
00777 transform_real_to_unit_cell_internal (const Triangulation<1,3>::cell_iterator &cell,
00778 const Point<3> &p,
00779 InternalData &mdata,
00780 Point<1> &p_unit) const;
00781
00782
00789 template <int dim, int spacedim=dim>
00790 struct StaticMappingQ1
00791 {
00792 static MappingQ1<dim, spacedim> mapping;
00793 };
00794
00795
00798
00799
00800 #ifndef DOXYGEN
00801
00802 template<int dim, int spacedim>
00803 inline
00804 double
00805 MappingQ1<dim,spacedim>::InternalData::shape (const unsigned int qpoint,
00806 const unsigned int shape_nr) const
00807 {
00808 Assert(qpoint*n_shape_functions + shape_nr < shape_values.size(),
00809 ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
00810 shape_values.size()));
00811 return shape_values [qpoint*n_shape_functions + shape_nr];
00812 }
00813
00814
00815
00816 template<int dim, int spacedim>
00817 inline
00818 double &
00819 MappingQ1<dim,spacedim>::InternalData::shape (const unsigned int qpoint,
00820 const unsigned int shape_nr)
00821 {
00822 Assert(qpoint*n_shape_functions + shape_nr < shape_values.size(),
00823 ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
00824 shape_values.size()));
00825 return shape_values [qpoint*n_shape_functions + shape_nr];
00826 }
00827
00828
00829 template<int dim, int spacedim>
00830 inline
00831 Tensor<1,dim>
00832 MappingQ1<dim,spacedim>::InternalData::derivative (const unsigned int qpoint,
00833 const unsigned int shape_nr) const
00834 {
00835 Assert(qpoint*n_shape_functions + shape_nr < shape_derivatives.size(),
00836 ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
00837 shape_derivatives.size()));
00838 return shape_derivatives [qpoint*n_shape_functions + shape_nr];
00839 }
00840
00841
00842
00843 template<int dim, int spacedim>
00844 inline
00845 Tensor<1,dim> &
00846 MappingQ1<dim,spacedim>::InternalData::derivative (const unsigned int qpoint,
00847 const unsigned int shape_nr)
00848 {
00849 Assert(qpoint*n_shape_functions + shape_nr < shape_derivatives.size(),
00850 ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
00851 shape_derivatives.size()));
00852 return shape_derivatives [qpoint*n_shape_functions + shape_nr];
00853 }
00854
00855
00856 template <int dim, int spacedim>
00857 inline
00858 Tensor<2,dim>
00859 MappingQ1<dim,spacedim>::InternalData::second_derivative (const unsigned int qpoint,
00860 const unsigned int shape_nr) const
00861 {
00862 Assert(qpoint*n_shape_functions + shape_nr < shape_second_derivatives.size(),
00863 ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
00864 shape_second_derivatives.size()));
00865 return shape_second_derivatives [qpoint*n_shape_functions + shape_nr];
00866 }
00867
00868
00869
00870 template <int dim, int spacedim>
00871 inline
00872 Tensor<2,dim> &
00873 MappingQ1<dim,spacedim>::InternalData::second_derivative (const unsigned int qpoint,
00874 const unsigned int shape_nr)
00875 {
00876 Assert(qpoint*n_shape_functions + shape_nr < shape_second_derivatives.size(),
00877 ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
00878 shape_second_derivatives.size()));
00879 return shape_second_derivatives [qpoint*n_shape_functions + shape_nr];
00880 }
00881
00882
00883
00884 template <int dim, int spacedim>
00885 inline
00886 bool
00887 MappingQ1<dim,spacedim>::preserves_vertex_locations () const
00888 {
00889 return true;
00890 }
00891
00892 #endif // DOXYGEN
00893
00894
00895
00896
00897 DEAL_II_NAMESPACE_CLOSE
00898
00899 #endif
00900