diff --git a/ntrfc/fluid/turbulence.py b/ntrfc/fluid/turbulence.py
index a71a73ece071ec17f210a7ba5b538b8d68ec0a40..9716f7809443dbb39d592c5045462844e6d92b9d 100644
--- a/ntrfc/fluid/turbulence.py
+++ b/ntrfc/fluid/turbulence.py
@@ -28,7 +28,7 @@ def calcFluc(velo):
velo = np.array(velo)
- return list(velo - mean_velo)
+ return velo - mean_velo
def calcRey(u, v, w):
@@ -36,12 +36,12 @@ def calcRey(u, v, w):
v_ = calcFluc(v)
w_ = calcFluc(w)
- u_u_ = list(np.multiply(u_, u_))
- v_v_ = list(np.multiply(v_, v_))
- w_w_ = list(np.multiply(w_, w_))
- u_v_ = list(np.multiply(u_, v_))
- u_w_ = list(np.multiply(u_, w_))
- v_w_ = list(np.multiply(v_, w_))
+ u_u_ = np.multiply(u_, u_)
+ v_v_ = np.multiply(v_, v_)
+ w_w_ = np.multiply(w_, w_)
+ u_v_ = np.multiply(u_, v_)
+ u_w_ = np.multiply(u_, w_)
+ v_w_ = np.multiply(v_, w_)
mean_u_u_ = np.mean(u_u_)
mean_v_v_ = np.mean(v_v_)
diff --git a/ntrfc/meshquality/nondimensionals.py b/ntrfc/meshquality/nondimensionals.py
index e1977279ad14d2a6f51162d192ef1c4755f8ffe6..85b2bc4b6229888952c38df75108bcb97e481648 100644
--- a/ntrfc/meshquality/nondimensionals.py
+++ b/ntrfc/meshquality/nondimensionals.py
@@ -48,29 +48,27 @@ def cellSpans(solutionMesh, calcFrom):
return cell_spans
-def getWalluTaus(solutionMesh, mu_0, rhofieldname, velfieldname):
+def get_wall_shear_stress(mesh, dynamic_viscosity, density_fieldname, velocity_fieldname):
"""
- Calculate the wall shear stress at various points on a surface in a CFD simulation.
+ Calculate the wall shear stress velocity at various points on a surface in a CFD simulation.
Parameters:
- solutionMesh (pv.PolyData): A mesh that represents the solution of the CFD simulation.
- surfaceMesh (pv.PolyData): A mesh that represents the surface on which the wall shear stress is to be calculated.
- mu_0 (float): The dynamic viscosity of the fluid.
- rhofieldname (str): The name of the field that contains the density of the fluid.
- velfieldname (str): The name of the field that contains the velocity of the fluid.
+ mesh (pv.PolyData): A mesh that represents the solution of the CFD simulation.
+ dynamic_viscosity (float): The dynamic viscosity of the fluid.
+ density_fieldname (str): The name of the field that contains the density of the fluid.
+ velocity_fieldname (str): The name of the field that contains the velocity of the fluid.
Returns:
- u_taus (np.ndarray): An array containing the wall shear stress at each point on the surface.
+ wall_shear_stress_velocity (np.ndarray): An array containing the velocity at which a fluid layer adjacent to the surface would need to move in order to experience the same shear stress as the actual fluid layer in contact with the surface.
"""
+ grad_velocity = mesh[f"grad_{velocity_fieldname}"].reshape(mesh.number_of_cells, 3, 3)
+ wall_normals = mesh["wallNormal"]
+ velocity_gradient_normal = vecAbs_list([np.dot(grad_velocity[i], wall_normals[i]) for i in range(mesh.number_of_cells)])
+ fluid_density = mesh[density_fieldname]
+ wall_shear_stress = velocity_gradient_normal * dynamic_viscosity
+ wall_shear_stress_velocity = (wall_shear_stress / fluid_density) ** 0.5
- graduwall = solutionMesh[f"grad_{velfieldname}"].reshape(solutionMesh.number_of_cells, 3, 3)
- wallnormals = solutionMesh["wallNormal"]
- dudy = vecAbs_list([np.dot(graduwall[i], wallnormals[i]) for i in range(solutionMesh.number_of_cells)])
- rhoWs = solutionMesh[rhofieldname]
- tauWs = dudy * mu_0
- u_taus = (tauWs / rhoWs) ** 0.5
-
- return u_taus
+ return wall_shear_stress_velocity
def constructWallMesh(surfaces):
@@ -81,8 +79,9 @@ def constructWallMesh(surfaces):
def compute_scalar_gradient(mesh, arrayname):
- mesh = mesh.compute_derivative(scalars=arrayname)
- mesh[f"grad_{arrayname}"] = mesh["gradient"].reshape((mesh.number_of_cells, 3, 3))
+ mesh = mesh.ptc()
+ mesh = mesh.compute_derivative(arrayname)
+ mesh[f"grad_{arrayname}"] = mesh["gradient"]
return mesh
@@ -91,7 +90,7 @@ def calc_dimensionless_yplus(volmesh, surfaces, use_velfield, use_rhofield, mu_0
surfacenormals_surface = surfaceMesh.extract_surface().compute_normals()
volmesh = compute_scalar_gradient(volmesh, use_velfield)
- walladjacentids = volmesh.find_containing_cell(surfacenormals_surface.points)
+ walladjacentids = volmesh.find_containing_cell(surfacenormals_surface.cell_centers().points)
volmesh_walladjacent = volmesh.extract_cells(walladjacentids)
volmesh_walladjacent["cellCenters"] = volmesh_walladjacent.cell_centers().points
volmesh_walladjacent["wallNormal"] = [
@@ -100,8 +99,8 @@ def calc_dimensionless_yplus(volmesh, surfaces, use_velfield, use_rhofield, mu_0
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["uTaus"] = get_wall_shear_stress(volmesh_walladjacent, mu_0, use_rhofield, use_velfield)
+ YPlus = np.abs(volmesh_walladjacent["cellCentersWallDistance"]) * volmesh_walladjacent["uTaus"] / mu_0
volmesh_walladjacent["yPlus"] = YPlus
return volmesh_walladjacent
@@ -123,7 +122,8 @@ def calc_dimensionless_gridspacing(volmesh, surfaces, use_velfield, use_rhofield
print("finding walladjacent cells")
volmesh = compute_scalar_gradient(volmesh, use_velfield)
- walladjacentids = volmesh.find_containing_cell(surfacenormals_surface.points)
+
+ walladjacentids = volmesh.find_containing_cell(surfacenormals_surface.cell_centers().points)
volmesh_walladjacent = volmesh.extract_cells(walladjacentids)
volmesh_walladjacent["cellCenters"] = volmesh_walladjacent.cell_centers().points
@@ -138,7 +138,7 @@ def calc_dimensionless_gridspacing(volmesh, surfaces, use_velfield, use_rhofield
volmesh_walladjacent["zSpan"] = np.array([i[2] for i in spanS]) # calculate cell span in span direction
print("calculating wall-shear and friction-velocity")
- volmesh_walladjacent["uTaus"] = getWalluTaus(volmesh_walladjacent, mu_0, use_rhofield, use_velfield)
+ volmesh_walladjacent["uTaus"] = get_wall_shear_stress(volmesh_walladjacent, mu_0, use_rhofield, use_velfield)
uTaus = volmesh_walladjacent["uTaus"]
volmesh_walladjacent["DeltaXPlus"] = volmesh_walladjacent[
"xSpan"] * uTaus / mu_0 # calculate cell span in flow direction
diff --git a/tests/meshquality/test_ntrfc_nondimensionals.py b/tests/meshquality/test_ntrfc_nondimensionals.py
index f68efc27ff54c1a6bab4992acdb9d0d7e7e9067c..8d2d577a0d54677d6dac8512d311d5fb1a15d6a5 100644
--- a/tests/meshquality/test_ntrfc_nondimensionals.py
+++ b/tests/meshquality/test_ntrfc_nondimensionals.py
@@ -213,7 +213,7 @@ def test_calc_dimensionless_gridspacing():
def test_compute_wallshearstress():
- from ntrfc.meshquality.nondimensionals import getWalluTaus
+ from ntrfc.meshquality.nondimensionals import get_wall_shear_stress
# 11**3 nodes are enough
nodes = 11
# needs to be something simple
@@ -271,7 +271,7 @@ def test_compute_wallshearstress():
surfacenormals_surface.point_data["Normals"][surfacenormals_surface.find_closest_point(i)]
for i in volmesh_walladjacent.points]
- utaus = getWalluTaus(volmesh_walladjacent, mu_0, "rho", "U")
+ utaus = get_wall_shear_stress(volmesh_walladjacent, mu_0, "rho", "U")
assert np.all(np.isclose(utaus, wallshearvelocity))