ENH: Adding volumeFraction to allow PairCollision models to take

account of nParticles to expand the effective parcel size for
collision.

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairCollision.C

@@ -213,7 +213,7 @@ void Foam::PairCollision<CloudType>::wallInteraction()
 
             const point& pos = p.position();
 
-            scalar r = p.d()/2;
+            scalar r = wallModel_->pREff(p);
 
             // real wallFace interactions
 

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairModel/PairModel.H

@@ -54,7 +54,7 @@ class PairModel
 {
     // Private data
 
-        //- The cloud dictionary
+        //- The CollisionModel dictionary
         const dictionary& dict_;
 
         //- Reference to the owner cloud class

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.C

@@ -45,13 +45,15 @@ void Foam::PairSpringSliderDashpot<CloudType>::findMinMaxProperties
 
         // Finding minimum diameter to avoid excessive arithmetic
 
-        RMin = min(p.d(), RMin);
+        scalar dEff = p.d()*cbrt(p.nParticle()*volumeFactor_);
+
+        RMin = min(dEff, RMin);
 
         rhoMax = max(p.rho(), rhoMax);
 
         UMagMax = max
         (
-            mag(p.U()) + mag(p.omega())*p.d()/2,
+            mag(p.U()) + mag(p.omega())*dEff/2,
             UMagMax
         );
     }
@@ -91,7 +93,8 @@ Foam::PairSpringSliderDashpot<CloudType>::PairSpringSliderDashpot
         (
             this->coeffDict().lookup("collisionResolutionSteps")
         )
-    )
+    ),
+    volumeFactor_(this->dict().lookupOrDefault("volumeFactor", 1.0))
 {
     scalar nu = this->owner().constProps().poissonsRatio();
 
@@ -155,7 +158,11 @@ void Foam::PairSpringSliderDashpot<CloudType>::evaluatePair
 {
     vector r_AB = (pA.position() - pB.position());
 
-    scalar normalOverlapMag = 0.5*(pA.d() + pB.d()) - mag(r_AB);
+    scalar dAEff = pA.d()*cbrt(pA.nParticle()*volumeFactor_);
+
+    scalar dBEff = pB.d()*cbrt(pB.nParticle()*volumeFactor_);
+
+    scalar normalOverlapMag = 0.5*(dAEff + dBEff) - mag(r_AB);
 
     if (normalOverlapMag > 0)
     {
@@ -166,7 +173,7 @@ void Foam::PairSpringSliderDashpot<CloudType>::evaluatePair
         vector U_AB = pA.U() - pB.U();
 
         // Effective radius
-        scalar R = 0.5*pA.d()*pB.d()/(pA.d() + pB.d());
+        scalar R = 0.5*dAEff*dBEff/(dAEff + dBEff);
 
         // Effective mass
         scalar M = pA.mass()*pB.mass()/(pA.mass() + pB.mass());
@@ -185,8 +192,8 @@ void Foam::PairSpringSliderDashpot<CloudType>::evaluatePair
 
         vector USlip_AB =
             U_AB - (U_AB & rHat_AB)*rHat_AB
-          + (pA.omega() ^ (pA.d()/2*-rHat_AB))
+          + (pA.omega() ^ (dAEff/2*-rHat_AB))
-          - (pB.omega() ^ (pB.d()/2*rHat_AB));
+          - (pB.omega() ^ (dBEff/2*rHat_AB));
 
         scalar deltaT = this->owner().mesh().time().deltaTValue();
 
@@ -241,8 +248,8 @@ void Foam::PairSpringSliderDashpot<CloudType>::evaluatePair
             pA.f() += fT_AB;
             pB.f() += -fT_AB;
 
-            pA.torque() += (pA.d()/2*-rHat_AB) ^ fT_AB;
+            pA.torque() += (dAEff/2*-rHat_AB) ^ fT_AB;
-            pB.torque() += (pB.d()/2*rHat_AB) ^ -fT_AB;
+            pB.torque() += (dBEff/2*rHat_AB) ^ -fT_AB;
         }
     }
 }

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H

@@ -74,6 +74,24 @@ class PairSpringSliderDashpot
         //  harmonic approximation of the collision period
         scalar collisionResolutionSteps_;
 
