CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Programming & Development

Diffusive flux in interMixingFoam

Register Blogs Community New Posts Updated Threads Search

Like Tree7Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   June 7, 2011, 12:12
Default Diffusive flux in interMixingFoam
  #1
Member
 
Join Date: Sep 2010
Location: Leipzig, Germany
Posts: 96
Rep Power: 16
oswald is on a distinguished road
Dear Foamers,

this is my first thread started, so here is a little introduction: I've been working with OF for almost one year now and my current version is 1.7.0. My main interests are multiphase flows and in the last time a little bit of mixing.

I've got a little question about the diffusive flux for alpha2->alpha3 in interMixingFoam. In alphaEqns.H there is the line

Code:
phiAlpha2 -= fvc::interpolate(Dc32)*mesh.magSf()*fvc::snGrad(alpha1)
Should it not be snGrad(alpha2) or snGrad(alpha3) at the end?

Sorry for bothering you if I'm wrong and thanks for your help in advance!

Best regards,
oswald
oswald is offline   Reply With Quote

Old   May 24, 2013, 10:12
Default bug in interMixingFoam?
  #2
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
Any news? I wonder about the same thing. Simulating channel flow with the three phases, alpha3 seems to dissolve slowly into alpha1 while alpha2 stays niceley seperated from alpha1 as it should! If anyone please could have a look befor I do a bug report?
vonboett is offline   Reply With Quote

Old   June 4, 2013, 11:46
Default
  #3
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
...and by the way, maybe someone could shed some light on:

// Add the diffusive flux for alpha3->alpha2
phiAlpha2 -= fvc::interpolate(Dc32)*mesh.magSf()*fvc::snGrad(al pha3);

were only Dc32 adds to diffusive flux but followed by:

// Solve for alpha2
fvScalarMatrix alpha2Eqn
(
fvm::ddt(alpha2)
+ fvc::div(phiAlpha2)
- fvm::laplacian(Dc23 + Dc32, alpha2)
);

Dc23 and Dc32 contribute. I am sorry if this is trivial, maybe I am a bit confused looking and comparing the different multiphase solvers. Inspired by twoLiquidMixingFoam I just
gave it a try to add the influence of turbulence via the turbulent schmidt number to the diffusion in interMixingFoam, and I noticed I have to account for it in the diffusive flux for alpha3->alpha 2 and in the alpha2 equation.
vonboett is offline   Reply With Quote

Old   June 26, 2013, 06:29
Default improved interMixingFoam in respect to mixing
  #4
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
So here is my adjusted alphaEqns.H and createFields.H with the following features:

-the two liquids are solved for and alpha1 is derived from 1 - alpha2 - alpha3

-the diffusive flux is added by *mesh.magSf()*fvc::snGrad(alpha2); and *mesh.magSf()*fvc::snGrad(alpha3); respectiveley, no alpha1 used here (still think this is a bug in interMixingFoam)

-eddy diffusion is included either by local turbulence viscosity and turbulent Schmidt number, or for laminar flow by local vorticity (grid dependent!) and mixing length derived in dependency to distance to boundary. So the reciprocal of the turbulent Schmidt number is introduced to transportProperties as well as a mixingLengthFactor which can be seen as an adjustable von Karman constant.

alphaEqns.H:

