ASL
0.1.7
Advanced Simulation Library
flowKDPGrowth.cc
Example: flowKDPGrowth
/*
* Advanced Simulation Library <http://asl.org.il>
*
* Copyright 2015 Avtech Scientific <http://avtechscientific.com>
*
*
* This file is part of Advanced Simulation Library (ASL).
*
* ASL is free software: you can redistribute it and/or modify it
* under the terms of the GNU Affero General Public License as
* published by the Free Software Foundation, version 3 of the License.
*
* ASL is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Affero General Public License for more details.
*
* You should have received a copy of the GNU Affero General Public License
* along with ASL. If not, see <http://www.gnu.org/licenses/>.
*
*/
#include <
utilities/aslParametersManager.h
>
#include <
math/aslTemplates.h
>
#include <
aslGeomInc.h
>
#include <
math/aslPositionFunction.h
>
#include <
aslDataInc.h
>
#include <
acl/aclGenerators.h
>
#include <
acl/aclMath/aclVectorOfElements.h
>
#include <
writers/aslVTKFormatWriters.h
>
#include <
num/aslLBGK.h
>
#include <
num/aslLBGKBC.h
>
#include <
num/aslBasicBC.h
>
#include <
num/aslCrystalGrowthBC.h
>
#include <
num/aslFDAdvectionDiffusion.h
>
#include <
utilities/aslTimer.h
>
using
asl::AVec
;
using
asl::makeAVec
;
asl::SPDistanceFunction
generateBath
(
asl::Block
& bl)
{
// double hBath(2.);
double
rBath(1.);
double
dx
(bl.
dx
);
asl::AVec<int>
size(bl.
getSize
());
asl::AVec<>
center(.5*
dx
*AVec<>(size));
auto
bath(-(
generateDFCylinderInf
(rBath,
makeAVec
(0.,0.,1.),
dx
*AVec<>(size)*.5) &
generateDFPlane
(
makeAVec
(0.,0.,1.), center*1.99) &
generateDFPlane
(
makeAVec
(0.,0.,-1.), center*0.)));
return
normalize
(bath,
dx
);
}
asl::SPDistanceFunction
generatePlatform
(
asl::Block
& bl)
{
double
rDisk(.9);
double
hDisk(0.1);
double
rAxis(0.05);
double
hAxis(.5);
double
wPillar(.2);
double
dPillar(.1);
double
dx
(bl.
dx
);
asl::AVec<int>
size(bl.
getSize
());
asl::AVec<>
center(.5*
dx
*AVec<>(size));
vector<asl::AVec<>> pillar1{
makeAVec
(wPillar*.5, dPillar*.5,0.),
makeAVec
(-wPillar*.5, dPillar*.5,0.),
makeAVec
(-wPillar*.5, -dPillar*.5,0.),
makeAVec
(wPillar*.5, -dPillar*.5,0.)};
vector<asl::AVec<>> pillar2{
makeAVec
(dPillar*.5, wPillar*.5,0.),
makeAVec
(-dPillar*.5, wPillar*.5,0.),
makeAVec
(-dPillar*.5, -wPillar*.5,0.),
makeAVec
(dPillar*.5, -wPillar*.5,0.)};
vector<asl::AVec<>> pillarC{
makeAVec
(center[0]+rDisk-dPillar*.5, center[1], 0.),
makeAVec
(center[0]-rDisk+dPillar*.5, center[1], 0.),
makeAVec
(center[0], center[1]+rDisk-dPillar*.5,0.),
makeAVec
(center[0], center[1]-rDisk+dPillar*.5,0.)};
vector<vector<asl::AVec<>>> pillarsPoints(4);
for
(
unsigned
int
i(0); i<4; ++i)
pillarsPoints[i].resize(4);
for
(
unsigned
int
i(0); i<4; ++i)
{
pillarsPoints[0][i] = pillar2[i] + pillarC[0];
pillarsPoints[1][i] = pillar2[i] + pillarC[1];
pillarsPoints[2][i] = pillar1[i] + pillarC[2];
pillarsPoints[3][i] = pillar1[i] + pillarC[3];
}
auto
diskBottom(
generateDFCylinder
(rDisk,
makeAVec
(0., 0., hDisk),
makeAVec
(center[0], center[1], .5*hDisk)));
auto
diskTop(
generateDFCylinder
(rDisk,
makeAVec
(0., 0., hDisk),
makeAVec
(center[0], center[1], -.5*hDisk - hAxis +
dx
*size[2])));
auto
axis(
generateDFCylinder
(rAxis,
makeAVec
(0., 0., hAxis+hDisk*.5),
makeAVec
(center[0], center[1], - .5*hAxis - hDisk*.25 +
dx
*size[2])));
auto
dfPillar1(
generateDFConvexPolygonPrism
(pillarsPoints[0]));
auto
dfPillar2(
generateDFConvexPolygonPrism
(pillarsPoints[1]));
auto
dfPillar3(
generateDFConvexPolygonPrism
(pillarsPoints[2]));
auto
dfPillar4(
generateDFConvexPolygonPrism
(pillarsPoints[3]));
auto
dfPillars((dfPillar1 | dfPillar2 | dfPillar3 | dfPillar4) &
generateDFPlane
(
makeAVec
(0.,0.,-1.),
makeAVec
(0.,0.,.5*hDisk)) &
generateDFPlane
(
makeAVec
(0.,0.,1.),
makeAVec
(0.,0.,-.5*hDisk - hAxis +
dx
*size[2])));
return
normalize
(diskBottom | diskTop | axis | dfPillars,
dx
);
}
asl::SPDistanceFunction
generateCrystal
(
asl::Block
& bl)
{
double
aCrystal(.5);
double
hCrystalBase(.5);
double
hCrystalPyramid(.5);
double
hDisk(0.1);
double
dx
(bl.
dx
);
asl::AVec<int>
size(bl.
getSize
());
asl::AVec<>
center(.5*
dx
*AVec<>(size));
auto
crystalB(
asl::generateDFConvexPolygonPrism
({center+
makeAVec
( aCrystal, aCrystal,0.),
center+
makeAVec
(-aCrystal, aCrystal,0.),
center+
makeAVec
(-aCrystal, -aCrystal,0.),
center+
makeAVec
( aCrystal, -aCrystal,0.)}) &
generateDFPlane
(
makeAVec
(0.,0.,-1.),
makeAVec
(0.,0., hDisk-.001)) &
generateDFPlane
(
makeAVec
(0.,0., 1.),
makeAVec
(0.,0., hDisk + hCrystalBase)));
auto
cCrPyrBase(
makeAVec
(center[0],center[1],hDisk+hCrystalBase-.01));
auto
crystalT(
asl::generateDFConvexPolygonPyramid
({cCrPyrBase+
makeAVec
( aCrystal, aCrystal,0.),
cCrPyrBase+
makeAVec
(-aCrystal, aCrystal,0.),
cCrPyrBase+
makeAVec
(-aCrystal, -aCrystal,0.),
cCrPyrBase+
makeAVec
( aCrystal, -aCrystal,0.)},
cCrPyrBase+
makeAVec
(0.,0.,hCrystalPyramid)));
return
normalize
(crystalB | crystalT,
dx
);
// return crystalB | crystalT;
}
double
getWRotation
(
double
t)
{
double
tPeriod(128);
double
wMax(6.*3.14*2./60.);
double
tPlato(tPeriod * .25);
double
tAcceleration(tPeriod * .1);
double
tStop(tPeriod * .05);
double
intPart;
double
tRel(modf(t/tPeriod, &intPart));
double
x(0);
if
(tRel<=tAcceleration)
x = tRel / tAcceleration;
if
(tRel>tAcceleration && tRel<=tAcceleration+tPlato)
x = 1.;
if
(tRel>tAcceleration+tPlato && tRel<=2.*tAcceleration+tPlato)
x = (2.*tAcceleration + tPlato - tRel) / tAcceleration;
if
(tRel>2.*tAcceleration+tPlato && tRel<=2.*tAcceleration+tPlato+tStop)
x = 0;
if
(tRel>2.*tAcceleration+tPlato+tStop && tRel<=3.*tAcceleration+tPlato+tStop)
x = -(tRel-2.*tAcceleration-tPlato-tStop) / tAcceleration;
if
(tRel>3.*tAcceleration+tPlato+tStop && tRel<=3.*tAcceleration+2.*tPlato+tStop)
x = -1.;
if
(tRel>3.*tAcceleration+2.*tPlato+tStop && tRel<=4.*tAcceleration+2.*tPlato+tStop)
x = -(4.*tAcceleration+2.*tPlato+tStop-tRel)/tAcceleration;
if
(tRel>4.*tAcceleration+2.*tPlato+tStop)
x = 0;
return
wMax*x;
// flux = -9.32e-5*(1.170 - c); c_0=0.326 ceq=0.267
}
typedef
float
FlT
;
//typedef double FlT;
typedef
asl::UValue<double>
Param
;
using
asl::AVec
;
using
asl::makeAVec
;
int
main
(
int
argc,
char
* argv[])
{
// Optionally add appParamsManager to be able to manipulate at least
// hardware parameters(platform/device) through command line/parameters file
asl::ApplicationParametersManager
appParamsManager(
"flowKDPGrowth"
,
"1.0"
);
appParamsManager.load(argc, argv);
Param
dx
(.02);
Param
dt(0.8e-2);
Param
nu(1e-2);
Param
difC(1e-2/300.);
// Param w(48.*3.14*2./60.);
// Angular velocity
Param
w(6.*3.14*2./60.);
Param
nuNum(nu.v()*dt.v()/
dx
.v()/
dx
.v());
Param
difCNum(difC.v()*dt.v()/
dx
.v()/
dx
.v());
// Angular velocity in one iteration
Param
wNum(w.v()*dt.v());
Param
c0(0.326);
AVec<int> size(
asl::makeAVec
(105.,105.,100.));
AVec<> gSize(
dx
.v()*AVec<>(size));
std::cout <<
"Data initialization..."
;
auto
templ(&
asl::d3q19
());
asl::Block
block(size,
dx
.v());
auto
bathMap(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
asl::initData
(bathMap,
generateBath
(block));
auto
platformCrysMap(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
asl::initData
(platformCrysMap,
generatePlatform
(block) |
generateCrystal
(block));
auto
bathPlatformMap(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
asl::initData
(bathPlatformMap,
generateBath
(block) |
generatePlatform
(block));
auto
bathPlatformCrystalMap(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
asl::initData
(bathPlatformCrystalMap,
generateBath
(block) |
generatePlatform
(block) |
generateCrystal
(block));
auto
crystalMap(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
asl::initData
(crystalMap,
generateCrystal
(block));
auto
cField(asl::generateDataContainerACL_SP<FlT>(block, 1, 1u));
asl::initData
(cField, c0.v());
std::cout <<
"Finished"
<< endl;
std::cout <<
"Numerics initialization..."
;
asl::SPLBGK
lbgk(
new
asl::LBGK
(block,
acl::generateVEConstant
(
FlT
(nuNum.v())),
templ));
// Set angular velocity in lbgk
lbgk->setOmega(
acl::generateVEConstant
(
makeAVec
(0.,0.,wNum.v())));
lbgk->init();
asl::SPLBGKUtilities
lbgkUtil(
new
asl::LBGKUtilities
(lbgk));
lbgkUtil->initF(
acl::generateVEConstant
(.0,.0,.0));
auto
nmDif(
asl::generateFDAdvectionDiffusion
(cField,
difCNum.v(),
lbgk->getVelocity(),
templ,
true
));
nmDif->init();
std::vector<asl::SPNumMethod> bc;
std::vector<asl::SPNumMethod> bcV;
std::vector<asl::SPNumMethod> bcDif;
// Position Function Angular Velocity Field
auto
vfBath(
asl::generatePFRotationField
(
makeAVec
(0.,0., wNum.v()/
dx
.v()), .5*gSize));
// Boundary condition
bc.push_back(
generateBCVelocity
(lbgk, vfBath, bathMap));
// Boundary condition for visualization
bcV.push_back(
generateBCVelocityVel
(lbgk, vfBath, bathMap));
bc.push_back(
generateBCNoSlip
(lbgk, platformCrysMap));
bcV.push_back(
generateBCNoSlipVel
(lbgk, platformCrysMap));
bcDif.push_back(
generateBCConstantGradient
(cField,
0.,
bathPlatformMap,
bathPlatformCrystalMap,
templ));
bcDif.push_back(
generateBCLinearGrowth2
(cField, 1.17,
-9.32e-6/difC.v()*
dx
.v(),
crystalMap,
bathPlatformCrystalMap,
templ));
// bcDif.push_back(generateBCConstantGradient2(cField, .1, crystalMap, bathPlatformCrystalMap, templ));
initAll
(bc);
initAll
(bcV);
initAll
(bcDif);
std::cout <<
"Finished"
<< endl;
std::cout <<
"Computing..."
;
asl::Timer
timer;
asl::Timer
timerBC;
asl::WriterVTKXML
writer(
"flowKDPGrowthRes"
);
writer.addScalars(
"mapBath"
, *bathMap);
writer.addScalars(
"mapPlatformCrys"
, *platformCrysMap);
writer.addScalars(
"mapBathPlatformCrystal"
, *bathPlatformCrystalMap);
writer.addScalars(
"mapCrys"
, *crystalMap);
writer.addScalars(
"rho"
, *lbgk->getRho());
writer.addScalars(
"c"
, *cField);
writer.addVector(
"v"
, *lbgk->getVelocity());
executeAll
(bcV);
executeAll
(bc);
executeAll
(bcDif);
writer.write();
timer.
start
();
timerBC.
reset
();
for
(
unsigned
int
i(0); i <= 8001 ; ++i)
{
lbgk->execute();
timerBC.
start
();
executeAll
(bcV);
executeAll
(bc);
timerBC.
stop
();
nmDif->execute();
timerBC.
start
();
executeAll
(bcDif);
timerBC.
stop
();
if
(!(i%2000))
{
cout << i << endl;
writer.write();
}
}
timer.
stop
();
cout <<
"Finished"
<< endl;
cout <<
"Computation statistic:"
<< endl;
cout <<
"Real Time = "
<< timer.
realTime
() <<
"; Processor Time = "
<< timer.
processorTime
() <<
"; Processor Load = "
<< timer.
processorLoad
() * 100 <<
"%"
<< endl;
return
0;
}
aclGenerators.h
aclVectorOfElements.h
aslBasicBC.h
aslCrystalGrowthBC.h
aslDataInc.h
aslFDAdvectionDiffusion.h
aslGeomInc.h
aslLBGK.h
aslLBGKBC.h
aslParametersManager.h
aslPositionFunction.h
aslTemplates.h
aslTimer.h
aslVTKFormatWriters.h
asl::AVec
Definition:
aslVectorsDynamicLength.h:40
asl::AVec::normalize
const AVec normalize(const AVec< T > &a)
Definition:
aslVectorsDynamicLengthOperations.h:152
asl::AVec::makeAVec
AVec< T > makeAVec(T a1)
Definition:
aslVectorsDynamicLength.h:176
asl::ApplicationParametersManager
Definition:
aslParametersManager.h:159
asl::Block
Definition:
aslBlocks.h:57
asl::Block::dx
double dx
Definition:
aslBlocks.h:66
asl::Block::getSize
const DV & getSize() const
Definition:
aslBlocks.h:208
asl::LBGK
Numerical method for fluid flow.
Definition:
aslLBGK.h:78
asl::LBGKUtilities
contains different kernels for preprocessing and posprocessing of data used by LBGK
Definition:
aslLBGK.h:138
asl::Timer
Definition:
aslTimer.h:34
asl::Timer::realTime
const double realTime() const
Definition:
aslTimer.h:45
asl::Timer::stop
void stop()
Definition:
aslTimer.h:44
asl::Timer::processorTime
const double processorTime() const
Definition:
aslTimer.h:46
asl::Timer::start
void start()
Definition:
aslTimer.h:43
asl::Timer::reset
void reset()
Definition:
aslTimer.h:48
asl::Timer::processorLoad
const double processorLoad() const
Definition:
aslTimer.h:47
asl::UValue< double >
asl::WriterVTKXML
Definition:
aslVTKFormatWriters.h:42
main
int main(int argc, char *argv[])
Definition:
flowKDPGrowth.cc:200
FlT
float FlT
Definition:
flowKDPGrowth.cc:192
generateCrystal
asl::SPDistanceFunction generateCrystal(asl::Block &bl)
Definition:
flowKDPGrowth.cc:130
generateBath
asl::SPDistanceFunction generateBath(asl::Block &bl)
Definition:
flowKDPGrowth.cc:47
Param
asl::UValue< double > Param
Definition:
flowKDPGrowth.cc:194
getWRotation
double getWRotation(double t)
Definition:
flowKDPGrowth.cc:159
generatePlatform
asl::SPDistanceFunction generatePlatform(asl::Block &bl)
Definition:
flowKDPGrowth.cc:68
asl::dx
acl::VectorOfElements dx(const TemplateVE &a)
differential operator
asl::generateBCConstantGradient
SPBCond generateBCConstantGradient(SPAbstractDataWithGhostNodes d, double v, const VectorTemplate *const t, const std::vector< SlicesNames > &sl)
Bondary condition that makes fixed gradient <>
asl::generateDFConvexPolygonPrism
SPDistanceFunction generateDFConvexPolygonPrism(std::vector< AVec< double >> points)
generates infinite prism with convex polygon at its base
asl::generateDFConvexPolygonPyramid
SPDistanceFunction generateDFConvexPolygonPyramid(std::vector< AVec< double >> points, AVec< double > a)
generates pyramid with convex polygon at its base and apex a
asl::generateDFCylinder
SPDistanceFunction generateDFCylinder(double r, const AVec< double > &l, const AVec< double > &c)
generates cylinder
asl::generateDFPlane
SPDistanceFunction generateDFPlane(const AVec< double > &n, const AVec< double > &p0)
asl::generateDFCylinderInf
SPDistanceFunction generateDFCylinderInf(double r, const AVec< double > &l, const AVec< double > &c)
generates infinite cylinder
asl::SPDistanceFunction
std::shared_ptr< DistanceFunction > SPDistanceFunction
Definition:
aslGeomInc.h:44
asl::generateFDAdvectionDiffusion
SPFDAdvectionDiffusion generateFDAdvectionDiffusion(SPDataWithGhostNodesACLData c, double diffustionCoeff, SPAbstractDataWithGhostNodes v, const VectorTemplate *vt, bool compressibilityCorrection=false)
asl::generatePFRotationField
SPPositionFunction generatePFRotationField(const AVec< double > &axis, const AVec< double > &c)
asl::d3q19
const VectorTemplate & d3q19()
Vector template.
asl::generateBCVelocityVel
SPNumMethod generateBCVelocityVel(SPLBGK nm, SPPositionFunction v, SPAbstractDataWithGhostNodes map)
asl::generateBCNoSlipVel
SPNumMethod generateBCNoSlipVel(SPLBGK nmU, SPAbstractDataWithGhostNodes map)
asl::generateBCVelocity
SPNumMethod generateBCVelocity(SPLBGK nm, SPPositionFunction v, SPAbstractDataWithGhostNodes map)
asl::generateBCNoSlip
SPBCond generateBCNoSlip(SPLBGK nm, const std::vector< SlicesNames > &sl)
acl::generateVEConstant
VectorOfElements generateVEConstant(T a)
Generates VectorOfElements with 1 Element acl::Constant with value a.
asl::SPLBGKUtilities
std::shared_ptr< LBGKUtilities > SPLBGKUtilities
Definition:
aslLBGK.h:161
asl::generateBCLinearGrowth2
SPNumMethod generateBCLinearGrowth2(SPAbstractDataWithGhostNodes d, double cEq, double beta, SPAbstractDataWithGhostNodes map, const VectorTemplate *const t)
asl::initAll
void initAll(std::vector< T * > &v)
Definition:
aslNumMethod.h:67
asl::makeAVec
AVec< T > makeAVec(T a1)
Definition:
aslVectorsDynamicLength.h:176
asl::SPLBGK
std::shared_ptr< LBGK > SPLBGK
Definition:
aslLBGK.h:133
asl::initData
void initData(SPAbstractData d, double a)
asl::executeAll
void executeAll(std::vector< T * > &v)
Definition:
aslNumMethod.h:55
Generated by
1.9.1