ASL  0.1.7
Advanced Simulation Library
multiphase_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 
30 #include <utilities/aslParametersManager.h>
31 #include <math/aslTemplates.h>
32 #include <aslGeomInc.h>
33 #include <aslDataInc.h>
34 #include <acl/aclGenerators.h>
35 #include <writers/aslVTKFormatWriters.h>
36 #include <num/aslLBGK.h>
37 #include <num/aslLBGKBC.h>
38 #include <utilities/aslTimer.h>
39 #include <num/aslFDMultiPhase.h>
40 #include <num/aslBasicBC.h>
41 
42 
43 typedef float FlT;
44 //typedef double FlT;
46 
47 using asl::AVec;
48 using asl::makeAVec;
49 
50 class Parameters
51 {
52  private:
53  void init();
54 
55  public:
57 
59 
62 
65 
68 
73 
77 
78 
79  void load(int argc, char * argv[]);
80  string getDir();
83 };
84 
85 
87  appParamsManager("multiphase_flow", "0.1"),
88  size(3),
89  dx(0.002, "dx", "space step"),
90  dt(1., "dt", "time step"),
91  tSimulation(2e-3, "simulation_time", "simulation time"),
92  tOutput(1e-4, "output_interval", "output interval"),
93  nu(4e-8, "nu", "viscosity"),
94  tubeL(0.5, "tubeL", "tube's length"),
95  tubeD(0.05, "tubeD", "tube's diameter"),
96  pumpL(0.025, "pumpL", "pump's length"),
97  pumpD(0.03, "pumpD", "pump's diameter"),
98  oilInVel(0.02, "oil_in_velocity", "flow velocity in the oil input"),
99  waterInVel(0.04, "water_in_velocity", "flow velocity in the water input"),
100  gasInVel(0.03, "gas_in_velocity", "flow velocity in the gas input")
101 {
102 }
103 
104 
105 void Parameters::load(int argc, char * argv[])
106 {
107  appParamsManager.load(argc, argv);
108 
109  init();
110 }
111 
112 
114 {
115  return appParamsManager.getDir();
116 }
117 
118 
120 {
121  nuNum = nu.v() * dt.v() / dx.v() / dx.v();
122  size[0] = tubeD.v() / dx.v() + 1;
123  size[1] = (tubeD.v() + 2 * pumpL.v()) / dx.v() + 1;
124  size[2] = tubeL.v() / dx.v() + 1;
125 }
126 
127 
128 void Parameters::init()
129 {
130  if (tubeD.v() < pumpD.v())
131  asl::errorMessage("Tube's diameter is smaller than pump's diameter");
132 
133  updateNumValues();
134 }
135 
136 // generate geometry of the mixer
138 {
139  asl::SPDistanceFunction mixerGeometry;
140  asl::AVec<double> orientation(asl::makeAVec(0., 0., 1.));
141  asl::AVec<double> center(asl::AVec<double>(params.size)*.5*params.dx.v());
142 
143  mixerGeometry = generateDFCylinderInf(params.tubeD.v() / 2., orientation, center);
144 
145  orientation[1] = 1.0;
146  orientation[2] = 0.0;
147  center[2]=params.pumpD.v() * 1.5;
148  mixerGeometry = mixerGeometry | generateDFCylinderInf(params.pumpD.v() / 2., orientation, center);
149 
150  return asl::normalize(-(mixerGeometry) | asl::generateDFInBlock(block, 0), params.dx.v());
151 }
152 
153 int main(int argc, char *argv[])
154 {
155  Parameters params;
156  params.load(argc, argv);
157 
158  std::cout << "Data initialization...";
159 
160  asl::Block block(params.size, params.dx.v());
161 
162  auto mpfMapMem(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
163  asl::initData(mpfMapMem, generateMixer(block, params));
164 
165  auto waterFrac(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
166  asl::initData(waterFrac, 0);
167 
168  std::cout << "Finished" << endl;
169 
170  std::cout << "Numerics initialization...";
171 
172  auto templ(&asl::d3q15());
173 
174  asl::SPLBGK lbgk(new asl::LBGK(block,
175  acl::generateVEConstant(FlT(params.nuNum.v())),
176  templ));
177 
178  lbgk->init();
179  asl::SPLBGKUtilities lbgkUtil(new asl::LBGKUtilities(lbgk));
180  lbgkUtil->initF(acl::generateVEConstant(.0, .0, .0));
181 
182  auto flowVel(lbgk->getVelocity());
183  auto nmWater(asl::generateFDMultiPhase(waterFrac, flowVel, templ, true));
184  nmWater->init();
185 
186  std::vector<asl::SPNumMethod> bc;
187  std::vector<asl::SPNumMethod> bcV;
188  std::vector<asl::SPNumMethod> bcDif;
189 
190  bc.push_back(generateBCNoSlip(lbgk, mpfMapMem));
191  bc.push_back(generateBCConstantPressure(lbgk,1.,{asl::ZE}));
192  bc.push_back(generateBCConstantPressureVelocity(lbgk, 1.,
193  makeAVec(0.,0.,params.oilInVel.v()),
194  {asl::Z0}));
195  bc.push_back(generateBCConstantPressureVelocity(lbgk, 1.,
196  makeAVec(0.,-params.waterInVel.v(),0.),
197  {asl::YE}));
198 
199  bcDif.push_back(generateBCNoSlipVel(lbgk, mpfMapMem));
200  bc.push_back(generateBCConstantGradient(waterFrac, 0., mpfMapMem, templ));
201  bc.push_back(generateBCConstantValue(waterFrac, 1., {asl::Y0, asl::YE}));
202  bc.push_back(generateBCConstantValue(waterFrac, 0., {asl::Z0, asl::ZE}));
203 
204  initAll(bc);
205  initAll(bcDif);
206  initAll(bcV);
207 
208  std::cout << "Finished" << endl;
209  std::cout << "Computing..." << endl;
210  asl::Timer timer;
211 
212  asl::WriterVTKXML writer(params.getDir() + "multiphase_flow");
213  writer.addScalars("map", *mpfMapMem);
214  writer.addScalars("water", *waterFrac);
215  writer.addScalars("rho", *lbgk->getRho());
216  writer.addVector("v", *flowVel);
217 
218  executeAll(bc);
219  executeAll(bcDif);
220  executeAll(bcV);
221 
222  writer.write();
223 
224  timer.start();
225  for (unsigned int i(1); i < 2001; ++i)
226  {
227  lbgk->execute();
228  executeAll(bcDif);
229  nmWater->execute();
230  executeAll(bc);
231 
232  if (!(i%200))
233  {
234  timer.stop();
235  cout << i << "/2000; time left (estimated): " << timer.estimatedRemainder(double(i)/2000.) << endl;
236  executeAll(bcV);
237  writer.write();
238  timer.start();
239  }
240  }
241  timer.stop();
242 
243  cout << "Finished" << endl;
244 
245  cout << "Computation statistic:" << endl;
246  cout << "Real Time = " << timer.realTime() << "; Processor Time = "
247  << timer.processorTime() << "; Processor Load = "
248  << timer.processorLoad() * 100 << "%" << endl;
249 
250  return 0;
251 }
asl::Parameter< double > pumpD
asl::Parameter< double > tubeL
asl::Parameter< double > tubeD
asl::Parameter< double > tOutput
asl::Parameter< double > pumpL
asl::Parameter< double > oilInVel
asl::UValue< double > nuNum
string getDir()
asl::Parameter< double > tSimulation
asl::Parameter< double > nu
asl::Parameter< double > waterInVel
asl::Block::DV size
asl::UValue< double > dt
asl::Parameter< double > gasInVel
asl::Parameter< double > dx
asl::ApplicationParametersManager appParamsManager
void updateNumValues()
void load(int argc, char *argv[])
AVec< T > makeAVec(T a1)
void load(int argc, char *argv[])
Numerical method for fluid flow.
Definition: aslLBGK.h:78
contains different kernels for preprocessing and posprocessing of data used by LBGK
Definition: aslLBGK.h:138
const T & v() const
std::string getDir()
const double realTime() const
Definition: aslTimer.h:45
void stop()
Definition: aslTimer.h:44
const double processorTime() const
Definition: aslTimer.h:46
void start()
Definition: aslTimer.h:43
const double processorLoad() const
Definition: aslTimer.h:47
const double estimatedRemainder(double completeness)
Returns estimated time till finishing current task based on its current completeness [0....
Definition: aslTimer.h:52
const T & v() const
Definition: aslUValue.h:43
void addVector(std::string name, AbstractData &data)
void addScalars(std::string name, AbstractData &data)
@ Y0
Definition: aslBCond.h:309
@ ZE
Definition: aslBCond.h:309
@ Z0
Definition: aslBCond.h:309
@ YE
Definition: aslBCond.h:309
acl::VectorOfElements dx(const TemplateVE &a)
differential operator
void errorMessage(cl_int status, const char *errorMessage)
Prints errorMessage and exits depending on the status.
SPBCond generateBCConstantGradient(SPAbstractDataWithGhostNodes d, double v, const VectorTemplate *const t, const std::vector< SlicesNames > &sl)
Bondary condition that makes fixed gradient <>
SPBCond generateBCConstantValue(SPAbstractDataWithGhostNodes d, double v, const std::vector< SlicesNames > &sl)
Bondary condition that puts fixed value in each point.
SPDistanceFunction generateDFInBlock(const Block &b, unsigned int nG)
generates map corresponding to external (ghost) part of the block
SPDistanceFunction normalize(SPDistanceFunction a, double dx)
SPDistanceFunction generateDFCylinderInf(double r, const AVec< double > &l, const AVec< double > &c)
generates infinite cylinder
std::shared_ptr< DistanceFunction > SPDistanceFunction
Definition: aslGeomInc.h:44
SPFDMultiPhase generateFDMultiPhase(SPDataWithGhostNodesACLData c, SPAbstractDataWithGhostNodes v, const VectorTemplate *vt, bool compressibilityCorrection=false)
const VectorTemplate & d3q15()
Vector template.
SPBCond generateBCConstantPressure(SPLBGK nm, double p, const std::vector< SlicesNames > &sl)
SPNumMethod generateBCNoSlipVel(SPLBGK nmU, SPAbstractDataWithGhostNodes map)
SPBCond generateBCConstantPressureVelocity(SPLBGK nm, double p, AVec<> v, const std::vector< SlicesNames > &sl)
SPBCond generateBCNoSlip(SPLBGK nm, const std::vector< SlicesNames > &sl)
VectorOfElements generateVEConstant(T a)
Generates VectorOfElements with 1 Element acl::Constant with value a.
int main(int argc, char *argv[])
float FlT
asl::UValue< double > Param
asl::SPDistanceFunction generateMixer(asl::Block &block, Parameters &params)
AVec< T > makeAVec(T a1)
void initAll(std::vector< T * > &v)
Definition: aslNumMethod.h:67
std::shared_ptr< LBGK > SPLBGK
Definition: aslLBGK.h:133
void initData(SPAbstractData d, double a)
std::shared_ptr< LBGKUtilities > SPLBGKUtilities
Definition: aslLBGK.h:161
void executeAll(std::vector< T * > &v)
Definition: aslNumMethod.h:55