diff --git a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNoSlipBCAlgorithm.cpp b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNoSlipBCAlgorithm.cpp index 13a7e517773ff32cfb93525e36622cc60e56ef5e..51fc2d5abdfe7cfa5bafb7ae21571c684cd26b02 100644 --- a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNoSlipBCAlgorithm.cpp +++ b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNoSlipBCAlgorithm.cpp @@ -74,36 +74,18 @@ void MultiphaseNoSlipBCAlgorithm::applyBC() D3Q27System::calcDensity(h, phi); - //LBMReal collFactorM = phi*collFactorL + (1-phi)*collFactorG; - //LBMReal collFactorM = collFactorL + (collFactorL - collFactorG)*(phi - phiH)/(phiH - phiL); - - //rho = phi + (1.0 - phi)*1.0/densityRatio; - //LBMReal rhoH = 1.0; - //LBMReal rhoL = 1.0/densityRatio; - //rho = rhoH + (rhoH - rhoL)*(phi - phiH)/(phiH - phiL); - calcMacrosFct(f, p1, vx1, vx2, vx3); - /*vx1/=(rho*c1o3); - vx2/=(rho*c1o3); - vx3/=(rho*c1o3);*/ - - //calcFeqFct(feq, rho, vx1, vx2, vx3); - //D3Q27System::calcMultiphaseFeq(feq, rho, p1, vx1, vx2, vx3); D3Q27System::calcMultiphaseFeqVB(feq, p1, vx1, vx2, vx3); D3Q27System::calcMultiphaseHeq(heq, phi, vx1, vx2, vx3); - //LBMReal collFactorM1 = 0.9; + for (int fdir = D3Q27System::FSTARTDIR; fdir<=D3Q27System::FENDDIR; fdir++) { if (bcPtr->hasNoSlipBoundaryFlag(fdir)) { //quadratic bounce back const int invDir = D3Q27System::INVDIR[fdir]; - //LBMReal q = bcPtr->getQ(invDir); - //LBMReal fReturn = ((1.0-q)/(1.0+q))*((f[invDir]-feq[invDir])/(1.0-collFactorM)+feq[invDir])+((q/(1.0+q))*(f[invDir]+f[fdir])); LBMReal fReturn = f[invDir]; distributions->setDistributionForDirection(fReturn, x1+D3Q27System::DX1[invDir], x2+D3Q27System::DX2[invDir], x3+D3Q27System::DX3[invDir], fdir); - - //LBMReal hReturn = ((1.0-q)/(1.0+q))*((h[invDir]-heq[invDir])/(1.0-collFactorPh)+heq[invDir])+((q/(1.0+q))*(h[invDir]+h[fdir])); LBMReal hReturn = h[invDir]; distributionsH->setDistributionForDirection(hReturn, x1+D3Q27System::DX1[invDir], x2+D3Q27System::DX2[invDir], x3+D3Q27System::DX3[invDir], fdir); } diff --git a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNonReflectingOutflowBCAlgorithm.cpp b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNonReflectingOutflowBCAlgorithm.cpp index 9ba171430686ba4f1873cd8be63fa3e06aa7aefb..d0a3b39c83c1d828801dd442fb387c4d4d33d589 100644 --- a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNonReflectingOutflowBCAlgorithm.cpp +++ b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseNonReflectingOutflowBCAlgorithm.cpp @@ -69,7 +69,6 @@ void MultiphaseNonReflectingOutflowBCAlgorithm::applyBC() LBMReal f[ENDF+1]; LBMReal ftemp[ENDF+1]; LBMReal h[D3Q27System::ENDF+1]; - //LBMReal heq[D3Q27System::ENDF+1]; LBMReal htemp[ENDF+1]; int nx1 = x1; @@ -91,25 +90,11 @@ void MultiphaseNonReflectingOutflowBCAlgorithm::applyBC() distributionsH->getDistribution(h, x1, x2, x3); distributionsH->getDistribution(htemp, nx1, nx2, nx3); - LBMReal phi, rho, p1, vx1, vx2, vx3; + LBMReal phi, p1, vx1, vx2, vx3; D3Q27System::calcDensity(h, phi); - //LBMReal collFactorM = phi*collFactorL + (1-phi)*collFactorG; - //LBMReal collFactorM = collFactorL + (collFactorL - collFactorG)*(phi - phiH)/(phiH - phiL); - - //rho = phi + (1.0 - phi)*1.0/densityRatio; - //LBMReal rhoH = 1.0; - //LBMReal rhoL = 1.0/densityRatio; - //rho = rhoH + (rhoH - rhoL)*(phi - phiH)/(phiH - phiL); - - - calcMacrosFct(f, p1, vx1, vx2, vx3); - /*vx1/=(rho*c1o3); - vx2/=(rho*c1o3); - vx3/=(rho*c1o3);*/ - switch (direction) { diff --git a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAdapter.cpp b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAdapter.cpp index b6692e3efdd5a22f5aa82d9d8c3b16d88854d6fd..7211bc3725b4a2607dc000c739f4bf8e98865013 100644 --- a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAdapter.cpp +++ b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAdapter.cpp @@ -236,24 +236,25 @@ void MultiphaseVelocityBCAdapter::init(const D3Q27Interactor* const& interactor, if( !this->isTimePeriodic() ) this->unsetTimeDependent(); else //bei peridoic die interavalle neu setzen: { - if( UbMath::equal(maxEndtime,BCFunction::INFCONST) ) - for(size_t pos=0; pos<vx1BCs.size(); ++pos) - { - vx1BCs[pos].setStartTime( vx1BCs[pos].getStartTime() + timeStep ); - vx1BCs[pos].setEndTime( vx1BCs[pos].getEndTime() + timeStep ); - } - if( UbMath::equal(maxEndtime,BCFunction::INFCONST) ) - for(size_t pos=0; pos<vx2BCs.size(); ++pos) - { - vx2BCs[pos].setStartTime( vx2BCs[pos].getStartTime() + timeStep ); - vx2BCs[pos].setEndTime( vx2BCs[pos].getEndTime() + timeStep ); - } - if( UbMath::equal(maxEndtime,BCFunction::INFCONST) ) - for(size_t pos=0; pos<vx3BCs.size(); ++pos) - { - vx3BCs[pos].setStartTime( vx3BCs[pos].getStartTime() + timeStep ); - vx3BCs[pos].setEndTime( vx3BCs[pos].getEndTime() + timeStep ); - } + if( UbMath::equal(maxEndtime,BCFunction::INFCONST) ) { + for (size_t pos = 0; pos < vx1BCs.size(); ++pos) { + vx1BCs[pos].setStartTime(vx1BCs[pos].getStartTime() + timeStep); + vx1BCs[pos].setEndTime(vx1BCs[pos].getEndTime() + timeStep); + } + } + if( UbMath::equal(maxEndtime,BCFunction::INFCONST) ) { + for (size_t pos = 0; pos < vx2BCs.size(); ++pos) { + vx2BCs[pos].setStartTime(vx2BCs[pos].getStartTime() + timeStep); + vx2BCs[pos].setEndTime(vx2BCs[pos].getEndTime() + timeStep); + } + } + + if( UbMath::equal(maxEndtime,BCFunction::INFCONST) ) { + for (size_t pos = 0; pos < vx3BCs.size(); ++pos) { + vx3BCs[pos].setStartTime(vx3BCs[pos].getStartTime() + timeStep); + vx3BCs[pos].setEndTime(vx3BCs[pos].getEndTime() + timeStep); + } + } this->init(interactor,time); } } diff --git a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAlgorithm.cpp b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAlgorithm.cpp index 5df0e09b619f2b0e7cd84cf8900751c0ebb71a23..80c1bc6518e5cfe27a4c1af9ca6bd2233d1b89d6 100644 --- a/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAlgorithm.cpp +++ b/src/cpu/VirtualFluidsCore/BoundaryConditions/MultiphaseVelocityBCAlgorithm.cpp @@ -71,27 +71,11 @@ void MultiphaseVelocityBCAlgorithm::applyBC() distributions->getDistributionInv(f, x1, x2, x3); distributionsH->getDistributionInv(h, x1, x2, x3); - LBMReal phi, rho, vx1, vx2, vx3, p1, phiBC; + LBMReal phi, vx1, vx2, vx3, p1, phiBC; D3Q27System::calcDensity(h, phi); - - //LBMReal collFactorM = phi*collFactorL + (1-phi)*collFactorG; - //LBMReal collFactorM = collFactorL + (collFactorL - collFactorG)*(phi - phiH)/(phiH - phiL); - - - - //rho = phi + (1.0 - phi)*1.0/densityRatio; - LBMReal rhoH = 1.0; - LBMReal rhoL = 1.0/densityRatio; - rho = rhoH + (rhoH - rhoL)*(phi - phiH)/(phiH - phiL); - calcMacrosFct(f, p1, vx1, vx2, vx3); - /*vx1/=(rho*c1o3); - vx2/=(rho*c1o3); - vx3/=(rho*c1o3);*/ - - //D3Q27System::calcMultiphaseFeq(feq, rho, p1, vx1, vx2, vx3); D3Q27System::calcMultiphaseFeqVB(feq, p1, vx1, vx2, vx3); D3Q27System::calcMultiphaseHeq(heq, phi, vx1, vx2, vx3); @@ -100,14 +84,14 @@ void MultiphaseVelocityBCAlgorithm::applyBC() int nx1 = x1; int nx2 = x2; int nx3 = x3; - int direction = -1; + //flag points in direction of fluid - if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::E)) { nx1 -= 1; direction = D3Q27System::E; } - else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::W)) { nx1 += 1; direction = D3Q27System::W; } - else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::N)) { nx2 -= 1; direction = D3Q27System::N; } - else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::S)) { nx2 += 1; direction = D3Q27System::S; } - else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::T)) { nx3 -= 1; direction = D3Q27System::T; } - else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::B)) { nx3 += 1; direction = D3Q27System::B; } + if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::E)) { nx1 -= 1; } + else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::W)) { nx1 += 1; } + else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::N)) { nx2 -= 1; } + else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::S)) { nx2 += 1; } + else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::T)) { nx3 -= 1; } + else if (bcPtr->hasVelocityBoundaryFlag(D3Q27System::B)) { nx3 += 1; } else UB_THROW(UbException(UB_EXARGS, "Danger...no orthogonal BC-Flag on velocity boundary...")); phiBC = bcPtr->getBoundaryPhaseField(); @@ -122,11 +106,6 @@ void MultiphaseVelocityBCAlgorithm::applyBC() distributionsH->setDistributionForDirection(hReturn, nx1, nx2, nx3, fdir); } } - - ////////////////////////////////// - - - for (int fdir = D3Q27System::FSTARTDIR; fdir<=D3Q27System::FENDDIR; fdir++) { @@ -135,12 +114,8 @@ void MultiphaseVelocityBCAlgorithm::applyBC() const int invDir = D3Q27System::INVDIR[fdir]; LBMReal q = bcPtr->getQ(invDir);// m+m q=0 stabiler LBMReal velocity = bcPtr->getBoundaryVelocity(invDir); - //LBMReal fReturn = ((1.0-q)/(1.0+q))*((f[invDir]-feq[invDir])/(1.0-collFactor)+feq[invDir])+((q*(f[invDir]+f[fdir])-velocity*rho*c1o3)/(1.0+q)); LBMReal fReturn = ((1.0-q)/(1.0+q))*((f[invDir]-feq[invDir])/(1.0-collFactor)+feq[invDir])+((q*(f[invDir]+f[fdir])-velocity)/(1.0+q)); distributions->setDistributionForDirection(fReturn, x1+D3Q27System::DX1[invDir], x2+D3Q27System::DX2[invDir], x3+D3Q27System::DX3[invDir], fdir); - - //LBMReal hReturn = ((1.0-q)/(1.0+q))*((h[invDir]-heq[invDir])/(1.0-collFactorM)+heq[invDir])+((q/(1.0+q))*(h[invDir]+h[fdir])); - //distributionsH->setDistributionForDirection(hReturn, x1+D3Q27System::DX1[invDir], x2+D3Q27System::DX2[invDir], x3+D3Q27System::DX3[invDir], fdir); } }