yade-dev team mailing list archive

[noreply@xxxxxxxxxxxxx]

------------------------------------------------------------
revno: 2765
committer: Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
branch nick: yade
timestamp: Fri 2011-02-25 20:11:02 +0100
message:
  - per-point imposed pressure mechanism
  - one fix in permeabilites along boundaries
modified:
  lib/triangulation/FlowBoundingSphere.cpp
  lib/triangulation/FlowBoundingSphere.hpp


--
lp:yade
https://code.launchpad.net/~yade-dev/yade/trunk

Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-dev/yade/trunk/+edit-subscription

=== modified file 'lib/triangulation/FlowBoundingSphere.cpp'
--- lib/triangulation/FlowBoundingSphere.cpp	2011-02-17 17:43:37 +0000
+++ lib/triangulation/FlowBoundingSphere.cpp	2011-02-25 19:11:02 +0000
@@ -194,7 +194,7 @@
         Compute_Permeability();
         clock.top("Compute_Permeability");
         /** END PERMEABILITY CALCULATION**/
-	
+
 	if(DEBUG_OUT) cerr << "TOTAL VOID VOLUME: " << Vporale <<endl;
 	if(DEBUG_OUT) cerr << "Porosity = " << V_porale_porosity / V_totale_porosity << endl;
 
@@ -639,6 +639,7 @@
         /*CALCULATION OF VORONOI CENTRES*/
 //        if ( !NewTes.Computed() ) NewTes.Compute();
         for (Finite_cells_iterator new_cell = NewTri.finite_cells_begin(); new_cell != cell_end; new_cell++) {
+		if (new_cell->info().Pcondition) continue;
                 old_cell = Tri.locate((Point) new_cell->info());
                 new_cell->info().p() = old_cell->info().p();
         }
@@ -874,12 +875,10 @@
 	//         if (DEBUG_OUT) std::cerr << "total facet surface "<< cell->info().facetSurfaces[j] << " with solid sectors : " << cell->info().facetSphereCrossSections[j][0] << " " << cell->info().facetSphereCrossSections[j][1] << " " << cell->info().facetSphereCrossSections[j][2] << " difference "<<sqrt(cell->info().facetSurfaces[j].squared_length())-cell->info().facetSphereCrossSections[j][0]-cell->info().facetSphereCrossSections[j][2]-cell->info().facetSphereCrossSections[j][1]<<endl;
 // cerr << "test14" << endl;
 	//handle symmetry (tested ok)
-	if (/*SLIP_ON_LATERALS &&*/ fictious_vertex>0)
-	{
+	if (SLIP_ON_LATERALS && fictious_vertex>0) {
 		//! Include a multiplier so that permeability will be K/2 or K/4 in symmetry conditions
 		Real mult= fictious_vertex==1 ? multSym1 : multSym2;
-		return Vpore/Ssolid*mult;
-	}
+		return Vpore/Ssolid*mult;}
 // 	cerr << "test15" << endl;
 	return Vpore/Ssolid;
 }
@@ -904,6 +903,12 @@
 			{(*it)->info().p() = bi.value;(*it)->info().Pcondition=true;}
                 }
         }
+        for (int n=0; n<imposedP.size();n++) {
+		Cell_handle cell=Tri.locate(imposedP[n].first);
+// 		cerr<<"cell found : "<<cell->vertex(0)->point()<<" "<<cell->vertex(1)->point()<<" "<<cell->vertex(2)->point()<<" "<<cell->vertex(3)->point()<<endl;
+// 		assert(cell);
+		cell->info().p()=imposedP[n].second;
+		cell->info().Pcondition=true;}
 }
 
 void FlowBoundingSphere::GaussSeidel()

=== modified file 'lib/triangulation/FlowBoundingSphere.hpp'
--- lib/triangulation/FlowBoundingSphere.hpp	2011-02-17 17:43:37 +0000
+++ lib/triangulation/FlowBoundingSphere.hpp	2011-02-25 19:11:02 +0000
@@ -43,6 +43,7 @@
 		bool meanK_LIMIT, meanK_STAT, distance_correction;
 		bool OUTPUT_BOUDARIES_RADII;
 		bool noCache;//flag for checking if cached values cell->unitForceVectors have been defined
+		vector<pair<Point,Real> > imposedP;
 
 		void initNewTri () {noCache=true; /*isLinearSystemSet=false; areCellsOrdered=false;*/}//set flags after retriangulation