{
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");

surfaceScalarField phir
(
IOobject
(
"phir",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
interface.cAlpha()*mag(phi/mesh.magSf())*interface.nHatf() //phase fraction * flux/cellsurface * cell face unit interface normal flux
);
// calculate vorticity to increase diffusion at regions of high vorticity
volScalarField vorticityMag
(
IOobject
(
"vorticityMag",
runTime.timeName(),
mesh,
IOobject::NO_READ
),
mag(fvc::curl(U))
);

for (int gCorr=0; gCorr<nAlphaCorr; gCorr++)
{
// Create the limiter to be used for all phase-fractions
scalarField allLambda(mesh.nFaces(), 1.0);

// Split the limiter into a surfaceScalarField
slicedSurfaceScalarField lambda
(
IOobject
(
"lambda",
mesh.time().timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh,
dimless,
allLambda,
false // Use slices for the couples
);


// Create the complete flux for alpha2
surfaceScalarField phiAlpha2
(
fvc::flux
(
phi,
alpha2,
alphaScheme
)
+ fvc::flux
(
-fvc::flux(phir, alpha1, alpharScheme),
alpha2,
alpharScheme
)
+ fvc::flux
(
-fvc::flux(-phir, alpha3, alpharScheme),
alpha2,
alpharScheme
)
);

// Create the bounded (upwind) flux for alpha2
surfaceScalarField phiAlpha2BD
(
upwind<scalar>(mesh, phi).flux(alpha2)
);

// Calculate the flux correction for alpha2
phiAlpha2 -= phiAlpha2BD;

// Calculate the limiter for alpha2
MULES::limiter
(
allLambda,
geometricOneField(),
alpha2,
phiAlpha2BD,
phiAlpha2,
zeroField(),
zeroField(),
1,
0,
3
);


// Construct the limited fluxes
phiAlpha2 = phiAlpha2BD + lambda*phiAlpha2;

// Create the complete flux for alpha3
surfaceScalarField phiAlpha3
(
fvc::flux
(
phi,
alpha3,
alphaScheme
)
+ fvc::flux
(
-fvc::flux(phir, alpha1, alpharScheme),
alpha3,
alpharScheme
)
);

// Create the bounded (upwind) flux for alpha3
surfaceScalarField phiAlpha3BD
(
upwind<scalar>(mesh, phi).flux(alpha3)
);

// Calculate the flux correction for alpha3
phiAlpha3 -= phiAlpha3BD;

// Further limit the limiter for alpha3
MULES::limiter
(
allLambda,
geometricOneField(),
alpha3,
phiAlpha3BD,
phiAlpha3,
zeroField(),
zeroField(),
1,
0,
3
);

// Construct the limited fluxes
phiAlpha3 = phiAlpha3BD + lambda*phiAlpha3;

volScalarField Dc23(D23*max(alpha3, scalar(0))*pos(alpha2));
volScalarField Dc32(D23*max(alpha2, scalar(0))*pos(alpha3));
/* calculate the turbulent mixing, orienting by P. Nielsen et al.:
"Turbulent diffusion of momentum and suspended particles: A finite-mixing-length theory" (2004) Physics of Fluids Vol. 16 Issue 7
but take the vorticity at the cell times the average cell size as a replacement for the settling velocity used in the paper. For mixing length:
By assuming that the mixing length is proportional to the distance to the bottom,
lm=κz (κ is von Karman constant, z is height above bed),Prandtl obtained the logarithmic velocity profile.
Assume the same but replace the von Karman constant with a model parameter, analog to the turbulent Schmidt number. See also Debris Flow: Mechanics, Prediction and Countermeasures (2007) by
T. Takahashi, chapter Models for Mechanics of Flow, p.83*/

volScalarField Vc23(mixingLengthFactor*vorticityMag*cellDiameter* wallDistance*max(alpha3, scalar(0))*pos(alpha2));
volScalarField Vc32(mixingLengthFactor*vorticityMag*cellDiameter* wallDistance*max(alpha2, scalar(0))*pos(alpha3));

Info<< "Phase-diffusion due to vorticity = "
<< " Min(vorticityMag*cellsize*wallDistance*mixingLengt hFactor) = " << (min(Vc32)).value()
<< " Max(vorticityMag*cellsize*wallDistance*mixingLengt hFactor) = " << (max(Vc32)).value()
<< endl;

// Add the diffusive flux for alpha3->alpha2
phiAlpha2 -= fvc::interpolate(Dc32 + Vc32 + alphatab*turbulence->nut() )*mesh.magSf()*fvc::snGrad(alpha2);
// Solve for alpha2
fvScalarMatrix alpha2Eqn
(
fvm::ddt(alpha2)
+ fvc::div(phiAlpha2)
+ fvm::SuSp(-fvc::div(phi), alpha2)//chegdan on October 16, 2010: http://www.cfd-online.com/Forums/ope...culations.html
- fvm::laplacian(Dc32 + Dc23 + Vc23 + Vc32 + alphatab*turbulence->nut() , alpha2)
);
alpha2Eqn.solve();
// Add the diffusive flux for alpha2->alpha3
phiAlpha3 -= fvc::interpolate(Dc23 + Vc23 + alphatab*turbulence->nut())*mesh.magSf()*fvc::snGrad(alpha3);
// Solve for alpha3
//alphatab is the reciprocal of the turbulent Schmidt number, introduce turbulent mixing see book by Fox "computational models for turbulent reacting flows"
fvScalarMatrix alpha3Eqn
(
fvm::ddt(alpha3)
+ fvc::div(phiAlpha3)
+ fvm::SuSp(-fvc::div(phi), alpha3)//chegdan on October 16, 2010: http://www.cfd-online.com/Forums/ope...culations.html
- fvm::laplacian(Dc23 + Dc32 + Vc23 + Vc32 + alphatab*turbulence->nut() , alpha3)
);
alpha3Eqn.solve();

// Construct the complete mass flux */
rhoPhi =
(phiAlpha3 + alpha3Eqn.flux())*(rho3 - rho1)
+ (phiAlpha2 + alpha2Eqn.flux())*(rho2 - rho1)
+ phi*rho1;

alpha1 = 1.0 - alpha2 - alpha3;
}

Info<< "Mud phase volume fraction = "
<< alpha2.weightedAverage(mesh.V()).value()
<< " Min(alpha2) = " << min(alpha2).value()
<< " Max(alpha2) = " << max(alpha2).value()
<< endl;

Info<< "Granular phase volume fraction = "
<< alpha3.weightedAverage(mesh.V()).value()
<< " Min(alpha3) = " << min(alpha3).value()
<< " Max(alpha3) = " << max(alpha3).value()
<< endl;

Info<< "Phase-diffusion due to turbulence = "
<< " Min(alphatab*turb->nut()) = " << (min(turbulence->nut())*alphatab).value()
<< " Max(alphatab*turb->nut()) = " << (max(turbulence->nut())*alphatab).value()
<< endl;
}







createFields.H

Info<< "Reading field p_rgh\n" << endl;
volScalarField p_rgh
(
IOobject
(
"p_rgh",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);

Info<< "Reading field alpha1\n" << endl;
volScalarField alpha1
(
IOobject
(
"alpha1",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);


Info<< "Reading field alpha2\n" << endl;
volScalarField alpha2
(
IOobject
(
"alpha2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);


Info<< "Reading field alpha3\n" << endl;
volScalarField alpha3
(
IOobject
(
"alpha3",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);

alpha3 == 1.0 - alpha1 - alpha2;


Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);

#include "createPhi.H"

threePhaseMixture threePhaseProperties(U, phi);

const dimensionedScalar& rho1 = threePhaseProperties.rho1();
const dimensionedScalar& rho2 = threePhaseProperties.rho2();
const dimensionedScalar& rho3 = threePhaseProperties.rho3();

dimensionedScalar D23(threePhaseProperties.lookup("D23"));

// Read the reciprocal of the turbulent Schmidt number
dimensionedScalar alphatab(threePhaseProperties.lookup("alphatab"));
// Read a constant in analogy to the von Karman constant, to gain mixing length in dependency to cell-distance to wall
dimensionedScalar mixingLengthFactor(threePhaseProperties.lookup("mi xingLengthFactor"));
// calculate an average cell size used to gain a mixing velocity in each cell from cell-vorticity, applied for eddy diffsuion
volScalarField cellDiameter(pow(0.1666*fvc::surfaceSum(mesh.magSf ()), 0.5));

// get cell distances to walls
volScalarField wallDistance(cellDiameter);
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(wallDistance, cellI)
{
// get the position for the particle
const vector& position = mesh.C()[cellI];

scalar distance = GREAT;
// loop over all patches
forAll(patches, patchI)
{
//if (isA<wallPolyPatch>(patches[patchI])) // not working
{
const polyPatch& patch = patches[patchI];
const vectorField& n = patch.faceNormals();
const vectorField& C = patch.faceCentres();
forAll(C, i)
{
vector normal = n[i]/mag(n[i]);
//scalar di = mag(position - C[i]); // distance to face centre
scalar di = mag( normal & ( position - C[i] ) );
distance = min(distance, di);
}
}
}
wallDistance[cellI] = distance;
// Info << "distance to wall is " << distance << endl;
}

// Need to store rho for ddt(rho, U)
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT
),
alpha1*rho1 + alpha2*rho2 + alpha3*rho3,
alpha1.boundaryField().types()
);
rho.oldTime();


// Mass flux
// Initialisation does not matter because rhoPhi is reset after the
// alpha solution before it is used in the U equation.
surfaceScalarField rhoPhi
(
IOobject
(
"rho*phi",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
rho1*phi
);


// Construct interface from alpha distribution
threePhaseInterfaceProperties interface(threePhaseProperties);


// Construct incompressible turbulence model
autoPtr<incompressible::turbulenceModel> turbulence
(
incompressible::turbulenceModel::New(U, phi, threePhaseProperties)
);


Info<< "Calculating field g.h\n" << endl;
volScalarField gh("gh", g & mesh.C());
surfaceScalarField ghf("ghf", g & mesh.Cf());

volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
p_rgh + rho*gh
);

label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell
(
p,
p_rgh,
mesh.solutionDict().subDict("PIMPLE"),
pRefCell,
pRefValue
);

if (p_rgh.needReference())
{
p += dimensionedScalar
(
"p",
p.dimensions(),
pRefValue - getRefCellValue(p, pRefCell)
);
p_rgh = p - rho*gh;
}

IObasicSourceList sources(mesh);





Attatched a comparison of the dambreak tutorial at t= 0.75 without turbulence but with mixingLengthFactor 0.4 increasing the diffusion due to vorticity. Anyone knows a testCase to test the solver, or maybe someone could check the code?
Attached Images
File Type: jpg debrisInterMixingFoamTut_T0,75.jpg (26.9 KB, 149 views)
File Type: jpg interMixingFoamTut_T0,75.jpg (26.4 KB, 131 views)
vonboett is offline   Reply With Quote

Old   July 2, 2013, 04:26
Default
  #5
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
Ok, in larger scale and longer runs, alpha2 and alpha3 diffuse into alpha1 when applying diffudive terms, and the massbalance (due to bounding?) is not quite right yet in the suggestion above because the liquid volume phase goes down continuously. But guess what, the regular interMixingFoam faces the same problem...
Any hints? I will try out what happends without bounding of alpha2 or alpha3.

Last edited by vonboett; July 11, 2013 at 12:04.
vonboett is offline   Reply With Quote

Old   September 26, 2013, 06:11
Default
  #6
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
Well, it works by now conserving the mass of the two liquids and keeping their concentrations bound between 0 and 1. But I had to introduce a solving step for the air phase. Will post the code here when tests are finished.
vonboett is offline   Reply With Quote

Old   February 5, 2014, 12:18
Default Showing both of the mixing fluids
  #7
New Member
 
anonymous
Join Date: Jan 2012
Location: Canada
Posts: 24
Rep Power: 14
Mahyar Javidi is on a distinguished road
Hi

I was wondering if anybody knows how to show both of the phases that are mixing with each other (in the interMixingFoam solver) in the paraview simultaneously.

Regards
Mahyar
Mahyar Javidi is offline   Reply With Quote

Old   February 6, 2014, 08:44
Default
  #8
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
Well I make a threshold for alpha1 (non-mixing phase) from 0 (or maybe -0.0001) to 0.5 and color it by alpha2, then I take a threshold of this threshold for alpha3 (full range from min to max) and set the representation to Wireframe and color it by alpha3
vonboett is offline   Reply With Quote

Old   February 9, 2014, 15:00
Default
  #9
New Member
 
anonymous
Join Date: Jan 2012
Location: Canada
Posts: 24
Rep Power: 14
Mahyar Javidi is on a distinguished road
Thanks for your reply. I appreciate it.

Mahyar
Mahyar Javidi is offline   Reply With Quote

Old   June 17, 2014, 06:13
Post about the solve of interMixingFoam
  #10
New Member
 
k zhang
Join Date: Jun 2012
Location: China
Posts: 9
Rep Power: 14
skykzhang is on a distinguished road
these days i am working in this solver. openFoam 1.6-ext. the solver interMixingFoam have no UEqn.H, so i try to look for it in the .dep file,

Quote:
interMixingFoam.dep: $(WM_PROJECT_DIR)/applications/solvers/multiphase/interFoam/UEqn.H
interMixingFoam.dep: $(WM_PROJECT_DIR)/applications/solvers/multiphase/interFoam/pEqn.H
$(OBJECTS_DIR)/interMixingFoam.o: $(EXE_DEP)
$(OBJECTS_DIR)/interMixingFoam.o:
@SOURCE_DIR=.
SOURCE=interMixingFoam.C ; $(Ctoo)
does it call to the UEqn.H from interFoam? interFoam is a two phase solver and interMixingFoam is three phase solver, the alpha in interFoam just one and in interMixing that are alpha1, alpha2 and alpha3. So can you tell me why the UEqn.H can call from interFoam??
skykzhang is offline   Reply With Quote

Old   June 17, 2014, 06:38
Default
  #11
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
Hi Zhang,

yes it calls UEqn.H from interFoam, thats OK since one solves only one phase-averaged Navier-Stokes equation neglecting drag between phases. I would like to know who maintains this solver. To my experience, it is not performing well in mass conservation or in respect to negative phase concentrations. Let me attatch my version of alphaEqns.H that gives higher priority to correct handling of the mixing fluid, shifting the weak parts to the air phase. This is for OF 2.2.x but a version for the newes OF exists. In my case, diffusion is not so important so I don't really care about the diffusive flux. (Sorry for not yet removing old comments from the code)
Attached Files
File Type: h alphaEqns.H (5.6 KB, 58 views)
Elham likes this.
vonboett is offline   Reply With Quote

Old   June 17, 2014, 08:27
Thumbs up
  #12
New Member
 
k zhang
Join Date: Jun 2012
Location: China
Posts: 9
Rep Power: 14
skykzhang is on a distinguished road
Quote:
Originally Posted by vonboett View Post
Hi Zhang,
thank u for your reply, vBoetticher. i have known the reason why it can call the UEqn.H from the interFoam.
and from your replies above, we can see you have done a lot of job about the solver.
the same to u , i have not know the guy who write the solver. from the origin of OF, i guess you can ask Hrvoje Jasak, maybe you can get the answer you need. good luck!

skykzhang
skykzhang is offline   Reply With Quote

Old   July 19, 2016, 06:43
Default
  #13
Member
 
Join Date: May 2016
Posts: 39
Rep Power: 10
dzordz is on a distinguished road
Quote:
Originally Posted by oswald View Post
Dear Foamers,

this is my first thread started, so here is a little introduction: I've been working with OF for almost one year now and my current version is 1.7.0. My main interests are multiphase flows and in the last time a little bit of mixing.

I've got a little question about the diffusive flux for alpha2->alpha3 in interMixingFoam. In alphaEqns.H there is the line

Code:
phiAlpha2 -= fvc::interpolate(Dc32)*mesh.magSf()*fvc::snGrad(alpha1)
Should it not be snGrad(alpha2) or snGrad(alpha3) at the end?

Sorry for bothering you if I'm wrong and thanks for your help in advance!

Best regards,
oswald
Hello,
I think I'm a little late to this party, but I hope it is still relevant. I have been playing around with interMixingFoam for a while now and since I had problems with spontaneous creation of a third phase plus not a sharp interface between third - diffusive and first-non diffusive phase (just like tutorial) I went digging into the source code.

So to answer your question. The way it is written with alpha1 is ok.
It is because you have a 3 phase system with condition alpha1+alpha2+alpha3=1.
So if you rewrite your code line you would get Dc32*grad(alpha1)=Dc32*grad(1-alpha2-alpha3)=-Dc32*grad(alpha2)-Dc32*grad(alpha3). Now the trick is in the way the laplacian part of differential equation is defined: laplacian(Dc23+Dc32, alpha2). this can be written in terms of divergence: div((Dc23+Dc32)*grad(alpha2)). When you combine this with the divergence term in diff. eq. and subtract the equal parts you get (here written only the diffusion part): div(Dc23*grad(alpha3) - Dc32*grad(alpha2). Exactly what you wanted - two diffusion fluxes.

Now if you would write down equations for alpha3 you would be left with div(Dc32*grad(alpha2) - Dc23*grad(alpha3)). Exactly the opposite fluxes. And now if you add the fluxes, the would sum to zero. Exactly as they should in the closed system.

Sorry for the long response. Hope it helps.
dzordz is offline   Reply With Quote

Old   July 28, 2016, 04:51
Default
  #14
Member
 
Join Date: May 2016
Posts: 39
Rep Power: 10
dzordz is on a distinguished road
Quote:
Originally Posted by dzordz View Post
... not a sharp interface between third - diffusive and first-non diffusive phase (just like tutorial) ...
...if anybody is still interested

After some days spent in the actual code I was able to find a problem. The problem is with the way mules is written into the solver. The complete flux for phase 1 is constructed and then the unbounded part of the flux is bounded by mules (this is done by calculating lambda coefficients). All good. Then the second flux is constructed and the unbounded part is again limited by mules. But it uses the same lambda coefficients (it basically overwrites it). This results in the non-sharp interface for the phase 3.

Solution: give second Mules different coefficients.

Cheers
dzordz is offline   Reply With Quote

Old   September 7, 2016, 05:33
Default
  #15
Member
 
Join Date: May 2016
Posts: 39
Rep Power: 10
dzordz is on a distinguished road
Quote:
Originally Posted by dzordz View Post
Solution: give second Mules different coefficients.
Better jet. Construct the limited flux for alpha 1 before you start creating the complete flux for alpha 2.
MrFyzzix likes this.
dzordz is offline   Reply With Quote

Old   April 5, 2018, 02:53
Default
  #16
Senior Member
 
Elham
Join Date: Oct 2009
Posts: 184
Rep Power: 17
Elham is on a distinguished road
Quote:
Originally Posted by dzordz View Post
Better jet. Construct the limited flux for alpha 1 before you start creating the complete flux for alpha 2.

Hi,

Do you know what is D23 or D32 in interMixingFoam?

Regards,

Elham
Elham is offline   Reply With Quote

Old   April 5, 2018, 03:46
Default
  #17
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
It is the diffusion constant that causes phase 2 to mix with phase 3 and phase 3 with phase 2. To my understanding, D23 should be equal to D32 to make sense.
Anyway, since there is no drift flux or any other momentum exchange implemented, be aware that the flux you create with D23 or D32 in the alpha equations does not account for the momentum transfer, so the thing is only valid for small D23 or D32. I once tried to link them to eddy diffusivity, but I stopped following that way due to the neglected momentum.
vonboett is offline   Reply With Quote

Old   April 5, 2018, 04:19
Default
  #18
Senior Member
 
Elham
Join Date: Oct 2009
Posts: 184
Rep Power: 17
Elham is on a distinguished road
Quote:
Originally Posted by vonboett View Post
It is the diffusion constant that causes phase 2 to mix with phase 3 and phase 3 with phase 2. To my understanding, D23 should be equal to D32 to make sense.
Anyway, since there is no drift flux or any other momentum exchange implemented, be aware that the flux you create with D23 or D32 in the alpha equations does not account for the momentum transfer, so the thing is only valid for small D23 or D32. I once tried to link them to eddy diffusivity, but I stopped following that way due to the neglected momentum.

Dear Albrecht,

So I am wondering if I put zero value for D23, there is no mixing between phases?
I have three phases and mixing is not a matter for me.
Another question is that if each phase continuity equation driven by the total continuity equation? I mean:

from the following eq:

ddt(rho) + div (rho.U) = 0

where

rho = alpha1*rho1 + alpha2*rho2+alpha3*rho3

Regards,

Elham
Elham is offline   Reply With Quote

Old   April 5, 2018, 11:52
Default
  #19
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
Yes, I would use D23 = D32 = 0 to block the diffusivity between the two miscible phases. InterMixingFoam uses the Volume-Of-Fluid method, where it keeps one phase unmixed, while a second phase can have two mixing components (like oil and water+alcohol, etc.). Between the unmixed phase an the other two phases, a surface is captured and here a surface drag term transfers momentum between the unmixing phase and the two other mixed phases. Due to the fluxes, all phases are advected in the alpha equations, but flux through the surface between the mixing and the non-mixing phase is corrected. Then, the whole distribution of phases is handled as if it would be one single phase with composition-weighted properties in each cell (density, viscosity etc) and then a single Navier-Stokes equation is solved for this mixture based on the PISO loop, where you make use of a volumetric continuity (incompressible fluids) in the pressure equation. This approach is valid for mixtures where you do not have drag between mixing phases and it is about 10 times faster then if you solve one Navier-Stokes equation set for each phase.
vonboett is offline   Reply With Quote

Old   April 5, 2018, 11:56
Default
  #20
Senior Member
 
Albrecht vBoetticher
Join Date: Aug 2010
Location: Zürich, Swizerland
Posts: 240
Rep Power: 17
vonboett is on a distinguished road
one note though: interMixingFoam is either keeping the overall mass right but allows negative concentrations or volumetric concentrations greater than one, or you have the concentrations limited to the physical range but may loose mass...
I tried to optimize this a bit in debrisInterMixing-2.3, you can download the code with the corresponding paper if you want to do the same.
vonboett is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Problem setting with chtmultiregionFoam Antonin OpenFOAM 10 April 24, 2012 10:50
Zero Diffusive flux Subramani Sockalingam FLUENT 1 May 27, 2009 11:36
Zero Diffusive flux Subramani Sockalingam Main CFD Forum 0 January 19, 2008 14:02
Help: diffusive flux BC at wall Quarkz Main CFD Forum 7 July 17, 2005 12:24
Replace periodic by inlet-outlet pair lego CFX 3 November 5, 2002 21:09


All times are GMT -4. The time now is 21:16.