Need to extract boundary layer thickness from 2D planar PIV data.

Question

Miguel Olvera 2024-7-5，19:52

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2134871-need-to-extract-boundary-layer-thickness-from-2d-planar-piv-data

回答： Shishir Reddy 2024-7-17，11:12

This is 2D Planar PIV data at a freestream velocity of ~ 13 m/s.

My code currently selects a velocity profile where max(profile_u) = 12.41 m/s. I want to use a velocity profile where max(profile_u) is closer to 13 m/s.

Scatter(u,y) will show all the velocity profiles I am working with.

clc;
close all;
clear all;
%% BASELINE
% Load the data
opts = detectImportOptions('FlatPlate_Baseline_avg0001(in).csv');
opts.VariableNamingRule = 'preserve';
data = readtable('FlatPlate_Baseline_avg0001(in).csv', opts);
% Extract variables
x = data.x / 1000; % Convert to meters
y = data.y / 1000; % Convert to meters
u = data.u;
v = data.v;
%%
% Filter data for x values between 0.05m to 0.08m
filter_indices = (x >= 0.05) & (x <= 0.08);
x_filtered = x(filter_indices);
y_filtered = y(filter_indices);
u_filtered = u(filter_indices);
v_filtered = v(filter_indices);
% Find unique x and y values within the filtered range
unique_x = unique(x_filtered);
unique_y = unique(y_filtered);
nu = 1.5e-5;
kappa = 0.41;
B = 5.0;
profile_spacing = 2 / 1000;
profile_x_positions = unique_x(unique_x >= 0.05);
wall_shear_stress = zeros(size(profile_x_positions));
u_tau_values = zeros(size(profile_x_positions));
delta_99_values = zeros(size(profile_x_positions));
delta_star_values = zeros(size(profile_x_positions));  % Displacement thickness
theta_values = zeros(size(profile_x_positions));      % Momentum thickness
for i = 1:length(profile_x_positions)
    profile_indices = find(abs(x - profile_x_positions(i)) < 1e-6);
    profile_y = y(profile_indices);
    profile_u = u(profile_indices);
    valid_indices = profile_y > 0 & profile_y * max(profile_u) / nu > 1e-6;
    profile_y = profile_y(valid_indices);
    profile_u = profile_u(valid_indices);
    if isempty(profile_y) || isempty(profile_u)
        fprintf('No valid data points at x = %f meters. Skipping...\n', profile_x_positions(i));
        continue;
    end
    % Estimate u_inf (freestream velocity) as the max value of the profile
    u_inf = max(profile_u);
    
    % f(eta) = u / u_inf
    f_eta = profile_u / u_inf;
    [f_eta_unique, unique_indices] = unique(f_eta);
    profile_y_unique = profile_y(unique_indices);
    % Find y_99 where f(eta) = 0.99
    y_99 = interp1(f_eta_unique, profile_y_unique, 0.99, 'linear', 'extrap');
    eta_99 = 0.99;
    % delta_99
    delta_99 = y_99 / eta_99;
    delta_99_values(i) = delta_99;
    % displacement thickness (delta_star)
    delta_star = abs(trapz(profile_y, (1 - f_eta)));
    delta_star_values(i) = delta_star;
    % momentum thickness (theta)
    theta = abs(trapz(profile_y, f_eta .* (1 - f_eta)));
    theta_values(i) = theta;
    Re_theta = (u_inf * theta)/nu;
    
    % Calculate wall shear stress and friction velocity
    fun_u_tau = @(u_tau) mean(profile_u - (u_tau / kappa) * log(profile_y * u_tau / nu + 1) - B * u_tau);
    try
        options = optimset('Display', 'off');
        u_tau = fzero(fun_u_tau, [0.01, 10], options);
        u_tau_values(i) = u_tau; 
    catch err
        fprintf('Cannot converge at x = %f meters: %s\n', profile_x_positions(i), err.message);
        continue;
    end
    rho = 1.225;
    tau_wall = rho * u_tau^2;
    wall_shear_stress(i) = tau_wall;
    % % Plot y/delta vs u/u_inf
    % figure;
    % y_over_delta = profile_y / delta_99;
    % u_over_u_inf = profile_u / u_inf;
    % plot(u_over_u_inf, y_over_delta, '-o');
    % title(sprintf('Velocity Profile at x = %.3f meters', profile_x_positions(i)));
    % xlabel('u / u_{\infty}');
    % ylabel('y / \delta');
    % grid on;
end
% Display the results
results_table = table(profile_x_positions, delta_99_values, delta_star_values, theta_values);
disp(results_table);
disp(['Momentum Thickness Reynolds Number (Re_theta): ', num2str(Re_theta)]);

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Shishir Reddy 2024-7-17，11:12

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2134871-need-to-extract-boundary-layer-thickness-from-2d-planar-piv-data#answer_1486772

在 MATLAB Online 中打开

Hi Miguel,

As per my understanding you would like to use a velocity profile where max(profile_u) is closer to 13 m/s for calculating boundary layer thickness.

Upon investigation, I see that the maximum velocity i.e., max(profile_u) does not match the expected freestream velocity of 13m/s because of noise or outliers in the provided data which could be the affecting factor for the maximum value.

This can be solved by normalizing the velocity values, such that the scaled value of max(profile_u) should be 13 m/s.

You can add the following line of code in your script to achieve scaling of the profile_u data.

profile_u = profile_u*(13/max(profile_u));

I hope this helps.