face_quad9.C
Go to the documentation of this file.00001 // The libMesh Finite Element Library. 00002 // Copyright (C) 2002-2012 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner 00003 00004 // This library is free software; you can redistribute it and/or 00005 // modify it under the terms of the GNU Lesser General Public 00006 // License as published by the Free Software Foundation; either 00007 // version 2.1 of the License, or (at your option) any later version. 00008 00009 // This library is distributed in the hope that it will be useful, 00010 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00011 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00012 // Lesser General Public License for more details. 00013 00014 // You should have received a copy of the GNU Lesser General Public 00015 // License along with this library; if not, write to the Free Software 00016 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00017 00018 // C++ includes 00019 00020 // Local includes 00021 #include "libmesh/side.h" 00022 #include "libmesh/edge_edge3.h" 00023 #include "libmesh/face_quad9.h" 00024 00025 namespace libMesh 00026 { 00027 00028 00029 00030 00031 // ------------------------------------------------------------ 00032 // Quad9 class static member initializations 00033 const unsigned int Quad9::side_nodes_map[4][3] = 00034 { 00035 {0, 1, 4}, // Side 0 00036 {1, 2, 5}, // Side 1 00037 {2, 3, 6}, // Side 2 00038 {3, 0, 7} // Side 3 00039 }; 00040 00041 00042 #ifdef LIBMESH_ENABLE_AMR 00043 00044 const float Quad9::_embedding_matrix[4][9][9] = 00045 { 00046 // embedding matrix for child 0 00047 { 00048 // 0 1 2 3 4 5 6 7 8 00049 { 1.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 0 00050 { 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 1 00051 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000 }, // 2 00052 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000 }, // 3 00053 { 0.375000, -0.125000, 0.00000, 0.00000, 0.750000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 4 00054 { 0.00000, 0.00000, 0.00000, 0.00000, 0.375000, 0.00000, -0.125000, 0.00000, 0.750000 }, // 5 00055 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, -0.125000, 0.00000, 0.375000, 0.750000 }, // 6 00056 { 0.375000, 0.00000, 0.00000, -0.125000, 0.00000, 0.00000, 0.00000, 0.750000, 0.00000 }, // 7 00057 { 0.140625, -0.0468750, 0.0156250, -0.0468750, 0.281250, -0.0937500, -0.0937500, 0.281250, 0.562500 } // 8 00058 }, 00059 00060 // embedding matrix for child 1 00061 { 00062 // 0 1 2 3 4 5 6 7 8 00063 { 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 0 00064 { 0.00000, 1.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 1 00065 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000, 0.00000 }, // 2 00066 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000 }, // 3 00067 { -0.125000, 0.375000, 0.00000, 0.00000, 0.750000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 4 00068 { 0.00000, 0.375000, -0.125000, 0.00000, 0.00000, 0.750000, 0.00000, 0.00000, 0.00000 }, // 5 00069 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.375000, 0.00000, -0.125000, 0.750000 }, // 6 00070 { 0.00000, 0.00000, 0.00000, 0.00000, 0.375000, 0.00000, -0.125000, 0.00000, 0.750000 }, // 7 00071 { -0.0468750, 0.140625, -0.0468750, 0.0156250, 0.281250, 0.281250, -0.0937500, -0.0937500, 0.562500 } // 8 00072 }, 00073 00074 // embedding matrix for child 2 00075 { 00076 // 0 1 2 3 4 5 6 7 8 00077 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000 }, // 0 00078 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000 }, // 1 00079 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000 }, // 2 00080 { 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 3 00081 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, -0.125000, 0.00000, 0.375000, 0.750000 }, // 4 00082 { 0.00000, 0.00000, 0.00000, 0.00000, -0.125000, 0.00000, 0.375000, 0.00000, 0.750000 }, // 5 00083 { 0.00000, 0.00000, -0.125000, 0.375000, 0.00000, 0.00000, 0.750000, 0.00000, 0.00000 }, // 6 00084 { -0.125000, 0.00000, 0.00000, 0.375000, 0.00000, 0.00000, 0.00000, 0.750000, 0.00000 }, // 7 00085 { -0.0468750, 0.0156250, -0.0468750, 0.140625, -0.0937500, -0.0937500, 0.281250, 0.281250, 0.562500 } // 8 00086 }, 00087 00088 // embedding matrix for child 3 00089 { 00090 // 0 1 2 3 4 5 6 7 8 00091 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000 }, // 0 00092 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000, 0.00000 }, // 1 00093 { 0.00000, 0.00000, 1.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000 }, // 2 00094 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 1.00000, 0.00000, 0.00000 }, // 3 00095 { 0.00000, 0.00000, 0.00000, 0.00000, 0.00000, 0.375000, 0.00000, -0.125000, 0.750000 }, // 4 00096 { 0.00000, -0.125000, 0.375000, 0.00000, 0.00000, 0.750000, 0.00000, 0.00000, 0.00000 }, // 5 00097 { 0.00000, 0.00000, 0.375000, -0.125000, 0.00000, 0.00000, 0.750000, 0.00000, 0.00000 }, // 6 00098 { 0.00000, 0.00000, 0.00000, 0.00000, -0.125000, 0.00000, 0.375000, 0.00000, 0.750000 }, // 7 00099 { 0.0156250, -0.0468750, 0.140625, -0.0468750, -0.0937500, 0.281250, 0.281250, -0.0937500, 0.562500 } // 8 00100 } 00101 }; 00102 00103 #endif 00104 00105 00106 00107 // ------------------------------------------------------------ 00108 // Quad9 class member functions 00109 00110 bool Quad9::is_vertex(const unsigned int i) const 00111 { 00112 if (i < 4) 00113 return true; 00114 return false; 00115 } 00116 00117 bool Quad9::is_edge(const unsigned int i) const 00118 { 00119 if (i < 4) 00120 return false; 00121 if (i > 7) 00122 return false; 00123 return true; 00124 } 00125 00126 bool Quad9::is_face(const unsigned int i) const 00127 { 00128 if (i > 7) 00129 return true; 00130 return false; 00131 } 00132 00133 bool Quad9::is_node_on_side(const unsigned int n, 00134 const unsigned int s) const 00135 { 00136 libmesh_assert_less (s, n_sides()); 00137 for (unsigned int i = 0; i != 3; ++i) 00138 if (side_nodes_map[s][i] == n) 00139 return true; 00140 return false; 00141 } 00142 00143 00144 00145 bool Quad9::has_affine_map() const 00146 { 00147 // make sure corners form a parallelogram 00148 Point v = this->point(1) - this->point(0); 00149 if (!v.relative_fuzzy_equals(this->point(2) - this->point(3))) 00150 return false; 00151 // make sure "horizontal" sides are straight 00152 v /= 2; 00153 if (!v.relative_fuzzy_equals(this->point(4) - this->point(0)) || 00154 !v.relative_fuzzy_equals(this->point(6) - this->point(3))) 00155 return false; 00156 // make sure "vertical" sides are straight 00157 // and the center node is centered 00158 v = (this->point(3) - this->point(0))/2; 00159 if (!v.relative_fuzzy_equals(this->point(7) - this->point(0)) || 00160 !v.relative_fuzzy_equals(this->point(5) - this->point(1)) || 00161 !v.relative_fuzzy_equals(this->point(8) - this->point(4))) 00162 return false; 00163 return true; 00164 } 00165 00166 00167 00168 dof_id_type Quad9::key (const unsigned int s) const 00169 { 00170 libmesh_assert_less (s, this->n_sides()); 00171 00172 switch (s) 00173 { 00174 case 0: 00175 00176 return 00177 this->compute_key (this->node(4)); 00178 00179 case 1: 00180 00181 return 00182 this->compute_key (this->node(5)); 00183 00184 case 2: 00185 00186 return 00187 this->compute_key (this->node(6)); 00188 00189 case 3: 00190 00191 return 00192 this->compute_key (this->node(7)); 00193 } 00194 00195 00196 // We will never get here... Look at the code above. 00197 libmesh_error(); 00198 return 0; 00199 } 00200 00201 00202 00203 AutoPtr<Elem> Quad9::build_side (const unsigned int i, 00204 bool proxy) const 00205 { 00206 libmesh_assert_less (i, this->n_sides()); 00207 00208 if (proxy) 00209 { 00210 AutoPtr<Elem> ap(new Side<Edge3,Quad9>(this,i)); 00211 return ap; 00212 } 00213 00214 else 00215 { 00216 Edge3* edge = new Edge3; 00217 00218 switch (i) 00219 { 00220 case 0: 00221 { 00222 edge->set_node(0) = this->get_node(0); 00223 edge->set_node(1) = this->get_node(1); 00224 edge->set_node(2) = this->get_node(4); 00225 00226 AutoPtr<Elem> ap(edge); return ap; 00227 } 00228 case 1: 00229 { 00230 edge->set_node(0) = this->get_node(1); 00231 edge->set_node(1) = this->get_node(2); 00232 edge->set_node(2) = this->get_node(5); 00233 00234 AutoPtr<Elem> ap(edge); return ap; 00235 } 00236 case 2: 00237 { 00238 edge->set_node(0) = this->get_node(2); 00239 edge->set_node(1) = this->get_node(3); 00240 edge->set_node(2) = this->get_node(6); 00241 00242 AutoPtr<Elem> ap(edge); return ap; 00243 } 00244 case 3: 00245 { 00246 edge->set_node(0) = this->get_node(3); 00247 edge->set_node(1) = this->get_node(0); 00248 edge->set_node(2) = this->get_node(7); 00249 00250 AutoPtr<Elem> ap(edge); return ap; 00251 } 00252 default: 00253 { 00254 libmesh_error(); 00255 } 00256 } 00257 } 00258 00259 // We will never get here... 00260 AutoPtr<Elem> ap(NULL); return ap; 00261 } 00262 00263 00264 00265 00266 00267 00268 00269 void Quad9::connectivity(const unsigned int sf, 00270 const IOPackage iop, 00271 std::vector<dof_id_type>& conn) const 00272 { 00273 libmesh_assert_less (sf, this->n_sub_elem()); 00274 libmesh_assert_not_equal_to (iop, INVALID_IO_PACKAGE); 00275 00276 conn.resize(4); 00277 00278 switch (iop) 00279 { 00280 case TECPLOT: 00281 { 00282 switch(sf) 00283 { 00284 case 0: 00285 // linear sub-quad 0 00286 conn[0] = this->node(0)+1; 00287 conn[1] = this->node(4)+1; 00288 conn[2] = this->node(8)+1; 00289 conn[3] = this->node(7)+1; 00290 return; 00291 00292 case 1: 00293 // linear sub-quad 1 00294 conn[0] = this->node(4)+1; 00295 conn[1] = this->node(1)+1; 00296 conn[2] = this->node(5)+1; 00297 conn[3] = this->node(8)+1; 00298 return; 00299 00300 case 2: 00301 // linear sub-quad 2 00302 conn[0] = this->node(7)+1; 00303 conn[1] = this->node(8)+1; 00304 conn[2] = this->node(6)+1; 00305 conn[3] = this->node(3)+1; 00306 return; 00307 00308 case 3: 00309 // linear sub-quad 3 00310 conn[0] = this->node(8)+1; 00311 conn[1] = this->node(5)+1; 00312 conn[2] = this->node(2)+1; 00313 conn[3] = this->node(6)+1; 00314 return; 00315 00316 default: 00317 libmesh_error(); 00318 } 00319 } 00320 00321 case VTK: 00322 { 00323 conn.resize(9); 00324 conn[0] = this->node(0); 00325 conn[1] = this->node(1); 00326 conn[2] = this->node(2); 00327 conn[3] = this->node(3); 00328 conn[4] = this->node(4); 00329 conn[5] = this->node(5); 00330 conn[6] = this->node(6); 00331 conn[7] = this->node(7); 00332 conn[8] = this->node(8); 00333 return; 00334 00335 /* 00336 switch(sf) 00337 { 00338 case 0: 00339 // linear sub-quad 0 00340 conn[0] = this->node(0); 00341 conn[1] = this->node(4); 00342 conn[2] = this->node(8); 00343 conn[3] = this->node(7); 00344 00345 return; 00346 00347 case 1: 00348 // linear sub-quad 1 00349 conn[0] = this->node(4); 00350 conn[1] = this->node(1); 00351 conn[2] = this->node(5); 00352 conn[3] = this->node(8); 00353 00354 return; 00355 00356 case 2: 00357 // linear sub-quad 2 00358 conn[0] = this->node(7); 00359 conn[1] = this->node(8); 00360 conn[2] = this->node(6); 00361 conn[3] = this->node(3); 00362 00363 return; 00364 00365 case 3: 00366 // linear sub-quad 3 00367 conn[0] = this->node(8); 00368 conn[1] = this->node(5); 00369 conn[2] = this->node(2); 00370 conn[3] = this->node(6); 00371 00372 return; 00373 00374 default: 00375 libmesh_error(); 00376 }*/ 00377 } 00378 00379 default: 00380 { 00381 libmesh_error(); 00382 } 00383 } 00384 00385 libmesh_error(); 00386 } 00387 00388 00389 00390 00391 unsigned int Quad9::n_second_order_adjacent_vertices (const unsigned int n) const 00392 { 00393 switch (n) 00394 { 00395 case 4: 00396 case 5: 00397 case 6: 00398 case 7: 00399 return 2; 00400 00401 case 8: 00402 return 4; 00403 00404 default: 00405 libmesh_error(); 00406 } 00407 00408 libmesh_error(); 00409 return libMesh::invalid_uint; 00410 } 00411 00412 00413 00414 unsigned short int Quad9::second_order_adjacent_vertex (const unsigned int n, 00415 const unsigned int v) const 00416 { 00417 libmesh_assert_greater_equal (n, this->n_vertices()); 00418 libmesh_assert_less (n, this->n_nodes()); 00419 00420 switch (n) 00421 { 00422 case 8: 00423 { 00424 libmesh_assert_less (v, 4); 00425 return static_cast<unsigned short int>(v); 00426 } 00427 00428 default: 00429 { 00430 libmesh_assert_less (v, 2); 00431 // use the matrix that we inherited from \p Quad 00432 return _second_order_adjacent_vertices[n-this->n_vertices()][v]; 00433 } 00434 } 00435 } 00436 00437 00438 00439 std::pair<unsigned short int, unsigned short int> 00440 Quad9::second_order_child_vertex (const unsigned int n) const 00441 { 00442 libmesh_assert_greater_equal (n, this->n_vertices()); 00443 libmesh_assert_less (n, this->n_nodes()); 00444 /* 00445 * the _second_order_vertex_child_* vectors are 00446 * stored in face_quad.C, since they are identical 00447 * for Quad8 and Quad9 (for the first 4 higher-order nodes) 00448 */ 00449 return std::pair<unsigned short int, unsigned short int> 00450 (_second_order_vertex_child_number[n], 00451 _second_order_vertex_child_index[n]); 00452 } 00453 00454 } // namespace libMesh 00455 00456
Site Created By: libMesh Developers
Last modified: February 05 2013 19:54:46 UTC
Hosted By: