pyNastran
0.5.0
pyNastran BDF Reader/Writer, OP2 Parser, and GUI
|
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 00026 from __future__ import division, print_function 00027 import sys 00028 from numpy import array 00029 from numpy.linalg import eigh 00030 00031 from .oes_objects import stressObject,strainObject 00032 00033 00034 class SolidStressObject(stressObject): 00035 """ 00036 # sCode=0 00037 N O N L I N E A R S T R E S S E S I N T E T R A H E D R O N S O L I D E L E M E N T S ( T E T R A ) 00038 ELEMENT GRID/ POINT STRESSES/ TOTAL STRAINS EQUIVALENT EFF. STRAIN EFF. CREEP 00039 ID GAUSS ID X Y Z XY YZ ZX STRESS PLAS/NLELAS STRAIN 00040 3 GRID CENTER 6.6667E+02 2.6667E+03 2.6667E+03 -1.3333E+03 2.6667E+03 -1.3333E+03 6.0000E+03 1.5000E-04 0.0 00041 1.6667E-05 6.6667E-05 6.6667E-05 -6.6667E-05 1.3333E-04 -6.6667E-05 00042 00043 # sCode=1 00044 S T R E S S E S I N H E X A H E D R O N S O L I D E L E M E N T S ( H E X A ) 00045 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN 00046 ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE VON MISES 00047 1 0GRID CS 8 GP 00048 CENTER X 4.499200E+02 XY -5.544791E+02 A 1.000000E+04 LX 0.00 0.69-0.72 -3.619779E+03 9.618462E+03 00049 Y 4.094179E+02 YZ 5.456968E-12 B -1.251798E+02 LY 0.00 0.72 0.69 00050 Z 1.000000E+04 ZX -4.547474E-13 C 9.845177E+02 LZ 1.00 0.00 0.00 00051 """ 00052 def __init__(self,dataCode,isSort1,iSubcase,dt=None): 00053 stressObject.__init__(self,dataCode,iSubcase) 00054 00055 self.eType = {} 00056 self.code = [self.formatCode,self.sortCode,self.sCode] 00057 00058 self.cid = {} # gridGauss 00059 self.oxx = {} 00060 self.oyy = {} 00061 self.ozz = {} 00062 self.txy = {} 00063 self.tyz = {} 00064 self.txz = {} 00065 self.o1 = {} 00066 self.o2 = {} 00067 self.o3 = {} 00068 self.ovmShear = {} 00069 00070 self.dt = dt 00071 if isSort1: 00072 if dt is not None: 00073 self.add = self.addSort1 00074 self.addNewEid = self.addNewEidSort1 00075 ### 00076 else: 00077 assert dt is not None 00078 self.add = self.addSort2 00079 self.addNewEid = self.addNewEidSort2 00080 ### 00081 00082 def addF06Data(self,data,transient): 00083 if transient is None: 00084 if not hasattr(self,'data'): 00085 self.data = [] 00086 self.data += data 00087 else: 00088 if not hasattr(self,'data'): 00089 self.data = {} 00090 if dt not in self.data: 00091 self.data[dt] = [] 00092 for line in data: 00093 self.data[dt] += data 00094 ### 00095 ### 00096 00097 def processF06Data(self): 00098 """ 00099 S T R E S S E S I N P E N T A H E D R O N S O L I D E L E M E N T S ( P E N T A ) 00100 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN 00101 ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE VON MISES 00102 2 0GRID CS 6 GP 00103 CENTER X -1.829319E+03 XY 7.883865E+01 A 1.033182E+04 LX-0.12 0.71 0.70 -2.115135E+03 1.232595E+04 00104 Y -1.825509E+03 YZ -1.415218E+03 B -2.080181E+03 LY-0.12 0.69-0.71 00105 Z 1.000023E+04 ZX -1.415218E+03 C -1.906232E+03 LZ 0.99 0.16 0.00 00106 """ 00107 eMap = {'CTETRA':5,'CPENTA':7,'CHEXA':9,'HEXA':9,'PENTA':7,'TETRA':5,} # +1 for the centroid 00108 if self.nonlinearFactor is None: 00109 ipack = [] 00110 i = 0 00111 n = 0 00112 while n<len(self.data): 00113 line = self.data[n] 00114 #print n,line 00115 00116 eType = line[0] 00117 eid = int(line[1]) 00118 #print "eType = ",eType 00119 nNodes = eMap[eType] 00120 #nodeID = None 00121 self.eType[eid] = eType 00122 self.oxx[eid] = {} 00123 self.oyy[eid] = {} 00124 self.ozz[eid] = {} 00125 self.txy[eid] = {} 00126 self.tyz[eid] = {} 00127 self.txz[eid] = {} 00128 self.o1[eid] = {} 00129 self.o2[eid] = {} 00130 self.o3[eid] = {} 00131 self.ovmShear[eid] = {} 00132 n += 1 00133 for j in xrange(nNodes): 00134 #print self.data[n] 00135 (blank,nodeID,x,oxx,xy,txy,a,o1,lx,d1,d2,d3,pressure,ovmShear) = self.data[n] 00136 (blank,blank, y,oyy,yz,tyz,b,o2,ly,d1,d2,d3,blank,blank) = self.data[n+1] 00137 (blank,blank, z,ozz,zx,txz,c,o3,lz,d1,d2,d3,blank,blank) = self.data[n+2] 00138 if nodeID.strip()=='CENTER': nodeID='C' 00139 else: nodeID = int(nodeID) 00140 00141 self.oxx[eid][nodeID] = float(oxx) 00142 self.oyy[eid][nodeID] = float(oyy) 00143 self.ozz[eid][nodeID] = float(ozz) 00144 self.txy[eid][nodeID] = float(txy) 00145 self.tyz[eid][nodeID] = float(tyz) 00146 self.txz[eid][nodeID] = float(txz) 00147 self.o1[eid][nodeID] = float(o1) 00148 self.o2[eid][nodeID] = float(o2) 00149 self.o3[eid][nodeID] = float(o3) 00150 self.ovmShear[eid][nodeID] = float(ovmShear) 00151 n+=3 00152 ### 00153 ### 00154 return 00155 raise NotImplementedError() 00156 00157 (dtName,dt) = transient 00158 self.dataCode['name'] = dtName 00159 if dt not in self.gridTypes: 00160 self.addNewTransient() 00161 00162 for line in data: 00163 (nodeID,gridType,t1,t2,t3,r1,r2,r3) = line 00164 self.gridTypes[dt][nodeID] = array([t1,t2,t3]) 00165 self.translations[dt][nodeID] = array([t1,t2,t3]) 00166 self.rotations[dt][nodeID] = array([r1,r2,r3]) 00167 ### 00168 del self.data 00169 00170 def deleteTransient(self,dt): 00171 del self.oxx[dt] 00172 del self.oyy[dt] 00173 del self.ozz[dt] 00174 del self.txy[dt] 00175 del self.tyz[dt] 00176 del self.txz[dt] 00177 del self.o1[dt] 00178 del self.o2[dt] 00179 del self.o3[dt] 00180 00181 def getTransients(self): 00182 k = self.oxx.keys() 00183 k.sort() 00184 return k 00185 00186 def addNewTransient(self,dt): 00187 """ 00188 initializes the transient variables 00189 """ 00190 self.oxx[dt] = {} 00191 self.oyy[dt] = {} 00192 self.ozz[dt] = {} 00193 self.txy[dt] = {} 00194 self.tyz[dt] = {} 00195 self.txz[dt] = {} 00196 self.o1[dt] = {} 00197 self.o2[dt] = {} 00198 self.o3[dt] = {} 00199 00200 #self.aCos[dt] = {} 00201 #self.bCos[dt] = {} 00202 #self.cCos[dt] = {} 00203 #self.pressure[dt] = {} 00204 self.ovmShear[dt] = {} 00205 00206 def addNewEid(self,eType,cid,dt,eid,nodeID,oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,aCos,bCos,cCos,pressure,ovm): 00207 #print "Solid Stress add..." 00208 #assert eid not in self.oxx 00209 assert cid >= 0 00210 assert eid >= 0 00211 self.eType[eid] = eType 00212 self.cid[eid] = cid 00213 self.oxx[eid] = {nodeID: oxx} 00214 self.oyy[eid] = {nodeID: oyy} 00215 self.ozz[eid] = {nodeID: ozz} 00216 self.txy[eid] = {nodeID: txy} 00217 self.tyz[eid] = {nodeID: tyz} 00218 self.txz[eid] = {nodeID: txz} 00219 self.o1[eid] = {nodeID: o1} 00220 self.o2[eid] = {nodeID: o2} 00221 self.o3[eid] = {nodeID: o3} 00222 #self.aCos[eid] = {nodeID: aCos} 00223 #self.bCos[eid] = {nodeID: bCos} 00224 #self.cCos[eid] = {nodeID: cCos} 00225 #self.pressure[eid] = {nodeID: pressure} 00226 self.ovmShear[eid] = {nodeID: ovm} 00227 msg = "*eid=%s nodeID=%s vm=%g" %(eid,nodeID,ovm) 00228 #print msg 00229 if nodeID==0: raise Exception(msg) 00230 00231 def addNewEidSort1(self,eType,cid,dt,eid,nodeID,oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,aCos,bCos,cCos,pressure,ovm): 00232 #print "Solid Stress add transient..." 00233 assert cid >= 0 00234 assert eid >= 0 00235 00236 #print "dt=%s eid=%s eType=%s" %(dt,eid,eType) 00237 if dt not in self.oxx: 00238 self.addNewTransient(dt) 00239 assert eid not in self.oxx[dt],self.oxx[dt] 00240 00241 self.eType[eid] = eType 00242 self.cid[eid] = cid 00243 self.oxx[dt][eid] = {nodeID: oxx} 00244 self.oyy[dt][eid] = {nodeID: oyy} 00245 self.ozz[dt][eid] = {nodeID: ozz} 00246 self.txy[dt][eid] = {nodeID: txy} 00247 self.tyz[dt][eid] = {nodeID: tyz} 00248 self.txz[dt][eid] = {nodeID: txz} 00249 00250 self.o1[dt][eid] = {nodeID: o1} 00251 self.o2[dt][eid] = {nodeID: o2} 00252 self.o3[dt][eid] = {nodeID: o3} 00253 00254 #self.aCos[dt][eid] = {nodeID: aCos} 00255 #self.bCos[dt][eid] = {nodeID: bCos} 00256 #self.cCos[dt][eid] = {nodeID: cCos} 00257 #self.pressure[dt][eid] = {nodeID: pressure} 00258 self.ovmShear[dt][eid] = {nodeID: ovm} 00259 msg = "*eid=%s nodeID=%s vm=%g" %(eid,nodeID,ovm) 00260 #print msg 00261 if nodeID==0: raise Exception(msg) 00262 00263 def add(self,dt,eid,nodeID,oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,aCos,bCos,cCos,pressure,ovm): 00264 #print "***add" 00265 msg = "eid=%s nodeID=%s vm=%g" %(eid,nodeID,ovm) 00266 #print msg 00267 #print self.oxx 00268 #print self.fiberDistance 00269 self.oxx[eid][nodeID] = oxx 00270 self.oyy[eid][nodeID] = oyy 00271 self.ozz[eid][nodeID] = ozz 00272 00273 self.txy[eid][nodeID] = txy 00274 self.tyz[eid][nodeID] = tyz 00275 self.txz[eid][nodeID] = txz 00276 00277 self.o1[eid][nodeID] = o1 00278 self.o2[eid][nodeID] = o2 00279 self.o3[eid][nodeID] = o3 00280 00281 #self.aCos[eid][nodeID] = aCos 00282 #self.bCos[eid][nodeID] = bCos 00283 #self.cCos[eid][nodeID] = cCos 00284 #self.pressure[eid][nodeID] = pressure 00285 self.ovmShear[eid][nodeID] = ovm 00286 if nodeID==0: raise Exception(msg) 00287 00288 def addSort1(self,dt,eid,nodeID,oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,aCos,bCos,cCos,pressure,ovm): 00289 #print "***add" 00290 msg = "eid=%s nodeID=%s vm=%g" %(eid,nodeID,ovm) 00291 #print msg 00292 #print self.oxx 00293 #print self.fiberDistance 00294 00295 #self.eType[eid] = eType 00296 #print "eid=%s nid=%s oxx=%s" %(eid,nodeID,oxx) 00297 self.oxx[dt][eid][nodeID] = oxx 00298 self.oyy[dt][eid][nodeID] = oyy 00299 self.ozz[dt][eid][nodeID] = ozz 00300 #print self.oxx 00301 self.txy[dt][eid][nodeID] = txy 00302 self.tyz[dt][eid][nodeID] = tyz 00303 self.txz[dt][eid][nodeID] = txz 00304 00305 self.o1[dt][eid][nodeID] = o1 00306 self.o2[dt][eid][nodeID] = o2 00307 self.o3[dt][eid][nodeID] = o3 00308 00309 #self.aCos[dt][eid][nodeID] = aCos 00310 #self.bCos[dt][eid][nodeID] = bCos 00311 #self.cCos[dt][eid][nodeID] = cCos 00312 #self.pressure[dt][eid][nodeID] = pressure 00313 self.ovmShear[dt][eid][nodeID] = ovm 00314 if nodeID==0: raise Exception(msg) 00315 00316 def getHeaders(self): 00317 headers = ['oxx','oyy','ozz','txy','tyz','txz'] 00318 if self.isVonMises(): 00319 headers.append('oVonMises') 00320 else: 00321 headers.append('oMaxShear') 00322 return headers 00323 00324 def pressure(self,o1,o2,o3): 00325 """ 00326 returns the hydrostatic pressure 00327 (o1+o2+o3)/-3. 00328 http://en.wikipedia.org/wiki/Stress_%28mechanics%29 00329 """ 00330 return (o1+o2+o3)/-3. 00331 00332 def directionalVectors(self,oxx,oyy,ozz,txy,tyz,txz): 00333 A = [[oxx,txy,txz], 00334 [txy,oyy,tyz], 00335 [txz,tyz,ozz]] 00336 (Lambda,v) = eigh(A) # a hermitian matrix is a symmetric-real matrix 00337 return v 00338 00339 def getF06_Header(self): 00340 if self.isVonMises(): 00341 vonMises = 'VON MISES' 00342 else: 00343 vonMises = 'MAX SHEAR' 00344 00345 tetraMsg = [' S T R E S S E S I N T E T R A H E D R O N S O L I D E L E M E N T S ( C T E T R A )\n', 00346 '0 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN \n', 00347 ' ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE %s \n' %(vonMises)] 00348 00349 pentaMsg = [' S T R E S S E S I N P E N T A H E D R O N S O L I D E L E M E N T S ( P E N T A )\n', 00350 '0 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN \n', 00351 ' ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE %s \n' %(vonMises)] 00352 00353 hexaMsg = [' S T R E S S E S I N H E X A H E D R O N S O L I D E L E M E N T S ( H E X A )\n', 00354 '0 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN \n', 00355 ' ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE %s \n' %(vonMises)] 00356 00357 #nNodes = {'CTETRA':4,'CPENTA':6,'CHEXA':8,'HEXA':8,'PENTA':6,'TETRA':4,} 00358 tetraEids = [] 00359 hexaEids = [] 00360 pentaEids = [] 00361 for eid,eType in sorted(self.eType.iteritems()): 00362 if eType=='CTETRA' or eType=='TETRA': 00363 tetraEids.append(eid) 00364 elif eType=='CPENTA' or eType=='PENTA': 00365 pentaEids.append(eid) 00366 elif eType=='CHEXA' or eType=='HEXA': 00367 hexaEids.append(eid) 00368 else: 00369 raise NotImplementedError('eType=|%r|' %(eType)) 00370 ### 00371 ### 00372 return (tetraMsg,pentaMsg,hexaMsg,tetraEids,hexaEids,pentaEids) 00373 00374 def writeF06(self,header,pageStamp,pageNum=1,f=None,isMagPhase=False): 00375 if self.nonlinearFactor is not None: 00376 return self.writeF06Transient(header,pageStamp,pageNum,f) 00377 msg = [] 00378 00379 (tetraMsg,pentaMsg,hexaMsg,tetraEids,hexaEids,pentaEids) = self.getF06_Header() 00380 #nNodes = {'CTETRA':4,'CPENTA':6,'CHEXA':8,'HEXA':8,'PENTA':6,'TETRA':4,} 00381 if tetraEids: 00382 msg += self.writeElement('CTETRA',4,tetraEids,header,tetraMsg,f) 00383 msg.append(pageStamp+str(pageNum)+'\n') 00384 pageNum+=1 00385 if hexaEids: 00386 msg += self.writeElement('CHEXA',8,hexaEids, header,hexaMsg,f) 00387 msg.append(pageStamp+str(pageNum)+'\n') 00388 pageNum+=1 00389 if pentaEids: 00390 msg += self.writeElement('CPENTA',6,pentaEids,header,pentaMsg,f) 00391 msg.append(pageStamp+str(pageNum)+'\n') 00392 pageNum+=1 00393 ### 00394 return (''.join(msg),pageNum-1) 00395 00396 def writeF06Transient(self,header,pageStamp,pageNum=1,f=None,isMagPhase=False): 00397 msg = [] 00398 (tetraMsg,pentaMsg,hexaMsg,tetraEids,hexaEids,pentaEids) = self.getF06_Header() 00399 dts = self.oxx.keys() 00400 for dt in dts: 00401 if tetraEids: 00402 msg += self.writeElementTransient('CTETRA',4,tetraEids,dt,header,tetraMsg,f) 00403 msg.append(pageStamp+str(pageNum)+'\n') 00404 pageNum+=1 00405 if hexaEids: 00406 msg += self.writeElementTransient('CHEXA',8,hexaEids, dt,header,hexaMsg,f) 00407 msg.append(pageStamp+str(pageNum)+'\n') 00408 pageNum+=1 00409 if pentaEids: 00410 msg += self.writeElementTransient('CPENTA',6,pentaEids,dt,header,pentaMsg,f) 00411 msg.append(pageStamp+str(pageNum)+'\n') 00412 pageNum+=1 00413 ### 00414 return (''.join(msg),pageNum-1) 00415 00416 def writeElement(self,eType,nNodes,eids,header,tetraMsg,f=None): 00417 msg = header+tetraMsg 00418 for eid in eids: 00419 #eType = self.eType[eid] 00420 00421 k = self.oxx[eid].keys() 00422 k.remove('C') 00423 k.sort() 00424 msg.append('0 %8s 0GRID CS %i GP\n' %(eid,nNodes) ) 00425 for nid in ['C']+k: 00426 oxx = self.oxx[eid][nid] 00427 oyy = self.oyy[eid][nid] 00428 ozz = self.ozz[eid][nid] 00429 txy = self.txy[eid][nid] 00430 tyz = self.tyz[eid][nid] 00431 txz = self.txz[eid][nid] 00432 00433 o1 = self.o1[eid][nid] 00434 o2 = self.o2[eid][nid] 00435 o3 = self.o3[eid][nid] 00436 ovm = self.ovmShear[eid][nid] 00437 p = (o1+o2+o3)/-3. 00438 00439 if nid=='C': nid='CENTER' 00440 #print "nid = |%r|" %(nid) 00441 A = [[oxx,txy,txz], 00442 [txy,oyy,tyz], 00443 [txz,tyz,ozz]] 00444 (Lambda,v) = eigh(A) # a hermitian matrix is a symmetric-real matrix 00445 00446 ([oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,p,ovm],isAllZeros) = self.writeFloats13E([oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,p,ovm]) 00447 msg.append('0 %8s X %13s XY %13s A %13s LX%5.2f%5.2f%5.2f %13s %-s\n' %(nid,oxx,txy,o1,v[0,1],v[0,2],v[0,0],p,ovm.strip())) 00448 msg.append(' %8s Y %13s YZ %13s B %13s LY%5.2f%5.2f%5.2f\n' %('', oyy,tyz,o2,v[1,1],v[1,2],v[1,0])) 00449 msg.append(' %8s Z %13s ZX %13s C %13s LZ%5.2f%5.2f%5.2f\n' %('', ozz,txz,o3,v[2,1],v[2,2],v[2,0])) 00450 ### 00451 ### 00452 if f is not None: 00453 f.write(''.join(msg)) 00454 ### 00455 return ''.join(msg) 00456 00457 def writeElementTransient(self,eType,nNodes,eids,dt,header,tetraMsg,f=None): 00458 dtLine = '%14s = %12.5E\n'%(self.dataCode['name'],dt) 00459 header[1] = dtLine 00460 msg = header+tetraMsg 00461 for eid in eids: 00462 #eType = self.eType[eid] 00463 00464 k = self.oxx[dt][eid].keys() 00465 k.remove('C') 00466 k.sort() 00467 msg.append('0 %8s 0GRID CS %i GP\n' %(eid,nNodes)) 00468 for nid in ['C']+k: 00469 oxx = self.oxx[dt][eid][nid] 00470 oyy = self.oyy[dt][eid][nid] 00471 ozz = self.ozz[dt][eid][nid] 00472 txy = self.txy[dt][eid][nid] 00473 tyz = self.tyz[dt][eid][nid] 00474 txz = self.txz[dt][eid][nid] 00475 00476 o1 = self.o1[dt][eid][nid] 00477 o2 = self.o2[dt][eid][nid] 00478 o3 = self.o3[dt][eid][nid] 00479 ovm = self.ovmShear[dt][eid][nid] 00480 p = (o1+o2+o3)/-3. 00481 00482 if nid=='C': nid='CENTER' 00483 #print "nid = |%r|" %(nid) 00484 A = [[oxx,txy,txz], 00485 [txy,oyy,tyz], 00486 [txz,tyz,ozz]] 00487 (Lambda,v) = eigh(A) # a hermitian matrix is a symmetric-real matrix 00488 00489 ([oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,p,ovm],isAllZeros) = self.writeFloats13E([oxx,oyy,ozz,txy,tyz,txz,o1,o2,o3,p,ovm]) 00490 msg.append('0 %8s X %13s XY %13s A %13s LX%5.2f%5.2f%5.2f %13s %-s\n' %(nid,oxx,txy,o1,v[0,1],v[0,2],v[0,0],p,ovm.strip())) 00491 msg.append(' %8s Y %13s YZ %13s B %13s LY%5.2f%5.2f%5.2f\n' %('', oyy,tyz,o2,v[1,1],v[1,2],v[1,0])) 00492 msg.append(' %8s Z %13s ZX %13s C %13s LZ%5.2f%5.2f%5.2f\n' %('', ozz,txz,o3,v[2,1],v[2,2],v[2,0])) 00493 ### 00494 ### 00495 if f is not None: 00496 f.write(''.join(msg)) 00497 msg = [''] 00498 ### 00499 return ''.join(msg) 00500 00501 def __repr__(self): 00502 #print "self.dt = ",self.dt 00503 #print "self.nf = ",self.nonlinearFactor 00504 if self.nonlinearFactor is not None: 00505 return self.__reprTransient__() 00506 00507 msg = '---SOLID STRESS---\n' 00508 headers = self.getHeaders() 00509 msg += '%-6s %6s %8s ' %('EID','eType','nodeID') 00510 for header in headers: 00511 msg += '%9s ' %(header) 00512 msg += '\n' 00513 for eid,oxxNodes in sorted(self.oxx.iteritems()): 00514 eType = self.eType[eid] 00515 for nid in sorted(oxxNodes): 00516 oxx = self.oxx[eid][nid] 00517 oyy = self.oyy[eid][nid] 00518 ozz = self.ozz[eid][nid] 00519 txy = self.txy[eid][nid] 00520 tyz = self.tyz[eid][nid] 00521 txz = self.txz[eid][nid] 00522 00523 #o1 = self.o1[eid][nid] 00524 #o2 = self.o2[eid][nid] 00525 #o3 = self.o3[eid][nid] 00526 ovm = self.ovmShear[eid][nid] 00527 msg += '%-6i %6s %8s ' %(eid,eType,nid) 00528 vals = [oxx,oyy,ozz,txy,tyz,txz,ovm] 00529 for val in vals: 00530 if abs(val)<1e-6: 00531 msg += '%9s ' %('0') 00532 else: 00533 msg += '%9i ' %(val) 00534 ### 00535 msg += '\n' 00536 #msg += "eid=%-4s eType=%-6s nid=%-2i oxx=%-5i oyy=%-5i ozz=%-5i txy=%-5i tyz=%-5i txz=%-5i ovm=%-5i\n" %(eid,eType,nid,oxx,oyy,ozz,txy,tyz,txz,ovm) 00537 ### 00538 ### 00539 return msg 00540 00541 def __reprTransient__(self): 00542 msg = '---SOLID STRESS---\n' 00543 msg += '%-6s %6s %8s ' %('EID','eType','nodeID') 00544 headers = self.getHeaders() 00545 for header in headers: 00546 msg += '%9s ' %(header) 00547 msg += '\n' 00548 00549 for dt,oxxs in sorted(self.oxx.iteritems()): 00550 msg += '%s = %g\n' %(self.dataCode['name'],dt) 00551 for eid,oxxNodes in sorted(oxxs.iteritems()): 00552 eType = self.eType[eid] 00553 for nid in sorted(oxxNodes): 00554 oxx = self.oxx[dt][eid][nid] 00555 oyy = self.oyy[dt][eid][nid] 00556 ozz = self.ozz[dt][eid][nid] 00557 txy = self.txy[dt][eid][nid] 00558 tyz = self.tyz[dt][eid][nid] 00559 txz = self.txz[dt][eid][nid] 00560 00561 #o1 = self.o1[dt][eid][nid] 00562 #o2 = self.o2[dt][eid][nid] 00563 #o3 = self.o3[dt][eid][nid] 00564 00565 ovm = self.ovmShear[dt][eid][nid] 00566 msg += '%-6i %6s %8s ' %(eid,eType,nid) 00567 vals = [oxx,oyy,ozz,txy,tyz,txz,ovm] 00568 for val in vals: 00569 if abs(val)<1e-6: 00570 msg += '%9s ' %('0') 00571 else: 00572 msg += '%9i ' %(val) 00573 ### 00574 msg += '\n' 00575 #msg += "eid=%-4s eType=%-6s nid=%-2i oxx=%-5i oyy=%-5i ozz=%-5i txy=%-5i tyz=%-5i txz=%-5i ovm=%-5i\n" %(eid,eType,nid,oxx,oyy,ozz,txy,tyz,txz,ovm) 00576 ### 00577 ### 00578 ### 00579 return msg 00580 00581 class SolidStrainObject(strainObject): 00582 """ 00583 # code=[1,0,11] 00584 S T R A I N S I N H E X A H E D R O N S O L I D E L E M E N T S ( H E X A ) 00585 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN 00586 ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE VON MISES 00587 1 0GRID CS 8 GP 00588 CENTER X 4.499200E+02 XY -5.544791E+02 A 1.000000E+04 LX 0.00 0.69-0.72 -3.619779E+03 9.618462E+03 00589 Y 4.094179E+02 YZ 5.456968E-12 B -1.251798E+02 LY 0.00 0.72 0.69 00590 Z 1.000000E+04 ZX -4.547474E-13 C 9.845177E+02 LZ 1.00 0.00 0.00 00591 00592 # code=[1,0,10] 00593 S T R A I N S I N T E T R A H E D R O N S O L I D E L E M E N T S ( C T E T R A ) 00594 CORNER ------CENTER AND CORNER POINT STRAINS--------- DIR. COSINES MEAN OCTAHEDRAL 00595 ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE SHEAR 00596 4 0GRID CS 4 GP 00597 CENTER X -2.288232E-04 XY 1.240506E-04 A 9.631978E-04 LX-0.10-0.71-0.70 -1.601805E-04 5.692614E-04 00598 Y -2.289814E-04 YZ -2.369997E-04 B -2.909276E-04 LY-0.10 0.71-0.70 00599 Z 9.383460E-04 ZX -2.369997E-04 C -1.917288E-04 LZ 0.99 0.00-0.15 00600 """ 00601 def __init__(self,dataCode,isSort1,iSubcase,dt=None): 00602 strainObject.__init__(self,dataCode,iSubcase) 00603 self.eType = {} 00604 00605 self.code = [self.formatCode,self.sortCode,self.sCode] 00606 00607 self.cid = {} 00608 self.exx = {} 00609 self.eyy = {} 00610 self.ezz = {} 00611 self.exy = {} 00612 self.eyz = {} 00613 self.exz = {} 00614 self.e1 = {} 00615 self.e2 = {} 00616 self.e3 = {} 00617 #self.aCos = {} 00618 #self.bCos = {} 00619 #self.cCos = {} 00620 #self.pressure = {} 00621 self.evmShear = {} 00622 00623 00624 self.nonlinearFactor = dt 00625 if isSort1: 00626 if dt is not None: 00627 self.add = self.addSort1 00628 self.addNewEid = self.addNewEidSort1 00629 ### 00630 else: 00631 assert dt is not None 00632 self.add = self.addSort2 00633 self.addNewEid = self.addNewEidSort2 00634 ### 00635 00636 def addF06Data(self,data,transient): 00637 if transient is None: 00638 if not hasattr(self,'data'): 00639 self.data = [] 00640 self.data += data 00641 else: 00642 if not hasattr(self,'data'): 00643 self.data = {} 00644 if dt not in self.data: 00645 self.data[dt] = [] 00646 for line in data: 00647 self.data[dt] += data 00648 ### 00649 ### 00650 00651 def processF06Data(self): 00652 """ 00653 S T R E S S E S I N P E N T A H E D R O N S O L I D E L E M E N T S ( P E N T A ) 00654 CORNER ------CENTER AND CORNER POINT STRESSES--------- DIR. COSINES MEAN 00655 ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE VON MISES 00656 2 0GRID CS 6 GP 00657 CENTER X -1.829319E+03 XY 7.883865E+01 A 1.033182E+04 LX-0.12 0.71 0.70 -2.115135E+03 1.232595E+04 00658 Y -1.825509E+03 YZ -1.415218E+03 B -2.080181E+03 LY-0.12 0.69-0.71 00659 Z 1.000023E+04 ZX -1.415218E+03 C -1.906232E+03 LZ 0.99 0.16 0.00 00660 """ 00661 eMap = {'CTETRA':5,'CPENTA':7,'CHEXA':9,'HEXA':9,'PENTA':7,'TETRA':5,} # +1 for the centroid 00662 if self.nonlinearFactor is None: 00663 pack = [] 00664 i = 0 00665 n = 0 00666 while n<len(self.data): 00667 line = self.data[n] 00668 #print n,line 00669 00670 eType = line[0] 00671 eid = int(line[1]) 00672 #print "eType = ",eType 00673 nNodes = eMap[eType] 00674 #nodeID = None 00675 self.eType[eid] = eType 00676 self.exx[eid] = {} 00677 self.eyy[eid] = {} 00678 self.ezz[eid] = {} 00679 self.exy[eid] = {} 00680 self.eyz[eid] = {} 00681 self.exz[eid] = {} 00682 self.e1[eid] = {} 00683 self.e2[eid] = {} 00684 self.e3[eid] = {} 00685 self.evmShear[eid] = {} 00686 n += 1 00687 for j in xrange(nNodes): 00688 #print self.data[n] 00689 (blank,nodeID,x,exx,xy,exy,a,e1,lx,d1,d2,d3,pressure,evmShear) = self.data[n] 00690 (blank,blank, y,eyy,yz,eyz,b,e2,ly,d1,d2,d3,blank,blank) = self.data[n+1] 00691 (blank,blank, z,ezz,zx,exz,c,e3,lz,d1,d2,d3,blank,blank) = self.data[n+2] 00692 if nodeID.strip()=='CENTER': nodeID='C' 00693 else: nodeID = int(nodeID) 00694 self.exx[eid][nodeID] = float(exx) 00695 self.eyy[eid][nodeID] = float(eyy) 00696 self.ezz[eid][nodeID] = float(ezz) 00697 self.exy[eid][nodeID] = float(exy) 00698 self.eyz[eid][nodeID] = float(eyz) 00699 self.exz[eid][nodeID] = float(exz) 00700 self.e1[eid][nodeID] = float(e1) 00701 self.e2[eid][nodeID] = float(e2) 00702 self.e3[eid][nodeID] = float(e3) 00703 self.evmShear[eid][nodeID] = float(evmShear) 00704 n+=3 00705 ### 00706 ### 00707 return 00708 raise NotImplementedError() 00709 (dtName,dt) = transient 00710 self.dataCode['name'] = dtName 00711 if dt not in self.gridTypes: 00712 self.addNewTransient() 00713 00714 for line in data: 00715 (nodeID,gridType,t1,t2,t3,r1,r2,r3) = line 00716 self.gridTypes[dt][nodeID] = array([t1,t2,t3]) 00717 self.translations[dt][nodeID] = array([t1,t2,t3]) 00718 self.rotations[dt][nodeID] = array([r1,r2,r3]) 00719 ### 00720 del self.data 00721 00722 def deleteTransient(self,dt): 00723 del self.exx[dt] 00724 del self.eyy[dt] 00725 del self.ezz[dt] 00726 del self.exy[dt] 00727 del self.eyz[dt] 00728 del self.exz[dt] 00729 del self.e1[dt] 00730 del self.e2[dt] 00731 del self.e3[dt] 00732 00733 def getTransients(self): 00734 k = self.exx.keys() 00735 k.sort() 00736 return k 00737 00738 def addNewTransient(self,dt): 00739 """ 00740 initializes the transient variables 00741 """ 00742 self.nonlinearFactor = dt 00743 self.exx[dt] = {} 00744 self.eyy[dt] = {} 00745 self.ezz[dt] = {} 00746 self.exy[dt] = {} 00747 self.eyz[dt] = {} 00748 self.exz[dt] = {} 00749 self.e1[dt] = {} 00750 self.e2[dt] = {} 00751 self.e3[dt] = {} 00752 00753 #self.aCos[dt] = {} 00754 #self.bCos[dt] = {} 00755 #self.cCos[dt] = {} 00756 #self.pressure[dt] = {} 00757 self.evmShear[dt] = {} 00758 ### 00759 00760 def addNewEid(self,eType,cid,dt,eid,nodeID,exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,aCos,bCos,cCos,pressure,evm): 00761 #print "Solid Strain add..." 00762 assert eid not in self.exx 00763 assert cid >= 0 00764 assert eid >= 0 00765 self.eType[eid] = eType 00766 self.cid[eid] = cid 00767 self.exx[eid] = {nodeID: exx} 00768 self.eyy[eid] = {nodeID: eyy} 00769 self.ezz[eid] = {nodeID: ezz} 00770 self.exy[eid] = {nodeID: exy} 00771 self.eyz[eid] = {nodeID: eyz} 00772 self.exz[eid] = {nodeID: exz} 00773 self.e1[eid] = {nodeID: e1} 00774 self.e2[eid] = {nodeID: e2} 00775 self.e3[eid] = {nodeID: e3} 00776 #self.aCos[eid] = {nodeID: aCos} 00777 #self.bCos[eid] = {nodeID: bCos} 00778 #self.cCos[eid] = {nodeID: cCos} 00779 #self.pressure[eid] = {nodeID: pressure} 00780 self.evmShear[eid] = {nodeID: evm} 00781 msg = "*eid=%s nodeID=%s evmShear=%g" %(eid,nodeID,evm) 00782 #print msg 00783 if nodeID==0: raise Exception(msg) 00784 00785 def addNewEidSort1(self,eType,cid,dt,eid,nodeID,exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,aCos,bCos,cCos,pressure,evm): 00786 #print "Solid Strain add..." 00787 assert cid >= 0 00788 assert eid >= 0 00789 00790 if dt not in self.exx: 00791 self.addNewTransient(dt) 00792 00793 assert eid not in self.exx[dt],self.exx[dt] 00794 self.eType[eid] = eType 00795 self.cid[eid] = cid 00796 self.exx[dt][eid] = {nodeID: exx} 00797 self.eyy[dt][eid] = {nodeID: eyy} 00798 self.ezz[dt][eid] = {nodeID: ezz} 00799 self.exy[dt][eid] = {nodeID: exy} 00800 self.eyz[dt][eid] = {nodeID: eyz} 00801 self.exz[dt][eid] = {nodeID: exz} 00802 self.e1[dt][eid] = {nodeID: e1} 00803 self.e2[dt][eid] = {nodeID: e2} 00804 self.e3[dt][eid] = {nodeID: e3} 00805 #self.aCos[dt][eid] = {nodeID: aCos} 00806 #self.bCos[dt][eid] = {nodeID: bCos} 00807 #self.cCos[dt][eid] = {nodeID: cCos} 00808 #self.pressure[dt][eid] = {nodeID: pressure} 00809 self.evmShear[dt][eid] = {nodeID: evm} 00810 msg = "*eid=%s nodeID=%s evmShear=%g" %(eid,nodeID,evm) 00811 #print msg 00812 if nodeID==0: raise Exception(msg) 00813 00814 def add(self,dt,eid,nodeID,exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,aCos,bCos,cCos,pressure,evm): 00815 #print "***add" 00816 msg = "eid=%s nodeID=%s vm=%g" %(eid,nodeID,evm) 00817 #print msg 00818 #print self.exx 00819 #print self.fiberDistance 00820 self.exx[eid][nodeID] = exx 00821 self.eyy[eid][nodeID] = eyy 00822 self.ezz[eid][nodeID] = ezz 00823 00824 self.exy[eid][nodeID] = exy 00825 self.eyz[eid][nodeID] = eyz 00826 self.exz[eid][nodeID] = exz 00827 self.e1[eid][nodeID] = e1 00828 self.e2[eid][nodeID] = e2 00829 self.e3[eid][nodeID] = e3 00830 00831 #self.aCos[eid][nodeID] = aCos 00832 #self.bCos[eid][nodeID] = bCos 00833 #self.cCos[eid][nodeID] = cCos 00834 #self.pressure[eid][nodeID] = pressure 00835 self.evmShear[eid][nodeID] = evm 00836 00837 if nodeID==0: raise Exception(msg) 00838 00839 def addSort1(self,dt,eid,nodeID,exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,aCos,bCos,cCos,pressure,evm): 00840 #print "***add" 00841 msg = "eid=%s nodeID=%s vm=%g" %(eid,nodeID,evm) 00842 #print msg 00843 #print self.exx 00844 #print self.fiberDistance 00845 00846 self.exx[dt][eid][nodeID] = exx 00847 self.eyy[dt][eid][nodeID] = eyy 00848 self.ezz[dt][eid][nodeID] = ezz 00849 00850 self.exy[dt][eid][nodeID] = exy 00851 self.eyz[dt][eid][nodeID] = eyz 00852 self.exz[dt][eid][nodeID] = exz 00853 00854 self.e1[dt][eid][nodeID] = e1 00855 self.e2[dt][eid][nodeID] = e2 00856 self.e3[dt][eid][nodeID] = e3 00857 00858 #self.aCos[dt][eid][nodeID] = aCos 00859 #self.bCos[dt][eid][nodeID] = bCos 00860 #self.cCos[dt][eid][nodeID] = cCos 00861 #self.pressure[dt][eid][nodeID] = pressure 00862 self.evmShear[dt][eid][nodeID] = evm 00863 00864 if nodeID==0: raise Exception(msg) 00865 00866 def getHeaders(self): 00867 headers = ['exx','eyy','ezz','exy','eyz','exz'] 00868 if self.isVonMises(): 00869 headers.append('evm') 00870 else: 00871 headers.append('maxShear') 00872 return headers 00873 00874 def ovm(self,o11,o22,o33,o12,o13,o23): 00875 """http://en.wikipedia.org/wiki/Von_Mises_yield_criterion""" 00876 ovm = 0.5 * ( (o11-o22)**2+(o22-o33)**2+(o11-o33)**2+6*(o23**2+o13**2+o12**2)) 00877 return ovm 00878 00879 #def ovmPlane(self,o11,o22,o12): 00880 #"""http://en.wikipedia.org/wiki/Von_Mises_yield_criterion""" 00881 #ovm = (o11**2+o22**2-o1*o2+3*o12**2)**0.5 00882 #return ovm 00883 00884 def octahedral(self,o11,o22,o33,o12,o13,o23): 00885 """http://en.wikipedia.org/wiki/Von_Mises_yield_criterion""" 00886 ovm = self.ovm(o11,o22,o33,o12,o13,o23) 00887 return ovm/3*2**0.5 00888 00889 def pressure(self,e1,e2,e3): 00890 """ 00891 returns the hydrostatic pressure 00892 (e1+e2+e3)/-3. 00893 http://en.wikipedia.org/wiki/Stress_%28mechanics%29 00894 """ 00895 return (e1+e2+e3)/-3. 00896 00897 def directionalVectors(self,exx,eyy,ezz,exy,eyz,exz): 00898 A = [[exx,exy,exz], 00899 [exy,eyy,eyz], 00900 [exz,eyz,ezz]] 00901 (Lambda,v) = eigh(A) # a hermitian matrix is a symmetric-real matrix 00902 return v 00903 00904 def getF06_Header(self): 00905 if self.isVonMises(): 00906 vonMises = 'VON MISES' 00907 else: 00908 vonMises = 'MAX SHEAR' 00909 00910 tetraMsg = [' S T R A I N S I N T E T R A H E D R O N S O L I D E L E M E N T S ( C T E T R A )\n', 00911 '0 CORNER ------CENTER AND CORNER POINT STRAINS--------- DIR. COSINES MEAN\n', 00912 ' ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE %s \n' %(vonMises)] 00913 00914 pentaMsg = [' S T R A I N S I N P E N T A H E D R O N S O L I D E L E M E N T S ( P E N T A )\n', 00915 '0 CORNER ------CENTER AND CORNER POINT STRAINS--------- DIR. COSINES MEAN\n', 00916 ' ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE %s \n' %(vonMises)] 00917 00918 hexaMsg = [' S T R A I N S I N H E X A H E D R O N S O L I D E L E M E N T S ( H E X A )\n', 00919 '0 CORNER ------CENTER AND CORNER POINT STRAINS--------- DIR. COSINES MEAN\n', 00920 ' ELEMENT-ID GRID-ID NORMAL SHEAR PRINCIPAL -A- -B- -C- PRESSURE %s \n' %(vonMises)] 00921 00922 00923 #nNodes = {'CTETRA':4,'CPENTA':6,'CHEXA':8,'HEXA':8,'PENTA':6,'TETRA':4,} 00924 tetraEids = [] 00925 hexaEids = [] 00926 pentaEids = [] 00927 for eid,eType in sorted(self.eType.iteritems()): 00928 if eType=='CTETRA' or eType=='TETRA': 00929 tetraEids.append(eid) 00930 elif eType=='CPENTA' or eType=='PENTA': 00931 pentaEids.append(eid) 00932 elif eType=='CHEXA' or eType=='HEXA': 00933 hexaEids.append(eid) 00934 else: 00935 raise NotImplementedError('eType=|%r|' %(eType)) 00936 ### 00937 ### 00938 return (tetraMsg,pentaMsg,hexaMsg,tetraEids,hexaEids,pentaEids) 00939 00940 def writeF06(self,header,pageStamp,pageNum=1,f=None,isMagPhase=False): 00941 if self.nonlinearFactor is not None: 00942 return self.writeF06Transient(header,pageStamp,pageNum,f) 00943 msg = [] 00944 00945 (tetraMsg,pentaMsg,hexaMsg,tetraEids,hexaEids,pentaEids) = self.getF06_Header() 00946 #nNodes = {'CTETRA':4,'CPENTA':6,'CHEXA':8,'HEXA':8,'PENTA':6,'TETRA':4,} 00947 if tetraEids: 00948 msg += self.writeElement('CTETRA',4,tetraEids,header,tetraMsg,f) 00949 msg.append(pageStamp+str(pageNum)+'\n') 00950 pageNum+=1 00951 if hexaEids: 00952 msg += self.writeElement('CHEXA',8,hexaEids,header,hexaMsg,f) 00953 msg.append(pageStamp+str(pageNum)+'\n') 00954 pageNum+=1 00955 if pentaEids: 00956 msg += self.writeElement('CPENTA',6,pentaEids,header,pentaMsg,f) 00957 msg.append(pageStamp+str(pageNum)+'\n') 00958 pageNum+=1 00959 ### 00960 return (''.join(msg),pageNum-1) 00961 00962 def writeF06Transient(self,header,pageStamp,pageNum=1,f=None,isMagPhase=False): 00963 msg = [] 00964 (tetraMsg,pentaMsg,hexaMsg,tetraEids,hexaEids,pentaEids) = self.getF06_Header() 00965 dts = self.exx.keys() 00966 for dt in dts: 00967 if tetraEids: 00968 msg += self.writeElementTransient('CTETRA',4,tetraEids,dt,header,tetraMsg,f) 00969 msg.append(pageStamp+str(pageNum)+'\n') 00970 pageNum+=1 00971 if hexaEids: 00972 msg += self.writeElementTransient('CHEXA',8,hexaEids, dt,header,hexaMsg,f) 00973 msg.append(pageStamp+str(pageNum)+'\n') 00974 pageNum+=1 00975 if pentaEids: 00976 msg += self.writeElementTransient('CPENTA',6,pentaEids,dt,header,pentaMsg,f) 00977 msg.append(pageStamp+str(pageNum)+'\n') 00978 pageNum+=1 00979 ### 00980 return (''.join(msg),pageNum-1) 00981 00982 def writeElement(self,eType,nNodes,eids,header,tetraMsg,f=None): 00983 msg = header+tetraMsg 00984 for eid in eids: 00985 #eType = self.eType[eid] 00986 00987 k = self.exx[eid].keys() 00988 k.remove('C') 00989 k.sort() 00990 msg.append('0 %8s 0GRID CS %i GP\n' %(eid,nNodes) ) 00991 for nid in ['C']+k: 00992 exx = self.exx[eid][nid] 00993 eyy = self.eyy[eid][nid] 00994 ezz = self.ezz[eid][nid] 00995 exy = self.exy[eid][nid] 00996 eyz = self.eyz[eid][nid] 00997 exz = self.exz[eid][nid] 00998 00999 e1 = self.e1[eid][nid] 01000 e2 = self.e2[eid][nid] 01001 e3 = self.e3[eid][nid] 01002 evm = self.evmShear[eid][nid] 01003 p = (e1+e2+e3)/-3. 01004 01005 if nid=='C': nid='CENTER' 01006 #print "nid = |%r|" %(nid) 01007 A = [[exx,exy,exz], 01008 [exy,eyy,eyz], 01009 [exz,eyz,ezz]] 01010 (Lambda,v) = eigh(A) # a hermitian matrix is a symmetric-real matrix 01011 01012 ([exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,p,evm],isAllZeros) = self.writeFloats13E([exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,p,evm]) 01013 msg.append('0 %8s X %13s XY %13s A %13s LX%5.2f%5.2f%5.2f %13s %-s\n' %(nid,exx,exy,e1,v[0,1],v[0,2],v[0,0],p,evm.strip())) 01014 msg.append(' %8s Y %13s YZ %13s B %13s LY%5.2f%5.2f%5.2f\n' %('', eyy,eyz,e2,v[1,1],v[1,2],v[1,0])) 01015 msg.append(' %8s Z %13s ZX %13s C %13s LZ%5.2f%5.2f%5.2f\n' %('', ezz,exz,e3,v[2,1],v[2,2],v[2,0])) 01016 ### 01017 ### 01018 if f is not None: 01019 f.write(''.join(msg)) 01020 msg = [''] 01021 ### 01022 return ''.join(msg) 01023 01024 def writeElementTransient(self,eType,nNodes,eids,dt,header,tetraMsg,pageStamp=1,f=None): 01025 dtLine = '%14s = %12.5E\n'%(self.dataCode['name'],dt) 01026 header[1] = dtLine 01027 msg = header+tetraMsg 01028 for eid in eids: 01029 #eType = self.eType[eid] 01030 01031 k = self.exx[dt][eid].keys() 01032 k.remove('C') 01033 k.sort() 01034 msgA = '0 %8s 0GRID CS %i GP\n' %(eid,nNodes) 01035 for nid in ['C']+k: 01036 exx = self.exx[dt][eid][nid] 01037 eyy = self.eyy[dt][eid][nid] 01038 ezz = self.ezz[dt][eid][nid] 01039 exy = self.exy[dt][eid][nid] 01040 eyz = self.eyz[dt][eid][nid] 01041 exz = self.exz[dt][eid][nid] 01042 01043 e1 = self.e1[dt][eid][nid] 01044 e2 = self.e2[dt][eid][nid] 01045 e3 = self.e3[dt][eid][nid] 01046 evm = self.evmShear[dt][eid][nid] 01047 p = (e1+e2+e3)/-3. 01048 01049 if nid=='C': nid='CENTER' 01050 #print "nid = |%r|" %(nid) 01051 A = [[exx,exy,exz], 01052 [exy,eyy,eyz], 01053 [exz,eyz,ezz]] 01054 (Lambda,v) = eigh(A) # a hermitian matrix is a symmetric-real matrix 01055 01056 ([exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,p,evm],isAllZeros) = self.writeFloats13E([exx,eyy,ezz,exy,eyz,exz,e1,e2,e3,p,evm]) 01057 msgA += '0 %8s X %13s XY %13s A %13s LX%5.2f%5.2f%5.2f %13s %-s\n' %(nid,exx,exy,e1,v[0,1],v[0,2],v[0,0],p,evm.strip()) 01058 msgA += ' %8s Y %13s YZ %13s B %13s LY%5.2f%5.2f%5.2f\n' %('', eyy,eyz,e2,v[1,1],v[1,2],v[1,0]) 01059 msgA += ' %8s Z %13s ZX %13s C %13s LZ%5.2f%5.2f%5.2f\n' %('', ezz,exz,e3,v[2,1],v[2,2],v[2,0]) 01060 ### 01061 ### 01062 if f is not None: 01063 f.write(''.join(msg)) 01064 msg = [''] 01065 ### 01066 return ''.join(msg) 01067 01068 def __repr__(self): 01069 if self.nonlinearFactor is not None: 01070 return self.__reprTransient__() 01071 01072 msg = '---SOLID STRESS---\n' 01073 headers = self.getHeaders() 01074 msg += '%-6s %6s %8s ' %('EID','eType','nodeID') 01075 for header in headers: 01076 msg += '%9s ' %(header) 01077 msg += '\n' 01078 for eid,oxxNodes in sorted(self.oxx.iteritems()): 01079 eType = self.eType[eid] 01080 for nid in sorted(oxxNodes): 01081 oxx = self.oxx[eid][nid] 01082 oyy = self.oyy[eid][nid] 01083 ozz = self.ozz[eid][nid] 01084 txy = self.txy[eid][nid] 01085 tyz = self.tyz[eid][nid] 01086 txz = self.txz[eid][nid] 01087 01088 #o1 = self.o1[eid][nid] 01089 #o2 = self.o2[eid][nid] 01090 #o3 = self.o3[eid][nid] 01091 ovm = self.ovmShear[eid][nid] 01092 msg += '%-6i %6s %8s ' %(eid,eType,nid) 01093 vals = [oxx,oyy,ozz,txy,tyz,txz,ovm] 01094 for val in vals: 01095 if abs(val)<1e-6: 01096 msg += '%9s ' %('0') 01097 else: 01098 msg += '%9i ' %(val) 01099 ### 01100 msg += '\n' 01101 #msg += "eid=%-4s eType=%-6s nid=%-2i oxx=%-5i oyy=%-5i ozz=%-5i txy=%-5i tyz=%-5i txz=%-5i ovm=%-5i\n" %(eid,eType,nid,oxx,oyy,ozz,txy,tyz,txz,ovm) 01102 ### 01103 ### 01104 return msg 01105 01106 def __repr__(self): 01107 if self.nonlinearFactor is not None: 01108 return self.__reprTransient__() 01109 01110 msg = '---SOLID STRAIN---\n' 01111 headers = self.getHeaders() 01112 msg += '%-6s %6s %8s ' %('EID','eType','nodeID') 01113 for header in headers: 01114 msg += '%10s ' %(header) 01115 msg += '\n' 01116 for eid,exxNodes in sorted(self.exx.iteritems()): 01117 eType = self.eType[eid] 01118 for nid in sorted(exxNodes): 01119 exx = self.exx[eid][nid] 01120 eyy = self.eyy[eid][nid] 01121 ezz = self.ezz[eid][nid] 01122 exy = self.exy[eid][nid] 01123 eyz = self.eyz[eid][nid] 01124 exz = self.exz[eid][nid] 01125 evm = self.evmShear[eid][nid] 01126 msg += '%-6i %6s %8s ' %(eid,eType,nid) 01127 vals = [exx,eyy,ezz,exy,eyz,exz,evm] 01128 for val in vals: 01129 if abs(val)<1e-6: 01130 msg += '%10s ' %('0') 01131 else: 01132 msg += '%10.3e ' %(val) 01133 ### 01134 msg += '\n' 01135 #msg += "eid=%-4s eType=%-6s nid=%-2i exx=%-5i eyy=%-5i ezz=%-5i exy=%-5i eyz=%-5i exz=%-5i evm=%-5i\n" %(eid,eType,nid,exx,eyy,ezz,exy,eyz,exz,evm) 01136 ### 01137 ### 01138 return msg 01139 01140 def __reprTransient__(self): 01141 #return '' 01142 msg = ['---SOLID STRAIN---\n'] 01143 headers = self.getHeaders() 01144 msg2 = '%-6s %6s %8s ' %('EID','eType','nodeID') 01145 for header in headers: 01146 msg2 += '%10s ' %(header) 01147 msg2 += '\n' 01148 msg.append(msg2) 01149 for dt,exxs in sorted(self.exx.iteritems()): 01150 msg.append('%s = %g\n' %(self.dataCode['name'],dt)) 01151 for eid,exxNodes in sorted(exxs.iteritems()): 01152 eType = self.eType[eid] 01153 msg2 = '' 01154 for nid in sorted(exxNodes): 01155 exx = self.exx[dt][eid][nid] 01156 eyy = self.eyy[dt][eid][nid] 01157 ezz = self.ezz[dt][eid][nid] 01158 exy = self.exy[dt][eid][nid] 01159 eyz = self.eyz[dt][eid][nid] 01160 exz = self.exz[dt][eid][nid] 01161 evm = self.evmShear[dt][eid][nid] 01162 msg2 += '%-6i %6s %8s ' %(eid,eType,nid) 01163 vals = [exx,eyy,ezz,exy,eyz,exz,evm] 01164 for val in vals: 01165 if abs(val)<1e-6: 01166 msg2 += '%10s ' %('0') 01167 else: 01168 msg2 += '%10.3e ' %(val) 01169 ### 01170 msg2 += '\n' 01171 msg.append(msg2) 01172 #msg += "eid=%-4s eType=%-6s nid=%-2i exx=%-5i eyy=%-5i ezz=%-5i exy=%-5i eyz=%-5i exz=%-5i evm=%-5i\n" %(eid,eType,nid,exx,eyy,ezz,exy,eyz,exz,evm) 01173 ### 01174 ### 01175 ### 01176 return ''.join(msg)