------------------------------------------------------------
revno: 3978
committer: bchareyre <bruno.chare...@grenoble-inp.fr>
timestamp: Wed 2016-11-23 15:35:05 +0100
message:
merging PhaseCluster code with r3975
modified:
pkg/pfv/TwoPhaseFlowEngine.cpp
pkg/pfv/TwoPhaseFlowEngine.hpp
pkg/pfv/UnsaturatedEngine.cpp
--
lp:yade
https://code.launchpad.net/~yade-pkg/yade/git-trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to
https://code.launchpad.net/~yade-pkg/yade/git-trunk/+edit-subscription
=== modified file 'pkg/pfv/TwoPhaseFlowEngine.cpp'
--- pkg/pfv/TwoPhaseFlowEngine.cpp 2016-11-23 10:43:14 +0000
+++ pkg/pfv/TwoPhaseFlowEngine.cpp 2016-11-23 14:35:05 +0000
@@ -17,6 +17,9 @@
#ifdef TWOPHASEFLOW
YADE_PLUGIN((TwoPhaseFlowEngineT));
YADE_PLUGIN((TwoPhaseFlowEngine));
+YADE_PLUGIN((PhaseCluster));
+
+PhaseCluster::~PhaseCluster(){}
void TwoPhaseFlowEngine::initialization()
{
@@ -618,24 +621,35 @@
{
RTriangulation& tri = solver->T[solver->currentTes].Triangulation();
FiniteCellsIterator cellEnd = tri.finite_cells_end();
+ clusters[0]->reset(); clusters[1]->reset();
for ( FiniteCellsIterator cell = tri.finite_cells_begin(); cell != cellEnd; cell++ ) {
- if (cell->info().isNWRes) cell->info().label=0;
- else if (cell->info().isWRes) cell->info().label=1;
+ if (cell->info().isNWRes) {cell->info().label=0; clusters[0]->pores.push_back(cell);}
+ else if (cell->info().isWRes) {cell->info().label=1; clusters[1]->pores.push_back(cell);
+ for (int facet = 0; facet < 4; facet ++) if (!cell->neighbor(facet)->info().isWRes) clusterGetFacet(clusters[1].get(),cell,facet);}
else if (cell->info().label>1) continue;
else cell->info().label=-1;
}
}
-void TwoPhaseFlowEngine:: updateCellLabel()
+void TwoPhaseFlowEngine::clusterGetFacet(PhaseCluster* cluster, CellHandle cell, int facet) {
+ cell->info().hasInterface = true;
+ cluster->interfacialArea += cell->info().poreThroatRadius[facet];//FIXME: define area correctly
+ if (cluster->entryRadius < cell->info().poreThroatRadius[facet]){
+ cluster->entryRadius = cell->info().poreThroatRadius[facet];
+ cluster->entryPore = cell->info().index;
+ }
+}
+
+void TwoPhaseFlowEngine::updateCellLabel()
{
- int currentLabel = getMaxCellLabel();
+ int currentLabel = getMaxCellLabel();//FIXME: A loop on cells for each new label?? is it serious??
updateReservoirLabel();
RTriangulation& tri = solver->T[solver->currentTes].Triangulation();
FiniteCellsIterator cellEnd = tri.finite_cells_end();
for ( FiniteCellsIterator cell = tri.finite_cells_begin(); cell != cellEnd; cell++ ) {
if (cell->info().label==-1) {
- clusters.push_back(shared_ptr<PhaseCluster>(new PhaseCluster));
- updateSingleCellLabelRecursion(cell,currentLabel+1,clusters.back()->pores);
+ clusters.push_back(shared_ptr<PhaseCluster> (new PhaseCluster()));
+ updateSingleCellLabelRecursion(cell,currentLabel+1,clusters.back().get());
currentLabel++;
}
}
@@ -652,10 +666,10 @@
return maxLabel;
}
-void TwoPhaseFlowEngine:: updateSingleCellLabelRecursion(CellHandle cell, int label, std::vector<CellHandle> *outputVector)
+void TwoPhaseFlowEngine::updateSingleCellLabelRecursion(CellHandle cell, int label, PhaseCluster* cluster)
{
cell->info().label=label;
- if (outputVector) outputVector->push_back(cell);
+ cluster->pores.push_back(cell);
for (int facet = 0; facet < 4; facet ++) {
CellHandle nCell = cell->neighbor(facet);
if (solver->T[solver->currentTes].Triangulation().is_infinite(nCell)) continue;
@@ -663,7 +677,10 @@
// if ( (nCell->info().isFictious) && (!isInvadeBoundary) ) continue;
//TODO:the following condition may relax to relate to nCell->info().hasInterface
if ( (nCell->info().saturation==cell->info().saturation) && (nCell->info().label!=cell->info().label) )
- updateSingleCellLabelRecursion(nCell,label,outputVector);
+ updateSingleCellLabelRecursion(nCell,label,cluster);
+ else {
+ clusterGetFacet(cluster,cell,facet);
+ }
}
}
=== modified file 'pkg/pfv/TwoPhaseFlowEngine.hpp'
--- pkg/pfv/TwoPhaseFlowEngine.hpp 2016-11-23 10:43:14 +0000
+++ pkg/pfv/TwoPhaseFlowEngine.hpp 2016-11-23 14:35:05 +0000
@@ -62,22 +62,23 @@
class PhaseCluster : public Serializable
{
double totalCellVolume;
+// CellHandle entryPoreHandle;
public :
virtual ~PhaseCluster();
- vector<TwoPhaseFlowEngineT::CellHandle>* pores;
+ vector<TwoPhaseFlowEngineT::CellHandle> pores;
TwoPhaseFlowEngineT::RTriangulation* tri;
+ void reset() {label=entryPore=-1;volume=entryRadius=interfacialArea=0; pores.clear();}
+
YADE_CLASS_BASE_DOC_ATTRS_INIT_CTOR_PY(PhaseCluster,Serializable,"Preliminary.",
((int,label,-1,,"Unique label of this cluster, should be reflected in pores of this cluster."))
((double,volume,0,,"cumulated volume of all pores."))
- ((double,entryPc,0,,"smallest entry capillary pressure."))
- ((int,entryPore,0,,"the pore of the cluster incident to the throat with smallest entry Pc."))
+ ((double,entryRadius,0,,"smallest entry capillary pressure."))
+ ((int,entryPore,-1,,"the pore of the cluster incident to the throat with smallest entry Pc."))
((double,interfacialArea,0,,"interfacial area of the cluster"))
,,,
)
};
REGISTER_SERIALIZABLE(PhaseCluster);
-YADE_PLUGIN((PhaseCluster));
-PhaseCluster::~PhaseCluster(){}
class TwoPhaseFlowEngine : public TwoPhaseFlowEngineT
{
@@ -126,7 +127,7 @@
void checkTrap(double pressure);
void updateReservoirLabel();
void updateCellLabel();
- void updateSingleCellLabelRecursion(CellHandle cell, int label, std::vector<CellHandle> *outputVector=NULL);
+ void updateSingleCellLabelRecursion(CellHandle cell, int label, PhaseCluster* cluster);
int getMaxCellLabel();
void invasion2();//without-trap
@@ -153,6 +154,11 @@
}
bool detectBridge(RTriangulation::Finite_edges_iterator& edge);
+ boost::python::list pyClusters() { boost::python::list ret;
+ for(vector<shared_ptr<PhaseCluster> >::iterator it=clusters.begin(); it!=clusters.end(); ++it) ret.append(*it);
+ return ret;}
+ void clusterGetFacet(PhaseCluster* cluster, CellHandle cell, int facet);//update cluster inetrfacial area and max entry radius wrt to a facet
+
//post-processing
void savePoreNetwork();
void saveVtk(const char* folder) {bool initT=solver->noCache; solver->noCache=false; solver->saveVtk(folder); solver->noCache=initT;}
@@ -210,7 +216,7 @@
((bool, isImbibitionActivated, false,, "Activates imbibition."))
- ,/*TwoPhaseFlowEngineT()*/,
+ ,/*clusters.resize(2);*//*TwoPhaseFlowEngineT()*/,
,
.def("getCellIsFictious",&TwoPhaseFlowEngine::cellIsFictious,"Check the connection between pore and boundary. If true, pore throat connects the boundary.")
.def("setCellIsNWRes",&TwoPhaseFlowEngine::setCellIsNWRes,"set status whether 'wetting reservoir' state")
@@ -239,6 +245,7 @@
.def("computeCapillaryForce",&TwoPhaseFlowEngine::computeCapillaryForce,"Compute capillary force. ")
.def("saveVtk",&TwoPhaseFlowEngine::saveVtk,(boost::python::arg("folder")="./VTK"),"Save pressure field in vtk format. Specify a folder name for output.")
.def("getPotentialPendularSpheresPair",&TwoPhaseFlowEngine::getPotentialPendularSpheresPair,"Get the list of sphere ID pairs of potential pendular liquid bridge.")
+ .def("getClusters",&TwoPhaseFlowEngine::pyClusters/*,(boost::python::arg("folder")="./VTK")*/,"Get the list of clusters.")
)
DECLARE_LOGGER;
=== modified file 'pkg/pfv/UnsaturatedEngine.cpp'
--- pkg/pfv/UnsaturatedEngine.cpp 2016-11-23 10:43:14 +0000
+++ pkg/pfv/UnsaturatedEngine.cpp 2016-11-23 14:35:05 +0000
@@ -37,27 +37,7 @@
double getSpecificInterfacialArea();
void printSomething();
-
- boost::python::list getPotentialPendularSpheresPair() {
- RTriangulation& Tri = solver->T[solver->currentTes].Triangulation();
- boost::python::list bridgeIds;
- FiniteEdgesIterator ed_it = Tri.finite_edges_begin();
- for ( ; ed_it!=Tri.finite_edges_end(); ed_it++ ) {
- if (detectBridge(ed_it)==true) {
- const VertexInfo& vi1=(ed_it->first)->vertex(ed_it->second)->info();
- const VertexInfo& vi2=(ed_it->first)->vertex(ed_it->third)->info();
- const int& id1 = vi1.id();
- const int& id2 = vi2.id();
- bridgeIds.append(boost::python::make_tuple(id1,id2));}}
- return bridgeIds;}
- bool detectBridge(RTriangulation::Finite_edges_iterator& edge);
-
- boost::python::list pyClusters() { boost::python::list ret;
- for(vector<shared_ptr<PhaseCluster> >::iterator it=clusters.begin(); it!=clusters.end(); ++it) ret.append(*it);
- return ret;}
-
virtual ~UnsaturatedEngine();
-
virtual void action();
YADE_CLASS_BASE_DOC_ATTRS_INIT_CTOR_PY(UnsaturatedEngine,TwoPhaseFlowEngine,"Preliminary version engine of a drainage model for unsaturated soils. Note:Air reservoir is on the top; water reservoir is on the bottom.(deprecated engine, use TwoPhaseFlowEngine instead)",
@@ -74,8 +54,6 @@
.def("getWindowsSaturation",&UnsaturatedEngine::getWindowsSaturation,(boost::python::arg("windowsID"),boost::python::arg("isSideBoundaryIncluded")), "get saturation of subdomain with windowsID. If isSideBoundaryIncluded=false (default), the pores of side boundary are excluded in saturation calculating; if isSideBoundaryIncluded=true (only in isInvadeBoundary=true drainage mode), the pores of side boundary are included in saturation calculating.")
.def("initializeCellWindowsID",&UnsaturatedEngine::initializeCellWindowsID,"Initialize cell windows index. A temporary function for comparison with experiments, will delete soon")
.def("printSomething",&UnsaturatedEngine::printSomething,"print debug.")
- .def("getPotentialPendularSpheresPair",&UnsaturatedEngine::getPotentialPendularSpheresPair,"Get the list of sphere ID pairs of potential pendular liquid bridge.")
- .def("getClusters",&UnsaturatedEngine::pyClusters/*,(boost::python::arg("folder")="./VTK")*/,"Get the list of clusters.")
)
DECLARE_LOGGER;
};
_______________________________________________
Mailing list: https://launchpad.net/~yade-dev
Post to : yade-dev@lists.launchpad.net
Unsubscribe : https://launchpad.net/~yade-dev
More help : https://help.launchpad.net/ListHelp
