From ca4412674ce40e794b1867d3c5bf51e679e6178d Mon Sep 17 00:00:00 2001
From: Konstantin Kutscher <kutscher@irmb.tu-bs.de>
Date: Fri, 4 May 2018 18:52:45 +0200
Subject: [PATCH] grid sub distance (q) calculation in D3Q27Interactor is
corrected
---
source/Applications/DLR-F16-Solid/f16.cpp | 100 ++----------------
.../Interactors/D3Q27Interactor.cpp | 26 +++--
source/VirtualFluidsCore/LBM/D3Q27System.cpp | 18 ++--
source/VirtualFluidsCore/LBM/D3Q27System.h | 33 ++++++
4 files changed, 67 insertions(+), 110 deletions(-)
diff --git a/source/Applications/DLR-F16-Solid/f16.cpp b/source/Applications/DLR-F16-Solid/f16.cpp
index ce6c22202..06681a118 100644
--- a/source/Applications/DLR-F16-Solid/f16.cpp
+++ b/source/Applications/DLR-F16-Solid/f16.cpp
@@ -126,6 +126,9 @@ void run(string configname)
const int baseLevel = 0;
+ //SPtr<GbObject3D> mic6(new GbCuboid3D(0.3, 0.015, -0.46+4.25*deltaXcoarse, 0.3+deltaXcoarse, 0.015+deltaXcoarse, -0.46+5.25*deltaXcoarse));
+ //if (myid==0) GbSystem3D::writeGeoObject(mic6.get(), pathOut+"/geo/mic6", WbWriterVtkXmlBinary::getInstance());
+
////////////////////////////////////////////////////////////////////////
//Grid
//////////////////////////////////////////////////////////////////////////
@@ -707,10 +710,14 @@ void run(string configname)
if (myid==0) GbSystem3D::writeGeoObject(mic4->getBoundingBox().get(), pathOut+"/geo/mic4", WbWriterVtkXmlBinary::getInstance());
SPtr<TimeseriesCoProcessor> tsp4(new TimeseriesCoProcessor(grid, stepMV, mic4, pathOut+"/mic/mic4", comm));
- SPtr<IntegrateValuesHelper> mic5(new IntegrateValuesHelper(grid, comm,0.3+deltaXfine, 0.015, 0.000517+0.00037+7.0*deltaXfine, 0.3+2.0*deltaXfine, 0.015+deltaXfine, 0.000517+0.00037+8.0*deltaXfine));
+ SPtr<IntegrateValuesHelper> mic5(new IntegrateValuesHelper(grid, comm, 0.3+deltaXfine, 0.015, 0.000517+0.00037+7.0*deltaXfine, 0.3+2.0*deltaXfine, 0.015+deltaXfine, 0.000517+0.00037+8.0*deltaXfine));
if (myid==0) GbSystem3D::writeGeoObject(mic5->getBoundingBox().get(), pathOut+"/geo/mic5", WbWriterVtkXmlBinary::getInstance());
SPtr<TimeseriesCoProcessor> tsp5(new TimeseriesCoProcessor(grid, stepMV, mic5, pathOut+"/mic/mic5", comm));
+ SPtr<IntegrateValuesHelper> mic6(new IntegrateValuesHelper(grid, comm, 0.3, 0.015, -0.46+4.25*deltaXcoarse, 0.3+deltaXcoarse, 0.015+deltaXcoarse, -0.46+5.25*deltaXcoarse));
+ if (myid==0) GbSystem3D::writeGeoObject(mic6->getBoundingBox().get(), pathOut+"/geo/mic6", WbWriterVtkXmlBinary::getInstance());
+ SPtr<TimeseriesCoProcessor> tsp6(new TimeseriesCoProcessor(grid, stepMV, mic6, pathOut+"/mic/mic6", comm));
+
omp_set_num_threads(numOfThreads);
SPtr<UbScheduler> stepGhostLayer(new UbScheduler(100));
SPtr<Calculator> calculator(new BasicCalculator(grid, stepGhostLayer, endTime));
@@ -723,6 +730,7 @@ void run(string configname)
calculator->addCoProcessor(tsp3);
calculator->addCoProcessor(tsp4);
calculator->addCoProcessor(tsp5);
+ calculator->addCoProcessor(tsp6);
//calculator->addCoProcessor(tav);
@@ -751,96 +759,6 @@ void run(string configname)
}
-//void test_run()
-//{
-// try
-// {
-//
-// SPtr<Communicator> comm = MPICommunicator::getInstance();.
-// int myid = comm->getProcessID();
-//
-// if (myid==0)
-// {
-// UBLOG(logINFO, "PID = "<<myid<<" Point 1");
-// UBLOG(logINFO, "PID = "<<myid<<" Total Physical Memory (RAM): "<<Utilities::getTotalPhysMem());
-// UBLOG(logINFO, "PID = "<<myid<<" Physical Memory currently used: "<<Utilities::getPhysMemUsed());
-// UBLOG(logINFO, "PID = "<<myid<<" Physical Memory currently used by current process: "<<Utilities::getPhysMemUsedByMe()/1073741824.0<<" GB");
-// }
-//
-// double g_minX1 = 0;
-// double g_minX2 = 0;
-// double g_minX3 = 0;
-//
-// double g_maxX1 = 5;
-// double g_maxX2 = 5;
-// double g_maxX3 = 5;
-//
-// int blockNx[3] ={ 5, 5, 5 };
-//
-// string pathOut = "d:/temp/DLR-F16-Solid-test";
-//
-// double deltaX = 1;
-// double rhoLB = 0.0;
-// double uLB = 0.0866025;
-// double nuLB = 0.001; //4.33013e-06;
-// SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter());
-// ////////////////////////////////////////////////////////////////////////
-// //Grid
-// //////////////////////////////////////////////////////////////////////////
-// SPtr<Grid3D> grid(new Grid3D(comm));
-// grid->setDeltaX(deltaX);
-// grid->setBlockNX(blockNx[0], blockNx[1], blockNx[2]);
-// grid->setPeriodicX1(false);
-// grid->setPeriodicX2(false);
-// grid->setPeriodicX3(false);
-//
-// SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3));
-// if (myid==0) GbSystem3D::writeGeoObject(gridCube.get(), pathOut+"/geo/gridCube", WbWriterVtkXmlASCII::getInstance());
-// GenBlocksGridVisitor genBlocks(gridCube);
-// grid->accept(genBlocks);
-//
-// WriteBlocksCoProcessor ppblocks(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathOut, WbWriterVtkXmlBinary::getInstance(), comm);
-// ppblocks.process(0);
-//
-// SPtr<LBMKernel> kernel = SPtr<LBMKernel>(new CompressibleCumulant4thOrderViscosityLBMKernel());
-// kernel->setNX(std::array<int, 3>{ {blockNx[0], blockNx[1], blockNx[2]}});
-// SPtr<BCProcessor> bcProc;
-// bcProc = SPtr<BCProcessor>(new BCProcessor());
-// kernel->setBCProcessor(bcProc);
-//
-// SetKernelBlockVisitor kernelVisitor(kernel, nuLB, 1e9, 12);
-// grid->accept(kernelVisitor);
-//
-// //initialization of distributions
-// InitDistributionsBlockVisitor initVisitor1;
-// initVisitor1.setVx1(0.001);
-// grid->accept(initVisitor1);
-//
-// SPtr<UbScheduler> stepSch(new UbScheduler(1));
-// SPtr<WriteMacroscopicQuantitiesCoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, stepSch, pathOut, WbWriterVtkXmlBinary::getInstance(), conv, comm));
-//
-// //omp_set_num_threads(numOfThreads);
-// SPtr<Calculator> calculator(new BasicCalculator(grid, stepSch, 2));
-// calculator->addCoProcessor(writeMQCoProcessor);
-//
-//
-// if (myid==0) UBLOG(logINFO, "Simulation-start");
-// calculator->calculate();
-// if (myid==0) UBLOG(logINFO, "Simulation-end");
-// }
-// catch (std::exception& e)
-// {
-// cerr<<e.what()<<endl<<flush;
-// }
-// catch (std::string& s)
-// {
-// cerr<<s<<endl;
-// }
-// catch (...)
-// {
-// cerr<<"unknown exception"<<endl;
-// }
-//}
int main(int argc, char* argv[])
{
diff --git a/source/VirtualFluidsCore/Interactors/D3Q27Interactor.cpp b/source/VirtualFluidsCore/Interactors/D3Q27Interactor.cpp
index ea0f08e11..820324c72 100644
--- a/source/VirtualFluidsCore/Interactors/D3Q27Interactor.cpp
+++ b/source/VirtualFluidsCore/Interactors/D3Q27Interactor.cpp
@@ -59,7 +59,9 @@ D3Q27Interactor::~D3Q27Interactor()
//////////////////////////////////////////////////////////////////////////
void D3Q27Interactor::initRayVectors()
{
- int fdir; double c1oS2 = UbMath::one_over_sqrt2;
+ int fdir;
+ double c1oS2 = UbMath::one_over_sqrt2;
+ double c1oS3 = UbMath::one_over_sqrt3;
fdir = D3Q27System::E; rayX1[fdir] = 1.0; rayX2[fdir] = 0.0; rayX3[fdir] = 0.0;
fdir = D3Q27System::W; rayX1[fdir] = -1.0; rayX2[fdir] = 0.0; rayX3[fdir] = 0.0;
fdir = D3Q27System::N; rayX1[fdir] = 0.0; rayX2[fdir] = 1.0; rayX3[fdir] = 0.0;
@@ -79,14 +81,14 @@ void D3Q27Interactor::initRayVectors()
fdir = D3Q27System::BN; rayX1[fdir] = 0.0; rayX2[fdir] = c1oS2; rayX3[fdir] = -c1oS2;
fdir = D3Q27System::TS; rayX1[fdir] = 0.0; rayX2[fdir] =-c1oS2; rayX3[fdir] = c1oS2;
- fdir = D3Q27System::TNW; rayX1[fdir] = -c1oS2; rayX2[fdir] = c1oS2; rayX3[fdir] = c1oS2;;
- fdir = D3Q27System::TNE; rayX1[fdir] = c1oS2; rayX2[fdir] = c1oS2; rayX3[fdir] = c1oS2;
- fdir = D3Q27System::TSW; rayX1[fdir] = -c1oS2; rayX2[fdir] =-c1oS2; rayX3[fdir] = c1oS2;
- fdir = D3Q27System::TSE; rayX1[fdir] = c1oS2; rayX2[fdir] =-c1oS2; rayX3[fdir] = c1oS2;
- fdir = D3Q27System::BNW; rayX1[fdir] = -c1oS2; rayX2[fdir] = c1oS2; rayX3[fdir] = -c1oS2;
- fdir = D3Q27System::BNE; rayX1[fdir] = c1oS2; rayX2[fdir] = c1oS2; rayX3[fdir] = -c1oS2;
- fdir = D3Q27System::BSW; rayX1[fdir] = -c1oS2; rayX2[fdir] =-c1oS2; rayX3[fdir] = -c1oS2;
- fdir = D3Q27System::BSE; rayX1[fdir] = c1oS2; rayX2[fdir] =-c1oS2; rayX3[fdir] = -c1oS2;
+ fdir = D3Q27System::TNW; rayX1[fdir] = -c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = c1oS3;
+ fdir = D3Q27System::TNE; rayX1[fdir] = c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = c1oS3;
+ fdir = D3Q27System::TSW; rayX1[fdir] = -c1oS3; rayX2[fdir] =-c1oS3; rayX3[fdir] = c1oS3;
+ fdir = D3Q27System::TSE; rayX1[fdir] = c1oS3; rayX2[fdir] =-c1oS3; rayX3[fdir] = c1oS3;
+ fdir = D3Q27System::BNW; rayX1[fdir] = -c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = -c1oS3;
+ fdir = D3Q27System::BNE; rayX1[fdir] = c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = -c1oS3;
+ fdir = D3Q27System::BSW; rayX1[fdir] = -c1oS3; rayX2[fdir] =-c1oS3; rayX3[fdir] = -c1oS3;
+ fdir = D3Q27System::BSE; rayX1[fdir] = c1oS3; rayX2[fdir] =-c1oS3; rayX3[fdir] = -c1oS3;
}
//////////////////////////////////////////////////////////////////////////
void D3Q27Interactor::initInteractor(const double& timeStep)
@@ -259,7 +261,11 @@ bool D3Q27Interactor::setDifferencesToGbObject3D(const SPtr<Block3D> block/*,con
vector<double> distNeigh(D3Q27System::FENDDIR+1, UbMath::sqrt2*deltaX1);
distNeigh[D3Q27System::E] = distNeigh[D3Q27System::W] = distNeigh[D3Q27System::N] = deltaX1;
distNeigh[D3Q27System::S] = distNeigh[D3Q27System::T] = distNeigh[D3Q27System::B] = deltaX1;
-
+ distNeigh[D3Q27System::NE] = distNeigh[D3Q27System::NW] = distNeigh[D3Q27System::SW] = distNeigh[D3Q27System::SE] = UbMath::sqrt2*deltaX1;
+ distNeigh[D3Q27System::TE] = distNeigh[D3Q27System::TN] = distNeigh[D3Q27System::TW] = distNeigh[D3Q27System::TS] = UbMath::sqrt2*deltaX1;
+ distNeigh[D3Q27System::BE] = distNeigh[D3Q27System::BN] = distNeigh[D3Q27System::BW] = distNeigh[D3Q27System::BS] = UbMath::sqrt2*deltaX1;
+ distNeigh[D3Q27System::TNE] = distNeigh[D3Q27System::TNW] = distNeigh[D3Q27System::TSE] = distNeigh[D3Q27System::TSW] = UbMath::sqrt3*deltaX1;
+ distNeigh[D3Q27System::BNE] = distNeigh[D3Q27System::BNW] = distNeigh[D3Q27System::BSE] = distNeigh[D3Q27System::BSW] = UbMath::sqrt3*deltaX1;
double q;
bool pointOnBoundary = false;
diff --git a/source/VirtualFluidsCore/LBM/D3Q27System.cpp b/source/VirtualFluidsCore/LBM/D3Q27System.cpp
index ef96b2dc9..c56164ab3 100644
--- a/source/VirtualFluidsCore/LBM/D3Q27System.cpp
+++ b/source/VirtualFluidsCore/LBM/D3Q27System.cpp
@@ -47,9 +47,9 @@ namespace D3Q27System
INV_BSW };
- // The x,y,z component for each normalized direction
- const double cNorm[3][ENDDIR] = {
- {
+ // The x,y,z component for each normalized direction
+ const double cNorm[3][ENDDIR] = {
+ {
double(DX1[0]), double(DX1[1]),
double(DX1[2]), double(DX1[3]),
double(DX1[4]), double(DX1[5]),
@@ -62,8 +62,8 @@ namespace D3Q27System
double(DX1[18]) / std::sqrt(double(3)), double(DX1[19]) / std::sqrt(double(3)),
double(DX1[20]) / std::sqrt(double(3)), double(DX1[21]) / std::sqrt(double(3)),
double(DX1[22]) / std::sqrt(double(3)), double(DX1[23]) / std::sqrt(double(3)),
- double(DX1[24]) / std::sqrt(double(3)), double(DX1[25]) / std::sqrt(double(3))
- },{
+ double(DX1[24]) / std::sqrt(double(3)), double(DX1[25]) / std::sqrt(double(3))
+ },{
double(DX2[0]), double(DX2[1]),
double(DX2[2]), double(DX2[3]),
double(DX2[4]), double(DX2[5]),
@@ -76,8 +76,8 @@ namespace D3Q27System
double(DX2[18]) / std::sqrt(double(3)), double(DX2[19]) / std::sqrt(double(3)),
double(DX2[20]) / std::sqrt(double(3)), double(DX2[21]) / std::sqrt(double(3)),
double(DX2[22]) / std::sqrt(double(3)), double(DX2[23]) / std::sqrt(double(3)),
- double(DX2[24]) / std::sqrt(double(3)), double(DX2[25]) / std::sqrt(double(3))
- },{
+ double(DX2[24]) / std::sqrt(double(3)), double(DX2[25]) / std::sqrt(double(3))
+ },{
double(DX3[0]), double(DX3[1]),
double(DX3[2]), double(DX3[3]),
double(DX3[4]), double(DX3[5]),
@@ -90,8 +90,8 @@ namespace D3Q27System
double(DX3[18]) / std::sqrt(double(3)), double(DX3[19]) / std::sqrt(double(3)),
double(DX3[20]) / std::sqrt(double(3)), double(DX3[21]) / std::sqrt(double(3)),
double(DX3[22]) / std::sqrt(double(3)), double(DX3[23]) / std::sqrt(double(3)),
- double(DX3[24]) / std::sqrt(double(3)), double(DX3[25]) / std::sqrt(double(3))
- }
+ double(DX3[24]) / std::sqrt(double(3)), double(DX3[25]) / std::sqrt(double(3))
+ }
};
}
diff --git a/source/VirtualFluidsCore/LBM/D3Q27System.h b/source/VirtualFluidsCore/LBM/D3Q27System.h
index 7e0d8f4f5..6ae3bf70e 100644
--- a/source/VirtualFluidsCore/LBM/D3Q27System.h
+++ b/source/VirtualFluidsCore/LBM/D3Q27System.h
@@ -729,6 +729,39 @@ namespace D3Q27System
distNeigh[TNE] = distNeigh[TNW] = distNeigh[TSE] = distNeigh[TSW] = sqrt(deltaX1*deltaX1+deltaX2*deltaX2+deltaX3*deltaX3);
distNeigh[BNE] = distNeigh[BNW] = distNeigh[BSE] = distNeigh[BSW] = sqrt(deltaX1*deltaX1+deltaX2*deltaX2+deltaX3*deltaX3);
}
+//////////////////////////////////////////////////////////////////////////
+ static inline void initRayVectors(double* const& rayX1, double* const& rayX2, double* const& rayX3)
+ {
+ int fdir;
+ double c1oS2 = UbMath::one_over_sqrt2;
+ double c1oS3 = UbMath::one_over_sqrt3;
+ fdir = E; rayX1[fdir] = 1.0; rayX2[fdir] = 0.0; rayX3[fdir] = 0.0;
+ fdir = W; rayX1[fdir] = -1.0; rayX2[fdir] = 0.0; rayX3[fdir] = 0.0;
+ fdir = N; rayX1[fdir] = 0.0; rayX2[fdir] = 1.0; rayX3[fdir] = 0.0;
+ fdir = S; rayX1[fdir] = 0.0; rayX2[fdir] = -1.0; rayX3[fdir] = 0.0;
+ fdir = T; rayX1[fdir] = 0.0; rayX2[fdir] = 0.0; rayX3[fdir] = 1.0;
+ fdir = B; rayX1[fdir] = 0.0; rayX2[fdir] = 0.0; rayX3[fdir] = -1.0;
+ fdir = NE; rayX1[fdir] = c1oS2; rayX2[fdir] = c1oS2; rayX3[fdir] = 0.0;
+ fdir = SW; rayX1[fdir] = -c1oS2; rayX2[fdir] = -c1oS2; rayX3[fdir] = 0.0;
+ fdir = SE; rayX1[fdir] = c1oS2; rayX2[fdir] = -c1oS2; rayX3[fdir] = 0.0;
+ fdir = NW; rayX1[fdir] = -c1oS2; rayX2[fdir] = c1oS2; rayX3[fdir] = 0.0;
+ fdir = TE; rayX1[fdir] = c1oS2; rayX2[fdir] = 0.0; rayX3[fdir] = c1oS2;
+ fdir = BW; rayX1[fdir] = -c1oS2; rayX2[fdir] = 0.0; rayX3[fdir] = -c1oS2;
+ fdir = BE; rayX1[fdir] = c1oS2; rayX2[fdir] = 0.0; rayX3[fdir] = -c1oS2;
+ fdir = TW; rayX1[fdir] = -c1oS2; rayX2[fdir] = 0.0; rayX3[fdir] = c1oS2;
+ fdir = TN; rayX1[fdir] = 0.0; rayX2[fdir] = c1oS2; rayX3[fdir] = c1oS2;
+ fdir = BS; rayX1[fdir] = 0.0; rayX2[fdir] =-c1oS2; rayX3[fdir] = -c1oS2;
+ fdir = BN; rayX1[fdir] = 0.0; rayX2[fdir] = c1oS2; rayX3[fdir] = -c1oS2;
+ fdir = TS; rayX1[fdir] = 0.0; rayX2[fdir] =-c1oS2; rayX3[fdir] = c1oS2;
+ fdir = TNE; rayX1[fdir] = c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = c1oS3;
+ fdir = TNW; rayX1[fdir] = -c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = c1oS3;
+ fdir = TSE; rayX1[fdir] = c1oS3; rayX2[fdir] = -c1oS3; rayX3[fdir] = c1oS3;
+ fdir = TSW; rayX1[fdir] = -c1oS3; rayX2[fdir] = -c1oS3; rayX3[fdir] = c1oS3;
+ fdir = BNE; rayX1[fdir] = c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = -c1oS3;
+ fdir = BNW; rayX1[fdir] = -c1oS3; rayX2[fdir] = c1oS3; rayX3[fdir] = -c1oS3;
+ fdir = BSE; rayX1[fdir] = c1oS3; rayX2[fdir] = -c1oS3; rayX3[fdir] = -c1oS3;
+ fdir = BSW; rayX1[fdir] = -c1oS3; rayX2[fdir] = -c1oS3; rayX3[fdir] = -c1oS3;
+ }
//////////////////////////////////////////////////////////////////////////
static inline LBMReal calcPress(const LBMReal* const f, LBMReal rho, LBMReal vx1, LBMReal vx2, LBMReal vx3)
{
--
GitLab