Skip to content

GitLab

Commit e14769f2 authored by Henry Weller's avatar Henry Weller
Browse files

lagrangian::ThermoSurfaceFilm: Updated adhesion->splash transition Weber number 

according to

        Bai et al, `Modelling of gasoline spray impingement', Atom. Sprays,
        vol 12, pp 1-27, 2002

Resolves bug-report https://bugs.openfoam.org/view.php?id=2478
parent 1fc82434
Hide whitespace changes
Inline Side-by-side
src/lagrangian/intermediate/submodels/Thermodynamic/SurfaceFilmModel/ThermoSurfaceFilm/ThermoSurfaceFilm.C
View file @ e14769f2
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -111,13 +111,13 @@ Foam::vector Foam::ThermoSurfaceFilm<CloudType>::splashDirection
const vector& nf
) const
{
// azimuthal angle [rad]
// Azimuthal angle [rad]
const scalar phiSi = twoPi*rndGen_.sample01<scalar>();
// ejection angle [rad]
// Ejection angle [rad]
const scalar thetaSi = pi/180.0*(rndGen_.sample01<scalar>()*(50 - 5) + 5);
// direction vector of new parcel
// Direction vector of new parcel
const scalar alpha = sin(thetaSi);
const scalar dcorr = cos(thetaSi);
const vector normal = alpha*(tanVec1*cos(phiSi) + tanVec2*sin(phiSi));
......@@ -226,7 +226,7 @@ void Foam::ThermoSurfaceFilm<CloudType>::drySplashInteraction
const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
const vector& nf = pp.faceNormals()[facei];
// local pressure
// Local pressure
const scalar pc = thermo_.thermo().p()[p.cell()];
// Retrieve parcel properties
......@@ -247,13 +247,13 @@ void Foam::ThermoSurfaceFilm<CloudType>::drySplashInteraction
// Critical Weber number
const scalar Wec = Adry_*pow(La, -0.183);
if (We < Wec) // adhesion - assume absorb
if (We < Wec) // Adhesion - assume absorb
{
absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
}
else // splash
else // Splash
{
// ratio of incident mass to splashing mass
// Ratio of incident mass to splashing mass
const scalar mRatio = 0.2 + 0.6*rndGen_.sample01<scalar>();
splashInteraction
(filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
......@@ -282,7 +282,7 @@ void Foam::ThermoSurfaceFilm<CloudType>::wetSplashInteraction
const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
const vector& nf = pp.faceNormals()[facei];
// local pressure
// Local pressure
const scalar pc = thermo_.thermo().p()[p.cell()];
// Retrieve parcel properties
......@@ -305,31 +305,31 @@ void Foam::ThermoSurfaceFilm<CloudType>::wetSplashInteraction
// Critical Weber number
const scalar Wec = Awet_*pow(La, -0.183);
if (We < 1) // adhesion - assume absorb
if (We < 2) // Adhesion - assume absorb
{
absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
}
else if ((We >= 1) && (We < 20)) // bounce
else if ((We >= 2) && (We < 20)) // Bounce
{
// incident angle of impingement
// Incident angle of impingement
const scalar theta = pi/2 - acos(U/mag(U) & nf);
// restitution coefficient
// Restitution coefficient
const scalar epsilon = 0.993 - theta*(1.76 - theta*(1.56 - theta*0.49));
// update parcel velocity
// Update parcel velocity
U = -epsilon*(Un) + 5/7*(Ut);
keepParticle = true;
return;
}
else if ((We >= 20) && (We < Wec)) // spread - assume absorb
else if ((We >= 20) && (We < Wec)) // Spread - assume absorb
{
absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
}
else // splash
else // Splash
{
// ratio of incident mass to splashing mass
// Ratio of incident mass to splashing mass
// splash mass can be > incident mass due to film entrainment
const scalar mRatio = 0.2 + 0.9*rndGen_.sample01<scalar>();
splashInteraction
......@@ -371,24 +371,24 @@ void Foam::ThermoSurfaceFilm<CloudType>::splashInteraction
const vector& posC = mesh.C()[p.cell()];
const vector& posCf = mesh.Cf().boundaryField()[pp.index()][facei];
// total mass of (all) splashed parcels
// Total mass of (all) splashed parcels
const scalar mSplash = m*mRatio;
// number of splashed particles per incoming particle
// Number of splashed particles per incoming particle
const scalar Ns = 5.0*(We/Wec - 1.0);
// average diameter of splashed particles
// Average diameter of splashed particles
const scalar dBarSplash = 1/cbrt(6.0)*cbrt(mRatio/Ns)*d + ROOTVSMALL;
// cumulative diameter splash distribution
// Cumulative diameter splash distribution
const scalar dMax = 0.9*cbrt(mRatio)*d;
const scalar dMin = 0.1*dMax;
const scalar K = exp(-dMin/dBarSplash) - exp(-dMax/dBarSplash);
// surface energy of secondary parcels [J]
// Surface energy of secondary parcels [J]
scalar ESigmaSec = 0;
// sample splash distribution to determine secondary parcel diameters
// Sample splash distribution to determine secondary parcel diameters
scalarList dNew(parcelsPerSplash_);
scalarList npNew(parcelsPerSplash_);
forAll(dNew, i)
......@@ -399,26 +399,26 @@ void Foam::ThermoSurfaceFilm<CloudType>::splashInteraction
ESigmaSec += npNew[i]*sigma*p.areaS(dNew[i]);
}
// incident kinetic energy [J]
// Incident kinetic energy [J]
const scalar EKIn = 0.5*m*magSqr(Urel);
// incident surface energy [J]
// Incident surface energy [J]
const scalar ESigmaIn = np*sigma*p.areaS(d);
// dissipative energy
// Dissipative energy
const scalar Ed = max(0.8*EKIn, np*Wec/12*pi*sigma*sqr(d));
// total energy [J]
// Total energy [J]
const scalar EKs = EKIn + ESigmaIn - ESigmaSec - Ed;
// switch to absorb if insufficient energy for splash
// Switch to absorb if insufficient energy for splash
if (EKs <= 0)
{
absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
return;
}
// helper variables to calculate magUns0
// Helper variables to calculate magUns0
const scalar logD = log(d);
const scalar coeff2 = log(dNew[0]) - logD + ROOTVSMALL;
scalar coeff1 = 0.0;
......@@ -427,7 +427,7 @@ void Foam::ThermoSurfaceFilm<CloudType>::splashInteraction
coeff1 += sqr(log(dNew[i]) - logD);
}
// magnitude of the normal velocity of the first splashed parcel
// Magnitude of the normal velocity of the first splashed parcel
const scalar magUns0 =
sqrt(2.0*parcelsPerSplash_*EKs/mSplash/(1.0 + coeff1/sqr(coeff2)));
......@@ -448,7 +448,7 @@ void Foam::ThermoSurfaceFilm<CloudType>::splashInteraction
pPtr->typeId() = splashParcelType_;
}
// perturb new parcels towards the owner cell centre
// Perturb new parcels towards the owner cell centre
pPtr->position() += 0.5*rndGen_.sample01<scalar>()*(posC - posCf);
pPtr->nParticle() = npNew[i];
......@@ -466,7 +466,7 @@ void Foam::ThermoSurfaceFilm<CloudType>::splashInteraction
nParcelsSplashed_++;
}
// transfer remaining part of parcel to film 0 - splashMass can be -ve
// Transfer remaining part of parcel to film 0 - splashMass can be -ve
// if entraining from the film
const scalar mDash = m - mSplash;
absorbInteraction(filmModel, p, pp, facei, mDash, keepParticle);
......@@ -613,11 +613,11 @@ bool Foam::ThermoSurfaceFilm<CloudType>::transferParcel
}
}
// transfer parcel/parcel interactions complete
// Transfer parcel/parcel interactions complete
return true;
}
// parcel not interacting with film
// Parcel not interacting with film
return false;
}
......
src/lagrangian/intermediate/submodels/Thermodynamic/SurfaceFilmModel/ThermoSurfaceFilm/ThermoSurfaceFilm.H
View file @ e14769f2
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -36,7 +36,7 @@ Description
Bai and Gosman, `Mathematical modelling of wall films formed by
impinging sprays', SAE 960626, 1996
Bai et al, `Modelling off gasoline spray impingement', Atom. Sprays,
Bai et al, `Modelling of gasoline spray impingement', Atom. Sprays,
vol 12, pp 1-27, 2002
......
Markdown is supported
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