CosmoBolognaLib
Free Software C++/Python libraries for cosmological calculations
Field3D.h
Go to the documentation of this file.
1 /*******************************************************************
2  * Copyright (C) 2010 by Federico Marulli and Alfonso Veropalumbo *
3  * federico.marulli3@unibo.it *
4  * *
5  * This program is free software; you can redistribute it and/or *
6  * modify it under the terms of the GNU General Public License as *
7  * published by the Free Software Foundation; either version 2 of *
8  * the License, or (at your option) any later version. *
9  * *
10  * This program is distributed in the hope that it will be useful, *
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
13  * GNU General Public License for more details. *
14  * *
15  * You should have received a copy of the GNU General Public *
16  * License along with this program; if not, write to the Free *
17  * Software Foundation, Inc., *
18  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
19  *******************************************************************/
20 
33 #ifndef __Field3D__
34 #define __Field3D__
35 
36 #include "Data.h"
37 
38 
39 namespace cbl {
40 
41  namespace data {
42 
51  class Field3D {
52 
53  protected:
54 
56  int m_nX;
57 
59  int m_nY;
60 
62  int m_nZ;
63 
65  int m_nZF;
66 
68  int m_nCells;
69 
72 
74  double m_deltaX;
75 
77  double m_deltaY;
78 
80  double m_deltaZ;
81 
83  double m_MinX;
84 
86  double m_MinY;
87 
89  double m_MinZ;
90 
92  double m_MaxX;
93 
95  double m_MaxY;
96 
98  double m_MaxZ;
99 
101  double m_Volume;
102 
104  std::vector<double> m_X;
105 
107  std::vector<double> m_Y;
108 
110  std::vector<double> m_Z;
111 
113  std::vector<double> m_kX;
114 
116  std::vector<double> m_kY;
117 
119  std::vector<double> m_kZ;
120 
130  long int inds_to_index (const int i, const int j, const int k) const { return k+m_nZ*(j+m_nY*i); }
131 
141  long int inds_to_index_Fourier (const int i, const int j, const int k) const { return k+m_nZF*(j+m_nY*i); }
142 
143 
144  public:
145 
150 
155  Field3D () = default;
156 
170  Field3D (const double deltaR, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
171 
187  Field3D (const int nx, const int ny, const int nz, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
188 
193  virtual ~Field3D () = default;
194 
196 
197 
202 
216  void set_parameters (const double deltaR, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
217 
233  void set_parameters (const int nx, const int ny, const int nz, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
234 
236 
237 
242 
247  int nx () const { return m_nX; }
248 
253  int ny() const {return m_nY;}
254 
259  int nz () const { return m_nZ; }
260 
265  int nzFourier () const { return m_nZF; }
266 
271  int nCells () const { return m_nCells; }
272 
277  int nCellsFourier () const { return m_nCells_Fourier; }
278 
283  double MinX () const { return m_MinX; }
284 
289  double MinY () const { return m_MinY; }
290 
295  double MinZ () const { return m_MinZ; }
296 
301  double MaxX () const { return m_MaxX; }
302 
307  double MaxY() const { return m_MaxY; }
308 
313  double MaxZ () const { return m_MaxZ; }
314 
319  double deltaX () const { return m_deltaX; }
320 
325  double deltaY () const { return m_deltaY; }
326 
331  double deltaZ () const { return m_deltaZ; }
332 
337  double Volume () const { return m_Volume; }
338 
346  double XX(const int i) const { return m_X[i];}
347 
355  double YY(const int i) const { return m_Y[i];}
356 
364  double ZZ(const int i) const { return m_Z[i];}
365 
374  double kX(const int i) const { return m_kX[i];}
375 
384  double kY(const int i) const { return m_kY[i];}
385 
394  double kZ(const int i) const { return m_kZ[i];}
395 
405  virtual double ScalarField (const int i, const int j, const int k) const
406  { (void)i; (void)j; (void)k; ErrorCBL("", "Scalarield", "Field3D.h"); return 0; }
407 
408 
416  virtual double ScalarField (const std::vector<double> pos) const
417  { (void)pos; ErrorCBL("", "Scalarield", "Field3D.h"); return 0; }
418 
424  virtual std::vector<double> ScalarField () const
425  { ErrorCBL("", "Scalarield", "Field3D.h"); std::vector<double> vv; return vv; }
426 
436  virtual std::vector<double> VectorField (const int i, const int j, const int k) const
437  { (void)i; (void)j; (void)k; ErrorCBL("", "VectorField", "Field3D.h"); return {0}; }
438 
446  virtual std::vector<double> VectorField (const std::vector<double> pos) const
447  { (void)pos; ErrorCBL("", "VectorField", "Field3D.h"); double vv = 0.; return {vv}; }
448 
458  virtual double ScalarField_FourierSpace_real (const int i, const int j, const int k) const
459  { (void)i; (void)j; (void)k; ErrorCBL("", "ScalarField_FourierSpace_real", "Field3D.h"); return 0; }
460 
470  virtual double ScalarField_FourierSpace_complex (const int i, const int j, const int k) const
471  { (void)i; (void)j; (void)k; ErrorCBL("", "ScalarField_FourierSpace_complex", "Field3D.h"); return 0; }
472 
478  virtual std::vector<double> ScalarField_FourierSpace_real () const
479  { ErrorCBL("", "Scalarield_FourierSpace_real", "Field3D.h"); std::vector<double> vv; return vv; }
480 
486  virtual std::vector<double> ScalarField_FourierSpace_complex () const
487  { ErrorCBL("", "Scalarield_FourierSpace_complex", "Field3D.h"); std::vector<double> vv; return vv; }
488 
498  virtual std::vector<double> VectorField_FourierSpace_real (const int i, const int j, const int k) const
499  { (void)i; (void)j; (void)k; ErrorCBL("", "VectorField_FourierSpace_real", "Field3D.h"); return {0}; }
500 
510  virtual std::vector<double> VectorField_FourierSpace_complex (const int i, const int j, const int k) const
511  { (void)i; (void)j; (void)k; ErrorCBL("", "VectorField_FourierSpace_complex", "Field3D.h"); return {0}; }
512 
514 
515 
520 
524  virtual void FourierTransformField ()
525  { ErrorCBL("", "FourierTransformField", "Field3D.h"); }
526 
531  { ErrorCBL("", "FourierAntiTransformField", "Field3D.h"); }
532 
537  virtual void GaussianConvolutionField (const double kernel_size)
538  { (void)kernel_size; ErrorCBL("", "GaussianConvolutionField", "Field3D.h"); }
539 
550  virtual void set_ScalarField (const double value, const int i, const int j, const int k, const bool add=0)
551  { (void)value; (void)i; (void)j; (void)k; (void)add; ErrorCBL("", "set_ScalarField", "Field3D.h"); }
552 
563  virtual void set_VectorField (const std::vector<double> value, const int i, const int j, const int k, const bool add=0)
564  { (void)value; (void)i; (void)j; (void)k; (void)add; ErrorCBL("", "set_Vectorield", "Field3D.h"); }
565 
577  virtual void set_ScalarField_FourierSpace_real (const double value, const int i, const int j, const int k, const bool add=0)
578  { (void)value; (void)i; (void)j; (void)k; (void)add; ErrorCBL("", "set_ScalarField_FourierSpace_real", "Field3D.h"); }
579 
591  virtual void set_ScalarField_FourierSpace_complex (const double value, const int i, const int j, const int k, const bool add=0)
592  { (void)value; (void)i; (void)j; (void)k; (void)add; ErrorCBL("", "set_ScalarField_FourierSpace_complex", "Field3D.h"); }
593 
606  virtual void set_VectorField_FourierSpace_real (const std::vector<double> value, const int i, const int j, const int k, const bool add=0)
607  { (void)value; (void)i; (void)j; (void)k; (void)add; ErrorCBL("", "set_VectorField_FourierSpace_real", "Field3D.h"); }
608 
621  virtual void set_VectorField_FourierSpace_complex (const std::vector<double> value, const int i, const int j, const int k, const bool add=0)
622  { (void)value; (void)i; (void)j; (void)k; (void)add; ErrorCBL("", "set_VectorField_FourierSpace_complex", "Field3D.h"); }
623 
627  virtual void reset()
628  { ErrorCBL("", "reset", "Field3D.h"); }
629 
631 
632  };
633 
634 
635  // ==================================================================================================================
636  // ==================================================================================================================
637  // ==================================================================================================================
638 
639 
648  class ScalarField3D : public Field3D {
649 
650  protected:
651 
653  double *m_field;
654 
656  fftw_complex *m_field_FourierSpace;
657 
658  public:
659 
664 
669  ScalarField3D () = default;
670 
684  ScalarField3D (const double deltaR, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
685 
701  ScalarField3D (const int nx, const int ny, const int nz, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
706  virtual ~ScalarField3D () = default;
707 
709 
710 
715 
728  void set_ScalarField (const double value, const int i, const int j, const int k, const bool add=0);
729 
744  void set_ScalarField_FourierSpace_real (const double value, const int i, const int j, const int k, const bool add=0);
745 
760  void set_ScalarField_FourierSpace_complex (const double value, const int i, const int j, const int k, const bool add=0);
761 
763 
764 
769 
779  double ScalarField (const int i, const int j, const int k) const;
780 
788  double ScalarField (const std::vector<double> pos) const;
789 
795  std::vector<double> ScalarField () const;
796 
807  double ScalarField_FourierSpace_real (const int i, const int j, const int k) const;
808 
819  double ScalarField_FourierSpace_complex (const int i, const int j, const int k) const;
820 
822 
823 
828 
832  void FourierTransformField ();
833 
839 
845  void GaussianConvolutionField (const double kernel_size);
846 
852  void reset();
853 
855  };
856 
857 
858  // ==================================================================================================================
859  // ==================================================================================================================
860  // ==================================================================================================================
861 
862 
871  class VectorField3D : public Field3D{
872  protected:
873 
875  std::vector<double *> m_field;
876 
878  std::vector<fftw_complex *> m_field_FourierSpace;
879 
880  public:
881 
886 
890  VectorField3D () = default;
891 
904  VectorField3D (const double deltaR, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
905 
920  VectorField3D (const int nx, const int ny, const int nz, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ);
925  ~VectorField3D () = default;
926 
928 
929 
934 
947  void set_VectorField (const std::vector<double> value, const int i, const int j, const int k, const bool add=0);
948 
961  void set_VectorField_FourierSpace_real (const std::vector<double> value, const int i, const int j, const int k, const bool add=0);
962 
975  void set_VectorField_FourierSpace_complex (const std::vector<double> value, const int i, const int j, const int k, const bool add=0);
976 
978 
979 
984 
994  std::vector<double> VectorField (const int i, const int j, const int k) const;
995 
1003  std::vector<double> VectorField (const std::vector<double> pos) const;
1004 
1015  std::vector<double> VectorField_FourierSpace_real (const int i, const int j, const int k) const;
1016 
1027  std::vector<double> VectorField_FourierSpace_complex (const int i, const int j, const int k) const;
1028 
1030 
1031 
1036 
1041  void FourierTransformField ();
1042 
1047  void FourierAntiTransformField ();
1048 
1054  void reset();
1055 
1057 
1058  };
1059 
1060  }
1061 }
1062 
1063 #endif
The class Data.
The class Field3D.
Definition: Field3D.h:51
double deltaX() const
get the private member m_deltaX
Definition: Field3D.h:319
double MaxY() const
get the private member m_MaxY
Definition: Field3D.h:307
int nCells() const
get the private member m_nCells
Definition: Field3D.h:271
virtual void FourierTransformField()
perform the Fourier transform on the field
Definition: Field3D.h:524
double kX(const int i) const
get the value of the X coordinates at the i-th cell, Fourier space
Definition: Field3D.h:374
double MinY() const
get the private member m_MinY
Definition: Field3D.h:289
double m_Volume
box volume
Definition: Field3D.h:101
int nzFourier() const
get the private member m_nZF
Definition: Field3D.h:265
int m_nCells
number of cells
Definition: Field3D.h:68
double m_MaxX
upper x bound
Definition: Field3D.h:92
double m_MinY
lower y bound
Definition: Field3D.h:86
virtual void set_ScalarField_FourierSpace_complex(const double value, const int i, const int j, const int k, const bool add=0)
set the value of the scalar field in Fourier space, complex part
Definition: Field3D.h:591
Field3D()=default
default constructor
virtual std::vector< double > VectorField_FourierSpace_complex(const int i, const int j, const int k) const
get the value of the vector field, Fourier space, complex part
Definition: Field3D.h:510
double kY(const int i) const
get the value of the Y coordinates at the i-th cell, Fourier space
Definition: Field3D.h:384
virtual void set_VectorField(const std::vector< double > value, const int i, const int j, const int k, const bool add=0)
set the value of the vectorr field
Definition: Field3D.h:563
std::vector< double > m_kX
coordinates of the cells along the kx-axis
Definition: Field3D.h:113
virtual void FourierAntiTransformField()
perform the anti-Fourier transform on the field
Definition: Field3D.h:530
double m_MinX
lower x bound
Definition: Field3D.h:83
virtual double ScalarField_FourierSpace_real(const int i, const int j, const int k) const
get the value of the scalar field, Fourier space, real part
Definition: Field3D.h:458
long int inds_to_index_Fourier(const int i, const int j, const int k) const
contract 3 indeces into one, Fourier space
Definition: Field3D.h:141
double XX(const int i) const
get the value of the X coordinates at the i-th cell
Definition: Field3D.h:346
double m_deltaY
Y cell size.
Definition: Field3D.h:77
virtual void reset()
set to 0 the fields
Definition: Field3D.h:627
double kZ(const int i) const
get the value of the Z coordinates at the i-th cell, Fourier space
Definition: Field3D.h:394
virtual double ScalarField(const int i, const int j, const int k) const
get the value of the scalar field
Definition: Field3D.h:405
double deltaZ() const
get the private member m_deltaZ
Definition: Field3D.h:331
int ny() const
get the private member m_nY
Definition: Field3D.h:253
int m_nX
number of cells along the x-axis
Definition: Field3D.h:56
double m_deltaX
X cell size.
Definition: Field3D.h:74
double m_MinZ
lower z bound
Definition: Field3D.h:89
double m_MaxY
upper y bound
Definition: Field3D.h:95
void set_parameters(const double deltaR, const double minX, const double maxX, const double minY, const double maxY, const double minZ, const double maxZ)
set the parameters
Definition: Field3D.cpp:63
int m_nZF
number of cells along the z-axis, Fourier space
Definition: Field3D.h:65
double MinX() const
get the private member m_MinX
Definition: Field3D.h:283
int nx() const
get the private member m_nX
Definition: Field3D.h:247
virtual void GaussianConvolutionField(const double kernel_size)
perform a smoothing of the field with a gaussian kernel
Definition: Field3D.h:537
double MaxZ() const
get the private member m_MaxZ
Definition: Field3D.h:313
virtual void set_VectorField_FourierSpace_complex(const std::vector< double > value, const int i, const int j, const int k, const bool add=0)
set the value of the vector field, Fourier space, complex part
Definition: Field3D.h:621
std::vector< double > m_kY
coordinates of the cells along the ky-axis
Definition: Field3D.h:116
std::vector< double > m_Y
coordinates of the cells along the y-axis
Definition: Field3D.h:107
virtual std::vector< double > ScalarField_FourierSpace_real() const
get the value of the scalar field, Fourier space, real part
Definition: Field3D.h:478
virtual void set_ScalarField(const double value, const int i, const int j, const int k, const bool add=0)
set the value of the scalar field
Definition: Field3D.h:550
std::vector< double > m_Z
coordinates of the cells along the z-axis
Definition: Field3D.h:110
virtual std::vector< double > ScalarField() const
get the value of the scalar field
Definition: Field3D.h:424
virtual ~Field3D()=default
default destructor
double ZZ(const int i) const
get the value of the Z coordinates at the i-th cell
Definition: Field3D.h:364
double deltaY() const
get the private member m_deltaY
Definition: Field3D.h:325
double m_MaxZ
upper z bound
Definition: Field3D.h:98
virtual double ScalarField_FourierSpace_complex(const int i, const int j, const int k) const
get the value of the scalar field, Fourier space, complex part
Definition: Field3D.h:470
long int inds_to_index(const int i, const int j, const int k) const
contract 3 indeces into one
Definition: Field3D.h:130
int m_nZ
number of cells along the z-axis
Definition: Field3D.h:62
virtual std::vector< double > ScalarField_FourierSpace_complex() const
get the value of the scalar field, Fourier space, complex part
Definition: Field3D.h:486
int m_nY
number of cells along the y-axis
Definition: Field3D.h:59
int nCellsFourier() const
get the private member m_nCells_Fourier
Definition: Field3D.h:277
double Volume() const
get the private member m_Volume
Definition: Field3D.h:337
virtual void set_ScalarField_FourierSpace_real(const double value, const int i, const int j, const int k, const bool add=0)
set the value of the scalar field in Fourier space, real part
Definition: Field3D.h:577
virtual std::vector< double > VectorField(const int i, const int j, const int k) const
get the value of the vector field
Definition: Field3D.h:436
int nz() const
get the private member m_nZ
Definition: Field3D.h:259
virtual std::vector< double > VectorField(const std::vector< double > pos) const
get the value of the vector field
Definition: Field3D.h:446
std::vector< double > m_X
coordinates of the cells along the x-axis
Definition: Field3D.h:104
virtual void set_VectorField_FourierSpace_real(const std::vector< double > value, const int i, const int j, const int k, const bool add=0)
set the value of the vector field, Fourier space, real part
Definition: Field3D.h:606
double MinZ() const
get the private member m_MinZ
Definition: Field3D.h:295
virtual double ScalarField(const std::vector< double > pos) const
get the value of the scalar field
Definition: Field3D.h:416
double m_deltaZ
Z cell size.
Definition: Field3D.h:80
std::vector< double > m_kZ
coordinates of the cells along the kz-axis
Definition: Field3D.h:119
int m_nCells_Fourier
number of cells, Fourier space
Definition: Field3D.h:71
virtual std::vector< double > VectorField_FourierSpace_real(const int i, const int j, const int k) const
get the value of the vector field, Fourier space, real part
Definition: Field3D.h:498
double YY(const int i) const
get the value of the Y coordinates at the i-th cell
Definition: Field3D.h:355
double MaxX() const
get the private member m_MaxX
Definition: Field3D.h:301
The class ScalarField3D.
Definition: Field3D.h:648
ScalarField3D()=default
default constructor
std::vector< double > ScalarField() const
get the values of the scalar field
Definition: Field3D.cpp:356
void set_ScalarField_FourierSpace_complex(const double value, const int i, const int j, const int k, const bool add=0)
set the value of the scalar field in Fourier space, real part
Definition: Field3D.cpp:307
void FourierTransformField()
perform the Fourier transform on the field
Definition: Field3D.cpp:217
void FourierAntiTransformField()
perform the anti-Fourier transform on the field
Definition: Field3D.cpp:241
double * m_field
scalar field
Definition: Field3D.h:653
void GaussianConvolutionField(const double kernel_size)
perform a smoothing of the field with a gaussian kernel
Definition: Field3D.cpp:257
void reset()
set to 0 the fields
Definition: Field3D.cpp:366
void set_ScalarField_FourierSpace_real(const double value, const int i, const int j, const int k, const bool add=0)
set the value of the scalar field in Fourier space, real part
Definition: Field3D.cpp:298
virtual ~ScalarField3D()=default
default destructor
void set_ScalarField(const double value, const int i, const int j, const int k, const bool add=0)
set the value of the scalar field
Definition: Field3D.cpp:289
fftw_complex * m_field_FourierSpace
fourier transform of the scalar field
Definition: Field3D.h:656
The class VectorField3D.
Definition: Field3D.h:871
void reset()
set to 0 the fields
Definition: Field3D.cpp:591
~VectorField3D()=default
default destructor
void FourierTransformField()
perform the anti-Fourier transform on the field
Definition: Field3D.cpp:449
void set_VectorField_FourierSpace_real(const std::vector< double > value, const int i, const int j, const int k, const bool add=0)
set the value of the vector field, Fourier space, real part
Definition: Field3D.cpp:529
void set_VectorField_FourierSpace_complex(const std::vector< double > value, const int i, const int j, const int k, const bool add=0)
set the value of the vector field, Fourier space, complex part
Definition: Field3D.cpp:540
std::vector< double * > m_field
vector field
Definition: Field3D.h:875
std::vector< fftw_complex * > m_field_FourierSpace
vector field in fourier space
Definition: Field3D.h:878
std::vector< double > VectorField_FourierSpace_complex(const int i, const int j, const int k) const
get the value of the vector field, Fourier space, complex part
Definition: Field3D.cpp:569
void set_VectorField(const std::vector< double > value, const int i, const int j, const int k, const bool add=0)
set the value of the vectorr field
Definition: Field3D.cpp:518
std::vector< double > VectorField_FourierSpace_real(const int i, const int j, const int k) const
get the value of the vector field, Fourier space, real part
Definition: Field3D.cpp:560
VectorField3D()=default
default constructor
std::vector< double > VectorField(const int i, const int j, const int k) const
get the value of the vector field
Definition: Field3D.cpp:551
void FourierAntiTransformField()
perform the anti-Fourier transform on the field
Definition: Field3D.cpp:490
The global namespace of the CosmoBolognaLib
Definition: CAMB.h:38
int ErrorCBL(const std::string msg, const std::string functionCBL, const std::string fileCBL, const cbl::glob::ExitCode exitCode=cbl::glob::ExitCode::_error_)
throw an exception: it is used for handling exceptions inside the CosmoBolognaLib
Definition: Kernel.h:780