yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #11420
[Branch ~yade-pkg/yade/git-trunk] Rev 3415: -add temp function for generating sample windows, calculating responding saturation...
------------------------------------------------------------
revno: 3415
committer: cyuan <chaoyuan2012@xxxxxxxxx>
timestamp: Fri 2014-04-04 19:12:06 +0200
message:
-add temp function for generating sample windows, calculating responding saturation...
modified:
lib/triangulation/def_types.h
pkg/dem/UnsaturatedEngine.cpp
pkg/dem/UnsaturatedEngine.hpp
--
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 'lib/triangulation/def_types.h'
--- lib/triangulation/def_types.h 2014-02-17 18:00:12 +0000
+++ lib/triangulation/def_types.h 2014-04-04 17:12:06 +0000
@@ -230,12 +230,14 @@
std::vector<double> poreRadius;//pore throat radius for drainage
double solidLine [4][4];//the length of intersecting line between sphere and facet. [i][j] is for sphere facet "i" and sphere facetVertices[i][j]. Last component for 1/sum_Lines in the facet.
double trapCapP;//for calculating the pressure of trapped phase, cell.Pressure_trapped = Pressure_air - trapCapP.
+ int windowsID;//a temp cell info for experiment comparison
UnsatCellInfo (void)
{
poreRadius.resize(4, 0);
isWaterReservoir = false; isAirReservoir = false; capillaryCellVolume = 0;
for (int k=0; k<4;k++) for (int l=0; l<3;l++) solidLine[k][l]=0;
trapCapP = 0;
+ windowsID = 0;
}
};
=== modified file 'pkg/dem/UnsaturatedEngine.cpp'
--- pkg/dem/UnsaturatedEngine.cpp 2014-03-29 23:32:02 +0000
+++ pkg/dem/UnsaturatedEngine.cpp 2014-04-04 17:12:06 +0000
@@ -111,10 +111,10 @@
Real UnsaturatedEngine::getSaturation(Solver& flow )
{
if (isPhaseTrapped) {
- return getSaturation1(solver);
+ return getSaturation1(flow);
}
else {
- return getSaturation2(solver);
+ return getSaturation2(flow);
}
}
@@ -1158,7 +1158,7 @@
template<class Solver>
void UnsaturatedEngine::clearImposedPressure ( Solver& flow ) { flow->imposedP.clear(); flow->IPCells.clear();}
-//----tempt function for Vahid Joekar-Niasar's data----
+//----temp function for Vahid Joekar-Niasar's data----
template<class Cellhandle >
double UnsaturatedEngine::getRadiusMin(Cellhandle cell, int j)
{
@@ -1379,7 +1379,92 @@
file.close();
}
-//----------end tempt function for Vahid Joekar-Niasar's data (clear later)---------------------
+//----------end temp function for Vahid Joekar-Niasar's data (clear later)---------------------
+//----------temp functions for comparison with experiment-----------------------
+template <class Solver>
+void UnsaturatedEngine::initializeCellWindowsID(Solver&flow)
+{
+ RTriangulation& tri = flow->T[flow->currentTes].Triangulation();
+ Finite_cells_iterator cell_end = tri.finite_cells_end();
+ for ( Finite_cells_iterator cell = tri.finite_cells_begin(); cell != cell_end; cell++ ) {
+ for (int i=1; i<(windowsNo+1); i++) {
+ if ( (cell->info()[1]>(flow->y_min+(i-1)*(flow->y_max-flow->y_min)/windowsNo) ) && (cell->info()[1] < (flow->y_min+i*(flow->y_max-flow->y_min)/windowsNo)) )
+ {cell->info().windowsID=i; break;}
+ }
+ }
+}
+template<class Solver>
+Real UnsaturatedEngine::getWindowsSaturation(Solver&flow, int windowsID)
+{
+ RTriangulation& tri = flow->T[flow->currentTes].Triangulation();
+ Finite_cells_iterator cell_end = tri.finite_cells_end();
+ for ( Finite_cells_iterator cell = tri.finite_cells_begin(); cell != cell_end; cell++ ) {
+ if (cell->info().windowsID==0) {cerr<<"Please initialize windowsID"<<endl;break;}
+ }
+ if (isPhaseTrapped) {
+ return getWindowsSaturation1(flow,windowsID);
+ }
+ else {
+ return getWindowsSaturation2(flow,windowsID);
+ }
+}
+template<class Solver>
+Real UnsaturatedEngine::getWindowsSaturation1(Solver&flow, int i)
+{
+ updatePressureReservoir(flow);
+ RTriangulation& tri = flow->T[flow->currentTes].Triangulation();
+ Real capillary_volume = 0.0; //total capillary volume
+ Real air_volume = 0.0; //air volume
+ Finite_cells_iterator cell_end = tri.finite_cells_end();
+
+ if (invadeBoundary==true) {
+ for ( Finite_cells_iterator cell = tri.finite_cells_begin(); cell != cell_end; cell++ ) {
+ if (tri.is_infinite(cell)) continue;
+ if (cell->info().Pcondition) continue;//NOTE:reservoirs cells should not be included in saturation
+ if (cell->info().windowsID != i) continue;
+ capillary_volume = capillary_volume + cell->info().capillaryCellVolume;
+ if (cell->info().isAirReservoir==true) {
+ air_volume = air_volume + cell->info().capillaryCellVolume;}}}
+ else {
+ for ( Finite_cells_iterator cell = tri.finite_cells_begin(); cell != cell_end; cell++ ) {
+ if (tri.is_infinite(cell)) continue;
+ if (cell->info().Pcondition) continue;
+ if (cell->info().isFictious) continue;
+ if (cell->info().windowsID != i) continue;
+ capillary_volume = capillary_volume + cell->info().capillaryCellVolume;
+ if (cell->info().isAirReservoir==true) {
+ air_volume = air_volume + cell->info().capillaryCellVolume;}}}
+ Real saturation = 1 - air_volume/capillary_volume;
+ return saturation;
+}
+template<class Solver>
+Real UnsaturatedEngine::getWindowsSaturation2(Solver& flow,int i)
+{
+ RTriangulation& tri = flow->T[flow->currentTes].Triangulation();
+ Real capillary_volume = 0.0;
+ Real water_volume = 0.0;
+ Finite_cells_iterator cell_end = tri.finite_cells_end();
+ if (invadeBoundary==true) {
+ for ( Finite_cells_iterator cell = tri.finite_cells_begin(); cell != cell_end; cell++ ) {
+ if (tri.is_infinite(cell)) continue;
+ if (cell->info().Pcondition) continue;
+ if (cell->info().windowsID != i) continue;
+ capillary_volume = capillary_volume + cell->info().capillaryCellVolume;
+ if (cell->info().p()==0) {
+ water_volume = water_volume + cell->info().capillaryCellVolume;}}}
+ else {
+ for ( Finite_cells_iterator cell = tri.finite_cells_begin(); cell != cell_end; cell++ ) {
+ if (tri.is_infinite(cell)) continue;
+ if (cell->info().Pcondition) continue;
+ if (cell->info().isFictious) continue;
+ if (cell->info().windowsID != i) continue;
+ capillary_volume = capillary_volume + cell->info().capillaryCellVolume;
+ if (cell->info().p()==0) {
+ water_volume = water_volume + cell->info().capillaryCellVolume;}}}
+ Real saturation = water_volume/capillary_volume;
+ return saturation;
+}
+//----------end temp functions for comparison with experiment-------------------
template <class Solver>
void UnsaturatedEngine::computeFacetPoreForcesWithCache(Solver& flow, bool onlyCache)
=== modified file 'pkg/dem/UnsaturatedEngine.hpp'
--- pkg/dem/UnsaturatedEngine.hpp 2014-03-29 15:55:08 +0000
+++ pkg/dem/UnsaturatedEngine.hpp 2014-04-04 17:12:06 +0000
@@ -63,7 +63,10 @@
TPL void invade (Solver& flow );
TPL Real getMinEntryValue (Solver& flow);
TPL Real getSaturation(Solver& flow);
-
+ TPL Real getWindowsSaturation(Solver&flow, int windowsID);
+ TPL Real getWindowsSaturation1(Solver&flow, int i);
+ TPL Real getWindowsSaturation2(Solver&flow, int i);
+
TPL void invadeSingleCell1(Cell_handle cell, double pressure, Solver& flow);
TPL void invade1 (Solver& flow );
TPL void checkTrap(Solver& flow, double pressure);//check trapped phase, define trapCapP.
@@ -98,7 +101,7 @@
TPL void savePoreBodyInfo(Solver& flow);
TPL void savePoreThroatInfo(Solver& flow);
TPL void debugTemp(Solver& flow);
-
+ TPL void initializeCellWindowsID(Solver&flow);
TPL void computeSolidLine(Solver& flow);
TPL void computeFacetPoreForcesWithCache(Solver& flow, bool onlyCache=false);
@@ -147,6 +150,7 @@
void _invade() {invade(solver);}
Real _getMinEntryValue() {return getMinEntryValue(solver);}
Real _getSaturation() {return getSaturation(solver);}
+ Real _getWindowsSaturation(int windowsID) {return getWindowsSaturation(solver,windowsID);}
void _saveListNodes() {saveListNodes(solver);}
void _saveListConnection() {saveListConnection(solver);}
void _saveLatticeNodeX(double x) {saveLatticeNodeX(solver,x);}
@@ -157,6 +161,7 @@
void _savePoreBodyInfo(){savePoreBodyInfo(solver);}
void _savePoreThroatInfo(){savePoreThroatInfo(solver);}
void _debugTemp(){debugTemp(solver);}
+ void _initializeCellWindowsID(){initializeCellWindowsID(solver);}
void _computeFacetPoreForcesWithCache(){computeFacetPoreForcesWithCache(solver);}
Vector3r _fluidForce(unsigned int id_sph) {return fluidForce(id_sph,solver);}
@@ -192,6 +197,7 @@
((bool, pressureForce, true,,"Compute the pressure field and associated fluid forces. WARNING: turning off means fluid flow is not computed at all."))
((bool, computeForceActivated, true,,"Activate capillary force computation. WARNING: turning off means capillary force is not computed at all, but the drainage can still work."))
((bool, invadeBoundary, false,,"Invade from boundaries."))
+ ((int, windowsNo, 10,, "Number of genrated windows/(zoomed samples)."))
,
/*deprec*/
,,
@@ -211,7 +217,8 @@
.def("getCell",&UnsaturatedEngine::_getCell,(python::arg("pos")),"get id of the cell containing (X,Y,Z).")
.def("testFunction",&UnsaturatedEngine::testFunction,"The playground for Chao's experiments.")
.def("buildTriangulation",&UnsaturatedEngine::_buildTriangulation,"Triangulate spheres of the current scene.")
- .def("getSaturation",&UnsaturatedEngine::_getSaturation,"get saturation.")
+ .def("getSaturation",&UnsaturatedEngine::_getSaturation,"get saturation.")
+ .def("getWindowsSaturation",&UnsaturatedEngine::_getWindowsSaturation,(python::arg("windowsID")), "get saturation of windowsID")
.def("getMinEntryValue",&UnsaturatedEngine::_getMinEntryValue,"get the minimum air entry pressure for the next invade step.")
.def("saveListNodes",&UnsaturatedEngine::_saveListNodes,"Save the list of nodes.")
.def("saveListConnection",&UnsaturatedEngine::_saveListConnection,"Save the connections between cells.")
@@ -223,6 +230,7 @@
.def("savePoreBodyInfo",&UnsaturatedEngine::_savePoreBodyInfo,"Save pore bodies positions/Voronoi centers and size/volume.")
.def("savePoreThroatInfo",&UnsaturatedEngine::_savePoreThroatInfo,"Save pore throat area, inscribed radius and perimeter.")
.def("debugTemp",&UnsaturatedEngine::_debugTemp,"debug temp file.")
+ .def("initializeCellWindowsID",&UnsaturatedEngine::_initializeCellWindowsID,"Initialize cell windows index. A temp function for comparison with experiments, will delete soon")
.def("invade",&UnsaturatedEngine::_invade,"Run the drainage invasion.")
.def("computeForce",&UnsaturatedEngine::_computeFacetPoreForcesWithCache,"Compute capillary force. ")
.def("fluidForce",&UnsaturatedEngine::_fluidForce,(python::arg("Id_sph")),"Return the fluid force on sphere Id_sph.")