Saras
Finite   Difference   Solver   for   Fluid   Dynamics   Simulations
poisson.h
Go to the documentation of this file.
1/********************************************************************************************************************************************
2 * Saras
3 *
4 * Copyright (C) 2019, Mahendra K. Verma
5 *
6 * All rights reserved.
7 *
8 * Redistribution and use in source and binary forms, with or without
9 * modification, are permitted provided that the following conditions are met:
10 * 1. Redistributions of source code must retain the above copyright
11 * notice, this list of conditions and the following disclaimer.
12 * 2. Redistributions in binary form must reproduce the above copyright
13 * notice, this list of conditions and the following disclaimer in the
14 * documentation and/or other materials provided with the distribution.
15 * 3. Neither the name of the copyright holder nor the
16 * names of its contributors may be used to endorse or promote products
17 * derived from this software without specific prior written permission.
18 *
19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
20 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
21 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
22 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR
23 * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
24 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
25 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
26 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
27 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
28 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 *
30 ********************************************************************************************************************************************
31 */
43#ifndef POISSON_H
44#define POISSON_H
45
46#include <blitz/array.h>
47#include <sys/time.h>
48#include <algorithm>
49#include <math.h>
50
51#include "plainsf.h"
52#include "grid.h"
53
54class poisson {
55 protected:
56 int vLevel, maxCount;
57
58 bool zeroBC;
59
60#ifdef TIME_RUN
61 real smothTimeComp;
62 real smothTimeTran;
63#endif
64
65 const grid &mesh;
66 const parser &inputParams;
67
68 blitz::Range all;
69
70 blitz::Array<blitz::Array<real, 3>, 1> pressureData;
71 blitz::Array<blitz::Array<real, 3>, 1> residualData;
72 blitz::Array<blitz::Array<real, 3>, 1> tmpDataArray;
73 blitz::Array<blitz::Array<real, 3>, 1> smoothedPres;
74
75 blitz::Array<blitz::RectDomain<3>, 1> stagFull;
76 blitz::Array<blitz::RectDomain<3>, 1> stagCore;
77 blitz::Array<int, 1> xEnd, yEnd, zEnd;
78
79 blitz::Array<int, 1> mgSizeArray;
80 blitz::Array<int, 1> strideValues;
81
82 blitz::TinyVector<int, 3> localSizeIndex;
83
84 blitz::Array<MPI_Request, 1> recvRequest;
85 blitz::Array<MPI_Status, 1> recvStatus;
86
87 blitz::Array<real, 1> hx, hy, hz;
88 blitz::Array<real, 1> hx2, hz2, hzhx;
89 blitz::Array<real, 1> hxhy, hyhz, hxhyhz;
90
91 blitz::Array<blitz::Array<real, 1>, 1> xixx, xix2;
92 blitz::Array<blitz::Array<real, 1>, 1> etyy, ety2;
93 blitz::Array<blitz::Array<real, 1>, 1> ztzz, ztz2;
94
95 blitz::Array<MPI_Datatype, 1> xMGArray;
96 blitz::Array<MPI_Datatype, 1> yMGArray;
97
98 blitz::Array<blitz::TinyVector<int, 3>, 1> mgSendLft, mgSendRgt;
99 blitz::Array<blitz::TinyVector<int, 3>, 1> mgRecvLft, mgRecvRgt;
100
101 blitz::Array<blitz::TinyVector<int, 3>, 1> mgSendFrn, mgSendBak;
102 blitz::Array<blitz::TinyVector<int, 3>, 1> mgRecvFrn, mgRecvBak;
103
104 static inline bool isOdd(int x) { return x % 2; };
105
106 void setLocalSizeIndex();
107 void initializeArrays();
108 void copyStaggrDerivs();
109 void setCoefficients();
110 void setStagBounds();
111
112 virtual void coarsen();
113 virtual void prolong();
114 virtual void computeResidual();
115 virtual void smooth(const int smoothCount);
116 virtual real computeError(const int normOrder);
117
118 virtual void solve() {};
119 virtual void imposeBC();
120 virtual void createMGSubArrays();
121 virtual void updatePads(blitz::Array<blitz::Array<real, 3>, 1> &data);
122
123 void vCycle();
124
125 public:
126 poisson(const grid &mesh, const parser &solParam);
127
128 void mgSolve(plainsf &inFn, const plainsf &rhs);
129
130 virtual real testProlong();
131 virtual real testTransfer();
132 virtual real testPeriodic();
133
134 virtual ~poisson();
135};
136
156class multigrid_d2: public poisson {
157 private:
158 void coarsen();
159 void prolong();
160 void computeResidual();
161 void smooth(const int smoothCount);
162 real computeError(const int normOrder);
163
164 void solve();
165
166 void imposeBC();
167 void initDirichlet();
168
169 void createMGSubArrays();
170 void updatePads(blitz::Array<blitz::Array<real, 3>, 1> &data);
171
172 blitz::Array<real, 1> xWall, zWall;
173
174 public:
175 multigrid_d2(const grid &mesh, const parser &solParam);
176
177 real testProlong();
178 real testTransfer();
179 real testPeriodic();
180
181 ~multigrid_d2() {};
182};
183
193class multigrid_d3: public poisson {
194 private:
195 void coarsen();
196 void prolong();
197 void computeResidual();
198 void smooth(const int smoothCount);
199 real computeError(const int normOrder);
200
201 void solve();
202
203 void imposeBC();
204 void initDirichlet();
205
206 void createMGSubArrays();
207 void updatePads(blitz::Array<blitz::Array<real, 3>, 1> &data);
208
209 blitz::Array<real, 2> xWall, yWall, zWall;
210
211 public:
212 multigrid_d3(const grid &mesh, const parser &solParam);
213
214 real testProlong();
215 real testTransfer();
216 real testPeriodic();
217
218 ~multigrid_d3() {};
219};
220
231#endif
Contains all the global variables related to the grid, its slices, limits, and grid derivatives used ...
Definition: grid.h:53
The derived class from poisson to perform multi-grid operations on a 2D grid.
Definition: poisson.h:156
real testTransfer()
Function to test whether strided data transfer is performing as expected.
Definition: poisson_d2.cc:510
real testPeriodic()
Function to test whether periodic BC is being implemented properly.
Definition: poisson_d2.cc:553
multigrid_d2(const grid &mesh, const parser &solParam)
Constructor of the multigrid_d2 class derived from the poisson class.
Definition: poisson_d2.cc:64
real testProlong()
Function to test whether prolongation operations are interpolating correctly.
Definition: poisson_d2.cc:481
The derived class from poisson to perform multi-grid operations on a 3D grid.
Definition: poisson.h:193
multigrid_d3(const grid &mesh, const parser &solParam)
Constructor of the multigrid_d3 class derived from the poisson class.
Definition: poisson_d3.cc:64
real testProlong()
Function to test whether prolongation operations are interpolating correctly.
Definition: poisson_d3.cc:689
real testTransfer()
Function to test whether strided data transfer is performing as expected.
Definition: poisson_d3.cc:722
real testPeriodic()
Function to test whether periodic BC is being implemented properly.
Definition: poisson_d3.cc:780
Contains all the global variables set by the user through the yaml file.
Definition: parser.h:63
Plain scalar field class to store simple scalar fields with no differentiation or interpolation.
Definition: plainsf.h:51
The base class poisson and its derived classes multigrid_d2 and multigrid_d3.
Definition: poisson.h:54
virtual real testPeriodic()
Function to test whether periodic BC is being implemented properly.
Definition: poisson.cc:643
void setLocalSizeIndex()
Function to calculate the local size indices of sub-domains after MPI domain decomposition.
Definition: poisson.cc:377
virtual void smooth(const int smoothCount)
Function to perform smoothing operation on the input array.
Definition: poisson.cc:550
virtual void imposeBC()
Function to impose the boundary conditions of Poisson solver at different levels of the V-cycle.
Definition: poisson.cc:578
virtual real testProlong()
Function to test whether prolongation operations are interpolating correctly.
Definition: poisson.cc:631
virtual void coarsen()
Function to coarsen the grid down the levels of the V-Cycle.
Definition: poisson.cc:509
void copyStaggrDerivs()
Function to copy the staggered grid derivatives from the grid class to local arrays.
Definition: poisson.cc:453
void setCoefficients()
Function to set the coefficients used for calculating laplacian and in smoothing.
Definition: poisson.cc:399
void initializeArrays()
Function to initialize the arrays used in multi-grid.
Definition: poisson.cc:293
virtual void computeResidual()
Function to compute the residual at the start of each V-Cycle.
Definition: poisson.cc:535
void mgSolve(plainsf &inFn, const plainsf &rhs)
The core, publicly accessible function of poisson to compute the solution for the Poisson equation.
Definition: poisson.cc:100
void setStagBounds()
Function to set the RectDomain variables for all future references throughout the poisson solver.
Definition: poisson.cc:326
virtual real testTransfer()
Function to test whether strided data transfer is performing as expected.
Definition: poisson.cc:619
virtual void prolong()
Function to perform prolongation on the array being solved.
Definition: poisson.cc:522
void vCycle()
Function to perform one loop of V-cycle.
Definition: poisson.cc:219
virtual void createMGSubArrays()
Function to create the MPI sub-array data types necessary to transfer data across sub-domains.
Definition: poisson.cc:606
poisson(const grid &mesh, const parser &solParam)
Constructor of the base poisson class.
Definition: poisson.cc:62
virtual void updatePads(blitz::Array< blitz::Array< real, 3 >, 1 > &data)
Function to update the pad points of the local sub-domains at different levels of the V-cycle.
Definition: poisson.cc:591
virtual real computeError(const int normOrder)
Function to compute the error at the end of each V-Cycle.
Definition: poisson.cc:563
Class declaration of grid.
Class declaration of plainsf - plain scalar field.