Skip to content
Snippets Groups Projects
Verified Commit 1f95c6a3 authored by Paul Dechamps's avatar Paul Dechamps :speech_balloon:
Browse files

(feat) Modify interpolators for new data structures

parent f0efec9d
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
blast/interfaces/interpolators/blMatchingInterpolator.py
+ 6
− 10
View file @ 1f95c6a3
......@@ -61,15 +61,11 @@ class MatchingInterpolator(Interpolator):
"""
if self.ndim == 2:
for ibody in range(len(iDict)):
for ireg in range(len(iDict[ibody])):
vDict[ibody][0][ireg].updateVars(iDict[ibody][ireg].V, iDict[ibody][ireg].M, iDict[ibody][ireg].Rho)
elif self.ndim == 3:
for ibody in range(len(iDict)):
for isec, ysec in enumerate(self._sections):
for ireg in range(2):
print(iDict[ibody][ireg].nodesCoord[iDict[ireg].nodesCoord[:,1] == ysec])
print(iDict[ibody][ireg].V[iDict[ireg].nodesCoord[:,1] == ysec])
vDict[ibody][isec][ireg].updateVars(iDict[ibody][ireg].V[iDict[ibody][ireg].nodesCoord[:,1] == ysec], iDict[ibody][ireg].Rho[iDict[ibody][ireg].nodesCoord[:,1] == ysec])
for i, ireg in enumerate(iDict[ibody]):
vreg = vDict[ibody][0][i]
vreg.updateVariables(ireg.getMach('all'), ireg.getVelocity('all'), ireg.getDensity('all'))
else:
raise RuntimeError('Incorrect number of dimensions', self.ndim)
def viscousToInviscid(self, iDict, vDict):
"""
......@@ -85,6 +81,6 @@ class MatchingInterpolator(Interpolator):
if self.ndim == 2:
for ibody in range(len(iDict)):
for iReg in range(len(iDict[ibody])):
iDict[ibody][iReg].blowingVel = vDict[ibody][0][iReg].blowingVel
iDict[ibody][iReg].setBlowingVelocity(vDict[ibody][0][iReg].getBlowingVelocity())
else:
raise RuntimeError('Incorrect number of dimensions', self.ndim)
\ No newline at end of file
blast/interfaces/interpolators/blRbfInterpolator.py
+ 145
− 130
View file @ 1f95c6a3
......@@ -34,157 +34,172 @@ class RbfInterpolator(Interpolator):
def inviscidToViscous(self, iDict, vDict):
self.tms['tot'].start()
for ibody in range(len(iDict)):
if iDict[ibody][0].sidesIdentified and iDict[ibody][0].type == 'wing':
for ibody, inviscid_regions in enumerate(iDict):
inviscid_airfoil = inviscid_regions[0]
inviscid_wake = inviscid_regions[1]
if inviscid_airfoil.sidesIdentified and inviscid_airfoil.getType() == 'wing':
self.tms['getIdict'].start()
xi_up = iDict[ibody][0].getUpperNodes()
Mi_up = iDict[ibody][0].getUpperMach()
rhoi_up = iDict[ibody][0].getUpperDensity()
vi_up = iDict[ibody][0].getUpperVelocity()
# lower side
xi_lw = iDict[ibody][0].getLowerNodes()
Mi_lw = iDict[ibody][0].getLowerMach()
rhoi_lw = iDict[ibody][0].getLowerDensity()
vi_lw = iDict[ibody][0].getLowerVelocity()
# wake
xi_wk = iDict[ibody][1].nodesCoord[:, :self.ndim]
Mi_wk = iDict[ibody][1].M
rhoi_wk = iDict[ibody][1].Rho
vi_wk = iDict[ibody][1].V[:,:self.ndim]
self.tms['getIdict'].stop()
# Extract inviscid data for interpolation
xi_up, Mi_up, rhoi_up, vi_up = (
inviscid_airfoil.getNodes('upper_fromLE')[:, :self.ndim],
inviscid_airfoil.getMach('upper_fromLE'),
inviscid_airfoil.getDensity('upper_fromLE'),
inviscid_airfoil.getVelocity('upper_fromLE')[:, :self.ndim],
)
xi_lw, Mi_lw, rhoi_lw, vi_lw = (
inviscid_airfoil.getNodes('lower_fromLE')[:, :self.ndim],
inviscid_airfoil.getMach('lower_fromLE'),
inviscid_airfoil.getDensity('lower_fromLE'),
inviscid_airfoil.getVelocity('lower_fromLE')[:, :self.ndim],
)
xi_wk, Mi_wk, rhoi_wk, vi_wk = (
inviscid_wake.getNodes('wake')[:, :self.ndim],
inviscid_wake.getMach('wake'),
inviscid_wake.getDensity('wake'),
inviscid_wake.getVelocity('wake')[:, :self.ndim],
)
self.tms['getIdict'].stop()
self.tms['getVdict'].start()
# Viscous side
xv_up = np.empty((0, self.ndim), dtype=float)
xv_lw = np.empty((0, self.ndim), dtype=float)
xv_wk = np.empty((0, self.ndim), dtype=float)
for iSec in range(len(vDict[ibody])):
xv_up = np.row_stack((xv_up, vDict[ibody][iSec][0].getUpperNodes()))
xv_lw = np.row_stack((xv_lw, vDict[ibody][iSec][0].getLowerNodes()))
xv_wk = np.row_stack((xv_wk, vDict[ibody][iSec][1].nodesCoord[:,:self.ndim]))
# Get viscous node
xv_up = np.vstack([sec[0].getNodes('upper_fromLE')[:, :self.ndim] for sec in vDict[ibody]])
xv_lw = np.vstack([sec[0].getNodes('lower_fromLE')[:, :self.ndim] for sec in vDict[ibody]])
xv_wk = np.vstack([sec[1].getNodes('wake')[:, :self.ndim] for sec in vDict[ibody]])
self.tms['getVdict'].stop()
# Interpolation
self.tms['rbf'].start()
Mv_up = self.__rbfToSection(xi_up, Mi_up, xv_up)
rhov_up = self.__rbfToSection(xi_up, rhoi_up, xv_up)
vv_up = np.zeros((xv_up.shape[0], self.ndim))
for idim in range(self.ndim):
vv_up[:,idim] = self.__rbfToSection(xi_up, vi_up[:,idim], xv_up)
Mv_lw = self.__rbfToSection(xi_lw, Mi_lw, xv_lw)
rhov_lw = self.__rbfToSection(xi_lw, rhoi_lw, xv_lw)
vv_lw = np.zeros((xv_lw.shape[0], self.ndim))
for idim in range(self.ndim):
vv_lw[:,idim] = self.__rbfToSection(xi_lw, vi_lw[:,idim], xv_lw)
Mv_wk = self.__rbfToSection(xi_wk, Mi_wk, xv_wk)
rhov_wk = self.__rbfToSection(xi_wk, rhoi_wk, xv_wk)
vv_wk = np.zeros((xv_wk.shape[0], self.ndim))
for idim in range(self.ndim):
vv_wk[:,idim] = self.__rbfToSection(xi_wk, vi_wk[:,idim], xv_wk)
self.tms['rbf'].stop()
# The viscous structure is sorted in selig format.
# Here, upper and lower side are 'opposite' (upper is from TE to LE
# while lower is from LE to TE). It is important that it always stays
# this way for this concatenation operation.
# Perform vectorized RBF interpolation
Mv_up, rhov_up = self.__rbfToSection(xi_up, Mi_up, xv_up), self.__rbfToSection(xi_up, rhoi_up, xv_up)
Mv_lw, rhov_lw = self.__rbfToSection(xi_lw, Mi_lw, xv_lw), self.__rbfToSection(xi_lw, rhoi_lw, xv_lw)
Mv_wk, rhov_wk = self.__rbfToSection(xi_wk, Mi_wk, xv_wk), self.__rbfToSection(xi_wk, rhoi_wk, xv_wk)
vv_up = np.column_stack([self.__rbfToSection(xi_up, vi_up[:, dim], xv_up) for dim in range(self.ndim)])
vv_lw = np.column_stack([self.__rbfToSection(xi_lw, vi_lw[:, dim], xv_lw) for dim in range(self.ndim)])
vv_wk = np.column_stack([self.__rbfToSection(xi_wk, vi_wk[:, dim], xv_wk) for dim in range(self.ndim)])
self.tms['rbf'].stop()
self.tms['setting'].start()
cpt_up = 0
cpt_lw = 0
cpt_wk = 0
from matplotlib import pyplot as plt
# Assign values to each viscous section
cpt_up, cpt_lw, cpt_wk = 0, 0, 0
for sec in vDict[ibody]:
M_airfoil = np.empty(sec[0].nPoints, dtype=float)
M_wake = np.empty(sec[1].nPoints, dtype=float)
rho_airfoil = np.empty(sec[0].nPoints, dtype=float)
rho_wake = np.empty(sec[1].nPoints, dtype=float)
v_airfoil = np.zeros((sec[0].nPoints, self.ndim), dtype=float)
v_wake = np.zeros((sec[1].nPoints, self.ndim), dtype=float)
# Mach number
M_airfoil = np.concatenate((Mv_up[cpt_up:cpt_up+sec[0].getUpperNodes().shape[0]],
Mv_lw[cpt_lw:cpt_lw+sec[0].getLowerNodes().shape[0]]))
M_wake = Mv_wk[cpt_wk:cpt_wk+sec[1].nPoints]
# Density
rho_airfoil = np.concatenate((rhov_up[cpt_up:cpt_up+sec[0].getUpperNodes().shape[0]],
rhov_lw[cpt_lw:cpt_lw+sec[0].getLowerNodes().shape[0]]))
rho_wake = rhov_wk[cpt_wk:cpt_wk+sec[1].nPoints]
# Velocity
for idim in range(self.ndim):
v_airfoil[:,idim] = np.concatenate((vv_up[cpt_up:cpt_up+sec[0].getUpperNodes().shape[0],idim],
vv_lw[cpt_lw:cpt_lw+sec[0].getLowerNodes().shape[0],idim]))
v_wake[:,idim] = vv_wk[cpt_wk:cpt_wk+sec[1].nPoints,idim]
cpt_up += sec[0].getUpperNodes().shape[0]
cpt_lw += sec[0].getLowerNodes().shape[0]
cpt_wk += sec[1].nPoints
sec[0].updateVars(v_airfoil, M_airfoil, rho_airfoil)
sec[1].updateVars(v_wake, M_wake, rho_wake)
n_up = sec[0].getNodes('upper_fromLE').shape[0]
n_lw = sec[0].getNodes('lower_fromLE').shape[0]
n_wk = sec[1].getnNodes()
# Concatenate Mach number and density
M_airfoil = np.concatenate((Mv_up[cpt_up:cpt_up + n_up], Mv_lw[cpt_lw:cpt_lw + n_lw]))
rho_airfoil = np.concatenate((rhov_up[cpt_up:cpt_up + n_up], rhov_lw[cpt_lw:cpt_lw + n_lw]))
M_wake, rho_wake = Mv_wk[cpt_wk:cpt_wk + n_wk], rhov_wk[cpt_wk:cpt_wk + n_wk]
# Concatenate velocity
v_airfoil = np.column_stack([
np.concatenate((vv_up[cpt_up:cpt_up + n_up, dim], vv_lw[cpt_lw:cpt_lw + n_lw, dim]))
for dim in range(self.ndim)
])
v_wake = vv_wk[cpt_wk:cpt_wk + n_wk]
# Update the variables
sec[0].updateVariables(M_airfoil, v_airfoil, rho_airfoil)
sec[1].updateVariables(M_wake, v_wake, rho_wake)
# Update counters
cpt_up += n_up
cpt_lw += n_lw
cpt_wk += n_wk
self.tms['setting'].stop()
else:
for iSec in range(len(vDict[ibody])):
for iReg, reg in enumerate(vDict[ibody][iSec]):
v = np.zeros((reg.nodesCoord.shape[0], 3))
M = np.zeros(reg.nodesCoord.shape[0])
rho = np.zeros(reg.nodesCoord.shape[0])
for iDim in range(self.ndim):
v[:,iDim] = self.__rbfToSection(iDict[ibody][iReg].nodesCoord[:,:(self.ndim if iDict[ibody][iReg].name == 'iWing' else 1)], iDict[ibody][iReg].V[:,iDim], reg.nodesCoord[:,:(self.ndim if 'vAirfoil' in reg.name else 1)])
M = self.__rbfToSection(iDict[ibody][iReg].nodesCoord[:,:(self.ndim if iDict[ibody][iReg].name == 'iWing' else 1)], iDict[ibody][iReg].M, reg.nodesCoord[:,:(self.ndim if 'vAirfoil' in reg.name else 1)])
rho = self.__rbfToSection(iDict[ibody][iReg].nodesCoord[:,:(self.ndim if iDict[ibody][iReg].name == 'iWing' else 1)], iDict[ibody][iReg].Rho, reg.nodesCoord[:,:(self.ndim if 'vAirfoil' in reg.name else 1)])
vDict[ibody][iSec][iReg].updateVars(v, M, rho)
# General case
for iSec, sections in enumerate(vDict[ibody]):
for iReg, reg in enumerate(sections):
n_nodes = reg.getnNodes()
slice_dim = self.ndim if 'vAirfoil' in reg._name else 1
xi = iDict[ibody][iReg].getNodes('all')[:, :slice_dim]
xv = reg.getNodes('all')[:, :slice_dim]
# Perform RBF interpolation
v = np.column_stack([
self.__rbfToSection(xi, iDict[ibody][iReg].getVelocity(reg='all')[:,dim], xv)
for dim in range(self.ndim)
])
M = self.__rbfToSection(xi, iDict[ibody][iReg].getMach('all'), xv)
rho = self.__rbfToSection(xi, iDict[ibody][iReg].getDensity('all'), xv)
# Update the variables
reg.updateVariables(M, v, rho)
self.tms['tot'].stop()
def viscousToInviscid(self, iDict, vDict):
for ibody in range(len(iDict)):
if iDict[ibody][0].sidesIdentified:
# Viscous side
xv_up = np.zeros((0, self.ndim))
xv_lw = np.zeros((0, self.ndim))
xv_wk = np.zeros((0, self.ndim))
blwv_up = np.zeros(0)
blwv_lw = np.zeros(0)
blwv_wk = np.zeros(0)
for iSec, sec in enumerate(vDict[ibody]):
xv_up = np.row_stack((xv_up, sec[0].getUpperElems()))
xv_lw = np.row_stack((xv_lw, sec[0].getLowerElems()))
xv_wk = np.row_stack((xv_wk, sec[1].elemsCoord[:,:self.ndim]))
blwv_up = np.concatenate((blwv_up, sec[0].getUpperBlowingVel()))
blwv_lw = np.concatenate((blwv_lw, sec[0].getLowerBlowingVel()))
blwv_wk = np.concatenate((blwv_wk, sec[1].blowingVel))
# Inviscid side
xi_up = iDict[ibody][0].getUpperElems()
xi_lw = iDict[ibody][0].getLowerElems()
xi_wk = iDict[ibody][1].elemsCoord[:,:self.ndim]
blwi_up = self.__rbfToSection(xv_up, blwv_up, xi_up[:,:self.ndim])
blwi_lw = self.__rbfToSection(xv_lw, blwv_lw, xi_lw[:,:self.ndim])
blwi_wk = self.__rbfToSection(xv_wk, blwv_wk, xi_wk[:,:self.ndim])
iDict[ibody][0].blowingVel = np.concatenate((blwi_up, blwi_lw))
iDict[ibody][1].blowingVel = blwi_wk
# Extract elements and blowing velocity for the viscous side
xv_up, xv_lw, xv_wk = [], [], []
blwv_up, blwv_lw, blwv_wk = [], [], []
for sec in vDict[ibody]:
xv_up.append(sec[0].getElms('upper_fromLE')[:, :self.ndim])
xv_lw.append(sec[0].getElms('lower_fromLE')[:, :self.ndim])
xv_wk.append(sec[1].getElms('wake')[:, :self.ndim])
blwv_up.append(sec[0].getBlowingVelocity('upper_fromLE'))
blwv_lw.append(sec[0].getBlowingVelocity('lower_fromLE'))
blwv_wk.append(sec[1].getBlowingVelocity('wake'))
xv_up, xv_lw, xv_wk = map(np.vstack, [xv_up, xv_lw, xv_wk])
blwv_up, blwv_lw, blwv_wk = map(np.concatenate, [blwv_up, blwv_lw, blwv_wk])
# Extract elements for the inviscid side
xi_up = iDict[ibody][0].getElms('upper_fromLE')[:, :self.ndim]
xi_lw = iDict[ibody][0].getElms('lower_fromLE')[:, :self.ndim]
xi_wk = iDict[ibody][1].getElms('wake')[:, :self.ndim]
# Perform RBF interpolation
blwi_up = self.__rbfToSection(xv_up, blwv_up, xi_up)
blwi_lw = self.__rbfToSection(xv_lw, blwv_lw, xi_lw)
blwi_wk = self.__rbfToSection(xv_wk, blwv_wk, xi_wk)
# Set interpolated blowing velocity for inviscid side
iDict[ibody][0].setBlowingVelocity(np.concatenate((blwi_up, blwi_lw)))
iDict[ibody][1].setBlowingVelocity(blwi_wk)
else:
if self.ndim == 2:
raise RuntimeError('2D not implemented')
for iReg, reg in enumerate(iDict[ibody]):
iDict[iReg].blowingVel = self.__rbfToSection(vDict[ibody][0][iReg].elemsCoord[:,:(self.ndim-1 if 'vAirfoil' in reg.name else 1)], vDict[0][iReg].blowingVel, reg.elemsCoord[:,:(self.ndim-1 if reg.name == 'iWing' else 1)])
elif self.ndim == 3:
for iReg, reg in enumerate(iDict[ibody]):
viscElemsCoord = np.zeros((0,3))
viscBlowing = np.zeros(0)
for iSec, sec in enumerate(vDict[ibody]):
viscElemsCoord = np.row_stack((viscElemsCoord, sec[iReg].elemsCoord[:,:3]))
viscBlowing = np.concatenate((viscBlowing, sec[iReg].blowingVel))
if len(self._sym) > 0:
# Collect elements and blowing velocity
viscElemsCoord, viscBlowing = [], []
for sec in vDict[ibody]:
viscElemsCoord.append(sec[iReg].getElms('all'))
viscBlowing.append(sec[iReg].getBlowingVelocity('all'))
# Apply symmetry transformations if needed
if self._sym:
for s in self._sym:
dummy = viscElemsCoord.copy()
dummy[:,1] = (dummy[:,1] - s)*(-1) + s
viscElemsCoord = np.row_stack((viscElemsCoord, dummy))
viscBlowing = np.concatenate((viscBlowing, viscBlowing))
reg.blowingVel = self.__rbfToSection(viscElemsCoord, viscBlowing, reg.elemsCoord[:,:3])
mirrored = np.copy(np.vstack(viscElemsCoord))
mirrored[:, 1] = (mirrored[:, 1] - s) * (-1) + s
viscElemsCoord.append(mirrored)
viscBlowing.append(np.concatenate(viscBlowing))
# Convert to numpy arrays efficiently
viscElemsCoord = np.vstack(viscElemsCoord)
viscBlowing = np.concatenate(viscBlowing)
# Apply interpolation and update
blw = self.__rbfToSection(viscElemsCoord, viscBlowing, reg.getElms('all'))
reg.setBlowingVelocity(blw)
else:
raise RuntimeError('Incorrect number of dimensions', self.cfg['nDim'])
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment