Skip to content
Snippets Groups Projects
Commit 7c5b41cb authored by many's avatar many
Browse files

add yplus computation

parent dd15aa70
No related branches found
No related tags found
1 merge request!25add yplus computation
Pipeline #8847 passed
Source diff could not be displayed: it is too large. Options to address this: view the blob.
...@@ -108,6 +108,7 @@ def gridSpacing(mu_0, volmesh): ...@@ -108,6 +108,7 @@ def gridSpacing(mu_0, volmesh):
return [Deltax, Deltay, Deltaz] return [Deltax, Deltay, Deltaz]
def constructWallMesh(surfaces): def constructWallMesh(surfaces):
wall = pv.UnstructuredGrid() wall = pv.UnstructuredGrid()
for surf in surfaces: for surf in surfaces:
...@@ -130,6 +131,25 @@ def compute_scalar_gradient(mesh,arrayname): ...@@ -130,6 +131,25 @@ def compute_scalar_gradient(mesh,arrayname):
return mesh return mesh
def calc_dimensionless_yplus(volmesh, surfaces, use_velfield, use_rhofield, mu_0):
surfaceMesh = constructWallMesh(surfaces)
surfacenormals_surface = surfaceMesh.extract_surface().compute_normals()
volmesh = compute_scalar_gradient(volmesh,use_velfield)
walladjacentids = volmesh.find_containing_cell(surfacenormals_surface.points)
volmesh_walladjacent = volmesh.extract_cells(walladjacentids)
volmesh_walladjacent["cellCenters"] = volmesh_walladjacent.cell_centers().points
volmesh_walladjacent["wallNormal"] = [
surfacenormals_surface.point_data["Normals"][surfacenormals_surface.find_closest_point(i)]
for i in volmesh_walladjacent.points]
distcompute = pv.PolyData(volmesh_walladjacent["cellCenters"]).compute_implicit_distance(
surfaceMesh.extract_surface())
volmesh_walladjacent["cellCentersWallDistance"]=distcompute["implicit_distance"]
volmesh_walladjacent["uTaus"] = getWalluTaus(volmesh_walladjacent, mu_0, use_rhofield, use_velfield)
YPlus = volmesh_walladjacent["cellCentersWallDistance"] * volmesh_walladjacent["uTaus"] / mu_0
volmesh_walladjacent["yPlus"]=YPlus
return volmesh_walladjacent
def calc_dimensionless_gridspacing(volmesh, surfaces, use_velfield, use_rhofield, mu_0): def calc_dimensionless_gridspacing(volmesh, surfaces, use_velfield, use_rhofield, mu_0):
""" """
:param volmesh: pyvista-vtk object :param volmesh: pyvista-vtk object
......
...@@ -72,3 +72,73 @@ def test_calc_dimensionless_gridspacing(): ...@@ -72,3 +72,73 @@ def test_calc_dimensionless_gridspacing():
runtest(height=2, length=1, width=1) runtest(height=2, length=1, width=1)
runtest(height=1, length=2, width=1) runtest(height=1, length=2, width=1)
runtest(height=1, length=1, width=2) runtest(height=1, length=1, width=2)
def test_calc_dimensionless_gridspacing():
"""
this method tests the nondimensional gridspacing postprocessing function calc_dimensionless_gridspacing
a volume mesh will be created. boundary meshes will be extracted of the volume mesh.
then a simple couette-velocity-field is defined
calc_dimensionless_gridspacing needs to compute accurately the Delta+-Values
"""
from ntrfc.meshquality.nondimensionals import calc_dimensionless_yplus
import pyvista as pv
import numpy as np
def runtest(height, length, width):
# 11**3 nodes are enough
nodes = 11
# needs to be something simple
mu_0 = 1 # dynamic viscosity
rho = 1
velocity = 2
gradient = velocity / height
# analytical solution
span_x = width / (nodes - 1)
span_y = height / (nodes - 1)
span_z = length / (nodes - 1)
wallshearstress = mu_0 * gradient
wallshearvelocity = np.sqrt(wallshearstress / rho)
deltayplus = wallshearvelocity * span_y/2 / mu_0
# define the mesh
xrng = np.arange(0, nodes, 1, dtype=np.float32)
yrng = np.arange(0, nodes, 1, dtype=np.float32)
zrng = np.arange(0, nodes, 1, dtype=np.float32)
x, y, z = np.meshgrid(xrng, yrng, zrng)
grid = pv.StructuredGrid(x, y, z)
# scale the mesh
grid.points /= nodes - 1
grid.points *= np.array([width, height, length])
# define velocityfield
bounds = grid.bounds
min_z = bounds[4]
grid["U"] = [gradient * (grid.cell_centers().points[::,1][i]-min_z) * np.array([0,0,1]) for i in range(grid.number_of_cells)]
grid["rho"] = np.ones(grid.number_of_cells)
# extract surface
surface = grid.extract_surface()
facecellids = surface.surface_indices()
upperwallids = []
lowerwallids = []
for faceid in facecellids:
cell = surface.extract_cells(faceid)
if all(cell.points[::, 1] == 0):
lowerwallids.append(faceid)
elif all(cell.points[::, 1] == height):
upperwallids.append(faceid)
lowerwall = surface.extract_cells(lowerwallids)
upperwall = surface.extract_cells(upperwallids)
lowerwall["U"] = [gradient * (lowerwall.cell_centers().points[::,1][i]-min_z) * np.array([0,0,1]) for i in range(lowerwall.number_of_cells)]
upperwall["U"] = [gradient * (upperwall.cell_centers().points[::,1][i]-min_z) * np.array([0,0,1]) for i in range(upperwall.number_of_cells)]
lowerwall["rho"] = np.ones(lowerwall.number_of_cells) * rho
upperwall["rho"] = np.ones(upperwall.number_of_cells) * rho
dimless_gridspacings = calc_dimensionless_yplus(grid, [lowerwall, upperwall], "U", "rho", mu_0)
assert all(np.isclose(deltayplus, dimless_gridspacings[
"yPlus"], rtol=0.05)), "calc_dimensionelss_gridspcing in x direction not accurate"
runtest(height=1, length=1, width=1)
runtest(height=2, length=1, width=1)
runtest(height=1, length=2, width=1)
runtest(height=1, length=1, width=2)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment