diff --git a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.C b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.C index 3882d0f1..cdc33354 100755 --- a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.C +++ b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.C @@ -70,12 +70,13 @@ IBVoidFraction::IBVoidFraction alphaLimited_(0), scaleUpVol_(readScalar(propsDict_.lookup("scaleUpVol"))) { - Info << "\n\n W A R N I N G - do not use in combination with differentialRegion model! \n\n" << endl; - //Info << "\n\n W A R N I N G - this model does not yet work properly! \n\n" << endl; - maxCellsPerParticle_=readLabel(propsDict_.lookup("maxCellsPerParticle")); + Info << "\n\n W A R N I N G - do not use in combination with differentialRegion model!\n\n" << endl; + maxCellsPerParticle_ = readLabel(propsDict_.lookup("maxCellsPerParticle")); - if(scaleUpVol_ < 1){ FatalError<< "scaleUpVol shloud be > 1."<< abort(FatalError); } - if(alphaMin_ > 1 || alphaMin_ < 0.01){ FatalError<< "alphaMin shloud be > 1 and < 0.01." << abort(FatalError); } + if (scaleUpVol_ < 1.0) + FatalError << "scaleUpVol shloud be > 1." << abort(FatalError); + if (alphaMin_ > 1.0 || alphaMin_ < 0.01) + FatalError << "alphaMin shloud be > 1 and < 0.01." << abort(FatalError); } @@ -93,25 +94,26 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction reAllocArrays(); - voidfractionNext_.ref()=1; + voidfractionNext_.ref() = 1.0; - for(int index=0; index< particleCloud_.numberOfParticles(); index++) + for (int index=0; index < particleCloud_.numberOfParticles(); index++) { //if(mask[index][0]) //{ //reset - for(int subcell=0;subcell= 0) { @@ -120,9 +122,9 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction //compute the voidfraction for the cell "particleCentreCellID vector cellCentrePosition = particleCloud_.mesh().C()[particleCenterCellID]; scalar fc = pointInParticle(index, positionCenter, cellCentrePosition); - vector minPeriodicParticlePos=positionCenter; + vector minPeriodicParticlePos = positionCenter; - if(particleCloud_.checkPeriodicCells()) //consider minimal distance to all periodic images of this particle + if (particleCloud_.checkPeriodicCells()) //consider minimal distance to all periodic images of this particle { fc = minPeriodicDistance(index,cellCentrePosition, positionCenter, globalBb, minPeriodicParticlePos); } @@ -133,7 +135,7 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction if (pointInParticle(index, minPeriodicParticlePos, coronaPoint) < 0.0) { - voidfractionNext_[particleCenterCellID] = 0; + voidfractionNext_[particleCenterCellID] = 0.0; } else { @@ -145,7 +147,7 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction vector vertexPosition = particleCloud_.mesh().points()[vertices[i]]; scalar fv = pointInParticle(index, positionCenter, vertexPosition); - if(particleCloud_.checkPeriodicCells()) //consider minimal distance to all periodic images of this particle + if (particleCloud_.checkPeriodicCells()) //consider minimal distance to all periodic images of this particle { fv = minPeriodicDistance(index, vertexPosition, positionCenter, globalBb, minPeriodicParticlePos); } @@ -173,12 +175,12 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction buildLabelHashSet(index, minPeriodicParticlePos, particleCenterCellID, hashSett, true); //Add cells of periodic particle images on same processor - if(particleCloud_.checkPeriodicCells()) + if (particleCloud_.checkPeriodicCells()) { int doPeriodicImage[3]; - for (int iDir=0; iDir<3; iDir++) + for (int iDir=0; iDir < 3; iDir++) { - doPeriodicImage[iDir]= 0; + doPeriodicImage[iDir] = 0; if ((minPeriodicParticlePos[iDir]+radius) > globalBb.max()[iDir]) { doPeriodicImage[iDir] = -1; @@ -209,10 +211,10 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction copyCounter++; } //y-direction - int currCopyCounter=copyCounter; + int currCopyCounter = copyCounter; if (doPeriodicImage[1] != 0) { - for(int yDirCop=0; yDirCop<=currCopyCounter; yDirCop++) + for (int yDirCop=0; yDirCop <= currCopyCounter; yDirCop++) { particlePosList.append( particlePosList[yDirCop] + vector( @@ -225,10 +227,10 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction } } //z-direction - currCopyCounter=copyCounter; + currCopyCounter = copyCounter; if (doPeriodicImage[2] != 0) { - for(int zDirCop=0; zDirCop<=currCopyCounter; zDirCop++) + for (int zDirCop=0; zDirCop <= currCopyCounter; zDirCop++) { particlePosList.append( particlePosList[zDirCop] + vector( @@ -244,7 +246,7 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction //add the nearest cell labels particleLabelList.append(particleCenterCellID); - for(int iPeriodicImage=1; iPeriodicImage<=copyCounter; iPeriodicImage++) + for (int iPeriodicImage=1; iPeriodicImage <= copyCounter; iPeriodicImage++) { label partCellId = particleCloud_.mesh().findNearestCell(particlePosList[iPeriodicImage]); particleLabelList.append(partCellId); @@ -258,17 +260,17 @@ void IBVoidFraction::setvoidFraction(double** const& mask,double**& voidfraction label hashSetLength = hashSett.size(); if (hashSetLength > maxCellsPerParticle_) { - FatalError<< "big particle algo found more cells ("<< hashSetLength - <<") than storage is prepared ("< 0) { cellsPerParticle_[index][0]=hashSetLength; hashSett.erase(particleCenterCellID); - for(label i=0;i= -eps && lambda_ <= 1.0+eps) lambda = lambda_; - } + } } if (lambda < 0.0) diff --git a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.H b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.H index 61eb32a2..e2072025 100755 --- a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.H +++ b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/IBVoidFraction/IBVoidFraction.H @@ -28,7 +28,7 @@ Description This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER). - void fraction model for the smooth representation of spheres with + void fraction model for the smooth representation of spheres with radius > cell length. contribution from Alice Hager @@ -64,7 +64,7 @@ private: const scalar alphaMin_; //NP min value of voidFraction - mutable bool alphaLimited_; + mutable bool alphaLimited_; const scalar scaleUpVol_; //NP scaling radius, keeping volume of particle diff --git a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.C b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.C index 6683147f..b5f13e8a 100644 --- a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.C +++ b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.C @@ -76,16 +76,20 @@ dividedVoidFraction::dividedVoidFraction { maxCellsPerParticle_ = numberOfMarkerPoints; - if(alphaMin_ > 1 || alphaMin_ < 0.01){ Warning << "alphaMin should be < 1 and > 0.01 !!!" << endl; } - if (propsDict_.found("interpolation")){ - interpolation_=true; + if (alphaMin_ > 1.0 || alphaMin_ < 0.01) + Warning << "alphaMin should be < 1 and > 0.01 !!!" << endl; + + if (propsDict_.found("interpolation")) + { + interpolation_ = true; Warning << "interpolation for dividedVoidFraction does not yet work correctly!" << endl; - Info << "Using interpolated voidfraction field - do not use this in combination with interpolation in drag model!"<< endl; + Info << "Using interpolated voidfraction field - do not use this in combination with interpolation in drag model!" << endl; } checkWeightNporosity(propsDict_); - if (propsDict_.found("verbose")) verbose_=true; + if (propsDict_.found("verbose")) + verbose_ = true; if (propsDict_.found("cfdemUseOnly")) { @@ -165,7 +169,7 @@ dividedVoidFraction::~dividedVoidFraction() void dividedVoidFraction::setvoidFraction(double** const& mask,double**& voidfractions,double**& particleWeights,double**& particleVolumes, double**& particleV) const { - if(cfdemUseOnly_) + if (cfdemUseOnly_) reAllocArrays(particleCloud_.numberOfParticles()); else reAllocArrays(); @@ -179,18 +183,18 @@ void dividedVoidFraction::setvoidFraction(double** const& mask,double**& voidfra scalar scaleRadius = pow(porosity(),1./3.); const boundBox& globalBb = particleCloud_.mesh().bounds(); - for(int index=0; index< particleCloud_.numberOfParticles(); index++) + for (int index=0; index < particleCloud_.numberOfParticles(); index++) { - if(!checkParticleType(index)) continue; //skip this particle if not correct type + if (!checkParticleType(index)) continue; //skip this particle if not correct type //if(mask[index][0]) //{ // reset - for(int subcell=0;subcell maxCellsPerParticle_ || cellsSet < 0) + if (cellsSet > maxCellsPerParticle_ || cellsSet < 0) { Info << "ERROR cellsSet =" << cellsSet << endl; } if (!procBoundaryCorrection_) { - // set source for particle center; source 1/nPts+weight of all subpoints that have not been found - scalar centreWeight = 1./nPoints*(nPoints-cellsSet); + // set source for particle center; source 1/nPts+weight of all subpoints that have not been found + scalar centreWeight = 1./nPoints*(nPoints-cellsSet); - // update voidfraction for each particle read - scalar newAlpha = voidfractionNext_[cellID]- volume*centreWeight/cellVol; - if(newAlpha > alphaMin_) voidfractionNext_[cellID] = newAlpha; - else - { - voidfractionNext_[cellID] = alphaMin_; - tooMuch_ += (alphaMin_-newAlpha) * cellVol; - } + // update voidfraction for each particle read + scalar newAlpha = voidfractionNext_[cellID]- volume*centreWeight/cellVol; + if (newAlpha > alphaMin_) + { + voidfractionNext_[cellID] = newAlpha; + } + else + { + voidfractionNext_[cellID] = alphaMin_; + tooMuch_ += (alphaMin_-newAlpha) * cellVol; + } - // store cellweight for each particle --- this should be done for subpoints as well!! - particleWeights[index][0] += centreWeight; + // store cellweight for each particle --- this should be done for subpoints as well!! + particleWeights[index][0] += centreWeight; - // store particleVolume for each particle - particleVolumes[index][0] += volume*centreWeight; - particleV[index][0] += volume*centreWeight; + // store particleVolume for each particle + particleVolumes[index][0] += volume*centreWeight; + particleV[index][0] += volume*centreWeight; } /*//OUTPUT @@ -286,7 +293,7 @@ void dividedVoidFraction::setvoidFraction(double** const& mask,double**& voidfra // bring voidfraction from Eulerian Field to particle array //interpolationCellPoint voidfractionInterpolator_(voidfractionNext_); //scalar voidfractionAtPos(0); - for(int index=0; index< particleCloud_.numberOfParticles(); index++) + for(int index=0; index < particleCloud_.numberOfParticles(); index++) { /*if(interpolation_) { @@ -318,14 +325,14 @@ void dividedVoidFraction::setvoidFraction(double** const& mask,double**& voidfra } else*/ { - for(int subcell=0;subcell= 0) + if (cellID >= 0) { voidfractions[index][subcell] = voidfractionNext_[cellID]; - } + } else { voidfractions[index][subcell] = -1.; diff --git a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.H b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.H index c23f12cb..9f4d91a4 100644 --- a/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.H +++ b/src/lagrangian/cfdemParticle/subModels/voidFractionModel/dividedVoidFraction/dividedVoidFraction.H @@ -67,7 +67,7 @@ private: const scalar alphaMin_; // min value of voidFraction - mutable bool alphaLimited_; + mutable bool alphaLimited_; mutable scalar tooMuch_; // particle volume which is lost due to voidFraction limitation