OpenFPM_pdata  4.1.0
Project that contain the implementation of distributed structures
 
Loading...
Searching...
No Matches
Vector 4 complex properties

Vector 4 complex properties

This example show how we can use complex properties in a vector

Initialization and vector creation

We first initialize the library and define useful constants

See also
Initialization
// initialize the library
openfpm_init(&argc,&argv);
// Here we define our domain a 2D box with internals from 0 to 1.0 for x and y
Box<2,float> domain({0.0,0.0},{1.0,1.0});
// Here we define the boundary conditions of our problem
size_t bc[2]={PERIODIC,PERIODIC};
// extended boundary around the domain, and the processor domain
Ghost<2,float> g(0.01);
// the scalar is the element at position 0 in the aggregate
constexpr int scalar = 0;
// the vector is the element at position 1 in the aggregate
constexpr int vector = 1;
// the tensor is the element at position 2 in the aggregate
constexpr int point = 2;
// A list1
constexpr int list = 3;
// A listA
constexpr int listA = 4;
// A list of list
constexpr int listlist = 5;
This class represent an N-dimensional box.
Definition Box.hpp:61

We also define a custom structure

// The custom structure
struct A
{
float p1;
int p2;
A() {};
A(float p1, int p2)
:p1(p1),p2(p2)
{}
};

After we initialize the library we can create a vector with complex properties with the following line

vector_dist<2,float, aggregate<float,
float[3],
vd(4096,domain,bc,g);
This class implement the point shape in an N-dimensional space.
Definition Point.hpp:28
Implementation of 1-D std::vector like structure.
Distributed vector.
aggregate of properties, from a list of object if create a struct that follow the OPENFPM native stru...

In this this particular case every particle carry a scalar, a vector in form of float[3], a Point, a list in form of vector of float and a list of custom structures, and a vector of vector. In general particles can have properties of arbitrary complexity.

Warning
For arbitrary complexity mean that we can use any openfpm data structure with and arbitrary nested complexity. For example a openfpm::vector<aggregate<grid_cpu<openfpm::vector<aggregate<double,double[3]>>>,openfpm::vector<float>> is valid
       particle
          *
        vector
         / \
        /   \
     grid    vector<float>
      /\
     /  \
double  double[3]

 * 

Our custom data-structure A is defined below. Note that this data-structure does not have pointers

// The custom structure
struct A
{
float p1;
int p2;
A() {};
A(float p1, int p2)
:p1(p1),p2(p2)
{}
};
Warning
custom data structure are allowed only if they does not have pointer. In case they have pointer we have to define how to serialize our data-structure
See also
Vector 4 property serialization

Assign values to properties

Assign values to properties does not changes, from the simple case. Consider now that each particle has a list, so when we can get the property listA for particle p and resize such list with vd.getProp<listA>(p).resize(...). We can add new elements at the end with vd.getProp<listA>(p).add(...) and get some element of this list with vd.getProp<listA>(p).get(i). More in general vd.getProp<listA>(p) return a reference to the openfpm::vector contained by the particle.

auto it = vd.getDomainIterator();
while (it.isNext())
{
auto p = it.get();
// we define x, assign a random position between 0.0 and 1.0
vd.getPos(p)[0] = (float)rand() / RAND_MAX;
// we define y, assign a random position between 0.0 and 1.0
vd.getPos(p)[1] = (float)rand() / RAND_MAX;
vd.getProp<scalar>(p) = 1.0;
vd.getProp<vector>(p)[0] = 1.0;
vd.getProp<vector>(p)[1] = 1.0;
vd.getProp<vector>(p)[2] = 1.0;
vd.getProp<point>(p).get(0) = 1.0;
vd.getProp<point>(p).get(1) = 1.0;
vd.getProp<point>(p).get(2) = 1.0;
size_t n_cp = (float)10.0 * rand()/RAND_MAX;
vd.getProp<listA>(p).resize(n_cp);
for (size_t i = 0 ; i < n_cp ; i++)
{
vd.getProp<list>(p).add(i + 10);
vd.getProp<list>(p).add(i + 20);
vd.getProp<list>(p).add(i + 30);
vd.getProp<listA>(p).get(i) = A(i+10.0,i+20.0);
}
vd.getProp<listlist>(p).resize(2);
vd.getProp<listlist>(p).get(0).resize(2);
vd.getProp<listlist>(p).get(1).resize(2);
vd.getProp<listlist>(p).get(0).get(0) = 1.0;
vd.getProp<listlist>(p).get(0).get(1) = 2.0;
vd.getProp<listlist>(p).get(1).get(0) = 3.0;
vd.getProp<listlist>(p).get(1).get(1) = 4.0;
// next particle
++it;
}

Mapping and ghost_get

Particles are redistributed across processors all properties are communicated but instead of using map we use map_list that we can use to select properties. A lot of time complex properties can be recomputed and communicate them is not a good idea. The same concept also apply for ghost_get. In general we choose which properties to communicate

See also
Mapping particles
Ghost
// Particles are redistribued across the processors but only the scalar,vector, and point properties
// are transfert
vd.map_list<scalar,vector,point,list,listA,listlist>();
// Synchronize the ghost
vd.ghost_get<scalar,vector,point,listA,listlist>();

Output and VTK visualization

Vector with complex properties can be still be visualized, because unknown properties are automatically excluded

