clear all; close all; clc;
U_inf = 100; % free stream velocity in m/s
L = 1; % plate length in m
nu = 1.5e-5; % kinematic viscosity of air in m^2/s
rho = 1.2; % air density in kg/m^3
Re = U_inf * L / nu; % Reynolds number
n = 1/7; % power law exponent for laminar flow
x = linspace(0.1, 1, 100); % distance from leading edge in m
delta = 5.0 * x / sqrt(Re); % boundary layer thickness
u = U_inf * (x ./ delta).^n % velocity profile
dudy = n * U_inf ./ delta .* (x ./ delta).^(n-1); % velocity gradient
tau = rho * u .* dudy % shear stress
R = (1./delta).*trapz(x.', tau .* u - tau.' * u) % Reynolds stress
plot(x, R);
xlabel('Distance from leading edge (m)');
ylabel('Reynolds stress (Pa)');
title('Reynolds stress profile');