How to use the Spatial Containers
Contents |
Using the Spatial Search Algorithm
Large-scale in Finite Element Methods, Discrete Element Method and Combined Finite-Discrete Element Method simulations involve contact of a large number of separate bodies. Processing contact interaction for all possible contacts would involve a total number of operations proportional to N^2, where N is the total number of separates bodies comprising the problem.
This would be very CPU intensive, and would limit the application that use and need to evaluate
the contact forces ,comprising a very small number (a few thousand) of separates bodies.
To reduce CPU requirements of processing contact interaction, it is necessary to eliminate couples of discrete elements that are far from each other and are not in contact. A set procedures designed to detect
bodies that are close to each other is usually called a contact detection
algorithm, or sometimes a contact search algorithm.
With this purpose was created in Kratos Program a contact search algorithm based on spatial decomposition,
so it can be used for any application and is unique in that it was implemented in a generic way.
Therefore it requires a contact search algorithm having the following characteristics:
- Robust,
- CPU efficient,
- RAM efficient
The space decomposition implemented in Kratos are:
- Trees : Kd-Trees, Octrees
- Bins
There are two kind of searches algorithms in Kratos:
- Points Searches: Used commonly for searching near points in a radius and can be used for
searches contacts for spheric and circle geometries. Use the Kd_tree,Octree and bins for spatial decomposition.
- Object Searches: The searches is based in the concept of objects. Is a generic class, thats mean, we can use either objects type even points. For now, use the bins for spatial decomposition. It could be interesting proof hybrid methods would be a
and future investigation lines.
Spatial decomposition based in bins. Main methods
Exist implemented in Kratos two versions of domain decomposition based in bins: dynamic and static bins. The difference between them is the dynamic allocation of memory. However both define the methods. The most common methods used are presented and described below:
- Bins based in points :
*SearchNearestPoint : Calculates the nearest point of a point.
*SearchInRadius : Get the list of point in a certain radius.
- Bins based in Objects :
* SearchObjects : Get the list of contacts of a certain objects even itself.
* SearchObjectsInner : The same as above, just without return itself.
* SearchConcts : Get the list of all pair contacts.
Using the search algorithm based in objects
Before start using this algorithm, is necessary to look some aspects:
- The dynamic or static bins based in object is a generic class. That means, a template
parameter is used. This template parameter is a user configure file.
- The user need to create and extra file (user configure file). This fail describes the types
and containers used by the bins. It has also the intersection method for contacts. As there are many geometric figures, so the intersection algorithm differs each other. Most of intersection algorithm can be found in computational geometric books and more of them is written c++ or c.
- Pointers and iterator: the bins works by pointers and references. So, in the configure file,
we need to know what is the class, a pointer to the class, and iterator of pointers class containers.
With this in hand, below describe how to create a configuration file and how to use the contact search algorithm.
The User Configure File
It say that is the main file created by the user, which defines the types, and operations containers necessary to make the spatial decomposition and the search for contacts. Is then the parameter that defines Bins.
This consists of four parts:
- Definition of the types and their containers.
- Function that calculates bounding box
- Function of intersection between objects
- Function of intercession between the object and cells.
An Example Configure File
Suppose that we have a set of discrete elements like a spheres. In our file called "discrete_configure_file.h" need to create the class, types and methods before mentioned. Is useful to create this file using the template of Kratos found in kratos/templates directory. Only we need to do is changes the default names and directives by the name of our file. Normally this file is stored in the custom_utilities directory of the application.
- Change the name of directives
#if !defined(KRATOS_DISCRETE_PARTICLE_CONFIGURE_INCLUDED) #define KRATOS_DISCRETE_PARTICLE_CONFIGURE_INCLUDED
- Change the name of class
///@}
///@name Kratos Classes
///@{
template <std::size_t TDimension>
class DiscreteParticleConfigure{
- Typing the typenames, objects and container of the particle
///@}
///@name Type Definitions
///@{
typedef Point<Dimension, double> PointType;
typedef std::vector<double>::iterator DistanceIteratorType;
typedef ModelPart::ElementsContainerType::ContainerType ContainerType;
typedef ContainerType::value_type PointerType;
typedef ContainerType::iterator IteratorType;
typedef ModelPart::ElementsContainerType::ContainerType ResultContainerType;
typedef ResultContainerType::value_type ResultPointerType;
typedef ResultContainerType::iterator ResultIteratorType;
typedef ContactPair<PointerType> ContactPairType;
typedef std::vector<ContactPairType> ContainerContactType;
typedef ContainerContactType::iterator IteratorContactType;
typedef ContainerContactType::value_type PointerContactType;
typedef std::vector<PointerType>::iterator PointerTypeIterator;
|
Another types and container different to the model part can be used. |
Note:
- Computing the bounding box
static inline void CalculateBoundingBox(const PointerType& rObject, PointType& rLowPoint, PointType& rHighPoint)
{
rHighPoint = rLowPoint = rObject->GetGeometry().GetPoint(0);
double radius = rObject->GetGeometry()(0)->GetSolutionStepValue(RADIUS);
for(std::size_t i = 0; i < 3; i++){
rLowPoint[i] += -radius;
rHighPoint[i] += radius;
}
}
|
It is necessary the name static inline void CalculateBoundingBox is well type, beacause the function is called by the bins. |
- The Intersection function between objects.
static inline bool Intersection(const PointerType& rObj_1, const PointerType& rObj_2){
array_1d<double, 3> rObj_2_to_rObj_1 = rObj_1->GetGeometry().GetPoint(0) - rObj_2->GetGeometry().GetPoint(0);
double distance_2 = inner_prod(rObj_2_to_rObj_1, rObj_2_to_rObj_1);
const double& radius_1 = rObj_1->GetGeometry()(0)->GetSolutionStepValue(RADIUS);
const double& radius_2 = rObj_2->GetGeometry()(0)->GetSolutionStepValue(RADIUS);
double radius_sum = radius_1 + radius_2;
bool intersect = (distance_2 - radius_sum * radius_sum) <= 0;
return intersect;
}
- Intersection cell and objects
static inline bool IntersectionBox(const PointerType& rObject, const PointType& rLowPoint, const PointType& rHighPoint){
array_1d<double, 3> center_of_particle = rObject->GetGeometry().GetPoint(0);
const double& radius = rObject->GetGeometry()(0)->GetSolutionStepValue(RADIUS);
bool intersect = (rLowPoint[0] - radius <= center_of_particle[0] && rLowPoint[1] - radius <= center_of_particle[1] && rLowPoint[2] - radius <= center_of_particle[2] &&
rHighPoint[0] + radius >= center_of_particle[0] && rHighPoint[1] + radius >= center_of_particle[1] && rHighPoint[2] + radius >= center_of_particle[2]);
return intersect;
}
|
Intersection cell and objects: If we do not have this function a simple return true is enougth. |
Calling the Constructor
Once the configure file is done, we are able to use Contact Search Algorithm, based in a domain decomposition. In the fail where it will be called , we need to add in the directives the corresponding files, those are the configure file and the bins. Maybe is required some redefinition of the typedef defined above in the configure file.
static const std::size_t space_dim = 2;
typedef DiscreteConfigure<space_dim> Configure;
typedef Configure::PointType PointType;
typedef PointType::CoordinatesArrayType CoordinatesArrayType;
typedef Configure::ContainerType ContainerType;
typedef Configure::PointerType PointerType;
typedef Configure::IteratorType IteratorType;
typedef Configure::ResultContainerType ResultContainerType;
typedef Configure::ResultPointerType ResultPointerType;
typedef Configure::ResultIteratorType ResultIteratorType;
typedef Configure::ContactPairType ContactPairType;
typedef Configure::ContainerContactType ContainerContactType;
typedef Configure::IteratorContactType IteratorContactType;
typedef Configure::PointerContactType PointerContactType;
typedef Configure::PointerTypeIterator PointerTypeIterator;
typedef ContainerContactType ContainerContactPair;
typedef IteratorContactType IteratorContainerContactPair;
typedef PointerContactType PointerContainerContactPair;
After that we call the constructor as follow:
- Using Dynamic Bins
IteratorType it_begin = mBoundaryElements.begin(); IteratorType it_end = mBoundaryElements.end(); BinsObjectDynamic<Configure> rBinsObjectDynamic(it_begin, it_end);
- Using Static Bins
BinsObjectStatic<Configure> Bins(it_begin, it_end); const std::size_t MaxNumberOfResults = 1000; std::size_t NumberOfResults = 0; ResultIteratorType begin;
|
Note that we can use the dynamic or static bins. The difference between them is the dynamic allocation of memory. Naturally, the static bins is faster because the we know size of containers. |
- Lastly, we can get the list of pairs of contacts
rBins.SearchContact(mPairContacts);
- However, is faster to do a loop and get the bodies that are in contact with a particular element.
For(IteratorType it =it_begin; it!=it_end; it++)
rBinsObjectDynamic.SearchObjectsInner(*it, Result);
where Result is the list of contact for a particular element. A complete configure file is found in /kratos/kratos/utilities/spatial_container_configure.h. Is a simple file, but the difficult part is the contact geometric calculus. However, most of intersection algorithms are in the literature and are written in a c++, matlab or c code.
References
- Munjiza. Ante, The Combined Finite-Discrete Element Method, John Wiley & Sons, Ltd. 2004
