v0.14.0
Loading...
Searching...
No Matches
Remodeling.hpp
Go to the documentation of this file.
1/** \file Remodeling.hpp
2* \ingroup bone_remodelling
3* \example Remodeling.hpp
4
5* \brief Implementation of element for bone remodeling
6
7*/
8
9/*
10 * This file is part of MoFEM.
11 * MoFEM is free software: you can redistribute it and/or modify it under
12 * the terms of the GNU Lesser General Public License as published by the
13 * Free Software Foundation, either version 3 of the License, or (at your
14 * option) any later version.
15 *
16 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19 * License for more details.
20 *
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
23
24#ifndef __REMODELING_HPP__
25#define __REMODELING_HPP__
26
27namespace BoneRemodeling {
28
29/**
30 * \brief Implementation of bone remodeling finite element
31 *
32 * Implementation base on paper \cite waffenschmidt2012anisotropic
33 * <http://biomechanics.stanford.edu/paper/IJSS12.pdf>
34 *
35 */
36struct Remodeling {
37
38 /**
39 * \brief Volume finite element
40 */
44 /** \brief Set integrate rule
45 */
46 int getRule(int order) { return 2 * (order - 1) + 1; };
47 };
48
49 /**
50 * \brief Not used at this stage. Could be used to do some calculations,
51 * before assembly of local elements.
52 */
53 struct FePrePostProcessRhs : public FEMethod {
54
55 MoFEMErrorCode preProcess() {
57 switch (ts_ctx) {
58 case CTX_TSSETIFUNCTION: {
59 snes_ctx = CTX_SNESSETFUNCTION;
60 snes_f = ts_F;
61 break;
62 }
63 default:
64 break;
65 }
67 }
68
69 MoFEMErrorCode postProcess() {
72 }
73 };
74
75 /**
76 * \brief Not used at this stage. Could be used to do some calculations,
77 * before assembly of local elements.
78 */
79 struct FePrePostProcessLhs : public FEMethod {
80
81 MoFEMErrorCode preProcess() {
82 //
84 switch (ts_ctx) {
85 case CTX_TSSETIJACOBIAN: {
86 snes_ctx = CTX_SNESSETJACOBIAN;
87 snes_B = ts_B;
88 break;
89 }
90 default:
91 break;
92 }
94 }
95
96 MoFEMErrorCode postProcess() {
99 }
100 };
101
102 /**
103 * Data structure for storing material parameters and evaluated values at
104 * integration points.
105 *
106 * @param oRder order of approximation
107 * @param lambda Lame parameter \f[ \lambda=\frac{E}{2(1+v)} \f]
108 * @param mu Lame parameter \f[ \mu=\frac{E\cdot\nu}{(1-2\nu)(1+\nu)} \f]
109 * @param c density evolution (growth) velocity [d/m^2] \f[ k_{p}^{\ast} \f]
110 * in Ellen Kuhl's paper
111 * @param m algorithmic exponent [-]
112 * @param n porosity exponent [-]
113 * @param rHo_ref reference density
114 * @param pSi_ref reference free energy
115 * @param R0 mass conduction coefficient
116 */
117 struct CommonData {
118
119 Mat A;
120 Vec F, D;
121 // Approximation
122 int oRder;
123 BitRefLevel bitLevel;
124
125 // DM
126 DMType dm_name; ///< dm (problem) name
127 DM dm; ///< Discretization manager
128 TS ts; ///< Time solver
129
130 boost::shared_ptr<Fe> feLhs; ///< FE to make left hand side
131 boost::shared_ptr<Fe> feRhs; ///< FE to make right hand side
132 boost::shared_ptr<FePrePostProcessRhs> preProcRhs;
133 boost::shared_ptr<FePrePostProcessLhs> preProcLhs;
134
135 boost::ptr_map<string, NeummanForcesSurface>
136 neumannForces; ///< Forces on surface
137 boost::ptr_map<string, NodalForce> nodalForces; ///< Nodal forces
138 boost::ptr_map<string, EdgeForce> edgeForces; ///< Forces on edges
139
140 boost::shared_ptr<NonlinearElasticElement> elasticPtr;
141 boost::shared_ptr<ElasticMaterials> elasticMaterialsPtr;
142
143 // Material parameters
144 double lambda; ///< Lame parameter
145 double mu; ///< Lame parameter
146 double c; ///< density evolution (growth) velocity [d/m^2]
147 double m; ///< algorithmic exponent [-]
148 double n; ///< porosity exponent [-]
149
150 double rHo_ref; ///< reference density
151 double rHo_max; ///< max density
152 double rHo_min; ///< min density
153 int b; ///< b exponent for bell function
154 double pSi_ref; ///< reference free energy
155 double R0; ///< mass conduction coefficient
156
157 double
158 cUrrent_psi; ///< current free energy for evaluating equilibrium state
159 double
160 cUrrent_mass; ///< current free energy for evaluating equilibrium state
161 PetscBool with_adol_c;
162 PetscBool is_atom_testing; ///< for atom tests
163 PetscBool less_post_proc; ///< reduce file size
165 Range tEts_all; // all blockset
166 Range tEts_block; // blockset with no remodelling
167
168 MoFEMErrorCode getParameters();
169 inline double getCFromDensity(const double &rho);
170 inline double getCFromDensityDiff(const double &rho);
171 // aDouble values
173 aC; ///< right Cauchy-Green deformation tensor \f[ C=F^{T}\cdot F \f]
174 adouble aPsi; ///< Elastic energy
175 ///< \f[\psi_{0}=\frac{\mu}{2}\left(\textrm{tr}(\mathbf{C})-3\right)-\mu\ln(J)+\frac{\lambda}{2}\ln^{2}(\ln
176 ///< J)^{2} \f]
177
178 const int tAg;
179 const int kEep;
180
181 CommonData() : tAg(1), kEep(0) {}
182
184 boost::shared_ptr<MatrixDouble>
185 gradDispPtr; ///< Ptr to gradient of displacements matrix container
186 boost::shared_ptr<VectorDouble>
187 rhoPtr; ///< Ptr to density matrix container
188 boost::shared_ptr<MatrixDouble> gradRhoPtr; ///< Gradient of density field
189 VectorDouble vecPsi; ///< Elastic energy density
190 MatrixDouble matS; ///< 2nd Piola stress
191 MatrixDouble matP; ///< 1st Piola stress
192 VectorDouble vecR; ///< Mass sorce
193 MatrixDouble matMaterialTangent; ///< Tangent matrix (derivative for F and
194 ///< 1st Piola stress)
195 MatrixDouble
196 matPushedMaterialTangent; ///< Pushed tangent matrix (derivative for F
197 ///< and 1st Piola stress)
198 MatrixDouble matGradientOfDeformation; ///< Gradient of deformation
199 MatrixDouble matInvF; ///< inverse of deformation gradient
200 VectorDouble vecDetF; ///< determinant of F
201 VectorDouble vecRhoDt; ///< Time derivative of density
202 };
204 };
205
208
209 Remodeling(MoFEM::Interface &m_field, CommonData &common_data)
210 : mField(m_field), commonData(common_data) {}
211
212 /**
213 * \brief Get parameters form line command or config file
214
215 * Read command line and config file to setup material and model parameters.
216
217 * @return Error code
218 */
219 MoFEMErrorCode getParameters();
220
221 /**
222 * \brief Set and add entities to approximation fields
223 * @return Error code
224 */
225 MoFEMErrorCode addFields();
226
227 /**
228 * \brief Set and add finite elements
229 * @return Error code
230 */
231 MoFEMErrorCode addElements();
232
233 /**
234 * \brief (Testing only) Set finite element to run mass transport problem only
235 * @return Error code
236 */
238
239 /**
240 * \brief (Testing only) Set finite element to run elastic problem only
241 * @return Error code
242 */
244
245 /**
246 * \brief Finite elements to calculate tractions
247 * @return Error code
248 */
249 MoFEMErrorCode addMomentumFluxes();
250
251 /**
252 * \brief Set problem and DM
253 * @return Error code
254 */
255 MoFEMErrorCode buildDM();
256
257 /**
258 * \brief Solve problem set up in DM
259 * @return Error code
260 */
261 MoFEMErrorCode solveDM();
262};
263
264/**
265* Calculate free energy
266*
267\f[\psi_{0}=\frac{\mu}{2}\left(\textrm{tr}(\mathbf{C})-3\right)-\mu\ln(J)+\frac{\lambda}{2}\ln^{2}(\ln
268J) \f]
269*
270*/
271template <class B1, class B2, class T>
272inline MoFEMErrorCode freeEnergy(Remodeling::CommonData &common_data, T &C,
273 B1 &psi) {
274
278 const double mu = common_data.mu;
279 const double lambda = common_data.lambda;
280
281 // Compressible Neo-Hookean
282 B2 det;
283 CHKERR determinantTensor3by3(C, det);
284 det = sqrt(det);
285 B2 traceC;
286 traceC = C(I, I);
287 B2 log_det = log(det);
288 psi = 0.5 * mu * (traceC - 3.0) - mu * log_det +
289 0.5 * lambda * log_det * log_det;
290
291 // St. Venant–Kirchhoff Material
292 // T E;
293 // E(I,J) = C(I,J);
294 // E(0,0) -= 1;
295 // E(1,1) -= 1;
296 // E(2,2) -= 1;
297 // E(I,J) *= 0.5;
298 // B2 traceE = E(I,I);
299 // psi = 0.5*lambda*traceE*traceE+mu*E(I,J)*E(I,J);
300
302}
303
304template <class B1, class T>
305inline MoFEMErrorCode recordFreeEnergy_dC(Remodeling::CommonData &common_data,
306 T &dC, B1 &psi) {
307
311 trace_on(common_data.tAg, common_data.kEep);
312 common_data.aC(I, J) <<= dC(I, J); // Set independent variables
313 CHKERR freeEnergy<adouble, adouble, FTensor::Tensor2_symmetric<adouble, 3>>(
314 common_data, common_data.aC, common_data.aPsi);
315 double psi_val;
316 common_data.aPsi >>= psi_val; // Set dependent variables
317 psi = psi_val;
318 trace_off();
320}
321/**
322 * \brief Evaluate density derivative with respect to time
323 in case of Backward Euler Method
324 \f[
325 \frac{\partial \rho}{\partial t} = \frac{\rho_{n+1} - \rho_{n}}{\Delta t}
326 \f]
327 The density within a finite element is approximated directly as
328 \f[
329 \rho = N r
330 \f]
331 */
334
336
338
339 MoFEMErrorCode doWork(int side, EntityType type,
340 DataForcesAndSourcesCore::EntData &data);
341};
342
343/**
344 * \brief Evaluate physical equations at integration points.
345
346 \f[
347 F_{ij}=u_{i,j}+\delta_{ij} \\
348 C_{ij}=F_{k,i}F_{k,j} \\
349 P_{ij}=F_{i,k}S_{k,j} \\
350 \psi_{0}=\frac{\mu}{2}\left[\textrm{tr}(C)-3\right]-\mu\ln(\sqrt{\det
351 C})+\frac{\lambda}{2}\ln^{2}(\sqrt{\det C}) \\
352 R=c\left[\left[\frac{\rho}{\rho_0}\right]^{-m}\Psi-\Psi_0\right] \\
353 \f]
354
355 */
358
360 VectorDouble vecC;
361
363
364 MoFEMErrorCode doWork(int side, EntityType type,
365 DataForcesAndSourcesCore::EntData &data);
366};
367
368/**
369 * \brief Used to post proc stresses, energy, source term.
370 *
371 * calculation of principal stresses. for e.g.
372 *
373 */
376
378 moab::Interface &postProcMesh;
379 std::vector<EntityHandle> &mapGaussPts;
380
381 OpPostProcStress(moab::Interface &post_proc_mesh,
382 std::vector<EntityHandle> &map_gauss_pts,
383 Remodeling::CommonData &common_data);
384
385 MoFEMErrorCode doWork(int side, EntityType type,
386 DataForcesAndSourcesCore::EntData &data);
387};
388
389/**
390 * Evaluate tangent material matrix
391
392 \f[
393 C_{IJ} = F_{iI}F_{iJ}\\
394 P_{iI} = F_{iI}S_{IJ}(\mathbf{C})\\
395 \delta \Psi = \delta v_{iI}P_{iI} = \delta v_{iJ} F_{iI}S_{IJ}\\
396 D_{iJkL} \delta u_{k,L} = S_{IJ}F_{iI,kL} \delta u_{k,L} +
397 F_{iI}S_{IJ,KL}F_{kL}\delta u_{k,L}\\
398 D_{iJkL} \delta u_{k,L} = S_{IJ}\delta_{ik}\delta_{IL} \delta u_{k,L} +
399 F_{iI}S_{IJ,KL}F_{kL}\delta u_{k,L}\\
400 D_{iJkL} \delta u_{k,L} = \left( S_{IJ}\delta_{ik}\delta_{IL} +
401 F_{iI}S_{IJ,KL}F_{kL} \right)\delta u_{k,L}\\ \f]
402
403 */
406
408
410
411 MoFEMErrorCode doWork(int side, EntityType type,
412 DataForcesAndSourcesCore::EntData &data);
413};
414
417
419 VectorDouble vecC;
420
422
423 MoFEMErrorCode doWork(int side, EntityType type,
424 DataForcesAndSourcesCore::EntData &data);
425};
426// TODO correct equations in documentation
427/**
428 * \brief Assemble residual for conservation of momentum (stresses)
429 *
430 \f[
431 R^{\phi}=\intop_{V}\left[\frac{\rho_{0}}{\rho_{0}^{\ast}}\right]^{n}
432 \mathbf{P}_{ij}\nabla N_{j}\,dV
433 \f]
434 */
437
439 VectorDouble nF; ///< Vector of the right hand side (internal forces)
440
442
443 MoFEMErrorCode doWork(int side, EntityType type,
444 DataForcesAndSourcesCore::EntData &data);
445};
446
447/**
448 * \brief Assemble residual for conservation of mass (density)
449 \f[
450 R^{\rho}=\intop_{V}N\frac{\partial\rho}{\partial t}-N\mathbf{R}-R_{0}\nabla
451 N\nabla\rho\,dV \f]
452 */
455
457 VectorDouble nF; ///< Vector of the right hand side (internal forces)
458
460
461 MoFEMErrorCode doWork(int side, EntityType type,
462 DataForcesAndSourcesCore::EntData &data);
463};
464
465/**
466 * \brief Diagonal block of tangent matrix \f$K_{\rho\rho}\f$
467 \f[
468 K^{\rho \rho}=\intop_{V} \Delta t N_i N_j - R_0 \nabla N_i \nabla N_j -c
469 \frac{n-m}{\rho_0}\left[\frac{\rho_{0}}{\rho_{0}^{\ast}}\right]^{n-m} \psi_0
470 N_i N_j \,dV \f]
471 */
474
476 MatrixDouble locK_rho_rho;
477 MatrixDouble transLocK_rho_rho;
478
480
481 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
482 EntityType col_type,
483 DataForcesAndSourcesCore::EntData &row_data,
484 DataForcesAndSourcesCore::EntData &col_data);
485};
486
487/**
488 * \brief Off diagonal block of tangent matrix \f$K_{\rho u}\f$
489 \f[
490 K_{\rho u}=-\intop_{V} c \left[\frac{\rho_{0}}{\rho_{0}^{\ast}}\right]^{n-m}
491 \nabla N_j P_{ij} N_i \,dV \f]
492 */
495
497 MatrixDouble locK_rho_F;
498
500 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
501 EntityType col_type,
502 DataForcesAndSourcesCore::EntData &row_data,
503 DataForcesAndSourcesCore::EntData &col_data);
504};
505
506/**
507 * \brief Off diagonal block of tangent matrix \f$K_{u u}\f$
508 \f[
509 K_{u u}=\left[\frac{\rho_{0}}{\rho_{0}^{\ast}}\right]^{n} \nabla N_j
510 D_{ijkl}\nabla N_l S_{jl} \f]
511 */
512template <bool ADOLC>
515
517 MatrixDouble locK_P_F;
518 MatrixDouble transLocK_P_F;
520
522
523 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
524 EntityType col_type,
525 DataForcesAndSourcesCore::EntData &row_data,
526 DataForcesAndSourcesCore::EntData &col_data);
527};
528
529/**
530 * \brief Off diagonal block of tangent matrix \f$K_{u \rho}\f$
531/f[
532 K_{u \rho}=\intop_{V} \left[\frac{n}{\rho_{0}}\right]
533\left[\frac{\rho_{0}}{\rho_{0}^{\ast}}\right]^{n} \nabla N_j P_{ij} N_i \,dV
534 /f]
535 */
538
540 MatrixDouble locK_P_RHO;
541
543
544 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
545 EntityType col_type,
546 DataForcesAndSourcesCore::EntData &row_data,
547 DataForcesAndSourcesCore::EntData &col_data);
548};
549
550/**
551 * Element used to post-process results at each time step
552 */
553struct MonitorPostProc : public FEMethod {
554
559 // DensityMapFe densityMaps;
560 PetscBool mass_postproc;
562 double rate;
563 bool iNit;
564 int pRT;
565
567 Remodeling::CommonData &common_data);
568
569 MoFEMErrorCode preProcess();
570 MoFEMErrorCode operator()();
571 MoFEMErrorCode postProcess();
572};
573
574/**
575 * Operator used to calculate potential energy at each time step, for
576 * postprocessing
577 */
580
584
586 Remodeling::CommonData &common_data,
587 Vec energy_vec, Vec mass_vec)
588 : MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator(
590 energVec(energy_vec), massVec(mass_vec), commonData(common_data) {}
591
592 MoFEMErrorCode doWork(int row_side, EntityType row_type,
593 DataForcesAndSourcesCore::EntData &row_data);
594};
595
596} // namespace BoneRemodeling
597
598#endif //__REMODELING_HPP__
static Index< 'J', 3 > J
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
constexpr int order
static double lambda
MoFEM::TsCtx * ts_ctx
Definition: level_set.cpp:1884
const double T
MoFEMErrorCode recordFreeEnergy_dC(Remodeling::CommonData &common_data, T &dC, B1 &psi)
Definition: Remodeling.hpp:305
MoFEMErrorCode freeEnergy(Remodeling::CommonData &common_data, T &C, B1 &psi)
Definition: Remodeling.hpp:272
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
constexpr IntegrationType I
double rho
Definition: plastic.cpp:191
constexpr auto field_name
Definition: Remodeling.hpp:553
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:556
int pRT
Definition: Remodeling.hpp:564
MoFEMErrorCode postProcess()
MoFEMErrorCode operator()()
MoFEMErrorCode preProcess()
MoFEM::Interface & mField
Definition: Remodeling.hpp:555
PetscBool mass_postproc
Definition: Remodeling.hpp:560
PetscBool equilibrium_flg
Definition: Remodeling.hpp:561
double rate
Definition: Remodeling.hpp:562
bool iNit
Definition: Remodeling.hpp:563
PostProcVolumeOnRefinedMesh postProcElastic
Definition: Remodeling.hpp:558
PostProcVolumeOnRefinedMesh postProc
Definition: Remodeling.hpp:557
Off diagonal block of tangent matrix .
Definition: Remodeling.hpp:494
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:496
Diagonal block of tangent matrix .
Definition: Remodeling.hpp:473
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:475
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
Definition: Remodeling.cpp:969
Assemble residual for conservation of mass (density)
Definition: Remodeling.hpp:454
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:889
VectorDouble nF
Vector of the right hand side (internal forces)
Definition: Remodeling.hpp:457
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:456
Off diagonal block of tangent matrix .
Definition: Remodeling.hpp:514
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
FTensor::Tensor2< double, 3, 3 > diffDiff
Definition: Remodeling.hpp:519
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:516
Off diagonal block of tangent matrix /f[ K_{u \rho}=\intop_{V} \left[\frac{n}{\rho_{0}}\right] \left...
Definition: Remodeling.hpp:537
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:539
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
Assemble residual for conservation of momentum (stresses)
Definition: Remodeling.hpp:436
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:438
VectorDouble nF
Vector of the right hand side (internal forces)
Definition: Remodeling.hpp:439
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:820
Evaluate physical equations at integration points.
Definition: Remodeling.hpp:357
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:359
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:294
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:732
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:407
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:396
Evaluate density derivative with respect to time in case of Backward Euler Method.
Definition: Remodeling.hpp:333
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:239
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:335
MoFEMErrorCode doWork(int row_side, EntityType row_type, DataForcesAndSourcesCore::EntData &row_data)
OpMassAndEnergyCalculation(const string &field_name, Remodeling::CommonData &common_data, Vec energy_vec, Vec mass_vec)
Definition: Remodeling.hpp:585
Used to post proc stresses, energy, source term.
Definition: Remodeling.hpp:375
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Definition: Remodeling.cpp:567
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:377
moab::Interface & postProcMesh
Definition: Remodeling.hpp:378
std::vector< EntityHandle > & mapGaussPts
Definition: Remodeling.hpp:379
VectorDouble vecPsi
Elastic energy density.
Definition: Remodeling.hpp:189
boost::shared_ptr< MatrixDouble > gradDispPtr
Ptr to gradient of displacements matrix container.
Definition: Remodeling.hpp:185
boost::shared_ptr< VectorDouble > rhoPtr
Ptr to density matrix container.
Definition: Remodeling.hpp:187
boost::shared_ptr< MatrixDouble > gradRhoPtr
Gradient of density field.
Definition: Remodeling.hpp:188
MatrixDouble matInvF
inverse of deformation gradient
Definition: Remodeling.hpp:199
VectorDouble vecRhoDt
Time derivative of density.
Definition: Remodeling.hpp:201
MatrixDouble matGradientOfDeformation
Gradient of deformation.
Definition: Remodeling.hpp:198
boost::shared_ptr< FePrePostProcessRhs > preProcRhs
Definition: Remodeling.hpp:132
double getCFromDensity(const double &rho)
Definition: Remodeling.cpp:218
double pSi_ref
reference free energy
Definition: Remodeling.hpp:154
boost::shared_ptr< Fe > feLhs
FE to make left hand side.
Definition: Remodeling.hpp:130
double n
porosity exponent [-]
Definition: Remodeling.hpp:148
FTensor::Tensor2_symmetric< adouble, 3 > aC
right Cauchy-Green deformation tensor
Definition: Remodeling.hpp:173
int b
b exponent for bell function
Definition: Remodeling.hpp:153
double getCFromDensityDiff(const double &rho)
Definition: Remodeling.cpp:225
double m
algorithmic exponent [-]
Definition: Remodeling.hpp:147
boost::shared_ptr< Fe > feRhs
FE to make right hand side.
Definition: Remodeling.hpp:131
PetscBool is_atom_testing
for atom tests
Definition: Remodeling.hpp:162
double c
density evolution (growth) velocity [d/m^2]
Definition: Remodeling.hpp:146
boost::ptr_map< string, NeummanForcesSurface > neumannForces
Forces on surface.
Definition: Remodeling.hpp:136
boost::ptr_map< string, EdgeForce > edgeForces
Forces on edges.
Definition: Remodeling.hpp:138
PetscBool less_post_proc
reduce file size
Definition: Remodeling.hpp:163
DM dm
Discretization manager.
Definition: Remodeling.hpp:127
boost::shared_ptr< ElasticMaterials > elasticMaterialsPtr
Definition: Remodeling.hpp:141
boost::shared_ptr< FePrePostProcessLhs > preProcLhs
Definition: Remodeling.hpp:133
boost::ptr_map< string, NodalForce > nodalForces
Nodal forces.
Definition: Remodeling.hpp:137
double R0
mass conduction coefficient
Definition: Remodeling.hpp:155
double cUrrent_psi
current free energy for evaluating equilibrium state
Definition: Remodeling.hpp:158
boost::shared_ptr< NonlinearElasticElement > elasticPtr
Definition: Remodeling.hpp:140
double cUrrent_mass
current free energy for evaluating equilibrium state
Definition: Remodeling.hpp:160
Volume finite element.
Definition: Remodeling.hpp:41
int getRule(int order)
Set integrate rule.
Definition: Remodeling.hpp:46
Fe(MoFEM::Interface &m_field)
Definition: Remodeling.hpp:42
Not used at this stage. Could be used to do some calculations, before assembly of local elements.
Definition: Remodeling.hpp:79
MoFEMErrorCode postProcess()
Definition: Remodeling.hpp:96
MoFEMErrorCode preProcess()
Definition: Remodeling.hpp:81
Not used at this stage. Could be used to do some calculations, before assembly of local elements.
Definition: Remodeling.hpp:53
MoFEMErrorCode postProcess()
Definition: Remodeling.hpp:69
MoFEMErrorCode preProcess()
Definition: Remodeling.hpp:55
Implementation of bone remodeling finite element.
Definition: Remodeling.hpp:36
MoFEMErrorCode addElementsTestingDensity()
(Testing only) Set finite element to run mass transport problem only
MoFEM::Interface & mField
Definition: Remodeling.hpp:206
MoFEMErrorCode addElements()
Set and add finite elements.
MoFEMErrorCode solveDM()
Solve problem set up in DM.
MoFEMErrorCode buildDM()
Set problem and DM.
MoFEMErrorCode getParameters()
Get parameters form line command or config file.
MoFEMErrorCode addElementsTestingElasticty()
(Testing only) Set finite element to run elastic problem only
Remodeling(MoFEM::Interface &m_field, CommonData &common_data)
Definition: Remodeling.hpp:209
MoFEMErrorCode addFields()
Set and add entities to approximation fields.
MoFEMErrorCode addMomentumFluxes()
Finite elements to calculate tractions.
Deprecated interface functions.
@ OPROW
operator doWork function is executed on FE rows
VolumeElementForcesAndSourcesCore(Interface &m_field, const EntityType type=MBTET)
Post processing.