Changeset 3518 in Sophya for trunk/Cosmo

Timestamp:

Sep 11, 2008, 12:37:39 PM (17 years ago)

Author:

cmv

Message:

possibilite de travailler en float cmv 11/09/2008

Location:

trunk/Cosmo/SimLSS

Files:

• Makefile (modified) (1 diff)
• cmvobserv3d.cc (modified) (1 diff)
• genefluct3d.cc (modified) (9 diffs)
• genefluct3d.h (modified) (4 diffs)

trunk/Cosmo/SimLSS/Makefile

@@ -12 +12 @@
 MYEXTINC = ${EXTLIBDIR}/Include
 #MYLIB = $(SOPHYAEXTSLBLIST) -L$(LIB) -lcmvsimbao -lfftw3 -lm
-MYLIB = $(SOPHYAEXTSLBLIST) -L$(LIB) -lcmvsimbao -lfftw3_threads -lfftw3 -lm
+MYLIB = $(SOPHYAEXTSLBLIST) -L$(LIB) -lcmvsimbao -lfftw3_threads -lfftw3f_threads -lfftw3 -lfftw3f -lm
 
 #--------------------------------------------------------------------------

trunk/Cosmo/SimLSS/cmvobserv3d.cc

@@ -379 +379 @@
  fluct3d.SetGrowthFactor(growth);
  fluct3d.LosComRedshift(0.001,-1);
- TArray< complex<r_8> >& pkgen = fluct3d.GetComplexArray();
- TArray<r_8>& rgen = fluct3d.GetRealArray();
+ TArray< complex<GEN3D_TYPE> >& pkgen = fluct3d.GetComplexArray();
+ TArray<GEN3D_TYPE>& rgen = fluct3d.GetRealArray();
  cout<<endl; fluct3d.Print();
  cout<<"\nMean number of galaxies per pixel = "<<ngal_by_mpc3*fluct3d.GetDVol()<<endl;

trunk/Cosmo/SimLSS/genefluct3d.cc

@@ -24 +24 @@
 #include "genefluct3d.h"
 
-#define WITH_FFTW_THREAD
+#if defined(GEN3D_FLOAT)
+#define GEN3D_FFTW_INIT_THREADS       fftwf_init_threads
+#define GEN3D_FFTW_CLEANUP_THREADS    fftwf_cleanup_threads
+#define GEN3D_FFTW_PLAN_WITH_NTHREADS fftwf_plan_with_nthreads
+#define GEN3D_FFTW_PLAN_DFT_R2C_3D    fftwf_plan_dft_r2c_3d
+#define GEN3D_FFTW_PLAN_DFT_C2R_3D    fftwf_plan_dft_c2r_3d
+#define GEN3D_FFTW_DESTROY_PLAN       fftwf_destroy_plan
+#define GEN3D_FFTW_EXECUTE            fftwf_execute
+#else
+#define GEN3D_FFTW_INIT_THREADS       fftw_init_threads
+#define GEN3D_FFTW_CLEANUP_THREADS    fftw_cleanup_threads
+#define GEN3D_FFTW_PLAN_WITH_NTHREADS fftw_plan_with_nthreads
+#define GEN3D_FFTW_PLAN_DFT_R2C_3D    fftw_plan_dft_r2c_3d
+#define GEN3D_FFTW_PLAN_DFT_C2R_3D    fftw_plan_dft_c2r_3d
+#define GEN3D_FFTW_DESTROY_PLAN       fftw_destroy_plan
+#define GEN3D_FFTW_EXECUTE            fftw_execute
+#endif
 
 #define MODULE2(_x_) ((double)((_x_).real()*(_x_).real() + (_x_).imag()*(_x_).imag()))
