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