59 const char param_file[] =
"param_file.petsc";
62 std::vector<boost::weak_ptr<RefEntity>> v_ref_entity;
63 std::vector<boost::weak_ptr<RefElement>> v_ref_element;
64 std::vector<boost::weak_ptr<Field>> v_field;
65 boost::weak_ptr<DofEntity> arc_dof;
66 boost::weak_ptr<MWLSApprox> mwls_approx;
67 std::vector<boost::weak_ptr<FEMethod>> vec_fe;
68 std::vector<boost::weak_ptr<NonlinearElasticElement>> v_elastic_ele_str;
72 PetscBool flg = PETSC_FALSE;
74 if (flg == PETSC_TRUE)
92 MeshsetsManager::brodcastMeshsets =
false;
96 auto core_log = logging::core::get();
98 LogManager::createSink(LogManager::getStrmWorld(),
"MEMCHECK"));
99 LogManager::setLog(
"MEMCHECK");
103 auto set_log_attr = [&](
auto c) {
106 set_log_attr(
"MEMCHECK");
107 MOFEM_LOG_C(
"MEMCHECK", Sev::inform,
"### Input parameter: -my_file %s",
115 <<
"File version " << cp.fileVersion.strVersion();
121 DMType dm_name =
"MOFEM";
125 auto add_bits = [&]() {
128 ->addToDatabaseBitRefLevelByType(MBTET,
BitRefLevel().set(),
134 std::vector<int> surface_ids;
136 std::vector<std::string> fe_surf_proj_list;
137 fe_surf_proj_list.push_back(
"SURFACE");
150 ->getAllEntitiesNotInDatabase(ents);
152 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets,
false);
153 for (
auto m : meshsets)
154 CHKERR moab.remove_entities(
m, ents);
155 CHKERR moab.delete_entities(ents);
157 auto first_number_string = [](std::string
const &str) {
158 std::size_t
const n = str.find_first_of(
"0123456789");
159 if (
n != std::string::npos) {
160 std::size_t
const m = str.find_first_not_of(
"0123456789",
n);
161 return str.substr(
n,
m != std::string::npos ?
m -
n :
m);
163 return std::string();
165 std::string str_number =
168 if (!str_number.empty())
169 cp.startStep = atoi(str_number.c_str()) + 1;
171 string cutting_surface_name =
173 boost::lexical_cast<std::string>(cp.startStep) +
".vtk";
177 cp.crackAccelerationFactor, cutting_surface_name,
VERBOSE,
188 string cutting_surface_name =
190 boost::lexical_cast<std::string>(cp.startStep) +
".vtk";
194 CHKERR cp.getInterface<
CPMeshCut>()->copySurface(cutting_surface_name);
196 CHKERR moab.get_entities_by_type(0, MBTET, vol,
false);
207 CHKERR moab.get_entities_by_type(0, MBTET, vol,
false);
215 CHKERR cp.createProblemDataStructures(surface_ids,
QUIET,
218 if (cp.propagateCrack == PETSC_TRUE) {
221 CHKERR cp.setMaterialPositionFromCoords();
222 CHKERR cp.setSpatialPositionFromCoords();
223 if (cp.addSingularity == PETSC_TRUE)
224 CHKERR cp.setSingularDofs(
"SPATIAL_POSITION");
233 CHKERR cp.createDMs(dm_elastic, dm_eigen_elastic, dm_material,
234 dm_crack_propagation, dm_material_forces,
235 dm_surface_projection, dm_crack_srf_area, surface_ids,
240 CHKERR cp.setSpatialPositionFromCoords();
241 if (cp.addSingularity == PETSC_TRUE) {
242 CHKERR cp.setSingularDofs(
"SPATIAL_POSITION");
245 cp.getLoadScale() = cp.getArcCtx()->getFieldData();
247 if (cp.solveEigenStressProblem) {
248 MOFEM_LOG(
"MEMCHECK", Sev::verbose) <<
"Solve Eigen ELASTIC problem";
250 dm_eigen_elastic, cp.getEigenArcCtx()->x0, 1);
251 CHKERR cp.postProcessDM(dm_eigen_elastic, -2,
"ELASTIC",
true);
254 MOFEM_LOG(
"MEMCHECK", Sev::verbose) <<
"Solve ELASTIC problem";
256 dm_elastic, cp.getArcCtx()->x0, cp.getLoadScale());
257 CHKERR cp.postProcessDM(dm_elastic, -1,
"ELASTIC",
true);
261 CHKERR cp.assembleMaterialForcesDM(PETSC_NULL);
262 CHKERR cp.assembleSmootherForcesDM(PETSC_NULL, surface_ids,
VERBOSE,
264 CHKERR cp.assembleCouplingForcesDM(PETSC_NULL,
VERBOSE,
false);
266 CHKERR cp.calculateElasticEnergy(dm_elastic);
268 dm_material_forces, dm_surface_projection, dm_crack_srf_area,
271 if (cp.doCutMesh == PETSC_TRUE) {
273 dm_crack_propagation, dm_elastic, dm_eigen_elastic, dm_material,
274 dm_material_forces, dm_surface_projection, dm_crack_srf_area,
278 dm_crack_propagation, dm_elastic, dm_material, dm_material_forces,
279 dm_surface_projection, dm_crack_srf_area, surface_ids);
285 CHKERR cp.setMaterialPositionFromCoords();
286 CHKERR cp.setSpatialPositionFromCoords();
287 if (cp.addSingularity == PETSC_TRUE)
288 CHKERR cp.setSingularDofs(
"SPATIAL_POSITION");
297 CHKERR cp.createDMs(dm_elastic, dm_eigen_elastic, dm_material,
298 dm_crack_propagation, dm_material_forces,
299 dm_surface_projection, dm_crack_srf_area, surface_ids,
305 CHKERR cp.createBcDM(dm_bc,
"BC_PROBLEM",
306 cp.mapBitLevel[
"spatial_domain"],
309 if (cp.solveEigenStressProblem) {
310 MOFEM_LOG(
"MEMCHECK", Sev::verbose) <<
"Solve Eigen ELASTIC problem";
312 dm_eigen_elastic, cp.getEigenArcCtx()->x0, 1);
313 CHKERR cp.postProcessDM(dm_eigen_elastic, -2,
"ELASTIC",
true);
317 MOFEM_LOG(
"MEMCHECK", Sev::verbose) <<
"Solve ELASTIC problem";
319 dm_elastic, cp.getArcCtx()->x0, cp.getLoadScale());
323 CHKERR cp.assembleMaterialForcesDM(PETSC_NULL);
325 CHKERR cp.assembleSmootherForcesDM(PETSC_NULL, surface_ids,
VERBOSE,
329 dm_material_forces, dm_surface_projection, dm_crack_srf_area,
331 CHKERR cp.calculateElasticEnergy(dm_elastic);
334 CHKERR cp.postProcessDM(dm_elastic, 0,
"ELASTIC",
true);
336 CHKERR cp.savePositionsOnCrackFrontDM(dm_elastic, PETSC_NULL, 2,
false);
337 CHKERR cp.writeCrackFont(cp.mapBitLevel[
"spatial_domain"], 0);
340 CHKERR cp.tetsingReleaseEnergyCalculation();
343 v_ref_entity.reserve(refined_ents_ptr->size());
344 for (
auto e : *refined_ents_ptr)
345 v_ref_entity.emplace_back(e);
348 v_ref_element.reserve(refined_finite_elements_ptr->size());
349 for (
auto e : *refined_finite_elements_ptr)
350 v_ref_element.emplace_back(e);
353 v_field.reserve(fields_ptr->size());
354 for (
auto e : *fields_ptr)
355 v_field.emplace_back(e);
357 arc_dof = cp.arcLengthDof;
358 mwls_approx = cp.mwlsApprox;
359 v_elastic_ele_str.emplace_back(cp.elasticFe);
360 v_elastic_ele_str.emplace_back(cp.materialFe);
361 vec_fe.emplace_back(cp.feLhs);
362 vec_fe.emplace_back(cp.feRhs);
363 vec_fe.emplace_back(cp.feMaterialRhs);
364 vec_fe.emplace_back(cp.feMaterialLhs);
365 vec_fe.emplace_back(cp.feEnergy);
366 vec_fe.emplace_back(cp.feCouplingElasticLhs);
367 vec_fe.emplace_back(cp.feCouplingMaterialLhs);
368 vec_fe.emplace_back(cp.feSmootherRhs);
369 vec_fe.emplace_back(cp.feSmootherLhs);
370 vec_fe.emplace_back(cp.feGriffithForceRhs);
371 vec_fe.emplace_back(cp.feGriffithConstrainsDelta);
372 vec_fe.emplace_back(cp.feGriffithConstrainsRhs);
373 vec_fe.emplace_back(cp.feGriffithConstrainsLhs);
374 vec_fe.emplace_back(cp.feSpringLhsPtr);
375 vec_fe.emplace_back(cp.feSpringRhsPtr);
376 vec_fe.emplace_back(cp.feRhsSimpleContact);
377 vec_fe.emplace_back(cp.feLhsSimpleContact);
378 vec_fe.emplace_back(cp.feMaterialAnaliticalTraction);
383 if (
auto a = arc_dof.lock()) {
385 <<
"cp.arcLengthDof.use_count() " <<
a.use_count();
389 for (
auto v : vec_fe)
390 if (
auto a =
v.lock()) {
391 MOFEM_LOG(
"MEMCHECK", Sev::error) <<
"fe " <<
a.use_count();
395 for (
auto v : v_elastic_ele_str)
396 if (
auto a =
v.lock()) {
398 <<
"elastic structure " <<
a.use_count();
402 if (
auto a = mwls_approx.lock()) {
404 <<
"cp.mwlsApprox.use_count() " <<
a.use_count();
408 for (
auto e : v_ref_entity) {
409 if (
auto a = e.lock()) {
410 MOFEM_LOG(
"MEMCHECK", Sev::error) <<
"v_ref_entity " <<
a.use_count();
415 for (
auto e : v_ref_element) {
416 if (
auto a = e.lock()) {
417 MOFEM_LOG(
"MEMCHECK", Sev::error) <<
"v_ref_element " <<
a.use_count();
422 for (
auto e : v_field) {
423 if (
auto a = e.lock()) {
424 MOFEM_LOG(
"MEMCHECK", Sev::error) <<
"v_field " <<
a.use_count();