#include <petscpf.h>
#include <petscdmswarm.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
#include <petsctime.h>
#include <petscsys.h>
#include <petscdmcomposite.h>
#include <petscsystypes.h>
#include "variables.h"
#include "ParticleSwarm.h"
#include "walkingsearch.h"
#include "grid.h"
#include "logging.h"
#include "io.h"
#include "interpolation.h"
#include "ParticleMotion.h"
#include "Boundaries.h"

Include dependency graph for wallfunction.h:

This graph shows which files directly or indirectly include this file:

Functions
void	noslip (UserCtx user, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts velocity_boundary, double normal_x, double normal_y, double normal_z)
	Applies no-slip wall boundary condition with linear interpolation.

void	freeslip (UserCtx user, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts velocity_boundary, double normal_x, double normal_y, double normal_z)
	Applies free-slip wall boundary condition.

double	E_coeff (double friction_velocity, double roughness_height, double kinematic_viscosity)
	Computes roughness-modified log-law coefficient E.

double	u_hydset_roughness (double kinematic_viscosity, double wall_distance, double friction_velocity, double roughness_height)
	Computes velocity from log-law for rough walls.

double	f_hydset (double kinematic_viscosity, double known_velocity, double wall_distance, double friction_velocity_guess, double roughness_height)
	Residual function for friction velocity equation (log-law with roughness)

double	df_hydset (double kinematic_viscosity, double known_velocity, double wall_distance, double friction_velocity_guess, double roughness_height)
	Numerical derivative of residual function.

double	find_utau_hydset (double kinematic_viscosity, double known_velocity, double wall_distance, double initial_guess, double roughness_height)
	Solves for friction velocity using Newton-Raphson iteration.

double	nu_t (double yplus)
	Computes turbulent eddy viscosity ratio (ν_t / ν)

double	integrate_1 (double kinematic_viscosity, double wall_distance, double friction_velocity, int integration_mode)
	Integrates eddy viscosity profile from wall to distance y.

double	taw (double kinematic_viscosity, double friction_velocity, double wall_distance, double velocity, double pressure_gradient_tangent)
	Computes wall shear stress with pressure gradient effects.

double	u_Cabot (double kinematic_viscosity, double wall_distance, double friction_velocity, double pressure_gradient_tangent, double wall_shear_stress)
	Computes velocity using Cabot wall function.

double	f_Cabot (double kinematic_viscosity, double velocity, double wall_distance, double friction_velocity_guess, double pressure_gradient_tangent, double pressure_gradient_normal)
	Residual function for Cabot wall function.

double	df_Cabot (double kinematic_viscosity, double velocity, double wall_distance, double friction_velocity_guess, double pressure_gradient_tangent, double pressure_gradient_normal)
	Numerical derivative for Cabot wall function.

void	find_utau_Cabot (double kinematic_viscosity, double velocity, double wall_distance, double initial_guess, double pressure_gradient_tangent, double pressure_gradient_normal, double friction_velocity, double wall_shear_velocity, double *wall_shear_normal)
	Solves for friction velocity using Cabot wall function.

double	u_Werner (double kinematic_viscosity, double wall_distance, double friction_velocity)
	Computes velocity using Werner-Wengle wall function.

double	f_Werner (double kinematic_viscosity, double velocity, double wall_distance, double friction_velocity)
	Residual function for Werner-Wengle iteration.

double	df_Werner (double kinematic_viscosity, double velocity, double wall_distance, double friction_velocity)
	Numerical derivative for Werner-Wengle iteration.

double	find_utau_Werner (double kinematic_viscosity, double velocity, double wall_distance, double initial_guess)
	Solves for friction velocity using Werner-Wengle wall function.

double	u_loglaw (double wall_distance, double friction_velocity, double roughness_length)
	Computes velocity using simple log-law (smooth wall with roughness offset)

double	find_utau_loglaw (double velocity, double wall_distance, double roughness_length)
	Solves for friction velocity using simple log-law (explicit formula)

void	wall_function (UserCtx user, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts velocity_boundary, PetscReal *friction_velocity, double normal_x, double normal_y, double normal_z)
	Applies standard wall function with Werner-Wengle model.

void	wall_function_loglaw (UserCtx user, double roughness_height, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts velocity_boundary, PetscReal *friction_velocity, double normal_x, double normal_y, double normal_z)
	Applies log-law wall function with roughness correction.

void	wall_function_Cabot (UserCtx user, double roughness_height, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts velocity_boundary, PetscReal *friction_velocity, double normal_x, double normal_y, double normal_z, double pressure_gradient_x, double pressure_gradient_y, double pressure_gradient_z, int iteration_count)
	Applies Cabot non-equilibrium wall function with pressure gradients.

Function Documentation

◆ noslip()

void noslip	(	UserCtx *	user,
		double	distance_reference,
		double	distance_boundary,
		Cmpnts	velocity_wall,
		Cmpnts	velocity_reference,
		Cmpnts *	velocity_boundary,
		double	normal_x,
		double	normal_y,
		double	normal_z
	)

Applies no-slip wall boundary condition with linear interpolation.

This function enforces a no-slip boundary condition (zero velocity at the wall) by linearly interpolating between the wall velocity (typically zero) and the velocity at a reference point in the flow.

MATHEMATICAL FORMULATION: For a point at distance sb from the wall, with a reference velocity Uc at distance sc from the wall: U_boundary = U_wall + (U_reference - U_wall) * (sb / sc)

PHYSICAL INTERPRETATION: This provides a first-order approximation assuming linear velocity variation in the near-wall region, which is valid in the viscous sublayer.

Parameters

[in]	user	Simulation context (unused but required for interface)
[in]	distance_reference	Wall-normal distance to reference point (sc)
[in]	distance_boundary	Wall-normal distance to boundary point (sb)
[in]	velocity_wall	Velocity at the wall (Ua), typically zero
[in]	velocity_reference	Velocity at reference point (Uc)
[out]	velocity_boundary	Computed velocity at boundary point (Ub)
[in]	normal_x	X-component of wall normal vector
[in]	normal_y	Y-component of wall normal vector
[in]	normal_z	Z-component of wall normal vector

Note: For moving walls, velocity_wall would be non-zero; The normal vector should point INTO the fluid domain

