software/Nstrophy
1
Fork 0
Code Issues Pull Requests Releases Activity

Values that seem to work

Browse Source
This commit is contained in:
Ian Jauslin 2018-02-19 16:00:07 +00:00
parent b845911397
commit cbf2b9854a
1 changed files with 47 additions and 16 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

63
src/main.c
View File

@ -72,15 +72,15 @@ int read_args(int argc, const char* argv[], ns_params* params, unsigned int* nst
// defaults
/*
params->h=6.103515625e-05;
params->K=16;
*nsteps=16777216;
params->nu=4.9632717887631524e-05;
*/
params->K=16;
params->h=1e-3/(2*params->K+1);
*nsteps=10000000;
params->nu=1./1024/(2*params->K+1);
*/
params->K=3;
params->h=0.0001220703125;
*nsteps=10000000;
params->nu=0.00048828125;
// loop over arguments
for(i=1;i<argc;i++){
@ -190,6 +190,7 @@ int enstrophy(ns_params params, unsigned int Nsteps){
unsigned int t;
int kx,ky;
fft_vects fft_vects;
double rescale;
// sizes
params.S=2*params.K+1;
@ -203,25 +204,55 @@ int enstrophy(ns_params params, unsigned int Nsteps){
tmp2=calloc(sizeof(_Complex double),params.S*params.S);
tmp3=calloc(sizeof(_Complex double),params.S*params.S);
srand(17);
// initial value
for(kx=-params.K;kx<=params.K;kx++){
for(ky=0;ky<=params.K;ky++){
u[KLOOKUP(0,ky,params.S)]=(-RAND_MAX*0.5+rand())*1.0/RAND_MAX+(-RAND_MAX*0.5+rand())*1.0/RAND_MAX*I;
}
for(kx=1;kx<=params.K;kx++){
for(ky=-params.K;ky<=params.K;ky++){
//u[KLOOKUP(kx,ky,params.S)]=kx*kx*ky*ky*exp(-(kx*kx+ky*ky));
if((kx==1 && ky==0) || (kx==-1 && ky==0)){
u[KLOOKUP(kx,ky,params.S)]=1;
}
else{
u[KLOOKUP(kx,ky,params.S)]=0;
}
u[KLOOKUP(kx,ky,params.S)]=(-RAND_MAX*0.5+rand())*1.0/RAND_MAX+(-RAND_MAX*0.5+rand())*1.0/RAND_MAX*I;
}
}
for(ky=-params.K;ky<=-1;ky++){
u[KLOOKUP(0,ky,params.S)]=conj(u[KLOOKUP(0,-ky,params.S)]);
}
for(kx=-params.K;kx<=-1;kx++){
for(ky=-params.K;ky<=params.K;ky++){
u[KLOOKUP(kx,ky,params.S)]=conj(u[KLOOKUP(-kx,-ky,params.S)]);
}
}
rescale=0;
for(kx=-params.K;kx<=params.K;kx++){
for(ky=-params.K;ky<=params.K;ky++){
rescale=rescale+((__real__ u[KLOOKUP(kx,ky,params.S)])*(__real__ u[KLOOKUP(kx,ky,params.S)])+(__imag__ u[KLOOKUP(kx,ky,params.S)])*(__imag__ u[KLOOKUP(kx,ky,params.S)]))*(kx*kx+ky*ky);
}
}
for(kx=-params.K;kx<=params.K;kx++){
for(ky=-params.K;ky<=params.K;ky++){
u[KLOOKUP(kx,ky,params.S)]=u[KLOOKUP(kx,ky,params.S)]*sqrt(155.1/rescale);
}
}
/*
for(kx=-params.K;kx<=params.K;kx++){
for(ky=-params.K;ky<=params.K;ky++){
printf("%d %d % .8e % .8e\n",kx,ky, __real__ u[KLOOKUP(kx,ky,params.S)], __imag__ u[KLOOKUP(kx,ky,params.S)]);
}
}
*/
// driving force
for(kx=-params.K;kx<=params.K;kx++){
for(ky=-params.K;ky<=params.K;ky++){
//params.g[KLOOKUP(kx,ky,params.S)]=sqrt(kx*kx*ky*ky)*exp(-(kx*kx+ky*ky));
if((kx==2 && ky==-1) || (kx==-2 && ky==1)){
params.g[KLOOKUP(kx,ky,params.S)]=1.0*params.K;
if(kx==2 && ky==-1){
params.g[KLOOKUP(kx,ky,params.S)]=0.5+sqrt(3)/2*I;
}
else if(kx==-2 && ky==1){
params.g[KLOOKUP(kx,ky,params.S)]=0.5-sqrt(3)/2*I;
}
else{
params.g[KLOOKUP(kx,ky,params.S)]=0;
@ -260,7 +291,7 @@ int enstrophy(ns_params params, unsigned int Nsteps){
avg=avg-(avg-alpha)/t;
}
if(t>0 && t%1000==0){
if(t>0 && t%8192==0){
fprintf(stderr,"%8d % .8e % .8e % .8e % .8e\n",t, __real__ avg, __imag__ avg, __real__ alpha, __imag__ alpha);
printf("%8d % .8e % .8e % .8e % .8e\n",t, __real__ avg, __imag__ avg, __real__ alpha, __imag__ alpha);
}