# Simulation of Faraday Waves with Matlab

15 views (last 30 days)
Thierry on 7 Jan 2023
Dear all,
I have a matlab code described here below which is running but which is incomplete. I would like to visualize the surface displacement of the fluid as a function of the x and y coordinates, but I get only a flat surface. Any help is welcomed !
% Define simulation parameters
rho = 1000; % fluid density [kg/m^3]
mu = 1e-3; % fluid viscosity [Pa*s]
g = 9.81; % acceleration due to gravity [m/s^2]
L = 1; % length of container [m]
W = 0.5; % width of container [m]
A = 0.1; % amplitude of oscillation [m]
omega = 2*pi; % frequency of oscillation [rad/s]
% Set up grid
Nx = 100;
Ny = 50;
x = linspace(0, L, Nx);
y = linspace(0, W, Ny);
[X, Y] = meshgrid(x, y);
% Define simulation parameters
rho = 1000; % fluid density [kg/m^3]
mu = 1e-3; % fluid viscosity [Pa*s]
g = 9.81; % acceleration due to gravity [m/s^2]
L = 1; % length of container [m]
W = 0.5; % width of container [m]
A = 0.1; % amplitude of oscillation [m]
omega = 2*pi; % frequency of oscillation [rad/s]
% Set up grid
Nx = 100;
Ny = 50;
x = linspace(0, L, Nx);
y = linspace(0, W, Ny);
[X, Y] = meshgrid(x, y);
% Set up time-stepping
dt = 0.001; % time step [s]
tmax = 10; % maximum time [s]
t = 0:dt:tmax; % time vector
% Set up initial conditions
u = zeros(Ny, Nx); % initial x velocity [m/s]
v = zeros(Ny, Nx); % initial y velocity [m/s]
eta = zeros(Ny, Nx); % initial displacement [m]
% Run simulation
for i = 2:length(t)
% Compute acceleration at current time step
[a, b] = acceleration(u, v, eta, rho, mu, g, A, omega, t(i));
% Update velocity and displacement using Euler method
u = u + a*dt;
v = v + b*dt;
eta = eta + v*dt;
end
% Visualize results
figure;
surf(X, Y, eta);
xlabel('x');
ylabel('y');
zlabel('displacement');
% Define function to compute acceleration
function [a, b] = acceleration(u, v, eta, rho, mu, g, A, omega, t)
% Compute acceleration using Navier-Stokes equations
a = -rho*g*eta - rho*A*omega^2*sin(omega*t);
b = -rho*g*eta - rho*A*omega^2*cos(omega*t);
end

Sulaymon Eshkabilov on 7 Jan 2023
Edited: Sulaymon Eshkabilov on 7 Jan 2023
Here is a corrected code of your exercise:
clearvars; clc
% Define simulation parameters
rho = 1000; % fluid density [kg/m^3]
mu = 1e-3; % fluid viscosity [Pa*s]
g = 9.81; % acceleration due to gravity [m/s^2]
L = 1; % length of container [m]
W = 0.5; % width of container [m]
A = 0.1; % amplitude of oscillation [m]
omega = 2*pi; % frequency of oscillation [rad/s]
% Set up grid
Nx = 100;
Ny = 50;
% Set up time-stepping
dt = 0.005; % time step [s]
tmax = 10; % maximum time [s]
t = 0:dt:(dt*(Nx-1)); % time vector
% Set up initial conditions
u = zeros(Ny, Nx); % initial x velocity [m/s]
v = zeros(Ny, Nx); % initial y velocity [m/s]
eta = zeros(Ny, Nx); % initial displacement [m]
% Run simulation
for ii = 2:Nx
for jj = 2:Ny
% Compute acceleration at current time step
a = -rho*g*eta(ii-1, jj-1) - rho*A*omega^2*sin(omega*t(ii));
b = -rho*g*eta(ii-1, jj-1) - rho*A*omega^2*cos(omega*t(ii));
% Update velocity and displacement using Euler method
u(ii,jj) = u(ii-1, jj-1) + a*dt;
v(ii,jj) = v(ii-1, jj-1) + b*dt;
eta(ii,jj) = eta(ii-1,jj-1) + v(ii, jj)*dt;
end
end
% Visualize results
figure;
surf(u, v, eta);
xlabel('u');
ylabel('v');
zlabel('displacement');

### Categories

Find more on Fluid Dynamics in Help Center and File Exchange

R2022b

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!