AFEPack
|
00001 00011 #ifndef __MovingMesh3D_h__ 00012 #define __MovingMesh3D_h__ 00013 00014 #include <lac/sparsity_pattern.h> 00015 #include <lac/sparse_matrix.h> 00016 #include "AMGSolver.h" 00017 #include "Geometry.h" 00018 00029 class MovingMesh3D : public Mesh<3,3> 00030 { 00031 public: 00032 typedef GeometryBM::bmark_t bound_t; 00033 struct Vertex : public afepack::Point<3> { 00034 int index; 00035 bound_t boundary_mark; 00036 }; 00037 struct Edge { 00038 int index; 00039 int vertex[2]; 00040 bound_t boundary_mark; 00041 }; 00042 struct Surface { 00043 int index; 00044 u_int n_edge; 00045 std::vector<int> edge; 00046 bound_t boundary_mark; 00047 double normal[3]; 00048 double logic_normal[3]; 00049 }; 00050 struct Domain { 00051 u_int n_vertex; 00052 u_int n_edge; 00053 u_int n_surf; 00054 std::vector<Vertex> physical_domain_vertex; 00055 std::vector<Vertex> logical_domain_vertex; 00056 std::vector<Edge> edge; 00057 std::vector<Surface> surface; 00058 }; 00059 00060 class Solver { 00061 protected: 00062 bool is_initialized; 00063 int n_project; 00064 std::vector<SparseMatrix<double> *> project_matrix; 00065 std::vector<SparseMatrix<double> *> project_matrix_r; 00066 std::vector<SparseMatrix<double> *> projected_matrix; 00067 std::vector<const std::vector<int> *> boundary_mark; 00068 u_int smooth_step; 00069 u_int min_order; 00070 const Domain * domain; 00071 public: 00072 Solver(); 00073 ~Solver(); 00074 public: 00075 u_int smoothStep() const {return smooth_step;} 00076 u_int& smoothStep() {return smooth_step;} 00077 u_int minimalOrder() const {return min_order;} 00078 u_int& minimalOrder() {return min_order;} 00079 00080 void reinit(SparseMatrix<double>&, 00081 const std::vector<int>&, 00082 const Domain&); 00083 void clear(); 00084 void solve(std::vector<Vector<double> >& x, 00085 const std::vector<Vector<double> >& r, 00086 u_int steps = 5) const; 00087 private: 00088 void Project(const SparseMatrix<double>& M, 00089 const std::vector<int>& bm, 00090 SparseMatrix<double> *& P, 00091 SparseMatrix<double> *& PMPT, 00092 SparseMatrix<double> *& Pt, 00093 std::vector<int> *& pbm); 00094 void Init(SparseMatrix<double>&, 00095 const std::vector<int>&, 00096 const Domain&); 00097 void GaussSidel(const SparseMatrix<double>& M, 00098 std::vector<Vector<double> >& x, 00099 const std::vector<Vector<double> >& r, 00100 const std::vector<int>& bm, 00101 const u_int& s) const; 00102 SparseMatrix<double> * getPMPT(const SparseMatrix<double> & P, 00103 const SparseMatrix<double> & M, 00104 const SparseMatrix<double> & Pt) const; 00105 void lazyPMPT(const SparseMatrix<double>& P, 00106 const SparseMatrix<double>& M, 00107 const SparseMatrix<double>& Pt, 00108 SparseMatrix<double>& A) const; 00109 }; 00110 00111 00112 private: 00113 Domain domain; 00114 std::vector<afepack::Point<3> > logical_node; 00115 double move_step_length; 00116 u_int n_move_step; 00117 std::vector<afepack::Point<3> > move_direction; 00118 std::vector<afepack::Point<3> > logical_move_direction; 00119 std::vector<float> mon; 00120 00121 u_int n_interior_node; 00122 u_int n_boundary_node; 00123 std::vector<int> index; 00124 std::vector<u_int> interior_node_index; 00125 std::vector<u_int> boundary_node_index; 00126 u_int n_boundary_constraint; 00127 std::vector<int> boundary_mark; 00128 00129 SparsityPattern spM; 00130 SparseMatrix<double> M; 00131 Solver solver; 00132 00133 static int primes[]; 00134 u_int solve_step; 00135 double tol; 00136 00137 public: 00138 MovingMesh3D(); 00139 virtual ~MovingMesh3D(); 00140 00141 public: 00142 const std::vector<float>& monitor() const {return mon;}; 00143 std::vector<float>& monitor() {return mon;}; 00144 const float& monitor(const u_int& i) const {return mon[i];}; 00145 float& monitor(const u_int& i) {return mon[i];} 00146 00147 const std::vector<afepack::Point<3> >& moveDirection() const {return move_direction;}; 00148 std::vector<afepack::Point<3> >& moveDirection() {return move_direction;}; 00149 const afepack::Point<3>& moveDirection(const u_int& i) const {return move_direction[i];}; 00150 afepack::Point<3>& moveDirection(const u_int& i) {return move_direction[i];}; 00151 00152 std::vector<double> moveDirection(const afepack::Point<3>& point, 00153 const int& element) const; 00154 std::vector<std::vector<double> > moveDirection(const std::vector<afepack::Point<3> >& point, 00155 const int& element) const; 00156 double moveDirectionDivergence(const u_int& element) const; 00157 00158 double moveStepLength() const {return move_step_length;}; 00159 double& moveStepLength() {return move_step_length;}; 00160 00161 u_int moveStep() const {return n_move_step;}; 00162 u_int& moveStep() {return n_move_step;}; 00163 00164 u_int solveStep() const {return solve_step;} 00165 u_int& solveStep() {return solve_step;} 00166 00167 double tolerence() const {return tol;} 00168 double& tolerence() {return tol;} 00169 00170 void moveMesh(); 00171 void outputPhysicalMesh(const std::string& file); 00172 void outputLogicalMesh(const std::string& file); 00173 virtual void getMonitor(); 00174 virtual void smoothMonitor(u_int step = 1); 00175 virtual void updateMesh(); 00176 virtual void updateSolution() = 0; 00177 virtual void outputSolution() = 0; 00178 virtual void getMoveStepLength(); 00179 void readDomain(const std::string& file); 00180 00181 private: 00182 void getLogicalMesh(); 00183 void getMoveDirection(); 00184 void readDummy(std::ifstream& is); 00185 void parseBoundary(); 00186 }; 00187 00188 #endif