v0.14.0
Functions | Variables
mass_calculation.cpp File Reference
#include <MoFEM.hpp>
#include <Projection10NodeCoordsOnField.hpp>
#include <metaImage.h>
#include <moab/Skinner.hpp>
#include <moab/AdaptiveKDTree.hpp>
#include <cholesky.hpp>
#include <DensityMaps.hpp>
#include <PostProcOnRefMesh.hpp>

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "\n"
 

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 42 of file mass_calculation.cpp.

42  {
43  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
44 
45  try {
46 
47  moab::Core mb_instance;
48  moab::Interface& moab = mb_instance;
49 
50  char mesh_file_name[255];
51  PetscBool flg_file = PETSC_FALSE;
52  char meta_file_name[255];
53  PetscBool meta_flg_file = PETSC_FALSE;
54  PetscInt order = 1;
55  PetscBool is_atom_test = PETSC_FALSE;
56  // char approximation_base[255];
57 
58  ierr = PetscOptionsBegin(PETSC_COMM_WORLD,"","Bone remodeling","none"); CHKERRQ(ierr);
59  ierr = PetscOptionsString(
60  "-my_file",
61  "mesh file name","",
62  "mesh.h5m",mesh_file_name,255,&flg_file
63  ); CHKERRQ(ierr);
64 
65  ierr = PetscOptionsString(
66  "-meta_image",
67  "meta image name","",
68  "file.mhd",meta_file_name,255,&meta_flg_file
69  ); CHKERRQ(ierr);
70 
71  // ierr = PetscOptionsString(
72  // "-base",
73  // "approximation base name","",
74  // "",approximation_base,255,PETSC_NULL
75  // ); CHKERRQ(ierr);
76 
77  ierr = PetscOptionsBool(
78  "-my_is_atom_test",
79  "is used with testing, exit with error when diverged","",
80  PETSC_FALSE,&is_atom_test,PETSC_NULL
81  ); CHKERRQ(ierr);
82 
83  ierr = PetscOptionsInt(
84  "-my_order",
85  "default approximation order","",
86  1,&order,PETSC_NULL
87  ); CHKERRQ(ierr);
88 
89  ierr = PetscOptionsEnd(); CHKERRQ(ierr);
90 
91  // #ifdef WITH_METAIO
92 
93  MetaImage meta_file(meta_file_name);
94  cout << "Quantity " << meta_file.Quantity() << endl;
95  const int *dim_size = meta_file.DimSize();
96  cout << "Image dimension: " << dim_size[0] << " " << dim_size[1] << " " << dim_size[2] << endl;
97  const double *elemet_size = meta_file.ElementSize();
98  cout << "Image element size: " << elemet_size[0] << " " << elemet_size[1] << " " << elemet_size[2] << endl;
99  const double *offset = meta_file.Offset();
100  cout << "Image offset: " << offset[0] << " " << offset[1] << " " << offset[2] << endl;
101  const double * transform_matrix = meta_file.TransformMatrix();
102  cout << "Image Transform Matrix:\n";
103  for(int rr = 0;rr<3;rr++) {
104  for(int cc = 0;cc<3;cc++) {
105  cout << transform_matrix[rr*3+cc] << " ";
106  }
107  cout << endl;
108  }
109 
110  //Read parameters from line command
111  if(flg_file != PETSC_TRUE) {
112  SETERRQ(PETSC_COMM_SELF,1,"*** ERROR -my_file (MESH FILE NEEDED)");
113  }
114  const char *option;
115  option = "";
116  rval = moab.load_file(mesh_file_name, 0, option); CHKERRQ_MOAB(rval);
117  ParallelComm* pcomm = ParallelComm::get_pcomm(&moab,MYPCOMM_INDEX);
118  if(pcomm == NULL) pcomm = new ParallelComm(&moab,PETSC_COMM_WORLD);
119 
120  MoFEM::Core core(moab);
121  MoFEM::Interface& m_field = core;
122  //meshset consisting all entities in mesh
123  EntityHandle root_set = moab.get_root_set();
124  // Seed all mesh entities to MoFEM database, those entities can be potentially used as finite elements
125  // or as entities which carry some approximation field.
126  BitRefLevel bit_level0;
127  bit_level0.set(0);
128  ierr = m_field.seed_ref_level_3D(root_set,bit_level0); CHKERRQ(ierr);
129  // Add field
130  ierr = m_field.add_field("RHO",H1,AINSWORTH_LEGENDRE_BASE,1); CHKERRQ(ierr);
131  ierr = m_field.add_field("MESH_NODE_POSITIONS",H1,AINSWORTH_LEGENDRE_BASE,3); CHKERRQ(ierr);
132 
133  // Add meshsets with material and boundary conditions
134  auto *meshsets_manager_ptr = m_field.getInterface<MeshsetsManager>();
135  ierr = meshsets_manager_ptr->setMeshsetFromFile(); CHKERRQ(ierr);
136 
137  PetscPrintf(PETSC_COMM_WORLD,"Read meshset. Added meshsets for bc.cfg\n");
138  for(_IT_CUBITMESHSETS_FOR_LOOP_(m_field,it)) {
139  PetscPrintf(
140  PETSC_COMM_WORLD,
141  "%s",static_cast<std::ostringstream&>(std::ostringstream().seekp(0) << *it << endl).str().c_str()
142  );
143  cerr << *it << endl;
144  }
145 
146  Range tets_all;
147  rval = m_field.get_moab().get_entities_by_type(0,MBTET,tets_all); CHKERRQ_MOAB(rval);
148 
149  for(_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(m_field,BLOCKSET,it)) {
150  string name = it->getName();
151  if (name.compare(0,14,"NO_REMODELLING") == 0) {
152  //nOremodellingBlock = true;
153  EntityHandle meshset = it->getMeshset();
154  Range tets_block;
155  rval = m_field.get_moab().get_entities_by_type(
156  meshset,MBTET,tets_block,true
157  ); CHKERRQ_MOAB(rval);
158  tets_all = subtract(tets_all,tets_block);
159  }
160  }
161 
162 
163  // if ( strcmp(approximation_base, "lobatto") == 0){
164  // ierr = m_field.add_field("RHO",L2,AINSWORTH_LOBATTO_BASE,1); CHKERRQ(ierr);
165  // } else if ( strcmp(approximation_base, "berstein") == 0){ //not implemented yet
166  // ierr = m_field.add_field("RHO",L2,AINSOWRTH_BERNSTEIN_BEZIER_BASE,1); CHKERRQ(ierr);
167  // }
168  // else {
169  // ierr = m_field.add_field("RHO",L2,AINSWORTH_LEGENDRE_BASE,1); CHKERRQ(ierr);
170  // }
171  // Add entities to field (root_mesh, i.e. on all mesh entities fields are approx.)
172  // On those entities and lower dimension entities approximation is spaned,
173  ierr = m_field.add_ents_to_field_by_type(tets_all,MBTET,"RHO"); CHKERRQ(ierr);
174  ierr = m_field.add_ents_to_field_by_type(root_set,MBTET,"MESH_NODE_POSITIONS"); CHKERRQ(ierr);
175 
176 
177  {
178  ierr = m_field.add_field("DISTANCE_FROM_SURFACE",H1,AINSWORTH_LEGENDRE_BASE,1); CHKERRQ(ierr);
179  Range vertices;
180  rval = moab.get_entities_by_type(root_set,MBVERTEX,vertices); CHKERRQ_MOAB(rval);
181  ierr = m_field.add_ents_to_field_by_type(vertices,MBTET,"DISTANCE_FROM_SURFACE"); CHKERRQ(ierr);
182  Range tets;
183  rval = moab.get_entities_by_type(root_set,MBTET,tets); CHKERRQ_MOAB(rval);
184  Range edges;
185  rval = moab.get_adjacencies(tets,1,true,edges,moab::Interface::UNION); CHKERRQ_MOAB(rval);
186  ierr = m_field.add_ents_to_field_by_type(edges,MBEDGE,"DISTANCE_FROM_SURFACE"); CHKERRQ(ierr);
187  }
188 
189  // Set approximation oder
190  ierr = m_field.set_field_order(root_set,MBVERTEX,"RHO",1); CHKERRQ(ierr);
191  ierr = m_field.set_field_order(root_set,MBEDGE,"RHO",order); CHKERRQ(ierr);
192  ierr = m_field.set_field_order(root_set,MBTRI,"RHO",order); CHKERRQ(ierr);
193  ierr = m_field.set_field_order(root_set,MBTET,"RHO",order); CHKERRQ(ierr);
194  ierr = m_field.set_field_order(root_set,MBEDGE,"DISTANCE_FROM_SURFACE",2); CHKERRQ(ierr);
195  ierr = m_field.set_field_order(root_set,MBVERTEX,"DISTANCE_FROM_SURFACE",1); CHKERRQ(ierr);
196 
197  // Apply 2nd order geometry approximation only on skin
198  {
199 
200  Skinner skin(&moab);
201  Range faces,tets;
202  rval = moab.get_entities_by_type(0,MBTET,tets); CHKERRQ_MOAB(rval);
203  rval = skin.find_skin(0,tets,false,faces); CHKERRQ_MOAB(rval);
204  Range edges;
205  rval = moab.get_adjacencies(
206  faces,1,false,edges,moab::Interface::UNION
207  ); CHKERRQ_MOAB(rval);
208  ierr = m_field.set_field_order(edges,"MESH_NODE_POSITIONS",2); CHKERRQ(ierr);
209  }
210  ierr = m_field.set_field_order(root_set,MBVERTEX,"MESH_NODE_POSITIONS",1); CHKERRQ(ierr);
211 
212  ierr = m_field.build_fields(); CHKERRQ(ierr);
213 
214  {
215  Projection10NodeCoordsOnField ent_method_material(m_field,"MESH_NODE_POSITIONS");
216  ierr = m_field.loop_dofs("MESH_NODE_POSITIONS",ent_method_material); CHKERRQ(ierr);
217  }
218 
219  SurfaceKDTree surface_distance(m_field);
220  {
221  Range tets;
222  rval = moab.get_entities_by_type(root_set,MBTET,tets); CHKERRQ_MOAB(rval);
223  ierr = surface_distance.takeASkin(tets); CHKERRQ(ierr);
224  ierr = surface_distance.buildTree(); CHKERRQ(ierr);
225  ierr = surface_distance.setDistanceFromSurface("DISTANCE_FROM_SURFACE"); CHKERRQ(ierr);
226  }
227 
228  // Add finite element (this defines element declaration comes later)
229  ierr = m_field.add_finite_element("DENSITY_MAP_ELEMENT"); CHKERRQ(ierr);
230  ierr = m_field.modify_finite_element_add_field_row("DENSITY_MAP_ELEMENT","RHO"); CHKERRQ(ierr);
231  ierr = m_field.modify_finite_element_add_field_col("DENSITY_MAP_ELEMENT","RHO"); CHKERRQ(ierr);
232  ierr = m_field.modify_finite_element_add_field_data("DENSITY_MAP_ELEMENT","RHO"); CHKERRQ(ierr);
233  ierr = m_field.modify_finite_element_add_field_data("DENSITY_MAP_ELEMENT","DISTANCE_FROM_SURFACE"); CHKERRQ(ierr);
234  ierr = m_field.modify_finite_element_add_field_data("DENSITY_MAP_ELEMENT","MESH_NODE_POSITIONS"); CHKERRQ(ierr);
235 
236  // Add entities to that element
237  ierr = m_field.add_ents_to_finite_element_by_type(tets_all,MBTET,"DENSITY_MAP_ELEMENT"); CHKERRQ(ierr);
238 
239  //build finite elements
240  ierr = m_field.build_finite_elements(); CHKERRQ(ierr);
241  //build adjacencies between elements and degrees of freedom
242  ierr = m_field.build_adjacencies(bit_level0); CHKERRQ(ierr);
243 
244  DM dm_rho;
245  DMType dm_name = "DM_RHO";
246  // Register DM problem
247  ierr = DMRegister_MoFEM(dm_name); CHKERRQ(ierr);
248  ierr = DMCreate(PETSC_COMM_WORLD,&dm_rho);CHKERRQ(ierr);
249  ierr = DMSetType(dm_rho,dm_name);CHKERRQ(ierr);
250  // Create DM instance
251  ierr = DMMoFEMCreateMoFEM(dm_rho,&m_field,dm_name,bit_level0); CHKERRQ(ierr);
252  // Configure DM form line command options (DM itself, solvers, pre-conditioners, ... )
253  ierr = DMSetFromOptions(dm_rho); CHKERRQ(ierr);
254  // Add elements to dm (only one here)
255  ierr = DMMoFEMAddElement(dm_rho,"DENSITY_MAP_ELEMENT"); CHKERRQ(ierr);
256  // Set up problem (number dofs, partition mesh, etc.)
257  //ierr = DMSetUp(dm_rho); CHKERRQ(ierr);
258 
259  {
260  ierr = DMSetUp(dm_rho); CHKERRQ(ierr);
261  // Do serial matrix
262  // const MoFEM::Problem *problem_ptr;
263  // ierr = DMMoFEMGetProblemPtr(dm_rho,&problem_ptr); CHKERRQ(ierr);
264  // ierr = m_field.build_problem(problem_ptr->getName()); CHKERRQ(ierr);
265  // // FIXME: That is not optimal partitioning
266  // ierr = m_field.partition_simple_problem(problem_ptr->getName()); CHKERRQ(ierr);
267  // ierr = m_field.partition_finite_elements(problem_ptr->getName());
268  // ierr = m_field.partition_ghost_dofs(problem_ptr->getName()); CHKERRQ(ierr);
269  }
270 
271  Vec volume_vec;
272  int volume_vec_ghost[] = { 0 };
273  ierr = VecCreateGhost(
274  PETSC_COMM_WORLD,(!m_field.get_comm_rank())?1:0,1,1,volume_vec_ghost,&volume_vec
275  ); CHKERRQ(ierr);
276  ierr = VecZeroEntries(volume_vec); CHKERRQ(ierr);
277  Vec mass_vec;
278  ierr = VecDuplicate(volume_vec,&mass_vec); CHKERRQ(ierr);
279  Vec mass_vec_approx;
280  ierr = VecDuplicate(volume_vec,&mass_vec_approx); CHKERRQ(ierr);
281 
282  DensityMapFe fe_volume(m_field);
283  auto &list_of_operators = fe_volume.getOpPtrVector();
284  list_of_operators.push_back(new OpVolumeCalculation("RHO",volume_vec));
285  DataFromMetaIO data_from_meta_io(meta_file);
286  ierr = data_from_meta_io.fromOptions(); CHKERRQ(ierr);
287  list_of_operators.push_back(new OpMassCalculation("RHO",mass_vec,data_from_meta_io));
288  ierr = DMoFEMLoopFiniteElements(dm_rho,"DENSITY_MAP_ELEMENT",&fe_volume); CHKERRQ(ierr);
289 
290  ierr = VecAssemblyBegin(volume_vec); CHKERRQ(ierr);
291  ierr = VecAssemblyEnd(volume_vec); CHKERRQ(ierr);
292  double volume;
293  ierr = VecSum(volume_vec,&volume); CHKERRQ(ierr);
294  ierr = PetscPrintf(PETSC_COMM_WORLD,"Volume %0.32g\n",volume); CHKERRQ(ierr);
295 
296  ierr = VecAssemblyBegin(mass_vec); CHKERRQ(ierr);
297  ierr = VecAssemblyEnd(mass_vec); CHKERRQ(ierr);
298  double mass;
299  ierr = VecSum(mass_vec,&mass); CHKERRQ(ierr);
300  ierr = PetscPrintf(PETSC_COMM_WORLD,"Mass %0.32g\n",mass); CHKERRQ(ierr);
301 
302  ierr = VecDestroy(&volume_vec); CHKERRQ(ierr);
303  ierr = VecDestroy(&mass_vec); CHKERRQ(ierr);
304 
305 
306  Mat A;
307  Vec F,D;
308  ierr = DMCreateMatrix(dm_rho,&A); CHKERRQ(ierr);
309  ierr = m_field.getInterface<VecManager>()->vecCreateGhost("DM_RHO",ROW,&F); CHKERRQ(ierr);
310  ierr = m_field.getInterface<VecManager>()->vecCreateGhost("DM_RHO",COL,&D); CHKERRQ(ierr);
311 
312  ierr = MatZeroEntries(A); CHKERRQ(ierr);
313  ierr = VecZeroEntries(F); CHKERRQ(ierr);
314 
315  CommonData common_data;
316  common_data.globA = A;
317  common_data.globF = F;
318 
319  common_data.distanceValue = boost::shared_ptr<VectorDouble>(new VectorDouble());
320  common_data.distanceGradient = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
321  boost::shared_ptr<VectorDouble> dist_value = common_data.distanceValue;
322  boost::shared_ptr<MatrixDouble> grad_value = common_data.distanceGradient;
323 
324  DensityMapFe fe_density_approximation(m_field);
325 
326  fe_density_approximation.getOpPtrVector().push_back
327  (new OpCalculateScalarFieldValues("DISTANCE_FROM_SURFACE",dist_value));
328  fe_density_approximation.getOpPtrVector().push_back
329  (new OpCalculateScalarFieldGradient<3>("DISTANCE_FROM_SURFACE",grad_value));
330  fe_density_approximation.getOpPtrVector().push_back
331  (new OpCalculateLhs("RHO","RHO",common_data,data_from_meta_io));
332  fe_density_approximation.getOpPtrVector().push_back
333  (new OpCalulatefRhoAtGaussPts("RHO",common_data,data_from_meta_io,surface_distance));
334  fe_density_approximation.getOpPtrVector().push_back
335  (new OpAssmbleRhs("RHO",common_data));
336  ierr = DMoFEMLoopFiniteElements(dm_rho,"DENSITY_MAP_ELEMENT",&fe_density_approximation); CHKERRQ(ierr);
337 
338  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
339  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
340  ierr = VecGhostUpdateBegin(F,ADD_VALUES,SCATTER_REVERSE);
341  ierr = VecGhostUpdateEnd(F,ADD_VALUES,SCATTER_REVERSE);
342  ierr = VecAssemblyBegin(F); CHKERRQ(ierr);
343  ierr = VecAssemblyEnd(F); CHKERRQ(ierr);
344 
345  // MatView(A,PETSC_VIEWER_DRAW_WORLD);
346  // std::string wait;
347  // std::cin >> wait;
348 
349  {
350  KSP ksp;
351  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp); CHKERRQ(ierr);
352  ierr = KSPSetOperators(ksp,A,A); CHKERRQ(ierr);
353  ierr = KSPSetFromOptions(ksp); CHKERRQ(ierr);
354  ierr = KSPSolve(ksp,F,D); CHKERRQ(ierr); CHKERRQ(ierr);
355  }
356 
357  ierr = VecGhostUpdateBegin(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
358  ierr = VecGhostUpdateEnd(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
359  ierr = m_field.getInterface<VecManager>()->setGlobalGhostVector("DM_RHO",COL,D,INSERT_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
360 
361  DensityMapFe fe_post_process_calulations(m_field);
362 
363  fe_post_process_calulations.getOpPtrVector().push_back(
364  new OpMassCalculationFromApprox("RHO",mass_vec_approx)
365  );
366  ierr = DMoFEMLoopFiniteElements(dm_rho,"DENSITY_MAP_ELEMENT",&fe_post_process_calulations); CHKERRQ(ierr);
367 
368  ierr = VecAssemblyBegin(mass_vec_approx); CHKERRQ(ierr);
369  ierr = VecAssemblyEnd(mass_vec_approx); CHKERRQ(ierr);
370  double mass_approx;
371  ierr = VecSum(mass_vec_approx,&mass_approx); CHKERRQ(ierr);
372  ierr = PetscPrintf(PETSC_COMM_WORLD,"Mass after approximation %0.32g\n",mass_approx); CHKERRQ(ierr);
373  if(is_atom_test) {
374  if( abs(mass_approx - volume) > 1.0e-2) {
375  SETERRQ(PETSC_COMM_SELF,MOFEM_ATOM_TEST_INVALID,"atom test diverged");
376  }
377  }
378  ierr = VecDestroy(&mass_vec_approx); CHKERRQ(ierr);
379 
380  PostProcVolumeOnRefinedMesh post_proc(m_field);
381  ierr = post_proc.generateReferenceElementMesh(); CHKERRQ(ierr);
382  ierr = post_proc.addFieldValuesPostProc("MESH_NODE_POSITIONS"); CHKERRQ(ierr);
383  ierr = post_proc.addFieldValuesPostProc("RHO"); CHKERRQ(ierr);
384  ierr = post_proc.addFieldValuesPostProc("DISTANCE_FROM_SURFACE"); CHKERRQ(ierr);
385  ierr = post_proc.addFieldValuesGradientPostProc("RHO"); CHKERRQ(ierr);
386  ierr = post_proc.addFieldValuesGradientPostProc("DISTANCE_FROM_SURFACE"); CHKERRQ(ierr);
387  ierr = DMoFEMLoopFiniteElements(dm_rho,"DENSITY_MAP_ELEMENT",&post_proc); CHKERRQ(ierr);
388  ierr = post_proc.writeFile("out.h5m"); CHKERRQ(ierr);
389 
390  ierr = VecDestroy(&F); CHKERRQ(ierr);
391  ierr = VecDestroy(&D); CHKERRQ(ierr);
392  ierr = MatDestroy(&A); CHKERRQ(ierr);
393  ierr = DMDestroy(&dm_rho);
394 
395  // Delete element which will be no longer used
396  ierr = m_field.delete_finite_element("DENSITY_MAP_ELEMENT",2); CHKERRQ(ierr);
397  if (m_field.get_comm_rank() == 0)
398  CHKERR moab.write_file("analysis_mesh.h5m");
399 
400 
401  } CATCH_ERRORS;
402 
403  MoFEM::Core::Finalize();
404 
405  return (0);
406 }

Variable Documentation

◆ help

char help[] = "\n"
static

Definition at line 34 of file mass_calculation.cpp.

MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Definition: UnknownInterface.hpp:93
MYPCOMM_INDEX
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215
MoFEM::CoreInterface::loop_dofs
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
EntityHandle
CHKERRQ_MOAB
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:454
MoFEM::CoreInterface::modify_finite_element_add_field_row
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
MoFEM::CoreInterface::get_comm_rank
virtual int get_comm_rank() const =0
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
BoneRemodeling::OpCalulatefRhoAtGaussPts
Assemble local vector containing density data.
Definition: DensityMaps.hpp:769
BoneRemodeling::OpMassCalculationFromApprox
Calculate mass after approximation.
Definition: DensityMaps.hpp:1050
MoFEM::Projection10NodeCoordsOnField
Projection of edge entities with one mid-node on hierarchical basis.
Definition: Projection10NodeCoordsOnField.hpp:24
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
MoFEM::CoreInterface::add_ents_to_field_by_type
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
MoFEM::DMMoFEMAddElement
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:501
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
BoneRemodeling::CommonData
Data structure for storing global and local matrix K and vector f.
Definition: DensityMaps.hpp:667
MoFEM::OpCalculateScalarFieldGradient
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Definition: UserDataOperators.hpp:1294
ROW
@ ROW
Definition: definitions.h:123
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1975
MoFEM::CoreInterface::add_ents_to_finite_element_by_type
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
MoFEM::Exceptions::rval
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:74
BoneRemodeling::SurfaceKDTree
Create KDTree on surface of the mesh and calculate distance.
Definition: DensityMaps.hpp:32
help
static char help[]
Definition: mass_calculation.cpp:34
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEM::CoreInterface::add_finite_element
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
BoneRemodeling::OpVolumeCalculation
Calculate volume of the model.
Definition: DensityMaps.hpp:636
MoFEM::CoreInterface::modify_finite_element_add_field_col
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
MoFEM::CoreInterface::build_finite_elements
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
BoneRemodeling::OpAssmbleRhs
Assemble RHS vector f.
Definition: DensityMaps.hpp:902
BoneRemodeling::CommonData::distanceValue
boost::shared_ptr< VectorDouble > distanceValue
Definition: DensityMaps.hpp:677
MoFEM::CoreInterface::add_field
virtual MoFEMErrorCode add_field(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
MoFEM::OpCalculateScalarFieldValues
Get value at integration points for scalar field.
Definition: UserDataOperators.hpp:82
COL
@ COL
Definition: definitions.h:123
MoFEM::CoreInterface::delete_finite_element
virtual MoFEMErrorCode delete_finite_element(const std::string name, int verb=DEFAULT_VERBOSITY)=0
delete finite element from mofem database
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:47
MoFEM::DMMoFEMCreateMoFEM
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMoFEM.cpp:118
MoFEM::VecManager
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:23
BoneRemodeling::CommonData::distanceGradient
boost::shared_ptr< MatrixDouble > distanceGradient
Definition: DensityMaps.hpp:678
mesh_file_name
char mesh_file_name[255]
Definition: mesh_smoothing.cpp:23
BoneRemodeling::CommonData::globF
Vec globF
Global vector.
Definition: DensityMaps.hpp:675
Range
BoneRemodeling::OpCalculateLhs
Assemble LHS matrix K.
Definition: DensityMaps.hpp:689
MoFEM::CoreTmp< 0 >::Initialize
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
PostProcVolumeOnRefinedMesh
Post processing.
Definition: PostProcOnRefMesh.hpp:955
_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
Definition: MeshsetsManager.hpp:71
BoneRemodeling::DataFromMetaIO
Load data from MetaImage file, translate grayscale values to densities.
Definition: DensityMaps.hpp:205
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1094
BoneRemodeling::OpMassCalculation
Calculate mass before approximation.
Definition: DensityMaps.hpp:957
BLOCKSET
@ BLOCKSET
Definition: definitions.h:148
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:66
MoFEM::MeshsetsManager
Interface for managing meshsets containing materials and boundary conditions.
Definition: MeshsetsManager.hpp:104
MoFEM::CoreInterface::build_fields
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
MoFEM::DeprecatedCoreInterface::seed_ref_level_3D
DEPRECATED MoFEMErrorCode seed_ref_level_3D(const EntityHandle meshset, const BitRefLevel &bit, int verb=-1)
seed 2D entities in the meshset and their adjacencies (only TETs adjacencies) in a particular BitRefL...
Definition: DeprecatedCoreInterface.cpp:18
MoFEM::CoreInterface::modify_finite_element_add_field_data
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_filed)=0
set finite element field data
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
ReactionDiffusionEquation::D
const double D
diffusivity
Definition: reaction_diffusion.cpp:20
_IT_CUBITMESHSETS_FOR_LOOP_
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
Definition: MeshsetsManager.hpp:34
BoneRemodeling::CommonData::globA
Mat globA
Global matrix.
Definition: DensityMaps.hpp:674
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MOFEM_ATOM_TEST_INVALID
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
MoFEM::CoreInterface::build_adjacencies
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
MoFEM::CoreInterface::set_field_order
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
MoFEM::DMoFEMLoopFiniteElements
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:590
BoneRemodeling::DensityMapFe
Definition: DensityMaps.hpp:625
F
@ F
Definition: free_surface.cpp:394
Generated by Doxygen 1.8.17 and hosted at University of Glasgow