yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #06205
[Branch ~yade-dev/yade/trunk] Rev 2558: - Reorganized debuggin'
------------------------------------------------------------
revno: 2558
committer: Emanuele Catalano <catalano@xxxxxxxxxxx
branch nick: yade
timestamp: Wed 2010-11-17 21:07:21 +0100
message:
- Reorganized debuggin'
- Found and corrected a error in boundary creator
- Changed names of .h files
removed:
lib/triangulation/FlowBoundingSphere.h
lib/triangulation/Network.h
added:
lib/triangulation/FlowBoundingSphere.hpp
lib/triangulation/Network.hpp
modified:
lib/triangulation/FlowBoundingSphere.cpp
lib/triangulation/Network.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 2010-11-17 15:19:09 +0000
+++ lib/triangulation/FlowBoundingSphere.cpp 2010-11-17 20:07:21 +0000
@@ -17,11 +17,11 @@
#include <assert.h>
#include <sys/stat.h>
#include <sys/types.h>
-#include "Network.h"
+#include "Network.hpp"
#include <omp.h>
// #define XVIEW
-#include "FlowBoundingSphere.h"//include after #define XVIEW
+#include "FlowBoundingSphere.hpp"//include after #define XVIEW
#ifdef XVIEW
// #include "Vue3D.h" //FIXME implicit dependencies will look for this class (out of tree) even ifndef XVIEW
@@ -468,7 +468,6 @@
// int id1 = cell->vertex(facetVertices[j][y])->info().id();
// int id2 = neighbour_cell->vertex(facetVertices[Tri.mirror_index(cell,j)][y])->info().id();
// cerr <<"id1/id2 : "<<id1<<" "<<id2<<endl;}
-
///TEST END
}
}
@@ -1031,6 +1030,7 @@
#endif
p_moy = sum_p/cell2;
dp_moy = sum_dp/cell2;
+
#ifdef GS_OPEN_MP
#pragma omp master
#endif
@@ -1130,7 +1130,6 @@
cout << "The outgoing average flow rate is = " << Q1 << " m^3/s " << endl;
cout << "The gradient of charge is = " << GradH << " [-]" << endl;
cout << "Darcy's velocity is = " << Vdarcy << " m/s" <<endl;
-
cout << "The permeability of the sample is = " << k << " m^2" <<endl;}
kFile << "the maximum superior pressure is = " << p_in_max << " the min is = " << p_in_min << endl;
=== removed file 'lib/triangulation/FlowBoundingSphere.h'
--- lib/triangulation/FlowBoundingSphere.h 2010-11-17 15:19:09 +0000
+++ lib/triangulation/FlowBoundingSphere.h 1970-01-01 00:00:00 +0000
@@ -1,183 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2010 by Emanuele Catalano <catalano@xxxxxxxxxxxxxxx> *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-#ifdef FLOW_ENGINE
-
-#ifndef _FLOWBOUNDINGSPHERE_H
-#define _FLOWBOUNDINGSPHERE_H
-
-#include "Operations.h"
-#include "Timer.h"
-#include "Tesselation.h"
-#include "basicVTKwritter.hpp"
-#include "Timer.h"
-#include "stdafx.h"
-#include "Empilement.h"
-#include "Network.h"
-
-#ifdef XVIEW
-#include "Vue3D.h" //FIXME implicit dependencies will look for this class (out of tree) even ifndef XVIEW
-#endif
-
-
-
-using namespace std;
-
-
-
-namespace CGT{
-
-class FlowBoundingSphere : public Network
-{
- public:
- virtual ~FlowBoundingSphere();
- FlowBoundingSphere();
-
-// int x_min_id, x_max_id, y_min_id, y_max_id, z_min_id, z_max_id;
-// int* boundsIds [6];
-// bool currentTes;
- bool SLIP_ON_LATERALS;
- double TOLERANCE;
- double RELAX;
- double ks; //Hydraulic Conductivity
- bool meanK_LIMIT, meanK_STAT, distance_correction;
-// bool DEBUG_OUT;
- bool OUTPUT_BOUDARIES_RADII;
- bool noCache;//flag for checking if cached values cell->unitForceVectors have been defined
-
- void initNewTri () {noCache=true; /*isLinearSystemSet=false; areCellsOrdered=false;*/}//set flags after retriangulation
-
- bool computeAllCells;//exececute computeHydraulicRadius for all facets and all spheres (double cpu time but needed for now in order to define crossSections correctly)
- double K_opt_factor;
- int Iterations;
-
- bool RAVERAGE;
-// Boundary boundaries [6];
- int walls_id[6];
-// short id_offset;
-// Boundary& boundary (int b) {return boundaries[b-id_offset];}
-
- void mplot ( char *filename);
- void Localize();
-
- void Compute_Permeability();
- void GaussSeidel ( );
-
-// void Compute_Forces ();
- void Fictious_cells ( );
-
-// Tesselation T [2];
-
-// double x_min, x_max, y_min, y_max, z_min, z_max, Rmoy;
-// Real Vsolid_tot, Vtotalissimo, Vporale, Ssolid_tot;
- double k_factor; //permeability moltiplicator
- std::string key; //to give to consolidation files a name with iteration number
- std::vector<double> Pressures; //for automatic write maximum pressures during consolidation
- bool tess_based_force; //allow the force computation method to be chosen from FlowEngine
- Real minPermLength; //min branch length for Poiseuille
-
- double P_SUP, P_INF, P_INS;
-
-// void AddBoundingPlanes ( Tesselation& Tes, double x_Min,double x_Max ,double y_Min,double y_Max,double z_Min,double z_Max );
-// void AddBoundingPlanes(bool yade);
-// void AddBoundingPlanes();
-// void AddBoundingPlanes(Real center[3], Real Extents[3], int id);
-
- Tesselation& Compute_Action ( );
- Tesselation& Compute_Action ( int argc, char *argv[ ], char *envp[ ] );
- Tesselation& LoadPositions(int argc, char *argv[ ], char *envp[ ]);
-// Vecteur external_force_single_fictious ( Cell_handle cell );
- void SpheresFileCreator ();
-// void Analytical_Consolidation ( );
- void DisplayStatistics();
-// void Boundary_Conditions ( RTriangulation& Tri );
- void Initialize_pressures ( double P_zero );
- /// Define forces using the same averaging volumes as for permeability
- void ComputeTetrahedralForces();
- /// Define forces spliting drag and buoyancy terms
- void ComputeFacetForces();
- void ComputeFacetForcesWithCache();
- void save_vtk_file ( );
- void MGPost ( );
-#ifdef XVIEW
- void Dessine_Triangulation ( Vue3D &Vue, RTriangulation &T );
- void Dessine_Short_Tesselation ( Vue3D &Vue, Tesselation &Tes );
-#endif
- double Permeameter ( double P_Inf, double P_Sup, double Section, double DeltaY, char *file );
- double Sample_Permeability( double& x_Min,double& x_Max ,double& y_Min,double& y_Max,double& z_Min,double& z_Max, string key);
- double Compute_HydraulicRadius (Cell_handle cell, int j );
- double PressureProfile ( char *filename, Real& time, int& intervals );
-
- double dotProduct ( Vecteur x, Vecteur y );
- double Compute_EffectiveRadius(Cell_handle cell, int j);
- double Compute_EquivalentRadius(Cell_handle cell, int j);
-
-// double crossProduct ( double x[3], double y[3] );
-
-// double surface_solid_facet ( Sphere ST1, Sphere ST2, Sphere ST3 );
-// Vecteur surface_double_fictious_facet ( Vertex_handle fSV1, Vertex_handle fSV2, Vertex_handle SV3 );
-// Vecteur surface_single_fictious_facet ( Vertex_handle fSV1, Vertex_handle SV2, Vertex_handle SV3 );
-
-// double surface_solid_double_fictious_facet ( Vertex_handle ST1, Vertex_handle ST2, Vertex_handle ST3 );
-
-// double surface_external_triple_fictious (Cell_handle cell, Boundary b );
-// double surface_external_triple_fictious ( Real position[3], Cell_handle cell, Boundary b );
-// double surface_external_double_fictious ( Cell_handle cell, Boundary b );
-
-// double surface_external_single_fictious ( Cell_handle cell, Boundary b );
-
- void GenerateVoxelFile ( );
-
- RTriangulation& Build_Triangulation ( Real x, Real y, Real z, Real radius, unsigned const id );
-
- void Build_Tessalation ( RTriangulation& Tri );
-
-// double spherical_triangle_area ( Sphere STA1, Sphere STA2, Sphere STA3, Point PTA1 );
-
-// double fast_spherical_triangle_area ( const Sphere& STA1, const Point& STA2, const Point& STA3, const Point& PTA1 );
-// Real solid_angle ( const Point& STA1, const Point& STA2, const Point& STA3, const Point& PTA1 );
-// double spherical_triangle_volume ( const Sphere& ST1, const Point& PT1, const Point& PT2, const Point& PT3 );
-// Real fast_solid_angle ( const Point& STA1, const Point& PTA1, const Point& PTA2, const Point& PTA3 );
-
- bool isInsideSphere ( double& x, double& y, double& z );
-
- void SliceField (char *filename);
- void ComsolField();
-
- void Interpolate ( Tesselation& Tes, Tesselation& NewTes );
-
-// double volume_single_fictious_pore ( Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1 );
- //Fast version, assign surface of facet for future forces calculations (pointing from PV2 to PV1)
-// double volume_single_fictious_pore ( const Vertex_handle& SV1, const Vertex_handle& SV2, const Vertex_handle& SV3, const Point& PV1, const Point& PV2, Vecteur& facetSurface);
-// double volume_double_fictious_pore ( Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1 );
- //Fast version, assign surface of facet for future forces calculations (pointing from PV2 to PV1)
-// double volume_double_fictious_pore (Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point& PV1, Point& PV2, Vecteur& facetSurface);
-
-// double PoreVolume (RTriangulation& Tri, Cell_handle cell);
- int Average_Cell_Velocity(int id_sphere, RTriangulation& Tri);
- void Average_Cell_Velocity();
- void Average_Grain_Velocity();
- void vtk_average_cell_velocity(RTriangulation &T, int id_sphere, int num_cells);
- void ApplySinusoidalPressure(RTriangulation& Tri, double Pressure, double load_intervals);
- void ApplySinusoidalPressure_Space_Time(RTriangulation& Tri, double Pressure, double load_intervals, double time, double dt);
-// double surface_external_triple_fictious(Real position[3], Cell_handle cell, Boundary b);
-// double surface_external_double_fictious(Cell_handle cell, Boundary b);
-// double surface_external_single_fictious(Cell_handle cell, Boundary b);
-
- ///TAUCS
- int SetLinearSystem();
- int SetLinearSystemFullGS();
- void VectorizedGaussSeidel ();
- int TaucsSolveTest();
- int PardisoSolveTest();
- ///END_TAUCS
-
-};
-
-} //namespace CGT
-#endif //FLOW_ENGINE
-
-#endif
=== added file 'lib/triangulation/FlowBoundingSphere.hpp'
--- lib/triangulation/FlowBoundingSphere.hpp 1970-01-01 00:00:00 +0000
+++ lib/triangulation/FlowBoundingSphere.hpp 2010-11-17 20:07:21 +0000
@@ -0,0 +1,183 @@
+/*************************************************************************
+* Copyright (C) 2010 by Emanuele Catalano <catalano@xxxxxxxxxxxxxxx> *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+#ifdef FLOW_ENGINE
+
+#ifndef _FLOWBOUNDINGSPHERE_H
+#define _FLOWBOUNDINGSPHERE_H
+
+#include "Operations.h"
+#include "Timer.h"
+#include "Tesselation.h"
+#include "basicVTKwritter.hpp"
+#include "Timer.h"
+#include "stdafx.h"
+#include "Empilement.h"
+#include "Network.hpp"
+
+#ifdef XVIEW
+#include "Vue3D.h" //FIXME implicit dependencies will look for this class (out of tree) even ifndef XVIEW
+#endif
+
+
+
+using namespace std;
+
+
+
+namespace CGT{
+
+class FlowBoundingSphere : public Network
+{
+ public:
+ virtual ~FlowBoundingSphere();
+ FlowBoundingSphere();
+
+// int x_min_id, x_max_id, y_min_id, y_max_id, z_min_id, z_max_id;
+// int* boundsIds [6];
+// bool currentTes;
+ bool SLIP_ON_LATERALS;
+ double TOLERANCE;
+ double RELAX;
+ double ks; //Hydraulic Conductivity
+ bool meanK_LIMIT, meanK_STAT, distance_correction;
+// bool DEBUG_OUT;
+ bool OUTPUT_BOUDARIES_RADII;
+ bool noCache;//flag for checking if cached values cell->unitForceVectors have been defined
+
+ void initNewTri () {noCache=true; /*isLinearSystemSet=false; areCellsOrdered=false;*/}//set flags after retriangulation
+
+ bool computeAllCells;//exececute computeHydraulicRadius for all facets and all spheres (double cpu time but needed for now in order to define crossSections correctly)
+ double K_opt_factor;
+ int Iterations;
+
+ bool RAVERAGE;
+// Boundary boundaries [6];
+ int walls_id[6];
+// short id_offset;
+// Boundary& boundary (int b) {return boundaries[b-id_offset];}
+
+ void mplot ( char *filename);
+ void Localize();
+
+ void Compute_Permeability();
+ void GaussSeidel ( );
+
+// void Compute_Forces ();
+ void Fictious_cells ( );
+
+// Tesselation T [2];
+
+// double x_min, x_max, y_min, y_max, z_min, z_max, Rmoy;
+// Real Vsolid_tot, Vtotalissimo, Vporale, Ssolid_tot;
+ double k_factor; //permeability moltiplicator
+ std::string key; //to give to consolidation files a name with iteration number
+ std::vector<double> Pressures; //for automatic write maximum pressures during consolidation
+ bool tess_based_force; //allow the force computation method to be chosen from FlowEngine
+ Real minPermLength; //min branch length for Poiseuille
+
+ double P_SUP, P_INF, P_INS;
+
+// void AddBoundingPlanes ( Tesselation& Tes, double x_Min,double x_Max ,double y_Min,double y_Max,double z_Min,double z_Max );
+// void AddBoundingPlanes(bool yade);
+// void AddBoundingPlanes();
+// void AddBoundingPlanes(Real center[3], Real Extents[3], int id);
+
+ Tesselation& Compute_Action ( );
+ Tesselation& Compute_Action ( int argc, char *argv[ ], char *envp[ ] );
+ Tesselation& LoadPositions(int argc, char *argv[ ], char *envp[ ]);
+// Vecteur external_force_single_fictious ( Cell_handle cell );
+ void SpheresFileCreator ();
+// void Analytical_Consolidation ( );
+ void DisplayStatistics();
+// void Boundary_Conditions ( RTriangulation& Tri );
+ void Initialize_pressures ( double P_zero );
+ /// Define forces using the same averaging volumes as for permeability
+ void ComputeTetrahedralForces();
+ /// Define forces spliting drag and buoyancy terms
+ void ComputeFacetForces();
+ void ComputeFacetForcesWithCache();
+ void save_vtk_file ( );
+ void MGPost ( );
+#ifdef XVIEW
+ void Dessine_Triangulation ( Vue3D &Vue, RTriangulation &T );
+ void Dessine_Short_Tesselation ( Vue3D &Vue, Tesselation &Tes );
+#endif
+ double Permeameter ( double P_Inf, double P_Sup, double Section, double DeltaY, char *file );
+ double Sample_Permeability( double& x_Min,double& x_Max ,double& y_Min,double& y_Max,double& z_Min,double& z_Max, string key);
+ double Compute_HydraulicRadius (Cell_handle cell, int j );
+ double PressureProfile ( char *filename, Real& time, int& intervals );
+
+ double dotProduct ( Vecteur x, Vecteur y );
+ double Compute_EffectiveRadius(Cell_handle cell, int j);
+ double Compute_EquivalentRadius(Cell_handle cell, int j);
+
+// double crossProduct ( double x[3], double y[3] );
+
+// double surface_solid_facet ( Sphere ST1, Sphere ST2, Sphere ST3 );
+// Vecteur surface_double_fictious_facet ( Vertex_handle fSV1, Vertex_handle fSV2, Vertex_handle SV3 );
+// Vecteur surface_single_fictious_facet ( Vertex_handle fSV1, Vertex_handle SV2, Vertex_handle SV3 );
+
+// double surface_solid_double_fictious_facet ( Vertex_handle ST1, Vertex_handle ST2, Vertex_handle ST3 );
+
+// double surface_external_triple_fictious (Cell_handle cell, Boundary b );
+// double surface_external_triple_fictious ( Real position[3], Cell_handle cell, Boundary b );
+// double surface_external_double_fictious ( Cell_handle cell, Boundary b );
+
+// double surface_external_single_fictious ( Cell_handle cell, Boundary b );
+
+ void GenerateVoxelFile ( );
+
+ RTriangulation& Build_Triangulation ( Real x, Real y, Real z, Real radius, unsigned const id );
+
+ void Build_Tessalation ( RTriangulation& Tri );
+
+// double spherical_triangle_area ( Sphere STA1, Sphere STA2, Sphere STA3, Point PTA1 );
+
+// double fast_spherical_triangle_area ( const Sphere& STA1, const Point& STA2, const Point& STA3, const Point& PTA1 );
+// Real solid_angle ( const Point& STA1, const Point& STA2, const Point& STA3, const Point& PTA1 );
+// double spherical_triangle_volume ( const Sphere& ST1, const Point& PT1, const Point& PT2, const Point& PT3 );
+// Real fast_solid_angle ( const Point& STA1, const Point& PTA1, const Point& PTA2, const Point& PTA3 );
+
+ bool isInsideSphere ( double& x, double& y, double& z );
+
+ void SliceField (char *filename);
+ void ComsolField();
+
+ void Interpolate ( Tesselation& Tes, Tesselation& NewTes );
+
+// double volume_single_fictious_pore ( Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1 );
+ //Fast version, assign surface of facet for future forces calculations (pointing from PV2 to PV1)
+// double volume_single_fictious_pore ( const Vertex_handle& SV1, const Vertex_handle& SV2, const Vertex_handle& SV3, const Point& PV1, const Point& PV2, Vecteur& facetSurface);
+// double volume_double_fictious_pore ( Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1 );
+ //Fast version, assign surface of facet for future forces calculations (pointing from PV2 to PV1)
+// double volume_double_fictious_pore (Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point& PV1, Point& PV2, Vecteur& facetSurface);
+
+// double PoreVolume (RTriangulation& Tri, Cell_handle cell);
+ int Average_Cell_Velocity(int id_sphere, RTriangulation& Tri);
+ void Average_Cell_Velocity();
+ void Average_Grain_Velocity();
+ void vtk_average_cell_velocity(RTriangulation &T, int id_sphere, int num_cells);
+ void ApplySinusoidalPressure(RTriangulation& Tri, double Pressure, double load_intervals);
+ void ApplySinusoidalPressure_Space_Time(RTriangulation& Tri, double Pressure, double load_intervals, double time, double dt);
+// double surface_external_triple_fictious(Real position[3], Cell_handle cell, Boundary b);
+// double surface_external_double_fictious(Cell_handle cell, Boundary b);
+// double surface_external_single_fictious(Cell_handle cell, Boundary b);
+
+ ///TAUCS
+ int SetLinearSystem();
+ int SetLinearSystemFullGS();
+ void VectorizedGaussSeidel ();
+ int TaucsSolveTest();
+ int PardisoSolveTest();
+ ///END_TAUCS
+
+};
+
+} //namespace CGT
+#endif //FLOW_ENGINE
+
+#endif
=== modified file 'lib/triangulation/Network.cpp'
--- lib/triangulation/Network.cpp 2010-11-17 15:19:09 +0000
+++ lib/triangulation/Network.cpp 2010-11-17 20:07:21 +0000
@@ -13,7 +13,7 @@
#include <sys/stat.h>
#include <sys/types.h>
-#include "Network.h"
+#include "Network.hpp"
#define FAST
@@ -565,4 +565,4 @@
// }
} //namespace CGT
-#endif //FLOW_ENGINE
\ No newline at end of file
+#endif //FLOW_ENGINE
=== removed file 'lib/triangulation/Network.h'
--- lib/triangulation/Network.h 2010-11-17 15:19:09 +0000
+++ lib/triangulation/Network.h 1970-01-01 00:00:00 +0000
@@ -1,89 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2010 by Emanuele Catalano <catalano@xxxxxxxxxxxxxxx> *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-#ifdef FLOW_ENGINE
-
-#ifndef _NETWORK_H
-#define _NETWORK_H
-
-#include "Operations.h"
-#include "Timer.h"
-#include "Tesselation.h"
-#include "basicVTKwritter.hpp"
-#include "Timer.h"
-#include "stdafx.h"
-#include "Empilement.h"
-
-
-
-namespace CGT{
-
-struct Boundary
-{
- Point p;
- Vecteur normal;
- int coordinate;
- bool flowCondition;//flowCondition=0, pressure is imposed // flowCondition=1, flow is imposed
- Real value;
-};
-
-class Network
-{
- public:
- virtual ~Network();
- Network();
-
- Tesselation T [2];
- bool currentTes;
- double x_min, x_max, y_min, y_max, z_min, z_max, Rmoy, SectionArea, Height, Vtotale;
- bool DEBUG_OUT;
- int nOfSpheres;
- int x_min_id, x_max_id, y_min_id, y_max_id, z_min_id, z_max_id;
- int* boundsIds [6];
- Point Corner_min;
- Point Corner_max;
- Real Vsolid_tot, Vtotalissimo, Vporale, Ssolid_tot;
- Boundary boundaries [6];
- Boundary& boundary (int b) {return boundaries[b-id_offset];}
- short id_offset;
- int vtk_infinite_vertices, vtk_infinite_cells, num_particles;
-
-// int F1, F2, Re1, Re2; //values between 0..3, refers to one cell's fictious(F)/real(Re) vertices
-// int facetRe1, facetRe2, facetRe3, facetF1, facetF2; //indices relative to the facet
- int fictious_vertex;
-// bool facet_detected;
-
-// void DisplayStatistics();
- void AddBoundingPlanes(bool yade);
- void AddBoundingPlanes();
- void Define_fictious_cells( );
- int Detect_facet_fictious_vertices (Cell_handle& cell, int& j);
- double Volume_Pore (Cell_handle cell);
- double Volume_Pore_VoronoiFraction ( Cell_handle& cell, int& j);
- double volume_single_fictious_pore(const Vertex_handle& SV1, const Vertex_handle& SV2, const Vertex_handle& SV3, const Point& PV1, const Point& PV2, Vecteur& facetSurface);
- double volume_double_fictious_pore(const Vertex_handle& SV1, const Vertex_handle& SV2, const Vertex_handle& SV3, const Point& PV1, const Point& PV2, Vecteur& facetSurface);
- double spherical_triangle_volume(const Sphere& ST1, const Point& PT1, const Point& PT2, const Point& PT3);
-
- Vecteur surface_single_fictious_facet(Vertex_handle fSV1, Vertex_handle SV2, Vertex_handle SV3);
-
- double fast_spherical_triangle_area(const Sphere& STA1, const Point& STA2, const Point& STA3, const Point& PTA1);
- Real fast_solid_angle(const Point& STA1, const Point& PTA1, const Point& PTA2, const Point& PTA3);
- double volume_double_fictious_pore(Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1);
- double volume_single_fictious_pore(Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1);
- double Surface_Solid_Pore( Cell_handle cell, int j, bool SLIP_ON_LATERALS);
- double spherical_triangle_area ( Sphere STA1, Sphere STA2, Sphere STA3, Point PTA1 );
-
- Vecteur surface_double_fictious_facet(Vertex_handle fSV1, Vertex_handle fSV2, Vertex_handle SV3);
-// Vecteur surface_single_fictious_facet(Vertex_handle fSV1, Vertex_handle SV2, Vertex_handle SV3);
- double surface_solid_double_fictious_facet(Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3);
- double surface_solid_facet(Sphere ST1, Sphere ST2, Sphere ST3);
-};
-
-} //namespaceCGT
-
-#endif
-
-#endif //FLOW_ENGINE
\ No newline at end of file
=== added file 'lib/triangulation/Network.hpp'
--- lib/triangulation/Network.hpp 1970-01-01 00:00:00 +0000
+++ lib/triangulation/Network.hpp 2010-11-17 20:07:21 +0000
@@ -0,0 +1,89 @@
+/*************************************************************************
+* Copyright (C) 2010 by Emanuele Catalano <catalano@xxxxxxxxxxxxxxx> *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+#ifdef FLOW_ENGINE
+
+#ifndef _NETWORK_H
+#define _NETWORK_H
+
+#include "Operations.h"
+#include "Timer.h"
+#include "Tesselation.h"
+#include "basicVTKwritter.hpp"
+#include "Timer.h"
+#include "stdafx.h"
+#include "Empilement.h"
+
+
+
+namespace CGT{
+
+struct Boundary
+{
+ Point p;
+ Vecteur normal;
+ int coordinate;
+ bool flowCondition;//flowCondition=0, pressure is imposed // flowCondition=1, flow is imposed
+ Real value;
+};
+
+class Network
+{
+ public:
+ virtual ~Network();
+ Network();
+
+ Tesselation T [2];
+ bool currentTes;
+ double x_min, x_max, y_min, y_max, z_min, z_max, Rmoy, SectionArea, Height, Vtotale;
+ bool DEBUG_OUT;
+ int nOfSpheres;
+ int x_min_id, x_max_id, y_min_id, y_max_id, z_min_id, z_max_id;
+ int* boundsIds [6];
+ Point Corner_min;
+ Point Corner_max;
+ Real Vsolid_tot, Vtotalissimo, Vporale, Ssolid_tot;
+ Boundary boundaries [6];
+ Boundary& boundary (int b) {return boundaries[b-id_offset];}
+ short id_offset;
+ int vtk_infinite_vertices, vtk_infinite_cells, num_particles;
+
+// int F1, F2, Re1, Re2; //values between 0..3, refers to one cell's fictious(F)/real(Re) vertices
+// int facetRe1, facetRe2, facetRe3, facetF1, facetF2; //indices relative to the facet
+ int fictious_vertex;
+// bool facet_detected;
+
+// void DisplayStatistics();
+ void AddBoundingPlanes(bool yade);
+ void AddBoundingPlanes();
+ void Define_fictious_cells( );
+ int Detect_facet_fictious_vertices (Cell_handle& cell, int& j);
+ double Volume_Pore (Cell_handle cell);
+ double Volume_Pore_VoronoiFraction ( Cell_handle& cell, int& j);
+ double volume_single_fictious_pore(const Vertex_handle& SV1, const Vertex_handle& SV2, const Vertex_handle& SV3, const Point& PV1, const Point& PV2, Vecteur& facetSurface);
+ double volume_double_fictious_pore(const Vertex_handle& SV1, const Vertex_handle& SV2, const Vertex_handle& SV3, const Point& PV1, const Point& PV2, Vecteur& facetSurface);
+ double spherical_triangle_volume(const Sphere& ST1, const Point& PT1, const Point& PT2, const Point& PT3);
+
+ Vecteur surface_single_fictious_facet(Vertex_handle fSV1, Vertex_handle SV2, Vertex_handle SV3);
+
+ double fast_spherical_triangle_area(const Sphere& STA1, const Point& STA2, const Point& STA3, const Point& PTA1);
+ Real fast_solid_angle(const Point& STA1, const Point& PTA1, const Point& PTA2, const Point& PTA3);
+ double volume_double_fictious_pore(Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1);
+ double volume_single_fictious_pore(Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3, Point PV1);
+ double Surface_Solid_Pore( Cell_handle cell, int j, bool SLIP_ON_LATERALS);
+ double spherical_triangle_area ( Sphere STA1, Sphere STA2, Sphere STA3, Point PTA1 );
+
+ Vecteur surface_double_fictious_facet(Vertex_handle fSV1, Vertex_handle fSV2, Vertex_handle SV3);
+// Vecteur surface_single_fictious_facet(Vertex_handle fSV1, Vertex_handle SV2, Vertex_handle SV3);
+ double surface_solid_double_fictious_facet(Vertex_handle SV1, Vertex_handle SV2, Vertex_handle SV3);
+ double surface_solid_facet(Sphere ST1, Sphere ST2, Sphere ST3);
+};
+
+} //namespaceCGT
+
+#endif
+
+#endif //FLOW_ENGINE
\ No newline at end of file
=== modified file 'pkg/dem/FlowEngine.cpp'
--- pkg/dem/FlowEngine.cpp 2010-11-17 15:19:09 +0000
+++ pkg/dem/FlowEngine.cpp 2010-11-17 20:07:21 +0000
@@ -50,10 +50,9 @@
if ( first ) Build_Triangulation( P_zero );
timingDeltas->checkpoint("Triangulating");
- eps_vol_max=0.f;
- std::ofstream eps_vol ("Epsilon_volume.txt", std::ios::app);
+
+ eps_vol_max=0.f;
UpdateVolumes ( );
- eps_vol << eps_vol_max << endl;
Eps_Vol_Cumulative += eps_vol_max;
if (Eps_Vol_Cumulative > ReTrg*EpsVolPercent_RTRG) {Update_Triangulation = true; ReTrg++;}
=== modified file 'pkg/dem/FlowEngine.hpp'
--- pkg/dem/FlowEngine.hpp 2010-11-17 15:19:09 +0000
+++ pkg/dem/FlowEngine.hpp 2010-11-17 20:07:21 +0000
@@ -9,7 +9,7 @@
#include<yade/core/PartialEngine.hpp>
#include<yade/pkg/dem/TriaxialCompressionEngine.hpp>
-#include<yade/lib/triangulation/FlowBoundingSphere.h>
+#include<yade/lib/triangulation/FlowBoundingSphere.hpp>
#include<yade/pkg/dem/TesselationWrapper.hpp>
class TesselationWrapper;
@@ -32,7 +32,6 @@
bool currentTes;
int id_offset;
// double IS;
- double eps_vol_max;
double Eps_Vol_Cumulative;
int ReTrg;
void Triangulate ();
@@ -71,7 +70,8 @@
((double,Tolerance,1e-06,,"Gauss-Seidel Tolerance"))
((double,Relax,1.9,,"Gauss-Seidel relaxation"))
((int,PermuteInterval,100000,,"Pore space re-triangulation period"))
- ((double,EpsVolPercent_RTRG,0.01,,"Percentuage of cumulate eps_vol at which retriangulation of pore space is performed"))
+ ((double, eps_vol_max, 0,,"Maximal absolute volumetric strain computed at each iteration"))
+ ((double, EpsVolPercent_RTRG,0.01,,"Percentuage of cumulate eps_vol at which retriangulation of pore space is performed"))
((bool,compute_K,true,,"Activates permeability measure within a granular sample"))
((bool,meanK_correction,true,,"Local permeabilities' correction through meanK threshold"))
((bool,meanK_opt,false,,"Local permeabilities' correction through an optimized threshold"))