pyNastran  0.5.0
pyNastran BDF Reader/Writer, OP2 Parser, and GUI
shell.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 #pylint: disable=C0103
00026 
00027 from __future__ import division, print_function
00028 import sys
00029 import copy
00030 from math import sin, cos, radians
00031 from itertools import izip
00032 
00033 from numpy import zeros, transpose, array
00034 
00035 from pyNastran.bdf.fieldWriter import set_blank_if_default
00036 from pyNastran.bdf.cards.baseCard import Property, Material
00037 
00038 class ShellProperty(Property):
00039     type = 'ShellProperty'
00040     def __init__(self, card=None, data=None):
00041         Property.__init__(self, card, data)
00042 
00043     def S(self):
00044         """"
00045         Calculates the compliance matrix for a lamina
00046         \f[ \large [Q] = [S]^{-1}  \f]
00047         @todo finish...if necessary...
00048         """
00049         raise NotImplementedError()
00050         return Q.inv()
00051 
00052     def ABDH(self):
00053         """
00054         tranforms load to strain/bending curvature taken at \f$ z=0 \f$
00055         \f[ \large  \left[ 
00056           \begin{array}{c}
00057               Nx  \\
00058               Ny  \\
00059               Nz  \\
00060               Mx  \\
00061               My  \\
00062               Mz  \\
00063           \end{array} \right] = 
00064 
00065           \left[ 
00066           \begin{array}{cc}
00067               A & B  \\
00068               B & D  \\
00069           \end{array} \right]
00070           
00071           \left[ 
00072           \begin{array}{c}
00073               \epsilon_{xx}  \\
00074               \epsilon_{yy}  \\
00075                 \gamma_{xy}  \\
00076                 \kappa_{xx}  \\
00077                 \kappa_{yy}  \\
00078                 \kappa_{xy}  \\
00079           \end{array} \right]
00080         \f] 
00081 
00082 
00083         [Nx] = [            ] [ e_xx0    ]
00084         [Ny] = [  [A]   [B] ] [ e_yy0    ]
00085         [Nz] = [            ] [ gamma_xy0]
00086         [Mx] = [            ] [ k_xx0    ]
00087         [My] = [  [B]   [D] ] [ k_yy0    ]
00088         [Mz] = [            ] [ k_xy0    ]
00089         
00090         \f[ \large  A_{ij} = \Sigma_{k=1}^N (\overline{Q_{ij}})_k \left( z_k  -z_{k-1}   \right) = \Sigma_{k=1}^N (Q_{ij})_k t_k                                                                                    \f]
00091         \f[ \large  B_{ij} = \Sigma_{k=1}^N (\overline{Q_{ij}})_k \left( z_k^2-z_{k-1}^2 \right) = \Sigma_{k=1}^N (Q_{ij})_k                           \left( \overline{z} t_k                      \right)         \f]
00092         \f[ \large  D_{ij} = \Sigma_{k=1}^N (\overline{Q_{ij}})_k \left( z_k^3-z_{k-1}^3 \right) = \Sigma_{k=1}^N (Q_{ij})_k                           \left( \overline{z}^2 t_k + \frac{t_k^3}{12} \right)         \f]
00093         \f[ \large  H_{ij} =                                                                       \Sigma_{k=1}^N (Q_{ij})_k \left( t_k -\frac{4}{t^2} \left( \overline{z}^2 t_k + \frac{t_k^3}{12} \right) \right) \f]
00094         
00095         p. 138 of "Introduction to Composite Material Design"
00096         """
00097         raise NotImplementedError()
00098         A = zeros(9,'d')
00099         B = copy.deepcopy(A)
00100         D = copy.deepcopy(A)
00101         H = copy.deepcopy(A)
00102 
00103         for (i, layer) in enumerate(self.layers):
00104             t = layer.t
00105             z0 = layer.z
00106             #z1 = z0+t
00107             zbar = (2*z0+t)/2. # (z1+z0)/2. 
00108             Qraw   = layer.Q() # needs E11, E22, G12, G13, nu12, nu21
00109             qlayer = self.Qall(Qout, layer.thetad)
00110             A += qlayer*t
00111             B += qlayer*t*z
00112             
00113             Dfactor = t*zbar*zbar+t**3/12.
00114             D += qlayer*Dfactor
00115             H += qlayer*(t-4./t**2 * Dfactor)
00116         B = 0.5*B
00117         D = D/3.
00118         H = H*5./4.
00119             
00120             
00121     def Qall(self, thetad):
00122         """
00123         Caculates the laminate tranformation  stiffness \f$ [Q]_{all} \f$
00124         \f[ \large  [Q]_{all} = [T]^{-1} [Q] [R][T][R]^{-1}  \f]
00125         \f[ \large  [Q]_{all} = [T]^{-1} [Q] [T]^{-T}        \f]
00126         
00127         p. 123 of "Introduction to Composite Material Design"
00128         """
00129         theta = radians(thetad)
00130         ct  = cos(theta)
00131         c2t = ct*ct
00132         c3t = ct*c2t
00133         c4t = c2t*c2t
00134 
00135         st  = sin(theta)
00136         s2t = st*st
00137         s3t = st*s2t
00138         s4t = s2t*s2t
00139         
00140         s2c2t = s2t*c2t
00141         #s4tpc4t = s4t+c4t
00142         
00143         Q11a = Q11*c4t + 2*(Q12+2*Q66)*s2c2t + Q22*s4t
00144         Q12a = (Q11+Q22-4*Q66)*s2c2t + Q12(s4t+c4t)
00145         Q22a = Q11*s4t + 2*(Q12+2*Q66)*s2c2t+Q22*c4t
00146         Q16a = (Q11-Q12-2*Q66)*st*c3t + (Q12-Q22+2*Q66)*s3t*ct
00147         Q26a = (Q11-Q12-2*Q66)*s3t*ct + (Q12-Q22+2*Q66)*st*c3t
00148         Q66a = (Q11+Q22-2*Q12-2*Q66)*s2c2t + Q66(s4t+c4t)
00149         Q44a = Q44*c2t + Q55*s2t
00150         Q55a = Q44*s2t + Q55*c2t
00151         Q45a = (Q55-Q44)*st*ct
00152         return array([Q11a, Q12a, Q22a, Q16a, Q26a, Q66a, Q44a, Q55a, Q45a])
00153 
00154     def Q(self):
00155         """"
00156         Calculates the stiffness matrix \f$ [Q] \f$ for a lamina
00157         @todo is this done?
00158         p. 114 "Introduction to Composite Material Design"
00159         """
00160         nu12 = self.nu12
00161         E1 = self.E1()
00162         E2 = self.E2()
00163         delta = 1-nu12*nu21
00164         Q11 = E1/delta
00165         Q12 = nu12*E2/delta
00166         Q22 = E2/delta
00167         Q66 = G12
00168         Q44 = G23
00169         Q55 = G13
00170         Qout = (Q11, Q22, Q12, Q44, Q55, Q66)
00171         return Qout
00172         
00173     def T(self, theta):
00174         """
00175         calculates the Transformation matricies \$ [T] \$ and  \f$ [T]^{-1} \f$
00176         @param self           the object pointer
00177         @param theta          in radians...
00178         @retval Tinv          the inverse transformation matrix
00179         @retval TinvTranspose the transposed inverse transformation matrix
00180         @todo document better
00181 
00182         tranformation matrix  \f$ [T] \f$
00183         \f[ \large  [T] = \left[ 
00184           \begin{array}{ccc}
00185               m^2 & n^2 &  2mn    \\
00186               n^2 & m^2 & -2mn    \\
00187               -mn & mn  & m^2-n^2 
00188           \end{array} \right]
00189         \f] 
00190 
00191                  [ m^2  n^2        2mn]
00192         [T]    = [ n^2  m^2       -2mn]   # transformation matrix
00193                  [ -mn   mn    m^2-n^2]
00194 
00195 
00196         inverse transformation matrix \f$ [T]^{-1} \f$
00197         \f[ \large  [T]^{-1} = \left[ 
00198           \begin{array}{ccc}
00199               m^2 & n^2 & -2mn    \\
00200               n^2 & m^2 &  2mn    \\
00201               mn  & -mn & m^2-n^2 
00202           \end{array} \right]
00203         \f] 
00204 
00205                  [ m^2  n^2       -2mn]
00206         [T]^-1 = [ n^2  m^2        2mn]   # inverse transform
00207                  [ mn   -mn    m^2-n^2]
00208 
00209         \f[ \large  \left[ 
00210             \begin{array}{c}
00211                  \sigma_{xx}   \\
00212                  \sigma_{yy}   \\
00213                  \sigma_{xy} 
00214             \end{array} \right]
00215 
00216              = [T]^{-1} [Q] [R][T]
00217           
00218             \left[ \begin{array}{c}
00219                  \epsilon_{xx} \\
00220                  \epsilon_{yy} \\
00221                  \frac{1}{2} \gamma_{xy}
00222             \end{array} \right]
00223         \f] 
00224         
00225         p.119 "Introduction to Composite Material Design"
00226         """
00227         n = cos(theta)
00228         m = sin(theta)
00229         Tinv  = zeros((6, 6), 'd')
00230         mm = m*m
00231         nn = n*n
00232         mn = m*n
00233         Tinv[0][0] = Tinv[1][1] = mm
00234         Tinv[1][0] = Tinv[0][1] = nn
00235         Tinv[0][2] = -2*mn
00236         Tinv[1][2] =  2*mn
00237 
00238         Tinv[2][0] =  mn
00239         Tinv[2][1] = -mn
00240         Tinv[2][2] =  mm-nn
00241         TinvT = transpose(Tinv)
00242         return (Tinv, TinvT)
00243 
00244 class PCOMP(ShellProperty):
00245     """
00246     PCOMP     701512   0.0+0 1.549-2                   0.0+0   0.0+0     SYM
00247               300704   3.7-2   0.0+0     YES  300704   3.7-2     45.     YES
00248               300704   3.7-2    -45.     YES  300704   3.7-2     90.     YES
00249               300705      .5   0.0+0     YES
00250     """
00251     type = 'PCOMP'
00252     def __init__(self, card=None, data=None): # not done, cleanup
00253         ShellProperty.__init__(self, card, data)
00254         
00255         if card:
00256             #fields = card.fields(1)
00257             ## Property ID
00258             self.pid  = card.field(1)
00259             ## Non-Structural Mass
00260             self.nsm  = card.field(3, 0.0)
00261             self.sb   = card.field(4, 0.0)
00262             self.ft   = card.field(5)
00263             self.TRef = card.field(6, 0.0)
00264             self.ge   = card.field(7, 0.0)
00265 
00266             ## symmetric flag - default = No Symmetry (NO)
00267             self.lam  = card.field(8)
00268             #print "lam = ",self.lam
00269 
00270             nPlyFields = card.nFields()-9 # -8 for the first 8 fields (1st line)
00271             #plyCards = card.fields(9)
00272 
00273             # counting plies
00274             nMajor    = nPlyFields//4
00275             nLeftover = nPlyFields % 4
00276             if nLeftover:
00277                 nMajor += 1
00278             nplies = nMajor
00279             #print "nplies = ",nplies
00280 
00281             iPly = 1
00282             plies = []
00283             midLast = None
00284             tLast = None
00285 
00286             ## supports single ply per line
00287             for i in xrange(9, 9+nplies*4, 4):
00288                 defaults = [midLast, tLast, 0.0, 'NO']
00289                 (mid, t, theta, sout) = card.fields(i, i+4, defaults)
00290                 ply = [mid, t, theta, sout]
00291                 
00292                 assert t > 0., 'thickness of PCOMP layer is invalid iLayer=%s t=%s ply=[mid,t,theta,sout]=%s' % (iPly, t, ply)
00293                 if ply != defaults: # if they're not all defaults...
00294                     plies.append(ply)
00295                 midLast = mid
00296                 tLast = t
00297                 iPly += 1
00298             #print "nplies = ",nplies
00299 
00300             ## list of plies
00301             self.plies = plies
00302 
00303             #self.plies = []
00304             #if self.lam == 'SYM':
00305             #    if nplies%2 == 1:  # 0th layer is the core layer
00306             #        plies[0][1] = plies[0][1]/2. # cut the thickness in half to make the ply have an even number of plies, is there a better way???
00307             #
00308             #    pliesLower = plies.reverse()
00309             #    self.plies = pliesLower+plies
00310             #    #print str(self)
00311             ###
00312             self.z0 = card.field(2, -0.5*self.Thickness())
00313         else:
00314             #print "len(data) = ",len(data)
00315             self.pid  = data[0]
00316             self.z0   = data[1]
00317             self.nsm  = data[2]
00318             self.sb   = data[3]
00319             self.ft   = data[4]
00320             self.TRef = data[5]
00321             self.ge   = data[6]
00322             self.lam  = data[7]
00323             Mid       = data[8]
00324             T         = data[9]
00325             Theta     = data[10]
00326             Sout      = data[11]
00327 
00328             self.plies = []
00329             #ply = [mid,t,theta,sout]
00330             for (mid, t, theta, sout) in izip(Mid, T, Theta, Sout):
00331                 if sout == 1:  ## @todo not sure  0=NO,1=YES (most likely)
00332                     sout = 'YES'
00333                 elif sout == 0:
00334                     sout = 'NO'
00335                 else:
00336                     raise RuntimeError('unsupported sout...needs debugging...'
00337                                        '\nPCOMP = %s' %(data))
00338                 self.plies.append([mid, t, theta, sout])
00339             ###
00340         ###
00341         #print self
00342         #print "nPlies = %s" %(self.nPlies())
00343 
00344     def hasCoreLayer(self):
00345         """is there a center layer (matters most for a symmetrical ply)"""
00346         return self.nPlies()%2==1 # True if has a core, False otherwise
00347 
00348     def nPlies(self):
00349         """
00350         how many plies are there?
00351         @note 
00352             returns half the number of actual plies if the ply is symmetrical (not considering the core)
00353         """
00354         return len(self.plies)
00355 
00356     def isSymmetrical(self):
00357         """is the laminate symmetrical laminate"""
00358         if self.lam == 'SYM':
00359             return True
00360         return False
00361 
00362     def isSameCard(self, prop, debug=False):
00363         if self.type != self.prop.type:
00364             return False
00365         fields2 = [prop.nsm, prop.sb, prop.ft, prop.TRef, prop.ge, prop.lam]
00366         fields1 = [self.nsm, self.sb, self.ft, self.TRef, self.ge, self.lam]
00367 
00368         for ply in self.plies:
00369             fields1 += ply
00370         for ply in prop.plies:
00371             fields2 += ply
00372 
00373         if debug:
00374             print("fields1=%s fields2=%s" % (fields1, fields2))
00375 
00376         for (field1, field2) in izip(fields1, fields2):
00377             if not self.isSame(field1, field2):
00378                 return False
00379             ###
00380         ###
00381         return True
00382 
00383     def crossReference(self, model):
00384         """
00385         links the material ID to the materials
00386         @param self the object pointer
00387         @param model a BDF object
00388         """
00389         for iPly in xrange(len(self.plies)):
00390             mid = self.plies[iPly][0]
00391             #print mid
00392             self.plies[iPly][0] = model.Material(mid)  # mid
00393         ###
00394 
00395     def Nsm(self):
00396         return self.nsm
00397 
00398     def Mid(self, iPly):
00399         """
00400         gets the material ID of the ith ply
00401         @param self the object pointer
00402         @param iPly the ply ID (starts from 0)
00403         """
00404         Mid = self.Material(iPly)
00405         if isinstance(Mid, int):
00406             return Mid
00407         return Mid.mid
00408 
00409     def Mids(self):
00410         """
00411         gets the material IDs of all the plies
00412         @param self the object pointer
00413         @retval mids the material IDs
00414         """
00415         mids = []
00416         for iPly in xrange(len(self.plies)):
00417             mids.append(self.Mid(iPly))
00418         return mids
00419 
00420     def Rho(self, iPly):
00421         """
00422         gets the density of the ith ply
00423         @param self the object pointer
00424         @param iPly the ply ID (starts from 0)
00425         """
00426         mid = self.Material(iPly)
00427         return mid.rho
00428 
00429     def Material(self, iPly):
00430         """
00431         gets the material of the ith ply (not the ID unless it's not
00432         cross-referenced)
00433         @param self the object pointer
00434         @param iPly the ply ID (starts from 0)
00435         """
00436         Mid = self.plies[iPly][0]
00437         return Mid
00438 
00439     def Thickness(self, iPly='all'):
00440         """
00441         gets the thickness of the ith ply
00442         @param self the object pointer
00443         @param iPly the string 'all' (default) or the mass per area of the ith
00444         ply
00445         """
00446         if iPly == 'all': # get all layers
00447             t = 0.
00448             for iply in xrange(len(self.plies)):
00449                 t += self.Thickness(iply)
00450             ###
00451             if self.isSymmetrical():
00452                 if self.hasCoreLayer():
00453                     # cut out the thickness of half the core layer
00454                     t -= self.Thickness(0)/2.
00455                 return t*2.
00456             return t
00457             ###
00458         ###
00459         else:
00460             t = self.plies[iPly][1]
00461             return t
00462         ###
00463 
00464     def Theta(self, iPly):
00465         """
00466         gets the ply angle of the ith ply (not the ID)
00467         @param self the object pointer
00468         @param iPly the ply ID (starts from 0)
00469         """
00470         Theta = self.plies[iPly][2]
00471         return Theta
00472 
00473     def sout(self, iPly):
00474         Sout = self.plies[iPly][3]
00475         return Sout
00476 
00477     def MassPerArea(self, iPly='all'):
00478         """
00479         mass = rho*A*t
00480         but area comes from the element
00481         mass/A = rho*t for the various layers
00482         the final mass calculation will be done later
00483         @param self the object pointer
00484         @param iPly the string 'all' (default) or the mass per area of the ith ply
00485         """
00486         if iPly == 'all': # get all layers
00487             massPerArea = 0.
00488             for iply in xrange(len(self.plies)):
00489                 massPerArea += self.MassPerArea(iply)
00490             ###
00491             if self.isSymmetrical():
00492                 if self.hasCoreLayer():
00493                     # cut out the thickness of half the core layer
00494                     massPerArea -= self.MassPerArea(0)/2.
00495                 return massPerArea*2.
00496             return massPerArea
00497         ###
00498         else:
00499             rho = self.Rho(iPly)
00500             t = self.plies[iPly][1]
00501             return rho*t+self.nsm
00502         ###
00503 
00504     def D(self):
00505         D = zeros([3, 3])
00506         #isSym = self.isSymmetrical()
00507         for iply in xrange(len(self.plies)):
00508             theta = self.Theta(iply)
00509             #t    = self.Thickness(iply)
00510             mat   = self.Material(iply)
00511             Di = mat.Dplate()
00512             transform = self.T(theta)
00513             D += Di*transform
00514         return D
00515 
00516     def rawFields(self):
00517         fields = ['PCOMP', self.pid, self.z0, self.nsm, self.sb, self.ft,
00518                   self.TRef, self.ge, self.lam,]
00519         #print "plies = ",self.plies
00520         for (iPly, ply) in enumerate(self.plies):
00521             (_mid, t, theta, sout) = ply
00522             mid = self.Mid(iPly)
00523             fields += [mid, t, theta, sout]
00524         return fields
00525 
00526     def reprFields(self):
00527         #print "t = ",self.Thickness()
00528         nsm  = set_blank_if_default(self.nsm,  0.0)
00529         sb   = set_blank_if_default(self.sb,   0.0)
00530         TRef = set_blank_if_default(self.TRef, 0.0)
00531         ge   = set_blank_if_default(self.ge,   0.0)
00532         z0   = set_blank_if_default(self.z0, -0.5*self.Thickness())
00533 
00534         fields = ['PCOMP', self.pid, z0, nsm, sb, self.ft, TRef, ge, self.lam,]
00535         #print "plies = ",self.plies
00536         for (iPly, ply) in enumerate(self.plies):
00537             (_mid, t, theta, sout) = ply
00538             mid = self.Mid(iPly)
00539             #theta = set_blank_if_default(theta,0.0)
00540             sout  = set_blank_if_default(sout,'NO')
00541             fields += [mid, t, theta, sout]
00542         return fields
00543 
00544 class PCOMPG(PCOMP):  ## @todo check for bugs in ply parser
00545     type = 'PCOMPG'
00546     def __init__(self, card=None, data=None):
00547         ShellProperty.__init__(self, card, data) ## @todo doesnt support data
00548         self.pid  = card.field(1)
00549         #self.z0 = 
00550         self.nsm  = card.field(3, 0.0)
00551         self.sb   = card.field(4, 0.0)
00552         self.ft   = card.field(5)
00553         self.TRef = card.field(6, 0.0)
00554         self.ge   = card.field(7, 0.0)
00555         self.lam  = card.field(8)
00556         fields = card.fields(9)
00557 
00558         T = 0. # thickness
00559         midLast = None
00560         tLast = None
00561         self.plies = []
00562         
00563         i=0
00564         n = 0
00565         while i < len(fields):
00566             gPlyID    = card.field(9+i)
00567             mid       = card.field(9+i+1, midLast)
00568             thickness = card.field(9+i+2, tLast) # can be blank 2nd time thru
00569             theta     = card.field(9+i+3, 0.0)
00570             sout      = card.field(9+i+4, 'NO')
00571             #print('n=%s gPlyID=%s mid=%s thickness=%s len=%s' %(n,gPlyID,mid,thickness,len(fields)))
00572 
00573             ply = [mid, thickness, theta, sout, gPlyID]
00574             #print("ply = %s" %(ply))
00575             self.plies.append(ply)
00576             #[mid,t,theta,sout] # PCOMP
00577             
00578             assert mid       is not None
00579             assert thickness is not None
00580             assert isinstance(mid,int),'mid=%s' %(mid)
00581             assert isinstance(thickness,float),'thickness=%s' %(thickness)
00582             midLast = mid
00583             tLast = thickness
00584             T += thickness
00585             i += 8
00586             n += 1
00587         ###
00588         self.z0 = card.field(2, -0.5*T)
00589 
00590     def GlobalPlyID(self, iPly):
00591         gPlyID = self.plies[iPly][4]
00592         return gPlyID
00593         
00594     def rawFields(self):
00595         fields = ['PCOMPG', self.pid, self.z0, self.nsm, self.sb, self.ft,
00596                   self.TRef, self.ge, self.lam,]
00597         for (iPly, ply) in enumerate(self.plies):
00598             (_mid, t, theta, sout, gPlyID) = ply
00599             mid = self.Mid(iPly)
00600             fields += [gPlyID, mid, t, theta, sout, None, None, None]
00601         return fields
00602 
00603     def reprFields(self):
00604         nsm  = set_blank_if_default(self.nsm,  0.0)
00605         sb   = set_blank_if_default(self.sb,   0.0)
00606         TRef = set_blank_if_default(self.TRef, 0.0)
00607         ge   = set_blank_if_default(self.ge,   0.0)
00608         z0 = set_blank_if_default(self.z0, -0.5*self.Thickness())
00609 
00610         fields = ['PCOMPG', self.pid, z0, nsm, sb, self.ft, TRef, ge, self.lam]
00611 
00612         for (iPly, ply) in enumerate(self.plies):
00613             (_mid, t, theta, sout, gPlyID) = ply
00614             mid = self.Mid(iPly)
00615             #theta = set_blank_if_default(theta,0.0)
00616             sout  = set_blank_if_default(sout, 'NO')
00617             fields += [gPlyID, mid, t, theta, sout, None, None, None]
00618         return fields
00619 
00620 class PSHEAR(ShellProperty):
00621     type = 'PSHEAR'
00622     def __init__(self, card=None, nPDAMP=0, data=None):
00623         """
00624         Defines the properties of a shear panel (CSHEAR entry).
00625         PSHEAR PID MID T NSM F1 F2
00626         """
00627         ShellProperty.__init__(self, card, data)
00628         if card:
00629             ## Property ID
00630             self.pid = card.field(1)
00631             ## Material ID
00632             self.mid = card.field(2)
00633             self.t   = card.field(3)
00634             self.nsm = card.field(4, 0.0)
00635             self.f1  = card.field(5, 0.0)
00636             self.f2  = card.field(6, 0.0)
00637         else:
00638             #(pid,mid,t,nsm,f1,f2) = out
00639             self.pid = data[0]
00640             self.mid = data[1]
00641             self.t   = data[2]
00642             self.nsm = data[3]
00643             self.f1  = data[4]
00644             self.f2  = data[5]
00645         ###
00646 
00647     def isSameCard(self, prop, debug=False):
00648         if self.type != prop.type:
00649             return False
00650         fields1 = self.rawFields()
00651         fields2 = prop.rawFields()
00652         if debug:
00653             print("fields1=%s fields2=%s" % (fields1, fields2))
00654         return self.isSameFields(fields1, fields2)
00655 
00656     def rawFields(self):
00657         fields = ['PSHEAR', self.pid, self.Mid(), self.t, self.nsm,
00658                   self.f1, self.f2]
00659         return fields
00660 
00661 class PSHELL(ShellProperty):
00662     """
00663     PSHELL PID MID1 T MID2 12I/T**3 MID3 TS/T NSM
00664     Z1 Z2 MID4
00665     PSHELL   41111   1      1.0000   1               1               0.02081"""
00666     type = 'PSHELL'
00667     def __init__(self, card=None, data=None):
00668         ShellProperty.__init__(self, card, data)
00669         if card:
00670             self.pid  = int(card.field(1))
00671             self.mid1 = card.field(2)
00672             self.t    = card.field(3)
00673 
00674             ## Material identification number for bending
00675             self.mid2 = card.field(4)
00676             ## \f$ \frac{12I}{t^3} \f$
00677             self.twelveIt3 = card.field(5, 1.0)  # poor name
00678             self.mid3  = card.field(6)
00679             self.tst   = card.field(7, 0.833333)
00680             self.nsm   = card.field(8, 0.0)
00681 
00682             tOver2 = self.t/2.
00683             self.z1    = card.field(9, -tOver2)
00684             self.z2    = card.field(10, tOver2)
00685             self.mid4  = card.field(11)
00686             #if self.mid2 is None:
00687             #    assert self.mid3 is None
00688             #else: # mid2 is defined
00689             #    #print "self.mid2 = ",self.mid2
00690             #    assert self.mid2 >= -1
00691             #    #assert self.mid3 >   0
00692 
00693             #if self.mid is not None and self.mid2 is not None:
00694             #    assert self.mid4==None
00695             ###
00696         else:
00697             self.pid       = data[0]
00698             self.mid1      = data[1]
00699             self.t         = data[2]
00700             self.mid2      = data[3]
00701             self.twelveIt3 = data[4]
00702             self.mid3      = data[5]
00703             self.tst       = data[6]
00704             self.nsm       = data[7]
00705             self.z1        = data[8]
00706             self.z2        = data[9]
00707             self.mid4      = data[10]
00708             #maxMid = max(self.mid,self.mid2,self.mid3,self.mid4)
00709         ###
00710 
00711         assert self.t > 0.0,('the thickness must be defined on the PSHELL'
00712                              ' card (Ti field not supported)')
00713 
00714     def mid(self):
00715         if isinstance(self.mid1, Material):
00716             return self.mid1
00717         return self.mid2
00718 
00719     def Mid(self):
00720         if isinstance(self.mid1, Material):
00721             return self.mid1.mid  # self.Mid1()
00722         return self.Mid2()
00723     
00724     def Mid1(self):
00725         if isinstance(self.mid1, Material):
00726             return self.mid1.mid
00727         return self.mid1
00728 
00729     def Mid2(self):
00730         if isinstance(self.mid2, Material):
00731             return self.mid2.mid
00732         return self.mid2
00733 
00734     def Mid3(self):
00735         if isinstance(self.mid3, Material):
00736             return self.mid3.mid
00737         return self.mid3
00738 
00739     def Mid4(self):
00740         if isinstance(self.mid4, Material):
00741             return self.mid4.mid
00742         return self.mid4
00743 
00744     def Thickness(self):
00745         return self.t
00746 
00747     def Rho(self):
00748         return self.mid().rho
00749 
00750     def Nsm(self):
00751         return self.nsm
00752 
00753     def MassPerArea(self):
00754         """
00755         calculates mass per area
00756         \f[ \large  \frac{mass}{A} = nsm + \rho t \f]
00757         """
00758         try:
00759             massPerArea = self.nsm + self.Rho()*self.t
00760         except:
00761             print("nsm=%s rho=%s t=%s" % (self.nsm, self.Rho(), self.t))
00762             raise
00763         return massPerArea
00764 
00765     def D(self):
00766         return self.mid().Dplate()
00767         
00768     def crossReference(self, model):
00769         if self.mid1:
00770             self.mid1 = model.Material(self.mid1)
00771         if self.mid2 and self.mid2 != -1:
00772             self.mid2 = model.Material(self.mid2)
00773         if self.mid3:
00774             self.mid3 = model.Material(self.mid3)
00775         if self.mid4:
00776             self.mid4 = model.Material(self.mid4)
00777 
00778     def writeCalculix(self, marker='markerDummyProp',
00779                       elementSet='ELsetDummyProp'):
00780         msg = '*SHELL SECTION,MATERIAL=M%s_%s,ELSET=%s,OFFSET=%s\n' % (marker, self.mid, elementSet, self.z1)
00781         msg += '** THICKNESS\n'
00782         msg += '%s\n\n' % (self.t)
00783         return msg
00784 
00785     def writeCodeAster(self):
00786         """
00787         http://www.caelinux.org/wiki/index.php/Contrib:KeesWouters/shell/static
00788         http://www.caelinux.org/wiki/index.php/Contrib:KeesWouters/platedynamics
00789         the angle_rep is a direction angle, use either angle(a,b) or
00790         vecteur(x,y,z)
00791         coque_ncou is the number of gauss nodes along the thickness, for
00792         linear analysis one node is sufficient.
00793         """
00794         msg = ''
00795         msg += "    COQUE=_F(GROUP_MA='P%s', # COQUE=PSHELL\n" % (self.pid)
00796         msg += "              EPAIS=%g, # EPAIS=thickness\n" % (self.t)
00797         msg += "              ANGL_REP=(0.,90.),     # ???\n"  ## @todo what is this?
00798         #msg += "              VECTEUR=(1.0,0.0,0.0,) #  local coordinate system\n"
00799         msg += "              EXCENTREMENT=%g,       # offset-Z1\n" % (self.z1)
00800         msg += "              COQUE_NCOU=1,          # Number of Integration Layers\n"
00801         msg += "              CARA=('NSM'), # ???\n"  ## @todo check
00802         msg += "              VALE=(%g),),\n"  % (self.nsm)
00803         return msg
00804 
00805     def rawFields(self):
00806         fields = ['PSHELL', self.pid, self.Mid1(), self.t, self.Mid2(),
00807                   self.twelveIt3, self.Mid3(), self.tst, self.nsm, self.z1,
00808                   self.z2, self.Mid4()]
00809         return fields
00810 
00811     def reprFields(self):
00812         twelveIt3 = set_blank_if_default(self.twelveIt3, 1.0)
00813         tst       = set_blank_if_default(self.tst, 0.833333)
00814         nsm       = set_blank_if_default(self.nsm, 0.0)
00815 
00816         tOver2 = self.t/2.
00817         z1 = set_blank_if_default(self.z1, -tOver2)
00818         z2 = set_blank_if_default(self.z2,  tOver2)
00819 
00820         fields = ['PSHELL', self.pid, self.Mid1(), self.t, self.Mid2(),
00821                   twelveIt3, self.Mid3(), tst, nsm, z1, z2, self.Mid4()]
00822         #print fields
00823         return fields
00824 
 All Classes Namespaces Files Functions Variables