OpenFPM_pdata  4.1.0
Project that contain the implementation of distributed structures
 
Loading...
Searching...
No Matches
staggered_dist_grid_util.hpp
1/*
2 * staggered_util.hpp
3 *
4 * Created on: Aug 19, 2015
5 * Author: i-bird
6 */
7
8#ifndef SRC_GRID_STAGGERED_DIST_GRID_UTIL_HPP_
9#define SRC_GRID_STAGGERED_DIST_GRID_UTIL_HPP_
10
11#include "util/common.hpp"
12#include "VTKWriter/VTKWriter.hpp"
13#include "util/convert.hpp"
14
15
23template<typename ele, typename vtk, bool has_attributes=has_attributes<ele>::value>
25{
33 vtk_write(vtk vtk_w, const std::string output, const size_t i)
34 {
35 vtk_w.write(output + "_" + ele::attributes::name[i] + ".vtk",ele::attributes::name[i]);
36 }
37};
38
46template<typename ele, typename vtk>
47struct vtk_write<ele,vtk,false>
48{
56 vtk_write(vtk vtk_w, const std::string output, const size_t i)
57 {
58 vtk_w.write(output + "_" + std::to_string(i) + ".vtk","attr" + std::to_string(i));
59 }
60};
61
62
66template<typename T>
67struct extends
68{
74 static inline size_t mul()
75 {
76 return 1;
77 }
78
84 static inline size_t dim()
85 {
86 return 0;
87 }
88};
89
91template<typename T,size_t N1>
92struct extends<T[N1]>
93{
99 static inline size_t mul()
100 {
101 return N1;
102 }
103
109 static inline size_t dim()
110 {
111 return 1;
112 }
113};
114
116template<typename T,size_t N1,size_t N2>
117struct extends<T[N1][N2]>
118{
124 static inline size_t mul()
125 {
126 return N1 * N2;
127 }
128
134 static inline size_t dim()
135 {
136 return 2;
137 }
138};
139
141template<typename T,size_t N1,size_t N2,size_t N3>
142struct extends<T[N1][N2][N3]>
143{
145 static inline size_t mul()
146 {
147 return N1 * N2 * N3;
148 }
149
155 static inline size_t dim()
156 {
157 return 3;
158 }
159};
160
162template<typename T,size_t N1,size_t N2,size_t N3,size_t N4>
163struct extends<T[N1][N2][N3][N4]>
164{
166 static inline size_t mul()
167 {
168 return N1 * N2 * N3 * N4;
169 }
170
176 static inline size_t dim()
177 {
178 return 4;
179 }
180};
181
183template<typename T,size_t N1,size_t N2,size_t N3,size_t N4,size_t N5>
184struct extends<T[N1][N2][N3][N4][N5]>
185{
187 static inline size_t mul()
188 {
189 return N1 * N2 * N3 * N4 * N5;
190 }
191
197 static inline size_t dim()
198 {
199 return 5;
200 }
201};
202
204template<typename T,size_t N1,size_t N2,size_t N3,size_t N4,size_t N5, size_t N6>
205struct extends<T[N1][N2][N3][N4][N5][N6]>
206{
208 static inline size_t mul()
209 {
210 return N1 * N2 * N3 * N4 * N5 * N6;
211 }
212
218 static inline size_t dim()
219 {
220 return 6;
221 }
222};
223
225template<typename T,size_t N1,size_t N2,size_t N3,size_t N4,size_t N5, size_t N6, size_t N7>
226struct extends<T[N1][N2][N3][N4][N5][N6][N7]>
227{
229 static inline size_t mul()
230 {
231 return N1 * N2 * N3 * N4 * N5 * N6 * N7;
232 }
233
239 static inline size_t dim()
240 {
241 return 7;
242 }
243};
244
246template<typename T,size_t N1,size_t N2,size_t N3,size_t N4,size_t N5, size_t N6, size_t N7, size_t N8>
247struct extends<T[N1][N2][N3][N4][N5][N6][N7][N8]>
248{
254 static inline size_t mul()
255 {
256 return N1 * N2 * N3 * N4 * N5 * N6 * N7 * N8;
257 }
258
264 static inline size_t dim()
265 {
266 return 8;
267 }
268};
269
271template<typename T,size_t N1,size_t N2,size_t N3,size_t N4,size_t N5, size_t N6, size_t N7, size_t N8, size_t N9>
272struct extends<T[N1][N2][N3][N4][N5][N6][N7][N8][N9]>
273{
279 static inline size_t mul()
280 {
281 return N1 * N2 * N3 * N4 * N5 * N6 * N7 * N8 * N9;
282 }
283
289 static inline size_t dim()
290 {
291 return 9;
292 }
293};
294
296template<typename T,size_t N1,size_t N2,size_t N3,size_t N4,size_t N5, size_t N6, size_t N7, size_t N8, size_t N9, size_t N10>
297struct extends<T[N1][N2][N3][N4][N5][N6][N7][N8][N9][N10]>
298{
304 static inline size_t mul()
305 {
306 return N1 * N2 * N3 * N4 * N5 * N6 * N7 * N8 * N9 * N10;
307 }
308
314 static inline size_t dim()
315 {
316 return 10;
317 }
318};
319
321
330template<typename T>
332{
344 template<unsigned int p_val, typename sg, typename v_g> static inline void write(sg & st_g, v_g & vg,size_t lg)
345 {
346 // Add a grid;
347 vg.add();
348 size_t k = vg.size() - 1;
349
350 // Get the source and destination grid
351 auto & g_src = st_g.get_loc_grid(lg);
352 auto & g_dst = vg.get(k);
353
354 // Set dimensions and memory
355 g_dst.resize(g_src.getGrid().getSize());
356
357 // copy
358
359 auto it = vg.get(k).getIterator();
360
361 while(it.isNext())
362 {
363 g_dst.template get<0>(it.get()) = g_src.template get<p_val>(it.get());
364
365 ++it;
366 }
367 }
368};
369
376template<typename T,size_t N1>
377struct write_stag<T[N1]>
378{
390 template<unsigned int p_val, typename sg, typename v_g> static inline void write(sg & st_g, v_g & vg,size_t lg)
391 {
392 for (size_t i = 0 ; i < N1 ; i++)
393 {
394 // Add a grid;
395 vg.add();
396 size_t k = vg.size() - 1;
397
398 // Get the source and destination grid
399 auto & g_src = st_g.get_loc_grid(lg);
400 auto & g_dst = vg.get(k);
401
402 // Set dimensions and memory
403 g_dst.resize(g_src.getGrid().getSize());
404
405 auto it = vg.get(k).getIterator();
406
407 while(it.isNext())
408 {
409 g_dst.template get<0>(it.get()) = g_src.template get<p_val>(it.get())[i];
410
411 ++it;
412 }
413 }
414 }
415};
416
418template<typename T,size_t N1,size_t N2>
419struct write_stag<T[N1][N2]>
420{
432 template<unsigned int p_val, typename sg, typename v_g> static inline void write(sg & st_g, v_g & vg,size_t lg)
433 {
434 for (size_t i = 0 ; i < N1 ; i++)
435 {
436 for (size_t j = 0 ; j < N2 ; j++)
437 {
438 // Add a grid;
439 vg.add();
440 size_t k = vg.size() - 1;
441
442 // Set dimensions and memory
443 vg.get(k).resize(st_g.get_loc_grid(lg).getGrid().getSize());
444
445 // copy
446 auto & g_src = st_g.get_loc_grid(lg);
447 auto & g_dst = vg.get(k);
448 auto it = vg.get(k).getIterator();
449
450 while(it.isNext())
451 {
452 g_dst.template get<0>(it.get()) = g_src.template get<p_val>(it.get())[i][j];
453
454 ++it;
455 }
456 }
457 }
458 }
459};
460
462
475template<unsigned int dim, typename v, bool has_pM = has_posMask<v>::value>
477{
479 // within the cell
480 openfpm::vector<comb<dim>> (& pos_prp)[boost::fusion::result_of::size<v>::type::value];
481
482public:
483
489 stag_set_position( openfpm::vector<comb<dim>> (& pos_prp)[boost::fusion::result_of::size<v>::type::value])
490 :pos_prp(pos_prp)
491 {}
492
498 template<typename T>
499 void operator()(T& t) const
500 {
501 // This is the type of the object we have to copy
502 typedef typename boost::mpl::at<v,typename boost::mpl::int_<T::value>>::type prop;
503
504 bool val = prop::stag_pos_mask[T::value];
505
506 if (val == false)
507 return;
508
509 // Dimension of the object
510 size_t dim_prp = extends<prop>::dim();
511
512 // It is a scalar
513 if (dim_prp == 0)
514 {
515 comb<dim> c;
516 c.zero();
517
518 // It stay in the center
519 pos_prp[T::value].add(c);
520 }
521 else if (dim_prp == 1)
522 {
523 // It stay on the object of dimension dim-1 (Negative part)
524 for (size_t i = 0 ; i < dim ; i++)
525 {
526 comb<dim> c;
527 c.zero();
528 c.value(i) = -1;
529
530 pos_prp[T::value].add(c);
531 }
532 }
533 else if (dim_prp == 2)
534 {
535 // Create an hypercube object
536 HyperCube<dim> hyp;
537
538 // Diagonal part live in
539 for (size_t i = 0 ; i < dim ; i++)
540 {
541 comb<dim> c1 = pos_prp[T::value-1].get(i);
542 for (size_t j = 0 ; j < dim ; j++)
543 {
544 comb<dim> c2;
545 c2.zero();
546 c2.value(i) = 1;
547
548 comb<dim> c_res = (c1 + c2) & 0x1;
549
550 pos_prp[T::value].add(c_res);
551 }
552 }
553 }
554 else if (dim_prp > 2)
555 {
556 std::cerr << __FILE__ << ":" << __LINE__ << " Tensor of rank bigger than 2 are not supported";
557 }
558 }
559};
560
562
573template<unsigned int dim, typename v>
574class stag_set_position<dim,v,false>
575{
576private:
577
579 // within the cell
580 openfpm::vector<comb<dim>> (& pos_prp)[boost::fusion::result_of::size<v>::type::value];
581
582
583public:
584
590 stag_set_position( openfpm::vector<comb<dim>> (& pos_prp)[boost::fusion::result_of::size<v>::type::value])
591 :pos_prp(pos_prp)
592 {}
593
599 template<typename T>
600 void operator()(T& t) const
601 {
602 // This is the type of the object we have to copy
603 typedef typename boost::mpl::at<v,typename boost::mpl::int_<T::value>>::type prop;
604
605 // Dimension of the object
606 size_t dim_prp = extends<prop>::dim();
607
608 // It is a scalar
609 if (dim_prp == 0)
610 {
611 comb<dim> c;
612 c.zero();
613
614 // It stay in the center
615 pos_prp[T::value].add(c);
616 }
617 else if (dim_prp == 1)
618 {
619 // It stay on the object of dimension dim-1 (Negative part)
620 for (size_t i = 0 ; i < dim ; i++)
621 {
622 comb<dim> c;
623 c.zero();
624 c.getComb()[i] = -1;
625
626 pos_prp[T::value].add(c);
627 }
628 }
629 else if (dim_prp == 2)
630 {
631 // Diagonal part live in
632 for (size_t i = 0 ; i < dim ; i++)
633 {
634 comb<dim> c1 = pos_prp[T::value-1].get(i);
635 for (size_t j = 0 ; j < dim ; j++)
636 {
637 comb<dim> c2;
638 c2.zero();
639 c2.getComb()[j] = 1;
640
641 comb<dim> c_res = (c2 + c1).flip();
642
643 pos_prp[T::value].add(c_res);
644 }
645 }
646 }
647 else if (dim_prp > 2)
648 {
649 std::cerr << __FILE__ << ":" << __LINE__ << " Tensor of rank bigger than 2 are not supported";
650 }
651 }
652};
653
660template<unsigned int dim, typename st_grid, typename St>
662{
663
664 size_t p_id;
665
666 // staggered grid to write
667 st_grid & st_g;
668
669public:
670
677 stag_create_and_add_grid(st_grid & st_g, size_t p_id)
678 :p_id(p_id),st_g(st_g)
679 {}
680
681 template<unsigned int p_val> void out_normal()
682 {
683 // property type
684 typedef typename boost::mpl::at< typename st_grid::value_type::type , typename boost::mpl::int_<p_val> >::type ele;
685
686 // create an openfpm format object from the property type
688
690
691 // Create a vector of grids
692
693 openfpm::vector< grid_cpu<dim, d_object > > vg(st_g.getN_loc_grid());
694
695 // for each domain grid
696 for (size_t i = 0 ; i < vg.size() ; i++)
697 {
698 // Set dimansions and memory
699 vg.get(i).resize(st_g.get_loc_grid(i).getGrid().getSize());
700
701 auto & g_src = st_g.get_loc_grid(i);
702 auto & g_dst = vg.get(i);
703
704 auto it = vg.get(i).getIterator();
705
706 while(it.isNext())
707 {
708 object_si_d< decltype(g_src.get_o(it.get())),decltype(g_dst.get_o(it.get())) ,OBJ_ENCAP,p_val>(g_src.get_o(it.get()),g_dst.get_o(it.get()));
709
710 ++it;
711 }
712
713 Point<dim,St> offset = st_g.getOffset(i);
714 Point<dim,St> spacing = st_g.getSpacing();
715 Box<dim,size_t> dom = st_g.getDomain(i);
716
717 vtk_w.add(g_dst,offset,spacing,dom);
718 }
719
720 vtk_w.write("vtk_grids_st_" + std::to_string(p_id) + "_" + std::to_string(p_val) + ".vtk",st_g.getPropNames(),"grids",file_type::BINARY);
721 }
722
723 template<unsigned int p_val> void out_staggered()
724 {
725 // property type
726 typedef typename boost::mpl::at< typename st_grid::value_type::type , typename boost::mpl::int_<p_val> >::type ele;
727
728 // Eliminate the extends
729 typedef typename std::remove_all_extents<ele>::type r_ele;
730
731 // create an openfpm format object from the property type
733
735
736 // Create a vector of grids
738 vg.reserve(st_g.getN_loc_grid() * extends<ele>::mul());
739
740 size_t k = 0;
741
742 // for each domain grid
743 for (size_t i = 0 ; i < st_g.getN_loc_grid() ; i++)
744 {
745 write_stag<ele>::template write<p_val, st_grid,openfpm::vector< grid_cpu<dim, d_object > > >(st_g,vg,i);
746
747 // for each component
748 for ( ; k < vg.size() ; k++)
749 {
750 Point<dim,St> offset = st_g.getOffset(i);
751 Point<dim,St> spacing = st_g.getSpacing();
752 Box<dim,size_t> dom = st_g.getDomain(i);
753
754 vtk_w.add(i,vg.get(k),offset,spacing,dom,st_g.c_prp[p_val].get(k));
755 }
756
757 k = vg.size();
758 }
759
760 vtk_write<typename st_grid::value_type,VTKWriter<boost::mpl::pair<grid_cpu<dim, d_object >,St>,VECTOR_ST_GRIDS>> v(vtk_w,"vtk_grids_st_" + std::to_string(p_id),p_val);
761 }
762
764 template<typename T>
765 void operator()(T& t)
766 {
767 if (st_g.is_staggered_prop(T::value) == false)
768 out_normal<T::value>();
769 else
770 out_staggered<T::value>();
771 }
772};
773
789template<typename Eqs_sys, typename it1_type, typename it2_type> bool checkIterator(const it1_type & it1, const it2_type & it2)
790{
791#ifdef SE_CLASS1
792
793 grid_key_dx<Eqs_sys::dims> it1_k = it1.getStop() - it1.getStart();
794 grid_key_dx<Eqs_sys::dims> it2_k = it2.getStop() - it2.getStart();
795
796 for (size_t i = 0 ; i < Eqs_sys::dims ; i++)
797 {
798 if (it1_k.get(i) != it2_k.get(i))
799 {
800 std::cerr << __FILE__ << ":" << __LINE__ << " error src iterator and destination iterator does not match in size\n";
801 return false;
802 }
803 }
804
805 return true;
806#else
807
808 return true;
809
810#endif
811}
812
825template<unsigned int dim, unsigned int n_prop, typename v_prp_id, typename v_prp_type>
827{
828/*#ifdef SE_CLASS3
829
830 // number of properties we are processing
831 typedef boost::mpl::size<v_prp_id> v_size_old;
832
833 typedef boost::mpl::int_<v_size_old::value-3> v_size;
834
835#else*/
836
837 // number of properties we are processing
838 typedef boost::mpl::size<v_prp_id> v_size;
839
840//#endif
841
842 // interpolation points for each property
843 openfpm::vector<std::vector<comb<dim>>> (& interp_pts)[v_size::value];
844
845 // staggered position for each property
846 const openfpm::vector<comb<dim>> (&stag_pos)[n_prop];
847
856 inline interp_points(openfpm::vector<std::vector<comb<dim>>> (& interp_pts)[v_size::value],const openfpm::vector<comb<dim>> (&stag_pos)[n_prop])
857 :interp_pts(interp_pts),stag_pos(stag_pos){};
858
860 template<typename T>
861 inline void operator()(T& t)
862 {
863 // This is the type of the object we have to copy
864 typedef typename boost::mpl::at_c<v_prp_type,T::value>::type prp_type;
865 typedef typename boost::mpl::at<v_prp_id,T>::type p_id;
866
867 interp_pts[T::value].resize(stag_pos[p_id::value].size());
868
869 for (size_t i = 0 ; i < stag_pos[p_id::value].size() ; i++)
870 {
871 // Create the interpolation points
872 interp_pts[T::value].get(i) = SubHyperCube<dim,dim - std::rank<prp_type>::value>::getCombinations_R(stag_pos[p_id::value].get(i),0);
873
874 // interp_point are -1,0,1, map the -1 to 0 and 1 to -1
875 for (size_t j = 0 ; j < interp_pts[T::value].get(i).size() ; j++)
876 {
877 for (size_t k = 0 ; k < dim ; k++)
878 interp_pts[T::value].get(i)[j].getComb()[k] = - ((interp_pts[T::value].get(i)[j].getComb()[k] == -1)?0:interp_pts[T::value].get(i)[j].getComb()[k]);
879 }
880 }
881 }
882};
883
884#endif /* SRC_GRID_STAGGERED_DIST_GRID_UTIL_HPP_ */
This class represent an N-dimensional box.
Definition Box.hpp:61
This class calculate elements of the hyper-cube.
Definition HyperCube.hpp:58
This class implement the point shape in an N-dimensional space.
Definition Point.hpp:28
This represent a sub-hyper-cube of an hyper-cube like a face or an edge of a cube.
grid_key_dx is the key to access any element in the grid
Definition grid_key.hpp:19
__device__ __host__ index_type get(index_type i) const
Get the i index.
Definition grid_key.hpp:503
Implementation of 1-D std::vector like structure.
It create separated grid for each properties to write them into a file.
stag_create_and_add_grid(st_grid &st_g, size_t p_id)
Constructor.
void operator()(T &t)
It call the copy function for each property.
stag_set_position(openfpm::vector< comb< dim > >(&pos_prp)[boost::fusion::result_of::size< v >::type::value])
vector containing the position of the properties in the cells (staggered properties are staggered)
this class is a functor for "for_each" algorithm
stag_set_position(openfpm::vector< comb< dim > >(&pos_prp)[boost::fusion::result_of::size< v >::type::value])
vector containing the position of the properties in the cells (staggered properties are staggered)
aggregate of properties, from a list of object if create a struct that follow the OPENFPM native stru...
Position of the element of dimension d in the hyper-cube of dimension dim.
Definition comb.hpp:35
signed char * getComb()
get the combination array pointer
Definition comb.hpp:260
signed char value(int i) const
get the index i of the combination
Definition comb.hpp:285
void zero()
Set all the elements to zero.
Definition comb.hpp:83
static size_t mul()
number of elements
static size_t mul()
number of elements
static size_t mul()
number of elements
static size_t mul()
number of elements
static size_t mul()
Matrix case return N1*N2 component.
static size_t mul()
Vector case return N1 component.
Classes to get the number of components of the properties.
static size_t dim()
Dimensionality.
static size_t mul()
Scalar case.
this class is a functor for "for_each" algorithm
interp_points(openfpm::vector< std::vector< comb< dim > > >(&interp_pts)[v_size::value], const openfpm::vector< comb< dim > >(&stag_pos)[n_prop])
constructor
void operator()(T &t)
It call the copy function for each property.
It copy the properties from one object to another.
This is a container to create a general object.
vtk_write(vtk vtk_w, const std::string output, const size_t i)
Add the grid with attributes name.
write a property that has attributes
vtk_write(vtk vtk_w, const std::string output, const size_t i)
Add the grid with attributes name.
static void write(sg &st_g, v_g &vg, size_t lg)
write the staggered grid
static void write(sg &st_g, v_g &vg, size_t lg)
write the staggered grid
Classes to copy each component into a grid and add to the VTKWriter the grid.
static void write(sg &st_g, v_g &vg, size_t lg)
write the staggered grid