diff --git a/src/fluid_solver.cpp b/src/fluid_solver.cpp
index 211a2550a03e4a05e0c8f52d7aea818d6e102c9b..9fff1689482e54c2327177465fac681b488ab787 100644
--- a/src/fluid_solver.cpp
+++ b/src/fluid_solver.cpp
@@ -51,12 +51,14 @@ void fluid_solver<rnumber>::impose_zero_modes()
\
template<> \
fluid_solver<R>::fluid_solver( \
+ const char *NAME, \
int nx, \
int ny, \
int nz, \
double DKX, \
double DKY, \
- double DKZ) : fluid_solver_base<R>(nx , ny , nz, \
+ double DKZ) : fluid_solver_base<R>(NAME, \
+ nx , ny , nz, \
DKX, DKY, DKZ) \
{ \
this->cvorticity = FFTW(alloc_complex)(this->cd->local_size);\
@@ -306,6 +308,52 @@ void fluid_solver<R>::step(double dt) \
this->force_divfree(this->cv[0]); \
\
this->iteration++; \
+} \
+ \
+template<> \
+int fluid_solver<R>::read(char field, char representation) \
+{ \
+ if ((field == 'v') && (representation == 'c')) \
+ return this->read_base("cvorticity", this->cvorticity); \
+ if ((field == 'v') && (representation == 'r')) \
+ { \
+ int read_result = this->read_base("rvorticity", this->rvorticity); \
+ if (!(read_result == EXIT_SUCCESS)) \
+ return read_result; \
+ else \
+ { \
+ FFTW(execute)(*((FFTW(plan)*)this->r2c_vorticity )); \
+ return EXIT_SUCCESS; \
+ } \
+ } \
+ if ((field == 'u') && (representation == 'c')) \
+ return this->read_base("cvelocity", this->cvelocity); \
+ if ((field == 'u') && (representation == 'r')) \
+ return this->read_base("rvelocity", this->rvelocity); \
+ return EXIT_FAILURE; \
+} \
+ \
+template<> \
+int fluid_solver<R>::write(char field, char representation) \
+{ \
+ if ((field == 'v') && (representation == 'c')) \
+ return this->write_base("cvorticity", this->cvorticity); \
+ if ((field == 'v') && (representation == 'r')) \
+ { \
+ FFTW(execute)(*((FFTW(plan)*)this->c2r_vorticity )); \
+ clip_zero_padding<R>(this->rd, this->rvorticity, 3); \
+ return this->write_base("rvorticity", this->rvorticity); \
+ } \
+ this->compute_velocity(this->cvorticity); \
+ if ((field == 'u') && (representation == 'c')) \
+ return this->write_base("cvelocity", this->cvelocity); \
+ if ((field == 'u') && (representation == 'r')) \
+ { \
+ FFTW(execute)(*((FFTW(plan)*)this->c2r_velocity )); \
+ clip_zero_padding<R>(this->rd, this->rvelocity, 3); \
+ return this->write_base("rvelocity", this->rvelocity); \
+ } \
+ return EXIT_FAILURE; \
}
/*****************************************************************************/
@@ -319,12 +367,12 @@ FLUID_SOLVER_DEFINITIONS(
fftwf_complex,
MPI_REAL4,
MPI_COMPLEX8)
-FLUID_SOLVER_DEFINITIONS(
- FFTW_MANGLE_DOUBLE,
- double,
- fftw_complex,
- MPI_REAL8,
- MPI_COMPLEX16)
+//FLUID_SOLVER_DEFINITIONS(
+// FFTW_MANGLE_DOUBLE,
+// double,
+// fftw_complex,
+// MPI_REAL8,
+// MPI_COMPLEX16)
/*****************************************************************************/
diff --git a/src/fluid_solver.hpp b/src/fluid_solver.hpp
index 7d5dc402ffd075c364da5e8590e4e67998d84c7d..f673f61e5898730e5ed12e9770c1452da8872f78 100644
--- a/src/fluid_solver.hpp
+++ b/src/fluid_solver.hpp
@@ -64,6 +64,7 @@ class fluid_solver:public fluid_solver_base<rnumber>
/* methods */
fluid_solver(
+ const char *NAME,
int nx,
int ny,
int nz,
@@ -77,6 +78,9 @@ class fluid_solver:public fluid_solver_base<rnumber>
void omega_nonlin(int src);
void step(double dt);
void impose_zero_modes(void);
+
+ int read(char field, char representation);
+ int write(char field, char representation);
};
#endif//FLUID_SOLVER
diff --git a/src/fluid_solver_base.cpp b/src/fluid_solver_base.cpp
index b38492b4ed4d94d327c01a9ad9de4e3bc87ff7da..ec82968d2ff55fe8ac87d7910e9ac1985d5f5286 100644
--- a/src/fluid_solver_base.cpp
+++ b/src/fluid_solver_base.cpp
@@ -21,6 +21,7 @@
#include <cassert>
#include <cmath>
+#include <cstring>
#include "fluid_solver_base.hpp"
#include "fftw_tools.hpp"
@@ -33,6 +34,7 @@
\
template<> \
fluid_solver_base<R>::fluid_solver_base( \
+ const char *NAME, \
int nx, \
int ny, \
int nz, \
@@ -40,6 +42,8 @@ fluid_solver_base<R>::fluid_solver_base( \
double DKY, \
double DKZ) \
{ \
+ strncpy(this->name, NAME, 256); \
+ this->name[255] = '\0'; \
this->iteration = 0; \
\
int ntmp[4]; \
@@ -259,6 +263,38 @@ void fluid_solver_base<R>::symmetrize(C *data, const int howmany) \
FFTW(free)(buffer); \
delete mpistatus; \
} \
+ \
+template<> \
+int fluid_solver_base<R>::read_base(const char *fname, R *data) \
+{ \
+ char full_name[512]; \
+ sprintf(full_name, "%s_%s_i%.5x", this->name, fname, this->iteration); \
+ return this->rd->read(full_name, (void*)data); \
+} \
+ \
+template<> \
+int fluid_solver_base<R>::read_base(const char *fname, C *data) \
+{ \
+ char full_name[512]; \
+ sprintf(full_name, "%s_%s_i%.5x", this->name, fname, this->iteration); \
+ return this->cd->read(full_name, (void*)data); \
+} \
+ \
+template<> \
+int fluid_solver_base<R>::write_base(const char *fname, R *data) \
+{ \
+ char full_name[512]; \
+ sprintf(full_name, "%s_%s_i%.5x", this->name, fname, this->iteration); \
+ return this->rd->write(full_name, (void*)data); \
+} \
+ \
+template<> \
+int fluid_solver_base<R>::write_base(const char *fname, C *data) \
+{ \
+ char full_name[512]; \
+ sprintf(full_name, "%s_%s_i%.5x", this->name, fname, this->iteration); \
+ return this->cd->write(full_name, (void*)data); \
+}
/*****************************************************************************/
@@ -272,12 +308,12 @@ FLUID_SOLVER_BASE_DEFINITIONS(
fftwf_complex,
MPI_REAL4,
MPI_COMPLEX8)
-FLUID_SOLVER_BASE_DEFINITIONS(
- FFTW_MANGLE_DOUBLE,
- double,
- fftw_complex,
- MPI_REAL8,
- MPI_COMPLEX16)
+//FLUID_SOLVER_BASE_DEFINITIONS(
+// FFTW_MANGLE_DOUBLE,
+// double,
+// fftw_complex,
+// MPI_REAL8,
+// MPI_COMPLEX16)
/*****************************************************************************/
diff --git a/src/fluid_solver_base.hpp b/src/fluid_solver_base.hpp
index 679a3342783ec294904b7eb2182c3090a4672378..aa08f275224f8e041d5d6973682290a998d3c2e5 100644
--- a/src/fluid_solver_base.hpp
+++ b/src/fluid_solver_base.hpp
@@ -46,6 +46,7 @@ class fluid_solver_base
ptrdiff_t normalization_factor;
/* simulation parameters */
+ char name[256];
int iteration;
/* physical parameters */
@@ -63,6 +64,7 @@ class fluid_solver_base
/* methods */
fluid_solver_base(
+ const char *NAME,
int nx,
int ny,
int nz,
@@ -75,6 +77,10 @@ class fluid_solver_base
void force_divfree(cnumber *a);
void symmetrize(cnumber *a, int howmany);
rnumber correl_vec(cnumber *a, cnumber *b);
+ int read_base(const char *fname, rnumber *data);
+ int read_base(const char *fname, cnumber *data);
+ int write_base(const char *fname, rnumber *data);
+ int write_base(const char *fname, cnumber *data);
};
diff --git a/test.py b/test.py
index 3f33eabecca8c81088687ebeca71927fd69fa589..95e11e49d85bd47b9f63366e78482a3574c75204 100755
--- a/test.py
+++ b/test.py
@@ -108,7 +108,7 @@ class stat_test(code):
//begincpp
fluid_solver<float> *fs;
char fname[512];
- fs = new fluid_solver<float>(nx, ny, nz);
+ fs = new fluid_solver<float>(simname, nx, ny, nz);
fs->nu = nu;
fs->iteration = iter0;
if (myrank == fs->cd->io_myrank)
@@ -116,28 +116,14 @@ class stat_test(code):
sprintf(fname, "%s_stats.bin", simname);
stat_file = fopen(fname, "wb");
}
-
- sprintf(fname, "%s_kvorticity_i%.5x", simname, fs->iteration);
- fs->cd->read(
- fname,
- (void*)fs->cvorticity);
+ fs->read('v', 'c');
fs->low_pass_Fourier(fs->cvorticity, 3, fs->kM);
fs->force_divfree(fs->cvorticity);
fs->symmetrize(fs->cvorticity, 3);
t = 0.0;
do_stats(fs);
- fftwf_execute(*((fftwf_plan*)fs->c2r_velocity));
- sprintf(fname, "%s_rvelocity_i%.5x", simname, fs->iteration);
- clip_zero_padding<float>(fs->rd, fs->rvelocity, 3);
- fs->rd->write(
- fname,
- (void*)fs->rvelocity);
- fftwf_execute(*((fftwf_plan*)fs->c2r_vorticity));
- sprintf(fname, "%s_rvorticity_i%.5x", simname, fs->iteration);
- clip_zero_padding<float>(fs->rd, fs->rvorticity, 3);
- fs->rd->write(
- fname,
- (void*)fs->rvorticity);
+ fs->write('u', 'r');
+ fs->write('v', 'r');
for (; fs->iteration < iter0 + niter_todo;)
{
fs->step(dt);
@@ -145,22 +131,9 @@ class stat_test(code):
do_stats(fs);
}
fclose(stat_file);
- sprintf(fname, "%s_kvorticity_i%.5x", simname, fs->iteration);
- fs->cd->write(
- fname,
- (void*)fs->cvorticity);
- fftwf_execute(*((fftwf_plan*)fs->c2r_vorticity));
- sprintf(fname, "%s_rvorticity_i%.5x", simname, fs->iteration);
- clip_zero_padding<float>(fs->rd, fs->rvorticity, 3);
- fs->rd->write(
- fname,
- (void*)fs->rvorticity);
- fftwf_execute(*((fftwf_plan*)fs->c2r_velocity));
- sprintf(fname, "%s_rvelocity_i%.5x", simname, fs->iteration);
- clip_zero_padding<float>(fs->rd, fs->rvelocity, 3);
- fs->rd->write(
- fname,
- (void*)fs->rvelocity);
+ fs->write('v', 'c');
+ fs->write('v', 'r');
+ fs->write('u', 'r');
delete fs;
//endcpp
"""
@@ -195,12 +168,12 @@ if __name__ == '__main__':
c.parameters['dt'] = 0.004
c.parameters['niter_todo'] = opt.nsteps
c.write_par(simname = 'test1')
- Kdata0.tofile("test1_kvorticity_i00000")
+ Kdata0.tofile("test1_cvorticity_i00000")
c.run(ncpu = opt.ncpu, simname = 'test1')
c.parameters['dt'] /= 2
c.parameters['niter_todo'] *= 2
c.write_par(simname = 'test2')
- Kdata0.tofile("test2_kvorticity_i00000")
+ Kdata0.tofile("test2_cvorticity_i00000")
c.run(ncpu = opt.ncpu, simname = 'test2')
Rdata = np.fromfile(
'test1_rvorticity_i00000',