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 # pylint: disable=C0103,R0902,R0904,R0914 00026 00027 from __future__ import division, print_function 00028 from itertools import izip 00029 from math import log 00030 #from math import (sin,sinh,cos,cosh,tan,tanh,sqrt,atan,atan2,acosh,acos,asin, 00031 # asinh,atanh) #,atanh2 # going to be used by DEQATN 00032 00033 from numpy import zeros #average 00034 from scipy.sparse import coo_matrix 00035 00036 from pyNastran.bdf.cards.baseCard import BaseCard 00037 00038 def ssq(*listA): 00039 """ 00040 sum of squares 00041 @note used for DEQATN 00042 """ 00043 out = 0. 00044 for x in listA: 00045 out += x*x 00046 return out 00047 00048 def sum2(*listA): 00049 """ 00050 sum of listA 00051 @note used for DEQATN 00052 """ 00053 return sum(listA) 00054 00055 def mod(x,y): 00056 """ 00057 x%y 00058 @note used for DEQATN 00059 """ 00060 return x%y 00061 00062 def logx(x,y): 00063 """ 00064 log base x of y 00065 @note used for DEQATN 00066 """ 00067 log(y,x) 00068 00069 def dim(x,y): 00070 """ 00071 @note used for DEQATN 00072 """ 00073 return x-min(x,y) 00074 00075 def db(p,pref): 00076 """ 00077 sound pressure in decibels 00078 would capitalize it, but you wouldnt be able to call the function... 00079 """ 00080 return 20.*log(p/pref) 00081 00082 class DEQATN(BaseCard):# needs work... 00083 type = 'DEQATN' 00084 def __init__(self,card=None,data=None): 00085 newCard = '' 00086 foundNone = False 00087 for field in card.card: 00088 if foundNone is False and field is not None: 00089 newCard += field+',' 00090 foundNone = True 00091 elif foundNone is True and field is not None: 00092 newCard += field 00093 00094 #if len(card.card)>1: 00095 #print "card.card = ",card.card 00096 #line0 = ','.join(card.card) 00097 #else: 00098 #line0 = ''.join(card.card) 00099 line0 = newCard 00100 self.eqID = line0[8:16] 00101 00102 assert len(self.eqID)==8,'len(eqID)==%s' %(len(self.eqID)) 00103 eq = line0[16:] 00104 eq = eq.replace(' ','').lower() 00105 (self.name,self.eq) = eq.split('=') 00106 #print("EQ = %s" %(self.eq)) 00107 00108 def evaluate(self,args): 00109 #eqLow = self.eq.lower() 00110 #eval(self.eq) 00111 pass 00112 00113 def __repr__(self): 00114 eq = self.name+'='+self.eq 00115 eqLine = eq[0:56] 00116 eq = eq[56:] 00117 fields = ['DEQATN ','%8s'%(self.eqID), eqLine] 00118 00119 if len(eq): 00120 eqLine = eq[0:72] 00121 eq = eq[72:] 00122 fields += [' '+eqLine] 00123 return ''.join(fields) 00124 00125 class NastranMatrix(BaseCard): 00126 """ 00127 Base class for the DMIG, DMIJ, DMIJI, DMIK matrices 00128 """ 00129 def __init__(self, card=None, data=None): 00130 self.name = card.field(1) 00131 #zero 00132 00133 ## 4-Lower Triangular; 5=Upper Triangular; 6=Symmetric; 8=Identity (m=nRows, n=m) 00134 self.ifo = int(card.field(3)) 00135 ## 1-Real, Single Precision; 2=Real,Double Precision; 3=Complex, Single; 4=Complex, Double 00136 self.tin = int(card.field(4)) 00137 ## 0-Set by cell precision 00138 self.tout = int(card.field(5,0)) 00139 00140 self.polar = card.field(6) 00141 self.ncol = card.field(8) 00142 00143 self.GCj = [] 00144 self.GCi = [] 00145 self.Real = [] 00146 if self.isComplex(): 00147 self.Complex = [] 00148 00149 def writeCodeAster(self): 00150 """ 00151 assume set 1 = MAAX1,MAAX2, etc. and 100/n % on each 00152 """ 00153 00154 # for real combination 00155 comm = 'K_Mtx_AB=COMB_MATR_ASSE(COMB_R=(\n' 00156 comm += ' _F(MATR_ASSE = K_Mtx_A,COEF_R = 1.),\n' 00157 comm += ' _F(MATR_ASSE = K_Mtx_B,COEF_R = 1.)));\n' 00158 00159 # for complex combination 00160 00161 comm += "K_Mtx_AB=COMB_MATR_ASSE(COMB_C=(\n" 00162 comm += "_F(MATR_ASSE=K_Mtx_A,COEF_C=('RI',0.7,0.3,),)\n" 00163 comm += "_F(MATR_ASSE=K_Mtx_B,COEF_C=('RI',0.7,0.3,),),),);\n" 00164 comm = 'K_Mtx=ASSE_MATRICE(MATR_ELEM=ElMtx_K,NUME_DDL=%s,);' 00165 return comm 00166 00167 def addColumn(self,card=None,data=None): 00168 #print "hi column" 00169 Gj = card.field(2) 00170 Cj = card.field(3) 00171 00172 nFields = card.nFields() 00173 nLoops = (nFields-5)//4 00174 minLoops = nLoops-1 00175 if minLoops<=0: 00176 minLoops = 1 00177 #assert nFields <= 8,'nFields=%s' %(nFields) 00178 00179 #print "minLoops = ",minLoops 00180 #print "nLoops = ",nLoops 00181 for i in xrange(minLoops): 00182 self.GCj.append((Gj,Cj)) 00183 00184 if self.isComplex(): 00185 for i in xrange(minLoops): 00186 n = 5+4*i 00187 Gi = card.field(n) 00188 Ci = card.field(n+1) 00189 self.GCi.append((Gi,Ci)) 00190 self.Real.append(card.field(n+2)) 00191 self.Complex.append(card.field(n+3)) 00192 else: 00193 for i in xrange(minLoops): 00194 n = 5+4*i 00195 Gi = card.field(n) 00196 Ci = card.field(n+1) 00197 self.GCi.append((Gi,Ci)) 00198 self.Real.append(card.field(n+2)) 00199 00200 assert len(self.GCj)==len(self.GCi),'(len(GCj)=%s len(GCi)=%s' %(len(self.GCj),len(self.GCi)) 00201 #if self.isComplex(): 00202 #self.Complex(card.field(v) 00203 00204 def getMatrix(self,isSparse=False): 00205 """ 00206 builds the Matrix 00207 @param self the object pointer 00208 @param isSparse should the matrix be returned as a sparse matrix (default=True). Slower for dense matrices. 00209 @retval M the matrix 00210 @retval rows dictionary of keys=rowID, values=(Grid,Component) for the matrix 00211 @retval cols dictionary of keys=columnID, values=(Grid,Component) for the matrix 00212 @warning isSparse WILL fail 00213 """ 00214 i=0 00215 rows = {} 00216 rowsReversed = {} 00217 for GCi in self.GCi: 00218 if GCi not in rows: 00219 rows[GCi] = i 00220 rowsReversed[i] = GCi 00221 i+=1 00222 #nRows = len(rows2) 00223 00224 j=0 00225 cols = {} 00226 colsReversed = {} 00227 for GCj in self.GCj: 00228 if GCj not in cols: 00229 cols[GCj] = j 00230 colsReversed[j] = GCj 00231 j+=1 00232 #nCols = len(cols2) 00233 00234 #A = ss.lil_matrix((3,3), dtype='d') # double precision 00235 00236 #rows=[]; cols=[]; data=[] 00237 #for i in xrange(3): 00238 # for j in xrange(3): 00239 # k = float((i+1)*(j+1)) 00240 # rows.append(i) 00241 # cols.append(j) 00242 # data.append(k) 00243 # A[i,j] = k 00244 # ### 00245 ### 00246 00247 #isSparse = False 00248 if isSparse: 00249 data = [] 00250 rows2 = [] 00251 cols2 = [] 00252 00253 nrows = 0 00254 ncols = 0 00255 if self.isComplex(): 00256 for (GCj, GCi, reali, complexi) in izip(self.GCj, self.GCi, self.Real, self.Complex): 00257 i = rows[GCi] 00258 j = cols[GCj] 00259 nrows = max(i, nrows) 00260 ncols = max(j, ncols) 00261 rows2.append(i) 00262 cols2.append(j) 00263 data.append(complex(reali, complexi)) 00264 else: 00265 for (GCj, GCi) in izip(self.GCj, self.GCi): 00266 i = rows[GCi] 00267 j = cols[GCj] 00268 nrows = max(i,nrows) 00269 ncols = max(j,ncols) 00270 rows2.append(i) 00271 cols2.append(j) 00272 data = self.Real 00273 #print("i=%s j=%s len(rows2)=%s len(cols2)=%s len(data)=%s" 00274 # % (i, j, len(self.GCi), len(self.GCj), len(data))) 00275 # ,dtype=Format 00276 #print(rows2) 00277 00278 #print("nrows=%s ncols=%s" %(nrows, ncols)) 00279 if self.ifo in [1,6]: 00280 nrows = max(nrows,ncols) 00281 ncols = nrows 00282 #print("nrows=%s ncols=%s" %(nrows, ncols)) 00283 00284 dType = self.getDType(self.tin) 00285 #A = coo_matrix( (entries,(rows,cols)),shape=(nrows,ncols),dtype=dType) # test 00286 M = coo_matrix( (data, (self.GCi, self.GCj)), shape=(nrows,ncols), dtype=dType) 00287 #M = coo_matrix( (data,(self.GCi,self.GCj)),shape=(i,j)) # old 00288 #M = coo_matrix( (data,(self.GCi,self.GCj)),shape=(nrows,ncols)) 00289 #print M.todense() 00290 #print M 00291 else: 00292 if self.isComplex(): 00293 M = zeros((i,j), dtype='complex128') 00294 for (GCj, GCi, reali, complexi) in izip(self.GCj, self.GCi, self.Real, self.Complex): 00295 i = rows[GCi] 00296 j = cols[GCj] 00297 M[i,j] = complex(reali,complexi) 00298 else: 00299 M = zeros((i,j), dtype='float64') 00300 for (GCj, GCi, reali) in izip(self.GCj, self.GCi, self.Real): 00301 i = rows[GCi] 00302 j = cols[GCj] 00303 M[i,j] = reali 00304 ### 00305 ### 00306 #print(M) 00307 return (M,rowsReversed,colsReversed) 00308 00309 def rename(self,newName): 00310 self.name = newName 00311 00312 def isReal(self): 00313 return not self.isComplex() 00314 00315 def isComplex(self): 00316 if self.tin in [3,4]: 00317 return True 00318 return False 00319 00320 def getDType(self,Type): 00321 if Type==1: 00322 dType = 'float32' 00323 elif Type==2: 00324 dType = 'float64' 00325 elif Type==3: 00326 dType = 'complex64' 00327 elif Type==4: 00328 dType = 'complex128' 00329 elif Type==0: 00330 if self.isComplex(): 00331 dType = 'complex128' 00332 else: 00333 dType = 'float64' 00334 else: 00335 raise RuntimeError("invalid option for matrix format") 00336 return dType 00337 00338 def __repr__(self): 00339 """ 00340 @todo support double precision 00341 """ 00342 msg = '\n$'+'-'*80 00343 msg += '\n$ %s Matrix %s\n' %(self.type, self.name) 00344 fields = [self.type, self.name, 0, self.ifo, self.tin, self.tout, self.polar, None, self.ncol] 00345 msg += self.printCard(fields) 00346 00347 if self.isComplex(): 00348 for (GCi, GCj, reali, imagi) in izip(self.GCi, self.GCj, self.Real, self.Complex): 00349 fields = [self.type, self.name, GCj[0], GCj[1], None, GCi[0], GCi[1], reali, imagi] 00350 msg += self.printCard(fields) 00351 else: 00352 for (GCi, GCj, reali) in izip(self.GCi, self.GCj, self.Real): 00353 fields = [self.type, self.name, GCj[0], GCj[1], None, GCi[0], GCi[1], reali, None] 00354 msg += self.printCard(fields) 00355 return msg 00356 00357 00358 class DMIG(NastranMatrix): 00359 """ 00360 Defines direct input matrices related to grid, extra, and/or scalar points. 00361 The matrix is defined by a single header entry and one or more column 00362 entries. A column entry is required for each column with nonzero elements. 00363 """ 00364 type = 'DMIG' 00365 def __init__(self, card=None, data=None): 00366 NastranMatrix.__init__(self, card, data) 00367 00368 class DMIJ(NastranMatrix): 00369 """ 00370 Direct Matrix Input at js-Set of the Aerodynamic Mesh 00371 Defines direct input matrices related to collation degrees-of-freedom 00372 (js-set) of aerodynamic mesh points for CAERO1, CAERO3, CAERO4 and CAERO5 00373 and for the slender body elements of CAERO2. These include W2GJ, FA2J and 00374 input pressures and downwashes associated with AEPRESS and AEDW entries. 00375 The matrix is described by a single header entry and one or more column 00376 entries. A column entry is required for each column with nonzero elements. 00377 For entering data for the interference elements of a CAERO2, use DMIJI 00378 or DMI. 00379 """ 00380 type = 'DMIJ' 00381 def __init__(self,card=None,data=None): 00382 NastranMatrix.__init__(self,card,data) 00383 00384 class DMIJI(NastranMatrix): 00385 """ 00386 Direct Matrix Input at js-Set of the Interference Body 00387 Defines direct input matrices related to collation degrees-of-freedom 00388 (js-set) of aerodynamic mesh points for the interference elements of CAERO2. 00389 These include W2GJ, FA2J and input pressures and downwashes associated with 00390 AEPRESS and AEDW entries. The matrix is described by a single header entry 00391 and one or more column entries. A column entry is required for each column 00392 with nonzero elements. For entering data for the slender elements of a 00393 CAERO2, or a CAERO1, 3, 4 or 5 use DMIJ or DMI. 00394 """ 00395 type = 'DMIJI' 00396 def __init__(self, card=None, data=None): 00397 NastranMatrix.__init__(self, card, data) 00398 00399 class DMIK(NastranMatrix): 00400 """ 00401 Direct Matrix Input at ks-Set of the Aerodynamic Mesh 00402 Defines direct input matrices related to physical (displacement) 00403 degrees-of-freedom (ks-set) of aerodynamic grid points. These include WKK, 00404 WTFACT and input forces associated with AEFORCE entries. The matrix is 00405 described by a single header entry and one or more column entries. A column 00406 entry is required for each column with nonzero elements. 00407 """ 00408 type = 'DMIK' 00409 def __init__(self,card=None,data=None): 00410 NastranMatrix.__init__(self,card,data) 00411 00412 class DMI(BaseCard): 00413 type = 'DMI' 00414 def __init__(self,card=None,data=None): 00415 self.name = card.field(1) 00416 #zero 00417 00418 ## Form of the matrix: 1=Square (not symmetric); 2=Rectangular; 00419 ## 3=Diagonal (m=nRows,n=1); 4=Lower Triangular; 5=Upper Triangular; 00420 ## 6=Symmetric; 8=Identity (m=nRows, n=m) 00421 self.form = int(card.field(3)) 00422 00423 ## 1-Real, Single Precision; 2=Real,Double Precision; 00424 ## 3=Complex, Single; 4=Complex, Double 00425 self.tin = int(card.field(4)) 00426 00427 ## 0-Set by cell precision 00428 self.tout = int(card.field(5, 0)) 00429 00430 self.nRows = int(card.field(7)) 00431 self.nCols = int(card.field(8)) 00432 00433 self.GCj = [] 00434 self.GCi = [] 00435 self.Real = [] 00436 00437 if self.isComplex(): 00438 self.Complex = [] 00439 00440 def addColumn(self, card=None, data=None): 00441 00442 if not self.isComplex(): # real 00443 self.readReal(card) 00444 00445 def readReal(self, card): 00446 ## column number 00447 j = card.field(2) 00448 00449 # counter 00450 i = 0 00451 fields = card.fields(3) 00452 00453 # Real, starts at A(i1,j), goes to A(i2,j) in a column 00454 while i<len(fields): 00455 i1 = fields[i] 00456 #print "i1 = ",i1 00457 if isinstance(i1,int): 00458 i+=1 00459 isDoneReadingFloats = False 00460 while not isDoneReadingFloats and i<len(fields): 00461 #print "i=%s len(fields)=%s" %(i,len(fields)) 00462 realValue = fields[i] 00463 if isinstance(realValue, int): 00464 isDoneReadingFloats = True 00465 elif isinstance(realValue, float): 00466 self.GCj.append(j) 00467 self.GCi.append(i1) 00468 self.Real.append(realValue) 00469 #print "i=%s j=%s value=%s" %(i1,j,realValue) 00470 i+=1 00471 else: 00472 #print "*i=%s j=%s value=%s type=%s" %(i1,j,realValue,type(realValue)) 00473 realValue = self.Real[-1] 00474 #print "*i=%s j=%s value=%s" %(i1,j,realValue) 00475 endI = fields[i+1] 00476 #print "*i=%s endI=%s j=%s value=%s" %(i1,endI,j,realValue) 00477 for ii in xrange(i1, endI+1): 00478 self.GCj.append(j) 00479 self.GCi.append(ii) 00480 self.Real.append(realValue) 00481 ### 00482 #print "i = ",3+i 00483 #print 'field i=',fields[i] 00484 i+=1 00485 isDoneReadingFloats = True 00486 ### 00487 ### 00488 00489 def readComplex(self,card): 00490 msg = 'complex matrices not supported in the DMI reader...' 00491 raise NotImplementedError(msg) 00492 ## column number 00493 j = card.field(2) 00494 00495 # counter 00496 i = 0 00497 fields = card.fields(3) 00498 00499 # Complex, starts at A(i1,j)+imag*A(i1,j), goes to A(i2,j) in a column 00500 while i<len(fields): 00501 i1 = fields[i] 00502 i+=1 00503 isDoneReadingFloats = False 00504 asdf 00505 while not isDoneReadingFloats and i<len(fields): 00506 #print("i=%s len(fields)=%s" %(i, len(fields))) 00507 realValue = fields[i] 00508 if isinstance(floatValue, int): 00509 isDoneReadingFloats = True 00510 elif isinstance(realValue, float): 00511 complexValue = fields[i+1] 00512 self.GCj.append(j) 00513 self.GCi.append(i1) 00514 self.Real.append(realValue) 00515 self.Complex.append(complexValue) 00516 i+=2 00517 else: 00518 asdf 00519 ### 00520 ### while floats 00521 ### while fields 00522 00523 def rename(self,newName): 00524 self.name = newName 00525 00526 def isReal(self): 00527 return not self.isComplex() 00528 00529 def isComplex(self): 00530 if self.tin in [3,4]: 00531 return True 00532 return False 00533 00534 def __repr__(self): 00535 """ 00536 @todo support shortened output format. There's a stupidly low 1000 00537 DMI cap, I assume this is entries and not matrices. 00538 @todo support double precision 00539 """ 00540 msg = '\n$'+'-'*80 00541 msg += '\n$ %s Matrix %s\n' %(self.type,self.name) 00542 fields = [self.type, self.name, 0, self.form, self.tin, self.tout, None, self.nRows, self.nCols] 00543 msg += self.printCard(fields) 00544 #msg += self.printCard(fields,size=16,isD=False) 00545 00546 if self.isComplex(): 00547 for (GCi, GCj, reali, imagi) in izip(self.GCi, self.GCj, self.Real, self.Complex): 00548 fields = [self.type, self.name, GCj, GCi, reali, imagi] 00549 msg += self.printCard(fields) 00550 else: 00551 for (GCi, GCj, reali) in izip(self.GCi, self.GCj, self.Real): 00552 fields = [self.type, self.name, GCj, GCi, reali] 00553 msg += self.printCard(fields) 00554 return msg 00555