MechSys: fem/debug.h Source File

00001 /*************************************************************************************
00002  * MechSys - A C++ library to simulate (Continuum) Mechanical Systems                *
00003  * Copyright (C) 2005 Dorival de Moraes Pedroso <dorival.pedroso at gmail.com>       *
00004  * Copyright (C) 2005 Raul Dario Durand Farfan  <raul.durand at gmail.com>           *
00005  *                                                                                   *
00006  * This file is part of MechSys.                                                     *
00007  *                                                                                   *
00008  * MechSys is free software; you can redistribute it and/or modify it under the      *
00009  * terms of the GNU General Public License as published by the Free Software         *
00010  * Foundation; either version 2 of the License, or (at your option) any later        *
00011  * version.                                                                          *
00012  *                                                                                   *
00013  * MechSys is distributed in the hope that it will be useful, but WITHOUT ANY        *
00014  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A   *
00015  * PARTICULAR PURPOSE. See the GNU General Public License for more details.          *
00016  *                                                                                   *
00017  * You should have received a copy of the GNU General Public License along with      *
00018  * MechSys; if not, write to the Free Software Foundation, Inc., 51 Franklin Street, *
00019  * Fifth Floor, Boston, MA 02110-1301, USA                                           *
00020  *************************************************************************************/
00021 
00022 #ifndef MECHSYS_FEM_DEBUG_H
00023 #define MECHSYS_FEM_DEBUG_H
00024 
00025 #ifdef HAVE_CONFIG_H
00026   #include "config.h"
00027 #else
00028   #ifndef REAL
00029     #define REAL double
00030   #endif
00031 #endif
00032 
00033 #include "util/numstreams.h"
00034 #include "fem/ele/tet10equilib.h"
00035 #include "linalg/vector.h"
00036 #include "linalg/matrix.h"
00037 
00038 using std::cout;
00039 using std::endl;
00040 using LinAlg::Vector;
00041 using LinAlg::Matrix;
00042 using Util::_6_3;
00043 
00044 namespace FEM
00045 {
00046 
00047 class Debug
00048 {
00049 public:
00050     void TestTet10ShapeAndDerivs();
00051 private:
00052 }; // class Debug
00053 
00054 
00056 
00057 
00058 inline void Debug::TestTet10ShapeAndDerivs() // {{{
00059 {
00061     Array<FEM::Node> N;
00062     N.resize(10);
00063     N[0].Initialize(0,  0, 0, 0);
00064     N[1].Initialize(1, 10, 0, 0);
00065     N[2].Initialize(2,  0,10, 0);
00066     N[3].Initialize(3,  0, 0,10);
00067     N[4].Initialize(4,  5, 0, 0);
00068     N[5].Initialize(5,  5, 5, 0);
00069     N[6].Initialize(6,  0, 5, 0);
00070     N[7].Initialize(7,  0, 0, 5);
00071     N[8].Initialize(8,  5, 0, 5);
00072     N[9].Initialize(9,  0, 5, 5);
00073 
00075     Tet10Equilib E;
00076 
00077     // Setup element connectivity and nodal vars
00078     for (int i=0; i<10; ++i)
00079         E.SetNode(i, &N[i]);
00080 
00081     cout << "\n/////////////////////////////////////////////////////////////////////////////////// Face shape\n";
00082     cout << "         s                 " << endl;
00083     cout << "         ^                 " << endl;
00084     cout << "         |                 " << endl;
00085     cout << "       2 @.                " << endl;
00086     cout << "         | '.              " << endl;
00087     cout << "         |   '. 4          " << endl;
00088     cout << "       5 @     @.          " << endl;
00089     cout << "         |       '.        " << endl;
00090     cout << "         |         '.      " << endl;
00091     cout << "         @-----@-----@-> r " << endl;
00092     cout << "         0     3     1     " << endl;
00093 
00094     Vector<REAL> fshape;
00095     REAL r,s;
00096     r=0;   s=0;   E._face_shape(r,s,fshape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<" : "; for (int i=0; i<6; ++i) { cout << _6_3()<<fshape(i); } cout << endl;
00097     r=1;   s=0;   E._face_shape(r,s,fshape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<" : "; for (int i=0; i<6; ++i) { cout << _6_3()<<fshape(i); } cout << endl;
00098     r=0;   s=1;   E._face_shape(r,s,fshape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<" : "; for (int i=0; i<6; ++i) { cout << _6_3()<<fshape(i); } cout << endl;
00099     r=0.5; s=0;   E._face_shape(r,s,fshape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<" : "; for (int i=0; i<6; ++i) { cout << _6_3()<<fshape(i); } cout << endl;
00100     r=0.5; s=0.5; E._face_shape(r,s,fshape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<" : "; for (int i=0; i<6; ++i) { cout << _6_3()<<fshape(i); } cout << endl;
00101     r=0;   s=0.5; E._face_shape(r,s,fshape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<" : "; for (int i=0; i<6; ++i) { cout << _6_3()<<fshape(i); } cout << endl;
00102 
00103     cout << "\n/////////////////////////////////////////////////////////////////////////////////// Element shape\n";
00104     cout << "                       t                                        " << endl;
00105     cout << "                       |                                        " << endl;
00106     cout << "                       |                                        " << endl;
00107     cout << "                       | 3                                      " << endl;
00108     cout << "                       @,                                       " << endl;
00109     cout << "                      /|`                                       " << endl;
00110     cout << "                      ||  `,                                    " << endl;
00111     cout << "                     / |    ',                                  " << endl;
00112     cout << "                     | |      \\                                " << endl;
00113     cout << "                    /  |       `.                               " << endl;
00114     cout << "                    |  |         `,  9                          " << endl;
00115     cout << "                   /   @ 7         `@                           " << endl;
00116     cout << "                   |   |             \\                         " << endl;
00117     cout << "                  /    |              `.                        " << endl;
00118     cout << "                  |    |                ',                      " << endl;
00119     cout << "               8 @     |                  \\                    " << endl;
00120     cout << "                 |     @.,,_       6       `.                   " << endl;
00121     cout << "                |     / 0   ``'-.,,@_        `.                 " << endl;
00122     cout << "                |    /              ``''-.,,_  ', 2             " << endl;
00123     cout << "               |    /                        ``'@.,,,           " << endl;
00124     cout << "               |   '                       ,.-``     ``''- s    " << endl;
00125     cout << "              |  ,@ 4                 _,-'`                     " << endl;
00126     cout << "              ' /                 ,.'`                          " << endl;
00127     cout << "             | /             _.@``                              " << endl;
00128     cout << "             '/          ,-'`   5                               " << endl;
00129     cout << "            |/      ,.-``                                       " << endl;
00130     cout << "            /  _,-``                                            " << endl;
00131     cout << "          .@ '`                                                 " << endl;
00132     cout << "         / 1                                                    " << endl;
00133     cout << "        /                                                       " << endl;
00134     cout << "       /                                                        " << endl;
00135     cout << "      r                                                         " << endl;
00136 
00137     Vector<REAL> shape;
00138     REAL t;
00139     r=0;   s=0;   t=0;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00140     r=1;   s=0;   t=0;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00141     r=0;   s=1;   t=0;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00142     r=0;   s=0;   t=1;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00143     r=0.5; s=0;   t=0;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00144     r=0.5; s=0.5; t=0;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00145     r=0;   s=0.5; t=0;   E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00146     r=0;   s=0;   t=0.5; E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00147     r=0.5; s=0;   t=0.5; E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00148     r=0;   s=0.5; t=0.5; E.Shape(r,s,t,shape); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : "; for (int i=0; i<10; ++i) { cout << _6_3()<<shape(i); } cout << endl;
00149 
00150     cout << "\n/////////////////////////////////////////////////////////////////////////////////// Derivs and Jacobian \n";
00151     Matrix<REAL> derivs;
00152     Matrix<REAL> J;
00153     r=0;   s=0;   t=0;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00154     r=1;   s=0;   t=0;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00155     r=0;   s=1;   t=0;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00156     r=0;   s=0;   t=1;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00157     r=0.5; s=0;   t=0;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00158     r=0.5; s=0.5; t=0;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00159     r=0;   s=0.5; t=0;   E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00160     r=0;   s=0;   t=0.5; E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00161     r=0.5; s=0;   t=0.5; E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00162     r=0;   s=0.5; t=0.5; E.Derivs(r,s,t,derivs); E.Jacobian(derivs,J); cout<<"r="<<_6_3()<<r<<", s="<<_6_3()<<s<<", t="<<_6_3()<<t<<" : \n"; cout<<"derivs=\n"<<derivs<<endl; cout<<"J=\n"<<J<<endl; cout<<"det_J="<<J.Det()<< endl;
00163 } // }}}
00164 
00165 }; // namespace FEM
00166 
00167 #endif // MECHSYS_FEM_DEBUG_H
00168 
00169 // vim:fdm=marker