Skip to content

Commit f8a7d50

Browse files
PhasicFlowPhasicFlow
authored
Merge pull request #50 from hamidrezanorouzi/main
diffusion model for porosity is added, base class porosity is modified
2 parents 8d1b8be + 7d68de1 commit f8a7d50

File tree

14 files changed

+345
-165
lines changed

14 files changed

+345
-165
lines changed

‎phasicFlowCoupling/Make/files‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@ couplingSystem/porosity/PIC.C
1212
couplingSystem/porosity/subDivision9.C
1313
couplingSystem/porosity/subDivision29.C
1414
couplingSystem/porosity/subDivision29Mod.C
15+
couplingSystem/porosity/diffusion.C
1516
couplingSystem/interaction/drag/drag.C
1617
couplingSystem/interaction/drag/ErgunWenYu.C
1718
couplingSystem/interaction/drag/DiFelice.C

‎phasicFlowCoupling/couplingSystem/couplingMesh.C‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -147,13 +147,13 @@ bool pFlow::coupling::couplingMesh::checkForDomainUpdate
147147
return true;
148148
}
149149

150-
if( abs(t-lastTimeUpdated_) < 0.98*fluidDt )
150+
if( std::abs(t-lastTimeUpdated_) < static_cast<real>(0.98*fluidDt) )
151151
{
152152
lastTimeUpdated_ = t;
153153
return true;
154154
}
155155

