Package solvcon :: Package io :: Module genesis
[hide private]
[frames] | no frames]

Source Code for Module solvcon.io.genesis

  1  # -*- coding: UTF-8 -*- 
  2  # 
  3  # Copyright (C) 2011 Yung-Yu Chen <yyc@solvcon.net>. 
  4  # 
  5  # This program is free software; you can redistribute it and/or modify 
  6  # it under the terms of the GNU General Public License as published by 
  7  # the Free Software Foundation; either version 2 of the License, or 
  8  # (at your option) any later version. 
  9  # 
 10  # This program is distributed in the hope that it will be useful, 
 11  # but WITHOUT ANY WARRANTY; without even the implied warranty of 
 12  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 
 13  # GNU General Public License for more details. 
 14  # 
 15  # You should have received a copy of the GNU General Public License along 
 16  # with this program; if not, write to the Free Software Foundation, Inc., 
 17  # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. 
 18   
 19  """ 
 20  Adapter to Genesis/ExodusII format. 
 21  """ 
 22   
 23  from .core import FormatIO 
 24  from .netcdf import NetCDF 
 25   
26 -class Genesis(NetCDF):
27 """ 28 Model of Genesis/ExodusII file, based on netCDF reader wrapper. The 29 Genesis/ExodusII file must meet the following criteria: (i) have at least 30 one block, (ii) each block needs a name, and (iii) each sideset (BC) needs 31 a name. 32 33 @ivar ndim: dimension (2 or 3). 34 @itype ndim: int 35 @ivar nnode: number of nodes. 36 @itype nnode: int 37 @ivar ncell: number of cells/elements. 38 @itype ncell: int 39 @ivar blks: list of tuples of (name, type_name, clnds) for each Genesis 40 block. 41 @itype blks: list 42 @ivar bcs: list of tuples of (name, elem, side) for each BC (Genesis 43 sideset). 44 @itype bcs: list 45 @ivar ndcrd: coordiate array. 46 @itype ndcrd: numpy.ndarray 47 """
48 - def __init__(self, *arg, **kw):
49 super(Genesis, self).__init__(*arg, **kw) 50 # shape data. 51 self.ndim = None 52 self.nnode = None 53 self.ncell = None 54 # blocks and BCs. 55 self.blks = None 56 self.bcs = None 57 # coordinate. 58 self.ndcrd = None 59 # mapper. 60 self.emap = None
61
62 - def get_attr_text(self, name, varname):
63 """ 64 Get the attribute text attached to an variable. 65 66 @param name: name of the attribute. 67 @type name: str 68 @param varname: name of the variable. 69 @type varname: str 70 @return: the text. 71 @rtype: str 72 """ 73 from ctypes import POINTER, c_int, c_char, byref, create_string_buffer 74 # get value ID. 75 varid = c_int() 76 retval = self.nc_inq_varid(self.ncid, varname, byref(varid)) 77 if retval != self.NC_NOERR: 78 raise IOError(self.nc_strerror(retval)) 79 # get text. 80 slen = self.get_dim('len_string') 81 buf = create_string_buffer(slen) 82 retval = self.nc_get_att_text(self.ncid, varid, name, buf) 83 if retval != self.NC_NOERR: 84 raise IOError(self.nc_strerror(retval)) 85 return buf.value
86
87 - def load(self):
88 """ 89 Load mesh data. 90 91 @return: nothing 92 """ 93 from ctypes import c_int, byref 94 from numpy import hstack 95 # meta data. 96 self.ndim = ndim = self.get_dim('num_dim') 97 self.nnode = nnode = self.get_dim('num_nodes') 98 self.ncell = self.get_dim('num_elem') 99 # blocks. 100 nblk = self.get_dim('num_el_blk') 101 slen = self.get_dim('len_string') 102 self.blks = self.get_lines('eb_names', (nblk, slen)) 103 for iblk in range(nblk): 104 ncell = self.get_dim('num_el_in_blk%d' % (iblk+1)) 105 clnnd = self.get_dim('num_nod_per_el%d' % (iblk+1)) 106 clnds = self.get_array('connect%d' % (iblk+1), 107 (ncell, clnnd), 'int32') 108 type_name = self.get_attr_text('elem_type', 'connect%d' % (iblk+1)) 109 self.blks[iblk] = (self.blks[iblk], type_name, clnds) 110 # BCs. 111 nbc = self.get_dim('num_side_sets') 112 self.bcs = self.get_lines('ss_names', (nbc, slen)) 113 for ibc in range(nbc): 114 nface = self.get_dim('num_side_ss%d'%(ibc+1)) 115 elem = self.get_array('elem_ss%d'%(ibc+1), (nface,), 'int32') 116 side = self.get_array('side_ss%d'%(ibc+1), (nface,), 'int32') 117 self.bcs[ibc] = (self.bcs[ibc], elem, side) 118 # coordinate. 119 large = c_int() 120 self.nc_get_att_int(self.ncid, self.NC_GLOBAL, 'file_size', 121 byref(large)) 122 if large.value: 123 vnames = ['coord%s' % ('xyz'[idm]) for idm in range(ndim)] 124 crds = [self.get_array(vn, nnode, 'float64') for vn in vnames] 125 self.ndcrd = hstack([crd.reshape((nnode,1)) for crd in crds]) 126 else: 127 self.ndcrd = self.get_array('coord', (ndim, nnode), 128 'float64').T.copy() 129 # mapper. 130 self.emap = self.get_array('elem_map', (self.ncell,), 'int32')
131
132 - def toblock(self, onlybcnames=None, bcname_mapper=None, fpdtype=None, 133 use_incenter=False):
134 """ 135 Convert Cubit/Genesis/ExodusII object to Block object. 136 137 @keyword onlybcnames: positively list wanted names of BCs. 138 @type onlybcnames: list 139 @keyword bcname_mapper: map name to bc type number. 140 @type bcname_mapper: dict 141 @keyword fpdtype: floating-point dtype. 142 @type fpdtype: str 143 @keyword use_incenter: use incenter when creating block. 144 @type use_incenter: bool 145 @return: Block object. 146 @rtype: solvcon.block.Block 147 """ 148 from ..block import Block 149 blk = Block(ndim=self.ndim, nnode=self.nnode, ncell=self.ncell, 150 fpdtype=fpdtype, use_incenter=use_incenter) 151 self._convert_interior_to(blk) 152 blk.build_interior() 153 self._convert_bc_to(blk, 154 onlynames=onlybcnames, name_mapper=bcname_mapper) 155 blk.build_boundary() 156 blk.build_ghost() 157 return blk
158 159 CLTPN_MAP = { 160 'SHELL4': 2, # CUBIT wants 2D quads in this way <shrug>. 161 'TRI3': 3, 162 'HEX8': 4, 163 'TETRA': 5, 164 'WEDGE': 6, 165 'PYRAMID': 7, 166 }
167 - def _convert_interior_to(self, blk):
168 """ 169 Convert interior connectivities to Block object. 170 171 @param blk: to-be-written Block object. 172 @type blk: solvcon.block.Block 173 @return: nothing 174 """ 175 from ..block import elemtype 176 # coordinate. 177 blk.ndcrd[:] = self.ndcrd[:] 178 # node definition. 179 ien = 0 180 for name, tname, clnds in self.blks: 181 ist = ien 182 ien += clnds.shape[0] 183 # type. 184 blk.cltpn[ist:ien] = self.CLTPN_MAP[tname] 185 # nodes. 186 nnd = elemtype[self.CLTPN_MAP[tname],2] 187 sclnds = blk.clnds[ist:ien] 188 sclnds[:,0] = nnd 189 sclnds[:,1:nnd+1] = clnds 190 sclnds[:,1:nnd+1] -= 1 191 if tname == 'PRISM': 192 arr = sclnds[:,2].copy() 193 sclnds[:,2] = sclnds[:,3] 194 sclnds[:,3] = arr 195 arr[:] = sclnds[:,5] 196 sclnds[:,5] = sclnds[:,6] 197 sclnds[:,6] = arr 198 # groups. 199 blk.grpnames = [it[0] for it in self.blks] 200 iblk = 0 201 ien = 0 202 for name, tname, clnds in self.blks: 203 ist = ien 204 ien += clnds.shape[0] 205 blk.clgrp[ist:ien] = iblk 206 iblk += 1
207
208 - def _convert_bc_to(self, blk, onlynames=None, name_mapper=None):
209 """ 210 Convert boundary condition information into Block object. 211 212 @param blk: to-be-written Block object. 213 @type blk: solvcon.block.Block 214 @keyword onlynames: positively list wanted names of BCs. 215 @type onlynames: list 216 @keyword name_mapper: map name to bc type and value dictionary; value 217 of the key can be a 2- or 3-tuple. If it is a 2-tuple (the usual 218 case), the first item is bc type and the second item is value array 219 dict. If it is a 3-tuple, the items are bc type, bc constructing 220 keywords, and value array dict. 221 @type name_mapper: dict 222 @return: nothing. 223 """ 224 # process all neutral bc objects. 225 for ibc in range(len(self.bcs)): 226 # extract boundary faces. 227 bc = self._tobc(ibc, blk) 228 # skip unwanted BCs. 229 if onlynames: 230 if bc.name not in onlynames: 231 continue 232 # recreate BC according to name mapping. 233 if name_mapper is not None: 234 # FIXME: this is a new treatment for bcmap dict. The old 235 # approach is a 2-tuple, but a 3-tuple should be better. The 236 # new approach has not been tested, but after fully tested, it 237 # should be adopted as default. 238 bpars = name_mapper.get(bc.name, (None, None, None)) 239 if len(bpars) == 3: 240 pass 241 elif len(bpars) == 2: 242 bpars = (bpars[0], {}, bpars[1]) 243 else: 244 raise ValueError('BC name must be mapped to a 2-/3-tuple ' 245 'for %s'%str(bc)) 246 bct, bkw, vdict = bpars 247 if bct is not None: 248 bc = bct(bc=bc, **bkw) 249 bc.feedValue(vdict) 250 # save to block object. 251 bc.sern = len(blk.bclist) 252 bc.blk = blk 253 blk.bclist.append(bc)
254 255 # define map for clfcs. 256 CLFCS_MAP = {} 257 CLFCS_MAP[2] = [0,0,1,2,3,4] # quadrilateral. 258 CLFCS_MAP[3] = [1,2,3] # triangle. 259 CLFCS_MAP[4] = [5,2,6,4,1,3] # hexahedron. 260 CLFCS_MAP[5] = [2,4,3,1] # tetrahedron. 261 CLFCS_MAP[6] = [4,5,3,1,2] # prism. 262 CLFCS_MAP[7] = [2,3,4,1,5] # pyramid.
263 - def _tobc(self, ibc, blk):
264 """ 265 Extract boundary condition information from self to become a BC object. 266 Only process element/cell type of boundary information. 267 268 @param ibc: index of the BC to be extracted. 269 @type ibc: int 270 @param blk: Block object for reference, nothing will be altered. 271 @type blk: solvcon.block.Block 272 @return: generic BC object. 273 @rtype: solvcon.boundcond.BC 274 """ 275 from numpy import empty 276 from ..boundcond import BC 277 clfcs_map = self.CLFCS_MAP 278 cltpn = blk.cltpn 279 clfcs = blk.clfcs 280 name, elem, side = self.bcs[ibc] 281 nbnd = elem.shape[0] 282 # extrace boundary face list. 283 facn = empty((nbnd,3), dtype='int32') 284 facn.fill(-1) 285 ibnd = 0 286 while ibnd < nbnd: 287 icl = elem[ibnd] - 1 288 tpn = cltpn[icl] 289 ifl = clfcs_map[tpn][side[ibnd]-1] 290 facn[ibnd,0] = clfcs[icl,ifl] 291 ibnd += 1 292 # craft BC object. 293 bc = BC(fpdtype=blk.fpdtype) 294 bc.name = name 295 slct = facn[:,0].argsort() # sort face list for bc object. 296 bc.facn = facn[slct] 297 # finish. 298 return bc
299
300 -class GenesisIO(FormatIO):
301 """ 302 Proxy to Cubit/Genesis/ExodusII file format. 303 """
304 - def load(self, stream, bcrej=None):
305 """ 306 Load block from stream with BC mapper applied. 307 308 @keyword stream: file name to be read. 309 @type stream: str 310 @keyword bcrej: names of the BC to reject. 311 @type bcrej: list 312 @return: the loaded block. 313 @rtype: solvcon.block.Block 314 """ 315 # load file into memory. 316 assert isinstance(stream, basestring) 317 gn = Genesis(stream) 318 gn.load() 319 gn.close_file() 320 # convert loaded neutral object into block object. 321 if bcrej: 322 onlybcnames = list() 323 for name, elem, side in gn.bcs: 324 if name not in bcrej: 325 onlybcnames.append(name) 326 else: 327 onlybcnames = None 328 blk = gn.toblock(onlybcnames=onlybcnames) 329 return blk
330