yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #07833
[Branch ~yade-dev/yade/trunk] Rev 2903: - FlowEngine does not depend on TriaxialCompressionEngine anymore
------------------------------------------------------------
revno: 2903
committer: Emanuele Catalano <catalano@xxxxxxxxxxx
branch nick: yade
timestamp: Mon 2011-08-22 17:37:27 +0200
message:
- FlowEngine does not depend on TriaxialCompressionEngine anymore
- Clean of useless variables
modified:
lib/triangulation/FlowBoundingSphere.cpp
pkg/dem/FlowEngine.cpp
pkg/dem/FlowEngine.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-07-18 13:53:54 +0000
+++ lib/triangulation/FlowBoundingSphere.cpp 2011-08-22 15:37:27 +0000
@@ -308,7 +308,6 @@
int num_cells = 0;
double facet_flow_rate = 0;
double volume_facet_translation = 0;
- std::ofstream oFile ( "Average_Cells_Velocity",std::ios::app );
Real tVel=0; Real tVol=0;
Finite_cells_iterator cell_end = Tri.finite_cells_end();
for ( Finite_cells_iterator cell = Tri.finite_cells_begin(); cell != cell_end; cell++ ) {
@@ -493,19 +492,18 @@
cell->info().isvisited=!ref;
}
if (DEBUG_OUT) {
- cout << "Facet scheme" <<endl;
+// cout << "Facet scheme" <<endl;
Vecteur TotalForce = nullVect;
for (Finite_vertices_iterator v = Tri.finite_vertices_begin(); v != Tri.finite_vertices_end(); ++v) {
if (!v->info().isFictious) {
TotalForce = TotalForce + v->info().forces;
- cout << "real_id = " << v->info().id() << " force = " << v->info().forces << endl;
+// cout << "real_id = " << v->info().id() << " force = " << v->info().forces << endl;
} else {
if (boundary(v->info().id()).flowCondition==1) TotalForce = TotalForce + v->info().forces;
- cout << "fictious_id = " << v->info().id() << " force = " << v->info().forces << endl;
+// cout << "fictious_id = " << v->info().id() << " force = " << v->info().forces << endl;
}
}
- cout << "TotalForce = "<< TotalForce << endl;
- }
+ cout << "TotalForce = "<< TotalForce << endl;}
}
void FlowBoundingSphere::ComputeFacetForcesWithCache()
@@ -569,20 +567,19 @@
noCache=false;//cache should always be defined after execution of this function
for (Finite_vertices_iterator v = Tri.finite_vertices_begin(); v != Tri.finite_vertices_end(); ++v) v->info().forces = 0*oldForces[v->info().id()]+1*v->info().forces;
if (DEBUG_OUT) {
- cout << "Facet cached scheme" <<endl;
+// cout << "Facet cached scheme" <<endl;
Vecteur TotalForce = nullVect;
for (Finite_vertices_iterator v = Tri.finite_vertices_begin(); v != Tri.finite_vertices_end(); ++v)
{
if (!v->info().isFictious) {
TotalForce = TotalForce + v->info().forces;
- if (DEBUG_OUT) cout << "real_id = " << v->info().id() << " force = " << v->info().forces << endl;
+// cout << "real_id = " << v->info().id() << " force = " << v->info().forces << endl;
} else {
if (boundary(v->info().id()).flowCondition==1) TotalForce = TotalForce + v->info().forces;
- if (DEBUG_OUT) cout << "fictious_id = " << v->info().id() << " force = " << v->info().forces << endl;
+// cout << "fictious_id = " << v->info().id() << " force = " << v->info().forces << endl;
}
}
- cout << "TotalForce = "<< TotalForce << endl;
- }
+ cout << "TotalForce = "<< TotalForce << endl;}
}
void FlowBoundingSphere::ComputeTetrahedralForces()
@@ -622,17 +619,17 @@
}
cell->info().isvisited=!ref;
}
- if (DEBUG_OUT) cout << "tetrahedral scheme" <<endl;
- Vecteur TotalForce = nullVect;
- for (Finite_vertices_iterator v = Tri.finite_vertices_begin(); v != Tri.finite_vertices_end(); ++v) {
- if (!v->info().isFictious) {
- TotalForce = TotalForce + v->info().forces;
- } else {
- if (boundary(v->info().id()).flowCondition==1) TotalForce = TotalForce + v->info().forces;
- if (DEBUG_OUT) cout << "fictious_id = " << v->info().id() << " force = " << v->info().forces << endl;
- }
- }
- cout << "TotalForce = "<< TotalForce << endl;
+// if (DEBUG_OUT) cout << "tetrahedral scheme" <<endl;
+// Vecteur TotalForce = nullVect;
+// for (Finite_vertices_iterator v = Tri.finite_vertices_begin(); v != Tri.finite_vertices_end(); ++v) {
+// if (!v->info().isFictious) {
+// TotalForce = TotalForce + v->info().forces;
+// } else {
+// if (boundary(v->info().id()).flowCondition==1) TotalForce = TotalForce + v->info().forces;
+// if (DEBUG_OUT) cout << "fictious_id = " << v->info().id() << " force = " << v->info().forces << endl;
+// }
+// }
+// if (DEBUG_OUT) cout << "TotalForce = "<< TotalForce << endl;
}
void FlowBoundingSphere::ApplySinusoidalPressure(RTriangulation& Tri, double Amplitude, double Average_Pressure, double load_intervals)
@@ -1198,7 +1195,7 @@
cout << "celle comunicanti in alto = " << cellQ1 << endl;}
double density = 1;
- double viscosity = 1.0;
+ double viscosity = VISCOSITY;
double gravity = 1;
double Vdarcy = Q1/Section;
double DeltaP = abs(P_Inf-P_Sup);
@@ -1227,8 +1224,8 @@
kFile << "The permeability of the sample is = " << k << " m^2" <<endl;
kFile << "The Darcy permeability of the sample is = " << k/9.869233e-13 << " darcys" << endl;
- cout << endl << "The hydraulic conductivity of the sample is = " << Ks << " m/s" << endl << endl;
- return Ks;
+ cout << endl << "The hydraulic conductivity of the sample is = " << k << " m^2" << endl << endl;
+ return k;
}
void FlowBoundingSphere::DisplayStatistics()
=== modified file 'pkg/dem/FlowEngine.cpp'
--- pkg/dem/FlowEngine.cpp 2011-07-18 13:53:54 +0000
+++ pkg/dem/FlowEngine.cpp 2011-08-22 15:37:27 +0000
@@ -43,18 +43,6 @@
ReTrg=1;
retriangulationLastIter=0;
}
- if (!triaxialCompressionEngine)
- {
- vector<shared_ptr<Engine> >::iterator itFirst = scene->engines.begin();
- vector<shared_ptr<Engine> >::iterator itLast = scene->engines.end();
- for (;itFirst!=itLast; ++itFirst) {
- if ((*itFirst)->getClassName() == "TriaxialCompressionEngine") {
-// cout << "stress controller engine found - FlowEngine" << endl;
- triaxialCompressionEngine = YADE_PTR_CAST<TriaxialCompressionEngine> (*itFirst);}}
- if (!triaxialCompressionEngine) cout << "stress controller engine NOT found" << endl;
- }
- currentStress = triaxialCompressionEngine->stress[triaxialCompressionEngine->wall_top][1];
- currentStrain = triaxialCompressionEngine->strain[1];
timingDeltas->start();
@@ -71,8 +59,6 @@
ReTrg++;
} else retriangulationLastIter++;
timingDeltas->checkpoint("Update_Volumes");
- ///Update boundary conditions
- BoundaryConditions();
///Compute flow and and forces here
flow->GaussSeidel();
@@ -229,7 +215,6 @@
flow->Vtotalissimo=0; flow->Vsolid_tot=0; flow->Vporale=0; flow->Ssolid_tot=0;
flow->SLIP_ON_LATERALS=slip_boundary;
- flow->key = triaxialCompressionEngine->Key;
flow->k_factor = permeability_factor;
flow->DEBUG_OUT = Debug;
flow->useSolver = useSolver;
@@ -304,22 +289,13 @@
flow->SectionArea = ( flow->x_max - flow->x_min ) * ( flow->z_max-flow->z_min );
flow->Vtotale = (flow->x_max-flow->x_min) * (flow->y_max-flow->y_min) * (flow->z_max-flow->z_min);
- if (triaxialCompressionEngine) {
- flow->y_min_id=triaxialCompressionEngine->wall_bottom_id;
- flow->y_max_id=triaxialCompressionEngine->wall_top_id;
- flow->x_max_id=triaxialCompressionEngine->wall_right_id;
- flow->x_min_id=triaxialCompressionEngine->wall_left_id;
- flow->z_min_id=triaxialCompressionEngine->wall_back_id;
- flow->z_max_id=triaxialCompressionEngine->wall_front_id;
- } else {
- flow->y_min_id=wallBottomId;
- flow->y_max_id=wallTopId;
- flow->x_max_id=wallRightId;
- flow->x_min_id=wallLeftId;
- flow->z_min_id=wallBackId;
- flow->z_max_id=wallFrontId;
- }
-
+ flow->y_min_id=wallBottomId;
+ flow->y_max_id=wallTopId;
+ flow->x_max_id=wallRightId;
+ flow->x_min_id=wallLeftId;
+ flow->z_min_id=wallBackId;
+ flow->z_max_id=wallFrontId;
+
flow->boundary ( flow->y_min_id ).useMaxMin = BOTTOM_Boundary_MaxMin;
flow->boundary ( flow->y_max_id ).useMaxMin = TOP_Boundary_MaxMin;
flow->boundary ( flow->x_max_id ).useMaxMin = RIGHT_Boundary_MaxMin;
@@ -336,8 +312,6 @@
flow->boundsIds[4]= &flow->z_min_id;
flow->boundsIds[5]= &flow->z_max_id;
- wall_thickness = triaxialCompressionEngine->thickness;
-
flow->Corner_min = CGT::Point(flow->x_min, flow->y_min, flow->z_min);
flow->Corner_max = CGT::Point(flow->x_max, flow->y_max, flow->z_max);
@@ -360,7 +334,7 @@
for (int i=0; i<6; i++)
{
- CGT::Vecteur Normal (triaxialCompressionEngine->normal[i].x(), triaxialCompressionEngine->normal[i].y(), triaxialCompressionEngine->normal[i].z());
+ CGT::Vecteur Normal (normal[i].x(), normal[i].y(), normal[i].z());
if (flow->boundary(*flow->boundsIds[i]).useMaxMin) flow->AddBoundingPlane (true, Normal, *flow->boundsIds[i]);
else {
const shared_ptr<Body>& wll = Body::byId ( *flow->boundsIds[i] , scene );
=== modified file 'pkg/dem/FlowEngine.hpp'
--- pkg/dem/FlowEngine.hpp 2011-07-18 13:53:54 +0000
+++ pkg/dem/FlowEngine.hpp 2011-08-22 15:37:27 +0000
@@ -25,16 +25,13 @@
class FlowEngine : public PartialEngine
{
private:
- shared_ptr<TriaxialCompressionEngine> triaxialCompressionEngine;
shared_ptr<FlowSolver> flow;
int retriangulationLastIter;
public :
- Vector3r gravity;
- int current_state;
- Real wall_thickness;
+ enum {wall_left=0, wall_right, wall_bottom, wall_top, wall_back, wall_front};
+ Vector3r normal [6];
bool currentTes;
int id_offset;
- double wall_up_y, wall_down_y;
double Eps_Vol_Cumulative;
int ReTrg;
void Triangulate ();
@@ -87,6 +84,7 @@
((double, Sinus_Average, 0,,"Pressure value (average) when sinusoidal pressure is applied"))
((bool, CachedForces, true,,"Des/Activate the cached forces calculation"))
((bool, Debug, false,,"Activate debug messages"))
+ ((double, wall_thickness,0.001,,"Walls thickness"))
((double,P_zero,0,,"Initial internal pressure for oedometer test"))
((double,Tolerance,1e-06,,"Gauss-Seidel Tolerance"))
((double,Relax,1.9,,"Gauss-Seidel relaxation"))
@@ -138,6 +136,12 @@
((bool, viscousShear, false,,"Compute viscous shear terms as developped by Donia Marzougui"))
,,
timingDeltas=shared_ptr<TimingDeltas>(new TimingDeltas);
+ normal[wall_bottom].y()=1;
+ normal[wall_top].y()=-1;
+ normal[wall_left].x()=1;
+ normal[wall_right].x()=-1;
+ normal[wall_front].z()=-1;
+ normal[wall_back].z()=1;
,
.def("imposeFlux",&FlowEngine::imposeFlux,(python::arg("pos"),python::arg("p")),"Impose incoming flux in boundary cell of location 'pos'.")
.def("imposePressure",&FlowEngine::imposePressure,(python::arg("pos"),python::arg("p")),"Impose pressure in cell of location 'pos'. The index of the condition is returned (for multiple imposed pressures at different points).")
