im following Dr. Lorena Barbara's 12 steps to Navier stokes equation (http://lorenabarba.com/blog/cfd-python-12-steps-to-navier-stokes/) structure but I'm not sure how to do the 2d linear convection. I was wondering if anyone here would be familiar on how to do it. Below is my the sample code:
#include <iostream>
#include <fstream>
#include <string>
#include <cstring>
#include <cstdlib>
#include <typeinfo>
#include <sstream>
#include <cmath>
void linspace_2d(double a, double b, double c, double d, double ** array){
double delta =(b-a)/(c-1);
for (int i=0; i<c; ++i){
for (int j=0; j<d; ++j){
array[i][j]= j*delta;
}
}
}
void convection_2d(const std::string& str, const int nx, const int ny, const int nt){
double c=1.; // speed
double dx=2.0/(nx-1.); // grid distance in x direction
double dy=2.0/(ny-1.); // grid distance in y direction
double sigma=0.2;
double dt=sigma*dx; // time step
double **space; // mesh grid
// Alocate memory
space = new double *[nx];
for (int i=0; i<nx; ++i){ space[i] = new double[ny]; }
double a=0, b=2;
linspace_2d(a, b, nx, ny, space); // function creates gives values to the mesh
// Initialize the array u and
double **u;
u = new double *[nx];
for (int i=0; i<nx; ++i){ u[i] = new double[ny]; }
// Set inital conditions
for (int i=0; i<nx; ++i){
for (int j=0; j<ny; ++j){
u[i][j] = 1.;
if ( ((double) i>=0.5/dx) && ((double) i<(1./dx+1.)) && ((double) j>=0.5/dy) && ((double) j<(1./dy+1.)) ){
u[i][j] = 2.;
}
}
// Iteration
for (int t=0; t<nt; ++t){
// Copy elements of array u into array un
double **un;
un = new double *[nx];
for (int x=0; x<nx; ++x){ un[x] = new double[ny]; }
for (int x=0; x<nx; ++x){
for (int y=0; y<ny; ++y){
un[x][y] = u[x][y];
}
}
// take timestep
for (int i=0; i<nx; ++i){
for (int j=0; j<ny; ++j){
u[i][j] = un[i][j] - (c*dt/dx*(un[i][j] - un[i-1][j])) - (c*dt/dy*(un[i][j]-un[i][j-1]));
}
}
}
Accessing
un[i-1][j]
andun[i][j-1]
is out of array bounds wheni == 0
andj == 0
respectively, hence the crash.Also you leak memory in large amounts allocating
un
every loop iteration and never deleting it.