ASL  0.1.7
Advanced Simulation Library
locomotive.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 
31 #include <math/aslTemplates.h>
32 #include <aslGeomInc.h>
34 #include <aslDataInc.h>
35 #include <acl/aclGenerators.h>
37 #include <num/aslLBGK.h>
38 #include <num/aslLBGKBC.h>
39 #include <utilities/aslTimer.h>
41 
42 
43 // typedef to switch to double precision
44 //typedef double FlT;
45 typedef float FlT;
46 
47 using asl::AVec;
48 using asl::makeAVec;
49 
50 // Generate geometry of the tunnel
52 {
53 
54  // Set length of the tunnel to the length (X size) of the block
55  double l(bl.getBPosition()[0] - bl.position[0] + bl.dx);
56  // Set radius of the tunnel to the ca. half of the block's height (Z size)
57  double rTunnel((bl.getBPosition()[2] - bl.position[2]) / 2.1);
58 
59  // Center of the tunnel (described as cylinder cut by a plane)
60  asl::AVec<> center(.5 * (bl.getBPosition() + bl.position));
61  center[1] = bl.position[1] + .25 * rTunnel;
62 
63  // Center of the ground plane (that cuts the cylinder)
64  asl::AVec<> centerG(center);
65  centerG[1] = bl.position[1];
66 
67  /* DF = DistanceFunction (part of the geometrical module of ASL)
68  1. Genarate cylinder
69  2. Generate ground plane
70  3. Conjunction of the cylinder and the plane ('&' - operator)
71  4. Space inversion ('-' - operator) */
72  auto tunnel(-(generateDFCylinder(rTunnel, makeAVec(l, 0., 0.), center) &
73  generateDFPlane(makeAVec(0., -1., 0.), centerG)));
74 
75  // Normalize DistanceFunction to the range [-1; 1]
76  return normalize(tunnel, bl.dx);
77 }
78 
79 
80 int main(int argc, char* argv[])
81 {
82  /* Convenience facility to manage simulation parameters (and also
83  hardware parameters defining platform and device for computations)
84  through command line and/or parameters file.
85  See `locomotive --help` for more information */
86  asl::ApplicationParametersManager appParamsManager("locomotive",
87  "1.0");
88 
89  /* Important: declare Parameters only after declaring
90  ApplicationParametersManager instance because each Parameter adds itself
91  to it automatically!
92  0.08 - default value; will be used if nothing else is provided during
93  runtime through command line or parameters file.
94  "dx" - option key; is used to specify this parameter through command line
95  and/or parameters file, like `locomotive --dx 0.05`
96  "space step" - option description; is used in the help output:
97  `locomotive -h` and as comment on parameters file generation:
98  `locomotive -g ./defaultParameters.ini`
99  "m" - parameter units; is used to complement the option description mentioned
100  above. Might be used for automatic unit conversion in future (to this end
101  it is recommended to use the notation of the Boost::Units library). */
102  asl::Parameter<FlT> dx(0.08, "dx", "space step", "m");
103  asl::Parameter<FlT> dt(1., "dt", "time step", "s");
104  asl::Parameter<FlT> nu(.001, "nu", "kinematic viscosity", "m^2/s");
105  asl::Parameter<unsigned int> iterations(10001, "iterations", "iterations number");
106  asl::Parameter<string> input("input", "path to the geometry input file");
107 
108  /* Load previously declared Parameters from command line and/or
109  parameters file. Use default values if neither is provided. */
110  appParamsManager.load(argc, argv);
111 
112  /* Set the size of the block to 40x10x15 m (in accordance with the
113  locomotive size read later on from the input file) */
114  AVec<int> size(makeAVec(40., 10., 15.) * (1. / dx.v()));
115 
116  /* Create block and shift it in accordance with the
117  position of the locomotive in the input file */
118  asl::Block bl(size, dx.v(), makeAVec(-30., 8.58, 1.53));
119 
120  // Define dimensionless viscosity value
121  FlT nuNum(nu.v() * dt.v() / dx.v() / dx.v());
122 
123  cout << "Data initialization... " << flush;
124 
125  // Read geometry of the locomotive from the file `locomotive.stl`
126  auto locomotive(asl::readSurface(input.v(), bl));
127 
128  // Create block for further use
129  asl::Block block(locomotive->getInternalBlock());
130 
131  // Generate memory data container for the tunnel
132  auto tunnelMap(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
133  // Place generated geometry of the tunnel into the tunnel data container
134  asl::initData(tunnelMap, generateTunnel(block));
135 
136  // Data container for air friction field
137  auto forceField(asl::generateDataContainerACL_SP<FlT>(block, 3, 1u));
138  // Initialization
139  asl::initData(forceField, makeAVec(0., 0., 0.));
140 
141  cout << "Finished" << endl;
142 
143  cout << "Numerics initialization... " << flush;
144 
145  // NOTE: the problem is considered in the reference frame related to the locomotive
146 
147  // Generate numerical method for air flow - LBGK (lattice Bhatnagar–Gross–Krook)
148  asl::SPLBGK lbgk(new asl::LBGKTurbulence(block,
149  acl::generateVEConstant(FlT(nu.v())),
150  &asl::d3q15()));
151  lbgk->init();
152  // Generate an instance for LBGK data initialization
153  asl::SPLBGKUtilities lbgkUtil(new asl::LBGKUtilities(lbgk));
154  // Initialize the LBGK internal data with the flow velocity of (0.1, 0, 0) in [lattice units]
155  lbgkUtil->initF(acl::generateVEConstant(.1, .0, .0));
156 
157  auto vfTunnel(asl::generatePFConstant(makeAVec(0.1, 0., 0.)));
158 
159  std::vector<asl::SPNumMethod> bc;
160  std::vector<asl::SPNumMethod> bcV;
161 
162  // Generate boundary conditions for the tunnel geometry. Constant velocity BC
163  bc.push_back(generateBCVelocity(lbgk, vfTunnel, tunnelMap));
164  // Generate boundary conditions for the tunnel geometry. Constant velocity BC
165  // This BC is used for visualization.
166  bcV.push_back(generateBCVelocityVel(lbgk, vfTunnel, tunnelMap));
167  bcV.push_back(generateBCNoSlipRho(lbgk, tunnelMap));
168 
169  // Generate boundary conditions for the locomotive geometry. Non-slip BC
170  bc.push_back(generateBCNoSlip(lbgk, locomotive));
171  bcV.push_back(generateBCNoSlipVel(lbgk, locomotive));
172  bcV.push_back(generateBCNoSlipRho(lbgk, locomotive));
173 
174  // Generate constant presure BC for in and out planes of the tunnel
175  bc.push_back(generateBCConstantPressureVelocity(lbgk, 1.,
176  makeAVec(0.1, 0., 0.),
177  {asl::X0, asl::XE}));
178 
179  // Initialization and building of all BC
180  initAll(bc);
181  initAll(bcV);
182 
183  // Generate a numerical method for computation of the air force field that acts on the locomotive
184  auto computeForce(generateComputeSurfaceForce(lbgk, forceField, locomotive));
185  computeForce->init();
186 
187  cout << "Finished" << endl;
188  cout << "Computing..." << endl;
189  asl::Timer timer;
190 
191  // Initialization of the output system
192  // Write the output to the directory containing the input parameters file (default "./")
193  asl::WriterVTKXML writer(appParamsManager.getDir() + "locomotive");
194  writer.addScalars("map", *locomotive);
195  writer.addScalars("tunnel", *tunnelMap);
196  writer.addScalars("rho", *lbgk->getRho());
197  writer.addVector("v", *lbgk->getVelocity());
198  writer.addVector("force", *forceField);
199 
200  // Execute all BC
201  executeAll(bc);
202  executeAll(bcV);
203  computeForce->execute();
204 
205  // First data output
206  writer.write();
207 
208  timer.start();
209  // Iteration loop
210  for (unsigned int i(1); i < iterations.v(); ++i)
211  {
212  // One iteration (timestep) of bulk numerical procedure
213  lbgk->execute();
214  // Execution of the BC procedures
215  executeAll(bc);
216  // Output and analysis scope
217  if (!(i%1000))
218  {
219  cout << i << endl;
220  // Execution of the visualization BC procedures
221  executeAll(bcV);
222  // Computation of the force field
223  computeForce->execute();
224  // Data writing
225  writer.write();
226  }
227  }
228  timer.stop();
229 
230  cout << "Finished" << endl;
231 
232  cout << "Computation statistic:" << endl;
233  cout << "Real Time = " << timer.realTime() << "; Processor Time = "
234  << timer.processorTime() << "; Processor Load = "
235  << timer.processorLoad() * 100 << "%" << endl;
236 
237  return 0;
238 }
int main(int argc, char *argv[])
Definition: locomotive.cc:80
const double realTime() const
Definition: aslTimer.h:45
std::shared_ptr< DistanceFunction > SPDistanceFunction
Definition: aslGeomInc.h:44
float FlT
const double processorTime() const
Definition: aslTimer.h:46
AVec< T > makeAVec(T a1)
const AVec normalize(const AVec< T > &a)
SPNumMethod generateBCNoSlipRho(SPLBGK nmU, SPAbstractDataWithGhostNodes map)
void initAll(std::vector< T * > &v)
Definition: aslNumMethod.h:67
AVec< T > makeAVec(T a1)
std::string getDir()
SPDataWithGhostNodesACLData readSurface(const string &fileName, double dx, acl::CommandQueue queue=acl::hardware.defaultQueue)
SPBCond generateBCConstantPressureVelocity(SPLBGK nm, double p, AVec<> v, const std::vector< SlicesNames > &sl)
double dx
Definition: aslBlocks.h:66
SPDistanceFunction generateDFPlane(const AVec< double > &n, const AVec< double > &p0)
std::shared_ptr< LBGKUtilities > SPLBGKUtilities
Definition: aslLBGK.h:161
std::shared_ptr< LBGK > SPLBGK
Definition: aslLBGK.h:133
const VectorTemplate & d3q15()
Vector template.
const T & v() const
cl_int flush(void)
Definition: cl.hpp:7042
acl::VectorOfElements dx(const TemplateVE &a)
differential operator
asl::SPDistanceFunction generateTunnel(asl::Block &bl)
Definition: locomotive.cc:51
SPNumMethod generateBCVelocityVel(SPLBGK nm, SPPositionFunction v, SPAbstractDataWithGhostNodes map)
void initData(SPAbstractData d, double a)
void executeAll(std::vector< T * > &v)
Definition: aslNumMethod.h:55
void addScalars(std::string name, AbstractData &data)
void stop()
Definition: aslTimer.h:44
void start()
Definition: aslTimer.h:43
VectorOfElements generateVEConstant(T a)
Generates VectorOfElements with 1 Element acl::Constant with value a.
SPNumMethod generateComputeSurfaceForce(SPLBGK nm, SPDataWithGhostNodesACLData fF, SPAbstractDataWithGhostNodes map)
float FlT
Definition: locomotive.cc:45
const V getBPosition() const
Definition: aslBlocks.h:214
SPDistanceFunction generateDFCylinder(double r, const AVec< double > &l, const AVec< double > &c)
generates cylinder
void load(int argc, char *argv[])
const double processorLoad() const
Definition: aslTimer.h:47
SPNumMethod generateBCNoSlipVel(SPLBGK nmU, SPAbstractDataWithGhostNodes map)
for velocity field
SPPositionFunction generatePFConstant(const AVec< double > &a)
contains different kernels for preprocessing and posprocessing of data used by LBGK
Definition: aslLBGK.h:137
SPNumMethod generateBCVelocity(SPLBGK nm, SPPositionFunction v, SPAbstractDataWithGhostNodes map)
SPBCond generateBCNoSlip(SPLBGK nm, const std::vector< SlicesNames > &sl)