diff --git a/Cassiopee/Connector/Connector/IBM.py b/Cassiopee/Connector/Connector/IBM.py index c4ffe0a13..42c82de1b 100644 --- a/Cassiopee/Connector/Connector/IBM.py +++ b/Cassiopee/Connector/Connector/IBM.py @@ -192,7 +192,7 @@ def prepareIBMData(t_case, t_out, tc_out, t_in=None, to=None, tbox=None, tinit=N if isinstance(t_case, str): tb = C.convertFile2PyTree(t_case) else: tb = Internal.copyTree(t_case) - + ## Note: cartesian = True is left as an input argument to avoid regressing during the non-regression test. ## In the near future the ref. values for the non-regression tests will be updated with cartesian=True. ## At this point, cartesian=True input argument can be deleted. @@ -212,7 +212,7 @@ def prepareIBMData(t_case, t_out, tc_out, t_in=None, to=None, tbox=None, tinit=N print("===========================================") print("Note: t_in is " + RED + "NOT" + END + " a " + RED + "CARTESIAN " + END + "grid") print("===========================================") - + ##[AJ] keep for now...will delete in the near future ##if isinstance(tc_out, str): ## if '/' in tc_out: fileoutpre = tc_out.split('/') @@ -303,7 +303,7 @@ def prepareIBMData(t_case, t_out, tc_out, t_in=None, to=None, tbox=None, tinit=N Cmpi.barrier() _redispatch__(t=t) if verbose: printTimeAndMemory__('blank by IBC bodies', time=python_time.time()-pt0) - + #=================== # STEP 4 : INTERP DATA CHIM #=================== @@ -368,214 +368,186 @@ def prepareIBMData(t_case, t_out, tc_out, t_in=None, to=None, tbox=None, tinit=N if tc2 is not None: return t, tc, tc2 else: return t, tc -#def prepareIBMDataExtrude(t_case, t_out, tc_out, t, to=None, -# depth=2, frontType=1, octreeMode=0, IBCType=1, -# verbose=True, check=False, balancing=False, distribute=False, twoFronts=False, cartesian=True, -# yplus=100., Lref=1., correctionMultiCorpsF42=False, blankingF42=False, wallAdaptF42=None, heightMaxF42=-1., -# tbox=None, extrusion='cart'): -# STILL IN DEV - DO NOT USE -# import Generator.IBM as G_IBM -# import time as python_time -# -# if isinstance(t_case, str): tb = C.convertFile2PyTree(t_case) -# else: tb = Internal.copyTree(t_case) -# -# tbF1 = None -# if tbox: -# tbF1 = Internal.getNodesFromNameAndType(tbox, '*KeepF1*', 'CGNSBase_t') -# tbox = None -# -# ##if isinstance(tc_out, str): fileoutpre = tc_out.split('/') -# ##else: fileoutpre = ['.','template.cgns'] -# -# refstate = Internal.getNodeFromName(tb, 'ReferenceState') -# flowEqn = Internal.getNodeFromName(tb, 'FlowEquationSet') -# -# Reynolds = Internal.getNodeFromName(tb, 'Reynolds') -# if Reynolds is not None: -# Reynolds = Internal.getValue(Reynolds) -# if Reynolds < 1.e5: frontType = 1 -# else: Reynolds = 1.e6 -# -# expand = 3 if frontType != 42 else 4 -# -# dimPb = Internal.getNodeFromName(tb, 'EquationDimension') -# if dimPb is None: raise ValueError('prepareIBMDataPara: EquationDimension is missing in input geometry tree.') -# dimPb = Internal.getValue(dimPb) -# if dimPb == 2: C._initVars(tb, 'CoordinateZ', 0.) -# -# model = Internal.getNodeFromName(tb, 'GoverningEquations') -# if model is None: raise ValueError('prepareIBMDataPara: GoverningEquations is missing in input geometry tree.') -# model = Internal.getValue(model) -# -# ibctypes = Internal.getNodesFromName(tb, 'ibctype') -# if ibctypes is None: raise ValueError('prepareIBMDataPara: ibc type is missing in input geometry tree.') -# ibctypes = list(set(Internal.getValue(ibc) for ibc in ibctypes)) -# -# if model == 'Euler': -# if any(ibc in ['Musker', 'MuskerMob', 'Mafzal', 'Log', 'TBLE', 'TBLE_FULL'] for ibc in ibctypes): -# raise ValueError("prepareIBMDataPara: governing equations (Euler) not consistent with ibc types %s"%(ibctypes)) -# -# #=================== -# # STEP 0 : GET FILAMENT BODIES Modified -# #=================== -# [filamentBases, isFilamentOnly, isOrthoProjectFirst, tb, tbFilament]=D_IBM.determineClosedSolidFilament__(tb) -# isWireModel=False -# for z in Internal.getZones(t_case): -# soldef = Internal.getNodeFromName(z,'.Solver#define') -# if soldef is not None: -# ibctype = Internal.getNodeFromName(soldef,'ibctype') -# if ibctype is not None: -# if Internal.getValue(ibctype) == "wiremodel": -# isWireModel=True -# break -# if isFilamentOnly: tbLocal=tbFilament -# else: tbLocal=Internal.merge([tb,tbFilament]) -# -# # Modification needed for Extrude -# # Keep F1 regions - for F1 & F42 synergy -# if tbF1: -# tbbBTmp = G.BB(tbF1) -# interDict_scale = X.getIntersectingDomains(tbbBTmp, t) -# for kk in interDict_scale: -# for kkk in interDict_scale[kk]: -# z=Internal.getNodeFromName(t, kkk) -# Internal._createUniqueChild(z, '.Solver#defineTMP', 'UserDefinedData_t') -# Internal._createUniqueChild(Internal.getNodeFromName1(z, '.Solver#defineTMP'), 'SaveF1', 'DataArray_t', value=1) -# node=Internal.getNodeFromName(t, kkk) -# -# ## New to extrude Prep -# if extrusion=='cyl': cartesian = False -# -# C._initVars(t, 'centers:cellN', 1.) -# X._applyBCOverlaps(t, depth=depth, loc='centers', val=2, cellNName='cellN') -# C._initVars(t,'{centers:cellNChim}={centers:cellN}') -# C._initVars(t, 'centers:cellN', 1.) -# -# #=================== -# # STEP 1 : BLANKING IBM Modified -# #=================== -# if verbose: pt0 = python_time.time(); printTimeAndMemory__('blank by IBC bodies', time=-1, functionName='setInterpDataIBMExtrude') -# _blankingIBM(t, tb, dimPb=dimPb, frontType=frontType, IBCType=IBCType, depth=depth, -# Reynolds=Reynolds, yplus=yplus, Lref=Lref, twoFronts=twoFronts, -# heightMaxF42=heightMaxF42, correctionMultiCorpsF42=correctionMultiCorpsF42, -# wallAdaptF42=wallAdaptF42, blankingF42=blankingF42, -# filamentBases=filamentBases, isFilamentOnly=isFilamentOnly, tbFilament=tbFilament, -# isWireModel=isWireModel) -# -# ##Modify needed for extrude & to replicate previous behavior -# ##Ghost kmin et kmax donneuse potentiel -# listvars_local =['cellNChim','cellNIBC'] -# for z in Internal.getZones(t): -# sol = Internal.getNodeFromName(z,'FlowSolution#Centers') -# for var in listvars_local: -# cellN = Internal.getNodeFromName(sol,var)[1] -# sh = numpy.shape(cellN) -# for k in [0,1, sh[2]-2, sh[2]-1]: -# for j in range(sh[1]): -# for i in range(sh[0]): -# if cellN[i,j,k] != 0: cellN[i,j,k] =1 -# C._initVars(t,'{centers:cellN}=maximum(0.,{centers:cellNChim})') # vaut -3, 0, 1, 2 initialement -# -# Cmpi.barrier() -# _redispatch__(t=t) -# if verbose: printTimeAndMemory__('blank by IBC bodies', time=python_time.time()-pt0, functionName='setInterpDataIBMExtrude') -# -# ##Same as prepareIBMDataPara -# #=================== -# # STEP 2 : INTERP DATA CHIM -# #=================== -# ## cellN mush be correct here -# ## _setInterpData uses cellN -# if verbose: pt0 = python_time.time(); printTimeAndMemory__('compute interpolation data (Abutting & Chimera)', time=-1, functionName='setInterpDataIBMExtrude') -# tc = C.node2Center(t) -# -# if Internal.getNodeFromType(t, "GridConnectivity1to1_t") is not None: -# Xmpi._setInterpData(t, tc, nature=1, loc='centers', storage='inverse', sameName=1, dim=dimPb, itype='abutting', order=2, cartesian=cartesian) -# Xmpi._setInterpData(t, tc, nature=1, loc='centers', storage='inverse', sameName=1, sameBase=1, dim=dimPb, itype='chimera', order=2, cartesian=cartesian) -# if verbose: printTimeAndMemory__('compute interpolation data (Abutting & Chimera)', time=python_time.time()-pt0, functionName='setInterpDataIBMExtrude') -# -# #=================== -# # STEP 3 : BUILD FRONT -# #=================== -# if verbose: pt0 = python_time.time(); printTimeAndMemory__('build IBM front', time=-1, functionName='setInterpDataIBMExtrude') -# t, tc, front, front2, frontWMM = buildFrontIBM(t, tc, dimPb=dimPb, frontType=frontType, -# cartesian=cartesian, twoFronts=twoFronts, check=check, -# isWireModel=isWireModel) -# if verbose: printTimeAndMemory__('build IBM front', time=python_time.time()-pt0, functionName='setInterpDataIBMExtrude') -# -# #=================== -# # STEP 4 : INTERP DATA IBM -# #=================== -# if verbose: pt0 = python_time.time(); printTimeAndMemory__('compute interpolation data (IBM)', time=-1, functionName='setInterpDataIBMExtrude') -# _setInterpDataIBM(t, tc, tb, front, front2=front2, dimPb=dimPb, frontType=frontType, IBCType=IBCType, depth=depth, -# Reynolds=Reynolds, yplus=yplus, Lref=Lref, -# cartesian=cartesian, twoFronts=twoFronts, check=check, -# isWireModel=isWireModel, tbFilament=tbFilament, isOrthoProjectFirst=isOrthoProjectFirst, isFilamentOnly=isFilamentOnly, -# frontWMM=frontWMM) -# if verbose: printTimeAndMemory__('compute interpolation data (IBM)', time=python_time.time()-pt0, functionName='setInterpDataIBMExtrude') -# -# -# #=================== -# # STEP 5 : INIT IBM Modified -# #=================== -# tc2 = Internal.copyTree(tc) if twoFronts or isWireModel else None -# -# if isWireModel: -# Internal._rmNodesByName(tc2, 'IBCD*') -# Internal._rmNodesByName(tc, '2_IBCD*') -# D_IBM._transformTc2(tc2) -# tc2 = Internal.rmNodesByName(tc2, 'ID*') -# NewIBCD= 141 -# _changeNameIBCD__(tc2,NewIBCD) -# tc = Internal.merge([tc,tc2]) -# tc2 = None -# -# ## New to extrude Prep -# if extrusion == 'cyl': -# T._cyl2Cart(t, (0,0,0),(1,0,0)) -# T._cyl2Cart(tc,(0,0,0),(1,0,0)) -# -# # modif info mesh in zonesubregion_t -# for z in Internal.getZones(tc): -# for zsr in Internal.getNodesFromType(z, "ZoneSubRegion_t"): -# zsrname = Internal.getName(zsr) -# zsrname = zsrname.split('_') -# if zsrname[0]=='IBCD': -# for var in ['C','W','I']: -# r = Internal.getNodeFromName(zsr,'CoordinateY_P'+var)[1] -# theta = Internal.getNodeFromName(zsr,'CoordinateZ_P'+var)[1] -# for l in range(numpy.size(r)): -# yy = r[l]*numpy.cos( theta[l] ) -# zz = r[l]*numpy.sin( theta[l] ) -# r[l]= yy; theta[l] = zz -# -# if distribute and Cmpi.size > 1: _redispatch__(t=t, tc=tc, tc2=tc2, twoFronts=twoFronts) -# -# ## New to extrude Prep -# vars = ['centers:TurbulentDistanceAllBC','centers:TurbulentDistanceWallBC', 'centers:cellNIBC_hole'] -# C._rmVars(t, vars) -# -# if isinstance(tc_out, str): -# tcp = Compressor.compressCartesian(tc) -# Cmpi.convertPyTree2File(tcp, tc_out, ignoreProcNodes=True) -# -# if twoFronts: -# tcp2 = Compressor.compressCartesian(tc2) -# tc2_out = tc_out.replace('tc', 'tc2') if 'tc' in tc_out else 'tc2.cgns' -# Cmpi.convertPyTree2File(tcp2, tc2_out, ignoreProcNodes=True) -# -# if isinstance(t_out, str): -# tp = Compressor.compressCartesian(t) -# Cmpi.convertPyTree2File(tp, t_out, ignoreProcNodes=True) -# -# _computeMeshInfo(t) -# -# if Cmpi.size > 1: Cmpi.barrier() -# if verbose: printTimeAndMemory__('initialize and clean', time=python_time.time()-pt0, functionName='setInterpDataIBMExtrude') -# -# if tc2 is not None: return t, tc, tc2 -# else: return t, tc +def prepareIBMDataExtrude(t_case, t_out, tc_out, t, to=None, + depth=2, frontType=1, octreeMode=0, IBCType=1, + verbose=True, check=False, twoFronts=False, cartesian=True, + yplus=100., Lref=1., correctionMultiCorpsF42=False, blankingF42=False, wallAdaptF42=None, heightMaxF42=-1., + tbox=None, extrusion='cart'): + import Generator.IBM as G_IBM + import time as python_time + + if isinstance(t_case, str): tb = C.convertFile2PyTree(t_case) + else: tb = Internal.copyTree(t_case) + + ##Will be changed in the near future - pending commit #165 + #cartesian = True + #if t: + # cartesian = G_IBM.checkCartesian(t, nghost=2) + # if cartesian: + # RED = "\033[1;31;40m" + # END = "\033[0m" + # print("===========================================") + # print("Note: t_in is a " + RED + "CARTESIAN " + END + "grid") + # print("===========================================") + # else: + # RED = "\033[1;31;40m" + # END = "\033[0m" + # print("===========================================") + # print("Note: t_in is " + RED + "NOT" + END + " a " + RED + "CARTESIAN " + END + "grid") + # print("===========================================") + + refstate = Internal.getNodeFromName(tb, 'ReferenceState') + flowEqn = Internal.getNodeFromName(tb, 'FlowEquationSet') + + Reynolds = Internal.getNodeFromName(tb, 'Reynolds') + if Reynolds is not None: + Reynolds = Internal.getValue(Reynolds) + if Reynolds < 1.e5: frontType = 1 + else: Reynolds = 1.e6 + + expand = 3 if frontType != 42 else 4 + + dimPb = Internal.getNodeFromName(tb, 'EquationDimension') + if dimPb is None: raise ValueError('prepareIBMDataPara: EquationDimension is missing in input geometry tree.') + dimPb = Internal.getValue(dimPb) + if dimPb == 2: C._initVars(tb, 'CoordinateZ', 0.) + + model = Internal.getNodeFromName(tb, 'GoverningEquations') + if model is None: raise ValueError('prepareIBMDataPara: GoverningEquations is missing in input geometry tree.') + model = Internal.getValue(model) + + ibctypes = Internal.getNodesFromName(tb, 'ibctype') + if ibctypes is None: raise ValueError('prepareIBMDataPara: ibc type is missing in input geometry tree.') + ibctypes = list(set(Internal.getValue(ibc) for ibc in ibctypes)) + + if model == 'Euler': + if any(ibc in ['Musker', 'MuskerMob', 'Mafzal', 'Log', 'TBLE', 'TBLE_FULL'] for ibc in ibctypes): + raise ValueError("prepareIBMDataPara: governing equations (Euler) not consistent with ibc types %s"%(ibctypes)) + + #=================== + # STEP 0 : GET FILAMENT BODIES + #=================== + tb, tbFilament = D_IBM.determineClosedSolidFilament__(tb) + isFilamentOnly, isWireModel = D_IBM.localWMMFlags__(tb, tbFilament) + + + if extrusion=='cyl': cartesian = False #__ + C._initVars(t, 'centers:cellN', 1.) # |modification needed for extrude + C._initVars(t, 'centers:cellNChim', 1.) # | + X._applyBCOverlaps(t, depth=depth, loc='centers', val=2, cellNName='cellNChim') # | + C._initVars(t, 'centers:cellN', 1.) #__ + + #=================== + # STEP 1 : GENERATE MESH + #=================== + ##SKIPPED - mesh is provided as an input + + #=================== + # STEP 2 : DIST2WALL + #=================== + ##SKIPPED - mesh is provided as an input & has TurbulentDistance already + + #=================== + # STEP 3 : BLANKING IBM + #=================== + if verbose: pt0 = python_time.time(); printTimeAndMemory__('blank by IBC bodies', time=-1, functionName='prepareIBMDataExtrude') + _blankingIBM(t, tb, dimPb=dimPb, frontType=frontType, IBCType=IBCType, depth=depth, + Reynolds=Reynolds, yplus=yplus, Lref=Lref, twoFronts=twoFronts, + heightMaxF42=heightMaxF42, correctionMultiCorpsF42=correctionMultiCorpsF42, + wallAdaptF42=wallAdaptF42, blankingF42=blankingF42, + tbFilament=tbFilament) + + ##set the kmin et kmax Ghost cells are potential donors #__ + listvars_local =['cellNChim','cellNIBC'] # | + for z in Internal.getZones(t): # | + sol = Internal.getNodeFromName(z,'FlowSolution#Centers') # | + for var in listvars_local: # | + cellN = Internal.getNodeFromName(sol,var)[1] # | Modification needed for extrude & + sh = numpy.shape(cellN) # | to replicate previous behavior + for k in [0,1, sh[2]-2, sh[2]-1]: # | + for j in range(sh[1]): # | + for i in range(sh[0]): # | + if cellN[i,j,k] >0: cellN[i,j,k] =1 # | + C._initVars(t,'{centers:cellN}=maximum(0.,{centers:cellNChim})')# vaut -3, 0, 1, 2 initialement #__ + Cmpi.barrier() + _redispatch__(t=t) + if verbose: printTimeAndMemory__('blank by IBC bodies', time=python_time.time()-pt0, functionName='prepareIBMDataExtrude') + #=================== + # STEP 4 : INTERP DATA CHIM + #=================== + ## REQUIREMENT:: cellN mush be correct here --> _setInterpData uses cellN + ## if cellN* is correct henceforth --> correct values at the end of prepareIBMDataExtrude + if verbose: pt0 = python_time.time(); printTimeAndMemory__('compute interpolation data (Abutting & Chimera)', time=-1, functionName='prepareIBMDataExtrude') + tc = C.node2Center(t) + if Internal.getNodeFromType(t, "GridConnectivity1to1_t") is not None: + Xmpi._setInterpData(t, tc, nature=1, loc='centers', storage='inverse', sameName=1, dim=dimPb, itype='abutting', order=2, cartesian=cartesian) + Xmpi._setInterpData(t, tc, nature=1, loc='centers', storage='inverse', sameName=1, sameBase=1, dim=dimPb, itype='chimera', order=2, cartesian=cartesian) + if verbose: printTimeAndMemory__('compute interpolation data (Abutting & Chimera)', time=python_time.time()-pt0, functionName='prepareIBMDataExtrude') + #=================== + # STEP 5 : BUILD FRONT + #=================== + if verbose: pt0 = python_time.time(); printTimeAndMemory__('build IBM front', time=-1, functionName='prepareIBMDataExtrude') + t, tc, front, front2, frontWMM = buildFrontIBM(t, tc, tb=tb, dimPb=dimPb, frontType=frontType, + cartesian=cartesian, twoFronts=twoFronts, check=check, + tbFilament=tbFilament) + if verbose: printTimeAndMemory__('build IBM front', time=python_time.time()-pt0, functionName='prepareIBMDataExtrude') + #=================== + # STEP 6 : INTERP DATA IBM + #=================== + if verbose: pt0 = python_time.time(); printTimeAndMemory__('compute interpolation data (IBM)', time=-1, functionName='prepareIBMDataExtrude') + _setInterpDataIBM(t, tc, tb, front, front2=front2, dimPb=dimPb, frontType=frontType, IBCType=IBCType, depth=depth, + Reynolds=Reynolds, yplus=yplus, Lref=Lref, + cartesian=cartesian, twoFronts=twoFronts, check=check, + tbFilament=tbFilament, frontWMM=frontWMM) + if verbose: printTimeAndMemory__('compute interpolation data (IBM)', time=python_time.time()-pt0, functionName='prepareIBMDataExtrude') + #=================== + # STEP 7 : INIT IBM + #=================== + if verbose: pt0 = python_time.time(); printTimeAndMemory__('initialize and clean', time=-1, functionName='prepareIBMDataExtrude') + tsave = Internal.copyTree(t) # Modification needed to by pass the initialization of t in the macro function initializeIBM + t = None # + t, tc, tc2 = initializeIBM(t, tc, tb, dimPb=dimPb, twoFronts=twoFronts, tbFilament=tbFilament) + t = Internal.copyTree(tsave) # Modification needed to by pass the initialization of t in the macro function initializeIBM + _redispatch__(t=t, tc=tc, tc2=tc2) + + if extrusion == 'cyl': #__ + T._cyl2Cart(t, (0,0,0),(1,0,0)) # | + T._cyl2Cart(tc,(0,0,0),(1,0,0)) # | + # modif info mesh in zonesubregion_t # | + for z in Internal.getZones(tc): # | + for zsr in Internal.getNodesFromType(z, "ZoneSubRegion_t"): # | + zsrname = Internal.getName(zsr) # | + zsrname = zsrname.split('_') # |Modification needed for extrude in cylindrical coordinates + if zsrname[0]=='IBCD': # | + for var in ['C','W','I']: # | + r = Internal.getNodeFromName(zsr,'CoordinateY_P'+var)[1] # | + theta = Internal.getNodeFromName(zsr,'CoordinateZ_P'+var)[1] # | + for l in range(numpy.size(r)): # | + yy = r[l]*numpy.cos( theta[l] ) # | + zz = r[l]*numpy.sin( theta[l] ) # | + r[l]= yy; theta[l] = zz # | + #__ + + if isinstance(tc_out, str): + tcp = Compressor.compressCartesian(tc) + Cmpi.convertPyTree2File(tcp, tc_out, ignoreProcNodes=True) + + if twoFronts: + tcp2 = Compressor.compressCartesian(tc2) + tc2_out = tc_out.replace('tc', 'tc2') if 'tc' in tc_out else 'tc2.cgns' + Cmpi.convertPyTree2File(tcp2, tc2_out, ignoreProcNodes=True) + + if isinstance(t_out, str): + tp = Compressor.compressCartesian(t) + Cmpi.convertPyTree2File(tp, t_out, ignoreProcNodes=True) + + _computeMeshInfo(t) + + if Cmpi.size > 1: Cmpi.barrier() + if verbose: printTimeAndMemory__('initialize and clean', time=python_time.time()-pt0, functionName='prepareIBMDataExtrude') + + if tc2 is not None: return t, tc, tc2 + else: return t, tc def _redispatch__(t=None, tc=None, tc2=None): import Distributor2.Mpi as D2mpi @@ -786,6 +758,11 @@ def _blankingIBM__(t, tb, dimPb=3, frontType=1, IBCType=1, depth=2, Reynolds=1.e tbFilament=None): isFilamentOnly, isWireModel = D_IBM.localWMMFlags__(tb, tbFilament) + + isSkipDist = False + minval = C.getMinValue(t, 'centers:TurbulentDistance'); + minval = Cmpi.allreduce(minval, op=Cmpi.MIN) + if minval<0: isSkipDist=True snear_min = 10e10 for z in Internal.getZones(tb): @@ -797,12 +774,10 @@ def _blankingIBM__(t, tb, dimPb=3, frontType=1, IBCType=1, depth=2, Reynolds=1.e if snearl is not None: snear_min = min(snear_min,snearl) snear_min = Cmpi.allreduce(snear_min, op=Cmpi.MIN) - if not isFilamentOnly: _blankByIBCBodies(t, tb, 'centers', dimPb) + if not isFilamentOnly: _blankByIBCBodies(t, tb, 'centers', dimPb) #cellN -> 0 outside domain ; 1 inside domain C._initVars(t, '{centers:cellNIBC}={centers:cellN}') - - if not isFilamentOnly: _signDistance(t) - + if not isFilamentOnly and not isSkipDist: _signDistance(t) C._initVars(t,'{centers:cellN}={centers:cellNIBC}') if tbFilament: @@ -823,7 +798,6 @@ def _blankingIBM__(t, tb, dimPb=3, frontType=1, IBCType=1, depth=2, Reynolds=1.e # _blankClosestTargetCells(t,cellNName='cellN', depth=depthL) else: raise ValueError('prepareIBMData: not valid IBCType. Check model.') - else: # F42: tracking of IB points using distance information # the classical algorithm (front 1) is first used to ensure a minimum of two rows of target points around the geometry @@ -993,8 +967,8 @@ def findIsoFront(cellNFront, Dist_1, Dist_2): C._rmVars(t,['centers:TurbulentDistanceFilament','centers:TurbulentDistanceFilamentWMM']) if tbFilament: C._rmVars(t,['centers:TurbulentDistanceSolid']) if isFilamentOnly: C._initVars(t,'{centers:cellNFilWMM}={centers:cellNFil}') - - if not isFilamentOnly: _removeBlankedGrids(t, loc='centers') + + if not isFilamentOnly: _removeBlankedGrids(t, loc='centers') return None def blankingIBM(t, tb, dimPb=3, frontType=1, IBCType=1, depth=2, Reynolds=1.e6, yplus=100, Lref=1., twoFronts=False, @@ -1034,7 +1008,6 @@ def _blankingIBM(t, tb, dimPb=3, frontType=1, IBCType=1, depth=2, Reynolds=1.e6, tbFilament=tbFilament) C._initVars(t, '{centers:cellNIBC}={centers:cellN}') - if IBCType == -1: C._initVars(t,'{centers:cellNDummy}=({centers:cellNIBC}>0.5)*({centers:cellNIBC}<1.5)') X._setHoleInterpolatedPoints(t,depth=1,dir=1,loc='centers',cellNName='cellNDummy',addGC=False) @@ -1047,11 +1020,10 @@ def _blankingIBM(t, tb, dimPb=3, frontType=1, IBCType=1, depth=2, Reynolds=1.e6, Internal.__FlowSolutionCenters__) else: C._initVars(t,'{centers:cellNFront}=logical_and({centers:cellNIBC}>0.5, {centers:cellNIBC}<1.5)') - if isWireModel: C._initVars(t,'{centers:cellNFrontFilWMM}={centers:cellNFilWMM}*({centers:cellNFilWMM}>0.5)+1*({centers:cellNFilWMM}<0.5)') C._initVars(t,'{centers:cellNFrontFilWMM}=logical_and({centers:cellNFrontFilWMM}>0.5, {centers:cellNFrontFilWMM}<1.5)') - + for z in Internal.getZones(t): if twoFronts: epsilon_dist = abs(C.getValue(z,'CoordinateX',1)-C.getValue(z,'CoordinateX',0)) @@ -1715,7 +1687,7 @@ def initializeIBM(t, tc, tb, tinit=None, tbCurvi=None, dimPb=3, twoFronts=False, tc2 = None _tcInitialize__(tc, tc2=tc2, ibctypes=ibctypes, isWireModel=isWireModel) - _tInitialize__(t, tinit=tinit, model=model, isWireModel=isWireModel) + if t: _tInitialize__(t, tinit=tinit, model=model, isWireModel=isWireModel) return t, tc, tc2 diff --git a/Cassiopee/Connector/test/IBMinterpolateExtrude2Dto3D_t1.py b/Cassiopee/Connector/test/IBMinterpolateExtrude2Dto3D_t1.py new file mode 100644 index 000000000..b61214fdd --- /dev/null +++ b/Cassiopee/Connector/test/IBMinterpolateExtrude2Dto3D_t1.py @@ -0,0 +1,114 @@ +# extrude 2D mesh to 3D with Cartesian approach & prepareIBMDataExtrude +import Apps.Fast.IBM as App +import Converter.PyTree as C +import Transform.PyTree as T +import Converter.Internal as Internal +import Generator.IBM as G_IBM +import Connector.IBM as X_IBM +import Connector.PyTree as X +import KCore.test as test +import sys, os +import numpy + + +LOCAL = test.getLocal() +bodySurfaceFile = '../../Apps/test/naca1DNS.cgns' + +# Prepare +vmin = 42 +dfars = 5 +snears = 1 +t, tc = X_IBM.prepareIBMData(bodySurfaceFile, None , None , + snears=snears , dfars=dfars , vmin=vmin, + check=False , frontType=1 , cartesian=False) + +#C.convertPyTree2File(t,LOCAL+'/t2D_checking.cgns') +#C.convertPyTree2File(tc,LOCAL+'/tc2D_checking.cgns') + +bodySurface = C.convertFile2PyTree(bodySurfaceFile) + +t2 = Internal.copyTree(t) +bodySurface2 = Internal.copyTree(bodySurface) + +## =========================================== +## Current Extrude & Interpolation 3D Approach +## =========================================== + +# Extrusion for 3D Mesh +extrusion = 'cart' +span = 1 +NPas = 4+1 #number of nodes +t3D, tb3D = G_IBM.extrudeCartesianZDir(t, bodySurface, extrusion=extrusion, NPas=NPas, span=span, dz=span/(NPas-1), isAutoPeriodic=True) + +for t in [t3D,tb3D]: + zmax = C.getMaxValue(t, 'CoordinateZ'); + zmin = C.getMinValue(t, 'CoordinateZ'); + zavg = (zmax+zmin)/2 + T._translate(t, (0,0,0-zavg)) + +#C.convertPyTree2File(t3D ,LOCAL+'/t3D_checking.cgns') +#C.convertPyTree2File(tb3D,LOCAL+'/tb3D_checking.cgns') + +# Interpolation 3D +t3D, tc3D = X_IBM.prepareIBMDataExtrude(tb3D, None, None, t3D, extrusion=extrusion) + +Internal._rmNodesByName(t3D, 'FlowEquationSet') +Internal._rmNodesByName(t3D, 'ReferenceState') + +Internal._rmNodesByName(tc3D, 'FlowEquationSet') +Internal._rmNodesByName(tc3D, 'ReferenceState') + +test.testT(t3D ,1) +test.testT(tc3D ,2) + +#C.convertPyTree2File(t3D ,LOCAL+'/t3D_checking.cgns') +#C.convertPyTree2File(tc3D,LOCAL+'/tc3D_checking.cgns') + +####Keep below - however, certain bug fixes in App.prepare1 are required +#### 1) recompute sign distance needs to be corrected (Nk[z[0]] += 2*ific-1 +c # -1) +#### 2) tb extrude needs to be corrected (if not self.isFilamentOnly: X_IBM._signDistance(t)) +## +#### ============================================ +#### Previous Extrude & Interpolation 3D Approach +#### ============================================ +## +##X._applyBCOverlaps(t2, depth=2, loc='centers', val=1, cellNName='cellN') +##C._initVars(t2, '{centers:TurbulentDistanceAllBC}={centers:TurbulentDistance}') +##C._initVars(t2, '{centers:cellNIBC_blank}=0*({centers:cellN}<1)+1*({centers:cellN}>0)') +##C._initVars(t2, '{centers:cellNIBC_hole}={centers:cellN}') +##X._applyBCOverlaps(t2, depth=2, loc='centers', val=2, cellNName='cellN') +## +#### Extrusion for 3D Mesh +##t3Dorig, tb3Dorig = App.extrudeCartesian(t2, bodySurface2, extrusion=extrusion, NPas=NPas, span=span, dz=span/(NPas-1), isAutoPeriodic=True) +##for t in [t3Dorig,tb3Dorig]: +## zmax = C.getMaxValue(t, 'CoordinateZ'); +## zmin = C.getMinValue(t, 'CoordinateZ'); +## zavg = (zmax+zmin)/2 +## T._translate(t, (0,0,0-zavg)) +## +##BCs = Internal.getNodesFromType(t3Dorig, "BC_t") +##for bc in BCs: +## if Internal.getValue(bc)=='BCautoperiod': +## nameSplit=bc[0].split(".") +## lenLocal=len(nameSplit) +## bc[0]='.'.join(nameSplit[0:lenLocal-1]) +## +### Interpolation 3D +##interpDataType = 0 +##order = 2 +##t3Dorig, tc3Dorig = App.prepare1(tb3Dorig, None, None, t_in=t3Dorig, extrusion=extrusion, interpDataType=interpDataType, order=order) +## +##C._rmVars(t3Dorig, ['centers:TurbulentDistanceAllBC','centers:cellNIBC_blank','centers:cellNIBC_hole']) +##Internal._rmNodesByName(t3Dorig, 'ReferenceState') +##Internal._rmNodesByName(t3Dorig, 'FlowEquationSet') +## +##Internal._rmNodesByName(tc3Dorig, 'ReferenceState') +##Internal._rmNodesByName(tc3Dorig, 'FlowEquationSet') +## +##C._rmVars(tc3Dorig, ['TurbulentDistanceAllBC','cellNIBC_blank','cellNIBC_hole']) +## +##test.testT(t3Dorig ,1) +##test.testT(tc3Dorig ,2) +## +###C.convertPyTree2File(t3Dorig , LOCAL+'/t3Dorig_checking.cgns') +###C.convertPyTree2File(tc3Dorig, LOCAL+'/tc3Dorig_checking.cgns') diff --git a/Cassiopee/Connector/test/ToBeCompletedIBMinterpolateExtrude2Dto3D.py b/Cassiopee/Connector/test/ToBeCompletedIBMinterpolateExtrude2Dto3D.py deleted file mode 100644 index 674fb5d11..000000000 --- a/Cassiopee/Connector/test/ToBeCompletedIBMinterpolateExtrude2Dto3D.py +++ /dev/null @@ -1,54 +0,0 @@ - -## DO NOT RUN THIS TEST CASE - IT IS IN PROGRESS AND WILL BE COMPLETED SHORTLY - -# extrude 2D mesh to 3D with Cartesian approach -import Converter.PyTree as C -import Converter.Mpi as Cmpi -import Transform.PyTree as T -import Converter.Internal as Internal -import Generator.IBM as G_IBM -import Geom.IBM as Geom -import Connector.IBM as X_IBM -import KCore.test as test -import sys - -LOCAL = test.getLocal() -bodySurfaceFile = LOCAL + '/naca0012.cgns' - -Lcharac = 0.03362355 -Lz = 0.2*Lcharac -dfar = 10.*Lcharac -vmin = 16; -snears = 1; - -###Generate 2D Mesh -t2D,tc2D=X_IBM.prepareIBMDataPara(bodySurfaceFile , None , None , - snears=snears , dfar=dfar , vmin=vmin ) -test.testT(t2D ,1) -test.testT(tc2D,2) -#C.convertPyTree2File(t2D,LOCAL+'/t2D_checking.cgns') - -####Extrusion for 3D Mesh -bodySurface = C.convertFile2PyTree(bodySurfaceFile) - -extrusion = 'cart' -span = Lz -NPas = 10+1 #number of nodes -t3D, tb3D = G_IBM.extrudeCartesianZDir(t2D, bodySurface, extrusion=extrusion, NPas=NPas, span=span, dz=span/(NPas-1), isAutoPeriodic=True) - -for t in [t3D,tb3D]: - zmax = C.getMaxValue(t, 'CoordinateZ'); - zmin = C.getMinValue(t, 'CoordinateZ'); - zavg = (zmax+zmin)/2 - T._translate(t, (0,0,0-zavg)) -test.testT(t3D ,3) -test.testT(tb3D ,4) - -#C.convertPyTree2File(t3D,LOCAL+'/t3D_checking.cgns') -#C.convertPyTree2File(tb3D,LOCAL+'/tb3D_checking.cgns') - -#####Interpolation 3D -t3D, tc3D = X_IBM.prepareIBMDataExtrude(tb3D, None, None, t3D, extrusion=extrusion) -test.testT(t3D ,5) -test.testT(tc3D ,6) -