Saturday 26 August 2023

Open - Source coded 3D Navier–Stokes equations in C++

     Here are 3D Navier–Stokes equations configured for lid-driven cavity flow. The syntax is C++. The results from this code are shown in Fig. 1. From the pressure-Poisson equation, I removed mixed derivative terms to improve solution stability. 😆

p[i][j][k] = ((p[i+1][j][k] + p[i-1][j][k]) * h * h + (p[i][j+1][k] + p[i][j-1][k]) * h * h + (p[i][j][k+1] + p[i][j][k-1]) * h * h) / (2 * (h * h + h * h + h * h)) - h * h * h * h * h * h / (2 * (h * h + h * h + h * h)) * (rho * (1 / dt * ((u(i+1, j, k) - u(i-1, j, k)) / (2 * h) + (v(i, j+1, k) - v(i, j-1, k)) / (2 * h) + (w(i, j, k+1) - w(i, j, k-1)) / (2 * h))));

p[0][j][k] = p[1][j][k];

p[num_i - 1][j][k] = p[num_i - 2][j][k];

p[i][0][k] = p[i][1][k];

p[i][num_j - 1][k] = p[i][num_j - 2][k];

p[i][j][0] = p[i][j][1];

p[i][j][num_k - 1] = 0.0;

u[i][j][k] = u[i][j][k] - u[i][j][k] * dt / (2 * h) * (u[i+1][j][k] - u[i-1][j][k]) - v[i][j][k] * dt / (2 * h) * (u[i][j+1][k] - u[i][j-1][k]) - w[i][j][k] * dt / (2 * h) * (u[i][j][k+1] - u[i][j][k-1]) - dt / (2 * rho * h) * (p[i+1][j][k] - p[i-1][j][k]) + nu * (dt / (h * h) * (u[i+1][j][k] - 2 * u[i][j][k] + u[i-1][j][k]) + dt / (h * h) * (u[i][j+1] [k] - 2 * u[i][j][k] + u[i][j-1][k]) + dt / (h * h) * (u[i][j][k+1] - 2 * u[i][j][k] + u[i][j][k-1]));

u[0][j][k] = 0.0;

u[num_i - 1][j][k] = 0.0;

u[i][0][k] = 0.0;

u[i][num_j - 1][k] = 0.0;

u[i][j][0] = 0.0;

u[i][j][num_k - 1] = -Uinf;

v[i][j][k] = v[i][j][k] - u[i][j][k] * dt / (2 * h) * (v[i+1][j][k] - v[i-1][j][k]) - v[i][j][k] * dt / (2 * h) * (v[i][j+1][k] - v[i][j-1][k]) - w[i][j][k] * dt / (2 * h) * (v[i][j][k+1] - v[i][j][k-1]) - dt / (2 * rho * h) * (p[i][j+1][k] - p[i][j-1][k]) + nu * (dt / (h * h) * (v[i+1][j][k] - 2 * v[i][j][k] + v[i-1][j][k]) + dt / (h * h) * (v[i][j+1][k] - 2 * v[i][j][k] + v[i][j-1][k]) + dt / (h * h) * (v[i][j][k+1] - 2 * v[i][j][k] + v[i][j][k-1]));

v[0][j][k] = 0.0;

v[num_i - 1][j][k] = 0.0;

v[i][0][k] = 0.0;

v[i][num_j - 1][k] = 0.0;

v[i][j][0] = 0.0;

v[i][j][num_k - 1] = 0.0;

w[i][j][k] = w[i][j][k] - u[i][j][k] * dt / (2 * h) * (w[i+1][j][k] - w[i-1][j][k]) - v[i][j][k] * dt / (2 * h) * (w[i][j+1][k] - w[i][j-1][k]) - w[i][j][k] * dt / (2 * h) * (w[i][j][k+1] - w[i][j][k-1]) - dt / (2 * rho * h) * (p[i][j][k+1] - p[i][j][k-1]) + nu * (dt / (h * h) * (w[i+1][j][k] - 2 * w[i][j][k] + w[i-1][j][k]) + dt / (h * h) * (w[i][j+1][k] - 2 * w[i][j][k] + w[i][j-1][k]) + dt / (h * h) * (w[i][j][k+1] - 2 * w[i][j][k] + w[i][j][k-1]));

w[0][j][k] = 0.0;

w[num_i - 1][j][k] = 0.0;

w[i][0][k] = 0.0;

w[i][num_j - 1][k] = 0.0;

w[i][j][0] = 0.0;

w[i][j][num_k - 1] = 0.0;


Fig. 1, Velocity and pressure iso-surfaces


     Of course constants need to be defined, such as grid spacing in space and time, density, kinematic viscosity. These equations have been validated, as you might have read and here already! Happy coding!

     If you want to hire me as your PhD student in the research projects related to turbo-machinery, aerodynamics, renewable energy, please reach out. Thank you very much for reading.

Monday 14 August 2023

Heaving Flat Plate Computational Fluid Dynamics (CFD) Simulation (In-House CFD Code)

     The adventure 🏕 🚵 I started a while ago to make my own CFD 🌬 code / software 💻 for my digital CV and to make a shiny new turbulence mode (one-day perhaps) is going along nicely. This post is about a 2D 10% flat plate undergoing forced heaving motion. Heaving motion is achieved by Eqn. 1.


hy = Ho*sin(2*π*fh*t)                                              Eqn. 1

     w.r.t. Eqn. 1 reduced frequency is defined as (fh*Ho/U∞)), fh is the frequency of oscillations, while t is the instantaneous time. Ho is the heaving amplitude and U∞ is the free stream velocity. hy is the position of the flat plate. The animation is shown in Fig. 1.

     The Strouhal number is 0.228 and Reynolds number is at 500. As we can see, the in-house CFD code works very well for this complex CFD simulation. Validation of this work will never be completed 😆. As soon, I will move on to the next project without completing this one. Anyway, discretized Navier-Stokes equations are available here, in both C++ and MATLAB formats if you want to validate this non sense yourself! Good Luck!

The animation from in-house CFD simulation

     If you want to hire me as your PhD student in the research projects related to turbo-machinery, aerodynamics, renewable energy, please reach out. Thank you very much for reading.