From a55065f4745f5f340eb6dffbd88fe2fb05a40625 Mon Sep 17 00:00:00 2001 From: Ian Jauslin Date: Wed, 25 May 2022 11:12:02 -0400 Subject: add L parameter --- src/main.c | 30 ++++++++++++++++++++---------- src/navier-stokes.c | 26 +++++++++++++++----------- src/navier-stokes.h | 15 ++++++++------- 3 files changed, 43 insertions(+), 28 deletions(-) (limited to 'src') diff --git a/src/main.c b/src/main.c index 27bfb7a..3d8867e 100644 --- a/src/main.c +++ b/src/main.c @@ -14,8 +14,8 @@ int print_usage(); // read command line arguments int read_args(int argc, const char* argv[], char** params, unsigned int* driving_force, unsigned int* command, unsigned int* nthreads); -int read_params(char* params, int* K1, int* K2, int* N1, int* N2, unsigned int* nsteps, double* nu, double* delta, unsigned int* print_freq); -int set_parameter(char* lhs, char* rhs, int* K1, int* K2, int* N1, int* N2, unsigned int* nsteps, double* nu, double* delta, unsigned int* print_freq, bool* setN1, bool* setN2); +int read_params(char* params, int* K1, int* K2, int* N1, int* N2, unsigned int* nsteps, double* nu, double* delta, double* L, unsigned int* print_freq); +int set_parameter(char* lhs, char* rhs, int* K1, int* K2, int* N1, int* N2, unsigned int* nsteps, double* nu, double* delta, double* L, unsigned int* print_freq, bool* setN1, bool* setN2); #define COMMAND_UK 1 @@ -35,7 +35,7 @@ int main ( int K1,K2; int N1,N2; unsigned int nsteps; - double nu,delta; + double nu,delta,L; _Complex double (*g)(int,int); int ret; unsigned int driving,command; @@ -51,7 +51,7 @@ int main ( return(-1); } // read params - ret=read_params(params, &K1, &K2, &N1, &N2, &nsteps, &nu, &delta, &print_freq); + ret=read_params(params, &K1, &K2, &N1, &N2, &nsteps, &nu, &delta, &L, &print_freq); if(ret<0){ return(-1); } @@ -71,16 +71,16 @@ int main ( // run command if (command==COMMAND_UK){ - uk(K1, K2, N1, N2, nsteps, nu, delta, g, print_freq, nthreads); + uk(K1, K2, N1, N2, nsteps, nu, delta, L, g, print_freq, nthreads); } else if (command==COMMAND_ENERGY){ - energy(K1, K2, N1, N2, nsteps, nu, delta, g, print_freq, nthreads); + energy(K1, K2, N1, N2, nsteps, nu, delta, L, g, print_freq, nthreads); } else if(command==COMMAND_ENSTROPHY){ - energy(K1, K2, N1, N2, nsteps, nu, delta, g, print_freq, nthreads); + energy(K1, K2, N1, N2, nsteps, nu, delta, L, g, print_freq, nthreads); } else if(command==COMMAND_QUIET){ - quiet(K1, K2, N1, N2, nsteps, nu, delta, g, nthreads); + quiet(K1, K2, N1, N2, nsteps, nu, delta, L, g, nthreads); } else if(command==0){ fprintf(stderr, "error: no command specified\n"); @@ -201,6 +201,7 @@ int read_params( unsigned int* nsteps, double* nu, double* delta, + double* L, unsigned int* print_freq ){ int ret; @@ -222,6 +223,7 @@ int read_params( *delta=0.0001220703125; //nu=2^-11 *nu=0.00048828125; + *L=2*M_PI; *nsteps=10000000; *print_freq=1000; @@ -244,7 +246,7 @@ int read_params( break; case ';': //set parameter - ret=set_parameter(buffer_lhs,buffer_rhs,K1,K2,N1,N2,nsteps,nu,delta,print_freq,&setN1,&setN2); + ret=set_parameter(buffer_lhs,buffer_rhs,K1,K2,N1,N2,nsteps,nu,delta,L,print_freq,&setN1,&setN2); if(ret<0){ return ret; } @@ -271,7 +273,7 @@ int read_params( // set last param if (*params!='\0'){ - ret=set_parameter(buffer_lhs,buffer_rhs,K1,K2,N1,N2,nsteps,nu,delta,print_freq,&setN1,&setN2); + ret=set_parameter(buffer_lhs,buffer_rhs,K1,K2,N1,N2,nsteps,nu,delta,L,print_freq,&setN1,&setN2); if(ret<0){ return ret; } @@ -305,6 +307,7 @@ int set_parameter( unsigned int* nsteps, double* nu, double* delta, + double* L, unsigned int* print_freq, bool* setN1, bool* setN2 @@ -380,6 +383,13 @@ int set_parameter( return(-1); } } + else if (strcmp(lhs,"L")==0){ + ret=sscanf(rhs,"%lf",L); + if(ret!=1){ + fprintf(stderr, "error: parameter 'L' should be a double\n got '%s'\n",rhs); + return(-1); + } + } else if (strcmp(lhs,"print_freq")==0){ ret=sscanf(rhs,"%u",print_freq); if(ret!=1){ diff --git a/src/navier-stokes.c b/src/navier-stokes.c index 99e41ab..89bfbbd 100644 --- a/src/navier-stokes.c +++ b/src/navier-stokes.c @@ -2,8 +2,6 @@ #include #include -#define M_PI 3.14159265358979323846 - // compute solution as a function of time int uk( int K1, @@ -13,6 +11,7 @@ int uk( unsigned int nsteps, double nu, double delta, + double L, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads @@ -44,7 +43,7 @@ int uk( // iterate for(t=0;t #include +#define M_PI 3.14159265358979323846 // arrays on which the ffts are performed typedef struct fft_vects { @@ -12,16 +13,16 @@ typedef struct fft_vects { } fft_vect; // compute u_k -int uk( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads); +int uk( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, double L, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads); // compute the energy as a function of time -int energy( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads); +int energy( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, double L, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads); // compute enstrophy -int enstrophy( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads); +int enstrophy( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, double L, _Complex double (*g)(int,int), unsigned int print_freq, unsigned int nthreads); // compute solution as a function of time, but do not print anything (useful for debugging) -int quiet( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, _Complex double (*g)(int,int), unsigned int nthreads); +int quiet( int K1, int K2, int N1, int N2, unsigned int nsteps, double nu, double delta, double L, _Complex double (*g)(int,int), unsigned int nthreads); // initialize vectors for computation @@ -33,15 +34,15 @@ int ns_free_tmps( _Complex double* u, _Complex double* tmp1, _Complex double *tm int ns_init_u( _Complex double* u, int K1, int K2); // next time step for Irreversible Navier-Stokes equation -int ins_step( _Complex double* u, int K1, int K2, int N1, int N2, double nu, double delta, _Complex double (*g)(int,int), fft_vect fft1, fft_vect fft2,fft_vect ifft, _Complex double* tmp1, _Complex double *tmp2, _Complex double *tmp3); +int ins_step( _Complex double* u, int K1, int K2, int N1, int N2, double nu, double delta, double L, _Complex double (*g)(int,int), fft_vect fft1, fft_vect fft2,fft_vect ifft, _Complex double* tmp1, _Complex double *tmp2, _Complex double *tmp3); // right side of Irreversible Navier-Stokes equation -int ins_rhs( _Complex double* out, _Complex double* u, int K1, int K2, int N1, int N2, double nu, _Complex double (*g)(int,int), fft_vect fft1, fft_vect fft2, fft_vect ifft); +int ins_rhs( _Complex double* out, _Complex double* u, int K1, int K2, int N1, int N2, double nu, double L, _Complex double (*g)(int,int), fft_vect fft1, fft_vect fft2, fft_vect ifft); // convolution term in right side of equation int ns_T( _Complex double* u, int K1, int K2, int N1, int N2, fft_vect fft1, fft_vect fft2, fft_vect ifft); -// convolution term in right side of equation +// convolution term in right side of equation (computed without fft) int ns_T_nofft( _Complex double* out, _Complex double* u, int K1, int K2, int N1, int N2); // compute alpha -- cgit v1.2.3-54-g00ecf