1/*
2 * staggered_grid.hpp
3 *
4 * Created on: Aug 19, 2015
5 * Author: i-bird
6 */
7
8#ifndef SRC_GRID_STAGGERED_DIST_GRID_HPP_
9#define SRC_GRID_STAGGERED_DIST_GRID_HPP_
10
11#include "Grid/grid_dist_id.hpp"
12#include "staggered_dist_grid_util.hpp"
13#include "VTKWriter/VTKWriter.hpp"
14#include "staggered_dist_grid_copy.hpp"
15#include <boost/mpl/vector_c.hpp>
16
17/*! \brief Implementation of the staggered grid
18 *
19 * \param dim Dimensionality of the staggered grid
20 * \param ele elements object on each dimensional objects, must be a stag_elements
21 *
22 * \verbatim
23
24 +--#--+--#--+--#--+--#--+--#--+--#--+
25 | | | | | | |
26 # * # * # * # * # * # * #
27 | | | | | | |
28 +--#--+--#--+--#--+--#--+--#--+--#--+
29 | | | | | | |
30 # * # * # * # * # * # * #
31 | | | | | | |
32 +--#--+--#--+--#--+--#--+--#--+--#--+
33 | | | | | | |
34 # * # * # * # * # * # * #
35 | | | | | | |
36 +--#--+--#--+--#--+--#--+--#--+--#--+
37 | | | | | | |
38 # * # * # * # * # * # * #
39 | | | | | | |
40 +--#--+--#--+--#--+--#--+--#--+--#--+
41 | | | | | | |
42 # * # * # * # * # * # * #
43 | | | | | | |
44 +--#--+--#--+--#--+--#--+--#--+--#--+
45
46\endverbatim
47
48 In the case of a 2D staggered grid we have 3 (in general dim+1 ) elements
49
50 + = vertex
51 # = edge
52 * = volume
53
54 ele = stag_ele<scalar<float>,Point_test<float>,scalar<float>>
55
56 It place a scalar on (*) an object Point_test<float> on (#) and an object scalar<float> on (+)
57
58 *
59 *
60 *
61 */
62template<unsigned int dim, typename St, typename T, typename Decomposition=CartDecomposition<dim,St>,typename Memory=HeapMemory , typename device_grid=grid_cpu<dim,T>>
63class staggered_grid_dist : public grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>
64{
65 //! position of the properties in the grid cell
66 openfpm::vector<comb<dim>> c_prp[T::max_prop_real];
67
68public:
69
70 //! Properties for each grid point
71 typedef T value_type;
72
73 //! Number of dimensions
74 static const unsigned int dims = dim;
75
76 /*! \brief This constructor is special, it construct an expanded grid that perfectly overlap with the previous
77 *
78 * The key-word here is "perfectly overlap". Using the default constructor you could create
79 * something similar, but because of rounding-off error it can happen that it is not perfectly overlapping
80 *
81 * \param g previous grid
82 * \param gh Ghost part in grid units
83 * \param ext extension of the grid (must be positive on every direction)
84 *
85 */
86 template<typename H> staggered_grid_dist(const grid_dist_id<dim,St,H,typename Decomposition::base_type,Memory,grid_cpu<dim,H>> & g,
87 const Ghost<dim,long int> & gh,
88 Box<dim,size_t> ext)
89 :grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>(g,gh,ext)
90 {
91 }
92
93 /*! \brief Constructor
94 *
95 * \param g_sz size of the staggered grid
96 * \param domain domain
97 * \param ghost part
98 *
99 *
100 */
101 staggered_grid_dist(const size_t (& g_sz)[dim],
102 const Box<dim,St> & domain,
103 const Ghost<dim,St> & ghost)
104 :grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>(g_sz,domain,ghost)
105 {}
106
107 /*! \brief set the staggered positions of the properties
108 *
109 * \tparam property p
110 *
111 * \param cmb a vector containing for each component the position in the cell-grid
112 *
113 */
114 template<unsigned int p> void setStagPosition(openfpm::vector<comb<dim>> & cmb)
115 {
116#ifdef SE_CLASS1
117 if (extends< typename boost::mpl::at<typename T::type,boost::mpl::int_<p> >::type >::mul() != cmb.size())
118 std::cerr << __FILE__ << ":" << __LINE__ << " error properties has " << extends< typename boost::mpl::at<typename T::type,boost::mpl::int_<p> >::type >::mul() << " components, but " << cmb.size() << "has been defined \n";
119#endif
120 c_prp[p] = cmb;
121 }
122
123 /*! \brief It set all the properties defined to be staggered on the default location
124 *
125 */
126 void setDefaultStagPosition()
127 {
128 // for each properties
129
130 stag_set_position<dim,typename T::type_real> ssp(c_prp);
131
132 boost::mpl::for_each_ref< boost::mpl::range_c<int,0,T::max_prop_real> >(ssp);
133 }
134
135 /*! \brief Copy the staggered grid into a normal one
136 *
137 *
138 * \tparam Grid_dst type of the destination Grid
139 * \tparam pos destination grid properties to fill
140 *
141 * \param g_dst destination grid
142 * \param pd padding of the grid compared to the destination grid
143 * \param start starting point
144 * \param stop stop point
145 *
146 */
147 template<typename Grid_dst ,unsigned int ... pos> bool to_normal(Grid_dst & g_dst, const Padding<dim> & pd, const long int (& start)[dim], const long int (& stop)[dim])
148 {
149 // interpolation points for each property
150 openfpm::vector<std::vector<comb<dim>>> interp_pos[sizeof...(pos)];
151
152 typedef boost::mpl::vector_c<unsigned int,pos ... > v_pos_type;
153
154 interp_points<dim,T::max_prop_real,v_pos_type,typename Grid_dst::value_type::type> itp(interp_pos,c_prp);
155 boost::mpl::for_each_ref<v_pos_type>(itp);
156
157 // shift the start and stop by the padding
158 grid_key_dx<dim> start_k = grid_key_dx<dim>(start);
159 grid_key_dx<dim> stop_k = grid_key_dx<dim>(stop);
160
161 // sub-grid iterator over the grid map
162 auto g_map_it = this->getSubDomainIterator(start_k,stop_k);
163
164 // Iterator over the destination grid
165 auto g_dst_it = g_dst.getDomainIterator();
166
167 // Check that the 2 iterator has the same size
168 checkIterator<Grid_dst,decltype(g_map_it),decltype(g_dst_it)>(g_map_it,g_dst_it);
169
170 while (g_map_it.isNext() == true)
171 {
172 typedef typename to_boost_vmpl<pos...>::type vid;
173 typedef boost::mpl::size<vid> v_size;
174
175 auto key_src = g_map_it.get();
176
177 // destination point
178 auto key_dst = g_dst_it.get();
179
180 // Transform this id into an id for the Eigen vector
181
182 interp_ele<vid,Grid_dst,typename std::remove_reference<decltype(*this)>::type,sizeof...(pos)> cp(key_dst,g_dst,*this,key_src,interp_pos);
183
184 // For each property in the target grid
185 boost::mpl::for_each_ref<boost::mpl::range_c<int,0,v_size::value>>(cp);
186
187 ++g_map_it;
188 ++g_dst_it;
189 }
190
191 return true;
192 }
193
194 /*! \brief Get the staggered positions
195 *
196 * \return for each property it contain a vector that specify where in the cell
197 * grid the properties live
198 *
199 */
200 const openfpm::vector<comb<dim>> (& getStagPositions()) [T::max_prop]
201 {
202 return c_prp;
203 }
204
205 /*! \brief Write a vtk file with the information of the staggered grid
206 *
207 * \param str vtk output file
208 *
209 */
210 void write(std::string str)
211 {
212 stag_create_and_add_grid<dim,staggered_grid_dist<dim,St,T,Decomposition,Memory,device_grid>,St> sgw(*this, this->getVC().getProcessUnitID());
213
214 boost::mpl::for_each_ref< boost::mpl::range_c<int,0,T::max_prop> >(sgw);
215 }
216
217 /*! \brief Return if the properties is a staggered property or not
218 *
219 * \param prp property to check
220 *
221 * \return true if the property is staggered
222 *
223 */
224 bool is_staggered_prop(size_t prp)
225 {
226 return c_prp[prp].size() != 0;
227 }
228
229 /*! \brief Return if the grid is staggered
230 *
231 * \return true
232 *
233 */
234 bool is_staggered()
235 {
236 return true;
237 }
238
239 friend class stag_create_and_add_grid<dim,staggered_grid_dist<dim,St,T,Decomposition,Memory,device_grid>,St>;
240};
241
242#endif /* SRC_GRID_STAGGERED_DIST_GRID_HPP_ */
243