DOLFIN
DOLFIN C++ interface
Loading...
Searching...
No Matches
Cell.h
1// Copyright (C) 2006-2015 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// Modified by Johan Hoffman 2006.
19// Modified by Andre Massing 2009.
20// Modified by Garth N. Wells 2010.
21// Modified by Jan Blechta 2013
22// Modified by Martin Alnaes, 2015
23
24#ifndef __CELL_H
25#define __CELL_H
26
27#include <memory>
28
29#include "CellType.h"
30#include "Mesh.h"
31#include "MeshEntity.h"
32#include "MeshEntityIteratorBase.h"
33#include "MeshFunction.h"
34#include <ufc.h>
35#include <dolfin/geometry/Point.h>
36
37namespace dolfin
38{
39
41
42 class Cell : public MeshEntity
43 {
44 public:
45
48
55 Cell(const Mesh& mesh, std::size_t index)
56 : MeshEntity(mesh, mesh.topology().dim(), index) {}
57
59 ~Cell() {}
60
63 { return _mesh->type().cell_type(); }
64
66 std::size_t num_vertices() const
67 { return _mesh->type().num_vertices(); }
68
73 std::size_t orientation() const
74 { return _mesh->type().orientation(*this); }
75
83 std::size_t orientation(const Point& up) const
84 { return _mesh->type().orientation(*this, up); }
85
98 double volume() const
99 { return _mesh->type().volume(*this); }
100
113 double h() const
114 { return _mesh->type().h(*this); }
115
128 double circumradius() const
129 { return _mesh->type().circumradius(*this); }
130
143 double inradius() const
144 {
145 // We would need facet areas
146 _mesh->init(_mesh->type().dim() - 1);
147
148 return _mesh->type().inradius(*this);
149 }
150
167 double radius_ratio() const
168 {
169 // We would need facet areas
170 _mesh->init(_mesh->type().dim() - 1);
171
172 return _mesh->type().radius_ratio(*this);
173 }
174
181 double squared_distance(const Point& point) const
182 { return _mesh->type().squared_distance(*this, point); }
183
190 double distance(const Point& point) const
191 {
192 return sqrt(squared_distance(point));
193 }
194
204 double normal(std::size_t facet, std::size_t i) const
205 { return _mesh->type().normal(*this, facet, i); }
206
214 Point normal(std::size_t facet) const
215 { return _mesh->type().normal(*this, facet); }
216
222 { return _mesh->type().cell_normal(*this); }
223
231 double facet_area(std::size_t facet) const
232 { return _mesh->type().facet_area(*this, facet); }
233
238 void order(const std::vector<std::int64_t>& local_to_global_vertex_indices)
239 { _mesh->type().order(*this, local_to_global_vertex_indices); }
240
248 bool ordered(const std::vector<std::int64_t>& local_to_global_vertex_indices) const
249 { return _mesh->type().ordered(*this, local_to_global_vertex_indices); }
250
259 bool contains(const Point& point) const;
260
268 bool collides(const Point& point) const;
269
277 bool collides(const MeshEntity& entity) const;
278
288 std::vector<Point>
289 intersection(const MeshEntity& entity) const;
290
291 // FIXME: This function is part of a UFC transition
293 void get_coordinate_dofs(std::vector<double>& coordinates) const
294 {
295 const MeshGeometry& geom = _mesh->geometry();
296 const std::size_t gdim = geom.dim();
297 const std::size_t geom_degree = geom.degree();
298 const std::size_t num_vertices = this->num_vertices();
299 const unsigned int* vertices = this->entities(0);
300
301 if (geom_degree == 1)
302 {
303 coordinates.resize(num_vertices*gdim);
304 for (std::size_t i = 0; i < num_vertices; ++i)
305 for (std::size_t j = 0; j < gdim; ++j)
306 coordinates[i*gdim + j] = geom.x(vertices[i])[j];
307 }
308 else if (geom_degree == 2)
309 {
310 const std::size_t tdim = _mesh->topology().dim();
311 const std::size_t num_edges = this->num_entities(1);
312 const unsigned int* edges = this->entities(1);
313
314 coordinates.resize((num_vertices + num_edges)*gdim);
315
316 for (std::size_t i = 0; i < num_vertices; ++i)
317 for (std::size_t j = 0; j < gdim; j++)
318 coordinates[i*gdim + j] = geom.x(vertices[i])[j];
319
320 for (std::size_t i = 0; i < num_edges; ++i)
321 {
322 const std::size_t entity_index
323 = (tdim == 1) ? index() : edges[i];
324 const std::size_t point_index
325 = geom.get_entity_index(1, 0, entity_index);
326 for (std::size_t j = 0; j < gdim; ++j)
327 coordinates[(i + num_vertices)*gdim + j] = geom.x(point_index)[j];
328 }
329 }
330 else
331 {
332 dolfin_error("Cell.h", "get coordinate_dofs", "Unsupported mesh degree");
333 }
334
335 }
336
337 // FIXME: This function is part of a UFC transition
339 void get_vertex_coordinates(std::vector<double>& coordinates) const
340 {
341 const std::size_t gdim = _mesh->geometry().dim();
342 const std::size_t num_vertices = this->num_vertices();
343 const unsigned int* vertices = this->entities(0);
344 coordinates.resize(num_vertices*gdim);
345 for (std::size_t i = 0; i < num_vertices; i++)
346 for (std::size_t j = 0; j < gdim; j++)
347 coordinates[i*gdim + j] = _mesh->geometry().x(vertices[i])[j];
348 }
349
350 // FIXME: This function is part of a UFC transition
352 void get_cell_data(ufc::cell& ufc_cell, int local_facet=-1) const
353 {
354 ufc_cell.geometric_dimension = _mesh->geometry().dim();
355 ufc_cell.local_facet = local_facet;
356 if (_mesh->cell_orientations().empty())
357 ufc_cell.orientation = -1;
358 else
359 {
360 dolfin_assert(index() < _mesh->cell_orientations().size());
361 ufc_cell.orientation = _mesh->cell_orientations()[index()];
362 }
363 ufc_cell.mesh_identifier = mesh_id();
364 ufc_cell.index = index();
365 }
366
367 // FIXME: This function is part of a UFC transition
369 void get_cell_topology(ufc::cell& ufc_cell) const
370 {
371 const MeshTopology& topology = _mesh->topology();
372
373 const std::size_t tdim = topology.dim();
374 ufc_cell.topological_dimension = tdim;
375 if (_mesh->cell_orientations().empty())
376 ufc_cell.orientation = -1;
377 else
378 {
379 dolfin_assert(index() < _mesh->cell_orientations().size());
380 ufc_cell.orientation = _mesh->cell_orientations()[index()];
381 }
382 ufc_cell.entity_indices.resize(tdim + 1);
383 for (std::size_t d = 0; d < tdim; d++)
384 {
385 ufc_cell.entity_indices[d].resize(num_entities(d));
386 if (topology.have_global_indices(d))
387 {
388 const std::vector<std::int64_t>& global_indices
389 = topology.global_indices(d);
390 for (std::size_t i = 0; i < num_entities(d); ++i)
391 ufc_cell.entity_indices[d][i] = global_indices[entities(d)[i]];
392 }
393 else
394 {
395 for (std::size_t i = 0; i < num_entities(d); ++i)
396 ufc_cell.entity_indices[d][i] = entities(d)[i];
397 }
398 }
399 ufc_cell.entity_indices[tdim].resize(1);
400 if (topology.have_global_indices(tdim))
401 ufc_cell.entity_indices[tdim][0] = global_index();
402 else
403 ufc_cell.entity_indices[tdim][0] = index();
404
405 // FIXME: Using the local cell index is inconsistent with UFC, but
406 // necessary to make DOLFIN run
407 // Local cell index
408 ufc_cell.index = ufc_cell.entity_indices[tdim][0];
409 }
410
411 };
412
415
416}
417
418#endif
virtual double volume(const MeshEntity &entity) const =0
Compute (generalized) volume of mesh entity.
virtual double facet_area(const Cell &cell, std::size_t facet) const =0
Compute the area/length of given facet with respect to the cell.
std::size_t num_vertices() const
Return number of vertices for cell.
Definition CellType.h:92
virtual double h(const MeshEntity &entity) const
Compute greatest distance between any two vertices.
Definition CellType.cpp:151
virtual double inradius(const Cell &cell) const
Compute inradius of cell.
Definition CellType.cpp:178
virtual Point cell_normal(const Cell &cell) const =0
Compute normal to given cell (viewed as embedded in 3D)
virtual double radius_ratio(const Cell &cell) const
Compute dim*inradius/circumradius for given cell.
Definition CellType.cpp:210
virtual void order(Cell &cell, const std::vector< std::int64_t > &local_to_global_vertex_indices) const =0
Order entities locally.
virtual std::size_t dim() const =0
Return topological dimension of cell.
bool ordered(const Cell &cell, const std::vector< std::int64_t > &local_to_global_vertex_indices) const
Check if entities are ordered.
Definition CellType.cpp:221
virtual double circumradius(const MeshEntity &entity) const =0
Compute circumradius of mesh entity.
Type cell_type() const
Return type of cell.
Definition CellType.h:72
virtual double squared_distance(const Cell &cell, const Point &point) const =0
Compute squared distance to given point.
virtual std::size_t orientation(const Cell &cell) const =0
Return orientation of the cell (assuming flat space)
Type
Enum for different cell types.
Definition CellType.h:51
virtual double normal(const Cell &cell, std::size_t facet, std::size_t i) const =0
Compute component i of normal of given facet with respect to the cell.
A Cell is a MeshEntity of topological codimension 0.
Definition Cell.h:43
double h() const
Definition Cell.h:113
double inradius() const
Definition Cell.h:143
double facet_area(std::size_t facet) const
Definition Cell.h:231
double distance(const Point &point) const
Definition Cell.h:190
double volume() const
Definition Cell.h:98
bool contains(const Point &point) const
Definition Cell.cpp:25
double radius_ratio() const
Definition Cell.h:167
void get_coordinate_dofs(std::vector< double > &coordinates) const
Get cell coordinate dofs (not vertex coordinates)
Definition Cell.h:293
void get_cell_topology(ufc::cell &ufc_cell) const
Fill UFC cell with topology data.
Definition Cell.h:369
bool ordered(const std::vector< std::int64_t > &local_to_global_vertex_indices) const
Definition Cell.h:248
void get_vertex_coordinates(std::vector< double > &coordinates) const
Get cell vertex coordinates (not coordinate dofs)
Definition Cell.h:339
~Cell()
Destructor.
Definition Cell.h:59
Cell(const Mesh &mesh, std::size_t index)
Definition Cell.h:55
double squared_distance(const Point &point) const
Definition Cell.h:181
Cell()
Create empty cell.
Definition Cell.h:47
Point normal(std::size_t facet) const
Definition Cell.h:214
std::vector< Point > intersection(const MeshEntity &entity) const
Definition Cell.cpp:41
Point cell_normal() const
Definition Cell.h:221
double circumradius() const
Definition Cell.h:128
double normal(std::size_t facet, std::size_t i) const
Definition Cell.h:204
bool collides(const Point &point) const
Definition Cell.cpp:30
CellType::Type type() const
Return type of cell.
Definition Cell.h:62
void get_cell_data(ufc::cell &ufc_cell, int local_facet=-1) const
Fill UFC cell with miscellaneous data.
Definition Cell.h:352
std::size_t orientation() const
Definition Cell.h:73
std::size_t orientation(const Point &up) const
Definition Cell.h:83
std::size_t num_vertices() const
Return number of vertices of cell.
Definition Cell.h:66
void order(const std::vector< std::int64_t > &local_to_global_vertex_indices)
Definition Cell.h:238
Base class for MeshEntityIterators.
Definition MeshEntityIteratorBase.h:37
Definition MeshEntity.h:43
const Mesh & mesh() const
Definition MeshEntity.h:99
const unsigned int * entities(std::size_t dim) const
Definition MeshEntity.h:163
std::size_t dim() const
Definition MeshEntity.h:106
std::size_t index() const
Definition MeshEntity.h:113
std::size_t num_entities(std::size_t dim) const
Definition MeshEntity.h:140
std::size_t mesh_id() const
Definition MeshEntity.h:175
std::int64_t global_index() const
Definition MeshEntity.h:122
MeshGeometry stores the geometry imposed on a mesh.
Definition MeshGeometry.h:40
std::size_t get_entity_index(std::size_t entity_dim, std::size_t order, std::size_t index) const
Get the index for an entity point in coordinates.
Definition MeshGeometry.h:153
std::size_t dim() const
Return Euclidean dimension of coordinate system.
Definition MeshGeometry.h:56
std::size_t degree() const
Return polynomial degree of coordinate field.
Definition MeshGeometry.h:60
double x(std::size_t n, std::size_t i) const
Return value of coordinate with local index n in direction i.
Definition MeshGeometry.h:99
Definition MeshTopology.h:53
std::size_t dim() const
Return topological dimension.
Definition MeshTopology.cpp:74
bool have_global_indices(std::size_t dim) const
Definition MeshTopology.h:121
const std::vector< std::int64_t > & global_indices(std::size_t d) const
Definition MeshTopology.h:113
Definition Mesh.h:84
std::size_t init(std::size_t dim) const
Definition Mesh.cpp:139
MeshTopology & topology()
Definition Mesh.h:220
CellType & type()
Definition Mesh.h:284
MeshGeometry & geometry()
Definition Mesh.h:234
const std::vector< int > & cell_orientations() const
Definition Mesh.cpp:437
Definition Point.h:41
Definition adapt.h:30
void dolfin_error(std::string location, std::string task, std::string reason,...)
Definition log.cpp:129
MeshEntityIteratorBase< Cell > CellIterator
A CellIterator is a MeshEntityIterator of topological codimension 0.
Definition Cell.h:414