Applies no-slip wall boundary condition with linear interpolation.

Local to this translation unit.

Definition at line 68 of file wallfunction.c.

{
    (void)user;
    (void)normal_x;
    (void)normal_y;
    (void)normal_z;
    // Compute velocity difference between reference point and wall
    double delta_u = velocity_reference.x - velocity_wall.x;
    double delta_v = velocity_reference.y - velocity_wall.y;
    double delta_w = velocity_reference.z - velocity_wall.z;
    
    // Linear interpolation factor
    double interpolation_factor = distance_boundary / distance_reference;
    
    // Apply linear interpolation
    (*velocity_boundary).x = interpolation_factor * delta_u;
    (*velocity_boundary).y = interpolation_factor * delta_v;
    (*velocity_boundary).z = interpolation_factor * delta_w;
    
    // Add wall velocity (shift to absolute frame)
    (*velocity_boundary).x += velocity_wall.x;
    (*velocity_boundary).y += velocity_wall.y;
    (*velocity_boundary).z += velocity_wall.z;
}

Here is the caller graph for this function:

◆ freeslip()

void freeslip	(	UserCtx *	user,
		double	distance_reference,
		double	distance_boundary,
		Cmpnts	velocity_wall,
		Cmpnts	velocity_reference,
		Cmpnts *	velocity_boundary,
		double	normal_x,
		double	normal_y,
		double	normal_z
	)

Applies free-slip wall boundary condition.

Free-slip conditions allow tangential flow but enforce zero normal velocity. This is appropriate for inviscid walls or symmetry planes where there is no shear stress but flow cannot penetrate the boundary.

DECOMPOSITION: Velocity is decomposed into normal and tangential components: U = U_n * n + U_t where U_n = U · n (normal component) U_t = U - U_n * n (tangential component)

BOUNDARY CONDITIONS:

Normal component: Interpolated from interior (∂U_n/∂n ≠ 0)
Tangential component: Extrapolated from interior (∂U_t/∂n = 0)

Parameters

[in]	user	Simulation context
[in]	distance_reference	Wall-normal distance to reference point
[in]	distance_boundary	Wall-normal distance to boundary point
[in]	velocity_wall	Velocity at the wall
[in]	velocity_reference	Velocity at reference point
[out]	velocity_boundary	Computed velocity at boundary point
[in]	normal_x	X-component of wall normal vector
[in]	normal_y	Y-component of wall normal vector
[in]	normal_z	Z-component of wall normal vector

Applies free-slip wall boundary condition.

Local to this translation unit.

Definition at line 101 of file wallfunction.c.

{
    (void)user;
    // Extract normal components of velocities
    double wall_normal_velocity = velocity_wall.x * normal_x + 
                                   velocity_wall.y * normal_y + 
                                   velocity_wall.z * normal_z;
    
    double reference_normal_velocity = velocity_reference.x * normal_x + 
                                        velocity_reference.y * normal_y + 
                                        velocity_reference.z * normal_z;
    
    // Interpolate normal velocity component
    double boundary_normal_velocity = wall_normal_velocity + 
        (reference_normal_velocity - wall_normal_velocity) * (distance_boundary / distance_reference);
    
    // Extract tangential velocity components (extrapolated from reference point)
    double tangential_u = velocity_reference.x - reference_normal_velocity * normal_x;
    double tangential_v = velocity_reference.y - reference_normal_velocity * normal_y;
    double tangential_w = velocity_reference.z - reference_normal_velocity * normal_z;
    
    // Reconstruct total velocity: U = U_t + U_n * n
    (*velocity_boundary).x = tangential_u + boundary_normal_velocity * normal_x;
    (*velocity_boundary).y = tangential_v + boundary_normal_velocity * normal_y;
    (*velocity_boundary).z = tangential_w + boundary_normal_velocity * normal_z;
}

Here is the caller graph for this function:

◆ E_coeff()

double E_coeff	(	double	friction_velocity,
		double	roughness_height,
		double	kinematic_viscosity
	)

Computes roughness-modified log-law coefficient E.

The coefficient E accounts for wall roughness effects on the log-law: u+ = (1/κ) ln(E * y+)

ROUGHNESS REGIMES:

Hydraulically smooth (ks+ < 2.25): E = exp(κ*B), no roughness effect
Transitional (2.25 < ks+ < 90): Smooth interpolation between regimes
Fully rough (ks+ > 90): E modified by roughness height where ks+ = u_τ * ks / ν (roughness Reynolds number)

Parameters

[in]	friction_velocity	Friction velocity u_τ
[in]	roughness_height	Equivalent sand grain roughness height ks
[in]	kinematic_viscosity	Kinematic viscosity ν

Returns: Roughness-modified log-law coefficient E

Note: This follows the Nikuradse sand grain roughness correlation

Computes roughness-modified log-law coefficient E.

Local to this translation unit.

Definition at line 140 of file wallfunction.c.

{
    // Compute roughness Reynolds number
    double roughness_reynolds = friction_velocity * roughness_height / kinematic_viscosity;
    
    double roughness_correction;
    
    if (roughness_reynolds <= ROUGHNESS_TRANSITION_YPLUS) {
        // Hydraulically smooth regime: no correction needed
        roughness_correction = 0.0;
    }
    else if (roughness_reynolds < FULLY_ROUGH_YPLUS) {
        // Transitional regime: smooth interpolation using sine function
        // This provides a gradual transition between smooth and rough wall behavior
        double max_correction = LOGLAW_B - 8.5 + (1.0 / KAPPA) * log(roughness_reynolds);
        
        roughness_correction = max_correction * 
            sin(0.4258 * (log(roughness_reynolds) - 0.811));
    }
    else {
        // Fully rough regime: maximum correction
        roughness_correction = LOGLAW_B - 8.5 + (1.0 / KAPPA) * log(roughness_reynolds);
    }
    
    // Return E = exp[κ(B - ΔB)]
    return exp(KAPPA * (LOGLAW_B - roughness_correction));
}

Here is the caller graph for this function:

◆ u_hydset_roughness()

double u_hydset_roughness	(	double	kinematic_viscosity,
		double	wall_distance,
		double	friction_velocity,
		double	roughness_height
	)

Computes velocity from log-law for rough walls.

Calculates the tangential velocity at a given wall distance using the roughness-modified log-law of the wall.

APPLICABLE REGIMES:

Viscous sublayer (y+ < 11.53): u+ = y+
Log-layer (11.53 < y+ < 300): u+ = (1/κ) ln(E * y+)
Outer layer (y+ > 300): Returns -1 (invalid, wake region)

Parameters

[in]	kinematic_viscosity	Kinematic viscosity ν
[in]	wall_distance	Distance from wall y
[in]	friction_velocity	Friction velocity u_τ
[in]	roughness_height	Equivalent sand grain roughness ks

Returns: Tangential velocity u_t, or -1 if outside valid range

Note: The upper limit y+ < 300 ensures we stay within the log-layer; For y+ > 300, wake effects become important and log-law is invalid

Computes velocity from log-law for rough walls.

Local to this translation unit.

Definition at line 172 of file wallfunction.c.

{
    // Normalized wall distance
    double yplus = friction_velocity * wall_distance / kinematic_viscosity;
    
    // Threshold values for regime transitions
    const double viscous_limit = VISCOUS_SUBLAYER_YPLUS;
    const double loglayer_limit = 300.0;
    
    double tangential_velocity;
    
    if (yplus <= viscous_limit) {
        // Viscous sublayer: linear velocity profile
        tangential_velocity = friction_velocity * yplus;
    }
    else if (yplus <= loglayer_limit) {
        // Log-layer: logarithmic velocity profile with roughness correction
        double E = E_coeff(friction_velocity, roughness_height, kinematic_viscosity);
        tangential_velocity = (friction_velocity / KAPPA) * log(E * yplus);
    }
    else {
        // Outside valid range (wake region or beyond)
        tangential_velocity = -1.0;
    }
    
    return tangential_velocity;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ f_hydset()

double f_hydset	(	double	kinematic_viscosity,
		double	known_velocity,
		double	wall_distance,
		double	friction_velocity_guess,
		double	roughness_height
	)

Residual function for friction velocity equation (log-law with roughness)

This function computes the residual for Newton-Raphson iteration: f(u_τ) = u_predicted(u_τ) - u_known where u_predicted comes from the log-law or linear law depending on y+.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	known_velocity	Known velocity at distance y
[in]	wall_distance	Distance from wall
[in]	friction_velocity_guess	Current guess for u_τ
[in]	roughness_height	Wall roughness height

Returns: Residual value (should be zero at solution)

Residual function for friction velocity equation (log-law with roughness)

Local to this translation unit.

Definition at line 209 of file wallfunction.c.

{
    double yplus = friction_velocity_guess * wall_distance / kinematic_viscosity;
    double residual;
    
    if (yplus <= VISCOUS_SUBLAYER_YPLUS) {
        // Viscous sublayer: u = u_τ * y+
        residual = friction_velocity_guess * yplus - known_velocity;
    }
    else {
        // Log-layer: u = (u_τ / κ) * ln(E * y+)
        double E = E_coeff(friction_velocity_guess, roughness_height, kinematic_viscosity);
        residual = friction_velocity_guess * (1.0 / KAPPA * log(E * yplus)) - known_velocity;
    }
    
    return residual;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ df_hydset()

double df_hydset	(	double	kinematic_viscosity,
		double	known_velocity,
		double	wall_distance,
		double	friction_velocity_guess,
		double	roughness_height
	)

Numerical derivative of residual function.

Computes df/du_τ using finite differences for Newton-Raphson iteration.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	known_velocity	Known velocity
[in]	wall_distance	Distance from wall
[in]	friction_velocity_guess	Current guess for u_τ
[in]	roughness_height	Wall roughness

Returns: Derivative df/du_τ

Numerical derivative of residual function.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: df_hydset()

Definition at line 235 of file wallfunction.c.

{
    const double perturbation = 1.e-7;
    
    double f_plus = f_hydset(kinematic_viscosity, known_velocity, wall_distance,
                            friction_velocity_guess + perturbation, roughness_height);
    
    double f_minus = f_hydset(kinematic_viscosity, known_velocity, wall_distance,
                             friction_velocity_guess - perturbation, roughness_height);
    
    return (f_plus - f_minus) / (2.0 * perturbation);
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ find_utau_hydset()

double find_utau_hydset	(	double	kinematic_viscosity,
		double	known_velocity,
		double	wall_distance,
		double	initial_guess,
		double	roughness_height
	)

Solves for friction velocity using Newton-Raphson iteration.

Given a known velocity at a known distance from the wall, this function iteratively solves for the friction velocity u_τ that satisfies the roughness-modified log-law or linear law.

ALGORITHM: Newton-Raphson: u_τ^(n+1) = u_τ^n - f(u_τ^n) / f'(u_τ^n) Convergence criterion: |u_τ^(n+1) - u_τ^n| < 1e-7

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	known_velocity	Velocity at reference point
[in]	wall_distance	Distance from wall to reference point
[in]	initial_guess	Initial guess for u_τ
[in]	roughness_height	Wall roughness height

Returns: Converged friction velocity u_τ

Note: Maximum 30 iterations; warns if convergence not achieved; Typical convergence in 3-5 iterations for good initial guess

Solves for friction velocity using Newton-Raphson iteration.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: find_utau_hydset()

Definition at line 256 of file wallfunction.c.

{
    double friction_velocity = initial_guess;
    double friction_velocity_old = initial_guess;
    
    const int max_iterations = 30;
    const double convergence_tolerance = 1.e-7;
    
    int iteration;
    for (iteration = 0; iteration < max_iterations; iteration++) {
        // Newton-Raphson update
        double residual = f_hydset(kinematic_viscosity, known_velocity, wall_distance,
                                  friction_velocity_old, roughness_height);
        double derivative = df_hydset(kinematic_viscosity, known_velocity, wall_distance,
                                     friction_velocity_old, roughness_height);
        
        friction_velocity = friction_velocity_old - residual / derivative;
        
        // Check convergence
        if (fabs(friction_velocity - friction_velocity_old) < convergence_tolerance) {
            break;
        }
        
        friction_velocity_old = friction_velocity;
    }
    
    // Warn if convergence not achieved
    if (iteration == max_iterations) {
        LOG_ALLOW(LOCAL, LOG_DEBUG, 
                 "WARNING: u_tau iteration did not converge after %d iterations. "
                 "Final difference: %le\n", 
                 iteration, fabs(friction_velocity - friction_velocity_old));
    }
    
    return friction_velocity;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ nu_t()

double nu_t ( double yplus )

Computes turbulent eddy viscosity ratio (ν_t / ν)

Uses the mixing length model with van Driest damping: ν_t / ν = κ * y+ * [1 - exp(-y+ / A+)]² where A+ ≈ 19 is the damping coefficient.

PHYSICAL INTERPRETATION:

Near wall (y+ → 0): ν_t → 0 (viscous effects dominate)
Far from wall (y+ >> A+): ν_t ~ κ * y+ (mixing length theory)
Damping function smoothly transitions between these limits

Parameters

[in] yplus Normalized wall distance y+ = y * u_τ / ν

Returns: Eddy viscosity ratio ν_t / ν

Note: This is the standard mixing length model with van Driest damping; Valid in the inner layer (y+ < 50-100)

Computes turbulent eddy viscosity ratio (ν_t / ν)

Local to this translation unit.

Definition at line 303 of file wallfunction.c.

{
    // van Driest damping function: [1 - exp(-y+ / A+)]²
    double damping_function = 1.0 - exp(-yplus / DAMPING_COEFFICIENT);
    
    // Mixing length model: ν_t / ν = κ * y+ * D²
    return KAPPA * yplus * damping_function * damping_function;
}

Here is the caller graph for this function:

◆ integrate_1()

double integrate_1	(	double	kinematic_viscosity,
		double	wall_distance,
		double	friction_velocity,
		int	integration_mode
	)

Integrates eddy viscosity profile from wall to distance y.

Computes integrals needed for non-equilibrium wall functions: If mode=0: ∫[0 to y] dy / (ν + ν_t) If mode=1: ∫[0 to y] y dy / (ν + ν_t)

These integrals appear in the solution of the momentum equation with pressure gradients in the near-wall region.

NUMERICAL METHOD:

Trapezoidal rule with 30 integration points
Higher-order correction at endpoints for accuracy

Parameters

[in]	kinematic_viscosity	Kinematic viscosity ν
[in]	wall_distance	Distance from wall y
[in]	friction_velocity	Friction velocity u_τ
[in]	integration_mode	0 for F integral, 1 for Fy integral

Returns: Value of the integral

Note: Used in Cabot wall function for pressure gradient effects

Integrates eddy viscosity profile from wall to distance y.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: integrate_1()

Definition at line 318 of file wallfunction.c.

{
    const int num_points = 30;
    double step_size = wall_distance / (num_points - 1);
    double integral_values[num_points];
    
    // Compute integrand at each point
    for (int i = 0; i < num_points; i++) {
        double current_distance = i * step_size;
        double yplus = current_distance * friction_velocity / kinematic_viscosity;
        double eddy_viscosity_ratio = nu_t(yplus);
        
        if (integration_mode == 0) {
            // Mode 0: integrand = 1 / [(1 + ν_t/ν) * ν]
            integral_values[i] = 1.0 / ((1.0 + eddy_viscosity_ratio) * kinematic_viscosity);
        }
        else {
            // Mode 1: integrand = y / [(1 + ν_t/ν) * ν]
            integral_values[i] = current_distance / 
                                ((1.0 + eddy_viscosity_ratio) * kinematic_viscosity);
        }
    }
    
    // Trapezoidal rule integration with corrected endpoints
    double integral_sum = 0.0;
    
    // Interior points: weight = (f[i-1] + 2*f[i] + f[i+1]) * dy/4
    for (int i = 1; i < num_points - 1; i++) {
        integral_sum += (integral_values[i-1] + 2.0 * integral_values[i] + 
                        integral_values[i+1]) * 0.25 * step_size;
    }
    
    // Endpoint correction for higher accuracy
    integral_sum += 0.1667 * step_size * 
                   (2.0 * integral_values[0] + integral_values[1] +
                    2.0 * integral_values[num_points-1] + integral_values[num_points-2]);
    
    return integral_sum;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ taw()

double taw	(	double	kinematic_viscosity,
		double	friction_velocity,
		double	wall_distance,
		double	velocity,
		double	pressure_gradient_tangent
	)

Computes wall shear stress with pressure gradient effects.

Solves the integrated momentum equation in the near-wall region: τ_w = (u - dp/dx * F_y) / F_1 where F_1 and F_y are integrals of the effective viscosity profile.

This accounts for non-equilibrium effects due to streamwise pressure gradients.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	friction_velocity	Friction velocity u_τ
[in]	wall_distance	Distance from wall
[in]	velocity	Velocity at wall_distance
[in]	pressure_gradient_tangent	Tangential pressure gradient dp/ds

Returns: Wall shear stress τ_w

Note: This is more accurate than equilibrium wall functions in separated flows

Computes wall shear stress with pressure gradient effects.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: taw()

Definition at line 369 of file wallfunction.c.

{
    // Compute integrated viscosity profiles
    double F1 = integrate_1(kinematic_viscosity, wall_distance, friction_velocity, 0);
    double Fy = integrate_1(kinematic_viscosity, wall_distance, friction_velocity, 1);
    
    // Apply momentum balance: τ_w * F1 = u - (dp/dx) * Fy
    return (1.0 / F1) * (velocity - pressure_gradient_tangent * Fy);
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ u_Cabot()

double u_Cabot	(	double	kinematic_viscosity,
		double	wall_distance,
		double	friction_velocity,
		double	pressure_gradient_tangent,
		double	wall_shear_stress
	)

Computes velocity using Cabot wall function.

Reconstructs velocity from wall shear stress and pressure gradient: u = τ_w * F1 + (dp/dx) * Fy

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	wall_distance	Distance from wall
[in]	friction_velocity	Friction velocity
[in]	pressure_gradient_tangent	Tangential pressure gradient
[in]	wall_shear_stress	Wall shear stress

Returns: Velocity at wall_distance

Computes velocity using Cabot wall function.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: u_Cabot()

Definition at line 386 of file wallfunction.c.

{
    double F1 = integrate_1(kinematic_viscosity, wall_distance, friction_velocity, 0);
    double Fy = integrate_1(kinematic_viscosity, wall_distance, friction_velocity, 1);
    
    return wall_shear_stress * F1 + pressure_gradient_tangent * Fy;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ f_Cabot()

double f_Cabot	(	double	kinematic_viscosity,
		double	velocity,
		double	wall_distance,
		double	friction_velocity_guess,
		double	pressure_gradient_tangent,
		double	pressure_gradient_normal
	)

Residual function for Cabot wall function.

Computes: f(u_τ) = u_τ - sqrt(|τ_w(u_τ)|) This enforces consistency between u_τ and τ_w.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	velocity	Velocity
[in]	wall_distance	Distance from wall
[in]	friction_velocity_guess	Guess for u_τ
[in]	pressure_gradient_tangent	Tangential pressure gradient
[in]	pressure_gradient_normal	Normal pressure gradient (currently unused)

Returns: Residual value

Residual function for Cabot wall function.

Local to this translation unit.

Definition at line 400 of file wallfunction.c.

{
    (void)pressure_gradient_normal;
    double wall_shear = taw(kinematic_viscosity, friction_velocity_guess,
                           wall_distance, velocity, pressure_gradient_tangent);
    
    return friction_velocity_guess - sqrt(fabs(wall_shear));
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ df_Cabot()

double df_Cabot	(	double	kinematic_viscosity,
		double	velocity,
		double	wall_distance,
		double	friction_velocity_guess,
		double	pressure_gradient_tangent,
		double	pressure_gradient_normal
	)

Numerical derivative for Cabot wall function.

Parameters

kinematic_viscosity	Parameter `kinematic_viscosity` passed to `df_Cabot()`.
velocity	Parameter `velocity` passed to `df_Cabot()`.
wall_distance	Parameter `wall_distance` passed to `df_Cabot()`.
friction_velocity_guess	Parameter `friction_velocity_guess` passed to `df_Cabot()`.
pressure_gradient_tangent	Parameter `pressure_gradient_tangent` passed to `df_Cabot()`.
pressure_gradient_normal	Parameter `pressure_gradient_normal` passed to `df_Cabot()`.

Returns: Scalar value produced by df_Cabot().

Numerical derivative for Cabot wall function.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: df_Cabot()

Definition at line 417 of file wallfunction.c.

{
    const double perturbation = 1.e-7;
    
    double f_plus = f_Cabot(kinematic_viscosity, velocity, wall_distance,
                           friction_velocity_guess + perturbation,
                           pressure_gradient_tangent, pressure_gradient_normal);
    
    double f_minus = f_Cabot(kinematic_viscosity, velocity, wall_distance,
                            friction_velocity_guess - perturbation,
                            pressure_gradient_tangent, pressure_gradient_normal);
    
    return (f_plus - f_minus) / (2.0 * perturbation);
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ find_utau_Cabot()

void find_utau_Cabot	(	double	kinematic_viscosity,
		double	velocity,
		double	wall_distance,
		double	initial_guess,
		double	pressure_gradient_tangent,
		double	pressure_gradient_normal,
		double *	friction_velocity,
		double *	wall_shear_velocity,
		double *	wall_shear_normal
	)

Solves for friction velocity using Cabot wall function.

This non-equilibrium wall function accounts for pressure gradient effects, making it more accurate in separated or strongly accelerating/decelerating flows.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	velocity	Velocity at reference point
[in]	wall_distance	Distance from wall
[in]	initial_guess	Initial guess for u_τ
[in]	pressure_gradient_tangent	Tangential pressure gradient
[in]	pressure_gradient_normal	Normal pressure gradient
[out]	friction_velocity	Converged friction velocity
[out]	wall_shear_velocity	Wall shear stress for velocity
[out]	wall_shear_normal	Wall shear stress for normal pressure gradient

Solves for friction velocity using Cabot wall function.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: find_utau_Cabot()

Definition at line 440 of file wallfunction.c.

{
    double current_guess = initial_guess;
    double new_guess;
    
    const int max_iterations = 30;
    const double convergence_tolerance = 1.e-10;
    
    int iteration;
    for (iteration = 0; iteration < max_iterations; iteration++) {
        double residual = f_Cabot(kinematic_viscosity, velocity, wall_distance,
                                 current_guess, pressure_gradient_tangent,
                                 pressure_gradient_normal);
        
        double derivative = df_Cabot(kinematic_viscosity, velocity, wall_distance,
                                    current_guess, pressure_gradient_tangent,
                                    pressure_gradient_normal);
        
        new_guess = fabs(current_guess - residual / derivative);
        
        if (fabs(current_guess - new_guess) < convergence_tolerance) {
            break;
        }
        
        current_guess = new_guess;
    }
    
    if (fabs(current_guess - new_guess) > 1.e-5 && iteration >= 29) {
        printf("\n!!!!!!!! Cabot Iteration Failed !!!!!!!!\n");
    }
    
    *friction_velocity = new_guess;
    *wall_shear_velocity = taw(kinematic_viscosity, new_guess, wall_distance,
                              velocity, pressure_gradient_tangent);
    *wall_shear_normal = taw(kinematic_viscosity, new_guess, wall_distance,
                            0.0, pressure_gradient_normal);
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ u_Werner()

double u_Werner	(	double	kinematic_viscosity,
		double	wall_distance,
		double	friction_velocity
	)

Computes velocity using Werner-Wengle wall function.

Algebraic wall function that provides explicit relation: u+ = y+ for y+ < 11.81 (viscous sublayer) u+ = A * (y+)^B for y+ ≥ 11.81 (power law) where A = 8.3, B = 1/7 are empirical constants.

ADVANTAGES:

Explicit (no iteration required)
Computationally cheap
Reasonable accuracy for equilibrium boundary layers

LIMITATIONS:

Less accurate than log-law for high Reynolds numbers
No pressure gradient effects

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	wall_distance	Distance from wall
[in]	friction_velocity	Friction velocity

Returns: Tangential velocity

Computes velocity using Werner-Wengle wall function.

Local to this translation unit.

Definition at line 489 of file wallfunction.c.

{
    double yplus = friction_velocity * wall_distance / kinematic_viscosity;
    
    // Werner-Wengle constants
    const double power_law_coefficient = 8.3;
    const double power_law_exponent = 1.0 / 7.0;
    
    double velocity;
    
    if (yplus <= VISCOUS_SUBLAYER_YPLUS) {
        // Viscous sublayer: u+ = y+
        velocity = yplus * friction_velocity;
    }
    else {
        // Power-law region: u+ = A * (y+)^B
        velocity = power_law_coefficient * pow(yplus, power_law_exponent) * friction_velocity;
    }
    
    return velocity;
}

Here is the caller graph for this function:

◆ f_Werner()

double f_Werner	(	double	kinematic_viscosity,
		double	velocity,
		double	wall_distance,
		double	friction_velocity
	)

Residual function for Werner-Wengle iteration.

Computes residual: f(u_τ) = u_τ² - g(u, y, ν) where g is derived from the velocity profile inversion.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	velocity	Known velocity
[in]	wall_distance	Distance from wall
[in]	friction_velocity	Guess for friction velocity

Returns: Residual value

Residual function for Werner-Wengle iteration.

Local to this translation unit.

Definition at line 516 of file wallfunction.c.

{
    const double power_law_coefficient = 8.3;
    const double power_law_exponent = 1.0 / 7.0;
    
    // Transition point between viscous and power-law regions
    double transition_distance = VISCOUS_SUBLAYER_YPLUS * kinematic_viscosity / friction_velocity;
    
    // Transition velocity
    double transition_velocity = kinematic_viscosity / (2.0 * transition_distance) *
        pow(power_law_coefficient, 2.0 / (1.0 - power_law_exponent));
    
    double residual;
    
    if (fabs(velocity) <= transition_velocity) {
        // Viscous sublayer regime
        residual = friction_velocity * friction_velocity - 
                  velocity / wall_distance * kinematic_viscosity;
    }
    else {
        // Power-law regime (more complex inversion formula)
        double term1 = 0.5 * (1.0 - power_law_exponent) *
                      pow(power_law_coefficient, (1.0 + power_law_exponent) / 
                          (1.0 - power_law_exponent)) *
                      pow(kinematic_viscosity / wall_distance, 1.0 + power_law_exponent);
        
        double term2 = (1.0 + power_law_exponent) / power_law_coefficient *
                      pow(kinematic_viscosity / wall_distance, power_law_exponent) *
                      fabs(velocity);
        
        residual = friction_velocity * friction_velocity -
                  velocity / fabs(velocity) * pow(term1 + term2, 2.0 / (1.0 + power_law_exponent));
    }
    
    return residual;
}

Here is the caller graph for this function:

◆ df_Werner()

double df_Werner	(	double	kinematic_viscosity,
		double	velocity,
		double	wall_distance,
		double	friction_velocity
	)

Numerical derivative for Werner-Wengle iteration.

Parameters

kinematic_viscosity	Parameter `kinematic_viscosity` passed to `df_Werner()`.
velocity	Parameter `velocity` passed to `df_Werner()`.
wall_distance	Parameter `wall_distance` passed to `df_Werner()`.
friction_velocity	Parameter `friction_velocity` passed to `df_Werner()`.

Returns: Scalar value produced by df_Werner().

Numerical derivative for Werner-Wengle iteration.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: df_Werner()

Definition at line 560 of file wallfunction.c.

{
    const double perturbation = 1.e-7;
    
    double f_plus = f_Werner(kinematic_viscosity, velocity, wall_distance,
                            friction_velocity + perturbation);
    
    double f_minus = f_Werner(kinematic_viscosity, velocity, wall_distance,
                             friction_velocity - perturbation);
    
    return (f_plus - f_minus) / (2.0 * perturbation);
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ find_utau_Werner()

double find_utau_Werner	(	double	kinematic_viscosity,
		double	velocity,
		double	wall_distance,
		double	initial_guess
	)

Solves for friction velocity using Werner-Wengle wall function.

Parameters

[in]	kinematic_viscosity	Kinematic viscosity
[in]	velocity	Velocity at reference point
[in]	wall_distance	Distance from wall
[in]	initial_guess	Initial guess for u_τ

Returns: Converged friction velocity

Solves for friction velocity using Werner-Wengle wall function.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: find_utau_Werner()

Definition at line 580 of file wallfunction.c.

{
    double current_guess = initial_guess;
    double new_guess;
    
    const int max_iterations = 20;
    const double convergence_tolerance = 1.e-7;
    
    int iteration;
    for (iteration = 0; iteration < max_iterations; iteration++) {
        double residual = f_Werner(kinematic_viscosity, velocity, wall_distance, current_guess);
        double derivative = df_Werner(kinematic_viscosity, velocity, wall_distance, current_guess);
        
        new_guess = current_guess - residual / derivative;
        
        if (fabs(current_guess - new_guess) < convergence_tolerance) {
            break;
        }
        
        current_guess = new_guess;
    }
    
    if (fabs(current_guess - new_guess) > 1.e-5 && iteration >= 19) {
        printf("\n!!!!!!!! Werner-Wengle Iteration Failed !!!!!!!!\n");
    }
    
    return new_guess;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ u_loglaw()

double u_loglaw	(	double	wall_distance,
		double	friction_velocity,
		double	roughness_length
	)

Computes velocity using simple log-law (smooth wall with roughness offset)

Simple logarithmic profile: u = (u_τ / κ) * ln((y + y0) / y0) where y0 is a roughness length scale.

Parameters

[in]	wall_distance	Distance from wall
[in]	friction_velocity	Friction velocity
[in]	roughness_length	Roughness length scale y0

Returns: Velocity

Note: This is a simplified model; use u_hydset_roughness for better accuracy

Computes velocity using simple log-law (smooth wall with roughness offset)

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: u_loglaw()

Definition at line 620 of file wallfunction.c.

{
    return friction_velocity * (1.0 / KAPPA) * log((roughness_length + wall_distance) / roughness_length);
}

Here is the caller graph for this function:

◆ find_utau_loglaw()

double find_utau_loglaw	(	double	velocity,
		double	wall_distance,
		double	roughness_length
	)

Solves for friction velocity using simple log-law (explicit formula)

Explicit inversion: u_τ = κ * u / ln((y + y0) / y0)

Parameters

[in]	velocity	Known velocity
[in]	wall_distance	Distance from wall
[in]	roughness_length	Roughness length scale

Returns: Friction velocity

Solves for friction velocity using simple log-law (explicit formula)

Local to this translation unit.

Definition at line 629 of file wallfunction.c.

{
    return KAPPA * velocity / log((wall_distance + roughness_length) / roughness_length);
}

Here is the caller graph for this function:

◆ wall_function()

void wall_function	(	UserCtx *	user,
		double	distance_reference,
		double	distance_boundary,
		Cmpnts	velocity_wall,
		Cmpnts	velocity_reference,
		Cmpnts *	velocity_boundary,
		PetscReal *	friction_velocity,
		double	normal_x,
		double	normal_y,
		double	normal_z
	)

Applies standard wall function with Werner-Wengle model.

This is a high-level interface that:

Decomposes velocity into normal/tangential components
Applies wall function to tangential velocity
Interpolates normal velocity
Reconstructs total velocity

Parameters

[in]	user	Simulation context
[in]	distance_reference	Distance to reference point
[in]	distance_boundary	Distance to boundary point
[in]	velocity_wall	Wall velocity
[in]	velocity_reference	Reference velocity
[out]	velocity_boundary	Output boundary velocity
[out]	friction_velocity	Output friction velocity
[in]	normal_x	X-component of wall normal
[in]	normal_y	Y-component of wall normal
[in]	normal_z	Z-component of wall normal

Applies standard wall function with Werner-Wengle model.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/wallfunction.h.

See also: wall_function()

Definition at line 648 of file wallfunction.c.

{
    SimCtx *simCtx = user->simCtx;
    double kinematic_viscosity = 1.0 / simCtx->ren;
    
    // Decompose velocity into components
    double delta_u = velocity_reference.x - velocity_wall.x;
    double delta_v = velocity_reference.y - velocity_wall.y;
    double delta_w = velocity_reference.z - velocity_wall.z;
    
    double normal_velocity = delta_u * normal_x + delta_v * normal_y + delta_w * normal_z;
    
    double tangential_u = delta_u - normal_velocity * normal_x;
    double tangential_v = delta_v - normal_velocity * normal_y;
    double tangential_w = delta_w - normal_velocity * normal_z;
    
    double tangential_magnitude = sqrt(tangential_u * tangential_u +
                                      tangential_v * tangential_v +
                                      tangential_w * tangential_w);
    
    // Apply Werner-Wengle wall function
    double utau_guess = 0.05;
    double tangential_modeled = u_Werner(kinematic_viscosity, distance_boundary, utau_guess);
    if (friction_velocity) {
        // Werner-Wengle path currently uses a fixed u_tau estimate.
        *friction_velocity = (PetscReal)utau_guess;
    }
    
    // Scale tangential components
    if (tangential_magnitude > 1.e-10) {
        tangential_u *= tangential_modeled / tangential_magnitude;
        tangential_v *= tangential_modeled / tangential_magnitude;
        tangential_w *= tangential_modeled / tangential_magnitude;
    }
    else {
        tangential_u = tangential_v = tangential_w = 0.0;
    }
    
    // Reconstruct velocity: U = U_t + U_n * n
    (*velocity_boundary).x = tangential_u + (distance_boundary / distance_reference) * normal_velocity * normal_x;
    (*velocity_boundary).y = tangential_v + (distance_boundary / distance_reference) * normal_velocity * normal_y;
    (*velocity_boundary).z = tangential_w + (distance_boundary / distance_reference) * normal_velocity * normal_z;
    
    // Add wall velocity
    (*velocity_boundary).x += velocity_wall.x;
    (*velocity_boundary).y += velocity_wall.y;
    (*velocity_boundary).z += velocity_wall.z;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ wall_function_loglaw()

void wall_function_loglaw	(	UserCtx *	user,
		double	roughness_height,
		double	distance_reference,
		double	distance_boundary,
		Cmpnts	velocity_wall,
		Cmpnts	velocity_reference,
		Cmpnts *	velocity_boundary,
		PetscReal *	friction_velocity,
		double	normal_x,
		double	normal_y,
		double	normal_z
	)

Applies log-law wall function with roughness correction.

This is the recommended wall function interface for most applications. It uses the roughness-corrected log-law which is accurate for both smooth and rough walls.

Parameters

[in]	user	Simulation context
[in]	roughness_height	Wall roughness height ks
[in]	distance_reference	Distance to reference point
[in]	distance_boundary	Distance to boundary point
[in]	velocity_wall	Wall velocity (typically zero)
[in]	velocity_reference	Reference velocity
[out]	velocity_boundary	Output boundary velocity
[out]	friction_velocity	Output friction velocity u_τ
[in]	normal_x	X-component of wall normal
[in]	normal_y	Y-component of wall normal
[in]	normal_z	Z-component of wall normal

Note: Sets velocity to reference value if y+ > 300 (outside log-layer)

Applies log-law wall function with roughness correction.

Local to this translation unit.

Definition at line 704 of file wallfunction.c.

{
    SimCtx *simCtx = user->simCtx;
    double kinematic_viscosity = 1.0 / simCtx->ren;
    
    // Decompose velocity
    double delta_u = velocity_reference.x - velocity_wall.x;
    double delta_v = velocity_reference.y - velocity_wall.y;
    double delta_w = velocity_reference.z - velocity_wall.z;
    
    double normal_velocity = delta_u * normal_x + delta_v * normal_y + delta_w * normal_z;
    
    double tangential_u = delta_u - normal_velocity * normal_x;
    double tangential_v = delta_v - normal_velocity * normal_y;
    double tangential_w = delta_w - normal_velocity * normal_z;
    
    double tangential_magnitude = sqrt(tangential_u * tangential_u +
                                      tangential_v * tangential_v +
                                      tangential_w * tangential_w);
    
    // Solve for friction velocity
    *friction_velocity = find_utau_hydset(kinematic_viscosity, tangential_magnitude,
                                         distance_reference, 0.001, roughness_height);
    
    // Compute wall function velocity
    double tangential_modeled = u_hydset_roughness(kinematic_viscosity, distance_boundary,
                                                   *friction_velocity, roughness_height);
    
    // Check if outside valid range (y+ > 300)
    if (tangential_modeled < 0.0) {
        tangential_modeled = tangential_magnitude;
    }
    
    // Scale tangential components
    if (tangential_magnitude > 1.e-10) {
        tangential_u *= tangential_modeled / tangential_magnitude;
        tangential_v *= tangential_modeled / tangential_magnitude;
        tangential_w *= tangential_modeled / tangential_magnitude;
    }
    else {
        tangential_u = tangential_v = tangential_w = 0.0;
    }
    
    // Reconstruct total velocity
    (*velocity_boundary).x = tangential_u + (distance_boundary / distance_reference) * normal_velocity * normal_x;
    (*velocity_boundary).y = tangential_v + (distance_boundary / distance_reference) * normal_velocity * normal_y;
    (*velocity_boundary).z = tangential_w + (distance_boundary / distance_reference) * normal_velocity * normal_z;
    
    // Add wall velocity
    (*velocity_boundary).x += velocity_wall.x;
    (*velocity_boundary).y += velocity_wall.y;
    (*velocity_boundary).z += velocity_wall.z;
}

Here is the call graph for this function:

Here is the caller graph for this function:

◆ wall_function_Cabot()

void wall_function_Cabot	(	UserCtx *	user,
		double	roughness_height,
		double	distance_reference,
		double	distance_boundary,
		Cmpnts	velocity_wall,
		Cmpnts	velocity_reference,
		Cmpnts *	velocity_boundary,
		PetscReal *	friction_velocity,
		double	normal_x,
		double	normal_y,
		double	normal_z,
		double	pressure_gradient_x,
		double	pressure_gradient_y,
		double	pressure_gradient_z,
		int	iteration_count
	)

Applies Cabot non-equilibrium wall function with pressure gradients.

This wall function accounts for pressure gradient effects and is most accurate in separated or strongly accelerating/decelerating flows.

Parameters

[in]	user	Simulation context
[in]	roughness_height	Wall roughness (currently unused in this version)
[in]	distance_reference	Distance to reference point
[in]	distance_boundary	Distance to boundary point
[in]	velocity_wall	Wall velocity
[in]	velocity_reference	Reference velocity
[out]	velocity_boundary	Output boundary velocity
[out]	friction_velocity	Output friction velocity
[in]	normal_x	X-component of wall normal
[in]	normal_y	Y-component of wall normal
[in]	normal_z	Z-component of wall normal
[in]	pressure_gradient_x	X-component of pressure gradient
[in]	pressure_gradient_y	Y-component of pressure gradient
[in]	pressure_gradient_z	Z-component of pressure gradient
[in]	iteration_count	Current iteration count

Note: Only updates friction velocity every 5 iterations for stability; Normal pressure gradient currently set to zero

Applies Cabot non-equilibrium wall function with pressure gradients.

Local to this translation unit.

Definition at line 766 of file wallfunction.c.

{
    (void)roughness_height;
    SimCtx *simCtx = user->simCtx;
    double kinematic_viscosity = 1.0 / simCtx->ren;
    
    // Decompose velocity
    double delta_u = velocity_reference.x - velocity_wall.x;
    double delta_v = velocity_reference.y - velocity_wall.y;
    double delta_w = velocity_reference.z - velocity_wall.z;
    
    double normal_velocity = delta_u * normal_x + delta_v * normal_y + delta_w * normal_z;
    
    double tangential_u = delta_u - normal_velocity * normal_x;
    double tangential_v = delta_v - normal_velocity * normal_y;
    double tangential_w = delta_w - normal_velocity * normal_z;
    
    double tangential_magnitude = sqrt(tangential_u * tangential_u +
                                      tangential_v * tangential_v +
                                      tangential_w * tangential_w);
    
    // Compute tangential pressure gradient
    double pressure_gradient_tangent = 0.0;
    if (tangential_magnitude > 1.e-10) {
        pressure_gradient_tangent = (pressure_gradient_x * tangential_u +
                                    pressure_gradient_y * tangential_v +
                                    pressure_gradient_z * tangential_w) / tangential_magnitude;
    }
    
    // Normal pressure gradient (currently set to zero)
    double pressure_gradient_normal = 0.0;
    
    // Solve for friction velocity (only on first iteration or periodically)
    if (iteration_count == 0 || iteration_count > 4) {
        double utau, wall_shear1, wall_shear2;
        find_utau_Cabot(kinematic_viscosity, tangential_magnitude, distance_reference,
                       0.01, pressure_gradient_tangent, pressure_gradient_normal,
                       &utau, &wall_shear1, &wall_shear2);
        *friction_velocity = utau;
    }
    
    // Compute wall function velocity with pressure gradient correction
    double tangential_modeled = u_Cabot(kinematic_viscosity, distance_boundary,
                                       *friction_velocity, pressure_gradient_tangent,
                                       taw(kinematic_viscosity, *friction_velocity,
                                           distance_reference, tangential_magnitude,
                                           pressure_gradient_tangent));
    
    // Scale tangential components
    if (tangential_magnitude > 1.e-10) {
        tangential_u *= tangential_modeled / tangential_magnitude;
        tangential_v *= tangential_modeled / tangential_magnitude;
        tangential_w *= tangential_modeled / tangential_magnitude;
    }
    else {
        tangential_u = tangential_v = tangential_w = 0.0;
    }
    
    // Reconstruct total velocity
    (*velocity_boundary).x = tangential_u + (distance_boundary / distance_reference) * normal_velocity * normal_x;
    (*velocity_boundary).y = tangential_v + (distance_boundary / distance_reference) * normal_velocity * normal_y;
    (*velocity_boundary).z = tangential_w + (distance_boundary / distance_reference) * normal_velocity * normal_z;
    
    // Add wall velocity
    (*velocity_boundary).x += velocity_wall.x;
    (*velocity_boundary).y += velocity_wall.y;
    (*velocity_boundary).z += velocity_wall.z;
}

Here is the call graph for this function:

Here is the caller graph for this function:

Functions

Function Documentation

◆ noslip()

◆ freeslip()

◆ E_coeff()

◆ u_hydset_roughness()

◆ f_hydset()

◆ df_hydset()

◆ find_utau_hydset()

◆ nu_t()

◆ integrate_1()

◆ taw()

◆ u_Cabot()

◆ f_Cabot()

◆ df_Cabot()

◆ find_utau_Cabot()

◆ u_Werner()

◆ f_Werner()

◆ df_Werner()

◆ find_utau_Werner()

◆ u_loglaw()

◆ find_utau_loglaw()

◆ wall_function()

◆ wall_function_loglaw()

◆ wall_function_Cabot()