Skip to content
Snippets Groups Projects
Commit ccb8cc57 authored by kutscher's avatar kutscher
Browse files

add flow past sphere and flow past cylinder benchmarks

parent 33f32301
No related branches found
No related tags found
1 merge request!18Feature/thixotropy
Hide whitespace changes
Inline Side-by-side
Showing
with 188 additions and 172 deletions
apps/cpu/Applications.cmake
+ 2
− 2
View file @ ccb8cc57
......@@ -11,8 +11,8 @@ add_subdirectory(${APPS_ROOT_CPU}/sphere)
#add_subdirectory(Applications/bananas2)
# add_subdirectory(Applications/plate)
# add_subdirectory(Applications/plate2)
##add_subdirectory(Applications/FlowAroundCylinder)
#add_subdirectory(${APPS_ROOT_CPU}/LaminarTubeFlow)
add_subdirectory(${APPS_ROOT_CPU}/FlowAroundCylinder)
add_subdirectory(${APPS_ROOT_CPU}/LaminarTubeFlow)
# add_subdirectory(Applications/LaminarTubeFlowConv)
#add_subdirectory(Applications/cylinderSt)
#add_subdirectory(Applications/mpichTest)
......
apps/cpu/FlowAroundCylinder/CMakeLists.txt
+ 2
− 17
View file @ ccb8cc57
......@@ -5,21 +5,6 @@ CMAKE_MINIMUM_REQUIRED(VERSION 2.8)
########################################################
PROJECT(cylinder)
INCLUDE(${APPS_ROOT}/IncludsList.cmake)
INCLUDE(${APPS_ROOT_CPU}/IncludsList.cmake)
#################################################################
### LOCAL FILES ###
#################################################################
FILE(GLOB SPECIFIC_FILES ${CMAKE_CURRENT_SOURCE_DIR}/*.h
${CMAKE_CURRENT_SOURCE_DIR}/*.cpp
${CMAKE_CURRENT_SOURCE_DIR}/*.hpp )
SET(ALL_SOURCES ${ALL_SOURCES} ${SPECIFIC_FILES})
SOURCE_GROUP(src FILES ${SPECIFIC_FILES})
SET(CAB_ADDITIONAL_LINK_LIBRARIES VirtualFluids)
#################################################################
### CREATE PROJECT ###
#################################################################
CREATE_CAB_PROJECT(cylinder BINARY)
vf_add_library(BUILDTYPE binary DEPENDS VirtualFluidsCore basics ${MPI_CXX_LIBRARIES} FILES cylinder.cpp )
apps/cpu/FlowAroundCylinder/cylinder.cfg
+ 4
− 4
View file @ ccb8cc57
......@@ -2,11 +2,11 @@ pathOut = d:/temp/cylinder_test
numOfThreads = 4
availMem = 15e9
refineLevel = 1
refineLevel = 0
blockNx = 25 41 41
uLB = 0.001
dx = 0.005
#dx = 0.01
#dx = 0.005
dx = 0.01
logToFile = false
......@@ -16,7 +16,7 @@ restartStep = 1000
cpStart = 1000
cpStep = 1000
outTime = 100
outTime = 10000
endTime = 100000
nupsStep = 100 100 10000000
\ No newline at end of file
apps/cpu/FlowAroundCylinder/cylinder.cpp
+ 11
− 8
View file @ ccb8cc57
......@@ -259,9 +259,9 @@ void run(string configname)
intHelper.setBC();
//domain decomposition
PQueuePartitioningGridVisitor pqPartVisitor(numOfThreads);
grid->accept(pqPartVisitor);
////domain decomposition
//PQueuePartitioningGridVisitor pqPartVisitor(numOfThreads);
//grid->accept(pqPartVisitor);
//initialization of distributions
InitDistributionsBlockVisitor initVisitor;
......@@ -297,10 +297,11 @@ void run(string configname)
SPtr<CoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, stepSch, pathOut, WbWriterVtkXmlBinary::getInstance(), conv, comm));
//double area = (2.0*radius*H)/(dx*dx);
//double v = 4.0*uLB/9.0;
//CalculateForcesCoProcessor fp(grid, stepSch, pathOut + "/results/forces.txt", comm, v, area);
//fp.addInteractor(cylinderInt);
double area = (2.0*radius*H)/(dx*dx);
double v = 4.0*uLB/9.0;
SPtr<UbScheduler> forceSch(new UbScheduler(100));
SPtr<CalculateForcesCoProcessor> fp = make_shared<CalculateForcesCoProcessor>(grid, forceSch, pathOut + "/results/forces.txt", comm, v, area);
fp->addInteractor(cylinderInt);
SPtr<UbScheduler> nupsSch(new UbScheduler(nupsStep[0], nupsStep[1], nupsStep[2]));
std::shared_ptr<CoProcessor> nupsCoProcessor(new NUPSCounterCoProcessor(grid, nupsSch, numOfThreads, comm));
......@@ -309,7 +310,9 @@ void run(string configname)
SPtr<UbScheduler> stepGhostLayer(new UbScheduler(1));
SPtr<Calculator> calculator(new BasicCalculator(grid, stepGhostLayer, endTime));
calculator->addCoProcessor(nupsCoProcessor);
calculator->addCoProcessor(writeMQCoProcessor);
calculator->addCoProcessor(fp);
calculator->addCoProcessor(writeMQCoProcessor);
if(myid == 0) UBLOG(logINFO,"Simulation-start");
calculator->calculate();
if(myid == 0) UBLOG(logINFO,"Simulation-end");
......
apps/cpu/HerschelBulkleySphere/hbsphere.cfg
+ 9
− 9
View file @ ccb8cc57
outputPath = d:/temp/Herschel-BulkleySphere
numOfThreads = 1
numOfThreads = 4
availMem = 8e9
logToFile = false
blocknx = 5 5 5
boundingBox = 40 40 40 #30*20=600**3=216000000
blocknx = 50 50 50
boundingBox = 150 50 50 #30*20=600**3=216000000
deltax = 1
refineLevel = 1
refineLevel = 0
radius = 5
velocity = 1e-3
n = 0.3
n = 0.9
Re = 1
Bn = 0.01
......@@ -20,8 +20,8 @@ Bn = 0.01
newStart = true
restartStep = 100000
cpStart = 10000
cpStep = 10000
cpStart = 100000
cpStep = 100000
outTime = 1
endTime = 1000
\ No newline at end of file
outTime = 10000
endTime = 1000000
\ No newline at end of file
apps/cpu/HerschelBulkleySphere/hbsphere.cpp
+ 31
− 21
View file @ ccb8cc57
......@@ -61,21 +61,21 @@ void bflow(string configname)
//bounding box
//double g_minX1 = 0;
//double g_minX2 = 0;
//double g_minX3 = 0;
double g_minX1 = 0;
double g_minX2 = 0;
double g_minX3 = 0;
//double g_maxX1 = boundingBox[0];
//double g_maxX2 = boundingBox[1];
//double g_maxX3 = boundingBox[2];
double g_maxX1 = boundingBox[0];
double g_maxX2 = boundingBox[1];
double g_maxX3 = boundingBox[2];
double g_minX1 = -boundingBox[0]/2.0;
double g_minX2 = -boundingBox[1] / 2.0;
double g_minX3 = -boundingBox[2]/2.0;
//double g_minX1 = -boundingBox[0]/2.0;
//double g_minX2 = -boundingBox[1] / 2.0;
//double g_minX3 = -boundingBox[2]/2.0;
double g_maxX1 = boundingBox[0]/2.0;
double g_maxX2 = boundingBox[1]/2.0;
double g_maxX3 = boundingBox[2]/2.0;
//double g_maxX1 = boundingBox[0]/2.0;
//double g_maxX2 = boundingBox[1]/2.0;
//double g_maxX3 = boundingBox[2]/2.0;
double blockLength = 3.0 * deltax;
......@@ -84,13 +84,20 @@ void bflow(string configname)
double Gamma = U / d;
double k = (U * d) / (Re * std::pow(Gamma, n - 1));
double tau0 = Bn * k * std::pow(Gamma, n);
//double k = 0.05; // (U * d) / (Re * std::pow(Gamma, n - 1));
//double tau0 = 3e-6; //Bn * k * std::pow(Gamma, n);
//double forcing = 8e-7;
double omegaMin = 1.0e-8;
SPtr<Thixotropy> thix = Thixotropy::getInstance();
thix->setPowerIndex(n);
thix->setViscosityParameter(k);
thix->setYieldStress(tau0);
thix->setOmegaMin(omegaMin);
SPtr<BCAdapter> noSlipBCAdapter(new NoSlipBCAdapter());
noSlipBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new HerschelBulkleyModelNoSlipBCAlgorithm()));
......@@ -117,6 +124,8 @@ void bflow(string configname)
SPtr<LBMKernel> kernel = SPtr<LBMKernel>(new HerschelBulkleyModelLBMKernel());
//SPtr<LBMKernel> kernel = SPtr<LBMKernel>(new BinghamModelLBMKernel());
kernel->setBCProcessor(bcProc);
//kernel->setForcingX1(forcing);
//kernel->setWithForcing(true);
SPtr<Grid3D> grid(new Grid3D(comm));
grid->setPeriodicX1(false);
......@@ -129,7 +138,7 @@ void bflow(string configname)
if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), outputPath + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance());
//sphere
SPtr<GbObject3D> sphere(new GbSphere3D(0, 0, 0, radius));
SPtr<GbObject3D> sphere(new GbSphere3D(75, 25, 25, radius));
GbSystem3D::writeGeoObject(sphere.get(), outputPath + "/geo/sphere", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> sphereInt(new D3Q27Interactor(sphere, grid, noSlipBCAdapter, Interactor3D::SOLID));
......@@ -147,13 +156,13 @@ void bflow(string configname)
//UBLOG(logINFO, "forcing = " << forcing);
UBLOG(logINFO, "rho = " << rhoLB);
UBLOG(logINFO, "U = " << U);
UBLOG(logINFO, "Re = " << Re);
UBLOG(logINFO, "Bn = " << Bn);
UBLOG(logINFO, "Re = " << (U * d) / (k * std::pow(Gamma, n - 1)));
UBLOG(logINFO, "Bn = " << tau0 / k * std::pow(Gamma, n));
UBLOG(logINFO, "k = " << k);
UBLOG(logINFO, "n = " << n);
UBLOG(logINFO, "tau0 = " << tau0);
UBLOG(logINFO, "deltax = " << deltax);
//UBLOG(logINFO, "number of levels = " << refineLevel + 1);
UBLOG(logINFO, "number of levels = " << refineLevel + 1);
UBLOG(logINFO, "number of threads = " << numOfThreads);
UBLOG(logINFO, "number of processes = " << comm->getNumberOfProcesses());
UBLOG(logINFO, "Preprozess - start");
......@@ -219,7 +228,7 @@ void bflow(string configname)
intHelper.addInteractor(wallYmaxInt);
intHelper.addInteractor(wallXminInt);
intHelper.addInteractor(wallXmaxInt);
//intHelper.addInteractor(sphereInt);
intHelper.addInteractor(sphereInt);
intHelper.selectBlocks();
if (myid == 0) UBLOG(logINFO, "deleteSolidBlocks - end");
//////////////////////////////////////
......@@ -281,6 +290,7 @@ void bflow(string configname)
//set connectors
InterpolationProcessorPtr iProcessor(new ThixotropyInterpolationProcessor());
static_pointer_cast<ThixotropyInterpolationProcessor>(iProcessor)->setOmegaMin(thix->getOmegaMin());
SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, k, iProcessor);
grid->accept(setConnsVisitor);
......@@ -297,10 +307,10 @@ void bflow(string configname)
SPtr<UbScheduler> visSch(new UbScheduler(outTime));
//SPtr<UbScheduler> visSch(new UbScheduler(10,1));
SPtr<WriteMacroscopicQuantitiesCoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, visSch, outputPath, WbWriterVtkXmlASCII::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm));
writeMQCoProcessor->process(0);
//writeMQCoProcessor->process(0);
double area = 4*UbMath::PI*radius*radius;
SPtr<UbScheduler> forceSch(new UbScheduler(1000));
double area = UbMath::PI*radius*radius;
SPtr<UbScheduler> forceSch(new UbScheduler(100));
SPtr<CalculateForcesCoProcessor> fp = make_shared<CalculateForcesCoProcessor>(grid, forceSch, outputPath + "/forces/forces.txt", comm, velocity, area);
fp->addInteractor(sphereInt);
......
apps/cpu/sphere/sphere.cpp
+ 88
− 102
View file @ ccb8cc57
......@@ -9,37 +9,36 @@ void run(string configname)
{
try
{
//Sleep(30000);
string machine = QUOTEME(CAB_MACHINE);
SPtr<Communicator> comm = MPICommunicator::getInstance();
int myid = comm->getProcessID();
int mybundle = comm->getBundleID();
int root = comm->getRoot();
ConfigurationFile config;
config.load(configname);
string pathname = config.getValue<string>("path");
double availMem = config.getValue<double>("memory");
string metafile = config.getValue<string>("metafile");
double outstep = config.getValue<double>("outstep");
double endstep = config.getValue<double>("endstep");
int numOfThreads = config.getValue<int>("threads");
const int refineLevel = config.getValue<int>("level");
bool test = config.getValue<bool>("test");
LBMReal radius = 10;
LBMReal uLB = 0.01;
//ConfigurationFile config;
//config.load(configname);
//string pathname = config.getValue<string>("path");
//double availMem = config.getValue<double>("memory");
//string metafile = config.getValue<string>("metafile");
//double outstep = config.getValue<double>("outstep");
//double endstep = config.getValue<double>("endstep");
//int numOfThreads = config.getValue<int>("threads");
//const int refineLevel = config.getValue<int>("level");
string outputPath = "d:/temp/sphereParab4";
double availMem = 8e9;
double outstep = 10000;
double endstep = 1e6;
int numOfThreads = 4;
omp_set_num_threads(numOfThreads);
int refineLevel = 0;
LBMReal radius = 5;
LBMReal uLB = 1e-3;
LBMReal Re = 1;
LBMReal rhoLB = 0.0;
//LBMReal nuLB = (uLB*2.0*radius)/Re;
//LBMReal nuLB = (uLB*L2)/Re;
LBMReal nuLB = 0.168666666667/100;
LBMReal nuLB = (uLB*2.0*radius)/Re;
double dp_LB = 1e-6;
double rhoLBinflow = dp_LB*3.0;
......@@ -47,11 +46,17 @@ void run(string configname)
SPtr<BCAdapter> noSlipBCAdapter(new NoSlipBCAdapter());
noSlipBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NoSlipBCAlgorithm()));
double H = 50;
//mu::Parser fct;
//fct.SetExpr("U");
//fct.DefineConst("U", uLB);
mu::Parser fct;
fct.SetExpr("U");
fct.SetExpr("16*U*x2*x3*(H-x2)*(H-x3)/H^4");
fct.DefineConst("U", uLB);
fct.DefineConst("H", H);
SPtr<BCAdapter> velBCAdapter(new VelocityBCAdapter(true, false, false, fct, 0, BCFunction::INFCONST));
velBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityBCAlgorithm()));
//velBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityBCAlgorithm()));
velBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new SimpleVelocityBCAlgorithm()));
SPtr<BCAdapter> denBCAdapter(new DensityBCAdapter(rhoLB));
denBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NonEqDensityBCAlgorithm()));
......@@ -61,39 +66,34 @@ void run(string configname)
bcVisitor.addBC(velBCAdapter);
bcVisitor.addBC(denBCAdapter);
SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter());
double dx = 1;
const int blocknx1 = 5;
const int blocknx2 = 5;
const int blocknx3 = 5;
const int gridNx1 = 8;//18;
const int gridNx2 = 8;// 11;
const int gridNx3 = 8;// 11;
//const int blocknx1 = 40;
//const int blocknx2 = 40;
//const int blocknx3 = 40;
//const int gridNx1 = 2;
//const int gridNx2 = 2;
//const int gridNx3 = 2;
const int blocknx1 = 50;
const int blocknx2 = 50;
const int blocknx3 = 50;
double L1 = gridNx1*blocknx1;
double L2, L3;
L2 = gridNx2*blocknx1;
L3 = gridNx3*blocknx1;
const int gridNx1 = 150;
const int gridNx2 = H;
const int gridNx3 = H;
double L1, L2, L3;
L1 = gridNx1;
L2 = gridNx2;
L3 = gridNx3;
SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter());
SPtr<Grid3D> grid(new Grid3D(comm));
grid->setDeltaX(dx);
grid->setBlockNX(blocknx1, blocknx2, blocknx3);
//sphere
//SPtr<GbObject3D> sphere(new GbSphere3D(L1 * 0.5, L2 * 0.5, L3 * 0.5, radius));
SPtr<GbObject3D> sphere(new GbSphere3D(75, 25, 25, radius));
GbSystem3D::writeGeoObject(sphere.get(), outputPath + "/geo/sphere", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> sphereInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(sphere, grid, noSlipBCAdapter, Interactor3D::SOLID));
if (grid->getTimeStep() == 0)
if (true)
{
const int baseLevel = 0;
......@@ -123,64 +123,48 @@ void run(string configname)
}
SPtr<GbObject3D> gridCube(new GbCuboid3D(d_minX1, d_minX2, d_minX3, d_maxX1, d_maxX2, d_maxX3));
if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), outputPath + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance());
GenBlocksGridVisitor genBlocks(gridCube);
grid->accept(genBlocks);
//sphere
//SPtr<GbObject3D> sphereRef(new GbSphere3D(L1/4.0, L2*0.5, L3*0.5, radius+1.0));
//GbSystem3D::writeGeoObject(sphereRef.get(),pathname + "/geo/sphereRef", WbWriterVtkXmlBinary::getInstance());
//sphere
SPtr<GbObject3D> sphere(new GbSphere3D(L1*0.5, L2*0.5, L3*0.5, radius));
//SPtr<GbObject3D> sphere(new GbSphere3D(L1/2.0-4.0, L2*0.5+4.0, L3*0.5+4.0, radius));
//SPtr<GbObject3D> sphere(new GbCuboid3D(L1/4.0-radius, L2/2.0-radius, L3/2.0-radius, L1/4.0+radius, L2/2.0+radius, L3/2.0+radius));
GbSystem3D::writeGeoObject(sphere.get(), pathname + "/geo/sphere", WbWriterVtkXmlBinary::getInstance());
double off = 0.0;
SPtr<GbObject3D> refCube(new GbCuboid3D(sphere->getX1Minimum() - off, sphere->getX2Minimum() - off, sphere->getX3Minimum(),
sphere->getX1Maximum() + off, sphere->getX2Maximum() + off, sphere->getX3Maximum()));
if (myid == 0) GbSystem3D::writeGeoObject(refCube.get(), pathname + "/geo/refCube", WbWriterVtkXmlBinary::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(refCube.get(), outputPath + "/geo/refCube", WbWriterVtkXmlBinary::getInstance());
if (refineLevel > 0)
{
if (myid == 0) UBLOG(logINFO, "Refinement - start");
RefineCrossAndInsideGbObjectHelper refineHelper(grid, refineLevel, comm);
refineHelper.addGbObject(sphere, refineLevel);
//refineHelper.addGbObject(refCube, refineLevel);
//refineHelper.addGbObject(sphere, refineLevel);
refineHelper.addGbObject(refCube, refineLevel);
refineHelper.refine();
if (myid == 0) UBLOG(logINFO, "Refinement - end");
}
//walls
GbCuboid3DPtr addWallYmin(new GbCuboid3D(d_minX1 - 4.0*blockLength, d_minX2 - 4.0*blockLength, d_minX3 - 4.0*blockLength, d_maxX1 + 4.0*blockLength, d_minX2, d_maxX3 + 4.0*blockLength));
if (myid == 0) GbSystem3D::writeGeoObject(addWallYmin.get(), pathname + "/geo/addWallYmin", WbWriterVtkXmlASCII::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(addWallYmin.get(), outputPath + "/geo/addWallYmin", WbWriterVtkXmlASCII::getInstance());
GbCuboid3DPtr addWallZmin(new GbCuboid3D(d_minX1 - 4.0*blockLength, d_minX2 - 4.0*blockLength, d_minX3 - 4.0*blockLength, d_maxX1 + 4.0*blockLength, d_maxX2 + 4.0*blockLength, d_minX3));
if (myid == 0) GbSystem3D::writeGeoObject(addWallZmin.get(), pathname + "/geo/addWallZmin", WbWriterVtkXmlASCII::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(addWallZmin.get(), outputPath + "/geo/addWallZmin", WbWriterVtkXmlASCII::getInstance());
GbCuboid3DPtr addWallYmax(new GbCuboid3D(d_minX1 - 4.0*blockLength, d_maxX2, d_minX3 - 4.0*blockLength, d_maxX1 + 4.0*blockLength, d_maxX2 + 4.0*blockLength, d_maxX3 + 4.0*blockLength));
if (myid == 0) GbSystem3D::writeGeoObject(addWallYmax.get(), pathname + "/geo/addWallYmax", WbWriterVtkXmlASCII::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(addWallYmax.get(), outputPath + "/geo/addWallYmax", WbWriterVtkXmlASCII::getInstance());
GbCuboid3DPtr addWallZmax(new GbCuboid3D(d_minX1 - 4.0*blockLength, d_minX2 - 4.0*blockLength, d_maxX3, d_maxX1 + 4.0*blockLength, d_maxX2 + 4.0*blockLength, d_maxX3 + 4.0*blockLength));
if (myid == 0) GbSystem3D::writeGeoObject(addWallZmax.get(), pathname + "/geo/addWallZmax", WbWriterVtkXmlASCII::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(addWallZmax.get(), outputPath + "/geo/addWallZmax", WbWriterVtkXmlASCII::getInstance());
//inflow
GbCuboid3DPtr geoInflow(new GbCuboid3D(d_minX1 - 4.0*blockLength, d_minX2 - 4.0*blockLength, d_minX3 - 4.0*blockLength, d_minX1, d_maxX2 + 4.0*blockLength, d_maxX3 + 4.0*blockLength));
if (myid == 0) GbSystem3D::writeGeoObject(geoInflow.get(), pathname + "/geo/geoInflow", WbWriterVtkXmlASCII::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(geoInflow.get(), outputPath + "/geo/geoInflow", WbWriterVtkXmlASCII::getInstance());
//outflow
GbCuboid3DPtr geoOutflow(new GbCuboid3D(d_maxX1, d_minX2 - 4.0*blockLength, d_minX3 - 4.0*blockLength, d_maxX1 + 4.0*blockLength, d_maxX2 + 4.0*blockLength, d_maxX3 + 4.0*blockLength));
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow.get(), pathname + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance());
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow.get(), outputPath + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance());
SPtr<CoProcessor> ppblocks(new WriteBlocksCoProcessor(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, WbWriterVtkXmlBinary::getInstance(), comm));
//sphere
SPtr<D3Q27Interactor> sphereInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(sphere, grid, noSlipBCAdapter, Interactor3D::SOLID));
SPtr<CoProcessor> ppblocks(new WriteBlocksCoProcessor(grid, SPtr<UbScheduler>(new UbScheduler(1)), outputPath, WbWriterVtkXmlBinary::getInstance(), comm));
//walls
SPtr<D3Q27Interactor> addWallYminInt(new D3Q27Interactor(addWallYmin, grid, noSlipBCAdapter, Interactor3D::SOLID));
......@@ -204,7 +188,7 @@ void run(string configname)
SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, D3Q27System::B));
InteractorsHelper intHelper(grid, metisVisitor);
//intHelper.addInteractor(sphereInt);
intHelper.addInteractor(sphereInt);
intHelper.addInteractor(addWallYminInt);
intHelper.addInteractor(addWallZminInt);
intHelper.addInteractor(addWallYmaxInt);
......@@ -217,11 +201,6 @@ void run(string configname)
//PQueuePartitioningGridVisitor pqPartVisitor(numOfThreads);
//grid->accept(pqPartVisitor);
////set connectors
//InterpolationProcessorPtr iProcessor(new IncompressibleOffsetInterpolationProcessor());
//SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
//grid->accept(setConnsVisitor);
ppblocks->process(0);
ppblocks.reset();
......@@ -241,8 +220,8 @@ void run(string configname)
UBLOG(logINFO, "Available memory per process = " << availMem << " bytes");
}
//SPtr<LBMKernel> kernel(new IncompressibleCumulantLBMKernel());
SPtr<LBMKernel> kernel(new CompressibleCumulantLBMKernel());
SPtr<LBMKernel> kernel(new IncompressibleCumulantLBMKernel());
//SPtr<LBMKernel> kernel(new CompressibleCumulantLBMKernel());
SPtr<BCProcessor> bcProcessor(new BCProcessor());
......@@ -270,14 +249,14 @@ void run(string configname)
//initialization of distributions
InitDistributionsBlockVisitor initVisitor;
//initVisitor.setVx1(fct);
initVisitor.setVx1(fct);
//initVisitor.setRho(fctRoh);
grid->accept(initVisitor);
//Postrozess
SPtr<UbScheduler> geoSch(new UbScheduler(1));
SPtr<CoProcessor> ppgeo(
new WriteBoundaryConditionsCoProcessor(grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm));
new WriteBoundaryConditionsCoProcessor(grid, geoSch, outputPath, WbWriterVtkXmlBinary::getInstance(), comm));
ppgeo->process(0);
ppgeo.reset();
......@@ -286,28 +265,34 @@ void run(string configname)
else
{
UBLOG(logINFO, "SetConnectors - start, id=" << myid);
}
//set connectors
//SPtr<InterpolationProcessor> iProcessor(new IncompressibleOffsetInterpolationProcessor());
SPtr<CompressibleOffsetMomentsInterpolationProcessor> iProcessor(new CompressibleOffsetMomentsInterpolationProcessor());
SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
//SPtr<ConnectorFactory> cFactory(new Block3DConnectorFactory());
//ConnectorBlockVisitor setConnsVisitor(comm, nuLB, iProcessor, cFactory);
grid->accept(setConnsVisitor);
UBLOG(logINFO, "SetConnectors - end, id=" << myid);
}
UBLOG(logINFO, "SetConnectors - start, id=" << myid);
//set connectors
SPtr<InterpolationProcessor> iProcessor(new IncompressibleOffsetInterpolationProcessor());
//SPtr<CompressibleOffsetMomentsInterpolationProcessor> iProcessor(new CompressibleOffsetMomentsInterpolationProcessor());
SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
grid->accept(setConnsVisitor);
UBLOG(logINFO, "SetConnectors - end, id=" << myid);
SPtr<UbScheduler> stepSch(new UbScheduler(outstep));
//stepSch->addSchedule(10000, 0, 1000000);
SPtr<WriteMacroscopicQuantitiesCoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, stepSch, pathname, WbWriterVtkXmlBinary::getInstance(), conv,comm));
SPtr<WriteMacroscopicQuantitiesCoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, stepSch, outputPath, WbWriterVtkXmlBinary::getInstance(), conv,comm));
SPtr<UbScheduler> nupsSch(new UbScheduler(10, 30, 100));
SPtr<CoProcessor> npr(new NUPSCounterCoProcessor(grid, nupsSch, numOfThreads, comm));
double area = UbMath::PI * radius * radius;
SPtr<UbScheduler> forceSch(new UbScheduler(1));
SPtr<CalculateForcesCoProcessor> fp = make_shared<CalculateForcesCoProcessor>(grid, forceSch, outputPath + "/forces/forces.txt", comm, uLB, area);
fp->addInteractor(sphereInt);
SPtr<UbScheduler> stepGhostLayer(new UbScheduler(1));
SPtr<Calculator> calculator(new BasicCalculator(grid, stepGhostLayer, endstep));
calculator->addCoProcessor(npr);
calculator->addCoProcessor(fp);
calculator->addCoProcessor(writeMQCoProcessor);
......@@ -336,14 +321,15 @@ int main(int argc, char* argv[])
{
if (argv != NULL)
{
if (argv[1] != NULL)
{
run(string(argv[1]));
}
else
{
cout << "Configuration file is missing!" << endl;
}
//if (argv[1] != NULL)
//{
// run(string(argv[1]));
//}
//else
//{
// cout << "Configuration file is missing!" << endl;
//}
run("sphere");
}
}
src/cpu/VirtualFluidsCore/Data/DataSet3D.h
+ 33
− 0
View file @ ccb8cc57
//=======================================================================================
// ____ ____ __ ______ __________ __ __ __ __
// \ \ | | | | | _ \ |___ ___| | | | | / \ | |
// \ \ | | | | | |_) | | | | | | | / \ | |
// \ \ | | | | | _ / | | | | | | / /\ \ | |
// \ \ | | | | | | \ \ | | | \__/ | / ____ \ | |____
// \ \ | | |__| |__| \__\ |__| \________/ /__/ \__\ |_______|
// \ \ | | ________________________________________________________________
// \ \ | | | ______________________________________________________________|
// \ \| | | | __ __ __ __ ______ _______
// \ | | |_____ | | | | | | | | | _ \ / _____)
// \ | | _____| | | | | | | | | | | \ \ \_______
// \ | | | | |_____ | \_/ | | | | |_/ / _____ |
// \ _____| |__| |________| \_______/ |__| |______/ (_______/
//
// This file is part of VirtualFluids. VirtualFluids is free software: you can
// redistribute it and/or modify it under the terms of the GNU General Public
// License as published by the Free Software Foundation, either version 3 of
// the License, or (at your option) any later version.
//
// VirtualFluids 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 General Public License
// for more details.
//
// You should have received a copy of the GNU General Public License along
// with VirtualFluids (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
//
//! \file DataSet3D.h
//! \ingroup Data
//! \author Konstantin Kutscher
//=======================================================================================
#ifndef DataSet3D_h
#define DataSet3D_h
......
src/cpu/VirtualFluidsCore/LBM/Thixotropy.cpp
+ 4
− 4
View file @ ccb8cc57
......@@ -18,7 +18,7 @@ void Thixotropy::setYieldStress(LBMReal yieldStress)
{
tau0 = yieldStress;
}
LBMReal Thixotropy::getYieldStress()
LBMReal Thixotropy::getYieldStress() const
{
return tau0;
}
......@@ -26,7 +26,7 @@ void Thixotropy::setViscosityParameter(LBMReal kParameter)
{
k = kParameter;
}
LBMReal Thixotropy::getViscosityParameter()
LBMReal Thixotropy::getViscosityParameter() const
{
return k;
}
......@@ -34,7 +34,7 @@ void Thixotropy::setPowerIndex(LBMReal index)
{
n = index;
}
LBMReal Thixotropy::getPowerIndex()
LBMReal Thixotropy::getPowerIndex() const
{
return n;
}
......@@ -43,7 +43,7 @@ void Thixotropy::setOmegaMin(LBMReal omega)
{
omegaMin = omega;
}
LBMReal Thixotropy::getOmegaMin()
LBMReal Thixotropy::getOmegaMin() const
{
return omegaMin;
}
......
src/cpu/VirtualFluidsCore/LBM/Thixotropy.h
+ 4
− 5
View file @ ccb8cc57
......@@ -12,16 +12,16 @@ public:
Thixotropy& operator=(Thixotropy const&) = delete;
static SPtr<Thixotropy> getInstance();
void setYieldStress(LBMReal tau0);
LBMReal getYieldStress();
LBMReal getYieldStress() const;
void setViscosityParameter(LBMReal k);
LBMReal getViscosityParameter();
LBMReal getViscosityParameter() const;
void setPowerIndex(LBMReal n);
LBMReal getPowerIndex();
LBMReal getPowerIndex() const;
void setOmegaMin(LBMReal omegaMin);
LBMReal getOmegaMin();
LBMReal getOmegaMin() const;
static LBMReal getBinghamCollFactor(LBMReal omegaInf, LBMReal shearRate, LBMReal drho);
static LBMReal getHerschelBulkleyCollFactor(LBMReal omegaInf, LBMReal shearRate, LBMReal drho);
......@@ -53,7 +53,6 @@ inline LBMReal Thixotropy::getHerschelBulkleyCollFactor(LBMReal omegaInf, LBMRea
LBMReal gammaDot = shearRate;
LBMReal omega = omegaInf;
LBMReal epsilon = 1;
LBMReal omegaMin = 0.7;
while (epsilon > 1e-10)
{
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment