v0.14.0
Public Member Functions | Public Attributes | List of all members
BoneRemodeling::SurfaceKDTree Struct Reference

Create KDTree on surface of the mesh and calculate distance. More...

#include <users_modules/bone_remodelling/src/DensityMaps.hpp>

Collaboration diagram for BoneRemodeling::SurfaceKDTree:
[legend]

Public Member Functions

 SurfaceKDTree (MoFEM::Interface &m_field)
 
PetscErrorCode takeASkin (Range &volume)
 Take a skin from a mesh. More...
 
PetscErrorCode buildTree ()
 Build tree. More...
 
PetscErrorCode findClosestPointToTheSurface (const double x, const double y, const double z, double &dx, double &dy, double &dz)
 Find point on surface which is closet to given point coordinates. More...
 
PetscErrorCode setDistanceFromSurface (const std::string field_name)
 Set distance to vertices on mesh. More...
 

Public Attributes

MoFEM::Interface & mField
 
EntityHandle kdTreeRootMeshset
 
AdaptiveKDTree kdTree
 
Range sKin
 

Detailed Description

Create KDTree on surface of the mesh and calculate distance.

Definition at line 32 of file DensityMaps.hpp.

Constructor & Destructor Documentation

◆ SurfaceKDTree()

BoneRemodeling::SurfaceKDTree::SurfaceKDTree ( MoFEM::Interface &  m_field)
inline

Definition at line 39 of file DensityMaps.hpp.

39  :
40  mField(m_field),
41  kdTree(&m_field.get_moab()) {
42  }

Member Function Documentation

◆ buildTree()

PetscErrorCode BoneRemodeling::SurfaceKDTree::buildTree ( )
inline

Build tree.

Returns
Error code

Definition at line 72 of file DensityMaps.hpp.

72  {
73  PetscFunctionBegin;
74  rval = kdTree.build_tree(sKin,&kdTreeRootMeshset); CHKERRQ_MOAB(rval);
75  PetscFunctionReturn(0);
76  }

◆ findClosestPointToTheSurface()

PetscErrorCode BoneRemodeling::SurfaceKDTree::findClosestPointToTheSurface ( const double  x,
const double  y,
const double  z,
double &  dx,
double &  dy,
double &  dz 
)
inline

Find point on surface which is closet to given point coordinates.

Parameters
x1st coordinate of point
y2nd coordinate of point
z3rd coordinate of point
dxcoordinate of distance vector
dycoordinate of distance vector
zcoordinate of distance vector
Returns
Error code

Definition at line 88 of file DensityMaps.hpp.

95  {
96  PetscFunctionBegin;
97  const double coords[] = {x,y,z};
98  double closest_point[3];
99  EntityHandle triangle;
100  rval = kdTree.closest_triangle(
101  kdTreeRootMeshset,
102  coords,
103  closest_point,
104  triangle
105  ); CHKERRQ_MOAB(rval);
106  cblas_daxpy(3,-1,coords,1,closest_point,1);
107  dx = closest_point[0];
108  dy = closest_point[1];
109  dz = closest_point[2];
110  PetscFunctionReturn(0);
111  }

◆ setDistanceFromSurface()

PetscErrorCode BoneRemodeling::SurfaceKDTree::setDistanceFromSurface ( const std::string  field_name)
inline

Set distance to vertices on mesh.

Parameters
field_nameName of field for distances
Returns
Error Code

Definition at line 119 of file DensityMaps.hpp.

