DOLFIN
DOLFIN C++ interface
Loading...
Searching...
No Matches
SubDomain.h
1// Copyright (C) 2007-2013 Anders Logg
2//
3// This file is part of DOLFIN.
4//
5// DOLFIN is free software: you can redistribute it and/or modify
6// it under the terms of the GNU Lesser General Public License as published by
7// the Free Software Foundation, either version 3 of the License, or
8// (at your option) any later version.
9//
10// DOLFIN is distributed in the hope that it will be useful,
11// but WITHOUT ANY WARRANTY; without even the implied warranty of
12// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13// GNU Lesser General Public License for more details.
14//
15// You should have received a copy of the GNU Lesser General Public License
16// along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17//
18// First added: 2007-04-10
19// Last changed: 2013-04-12
20
21#ifndef __SUB_DOMAIN_H
22#define __SUB_DOMAIN_H
23
24#include <cstddef>
25#include <map>
26#include <dolfin/common/constants.h>
27#include <Eigen/Dense>
28
29namespace dolfin
30{
31
32 // Forward declarations
33 class Mesh;
34 template <typename T> class MeshFunction;
35 template <typename T> class MeshValueCollection;
36 template <typename T> class Array;
37
41
43 {
44 public:
45
51 SubDomain(const double map_tol=1.0e-10);
52
54 virtual ~SubDomain();
55
65 virtual bool inside(const Array<double>& x, bool on_boundary) const;
66
76 virtual bool inside(Eigen::Ref<const Eigen::VectorXd> x, bool on_boundary) const;
77
85 virtual void map(const Array<double>& x, Array<double>& y) const;
86
87
95 virtual void map(Eigen::Ref<const Eigen::VectorXd> x, Eigen::Ref<Eigen::VectorXd> y) const;
96
97
102 virtual void snap(Array<double>& x) const;
103
108 virtual void snap(Eigen::Ref<Eigen::VectorXd> x) const;
109
110 //--- Marking of Mesh ---
111
120 void mark_cells(Mesh& mesh,
121 std::size_t sub_domain,
122 bool check_midpoint=true) const;
123
132 void mark_facets(Mesh& mesh,
133 std::size_t sub_domain,
134 bool check_midpoint=true) const;
135
147 void mark(Mesh& mesh,
148 std::size_t dim,
149 std::size_t sub_domain,
150 bool check_midpoint=true) const;
151
152 //--- Marking of MeshFunction ---
153
162 void mark(MeshFunction<std::size_t>& sub_domains,
163 std::size_t sub_domain,
164 bool check_midpoint=true) const;
165
174 void mark(MeshFunction<int>& sub_domains,
175 int sub_domain,
176 bool check_midpoint=true) const;
177
186 void mark(MeshFunction<double>& sub_domains,
187 double sub_domain,
188 bool check_midpoint=true) const;
189
198 void mark(MeshFunction<bool>& sub_domains,
199 bool sub_domain,
200 bool check_midpoint=true) const;
201
202 //--- Marking of MeshValueCollection ---
203
214 void mark(MeshValueCollection<std::size_t>& sub_domains,
215 std::size_t sub_domain,
216 const Mesh& mesh,
217 bool check_midpoint=true) const;
218
229 void mark(MeshValueCollection<int>& sub_domains,
230 int sub_domain,
231 const Mesh& mesh,
232 bool check_midpoint=true) const;
233
244 void mark(MeshValueCollection<double>& sub_domains,
245 double sub_domain,
246 const Mesh& mesh,
247 bool check_midpoint=true) const;
248
259 void mark(MeshValueCollection<bool>& sub_domains,
260 bool sub_domain,
261 const Mesh& mesh,
262 bool check_midpoint=true) const;
263
268 std::size_t geometric_dimension() const;
269
274 virtual void set_property(std::string name, double value);
275
280 virtual double get_property(std::string name) const;
281
286 const double map_tolerance;
287
288 private:
289
292 template<typename S, typename T>
293 void apply_markers(S& sub_domains,
294 T sub_domain,
295 const Mesh& mesh,
296 bool check_midpoint) const;
297
298 template<typename T>
299 void apply_markers(std::map<std::size_t, std::size_t>& sub_domains,
300 std::size_t dim,
301 T sub_domain,
302 const Mesh& mesh,
303 bool check_midpoint) const;
304
305 // Friends
306 friend class DirichletBC;
307 friend class PeriodicBC;
308
309 // Geometric dimension, needed for SWIG interface, will be set before
310 // calls to inside() and map()
311 mutable std::size_t _geometric_dimension;
312
313 };
314
315}
316
317#endif
Definition SubDomain.h:36
Interface for setting (strong) Dirichlet boundary conditions.
Definition DirichletBC.h:125
Definition RegularCutRefinement.h:31
Definition SubDomain.h:35
Definition Mesh.h:84
Definition SubDomain.h:43
virtual double get_property(std::string name) const
Definition SubDomain.cpp:401
SubDomain(const double map_tol=1.0e-10)
Definition SubDomain.cpp:40
void mark(Mesh &mesh, std::size_t dim, std::size_t sub_domain, bool check_midpoint=true) const
Definition SubDomain.cpp:107
virtual ~SubDomain()
Destructor.
Definition SubDomain.cpp:46
virtual bool inside(const Array< double > &x, bool on_boundary) const
Definition SubDomain.cpp:51
virtual void set_property(std::string name, double value)
Definition SubDomain.cpp:394
std::size_t geometric_dimension() const
Definition SubDomain.cpp:178
const double map_tolerance
Definition SubDomain.h:286
void mark_facets(Mesh &mesh, std::size_t sub_domain, bool check_midpoint=true) const
Definition SubDomain.cpp:100
virtual void snap(Array< double > &x) const
Definition SubDomain.cpp:80
void mark_cells(Mesh &mesh, std::size_t sub_domain, bool check_midpoint=true) const
Definition SubDomain.cpp:93
virtual void map(const Array< double > &x, Array< double > &y) const
Definition SubDomain.cpp:65
Definition adapt.h:30