@@ -64 +80 @@
 {
  Nx_ = Ny_ = Nz_ = 0;
- is_set_fftw_plan = false;
+ is_set_fft_plan = false;
  nthread_ = 0;
  lp_ = 0;
@@ -134 +150 @@
    T_.ReSize(3,SzK_);
    array_allocated_ = true;
-   if(lp_>1) cout<<"  allocating: "<<T_.Size()*sizeof(complex<r_8>)/1.e6<<" Mo"<<endl;
+   if(lp_>1) cout<<"  allocating: "<<T_.Size()*sizeof(complex<GEN3D_TYPE>)/1.e6<<" Mo"<<endl;
  } catch (...) {
    cout<<"GeneFluct3D::setalloc_Error: Problem allocating T_"<<endl;
@@ -147 +163 @@
    else        R_ = ArrCastC2R(T_);
  // On remplit les pointeurs
- fdata_ = (fftw_complex *) (&T_(0,0,0));
- data_ = (double *) (&R_(0,0,0));
+ fdata_ = (GEN3D_FFTW_COMPLEX *) (&T_(0,0,0));
+ data_ = (GEN3D_TYPE *) (&R_(0,0,0));
 }
 
 void GeneFluct3D::init_fftw(void)
 {
- if( is_set_fftw_plan ) delete_fftw();
+ if( is_set_fft_plan ) delete_fftw();
 
  // --- Initialisation de fftw3 (attention data est sur-ecrit a l'init)
@@ -160 +176 @@
  if(nthread_>0) {
    cout<<"...Computing with "<<nthread_<<" threads"<<endl;
-   fftw_init_threads();
-   fftw_plan_with_nthreads(nthread_);
+   GEN3D_FFTW_INIT_THREADS();
+   GEN3D_FFTW_PLAN_WITH_NTHREADS(nthread_);
  }
 #endif
  if(lp_>1) cout<<"...forward plan"<<endl;
- pf_ = fftw_plan_dft_r2c_3d(Nx_,Ny_,Nz_,data_,fdata_,FFTW_ESTIMATE);
+ pf_ = GEN3D_FFTW_PLAN_DFT_R2C_3D(Nx_,Ny_,Nz_,data_,fdata_,FFTW_ESTIMATE);
  if(lp_>1) cout<<"...backward plan"<<endl;
- pb_ = fftw_plan_dft_c2r_3d(Nx_,Ny_,Nz_,fdata_,data_,FFTW_ESTIMATE);
- is_set_fftw_plan = true;
+ pb_ = GEN3D_FFTW_PLAN_DFT_C2R_3D(Nx_,Ny_,Nz_,fdata_,data_,FFTW_ESTIMATE);
+ is_set_fft_plan = true;
 }
 
 void GeneFluct3D::delete_fftw(void)
 {
- if( !is_set_fftw_plan ) return;
- fftw_destroy_plan(pf_);
- fftw_destroy_plan(pb_);
+ if( !is_set_fft_plan ) return;
+ GEN3D_FFTW_DESTROY_PLAN(pf_);
+ GEN3D_FFTW_DESTROY_PLAN(pb_);
 #ifdef WITH_FFTW_THREAD
- if(nthread_>0) fftw_cleanup_threads();
+ if(nthread_>0) GEN3D_FFTW_CLEANUP_THREADS();
 #endif
- is_set_fftw_plan = false;
+ is_set_fft_plan = false;
 }
 
@@ -566 +582 @@
  // --- realisation d'un tableau de tirage gaussiens
  if(lp_>0) cout<<"...before fft real ---"<<endl;
- fftw_execute(pf_);
+ GEN3D_FFTW_EXECUTE(pf_);
 
  // --- On remplit avec une realisation
  if(lp_>0) cout<<"...before Fourier realization filling"<<endl;
- T_(0,0,0) = complex<r_8>(0.);  // on coupe le continue et on l'initialise
+ T_(0,0,0) = complex<GEN3D_TYPE>(0.);  // on coupe le continue et on l'initialise
  long lmod = Nx_/10; if(lmod<1) lmod=1;
  for(long i=0;i<Nx_;i++) {
@@ -609 +625 @@
 {
  // --- RaZ du tableau
- T_ = complex<r_8>(0.);
+ T_ = complex<GEN3D_TYPE>(0.);
 
  // --- On remplit avec une realisation
@@ -830 +846 @@
 
  // On fait la FFT
- fftw_execute(pb_);
+ GEN3D_FFTW_EXECUTE(pb_);
 }
 
@@ -840 +856 @@
 
  // On fait la FFT
- fftw_execute(pf_);
+ GEN3D_FFTW_EXECUTE(pf_);
  // On corrige du pb de la normalisation de FFTW3
  double v = (double)NRtot_;

trunk/Cosmo/SimLSS/genefluct3d.h

@@ -18 +18 @@
 #include "pkspectrum.h"
 
+#define WITH_FFTW_THREAD
+//#define GEN3D_FLOAT
+
+#if defined(GEN3D_FLOAT)
+#define GEN3D_TYPE r_4
+#define GEN3D_FFTW_PLAN fftwf_plan
+#define GEN3D_FFTW_COMPLEX fftwf_complex
+#else
+#define GEN3D_TYPE r_8
+#define GEN3D_FFTW_PLAN fftw_plan
+#define GEN3D_FFTW_COMPLEX fftw_complex
+#endif
 
 namespace SOPHYA {
@@ -41 +53 @@
   long LosComRedshift(double zinc=0.001,long npoints=-1);
 
-  TArray< complex<r_8> >& GetComplexArray(void) {return T_;}
-  fftw_complex* GetComplexPointer(void) {return fdata_;}
-  TArray<r_8>& GetRealArray(void) {return R_;}
-  r_8* GetRealPointer(void) {return data_;}
+  TArray< complex<GEN3D_TYPE> >& GetComplexArray(void) {return T_;}
+  GEN3D_FFTW_COMPLEX * GetComplexPointer(void) {return fdata_;}
+  TArray<GEN3D_TYPE>& GetRealArray(void) {return R_;}
+  GEN3D_TYPE* GetRealPointer(void) {return data_;}
 
   // Pour adressage data_[ip]
@@ -157 +169 @@
 
   // la gestion de la FFT
-  bool is_set_fftw_plan;
-  fftw_plan pf_,pb_;
+  bool is_set_fft_plan;
+  GEN3D_FFTW_PLAN pf_,pb_;
   unsigned short nthread_;
   int lp_;
@@ -164 +176 @@
   // le stockage du Cube de donnees et les pointeurs
   bool array_allocated_;  // true if array has been allocated
-  TArray< complex<r_8> > T_;
-  fftw_complex *fdata_;
-  TArray<r_8> R_;
-  double *data_;
+  TArray< complex<GEN3D_TYPE> > T_;
+  GEN3D_FFTW_COMPLEX *fdata_;
+  TArray<GEN3D_TYPE> R_;
+  GEN3D_TYPE *data_;
 
   // l'observateur