156-
if( abs(t-(lastTimeUpdated_+domainUpdateInterval_)) < 0.98*fluidDt)
156+
if( std::abs(t-(lastTimeUpdated_+domainUpdateInterval_)) < static_cast<real>(0.98*fluidDt))
157157
{
158158
lastTimeUpdated_ = t;
159159
return true;

‎phasicFlowCoupling/couplingSystem/porosity/PIC.C‎

Lines changed: 3 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -48,24 +48,20 @@ bool pFlow::coupling::PIC::internalFieldUpdate()
4848

4949
auto& solidVol = solidVoldTmp.ref();
5050

51-
numInMesh_ = 0;
5251
size_t numPar = centerMass_.size();
5352

54-
#pragma omp parallel for reduction(+:numInMesh_)
53+
#pragma omp parallel for
5554
for(size_t i=0; i<numPar; i++)
5655
{
57-
58-
auto cellId = cMesh_.findCellTree(centerMass_[i], parCellIndex_[i]);
56+
const auto cellId = parCellIndex_[i];
5957
if( cellId >= 0 )
6058
{
6159
#pragma omp atomic
6260
solidVol[cellId] +=
6361
static_cast<real>(3.14159265358979/6)*
6462
pFlow::pow(particleDiameter_[i], static_cast<real>(3.0));
65-
numInMesh_++;
63+
6664
}
67-
parCellIndex_[i] = cellId;
68-
6965
}
7066

7167
this->ref() = Foam::max(

‎phasicFlowCoupling/couplingSystem/porosity/PIC.hpp‎

Lines changed: 0 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -71,10 +71,6 @@ class PIC
7171

7272
bool internalFieldUpdate() override;
7373

74-
int32 numInMesh()const override
75-
{
76-
return numInMesh_;
77-
}
7874

7975
};
8076

‎phasicFlowCoupling/couplingSystem/porosity/diffusion.C‎

Lines changed: 135 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,135 @@
1+
/*------------------------------- phasicFlow ---------------------------------
2+
O C enter of
3+
O O E ngineering and
4+
O O M ultiscale modeling of
5+
OOOOOOO F luid flow
6+
------------------------------------------------------------------------------
7+
Copyright (C): www.cemf.ir
8+
email: hamid.r.norouzi AT gmail.com
9+
------------------------------------------------------------------------------
10+
Licence:
11+
This file is part of phasicFlow code. It is a free software for simulating
12+
granular and multiphase flows. You can redistribute it and/or modify it under
13+
the terms of GNU General Public License v3 or any other later versions.
14+
15+
phasicFlow is distributed to help others in their research in the field of
16+
granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the
17+
implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
18+
19+
-----------------------------------------------------------------------------*/
20+
21+
// from OpenFOAM
22+
#include "fvCFD.H"
23+
24+
25+
#include "diffusion.hpp"
26+
27+
28+
Foam::tmp<Foam::fvMatrix<Foam::scalar>> pFlow::coupling::diffusion::fvmDdt
29+
(
30+
const Foam::volScalarField& sField
31+
)
32+
{
33+
Foam::tmp<fvMatrix<Foam::scalar>> tfvm
34+
(
35+
new Foam::fvMatrix<Foam::scalar>
36+
(
37+
sField,
38+
sField.dimensions()*Foam::dimVol/Foam::dimTime
39+
)
40+
);
41+
42+
Foam::fvMatrix<Foam::scalar>& fvm = tfvm.ref();
43+
44+
const Foam::scalar rDeltaT = 1.0/dt_.value();
45+
46+
fvm.diag() = rDeltaT*sField.mesh().Vsc();
47+
48+
if (sField.mesh().moving())
49+
{
50+
fvm.source() = rDeltaT*sField.oldTime().primitiveField()*sField.mesh().Vsc0();
51+
}
52+
else
53+
{
54+
fvm.source() = rDeltaT*sField.oldTime().primitiveField()*sField.mesh().Vsc();
55+
}
56+
57+
return tfvm;
58+
}
59+
60+
pFlow::coupling::diffusion::diffusion(
61+
Foam::dictionary dict,
62+
couplingMesh& cMesh,
63+
MPI::centerMassField& centerMass,
64+
MPI::realProcCMField& parDiam)
65+
:
66+
PIC(dict, cMesh, centerMass, parDiam),
67+
nSteps_(Foam::max(1,dict.lookup<Foam::label>("nSteps"))),
68+
intTime_(dict.lookup<Foam::scalar>("intTime")),
69+
dt_("dt", Foam::dimTime, intTime_/nSteps_),
70+
picSolDict_("picSolDict")
71+
{
72+
73+
picSolDict_.add("relTol", 0);
74+
picSolDict_.add("tolerance", 1.0e-8);
75+
picSolDict_.add("solver", "smoothSolver");
76+
picSolDict_.add("smoother", "symGaussSeidel");
77+
}
78+
79+
80+
bool pFlow::coupling::diffusion::internalFieldUpdate()
81+
{
82+
83+
auto solidVoldTmp = Foam::volScalarField::Internal::New(
84+
"solidVol",
85+
this->mesh(),
86+
Foam::dimensioned("solidVol", Foam::dimVol, Foam::scalar(0))
87+
);
88+
89+
auto& solidVol = solidVoldTmp.ref();
90+
91+
size_t numPar = centerMass_.size();
92+
93+
#pragma omp parallel for
94+
for(size_t i=0; i<numPar; i++)
95+
{
96+
const auto cellId = parCellIndex_[i];
97+
if( cellId >= 0 )
98+
{
99+
#pragma omp atomic
100+
solidVol[cellId] +=
101+
static_cast<real>(3.14159265358979/6)*
102+
pFlow::pow(particleDiameter_[i], static_cast<real>(3.0));
103+
104+
}
105+
}
106+
107+
auto picAlphaTmp = Foam::volScalarField::New(
108+
"picAlpha",
109+
this->mesh(),
110+
Foam::dimensioned("picAlpha", Foam::dimless, Foam::scalar(0)),
111+
"zeroGradient"
112+
);
113+
114+
volScalarField& picAlpha = picAlphaTmp.ref();
115+
116+
117+
picAlpha.ref() = Foam::max( 1 - solidVol/this->mesh().V(), 0.0);
118+
picAlpha.correctBoundaryConditions();
119+
120+
121+
// start of Time loop
122+
for(Foam::label i=0; i<nSteps_; i++)
123+
{
124+
picAlpha.storeOldTime();
125+
fvScalarMatrix alphaEq
126+
(
127+
fvmDdt(picAlpha) - fvm::laplacian(picAlpha)
128+
);
129+
alphaEq.solve(picSolDict_);
130+
}
131+
132+
this->ref() = picAlpha.internalField();
133+
134+
return true;
135+
}

‎phasicFlowCoupling/couplingSystem/porosity/diffusion.hpp‎

Lines changed: 91 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,91 @@
1+
/*------------------------------- phasicFlow ---------------------------------
2+
O C enter of
3+
O O E ngineering and
4+
O O M ultiscale modeling of
5+
OOOOOOO F luid flow
6+
------------------------------------------------------------------------------
7+
Copyright (C): www.cemf.ir
8+
email: hamid.r.norouzi AT gmail.com
9+
------------------------------------------------------------------------------
10+
Licence:
11+
This file is part of phasicFlow code. It is a free software for simulating
12+
granular and multiphase flows. You can redistribute it and/or modify it under
13+
the terms of GNU General Public License v3 or any other later versions.
14+
15+
phasicFlow is distributed to help others in their research in the field of
16+
granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the
17+
implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
18+
19+
-----------------------------------------------------------------------------*/
20+
21+
#ifndef __diffusion_hpp__
22+
#define __diffusion_hpp__
23+
24+
#include "virtualConstructor.hpp"
25+
26+
// from phasicFlow-coupling
27+
#include "PIC.hpp"
28+
29+
30+
namespace pFlow::coupling
31+
{
32+
33+
/**
34+
* Particle In Cell (diffusion) model for porosity calculation
35+
*
36+
* This model only considers the particle center and if the particle center
37+
* resides inside a cell, it is assumed that the whole volume of the particle
38+
* is located in that cell.
39+
*
40+
*/
41+
class diffusion
42+
:
43+
public PIC
44+
{
45+
private:
46+
47+
Foam::label nSteps_;
48+
49+
Foam::scalar intTime_;
50+
51+
Foam::dimensionedScalar dt_;
52+
53+
Foam::dictionary picSolDict_;
54+
55+
Foam::tmp<Foam::fvMatrix<Foam::scalar>> fvmDdt
56+
(
57+
const Foam::volScalarField& sField
58+
);
59+
60+
public:
61+
62+
/// Type info
63+
TypeInfo("diffusion");
64+
65+
/// Construc from dictionary
66+
diffusion(
67+
Foam::dictionary dict,
68+
couplingMesh& cMesh,
69+
MPI::centerMassField& centerMass,
70+
MPI::realProcCMField& parDiam);
71+
72+
/// Destructor
73+
virtual ~diffusion() = default;
74+
75+
/// Add this constructor to the list of virtual constructors
76+
add_vCtor
77+
(
78+
porosity,
79+
diffusion,
80+
dictionary
81+
);
82+
83+
bool internalFieldUpdate() override;
84+
85+
86+
};
87+
88+
} // pFlow::coupling
89+
90+
91+
#endif

‎phasicFlowCoupling/couplingSystem/porosity/porosity.C‎

Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -24,6 +24,23 @@ Licence:
2424
#include "streams.hpp"
2525
#include "Timer.hpp"
2626

27+
28+
void pFlow::coupling::porosity::mapCenters()
29+
{
30+
numInMesh_ = 0;
31+
size_t numPar = centerMass_.size();
32+
33+
34+
#pragma omp parallel for reduction(+:numInMesh_)
35+
for(size_t i = 0; i<numPar; i++)
36+
{
37+
auto cellId = cMesh_.findCellTree(centerMass_[i], parCellIndex_[i]);
38+
parCellIndex_[i] = cellId;
39+
if( cellId >= 0 ) numInMesh_++;
40+
}
41+
}
42+
43+
2744
pFlow::coupling::porosity::porosity(
2845
Foam::dictionary dict,
2946
couplingMesh& cMesh,
@@ -61,6 +78,7 @@ void pFlow::coupling::porosity::calculatePorosity()
6178
{
6279
Timer t;
6380
t.start();
81+
mapCenters();
6482
this->internalFieldUpdate();
6583
t.end();
6684
output<<"mapping execution time " << t.lastTime()<<endl;

‎phasicFlowCoupling/couplingSystem/porosity/porosity.hpp‎

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -62,6 +62,10 @@ class porosity
6262
/// cell indices of particles in this processor
6363
MPI::procCMField<Foam::label> parCellIndex_;
6464

65+
int32 numInMesh_ = 0;
66+
67+
void mapCenters();
68+
6569
public:
6670

6771
/// Type info
@@ -150,8 +154,10 @@ class porosity
150154
bool internalFieldUpdate() = 0;
151155

152156
/// Return number of center mass points found in this mesh (processor)
153-
virtual
154-
int32 numInMesh()const = 0;
157+
int32 numInMesh()const
158+
{
159+
return numInMesh_;
160+
}
155161

156162
/// Report (output) number of center mass points found in all processors
157163
/// It is effective only in master processor

0 commit comments

Comments
 (0)