119  {
120  PetscFunctionBegin;
121  auto field_ents = mField.get_field_ents();
122  auto field_bit_number = mField.get_field_bit_number(field_name);
123  auto it = field_ents->lower_bound(
124  FieldEntity::getLoBitNumberUId(field_bit_number));
125  auto hi_it = field_ents->upper_bound(
126  FieldEntity::getHiBitNumberUId(field_bit_number));
127  // Get values at nodes first
128  for(;it!=hi_it;it++) {
129  EntityType type = it->get()->getEntType();
130  if(type!=MBVERTEX) {
131  continue;
132  }
133  EntityHandle ent = it->get()->getEnt();
134  double coords[3];
135  rval = mField.get_moab().get_coords(&ent,1,coords); CHKERRQ_MOAB(rval);
136  double delta[3];
137  ierr = findClosestPointToTheSurface(
138  coords[0],coords[1],coords[2],delta[0],delta[1],delta[2]
139  ); CHKERRQ(ierr);
140  // Get vector of DOFs on vertex
141  VectorAdaptor dofs = it->get()->getEntFieldData();
142  // Set values
143  dofs[0] = cblas_dnrm2(3,delta,1);
144  }
145  // Get values at edges and project those values on approximation base
146  it = field_ents->lower_bound(
147  FieldEntity::getLoBitNumberUId(field_bit_number));
148  for(;it!=hi_it;it++) {
149  EntityType type = it->get()->getEntType();
150  if(type!=MBEDGE) {
151  continue;
152  }
153  const EntityHandle *conn;
154  int num_ndodes;
155  EntityHandle ent = it->get()->getEnt();
156  rval = mField.get_moab().get_connectivity(ent,conn,num_ndodes,true); CHKERRQ_MOAB(rval);
157  double coords[3*num_ndodes];
158  rval = mField.get_moab().get_coords(conn,num_ndodes,coords); CHKERRQ_MOAB(rval);
159  cblas_daxpy(3,1,&coords[3],1,coords,1);
160  cblas_dscal(3,0.5,coords,1);
161  double delta[3];
162  ierr = findClosestPointToTheSurface(
163  coords[0],coords[1],coords[2],delta[0],delta[1],delta[2]
164  ); CHKERRQ(ierr);
165 
166  auto conn_it0 = field_ents->find(
167  FieldEntity::getLoLocalEntityBitNumber(field_bit_number, conn[0]));
168  auto conn_it1 = field_ents->find(
169  FieldEntity::getLoLocalEntityBitNumber(field_bit_number, conn[1]));
170 
171  if(conn_it0==field_ents->end()) {
172  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"entity not found");
173  }
174  if (conn_it1 == field_ents->end()) {
175  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"entity not found");
176  }
177 
178  VectorAdaptor vec_conn0 = conn_it0->get()->getEntFieldData();
179  VectorAdaptor vec_conn1 = conn_it1->get()->getEntFieldData();
180  // cerr << vec_conn0 << " " << vec_conn1 << endl;
181  double ave_delta = 0.5*(vec_conn0[0]+vec_conn1[0]);
182  double edge_shape_function_val = 0.25;
183  if(it->get()->getApproxBase()==AINSWORTH_LEGENDRE_BASE) {
184  } else if(it->get()->getApproxBase()==AINSWORTH_LOBATTO_BASE) {
185  edge_shape_function_val *= LOBATTO_PHI0(0);
186  } else {
187  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Base not implemented");
188  }
189  // Get vector of DOFs on vertex
190  VectorAdaptor dofs = it->get()->getEntFieldData();
191  // cerr << mid_val_dx << " " << ave_dx << endl;
192  // Set values
193  dofs[0] = (cblas_dnrm2(3,delta,1)-ave_delta)/edge_shape_function_val;
194  }
195  PetscFunctionReturn(0);
196  }

◆ takeASkin()

PetscErrorCode BoneRemodeling::SurfaceKDTree::takeASkin ( Range &  volume)
inline

Take a skin from a mesh.

Parameters
volumeVolume elements in the mesh
Returns
Error code

Definition at line 55 of file DensityMaps.hpp.

55  {
56  PetscFunctionBegin;
57  Skinner skin(&mField.get_moab());
58  rval = skin.find_skin(0,volume,false,sKin); CHKERRQ_MOAB(rval);
59  // Only for debugging
60  // EntityHandle meshset;
61  // rval = mField.get_moab().create_meshset(MESHSET_SET,meshset);
62  // rval = mField.get_moab().add_entities(meshset,sKin); CHKERRQ(rval);
63  // rval = mField.get_moab().write_file("skin.vtk","VTK","",&meshset,1); CHKERRQ(rval);
64  PetscFunctionReturn(0);
65  }

Member Data Documentation

◆ kdTree

AdaptiveKDTree BoneRemodeling::SurfaceKDTree::kdTree

Definition at line 36 of file DensityMaps.hpp.

◆ kdTreeRootMeshset

EntityHandle BoneRemodeling::SurfaceKDTree::kdTreeRootMeshset

Definition at line 35 of file DensityMaps.hpp.

◆ mField

MoFEM::Interface& BoneRemodeling::SurfaceKDTree::mField

Definition at line 34 of file DensityMaps.hpp.

◆ sKin

Range BoneRemodeling::SurfaceKDTree::sKin

Definition at line 47 of file DensityMaps.hpp.

The documentation for this struct was generated from the following file:
EntityHandle
CHKERRQ_MOAB
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:467
MoFEM::CoreInterface::get_field_bit_number
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
BoneRemodeling::SurfaceKDTree::kdTreeRootMeshset
EntityHandle kdTreeRootMeshset
Definition: DensityMaps.hpp:35
MoFEM::Exceptions::rval
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:74
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
delta
static constexpr double delta
Definition: prism_elements_from_surface.cpp:18
convert.type
type
Definition: convert.py:64
BoneRemodeling::SurfaceKDTree::kdTree
AdaptiveKDTree kdTree
Definition: DensityMaps.hpp:36
MoFEM::CoreInterface::get_field_ents
virtual const FieldEntity_multiIndex * get_field_ents() const =0
Get the field ents object.
AINSWORTH_LOBATTO_BASE
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
BoneRemodeling::SurfaceKDTree::mField
MoFEM::Interface & mField
Definition: DensityMaps.hpp:34
BoneRemodeling::SurfaceKDTree::sKin
Range sKin
Definition: DensityMaps.hpp:47
MoFEM::Types::VectorAdaptor
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:115
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
LOBATTO_PHI0
#define LOBATTO_PHI0(x)
Definitions taken from Hermes2d code.
Definition: base_functions.h:126
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
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
BoneRemodeling::SurfaceKDTree::findClosestPointToTheSurface
PetscErrorCode findClosestPointToTheSurface(const double x, const double y, const double z, double &dx, double &dy, double &dz)
Find point on surface which is closet to given point coordinates.
Definition: DensityMaps.hpp:88
Generated by Doxygen 1.8.17 and hosted at University of Glasgow