ASL  0.1.6
Advanced Simulation Library
multicomponent_flow.cc
Go to the documentation of this file.
1 /*
2  * Advanced Simulation Library <http://asl.org.il>
3  *
4  * Copyright 2015 Avtech Scientific <http://avtechscientific.com>
5  *
6  *
7  * This file is part of Advanced Simulation Library (ASL).
8  *
9  * ASL is free software: you can redistribute it and/or modify it
10  * under the terms of the GNU Affero General Public License as
11  * published by the Free Software Foundation, version 3 of the License.
12  *
13  * ASL is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  * GNU Affero General Public License for more details.
17  *
18  * You should have received a copy of the GNU Affero General Public License
19  * along with ASL. If not, see <http://www.gnu.org/licenses/>.
20  *
21  */
22 
23 
28 #include <utilities/aslParametersManager.h>
29 #include <math/aslTemplates.h>
30 #include <aslGeomInc.h>
31 #include <aslDataInc.h>
32 #include <acl/aclGenerators.h>
33 #include <writers/aslVTKFormatWriters.h>
34 #include <num/aslLBGK.h>
35 #include <num/aslLBGKBC.h>
36 #include <utilities/aslTimer.h>
37 #include <num/aslFDAdvectionDiffusion.h>
38 #include <num/aslBasicBC.h>
39 
40 // typedef to switch to double precision
41 //typedef double FlT;
42 typedef float FlT;
43 
44 using asl::AVec;
45 using asl::makeAVec;
46 
47 class Parameters
48 {
49  private:
50  void init();
51 
52  public:
54 
58 
61 
64 
69 
73 
74 
75  void load(int argc, char * argv[]);
76  Parameters();
77  void updateNumValues();
78 };
79 
80 
82  appParamsManager("multicomponent_flow", "0.1"),
83  size(3),
84  dx(0.0005, "dx", "space step"),
85  dt(1., "dt", "time step"),
86  tSimulation(2e-3, "simulation_time", "simulation time"),
87  tOutput(1e-4, "output_interval", "output interval"),
88  nu(4e-8/1.6, "nu", "viscosity"),
89  tubeL(0.25, "tubeL", "tube's length"),
90  tubeD(0.05, "tubeD", "tube's diameter"),
91  pumpL(0.025, "pumpL", "pump's length"),
92  pumpD(0.03, "pumpD", "pump's diameter"),
93  component1InVel(0.16, "component1_in_velocity", "flow velocity in the component1 input"),
94  component2InVel(0.08, "component2_in_velocity", "flow velocity in the component2 input"),
95  component3InVel(0.1, "component3_in_velocity", "flow velocity in the component3 input")
96 {
97 }
98 
99 
100 void Parameters::load(int argc, char * argv[])
101 {
102  appParamsManager.load(argc, argv);
103  init();
104 }
105 
106 
108 {
109  nuNum = nu.v() * dt.v() / dx.v() / dx.v();
110  size[0] = tubeD.v() / dx.v() + 1;
111  size[1] = (tubeD.v() + 2 * pumpL.v()) / dx.v() + 1;
112  size[2] = tubeL.v() / dx.v() + 1;
113 }
114 
115 
116 void Parameters::init()
117 {
118  if (tubeD.v() < pumpD.v())
119  asl::errorMessage("Tube's diameter is smaller than pump's diameter");
120 
121  updateNumValues();
122 }
123 
124 // Generate geometry of the mixer (cross-coupled pipes)
126 {
127  asl::SPDistanceFunction mixerGeometry;
128  asl::AVec<double> orientation(asl::makeAVec(0., 0., 1.));
129  asl::AVec<double> center(asl::AVec<double>(params.size) * .5 * params.dx.v());
130 
131  mixerGeometry = generateDFCylinderInf(params.tubeD.v() / 2., orientation,
132  center);
133 
134  orientation[1] = 1.0;
135  orientation[2] = 0.0;
136  center[2] = params.pumpD.v() * 1.5;
137  mixerGeometry = mixerGeometry | generateDFCylinderInf(params.pumpD.v() / 2.,
138  orientation, center);
139 
140  return asl::normalize(-(mixerGeometry) | asl::generateDFInBlock(block, 0),
141  params.dx.v());
142 }
143 
144 int main(int argc, char *argv[])
145 {
146  Parameters params;
147  params.load(argc, argv);
148 
149  cout << "Data initialization..." << endl;
150 
151  asl::Block block(params.size, params.dx.v());
152 
153  auto mcfMapMem(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
154  asl::initData(mcfMapMem, generateMixer(block, params));
155 
156  auto component1Frac(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
157  asl::initData(component1Frac, 0);
158  auto component3Frac(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
159  asl::initData(component3Frac, 0);
160 
161 
162  cout << "Finished" << endl;
163 
164  cout << "Numerics initialization..." << endl;
165 
166  auto templ(&asl::d3q15());
167 
168  asl::SPLBGK lbgk(new asl::LBGK(block,
169  acl::generateVEConstant(FlT(params.nuNum.v())),
170  templ));
171 
172  lbgk->init();
173  asl::SPLBGKUtilities lbgkUtil(new asl::LBGKUtilities(lbgk));
174  lbgkUtil->initF(acl::generateVEConstant(.0, .0, .0));
175 
176  auto flowVel(lbgk->getVelocity());
177  auto nmcomponent1(asl::generateFDAdvectionDiffusion(component1Frac, 0.01,
178  flowVel, templ, true));
179  nmcomponent1->init();
180  auto nmcomponent3(asl::generateFDAdvectionDiffusion(component3Frac, 0.01,
181  flowVel, templ));
182  nmcomponent3->init();
183 
184  std::vector<asl::SPNumMethod> bc;
185  std::vector<asl::SPNumMethod> bcV;
186  std::vector<asl::SPNumMethod> bcDif;
187 
188  bc.push_back(generateBCNoSlip(lbgk, mcfMapMem));
189  bc.push_back(generateBCConstantPressure(lbgk, 1., {asl::ZE}));
190  bc.push_back(generateBCConstantPressureVelocity(lbgk, 1.,
191  makeAVec(0., 0., params.component2InVel.v()),
192  {asl::Z0}));
193  bc.push_back(generateBCConstantPressureVelocity(lbgk, 1.,
194  makeAVec(0., -params.component1InVel.v(), 0.),
195  {asl::YE}));
196  bc.push_back(generateBCConstantPressureVelocity(lbgk, 1.,
197  makeAVec(0., params.component3InVel.v(), 0.),
198  {asl::Y0}));
199 
200  bcDif.push_back(generateBCNoSlipVel(lbgk, mcfMapMem));
201  bc.push_back(generateBCConstantGradient(component1Frac, 0., mcfMapMem, templ));
202  bc.push_back(generateBCConstantGradient(component3Frac, 0., mcfMapMem, templ));
203  bc.push_back(generateBCConstantValue(component1Frac, 1., {asl::YE}));
204  bc.push_back(generateBCConstantValue(component3Frac, 0., {asl::YE, asl::Z0, asl::ZE}));
205  bc.push_back(generateBCConstantValue(component1Frac, 0., {asl::Y0, asl::Z0, asl::ZE}));
206  bc.push_back(generateBCConstantValue(component3Frac, 1., {asl::Y0}));
207 // bc.push_back(generateBCConstantGradient(component1Frac, 0.,templ, {asl::ZE}));
208 // bc.push_back(generateBCConstantGradient(component3Frac, 0.,templ, {asl::ZE}));
209 
210  initAll(bc);
211  initAll(bcDif);
212  initAll(bcV);
213 
214  cout << "Finished" << endl;
215  cout << "Computing..." << endl;
216  asl::Timer timer;
217 
218  asl::WriterVTKXML writer("multicomponent_flow");
219  writer.addScalars("map", *mcfMapMem);
220  writer.addScalars("component1", *component1Frac);
221  writer.addScalars("component3", *component3Frac);
222  writer.addScalars("rho", *lbgk->getRho());
223  writer.addVector("v", *flowVel);
224 
225  executeAll(bc);
226  executeAll(bcDif);
227  executeAll(bcV);
228 
229  writer.write();
230 
231  timer.start();
232  for (unsigned int i(1); i < 10001; ++i)
233  {
234  lbgk->execute();
235  executeAll(bcDif);
236  nmcomponent1->execute();
237  nmcomponent3->execute();
238  executeAll(bc);
239 
240  if (!(i%100))
241  {
242  timer.stop();
243  cout << i << "/10000; time left (expected): " << timer.getLeftTime(double(i)/10000.) << endl;
244  executeAll(bcV);
245  writer.write();
246  timer.resume();
247  }
248  }
249  timer.stop();
250 
251  cout << "Finished" << endl;
252 
253  cout << "Computation statistic:" << endl;
254  cout << "time = " << timer.getTime() << "; clockTime = "
255  << timer.getClockTime() << "; load = "
256  << timer.getProcessorLoad() * 100 << "%" << endl;
257 
258  return 0;
259 }
const double getTime() const
Definition: aslTimer.h:46
asl::Block::DV size
const double getProcessorLoad() const
Definition: aslTimer.h:48
float FlT
asl::ApplicationParametersManager appParamsManager
asl::Parameter< double > nu
Numerical method for fluid flow.
Definition: aslLBGK.h:77
void resume()
Definition: aslTimer.h:44
const double getClockTime() const
Definition: aslTimer.h:47
SPDistanceFunction normalize(SPDistanceFunction a, double dx)
void errorMessage(cl_int status, const char *errorMessage)
Prints errorMessage and exits depending on the status.
void addVector(std::string name, AbstractData &data)
const T & v() const
Definition: aslUValue.h:43
asl::Parameter< double > tubeL
asl::Parameter< double > pumpD
void initAll(std::vector< T * > &v)
Definition: aslNumMethod.h:67
SPBCond generateBCConstantPressureVelocity(SPLBGK nm, double p, AVec<> v, const std::vector< SlicesNames > &sl)
SPDistanceFunction generateDFCylinderInf(double r, const AVec< double > &l, const AVec< double > &c)
generates infinite cylinder
asl::UValue< double > dt
const VectorTemplate & d3q15()
Vector template.
const T & v() const
SPFDAdvectionDiffusion generateFDAdvectionDiffusion(SPDataWithGhostNodesACLData c, double diffustionCoeff, SPAbstractDataWithGhostNodes v, const VectorTemplate *vt, bool compressibilityCorrection=false)
asl::Parameter< double > tSimulation
AVec< T > makeAVec(T a1)
asl::Parameter< double > tOutput
void load(int argc, char *argv[])
asl::Parameter< double > component1InVel
asl::Parameter< double > pumpL
void initData(SPAbstractData d, double a)
asl::Parameter< double > component2InVel
asl::UValue< double > nuNum
const double getLeftTime(double fractTime)
Definition: aslTimer.h:51
SPBCond generateBCConstantGradient(SPAbstractDataWithGhostNodes d, double v, const VectorTemplate *const t, const std::vector< SlicesNames > &sl)
Bondary condition that makes fixed gradient.
void executeAll(std::vector< T * > &v)
Definition: aslNumMethod.h:55
void addScalars(std::string name, AbstractData &data)
SPBCond generateBCConstantPressure(SPLBGK nm, double p, const std::vector< SlicesNames > &sl)
asl::Parameter< double > component3InVel
void stop()
Definition: aslTimer.h:45
asl::Parameter< double > dx
void start()
Definition: aslTimer.h:43
VectorOfElements generateVEConstant(T a)
Generates VectorOfElements with 1 Element acl::Constant with value a.
asl::Parameter< double > tubeD
SPBCond generateBCConstantValue(SPAbstractDataWithGhostNodes d, double v, const std::vector< SlicesNames > &sl)
Bondary condition that puts fixed value in each point.
void updateNumValues()
asl::SPDistanceFunction generateMixer(asl::Block &block, Parameters &params)
asl::Parameter< double > flowVel
void load(int argc, char *argv[])
SPDistanceFunction generateDFInBlock(const Block &b, unsigned int nG)
generates map corresponding to external (ghost) part of the block
SPNumMethod generateBCNoSlipVel(SPLBGK nmU, SPAbstractDataWithGhostNodes map)
for velocity field
contains different kernels for preprocessing and posprocessing of data used by LBGK ...
Definition: aslLBGK.h:137
int main(int argc, char *argv[])
SPBCond generateBCNoSlip(SPLBGK nm, const std::vector< SlicesNames > &sl)