Line data Source code
1 : ///////////////////////////////////////////////////////////////////////////////
2 : // Copyright 2019 Lakshman Anumolu, Raunak Bardia.
3 : //
4 : // Redistribution and use in source and binary forms, with or without
5 : // modification, are permitted provided that the following conditions are
6 : // met:
7 : //
8 : // 1. Redistributions of source code must retain the above copyright notice,
9 : // this list of conditions and the following disclaimer.
10 : //
11 : // 2. Redistributions in binary form must reproduce the above copyright notice,
12 : // this list of conditions and the following disclaimer in the documentation
13 : // and/or other materials provided with the distribution.
14 : //
15 : // 3. Neither the name of the copyright holder nor the names of its contributors
16 : // may be used to endorse or promote products derived from this software without
17 : // specific prior written permission.
18 : //
19 : // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
20 : // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
21 : // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
22 : // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
23 : // HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
24 : // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
25 : // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
26 : // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
27 : // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
28 : // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
29 : // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30 : ///////////////////////////////////////////////////////////////////////////////
31 :
32 : #include "gals/utilities/grid.h"
33 :
34 : #include <limits.h>
35 : #include <math.h>
36 : #include <fstream>
37 :
38 : #include "gals/input-parser.h"
39 :
40 : template <typename T, int DIM>
41 63 : GALS::CPU::Grid<T, DIM>::Grid(int nx, int ny, int nz)
42 63 : : m_dimension(DIM), m_nx(nx), m_ny(ny), m_nz(nz), m_pad(1), m_total_cells(nx * ny * nz)
43 : {
44 63 : m_mask[0] = 1, m_mask[1] = 0, m_mask[2] = 0;
45 :
46 32 : if (dim == 2)
47 : m_mask[1] = 1;
48 : else if (dim == 3)
49 : m_mask[1] = 1, m_mask[2] = 1;
50 :
51 63 : m_box_min = Vec3<T>(INT_MIN, INT_MIN, INT_MIN);
52 63 : m_box_max = Vec3<T>(INT_MAX, INT_MAX, INT_MAX);
53 :
54 63 : m_dx = GALS::CPU::Vec3<T>(INT_MAX, INT_MAX, INT_MAX), m_one_over_dx = GALS::CPU::Vec3<T>(INT_MIN, INT_MIN, INT_MIN);
55 63 : }
56 :
57 : template <typename T, int DIM>
58 3 : GALS::CPU::Grid<T, DIM>::Grid(int nx, int ny) : Grid(nx, ny, 1)
59 : {
60 3 : }
61 :
62 : template <typename T, int DIM>
63 2 : GALS::CPU::Grid<T, DIM>::Grid(int nx) : Grid(nx, 1, 1)
64 : {
65 2 : }
66 :
67 : template <typename T, int DIM>
68 63 : GALS::CPU::Grid<T, DIM>::~Grid()
69 : {
70 63 : m_grid.clear();
71 : m_grid.shrink_to_fit();
72 63 : }
73 :
74 : template <typename T, int DIM>
75 147457 : GALS::CPU::Vec3<T>& GALS::CPU::Grid<T, DIM>::x(const int i, const int j, const int k)
76 : {
77 294914 : return m_grid[this->index(i, j, k)];
78 : }
79 :
80 : template <typename T, int DIM>
81 3 : const int GALS::CPU::Grid<T, DIM>::dimension() const
82 : {
83 3 : return m_dimension;
84 : }
85 :
86 : template <typename T, int DIM>
87 280 : const int GALS::CPU::Grid<T, DIM>::size() const
88 : {
89 560 : return m_grid.size();
90 : }
91 :
92 : template <typename T, int DIM>
93 40 : const int* GALS::CPU::Grid<T, DIM>::getMask() const
94 : {
95 40 : return m_mask;
96 : }
97 :
98 : template <typename T, int DIM>
99 16590 : const GALS::CPU::Vec3<int> GALS::CPU::Grid<T, DIM>::numCells() const
100 : {
101 16590 : return Vec3<int>(m_nx, m_ny, m_nz);
102 : }
103 :
104 : template <typename T, int DIM>
105 43 : const size_t GALS::CPU::Grid<T, DIM>::totalCells() const
106 : {
107 43 : return m_total_cells;
108 : }
109 :
110 : template <typename T, int DIM>
111 310 : const int GALS::CPU::Grid<T, DIM>::getPadding() const
112 : {
113 310 : return m_pad;
114 : }
115 :
116 : template <typename T, int DIM>
117 6299686 : const std::size_t GALS::CPU::Grid<T, DIM>::index(const int i, const int j, const int k) const
118 : {
119 18899058 : const std::size_t idx = ((k + m_pad) * (m_nx + 2 * m_pad) * (m_ny + 2 * m_pad) * m_mask[2]) +
120 12599372 : ((j + m_pad) * (m_nx + 2 * m_pad) * m_mask[1]) + ((i + m_pad) * m_mask[0]);
121 :
122 6299686 : return idx;
123 : }
124 :
125 : template <typename T, int DIM>
126 682536 : const std::size_t GALS::CPU::Grid<T, DIM>::index(const Vec3<int> node_id) const
127 : {
128 682536 : return this->index(node_id[0], node_id[1], node_id[2]);
129 : }
130 :
131 : template <typename T, int DIM>
132 5870 : const GALS::CPU::Vec3<int> GALS::CPU::Grid<T, DIM>::baseNodeId(const Vec3<T>& x) const
133 : {
134 5870 : Vec3<int> base_node_id(INT_MAX, INT_MAX, INT_MAX);
135 :
136 38998 : for (int axis = 0; axis < DIM; ++axis) {
137 16564 : base_node_id[axis] = floor(((x[axis] - m_box_min[axis]) * m_one_over_dx[axis]) - 0.5);
138 : }
139 :
140 5870 : return base_node_id;
141 : }
142 :
143 : template <typename T, int DIM>
144 5789 : const GALS::CPU::Vec3<T>& GALS::CPU::Grid<T, DIM>::dX() const
145 : {
146 5789 : return m_dx;
147 : }
148 :
149 : template <typename T, int DIM>
150 5872 : const GALS::CPU::Vec3<T>& GALS::CPU::Grid<T, DIM>::oneOverDX() const
151 : {
152 5872 : return m_one_over_dx;
153 : }
154 :
155 : template <typename T, int DIM>
156 2556764 : const GALS::CPU::Vec3<T>& GALS::CPU::Grid<T, DIM>::operator()(const int i, const int j, const int k) const
157 : {
158 5113528 : return m_grid[this->index(i, j, k)];
159 : }
160 :
161 : template <typename T, int DIM>
162 6087 : const GALS::CPU::Vec3<T>& GALS::CPU::Grid<T, DIM>::operator()(const Vec3<int> node_id) const
163 : {
164 12174 : return m_grid[this->index(node_id)];
165 : }
166 :
167 : template <typename T, int DIM>
168 9 : void GALS::CPU::Grid<T, DIM>::setPadding(const int pad)
169 : {
170 9 : m_pad = pad;
171 9 : }
172 :
173 : template <typename T, int DIM>
174 57 : void GALS::CPU::Grid<T, DIM>::generate(T x_min, T x_max, T y_min, T y_max, T z_min, T z_max)
175 : {
176 114 : if (m_grid.size()) m_grid.clear(), m_grid.shrink_to_fit();
177 :
178 57 : m_box_min[0] = x_min, m_box_min[1] = y_min, m_box_min[2] = z_min;
179 57 : m_box_max[0] = x_max, m_box_max[1] = y_max, m_box_max[2] = z_max;
180 :
181 57 : m_grid.resize((m_nz + 2 * m_pad * m_mask[2]) * (m_ny + 2 * m_pad * m_mask[1]) * (m_nx + 2 * m_pad * m_mask[0]));
182 :
183 57 : std::vector<T> domain_min({x_min, y_min, z_min}), domain_min_new(3);
184 :
185 57 : m_dx[0] = (x_max - x_min) / m_nx;
186 57 : m_dx[1] = (y_max - y_min) / m_ny;
187 57 : m_dx[2] = (z_max - z_min) / m_nz;
188 :
189 152 : for (int i = 0; i < DIM; ++i) m_one_over_dx[i] = static_cast<T>(1.) / m_dx[i];
190 :
191 342 : for (int i = 0; i < DIM; ++i) domain_min_new[i] = domain_min[i] + (m_dx[i] * 0.5) - (m_dx[i] * m_mask[i]);
192 :
193 114 : Vec3<T> elem;
194 :
195 57 : int i_min = -m_pad * m_mask[0], j_min = -m_pad * m_mask[1], k_min = -m_pad * m_mask[2];
196 :
197 114 : GALS::CPU::Vec3<int> vec;
198 :
199 2955 : for (int i = i_min; i < m_nx + m_pad * m_mask[0]; ++i) {
200 86507 : for (int j = j_min; j < m_ny + m_pad * m_mask[1]; ++j) {
201 337003 : for (int k = k_min; k < m_nz + m_pad * m_mask[2]; ++k) {
202 147237 : vec[0] = i - i_min;
203 147237 : vec[1] = j - j_min;
204 147237 : vec[2] = k - k_min;
205 :
206 951379 : for (int axis = 0; axis < DIM; ++axis) elem[axis] = domain_min_new[axis] + vec[axis] * m_dx[axis];
207 :
208 147237 : this->x(i, j, k) = elem;
209 : }
210 : }
211 : }
212 :
213 : // Update total cells.
214 57 : m_total_cells = 1;
215 152 : for (int d = 0; d < DIM; ++d) m_total_cells *= this->numCells()[d];
216 57 : }
217 :
218 : template class GALS::CPU::Grid<double, 1>;
219 : template class GALS::CPU::Grid<double, 2>;
220 3 : template class GALS::CPU::Grid<double, 3>;
|