interpolation of a sinusoidal tunnel data

hallo eveyone, i have fluent data of a flow instive a wavy tunnel as present here. I export the data as ACII file and i want to crate matrixs so i can work on the data preaty easely. the problem is that all the scatterd interpolation mathods do not work smothly on a wavy wall, wich asagnifetly decline. i need bater mathod to work with that interpolat on the wall in between those data points and not in the normal x y deraction? any segestions ??? may data is simpal scatterd data of x,y,and u

[Matlab Pro]

Hi @ran -

please share that raw data.

+ you have 1000 of Englishe Spelling mistakes, it is hard to understand what exactly you want.

Please try to correct it and to submit a clear question

[Mathieu NOE]

pls share your data and yur code if you have one

[ran]

Sorry, but the main point is that I am attempting to interpolate data from Fluent, which are on a sinusoidal wall (sinusoidal tunnel). The code and data are straightforward and I hope they are helpful. those waves are not origin in the fluent results and artifact from interpolation, you can see some "steps on the wall wich i asumed creat those

% Let the user choose the folder containing the text files
folder_path = uigetdir('Select the folder containing the text files');

% Check if the user canceled folder selection
if folder_path == 0
    disp('Folder selection canceled. Exiting...');
    return;
end

% Get a list of all the text files in the folder
file_list = dir(fullfile(folder_path, '*.txt'));

% Initialize cell arrays to store data for each variable
x_cell = cell(length(file_list), 1);
y_cell = cell(length(file_list), 1);
u_cell = cell(length(file_list), 1);
v_cell = cell(length(file_list), 1);
T_cell = cell(length(file_list), 1);
P_cell = cell(length(file_list), 1);
Tw_cell = cell(length(file_list), 1);
% rho_cell = cell(length(file_list), 1);
% k_cell = cell(length(file_list), 1);
% dissRate_cell = cell(length(file_list), 1);
% Visxosity_cell = cell(length(file_list), 1);
% Cp_cell = cell(length(file_list), 1);
% thermal_cell = cell(length(file_list), 1);

% Loop through each file and read its contents
for i = 1:length(file_list)
    % Get the file name
    file_name = fullfile(folder_path, file_list(i).name);
    
    % Read the data from the file
    data = readtable(file_name);
    
    % Store the data in the cell arrays with the file name as the cell name
    [~, file_name_only, ~] = fileparts(file_name);
    x_cell{i} = struct('file_name', file_name_only, 'data', data.x_coordinate);
    y_cell{i} = struct('file_name', file_name_only, 'data', data.y_coordinate);
    u_cell{i} = struct('file_name', file_name_only, 'data', data.x_velocity);
    v_cell{i} = struct('file_name', file_name_only, 'data', data.y_velocity);
    T_cell{i} = struct('file_name', file_name_only, 'data', data.temperature);
    P_cell{i} = struct('file_name', file_name_only, 'data', data.pressure);
    Tw_cell{i} = struct('file_name', file_name_only, 'data', data.wall_temperature);
    % rho_cell{i} = struct('file_name', file_name_only, 'data', data.density);
    % k_cell{i} = struct('file_name', file_name_only, 'data', data.turb_kinetic_energy);
    % dissRate_cell{i} = struct('file_name', file_name_only, 'data', data.turb_diss_rate);
    % Visxosity_cell{i} = struct('file_name', file_name_only, 'data', data.viscosity_lam);
    % Cp_cell{i} = struct('file_name', file_name_only, 'data', data.specific_heat_cp);
    % thermal_cell{i} = struct('file_name', file_name_only, 'data', data.thermal_conductivity_lam);
end
%%
% Initialize structures to store wall coordinates and number of y coordinates
wall_coords = struct('x', {}, 'y', {});
num_y_coordinates = struct('file_name', {}, 'num_y', {});
y_coordinates_end_tunnel = struct('file_name', {}, 'num_y', {});

