Chaste Commit::30a3e656d4b131f8c595cc6eb2becd297337570f
AbstractBoxDomainPdeModifier.cpp
1/*
2
3Copyright (c) 2005-2025, University of Oxford.
4All rights reserved.
5
6University of Oxford means the Chancellor, Masters and Scholars of the
7University of Oxford, having an administrative office at Wellington
8Square, Oxford OX1 2JD, UK.
9
10This file is part of Chaste.
11
12Redistribution and use in source and binary forms, with or without
13modification, are permitted provided that the following conditions are met:
14 * Redistributions of source code must retain the above copyright notice,
15 this list of conditions and the following disclaimer.
16 * Redistributions in binary form must reproduce the above copyright notice,
17 this list of conditions and the following disclaimer in the documentation
18 and/or other materials provided with the distribution.
19 * Neither the name of the University of Oxford nor the names of its
20 contributors may be used to endorse or promote products derived from this
21 software without specific prior written permission.
22
23THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
24AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
25IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
26ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
27LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
28CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
29GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
32OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33
34*/
35
36#include "AbstractBoxDomainPdeModifier.hpp"
37#include "ReplicatableVector.hpp"
38#include "LinearBasisFunction.hpp"
39
40template<unsigned DIM>
42 boost::shared_ptr<AbstractBoundaryCondition<DIM> > pBoundaryCondition,
43 bool isNeumannBoundaryCondition,
44 boost::shared_ptr<ChasteCuboid<DIM> > pMeshCuboid,
45 double stepSize,
46 Vec solution)
47 : AbstractPdeModifier<DIM>(pPde,
48 pBoundaryCondition,
49 isNeumannBoundaryCondition,
50 solution),
51 mpMeshCuboid(pMeshCuboid),
52 mStepSize(stepSize),
53 mSetBcsOnBoxBoundary(true),
54 mSetBcsOnBoundingSphere(false),
55 mUseVoronoiCellsForInterpolation(false),
56 mTypicalCellRadius(0.5) // defaults to 0.5
57{
58 if (pMeshCuboid)
59 {
60 // We only need to generate mpFeMesh once, as it does not vary with time
62 this->mDeleteFeMesh = true;
63 // Initialise the boundary nodes
64 this->mIsDirichletBoundaryNode = std::vector<double>(this->mpFeMesh->GetNumNodes(), 0.0);
65 }
66}
67
68template<unsigned DIM>
70{
71 return mStepSize;
72}
73
74template<unsigned DIM>
76{
77 mSetBcsOnBoxBoundary = setBcsOnBoxBoundary;
78}
79
80template<unsigned DIM>
82{
83 return mSetBcsOnBoxBoundary;
84}
85
86template<unsigned DIM>
87void AbstractBoxDomainPdeModifier<DIM>::SetBcsOnBoundingSphere(const bool setBcsOnBoundingSphere)
88{
89 mSetBcsOnBoundingSphere = setBcsOnBoundingSphere;
90}
91
92template<unsigned DIM>
94{
95 return mSetBcsOnBoundingSphere;
96}
97
98template<unsigned DIM>
99void AbstractBoxDomainPdeModifier<DIM>::SetUseVoronoiCellsForInterpolation(const bool useVoronoiCellsForInterpolation)
100{
101 mUseVoronoiCellsForInterpolation = useVoronoiCellsForInterpolation;
102}
103
104template<unsigned DIM>
106{
107 return mUseVoronoiCellsForInterpolation;
108}
109
110template<unsigned DIM>
112{
113 assert(mTypicalCellRadius>=0.0);
114 mTypicalCellRadius = typicalCellRadius;
115}
116
117template<unsigned DIM>
119{
120 return mTypicalCellRadius;
121}
122
123template<unsigned DIM>
125 std::shared_ptr<BoundaryConditionsContainer<DIM,DIM,1> > pBcc)
126{
127 if (!this->mSetBcsOnBoxBoundary)
128 {
129 // Here we approximate the cell population by the bounding sphere and apply the boundary conditions outside the sphere.
130 if (this->mSetBcsOnBoundingSphere)
131 {
132 // First find the centre of the tissue by choosing the midpoint of the extrema.
133 c_vector<double, DIM> tissue_maxima = zero_vector<double>(DIM);
134 c_vector<double, DIM> tissue_minima = zero_vector<double>(DIM);
135 for (unsigned i = 0; i < DIM; i++)
136 {
137 tissue_maxima[i] = -DBL_MAX;
138 tissue_minima[i] = DBL_MAX;
139 }
140
141
142 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
143 cell_iter != rCellPopulation.End();
144 ++cell_iter)
145 {
146 const ChastePoint<DIM>& r_position_of_cell = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
147
148 for (unsigned i = 0; i < DIM; i++)
149 {
150 if (r_position_of_cell[i] > tissue_maxima[i])
151 {
152 tissue_maxima[i] = r_position_of_cell[i];
153 }
154 if (r_position_of_cell[i] < tissue_minima[i])
155 {
156 tissue_minima[i] = r_position_of_cell[i];
157 }
158 }
159 }
160
161 c_vector<double, DIM> tissue_centre = 0.5*(tissue_maxima + tissue_minima);
162
163
164 double tissue_radius = 0.0;
165 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
166 cell_iter != rCellPopulation.End();
167 ++cell_iter)
168 {
169 const ChastePoint<DIM>& r_position_of_cell = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
170
171 double radius = norm_2(tissue_centre - r_position_of_cell.rGetLocation());
172
173 if (tissue_radius < radius)
174 {
175 tissue_radius = radius;
176 }
177 }
178
179 // Apply boundary condition to the nodes outside the tissue_radius
180 if (this->IsNeumannBoundaryCondition())
181 {
183 }
184 else
185 {
186 // Impose any Dirichlet boundary conditions
187 for (unsigned i=0; i<this->mpFeMesh->GetNumNodes(); i++)
188 {
189 double radius = norm_2(tissue_centre - this->mpFeMesh->GetNode(i)->rGetLocation());
190
191 if (radius > tissue_radius)
192 {
193 pBcc->AddDirichletBoundaryCondition(this->mpFeMesh->GetNode(i), this->mpBoundaryCondition.get(), 0, false);
194 this->mIsDirichletBoundaryNode[i] = 1.0;
195 }
196 }
197 }
198 }
199 else // Set pde nodes as boundary node if elements don't contain cells or nodes aren't within 0.5CD of a cell centre
200 {
201 // Get the set of coarse element indices that contain cells
202 std::set<unsigned> coarse_element_indices_in_map;
203 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
204 cell_iter != rCellPopulation.End();
205 ++cell_iter)
206 {
207 coarse_element_indices_in_map.insert(this->mCellPdeElementMap[*cell_iter]);
208 }
209
210 // Find the node indices associated with elements whose indices are NOT in the set coarse_element_indices_in_map
211 std::set<unsigned> coarse_mesh_boundary_node_indices;
212 for (unsigned i=0; i<this->mpFeMesh->GetNumElements(); i++)
213 {
214 if (coarse_element_indices_in_map.find(i) == coarse_element_indices_in_map.end())
215 {
216 Element<DIM,DIM>* p_element = this->mpFeMesh->GetElement(i);
217 for (unsigned j=0; j<DIM+1; j++)
218 {
219 unsigned node_index = p_element->GetNodeGlobalIndex(j);
220 coarse_mesh_boundary_node_indices.insert(node_index);
221 }
222 }
223 }
224
225 // Also remove nodes that are within the typical cell radius from the centre of a cell.
226 std::set<unsigned> nearby_node_indices;
227 for (unsigned int coarse_mesh_boundary_node_idx : coarse_mesh_boundary_node_indices)
228 {
229 bool remove_node = false;
230
231 c_vector<double,DIM> node_location = this->mpFeMesh->GetNode(coarse_mesh_boundary_node_idx)->rGetLocation();
232
233 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
234 cell_iter != rCellPopulation.End();
235 ++cell_iter)
236 {
237 const ChastePoint<DIM>& r_position_of_cell = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
238
239 double separation = norm_2(node_location - r_position_of_cell.rGetLocation());
240
241 if (separation <= mTypicalCellRadius)
242 {
243 remove_node = true;
244 break;
245 }
246 }
247
248 if (remove_node)
249 {
250 // Node near cell so set it to be removed from boundary set
251 nearby_node_indices.insert(coarse_mesh_boundary_node_idx);
252 }
253 }
254
255 // Remove nodes that are near cells from boundary set
256 for (unsigned int nearby_node_idx : nearby_node_indices)
257 {
258 coarse_mesh_boundary_node_indices.erase(nearby_node_idx);
259 }
260
261
262 // Apply boundary condition to the nodes in the set coarse_mesh_boundary_node_indices
263 if (this->IsNeumannBoundaryCondition())
264 {
266 }
267 else
268 {
269 // Impose any Dirichlet boundary conditions
270 for (unsigned int coarse_mesh_boundary_node_idx : coarse_mesh_boundary_node_indices)
271 {
272 pBcc->AddDirichletBoundaryCondition(this->mpFeMesh->GetNode(coarse_mesh_boundary_node_idx), this->mpBoundaryCondition.get(), 0, false);
273 this->mIsDirichletBoundaryNode[coarse_mesh_boundary_node_idx] = 1.0;
274 }
275 }
276 }
277 }
278 else // Apply BC at boundary of box domain FE mesh
279 {
280 if (this->IsNeumannBoundaryCondition())
281 {
282 // Impose any Neumann boundary conditions
283 for (auto elem_iter = this->mpFeMesh->GetBoundaryElementIteratorBegin();
284 elem_iter != this->mpFeMesh->GetBoundaryElementIteratorEnd();
285 ++elem_iter)
286 {
287 pBcc->AddNeumannBoundaryCondition(*elem_iter, this->mpBoundaryCondition.get());
288 }
289 }
290 else
291 {
292 // Impose any Dirichlet boundary conditions
293 for (auto node_iter = this->mpFeMesh->GetBoundaryNodeIteratorBegin();
294 node_iter != this->mpFeMesh->GetBoundaryNodeIteratorEnd();
295 ++node_iter)
296 {
297 pBcc->AddDirichletBoundaryCondition(*node_iter, this->mpBoundaryCondition.get());
298 this->mIsDirichletBoundaryNode[(*node_iter)->GetIndex()] = 1.0;
299 }
300 }
301 }
302}
303
304
305
306template<unsigned DIM>
308{
309 AbstractPdeModifier<DIM>::SetupSolve(rCellPopulation, outputDirectory);
310
311 InitialiseCellPdeElementMap(rCellPopulation);
312}
313
314template<unsigned DIM>
315void AbstractBoxDomainPdeModifier<DIM>::GenerateFeMesh(boost::shared_ptr<ChasteCuboid<DIM> > pMeshCuboid, double stepSize)
316{
317 // Create a regular coarse tetrahedral mesh
318 this->mpFeMesh = new TetrahedralMesh<DIM,DIM>();
319
320 GenerateAndReturnFeMesh(pMeshCuboid, stepSize, this->mpFeMesh);
321}
322
323template<unsigned DIM>
325{
326 // Create a regular coarse tetrahedral mesh
327 switch (DIM)
328 {
329 case 1:
330 pMesh->ConstructRegularSlabMesh(stepSize, pMeshCuboid->GetWidth(0));
331 break;
332 case 2:
333 pMesh->ConstructRegularSlabMesh(stepSize, pMeshCuboid->GetWidth(0), pMeshCuboid->GetWidth(1));
334 break;
335 case 3:
336 pMesh->ConstructRegularSlabMesh(stepSize, pMeshCuboid->GetWidth(0), pMeshCuboid->GetWidth(1), pMeshCuboid->GetWidth(2));
337 break;
338 default:
340 }
341
342 // Get centroid of meshCuboid
343 ChastePoint<DIM> upper = pMeshCuboid->rGetUpperCorner();
344 ChastePoint<DIM> lower = pMeshCuboid->rGetLowerCorner();
345 c_vector<double,DIM> centre_of_cuboid = 0.5*(upper.rGetLocation() + lower.rGetLocation());
346
347 // Find the centre of the PDE mesh
348 c_vector<double,DIM> centre_of_coarse_mesh = zero_vector<double>(DIM);
349 for (unsigned i=0; i<pMesh->GetNumNodes(); i++)
350 {
351 centre_of_coarse_mesh += pMesh->GetNode(i)->rGetLocation();
352 }
353 centre_of_coarse_mesh /= pMesh->GetNumNodes();
354
355 // Now move the mesh to the correct location
356 pMesh->Translate(centre_of_cuboid - centre_of_coarse_mesh);
357}
358
359template<unsigned DIM>
361{
362 // Store the PDE solution in an accessible form
363 ReplicatableVector solution_repl(this->mSolution);
364
365 if (mUseVoronoiCellsForInterpolation)
366 {
367 unsigned num_nodes = rCellPopulation.GetNumNodes();
368
369 std::vector<double> cell_data(num_nodes, -1);
370 std::vector<unsigned> num_cells(num_nodes, -1);
371
372 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
373 cell_iter != rCellPopulation.End();
374 ++cell_iter)
375 {
376 unsigned cell_location_index = rCellPopulation.GetLocationIndexUsingCell(*cell_iter);
377 cell_data[cell_location_index]=0.0;
378 num_cells[cell_location_index]=0;
379 }
380
381 // Loop over nodes of the finite element mesh and work out which voronoi region the node is in.
382 for (typename TetrahedralMesh<DIM,DIM>::NodeIterator node_iter = this->mpFeMesh->GetNodeIteratorBegin();
383 node_iter != this->mpFeMesh->GetNodeIteratorEnd();
384 ++node_iter)
385 {
386 unsigned node_index = node_iter->GetIndex();
387
388 c_vector<double,DIM> node_location = node_iter->rGetLocation();
389
390 double closest_separation = DBL_MAX;
391 unsigned nearest_cell = UNSIGNED_UNSET;
392
393 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
394 cell_iter != rCellPopulation.End();
395 ++cell_iter)
396 {
397 c_vector<double, DIM> cell_location = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
398
399 double separation = norm_2(node_location - cell_location);
400
401 if (separation < closest_separation)
402 {
403 closest_separation = separation;
404 nearest_cell = rCellPopulation.GetLocationIndexUsingCell(*cell_iter);
405 }
406 }
407 assert(closest_separation<DBL_MAX);
408
409 cell_data[nearest_cell] = cell_data[nearest_cell] + solution_repl[node_index];
410 num_cells[nearest_cell] = num_cells[nearest_cell] + 1;
411 }
412
413 // Now calculate the solution in the cell by averaging over all nodes in the voronoi region.
414 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
415 cell_iter != rCellPopulation.End();
416 ++cell_iter)
417 {
418 unsigned cell_location_index = rCellPopulation.GetLocationIndexUsingCell(*cell_iter);
419
420
421 if (num_cells[cell_location_index]==0)
422 {
423 EXCEPTION("One or more of the cells doesnt contain any pde nodes so cant use Voronoi CellData calculation in the ");
424 }
425
426 double solution_at_cell = cell_data[cell_location_index]/num_cells[cell_location_index];
427
428 cell_iter->GetCellData()->SetItem(this->mDependentVariableName, solution_at_cell);
429 }
430
431 if (this->mOutputGradient)
432 {
433 // This is not implemented yet
435 }
436 }
437 else // Interpolate solutions
438 {
439 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
440 cell_iter != rCellPopulation.End();
441 ++cell_iter)
442 {
443 // The cells are not nodes of the mesh, so we must interpolate
444 double solution_at_cell = 0.0;
445
446 // Find the element in the FE mesh that contains this cell. CellElementMap has been updated so use this.
447 unsigned elem_index = mCellPdeElementMap[*cell_iter];
448 Element<DIM,DIM>* p_element = this->mpFeMesh->GetElement(elem_index);
449
450 const ChastePoint<DIM>& node_location = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
451
452 c_vector<double,DIM+1> weights = p_element->CalculateInterpolationWeights(node_location);
453
454 for (unsigned i=0; i<DIM+1; i++)
455 {
456 double nodal_value = solution_repl[p_element->GetNodeGlobalIndex(i)];
457 solution_at_cell += nodal_value * weights(i);
458 }
459
460 cell_iter->GetCellData()->SetItem(this->mDependentVariableName, solution_at_cell);
461
462 if (this->mOutputGradient)
463 {
464 // Now calculate the gradient of the solution and store this in CellVecData
465 c_vector<double, DIM> solution_gradient = zero_vector<double>(DIM);
466
467 // Calculate the basis functions at any point (e.g. zero) in the element
468 c_matrix<double, DIM, DIM> jacobian, inverse_jacobian;
469 double jacobian_det;
470 this->mpFeMesh->GetInverseJacobianForElement(elem_index, jacobian, jacobian_det, inverse_jacobian);
471 const ChastePoint<DIM> zero_point;
472 c_matrix<double, DIM, DIM+1> grad_phi;
473 LinearBasisFunction<DIM>::ComputeTransformedBasisFunctionDerivatives(zero_point, inverse_jacobian, grad_phi);
474
475 for (unsigned node_index=0; node_index<DIM+1; node_index++)
476 {
477 double nodal_value = solution_repl[p_element->GetNodeGlobalIndex(node_index)];
478
479 for (unsigned j=0; j<DIM; j++)
480 {
481 solution_gradient(j) += nodal_value* grad_phi(j, node_index);
482 }
483 }
484
485 switch (DIM)
486 {
487 case 1:
488 cell_iter->GetCellData()->SetItem(this->mDependentVariableName+"_grad_x", solution_gradient(0));
489 break;
490 case 2:
491 cell_iter->GetCellData()->SetItem(this->mDependentVariableName+"_grad_x", solution_gradient(0));
492 cell_iter->GetCellData()->SetItem(this->mDependentVariableName+"_grad_y", solution_gradient(1));
493 break;
494 case 3:
495 cell_iter->GetCellData()->SetItem(this->mDependentVariableName+"_grad_x", solution_gradient(0));
496 cell_iter->GetCellData()->SetItem(this->mDependentVariableName+"_grad_y", solution_gradient(1));
497 cell_iter->GetCellData()->SetItem(this->mDependentVariableName+"_grad_z", solution_gradient(2));
498 break;
499 default:
501 }
502 }
503 }
504 }
505}
506
507template<unsigned DIM>
509{
510 mCellPdeElementMap.clear();
511
512 // Find the element of mpFeMesh that contains each cell and populate mCellPdeElementMap
513 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
514 cell_iter != rCellPopulation.End();
515 ++cell_iter)
516 {
517 const ChastePoint<DIM>& r_position_of_cell = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
518 unsigned elem_index = this->mpFeMesh->GetContainingElementIndex(r_position_of_cell);
519 mCellPdeElementMap[*cell_iter] = elem_index;
520 }
521}
522
523template<unsigned DIM>
525{
526 // Find the element of mpCoarsePdeMesh that contains each cell and populate mCellPdeElementMap
527 for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = rCellPopulation.Begin();
528 cell_iter != rCellPopulation.End();
529 ++cell_iter)
530 {
531 const ChastePoint<DIM>& r_position_of_cell = rCellPopulation.GetLocationOfCellCentre(*cell_iter);
532 unsigned elem_index = this->mpFeMesh->GetContainingElementIndexWithInitialGuess(r_position_of_cell, mCellPdeElementMap[*cell_iter]);
533 mCellPdeElementMap[*cell_iter] = elem_index;
534 }
535}
536
537template<unsigned DIM>
539{
540 // No parameters to output, so just call method on direct parent class
542}
543
544// Explicit instantiation
#define EXCEPTION(message)
#define NEVER_REACHED
const unsigned UNSIGNED_UNSET
Definition Exception.hpp:53
AbstractBoxDomainPdeModifier(boost::shared_ptr< AbstractLinearPde< DIM, DIM > > pPde=boost::shared_ptr< AbstractLinearPde< DIM, DIM > >(), boost::shared_ptr< AbstractBoundaryCondition< DIM > > pBoundaryCondition=boost::shared_ptr< AbstractBoundaryCondition< DIM > >(), bool isNeumannBoundaryCondition=true, boost::shared_ptr< ChasteCuboid< DIM > > pMeshCuboid=boost::shared_ptr< ChasteCuboid< DIM > >(), double stepSize=1.0, Vec solution=nullptr)
void OutputSimulationModifierParameters(out_stream &rParamsFile) override
void SetTypicalCellRadius(double typicalCellRadius)
void GenerateFeMesh(boost::shared_ptr< ChasteCuboid< DIM > > pMeshCuboid, double stepSize)
void GenerateAndReturnFeMesh(boost::shared_ptr< ChasteCuboid< DIM > > pMeshCuboid, double stepSize, TetrahedralMesh< DIM, DIM > *pMesh)
void SetupSolve(AbstractCellPopulation< DIM, DIM > &rCellPopulation, std::string outputDirectory) override
boost::shared_ptr< ChasteCuboid< DIM > > mpMeshCuboid
void ConstructBoundaryConditionsContainerHelper(AbstractCellPopulation< DIM, DIM > &rCellPopulation, std::shared_ptr< BoundaryConditionsContainer< DIM, DIM, 1 > > pBcc)
void SetUseVoronoiCellsForInterpolation(bool useVoronoiCellsForInterpolation)
void UpdateCellData(AbstractCellPopulation< DIM, DIM > &rCellPopulation)
void SetBcsOnBoxBoundary(bool setBcsOnBoxBoundary)
void UpdateCellPdeElementMap(AbstractCellPopulation< DIM, DIM > &rCellPopulation)
void InitialiseCellPdeElementMap(AbstractCellPopulation< DIM, DIM > &rCellPopulation)
void SetBcsOnBoundingSphere(bool setBcsOnBoundingSphere)
virtual c_vector< double, SPACE_DIM > GetLocationOfCellCentre(CellPtr pCell)=0
unsigned GetLocationIndexUsingCell(CellPtr pCell)
virtual unsigned GetNumNodes()=0
unsigned GetNodeGlobalIndex(unsigned localIndex) const
virtual void Translate(const c_vector< double, SPACE_DIM > &rDisplacement)
virtual unsigned GetNumNodes() const
Node< SPACE_DIM > * GetNode(unsigned index) const
void SetupSolve(AbstractCellPopulation< DIM, DIM > &rCellPopulation, std::string outputDirectory) override
TetrahedralMesh< DIM, DIM > * mpFeMesh
void OutputSimulationModifierParameters(out_stream &rParamsFile) override
std::vector< double > mIsDirichletBoundaryNode
void ConstructRegularSlabMesh(double spaceStep, double width, double height=0, double depth=0)
c_vector< double, DIM > & rGetLocation()
c_vector< double, SPACE_DIM+1 > CalculateInterpolationWeights(const ChastePoint< SPACE_DIM > &rTestPoint)
Definition Element.cpp:224
static void ComputeTransformedBasisFunctionDerivatives(const ChastePoint< ELEMENT_DIM > &rPoint, const c_matrix< double, ELEMENT_DIM, ELEMENT_DIM > &rInverseJacobian, c_matrix< double, ELEMENT_DIM, ELEMENT_DIM+1 > &rReturnValue)