4
4
// Modified by Kristian Oelgaard, 2007
6
6
// First added: 2007-04-10
7
// Last changed: 2008-05-18
7
// Last changed: 2008-05-21
9
9
#include <dolfin/common/constants.h>
10
10
#include <dolfin/log/log.h>
11
11
#include <dolfin/mesh/Mesh.h>
12
#include <dolfin/mesh/MeshData.h>
12
13
#include <dolfin/mesh/Vertex.h>
13
14
#include <dolfin/mesh/Cell.h>
14
15
#include <dolfin/mesh/Facet.h>
32
33
: BoundaryCondition(), g(g), _mesh(mesh),
33
34
sub_domains(0), sub_domain(0), sub_domains_local(false),
34
35
method(method), user_sub_domain(&sub_domain)
37
37
// Initialize sub domain markers
38
initFromSubDomain(sub_domain);
40
40
//-----------------------------------------------------------------------------
41
41
DirichletBC::DirichletBC(Function& g,
51
51
//-----------------------------------------------------------------------------
52
52
DirichletBC::DirichletBC(Function& g,
56
: BoundaryCondition(), g(g), _mesh(mesh),
57
sub_domains(0), sub_domain(sub_domain), sub_domains_local(false),
58
method(method), user_sub_domain(0)
60
// Initialize sub domain markers
63
//-----------------------------------------------------------------------------
64
DirichletBC::DirichletBC(Function& g,
54
66
SubDomain& sub_domain,
55
67
const SubSystem& sub_system,
58
70
sub_domains(0), sub_domain(0), sub_domains_local(false),
59
71
sub_system(sub_system), method(method), user_sub_domain(&sub_domain)
61
// Set sub domain markers
73
// Initialize sub domain markers
74
initFromSubDomain(sub_domain);
64
76
//-----------------------------------------------------------------------------
65
77
DirichletBC::DirichletBC(Function& g,
76
88
//-----------------------------------------------------------------------------
77
89
DirichletBC::DirichletBC(Function& g,
92
const SubSystem& sub_system,
80
94
: BoundaryCondition(), g(g), _mesh(mesh),
81
sub_domains(0), sub_domain(0), sub_domains_local(false),
82
method(method), user_sub_domain(0)
95
sub_domains(), sub_domain(sub_domain), sub_domains_local(false),
96
sub_system(sub_system), method(method), user_sub_domain(0)
84
// Create sub domain for entire boundary
85
class EntireBoundary : public SubDomain
88
bool inside(const real* x, bool on_boundary) const
94
EntireBoundary sub_domain;
96
98
// Initialize sub domain markers
99
101
//-----------------------------------------------------------------------------
100
102
DirichletBC::~DirichletBC()
174
176
Progress p(s, _mesh.size(D - 1));
175
177
for (FacetIterator facet(_mesh); !facet.end(); ++facet)
177
// Skip facets not inside the sub domain
179
// Skip facets not inside the sub domain
178
180
if ((*sub_domains)(*facet) != sub_domain)
238
240
//-----------------------------------------------------------------------------
239
void DirichletBC::init(SubDomain& sub_domain)
241
void DirichletBC::initFromSubDomain(SubDomain& sub_domain)
241
243
cout << "Creating sub domain markers for boundary condition." << endl;
252
254
sub_domain.mark(*sub_domains, 0);
254
256
//-----------------------------------------------------------------------------
257
void DirichletBC::initFromMesh()
259
cout << "Creating sub domain markers for boundary condition." << endl;
262
Array<uint>* facet_cells = _mesh.data().array("boundary facet cells");
263
Array<uint>* facet_numbers = _mesh.data().array("boundary facet numbers");
264
Array<uint>* indicators = _mesh.data().array("boundary indicators");
268
error("Mesh data \"boundary facet cells\" not available.");
270
error("Mesh data \"boundary facet numbers\" not available.");
272
error("Mesh data \"boundary indicators\" not available.");
275
const uint size = facet_cells->size();
276
dolfin_assert(size == facet_numbers->size());
277
dolfin_assert(size == indicators->size());
279
// Create mesh function for sub domain markers on facets
280
const uint dim = _mesh.topology().dim();
282
sub_domains = new MeshFunction<uint>(_mesh, dim - 1);
283
sub_domains_local = true;
285
// Compute the maximum boundary indicator
287
for (uint i = 0; i < size; i++)
288
maxid = std::max(maxid, (*indicators)[i]);
289
cout << "maxid = " << maxid << endl;
291
// Mark everything with the maximum value + 1
292
(*sub_domains) = maxid + 1;
294
// Mark facets according to data
295
for (uint i = 0; i < size; i++)
297
// Get cell incident with facet
298
Cell cell(_mesh, (*facet_cells)[i]);
301
uint facet = cell.entities(dim - 1)[(*facet_numbers)[i]];
304
sub_domains->set(facet, (*indicators)[i]);
307
//-----------------------------------------------------------------------------
255
308
void DirichletBC::computeBCTopological(std::map<uint, real>& boundary_values,
257
310
BoundaryCondition::LocalData& data)
470
523
std::map<uint, real>::const_iterator boundary_value;
472
525
for (boundary_value = boundary_values.begin(); boundary_value != boundary_values.end(); ++boundary_value)
474
dofs[i++] = boundary_value->first;
526
dofs[i++] = boundary_value->first;
477
528
// Modify linear system (A_ii = 1)
478
529
A.zero(boundary_values.size(), dofs);