for i = 1:length(file_list)
    % Find wall coordinates
    xwall = x_cell{i}.data(Tw_cell{i}.data > 0);
    ywall = y_cell{i}.data(Tw_cell{i}.data > 0);

    % Find the largest wall X-coordinate (end of the tunnel)
    largest_wall_X = max(xwall);

    % Find the index of the largest wall X-coordinate
    largest_wall_X_idx = find(xwall,1, 'last');   
   

        % Determine the distance between the largest wall X-coordinate and the preceding one
        distance_between = max(gradient(sort(xwall)));

        % Define the search range for X-coordinates
        EndOfTunnel = find(largest_wall_X - 0.7 * distance_between < x_cell{i}.data & x_cell{i}.data <= largest_wall_X);

        % Extract the corresponding y-coordinates
        y_coordinates_end_tunnel(i).file_name = x_cell{i}.file_name;
        y_coordinates_end_tunnel(i).num_y = y_cell{i}.data(EndOfTunnel);

        % Store wall coordinates
        wall_coords(i).x = xwall;
        wall_coords(i).y = ywall;

        % Count the number of y-coordinates
        num_y_coordinates(i).file_name = x_cell{i}.file_name;
        num_y_coordinates(i).num_y = length(y_coordinates_end_tunnel(i).num_y);
    
end

%%
% Initialize structures to store interpolated data
interp_u_cell = cell(length(file_list), 1);
interp_v_cell = cell(length(file_list), 1);
interp_T_cell = cell(length(file_list), 1);
interp_P_cell = cell(length(file_list), 1);
interp_Tw_cell = cell(length(file_list), 1);
interp_x_cell = cell(length(file_list), 1);
interp_y_cell = cell(length(file_list), 1);

for i = 1:length(file_list)
    % Find wall coordinates
    xwall = x_cell{i}.data(Tw_cell{i}.data > 0);
    ywall = y_cell{i}.data(Tw_cell{i}.data > 0);
    
    
    % Find the largest wall X-coordinate (end of the tunnel)
    largest_wall_X = max(xwall);

    % Calculate the distance between the maximum and minimum wall X-coordinates
    distance_total = max(xwall) - min(xwall);
    
    % Define the wavelength
    wavelength = 28.8/1000;  % m

    % Calculate the number of sections based on the wavelength
    num_sections = floor(distance_total / wavelength);

    % Initialize structures to store interpolated data for each section
    interp_u = cell(num_sections, 1);
    interp_v = cell(num_sections, 1);
    interp_T = cell(num_sections, 1);
    interp_P = cell(num_sections, 1);
    interp_Tw = cell(num_sections, 1);
    interp_x= cell(num_sections, 1);
    interp_y = cell(num_sections, 1);

    % Initialize starting point for creating sections
    current_x = largest_wall_X;

    % Loop to create sections and interpolate data within each section
    for j = 1:num_sections
        % Find the end point of the current section
        end_x = current_x - wavelength;

        % Find the indices of x coordinates within the current section
        section_indices = find(x_cell{i}.data >= end_x & x_cell{i}.data < current_x);
        wall_length = length(xwall(xwall >= end_x & xwall < current_x));
        
         % Extract the x and y coordinates within the current section
        section_x = x_cell{i}.data(section_indices);
        section_y = y_cell{i}.data(section_indices);
        section_u = u_cell{i}.data(section_indices);
        section_v = v_cell{i}.data(section_indices);
        section_T = T_cell{i}.data(section_indices);
        section_P = P_cell{i}.data(section_indices);
        section_Tw = Tw_cell{i}.data(section_indices);
        
      
        % Create linearly spaced vectors for interpolation
        % interp_X = linspace(min(section_x), max(section_x), wall_length * 5);
        % interp_Y = linspace(min(section_y), max(section_y), num_y_coordinates(i).num_y * 3);
        % [interp_X_grid, interp_Y_grid] = meshgrid(interp_X, interp_Y);

        interp_X = linspace(min(section_x), max(section_x), 2500);
        interp_Y = linspace(min(section_y), max(section_y),2500*0.694444444);

        % Construct meshgrid
        [interp_x{j}, interp_y{j}] = meshgrid(interp_X, interp_Y);
        % Initialize matrix to hold combined upper and lower wall boundaries
    

            % Interpolation using scatteredInterpolant on normalized coordinates
            F_u = scatteredInterpolant(section_x, section_y, section_u, 'natural', 'none');
            F_v = scatteredInterpolant(section_x, section_y, section_v, 'natural', 'none');
            F_T = scatteredInterpolant(section_x, section_y, section_T, 'natural', 'none');
            F_P = scatteredInterpolant(section_x, section_y, section_P, 'natural', 'none');
            
            
            % Interpolate and re-normalize the data
            interp_u{j} = F_u(interp_x{i}, interp_y{i});
            interp_v{j} = F_v(interp_x{i}, interp_y{i});
            interp_T{j} = F_T(interp_x{i}, interp_y{i});
            interp_P{j} = F_P(interp_x{i}, interp_y{i});
            
            % Update current_x for the next section
        current_x = end_x;
    end

 
