Implementation of element for bone remodeling.
#ifndef __REMODELING_HPP__
#define __REMODELING_HPP__
struct Remodeling {
};
struct FePrePostProcessRhs :
public FEMethod {
case CTX_TSSETIFUNCTION: {
snes_ctx = CTX_SNESSETFUNCTION;
snes_f = ts_F;
break;
}
default:
break;
}
}
}
};
struct FePrePostProcessLhs :
public FEMethod {
case CTX_TSSETIJACOBIAN: {
snes_ctx = CTX_SNESSETJACOBIAN;
snes_B = ts_B;
break;
}
default:
break;
}
}
}
};
boost::shared_ptr<Fe>
feLhs;
boost::shared_ptr<Fe>
feRhs;
boost::shared_ptr<FePrePostProcessRhs>
preProcRhs;
boost::shared_ptr<FePrePostProcessLhs>
preProcLhs;
boost::ptr_map<string, NeummanForcesSurface>
boost::shared_ptr<NonlinearElasticElement>
elasticPtr;
double
double
struct DataContainers {
boost::shared_ptr<MatrixDouble>
boost::shared_ptr<VectorDouble>
MatrixDouble
};
};
};
template <class B1, class B2, class T>
inline MoFEMErrorCode
freeEnergy(Remodeling::CommonData &common_data,
T &C,
B1 &psi) {
const double mu = common_data.mu;
const double lambda = common_data.lambda;
B2 det;
CHKERR determinantTensor3by3(C, det);
det = sqrt(det);
B2 traceC;
B2 log_det = log(det);
psi = 0.5 *
mu * (traceC - 3.0) -
mu * log_det +
0.5 *
lambda * log_det * log_det;
}
template <class B1, class T>
trace_on(common_data.tAg, common_data.kEep);
common_data.aC(
I,
J) <<= dC(
I,
J);
CHKERR freeEnergy<adouble, adouble, FTensor::Tensor2_symmetric<adouble, 3>>(
common_data, common_data.aC, common_data.aPsi);
double psi_val;
common_data.aPsi >>= psi_val;
psi = psi_val;
trace_off();
}
struct OpGetRhoTimeDirevative
OpGetRhoTimeDirevative(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &data);
};
DataForcesAndSourcesCore::EntData &data);
};
std::vector<EntityHandle> &map_gauss_pts,
Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &data);
};
struct OpCalculateStressTangent
OpCalculateStressTangent(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &data);
};
struct OpCalculateStressTangentWithAdolc
OpCalculateStressTangentWithAdolc(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &data);
};
struct OpAssmbleStressRhs
OpAssmbleStressRhs(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &data);
};
struct OpAssmbleRhoRhs
OpAssmbleRhoRhs(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &data);
};
struct OpAssmbleRhoLhs_dRho
OpAssmbleRhoLhs_dRho(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &row_data,
DataForcesAndSourcesCore::EntData &col_data);
};
struct OpAssmbleRhoLhs_dF
OpAssmbleRhoLhs_dF(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &row_data,
DataForcesAndSourcesCore::EntData &col_data);
};
template <bool ADOLC>
struct OpAssmbleStressLhs_dF
OpAssmbleStressLhs_dF(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &row_data,
DataForcesAndSourcesCore::EntData &col_data);
};
struct OpAssmbleStressLhs_dRho
OpAssmbleStressLhs_dRho(Remodeling::CommonData &common_data);
DataForcesAndSourcesCore::EntData &row_data,
DataForcesAndSourcesCore::EntData &col_data);
};
Remodeling::CommonData &common_data);
};
struct OpMassAndEnergyCalculation
OpMassAndEnergyCalculation(
const string &
field_name,
Remodeling::CommonData &common_data,
Vec energy_vec, Vec mass_vec)
DataForcesAndSourcesCore::EntData &row_data);
};
}
#endif
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
MoFEMErrorCode recordFreeEnergy_dC(Remodeling::CommonData &common_data, T &dC, B1 &psi)
MoFEMErrorCode freeEnergy(Remodeling::CommonData &common_data, T &C, B1 &psi)
constexpr IntegrationType I
constexpr auto field_name
Remodeling::CommonData & commonData
MoFEMErrorCode postProcess()
MoFEMErrorCode operator()()
MoFEMErrorCode preProcess()
MoFEM::Interface & mField
PetscBool equilibrium_flg
PostProcVolumeOnRefinedMesh postProcElastic
PostProcVolumeOnRefinedMesh postProc
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
MatrixDouble transLocK_rho_rho
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
MatrixDouble locK_rho_rho
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
VectorDouble nF
Vector of the right hand side (internal forces)
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
MatrixDouble transLocK_P_F
FTensor::Tensor2< double, 3, 3 > diffDiff
Remodeling::CommonData & commonData
Remodeling::CommonData & commonData
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
VectorDouble nF
Vector of the right hand side (internal forces)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int row_side, EntityType row_type, DataForcesAndSourcesCore::EntData &row_data)
Remodeling::CommonData & commonData
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Remodeling::CommonData & commonData
moab::Interface & postProcMesh
std::vector< EntityHandle > & mapGaussPts
VectorDouble vecDetF
determinant of F
VectorDouble vecPsi
Elastic energy density.
VectorDouble vecR
Mass sorce.
boost::shared_ptr< MatrixDouble > gradDispPtr
Ptr to gradient of displacements matrix container.
boost::shared_ptr< VectorDouble > rhoPtr
Ptr to density matrix container.
MatrixDouble matPushedMaterialTangent
MatrixDouble matP
1st Piola stress
MatrixDouble matS
2nd Piola stress
boost::shared_ptr< MatrixDouble > gradRhoPtr
Gradient of density field.
MatrixDouble matMaterialTangent
MatrixDouble matInvF
inverse of deformation gradient
VectorDouble vecRhoDt
Time derivative of density.
MatrixDouble matGradientOfDeformation
Gradient of deformation.
boost::shared_ptr< FePrePostProcessRhs > preProcRhs
double getCFromDensity(const double &rho)
double pSi_ref
reference free energy
boost::shared_ptr< Fe > feLhs
FE to make left hand side.
double n
porosity exponent [-]
FTensor::Tensor2_symmetric< adouble, 3 > aC
right Cauchy-Green deformation tensor
int b
b exponent for bell function
double getCFromDensityDiff(const double &rho)
double m
algorithmic exponent [-]
DMType dm_name
dm (problem) name
boost::shared_ptr< Fe > feRhs
FE to make right hand side.
PetscBool is_atom_testing
for atom tests
double c
density evolution (growth) velocity [d/m^2]
boost::ptr_map< string, NeummanForcesSurface > neumannForces
Forces on surface.
MoFEMErrorCode getParameters()
boost::ptr_map< string, EdgeForce > edgeForces
Forces on edges.
double rHo_max
max density
PetscBool less_post_proc
reduce file size
DM dm
Discretization manager.
boost::shared_ptr< ElasticMaterials > elasticMaterialsPtr
boost::shared_ptr< FePrePostProcessLhs > preProcLhs
boost::ptr_map< string, NodalForce > nodalForces
Nodal forces.
double rHo_min
min density
double rHo_ref
reference density
double lambda
Lame parameter.
double R0
mass conduction coefficient
double cUrrent_psi
current free energy for evaluating equilibrium state
boost::shared_ptr< NonlinearElasticElement > elasticPtr
double cUrrent_mass
current free energy for evaluating equilibrium state
int getRule(int order)
Set integrate rule.
MoFEMErrorCode postProcess()
MoFEMErrorCode preProcess()
MoFEMErrorCode postProcess()
MoFEMErrorCode preProcess()
MoFEMErrorCode addElementsTestingDensity()
(Testing only) Set finite element to run mass transport problem only
MoFEM::Interface & mField
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
MoFEMErrorCode addFields()
Set and add entities to approximation fields.
MoFEMErrorCode addMomentumFluxes()
Finite elements to calculate tractions.
Deprecated interface functions.
Volume finite element base.