pyNastran  0.5.0
pyNastran BDF Reader/Writer, OP2 Parser, and GUI
massMatrix.py
Go to the documentation of this file.
00001 ## GNU Lesser General Public License
00002 ## 
00003 ## Program pyNastran - a python interface to NASTRAN files
00004 ## Copyright (C) 2011-2012  Steven Doyle, Al Danial
00005 ## 
00006 ## Authors and copyright holders of pyNastran
00007 ## Steven Doyle <mesheb82@gmail.com>
00008 ## Al Danial    <al.danial@gmail.com>
00009 ## 
00010 ## This file is part of pyNastran.
00011 ## 
00012 ## pyNastran is free software: you can redistribute it and/or modify
00013 ## it under the terms of the GNU Lesser General Public License as published by
00014 ## the Free Software Foundation, either version 3 of the License, or
00015 ## (at your option) any later version.
00016 ## 
00017 ## pyNastran is distributed in the hope that it will be useful,
00018 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
00019 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00020 ## GNU General Public License for more details.
00021 ## 
00022 ## You should have received a copy of the GNU Lesser General Public License
00023 ## along with pyNastran.  If not, see <http://www.gnu.org/licenses/>.
00024 ## 
00025 from sympy import Symbol, Matrix, integrate
00026 
00027 z = Symbol('z')
00028 n = Symbol('n')
00029 b = Symbol('b')
00030 
00031 p = Symbol('p')
00032 A = Symbol('A')
00033 L = Symbol('L')
00034 t = Symbol('t')
00035 v = Symbol('V')
00036 
00037 def bar():
00038     N1 = (1-z)/2
00039     N2 = (1+z)/2
00040     NT = Matrix([N1,N2])
00041     pdV = p*A*L/2
00042     M = makeM(pdV,NT)
00043     print "Mbar = \n",M
00044 
00045 def truss():
00046     N1 = (1-z)/2
00047     N2 = (1+z)/2
00048     NT = Matrix([[N1,0,N2,0],
00049                  [0,N1,0,N2]])
00050     pdV = p*A*L/2
00051     M = makeM(pdV,NT)
00052     print "Mtruss = \n",M
00053 
00054 def quad():
00055     N1 = (1-z)*(1-n)/4
00056     N2 = (1+z)*(1-n)/4
00057     N3 = (1+z)*(1+n)/4
00058     N4 = (1-z)*(1+n)/4
00059 
00060     N = Matrix([[N1,0,N2,0,N3,0,N4,0],
00061                 [0,N1,0,N2,0,N3,0,N4]])
00062     NT = N.transpose()
00063     pdV = p*A*t/4 # 4 refers to number of nodes??
00064     Jacobian = Matrix([ [0,0,0,0],  # not done
00065                         [0,0,0,0],
00066                         [0,0,0,0],
00067                         [0,0,0,0]])
00068     factorI = Jacobian
00069     M = makeM(pdV,NT,factorI,levels=2)
00070     print "Mquad = \n",M
00071 
00072 def tet4():
00073     N1 = z
00074     N2 = n
00075     N3 = b
00076     N4 = 1-z-n-b
00077     
00078     X  = Matrix([[1,x,y,z]])
00079     X2 = Matrix([[1,x1,y1,z1],
00080                  [1,x2,y2,z2],
00081                  [1,x3,y3,z3],
00082                  [1,x4,y4,z4]])
00083     N = X*X2.inv()
00084     N1 = N[0,0]
00085     N2 = N[0,1]
00086     N3 = N[0,2]
00087     N4 = N[0,3]
00088     
00089     N = Matrix([[N1,0,0,N2,0,0,N3,0,0,N4,0,0],
00090                 [0,N1,0,0,N2,0,0,N3,0,0,N4,0],
00091                 [0,0,N1,0,0,N2,0,0,N3,0,0,N4]])
00092     NT = N.transpose()
00093     #pdV = p*v
00094     pdV = 3
00095     factorI = 1
00096     M = makeM(pdV,NT,factorI,levels=3)
00097     print "Mtet = \n",M
00098 
00099 def makeM(pdV,NT,factorI=1,levels=1):
00100     N = NT.transpose()
00101 
00102     #print "N = \n",N
00103 
00104     print "size(NT) = ",NT.shape
00105     print "size(N) = ",N.shape
00106 
00107     NtN = NT*N
00108     B = []
00109     print "NtN = \n",NtN
00110     print "size(NtN) = ",NtN.shape
00111 
00112     M = integrate(NtN*factorI,z)
00113     Mp1 = M.subs(z,1)
00114     Mm1 = M.subs(z,-1)
00115     M2 = Mp1-Mm1
00116     if levels>=2:
00117         M3 = integrate(M2,n)
00118         M3p1 = M3.subs(n,1)
00119         M3m1 = M3.subs(n,-1)
00120         M4 = M3p1-M3m1
00121         M2 = M4
00122         if levels>=3:
00123             print "M4 = ",M4
00124             M5 = integrate(M4,b)
00125             M5p1 = M5.subs(b,1)
00126             M5m1 = M5.subs(b,-1)
00127             M6 = M5p1-M5m1
00128             M2 = M6
00129             print "M6 = ",M6
00130         ###
00131     ###
00132     print "pdV = ",pdV
00133     MM = pdV*M2
00134     MM.simplify()
00135     return MM
00136 
00137 #bar()
00138 #truss()
00139 tet4()
 All Classes Namespaces Files Functions Variables