Line data Source code
1 : //
2 : // Class FFTTruncatedGreenPeriodicPoissonSolver
3 : // Poisson solver for periodic boundaries, based on FFTs.
4 : // Solves laplace(phi) = -rho, and E = -grad(phi).
5 : //
6 : // Uses a convolution with a Green's function given by:
7 : // G(r) = forceConstant * erf(alpha * r) / r,
8 : // alpha = controls long-range interaction.
9 : //
10 : //
11 :
12 : #ifndef IPPL_FFT_TRUNCATED_GREEN_PERIODIC_POISSON_SOLVER_H_SOLVER_H_
13 : #define IPPL_FFT_TRUNCATED_GREEN_PERIODIC_POISSON_SOLVER_H_SOLVER_H_
14 :
15 : #include "Types/Vector.h"
16 :
17 : #include "Field/Field.h"
18 :
19 : #include "FFT/FFT.h"
20 : #include "FieldLayout/FieldLayout.h"
21 : #include "Meshes/UniformCartesian.h"
22 : #include "Poisson.h"
23 :
24 : namespace ippl {
25 : template <typename FieldLHS, typename FieldRHS>
26 : class FFTTruncatedGreenPeriodicPoissonSolver : public Poisson<FieldLHS, FieldRHS> {
27 : constexpr static unsigned Dim = FieldLHS::dim;
28 : using Trhs = typename FieldRHS::value_type;
29 : using mesh_type = typename FieldRHS::Mesh_t;
30 :
31 : public:
32 : // type of output
33 : using Base = Poisson<FieldLHS, FieldRHS>;
34 :
35 : // types for LHS and RHS
36 : using typename Base::lhs_type, typename Base::rhs_type;
37 :
38 : // define a type for the 3 dimensional real to complex Fourier transform
39 : typedef FFT<RCTransform, FieldRHS> FFT_t;
40 :
41 : // define a type for a 3 dimensional field (e.g. charge density field)
42 : // define a type of Field with integers to be used for the helper Green's function
43 : // also define a type for the Fourier transformed complex valued fields
44 : typedef FieldRHS Field_t;
45 : typedef Field<int, Dim, mesh_type, typename FieldLHS::Centering_t> IField_t;
46 : typedef typename FFT_t::ComplexField CxField_t;
47 : typedef Vector<Trhs, Dim> Vector_t;
48 :
49 : // define type for field layout
50 : typedef FieldLayout<Dim> FieldLayout_t;
51 :
52 : // constructor and destructor
53 : FFTTruncatedGreenPeriodicPoissonSolver();
54 : FFTTruncatedGreenPeriodicPoissonSolver(rhs_type& rhs, ParameterList& params);
55 : FFTTruncatedGreenPeriodicPoissonSolver(lhs_type& lhs, rhs_type& rhs, ParameterList& params);
56 0 : ~FFTTruncatedGreenPeriodicPoissonSolver() = default;
57 :
58 : // override the setRhs function of the Solver class
59 : // since we need to call initializeFields()
60 : void setRhs(rhs_type& rhs) override;
61 :
62 : // solve the Poisson equation
63 : // more specifically, compute the scalar potential given a density field rho
64 : void solve() override;
65 :
66 : // function called in the constructor to initialize the fields
67 : void initializeFields();
68 :
69 : // compute standard Green's function
70 : void greensFunction();
71 :
72 :
73 : private:
74 : Field_t grn_m; // the Green's function
75 :
76 : CxField_t rhotr_m;
77 : CxField_t grntr_m;
78 : CxField_t tempFieldComplex_m;
79 :
80 : // fields that facilitate the calculation in greensFunction()
81 : IField_t grnIField_m[Dim];
82 :
83 : // the FFT object
84 : std::unique_ptr<FFT_t> fft_m;
85 :
86 : // mesh and layout objects for rho_m (RHS)
87 : mesh_type* mesh_mp;
88 : FieldLayout_t* layout_mp;
89 :
90 : // mesh and layout objects for the Fourier transformed Complex fields
91 : std::unique_ptr<mesh_type> meshComplex_m;
92 : std::unique_ptr<FieldLayout_t> layoutComplex_m;
93 :
94 : // domains for the various fields
95 : NDIndex<Dim> domain_m; // physical domain
96 : NDIndex<Dim> domainComplex_m; // Fourier domain
97 :
98 : // mesh spacing and mesh size
99 : Vector_t hr_m;
100 : Vector<int, Dim> nr_m;
101 :
102 : protected:
103 0 : void setDefaultParameters() override {
104 : using heffteBackend = typename FFT_t::heffteBackend;
105 0 : heffte::plan_options opts = heffte::default_options<heffteBackend>();
106 0 : this->params_m.add("use_pencils", opts.use_pencils);
107 0 : this->params_m.add("use_reorder", opts.use_reorder);
108 0 : this->params_m.add("use_gpu_aware", opts.use_gpu_aware);
109 0 : this->params_m.add("r2c_direction", 0);
110 0 : this->params_m.template add<Trhs>("alpha", 1);
111 0 : this->params_m.template add<Trhs>("force_constant", 1);
112 :
113 0 : switch (opts.algorithm) {
114 0 : case heffte::reshape_algorithm::alltoall:
115 0 : this->params_m.add("comm", a2a);
116 0 : break;
117 0 : case heffte::reshape_algorithm::alltoallv:
118 0 : this->params_m.add("comm", a2av);
119 0 : break;
120 0 : case heffte::reshape_algorithm::p2p:
121 0 : this->params_m.add("comm", p2p);
122 0 : break;
123 0 : case heffte::reshape_algorithm::p2p_plined:
124 0 : this->params_m.add("comm", p2p_pl);
125 0 : break;
126 0 : default:
127 0 : throw IpplException("FFTTruncatedGreenPeriodicPoissonSolver::setDefaultParameters",
128 : "Unrecognized heffte communication type");
129 : }
130 0 : }
131 : };
132 : } // namespace ippl
133 :
134 : #include "PoissonSolvers/FFTTruncatedGreenPeriodicPoissonSolver.hpp"
135 : #endif // IPPL_FFT_TRUNCATED_GREEN_PERIODIC_POISSON_SOLVER_H_SOLVER_H_
|
