yade-dev team mailing list archive

[noreply@xxxxxxxxxxxxx]

------------------------------------------------------------
revno: 4024
committer: bchareyre <bruno.chareyre@xxxxxxxxxxxxxxx>
timestamp: Tue 2017-03-14 18:54:25 +0100
message:
  new data 'crackArea' added to DFNFlow cells and updated automatically - first step toward solving the artificial compressibility of fluids in cracks
modified:
  pkg/pfv/DFNFlow.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/DFNFlow.cpp'
--- pkg/pfv/DFNFlow.cpp	2017-03-07 17:06:18 +0000
+++ pkg/pfv/DFNFlow.cpp	2017-03-14 17:54:25 +0000
@@ -14,7 +14,7 @@
 //it will save compilation time for everyone else
 //when you want it compiled, you can pass -DDFNFLOW to cmake, or just uncomment the following line
 
-// #define DFNFLOW
+#define DFNFLOW
 
 #ifdef DFNFLOW
 #include "FlowEngine_DFNFlowEngineT.hpp"
@@ -24,10 +24,11 @@
 	public:
 	Real anotherVariable;
 	bool crack;
+	Real crackArea;// the volume of cracks
 // 	bool preExistingJoint;
 // 	void anotherFunction() {};
 // 	DFNCellInfo() : FlowCellInfo(),crack(false)  {}
-	DFNCellInfo() : crack(false)  {}
+	DFNCellInfo() : crack(false), crackArea(0) {}
 };
 
 class DFNVertexInfo : public FlowVertexInfo_DFNFlowEngineT {
@@ -127,19 +128,20 @@
 {
 	public :
 	void trickPermeability();
-	void trickPermeability (RTriangulation::Facet_circulator& facet,Real aperture, Real residualAperture);
+	void trickPermeability (RTriangulation::Facet_circulator& facet,Real aperture, Real residualAperture, RTriangulation::Finite_edges_iterator& edge);
 	void trickPermeability (RTriangulation::Finite_edges_iterator& edge,Real aperture, Real residualAperture);
 	void setPositionsBuffer(bool current);
 // 	void computeTotalFractureArea(Real totalFracureArea,bool printFractureTotalArea);/// Trying to get fracture's surface
 	Real totalFracureArea; /// Trying to get fracture's surface
+	CELL_SCALAR_GETTER(double,.crackArea,crackArea)
 
 	YADE_CLASS_BASE_DOC_ATTRS_INIT_CTOR_PY(DFNFlowEngine,DFNFlowEngineT,"This is an enhancement of the FlowEngine for intact and fractured rocks that takes into acount pre-existing discontinuities and bond breakage between particles. The local conductivity around the broken link is calculated according to parallel plates model",
 	((Real, jointResidual, 0,,"calibration parameter for residual aperture of joints"))
 	((bool, updatePositions, false,,"update particles positions when rebuilding the mesh (experimental)"))
  	((bool, printFractureTotalArea, 0,,"The final fracture area computed through the network")) /// Trying to get fracture's surface
-	,
-	,
-	,
+	((bool, calcCrackArea, true,,"The amount of crack per pore () is updated if calcCrackArea=True")) /// Trying to get fracture's surface
+	,,,
+	.def("getCrackArea",&DFNFlowEngine::crackArea,(boost::python::arg("id")),"get the cracked area within cell 'id'.")
 // 	.def("computeTotalFractureArea",&DFNFlowEngineT::computeTotalFractureArea," Compute and print the total fracture area of the network") /// Trying to get fracture's surface
 // 	.def("trickPermeability",&DFNFlowEngineT::trickPermeability,"measure the mean trickPermeability in the period")
 	)
@@ -168,7 +170,7 @@
 	}
 }
 
-void DFNFlowEngine::trickPermeability(RTriangulation::Facet_circulator& facet, Real aperture, Real residualAperture)
+void DFNFlowEngine::trickPermeability(RTriangulation::Facet_circulator& facet, Real aperture, Real residualAperture, RTriangulation::Finite_edges_iterator& ed_it)
 {
 	const RTriangulation::Facet& currentFacet = *facet; // seems verbose but facet->first was declaring a junk cell and crashing program (https://bugs.launchpad.net/yade/+bug/1666339)
 	const RTriangulation& Tri = solver->T[solver->currentTes].Triangulation();
@@ -188,15 +190,25 @@
 	totalFracureArea += networkFractureArea; /// Trying to get fracture's surface 
 // 	cout <<" ------------------ The total surface area up to here is --------------------" << totalFracureArea << endl;
 // 	printFractureTotalArea = totalFracureArea; /// Trying to get fracture's surface 
+	if (calcCrackArea) {
+			CVector edge = ed_it->first->vertex(ed_it->second)->point().point() - ed_it->first->vertex(ed_it->third)->point().point();
+			CVector unitV = edge*(1./sqrt(edge.squared_length()));
+			Point p3 = ed_it->first->vertex(ed_it->third)->point().point() + unitV*(cell1->info() - ed_it->first->vertex(ed_it->third)->point().point())*unitV;
+			Real halfCrackArea = 0.25*sqrt(std::abs(cross_product(CellCentre1-p3,CellCentre2-p3).squared_length()));//
+			cell1->info().crackArea += halfCrackArea;
+			cell2->info().crackArea += halfCrackArea;
+		}
 }
 
 void DFNFlowEngine::trickPermeability(RTriangulation::Finite_edges_iterator& edge, Real aperture, Real residualAperture)
 {
 	const RTriangulation& Tri = solver->T[solver->currentTes].Triangulation();
+	
+	
 	RTriangulation::Facet_circulator facet1 = Tri.incident_facets(*edge);
 	RTriangulation::Facet_circulator facet0=facet1++;
-	trickPermeability(facet0, aperture,residualAperture);
-	while ( facet1!=facet0 ) {trickPermeability(facet1, aperture, residualAperture); facet1++;}
+	trickPermeability(facet0, aperture,residualAperture, edge);
+	while ( facet1!=facet0 ) trickPermeability(facet1, aperture, residualAperture, edge); facet1++;
 	/// Needs the fracture surface for this edge?
 // 	double edgeArea = solver->T[solver->currentTes].computeVFacetArea(edge); cout<<"edge area="<<edgeArea<<endl;
 }