See also
Visualization, write VTK files
vd.write("particles");

Print 4 particles in the ghost area

Here we print that the first 4 particles to show that the list of A and the list of list are filled and the ghosts contain the correct information

size_t fg = vd.size_local();
Vcluster<> & v_cl = create_vcluster();
if (v_cl.getProcessUnitID() == 0)
{
for ( ; fg < vd.size_local()+4 ; fg++)
{
std::cout << "List of A" << std::endl;
for (size_t i = 0 ; i < vd.getProp<listA>(fg).size() ; i++)
std::cout << "Element: " << i << " p1=" << vd.getProp<listA>(fg).get(i).p1 << " p2=" << vd.getProp<listA>(fg).get(i).p2 << std::endl;
std::cout << "List of list" << std::endl;
for (size_t i = 0 ; i < vd.getProp<listlist>(fg).size() ; i++)
{
for (size_t j = 0 ; j < vd.getProp<listlist>(fg).get(i).size() ; j++)
std::cout << "Element: " << i << " " << j << " " << vd.getProp<listlist>(fg).get(i).get(j) << std::endl;
}
}
}
size_t getProcessUnitID()
Get the process unit id.
Implementation of VCluster class.
Definition VCluster.hpp:59

Finalize

At the very end of the program we have always to de-initialize the library

openfpm_finalize();

Full code

#include "Vector/vector_dist.hpp"
int main(int argc, char* argv[])
{
// initialize the library
openfpm_init(&argc,&argv);
// Here we define our domain a 2D box with internals from 0 to 1.0 for x and y
Box<2,float> domain({0.0,0.0},{1.0,1.0});
// Here we define the boundary conditions of our problem
size_t bc[2]={PERIODIC,PERIODIC};
// extended boundary around the domain, and the processor domain
Ghost<2,float> g(0.01);
// the scalar is the element at position 0 in the aggregate
constexpr int scalar = 0;
// the vector is the element at position 1 in the aggregate
constexpr int vector = 1;
// the tensor is the element at position 2 in the aggregate
constexpr int point = 2;
// A list1
constexpr int list = 3;
// A listA
constexpr int listA = 4;
// A list of list
constexpr int listlist = 5;
// The custom structure
struct A
{
float p1;
int p2;
A() {};
A(float p1, int p2)
:p1(p1),p2(p2)
{}
};
vector_dist<2,float, aggregate<float,
float[3],
vd(4096,domain,bc,g);
auto it = vd.getDomainIterator();
while (it.isNext())
{
auto p = it.get();
// we define x, assign a random position between 0.0 and 1.0
vd.getPos(p)[0] = (float)rand() / RAND_MAX;
// we define y, assign a random position between 0.0 and 1.0
vd.getPos(p)[1] = (float)rand() / RAND_MAX;
vd.getProp<scalar>(p) = 1.0;
vd.getProp<vector>(p)[0] = 1.0;
vd.getProp<vector>(p)[1] = 1.0;
vd.getProp<vector>(p)[2] = 1.0;
vd.getProp<point>(p).get(0) = 1.0;
vd.getProp<point>(p).get(1) = 1.0;
vd.getProp<point>(p).get(2) = 1.0;
size_t n_cp = (float)10.0 * rand()/RAND_MAX;
vd.getProp<listA>(p).resize(n_cp);
for (size_t i = 0 ; i < n_cp ; i++)
{
vd.getProp<list>(p).add(i + 10);
vd.getProp<list>(p).add(i + 20);
vd.getProp<list>(p).add(i + 30);
vd.getProp<listA>(p).get(i) = A(i+10.0,i+20.0);
}
vd.getProp<listlist>(p).resize(2);
vd.getProp<listlist>(p).get(0).resize(2);
vd.getProp<listlist>(p).get(1).resize(2);
vd.getProp<listlist>(p).get(0).get(0) = 1.0;
vd.getProp<listlist>(p).get(0).get(1) = 2.0;
vd.getProp<listlist>(p).get(1).get(0) = 3.0;
vd.getProp<listlist>(p).get(1).get(1) = 4.0;
// next particle
++it;
}
// Particles are redistribued across the processors but only the scalar,vector, and point properties
// are transfert
vd.map_list<scalar,vector,point,list,listA,listlist>();
// Synchronize the ghost
vd.ghost_get<scalar,vector,point,listA,listlist>();
vd.write("particles");
size_t fg = vd.size_local();
Vcluster<> & v_cl = create_vcluster();
if (v_cl.getProcessUnitID() == 0)
{
for ( ; fg < vd.size_local()+4 ; fg++)
{
std::cout << "List of A" << std::endl;
for (size_t i = 0 ; i < vd.getProp<listA>(fg).size() ; i++)
std::cout << "Element: " << i << " p1=" << vd.getProp<listA>(fg).get(i).p1 << " p2=" << vd.getProp<listA>(fg).get(i).p2 << std::endl;
std::cout << "List of list" << std::endl;
for (size_t i = 0 ; i < vd.getProp<listlist>(fg).size() ; i++)
{
for (size_t j = 0 ; j < vd.getProp<listlist>(fg).get(i).size() ; j++)
std::cout << "Element: " << i << " " << j << " " << vd.getProp<listlist>(fg).get(i).get(j) << std::endl;
}
}
}
openfpm_finalize();
}