+        //- Volume factor for determining the equivalent size of a
+        //  parcel where nParticles is not 1.  The equivalent size of
+        //  the parcel is
+        //      parcelEquivVolume = volumeFactor*nParticles*p.volume()
+        //  so
+        //      parcelEquivD = cbrt(volumeFactor*nParticles)*p.d()
+        //  + When volumeFactor = 1, the particles are compressed
+        //    together so that the equivalent volume of the parcel is
+        //    the sum of the constituent particles
+        //  + When volumeFactor = 3*sqrt(2)/pi, the particles are
+        //    close packed, but uncompressed.
+        //  + When volumeFactor > 3*sqrt(2)/pi, the particles loosely
+        //    grouped.
+        // 3*sqrt(2)/pi = 1.350474 is the volume factor for close
+        // packing, i.e pi/(3*sqrt(2)) is the maximum close packing
+        // factor
+        scalar volumeFactor_;
+
 
     // Private Member Functions
 
@@ -104,6 +122,12 @@ public:
 
     // Member Functions
 
+        //- Return the volumeFactor
+        inline scalar volumeFactor() const
+        {
+            return volumeFactor_;
+        }
+
         //- Whether the PairModel has a timestep limit that will
         //  require subCycling
         virtual bool controlsTimestep() const;

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallModel/WallModel.H

@@ -54,7 +54,7 @@ class WallModel
 {
     // Private data
 
-        //- The cloud dictionary
+        //- The CollisionModel dictionary
         const dictionary& dict_;
 
         //- Reference to the owner cloud class
@@ -120,6 +120,9 @@ public:
 
     // Member Functions
 
+        //- Return the effective radius for a particle for the model
+        virtual scalar pREff(const typename CloudType::parcelType& p) const = 0;
+
         //- Whether the WallModel has a timestep limit that will
         //  require subCycling
         virtual bool controlsTimestep() const = 0;

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.C

@@ -44,13 +44,16 @@ void Foam::WallSpringSliderDashpot<CloudType>::findMinMaxProperties
         const typename CloudType::parcelType& p = iter();
 
         // Finding minimum diameter to avoid excessive arithmetic
-        rMin = min(p.d(), rMin);
+
+        scalar dEff = p.d()*cbrt(p.nParticle()*volumeFactor_);
+
+        rMin = min(dEff, rMin);
 
         rhoMax = max(p.rho(), rhoMax);
 
         UMagMax = max
         (
-            mag(p.U()) + mag(p.omega())*p.d()/2,
+            mag(p.U()) + mag(p.omega())*dEff/2,
             UMagMax
         );
     }
@@ -69,18 +72,17 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
     const WallSiteData<vector>& data,
     scalar pNu,
     scalar pE,
+    scalar pREff,
     scalar Estar,
     scalar kN,
     scalar Gstar
 ) const
 {
-    scalar pR = p.d()/2;
-
     vector r_PW = p.position() - site;
 
     vector U_PW = p.U() - data.wallData();
 
-    scalar normalOverlapMag = pR - mag(r_PW);
+    scalar normalOverlapMag = pREff - mag(r_PW);
 
     vector rHat_PW = r_PW/(mag(r_PW) + VSMALL);
 
@@ -94,7 +96,7 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
 
     vector USlip_PW =
         U_PW - (U_PW & rHat_PW)*rHat_PW
-      + (p.omega() ^ (pR*-rHat_PW));
+      + (p.omega() ^ (pREff*-rHat_PW));
 
     scalar deltaT = this->owner().mesh().time().deltaTValue();
 
@@ -108,7 +110,7 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
 
     if (tangentialOverlapMag > VSMALL)
     {
-        scalar kT = 8.0*sqrt(pR*normalOverlapMag)*Gstar;
+        scalar kT = 8.0*sqrt(pREff*normalOverlapMag)*Gstar;
 
         scalar etaT = etaN;
 
@@ -134,7 +136,7 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
 
         p.f() += fT_PW;
 
-        p.torque() += (pR*-rHat_PW) ^ fT_PW;
+        p.torque() += (pREff*-rHat_PW) ^ fT_PW;
     }
 }
 
@@ -160,7 +162,8 @@ Foam::WallSpringSliderDashpot<CloudType>::WallSpringSliderDashpot
         (
             this->coeffDict().lookup("collisionResolutionSteps")
         )
-    )
+    ),
+    volumeFactor_(this->dict().lookupOrDefault("volumeFactor", 1.0))
 {}
 
 
