From d3ccfa93e6642cbb2c2f77ba47cf4745a4864968 Mon Sep 17 00:00:00 2001 From: kanschat Date: Mon, 5 Jan 2009 16:52:48 +0000 Subject: [PATCH] introduce function transform() into Mappings git-svn-id: https://svn.dealii.org/trunk@18075 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/mapping.h | 60 ++++++++++++++++++- .../deal.II/include/fe/mapping_cartesian.h | 14 ++++- deal.II/deal.II/include/fe/mapping_q.h | 14 ++++- deal.II/deal.II/include/fe/mapping_q1.h | 14 ++++- .../deal.II/source/fe/mapping_cartesian.cc | 46 +++++++++++++- deal.II/deal.II/source/fe/mapping_q.cc | 46 +++++++++++++- deal.II/deal.II/source/fe/mapping_q1.cc | 45 +++++++++++++- 7 files changed, 232 insertions(+), 7 deletions(-) diff --git a/deal.II/deal.II/include/fe/mapping.h b/deal.II/deal.II/include/fe/mapping.h index 6ce98fb94a..fe9ad1f52e 100644 --- a/deal.II/deal.II/include/fe/mapping.h +++ b/deal.II/deal.II/include/fe/mapping.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors +// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -34,6 +34,37 @@ template class FESubfaceValues; //TODO: Offset in transform functions should be replaced by initializing VectorSlice correctly + /** + * The transformation type used + * for the Mapping::transform() functions. + * + * Special finite elements may + * need special Mapping from the + * reference cell to the actual + * mesh cell. In order to be most + * flexible, this enum provides + * an extensible interface for + * arbitrary + * transformations. Nevertheless, + * these must be implemented in + * the transform() functions of + * inheriting classes in order to + * work. + */ +enum MappingType +{ +/// Covariant mapping + mapping_covariant = 0x0001, +/// Contravariant mapping + mapping_contravariant = 0x0002, +/// The Piola transform usually used for Hdiv elements + mapping_piola = 0x0003, +/// The mapping used for Raviart-Thomas elements + mapping_raviart_thomas = mapping_piola, +/// The mapping used for BDM elements + mapping_bdm = mapping_piola +}; + /** * Abstract base class for mapping classes. @@ -256,6 +287,31 @@ class Mapping : public Subscriptor bool first_cell; }; + /** + * Transform a field of + * vectors accorsing to + * the selected MappingType. + */ + virtual + void + transform (const VectorSlice > > input, + VectorSlice > > output, + const InternalDataBase &internal, + const MappingType type) const = 0; + + + /** + * Transform a field of + * rank two tensors accorsing to + * the selected MappingType. + */ + virtual + void + transform (const VectorSlice > > input, + VectorSlice > > output, + const InternalDataBase &internal, + const MappingType type) const = 0; + /** * Transform a field of covariant * vectors. The covariant @@ -307,6 +363,8 @@ class Mapping : public Subscriptor const unsigned int offset, VectorSlice > > output, const InternalDataBase &internal) const = 0; + +//TODO: The argument list here seems wrong /** * Transform a field of diff --git a/deal.II/deal.II/include/fe/mapping_cartesian.h b/deal.II/deal.II/include/fe/mapping_cartesian.h index 6656512626..34c14b6b2a 100644 --- a/deal.II/deal.II/include/fe/mapping_cartesian.h +++ b/deal.II/deal.II/include/fe/mapping_cartesian.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -92,6 +92,18 @@ class MappingCartesian : public Mapping std::vector > &normal_vectors, std::vector &cell_JxW_values) const ; + virtual void + transform (const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType type) const; + + virtual void + transform (const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType type) const; + virtual void transform_covariant (const VectorSlice > > input, diff --git a/deal.II/deal.II/include/fe/mapping_q.h b/deal.II/deal.II/include/fe/mapping_q.h index 0ae4531ef9..6995b18a83 100644 --- a/deal.II/deal.II/include/fe/mapping_q.h +++ b/deal.II/deal.II/include/fe/mapping_q.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -116,6 +116,18 @@ class MappingQ : public MappingQ1 const typename Triangulation::cell_iterator &cell, const Point &p) const; + virtual void + transform (const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType type) const; + + virtual void + transform (const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType type) const; + /** * Implementation of the * interface in Mapping. diff --git a/deal.II/deal.II/include/fe/mapping_q1.h b/deal.II/deal.II/include/fe/mapping_q1.h index 7935b2854c..73c8b50321 100644 --- a/deal.II/deal.II/include/fe/mapping_q1.h +++ b/deal.II/deal.II/include/fe/mapping_q1.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -77,6 +77,18 @@ class MappingQ1 : public Mapping const typename Triangulation::cell_iterator &cell, const Point &p) const; + virtual void + transform (const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType type) const; + + virtual void + transform (const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType type) const; + virtual void transform_covariant (const VectorSlice > > input, const unsigned int offset, diff --git a/deal.II/deal.II/source/fe/mapping_cartesian.cc b/deal.II/deal.II/source/fe/mapping_cartesian.cc index 604f819c03..ffb5c9ab23 100644 --- a/deal.II/deal.II/source/fe/mapping_cartesian.cc +++ b/deal.II/deal.II/source/fe/mapping_cartesian.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -504,6 +504,50 @@ MappingCartesian<1>::fill_fe_subface_values (const Triangulation<1>::cell_iterat #endif +template +void +MappingCartesian::transform ( + const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType mapping_type) const +{ + switch (mapping_type) + { + case mapping_covariant: + transform_covariant(input, 0, output, internal); + return; +// case mapping_contravariant: +// transform_contravariant(input, 0, output, internal); +// return; + default: + Assert(false, ExcNotImplemented()); + } +} + +template +void +MappingCartesian::transform ( + const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType mapping_type) const +{ + switch (mapping_type) + { + case mapping_covariant: + transform_covariant(input, 0, output, internal); + return; +// case mapping_contravariant: +// transform_contravariant(input, 0, output, internal); +// return; + default: + Assert(false, ExcNotImplemented()); + } +} + + + template void MappingCartesian::transform_covariant ( diff --git a/deal.II/deal.II/source/fe/mapping_q.cc b/deal.II/deal.II/source/fe/mapping_q.cc index 880139da6a..48d14df0e2 100644 --- a/deal.II/deal.II/source/fe/mapping_q.cc +++ b/deal.II/deal.II/source/fe/mapping_q.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -1201,6 +1201,50 @@ add_quad_support_points(const typename Triangulation::cell_iterato } +template +void +MappingQ::transform ( + const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType mapping_type) const +{ + switch (mapping_type) + { + case mapping_covariant: + transform_covariant(input, 0, output, internal); + return; +// case mapping_contravariant: +// transform_contravariant(input, 0, output, internal); +// return; + default: + Assert(false, ExcNotImplemented()); + } +} + + +template +void +MappingQ::transform ( + const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType mapping_type) const +{ + switch (mapping_type) + { + case mapping_covariant: + transform_covariant(input, 0, output, internal); + return; +// case mapping_contravariant: +// transform_contravariant(input, 0, output, internal); +// return; + default: + Assert(false, ExcNotImplemented()); + } +} + + template void diff --git a/deal.II/deal.II/source/fe/mapping_q1.cc b/deal.II/deal.II/source/fe/mapping_q1.cc index e116a9311d..17f59a091e 100644 --- a/deal.II/deal.II/source/fe/mapping_q1.cc +++ b/deal.II/deal.II/source/fe/mapping_q1.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -1165,6 +1165,49 @@ MappingQ1<1,2>::fill_fe_subface_values (const Triangulation<1,2>::cell_iterator #endif +template +void +MappingQ1::transform ( + const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType mapping_type) const +{ + switch (mapping_type) + { + case mapping_covariant: + transform_covariant(input, 0, output, internal); + return; +// case mapping_contravariant: +// transform_contravariant(input, 0, output, internal); +// return; + default: + Assert(false, ExcNotImplemented()); + } +} + +template +void +MappingQ1::transform ( + const VectorSlice > > input, + VectorSlice > > output, + const typename Mapping::InternalDataBase &internal, + const MappingType mapping_type) const +{ + switch (mapping_type) + { + case mapping_covariant: + transform_covariant(input, 0, output, internal); + return; +// case mapping_contravariant: +// transform_contravariant(input, 0, output, internal); +// return; + default: + Assert(false, ExcNotImplemented()); + } +} + + template void MappingQ1::transform_covariant ( -- 2.39.5