% Update cell arrays to contain the combined matrices
interp_u_cell{i} = struct('file_name', x_cell{i}.file_name, 'data', interp_u);
interp_v_cell{i} = struct('file_name', x_cell{i}.file_name, 'data', interp_v);
interp_T_cell{i} = struct('file_name', x_cell{i}.file_name, 'data', interp_T);
interp_P_cell{i} = struct('file_name', x_cell{i}.file_name, 'data', interp_P);
interp_x_cell{i} = struct('file_name', x_cell{i}.file_name, 'data', interp_x);
interp_y_cell{i} = struct('file_name', x_cell{i}.file_name, 'data', interp_y);
end

folder_path = uigetdir('Select the folder containing the text files');
% Define the file name without leading space characters
% Define the file name without leading space characters
filename = 'StaedyFlowM028Re8000.mat';

% Concatenate the folder path and file name to create the full path
full_path = fullfile(folder_path, filename);

% Save the structures to the file using MAT-file version 7.3
save(full_path, 'interp_u_cell', 'interp_v_cell', 'interp_T_cell', 'interp_P_cell', 'interp_Tw_cell', 'interp_x_cell', 'interp_y_cell', '-v7.3');

Answer 1

[Mathieu NOE]

hello again

so I looked a bit some of your data , and my first question would be , why not first plot your (scattered) data with quiver or a colored quiver function, before doing some meshgrid & interpolation

here I used this Fex submission : Quiver - magnitude-dependent color in 2D and 3D - File Exchange - MATLAB Central (mathworks.com)

and you can see that there are no artifacts in the plot rendering.

then , if you really want to do meshgrid & interpolation, I noticed than 500 or 2500 points of interpolation, nor the method used had a great effect on reducing (to zero) the little "steps" close to the boundary; the only fix I could suggest here is to add a bit of 2D smoothing

here I used this Fex submission, but there are other alternatives if you have your own prefered function for this task : smooth2a - File Exchange - MATLAB Central (mathworks.com)

Quiver plot :

Filled contour plot - before smoothing :

Filled contour plot - after smoothing :

load('matlab.mat')


%   Name                Size             Bytes  Class     Attributes
% 
%   section_P       37502x1             300016  double              
%   section_T       37502x1             300016  double              
%   section_Tw      37502x1             300016  double              
%   section_u       37502x1             300016  double              
%   section_v       37502x1             300016  double              
%   section_x       37502x1             300016  double              
%   section_y       37502x1             300016  double  

figure(1) % quiver plot (with colors)
quiverC2D(section_x,section_y,section_u,section_v,1); 
% Fex : https://fr.mathworks.com/matlabcentral/fileexchange/58527-quiver-magnitude-dependent-color-in-2d-and-3d?s_tid=ta_fx_results
colorbar('vert');

N = 1000; % 2500
interp_X = linspace(min(section_x), max(section_x), N);
interp_Y = linspace(min(section_y), max(section_y),round(N*0.694444444)); % "round" added
% Construct meshgrid
[interp_x, interp_y] = meshgrid(interp_X, interp_Y);

%  velocity field
velo = sqrt(section_u.^2+section_v.^2).*sign(section_u); 
% NB sign of section_u is used to show forward or backward oriented flow

% Interpolation using scatteredInterpolant on normalized coordinates
method = 'natural';
F_velo = scatteredInterpolant(section_x, section_y, velo, method, 'none');
Z = F_velo(interp_x, interp_y);

figure(2) % contourf plot 
grid on
[hq,hqlbl] = contourf(interp_X,interp_Y,Z);
clabel(hq,hqlbl)
colorbar('vert');

% contourf plot with smoothing
ZS = smooth2a(Z,8,8);

figure(3) % contourf plot 
grid on
[hq,hqlbl] = contourf(interp_X,interp_Y,ZS);
clabel(hq,hqlbl)
colorbar('vert');

[ran]

Thanks, of course I checked beforehand to see if the data was corrupted when I transferred it from Fluent. These results indicate that the interpolation caused the artefact, which puzzled me as I couldn't find a good solution for that problem anywhere online. Thanks for the solution; I was looking for a good smoothing method to resolve the issue.

[ran]

i thx u allot my frined but i found a new problem, as u see in the pic the intarpulation is problematic even more beacuse it intarpulat from the inlet outlet to the wall. need diper sulution

[Mathieu NOE]

how did you generate this image ?

Answer 2

[ran]

using delaunay

[ran]

i solve it using alphaShape