@@ -173,6 +176,15 @@ Foam::WallSpringSliderDashpot<CloudType>::~WallSpringSliderDashpot()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class CloudType>
+Foam::scalar Foam::WallSpringSliderDashpot<CloudType>::pREff
+(
+    const typename CloudType::parcelType& p
+) const
+{
+    return p.d()/2*cbrt(p.nParticle()*volumeFactor_);
+}
+
 template<class CloudType>
 bool Foam::WallSpringSliderDashpot<CloudType>::controlsTimestep() const
 {
@@ -225,9 +237,11 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
 
     scalar pE = this->owner().constProps().youngsModulus();
 
+    scalar pREff = this->pREff(p);
+
     scalar Estar = 1/((1 - sqr(pNu))/pE + (1 - sqr(nu_))/E_);
 
-    scalar kN = (4.0/3.0)*sqrt(p.d()/2)*Estar;
+    scalar kN = (4.0/3.0)*sqrt(pREff)*Estar;
 
     scalar GStar = 1/(2*((2 + pNu - sqr(pNu))/pE + (2 + nu_ - sqr(nu_))/E_));
 
@@ -240,6 +254,7 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
             flatSiteData[siteI],
             pNu,
             pE,
+            pREff,
             Estar,
             kN,
             GStar
@@ -257,6 +272,7 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
             sharpSiteData[siteI],
             pNu,
             pE,
+            pREff,
             Estar,
             kN,
             GStar

src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.H

@@ -69,6 +69,24 @@ class WallSpringSliderDashpot
         //  harmonic approximation of the collision period
         scalar collisionResolutionSteps_;
 
+        //- Volume factor for determining the equivalent size of a
+        //  parcel where nParticles is not 1.  The equivalent size of
+        //  the parcel is
+        //      parcelEquivVolume = volumeFactor*nParticles*p.volume()
+        //  so
+        //      parcelEquivD = cbrt(volumeFactor*nParticles)*p.d()
+        //  + When volumeFactor = 1, the particles are compressed
+        //    together so that the equivalent volume of the parcel is
+        //    the sum of the constituent particles
+        //  + When volumeFactor = 3*sqrt(2)/pi, the particles are
+        //    close packed, but uncompressed.
+        //  + When volumeFactor > 3*sqrt(2)/pi, the particles loosely
+        //    grouped.
+        // 3*sqrt(2)/pi = 1.350474 is the volume factor for close
+        // packing, i.e pi/(3*sqrt(2)) is the maximum close packing
+        // factor
+        scalar volumeFactor_;
+
 
     // Private Member Functions
 
@@ -89,6 +107,7 @@ class WallSpringSliderDashpot
             const WallSiteData<vector>& data,
             scalar pNu,
             scalar pE,
+            scalar pREff,
             scalar Estar,
             scalar kN,
             scalar Gstar
@@ -113,6 +132,15 @@ public:
 
     // Member Functions
 
+        //- Return the volumeFactor
+        inline scalar volumeFactor() const
+        {
+            return volumeFactor_;
+        }
+
+        //- Return the effective radius for a particle for the model
+        virtual scalar pREff(const typename CloudType::parcelType& p) const;
+
         //- Whether the WallModel has a timestep limit that will
         //  require subCycling
         virtual bool controlsTimestep() const;

src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C

@@ -217,7 +217,7 @@ Foam::scalar Foam::InjectionModel<CloudType>::setNumberOfParticles
         }
         case pbFixed:
         {
-            nP = nParticlesFixed_;
+            nP = nParticleFixed_;
             break;
         }
         default:
@@ -290,7 +290,7 @@ Foam::InjectionModel<CloudType>::InjectionModel(CloudType& owner)
     nInjections_(0),
     parcelsAddedTotal_(0),
     parcelBasis_(pbNumber),
-    nParticlesFixed_(0.0),
+    nParticleFixed_(0.0),
     time0_(0.0),
     timeStep0_(0.0)
 {
@@ -316,7 +316,7 @@ Foam::InjectionModel<CloudType>::InjectionModel
     nInjections_(0),
     parcelsAddedTotal_(0),
     parcelBasis_(pbNumber),
-    nParticlesFixed_(0.0),
+    nParticleFixed_(0.0),
     time0_(owner.db().time().value()),
     timeStep0_(0.0)
 {
@@ -340,11 +340,11 @@ Foam::InjectionModel<CloudType>::InjectionModel
     {
         parcelBasis_ = pbFixed;
 
-        Info<< "    Choosing nParticles to be a fixed value, massTotal "
+        Info<< "    Choosing nParticle to be a fixed value, massTotal "
             << "variable now does not determine anything."
             << endl;
 
-        nParticlesFixed_ = readScalar(coeffDict_.lookup("nParticles"));
+        nParticleFixed_ = readScalar(coeffDict_.lookup("nParticle"));
     }
     else
     {

src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H

@@ -138,9 +138,9 @@ protected:
             //- Parcel basis enumeration
             parcelBasis parcelBasis_;
 
-            //- nParticles to assign to parcels when the 'fixed' basis
+            //- nParticle to assign to parcels when the 'fixed' basis
             //  is selected
-            scalar nParticlesFixed_;
+            scalar nParticleFixed_;
 
             //- Continuous phase time at start of injection time step [s]
             scalar time0_;