WarpX
IntegratedGreenFunctionSolver.H
Go to the documentation of this file.
1 /* Copyright 2019-2024 Arianna Formenti, Remi Lehe
2  *
3  * This file is part of ABLASTR.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 #ifndef ABLASTR_IGF_SOLVER_H
8 #define ABLASTR_IGF_SOLVER_H
9 
10 #include <AMReX_BoxArray.H>
11 #include <AMReX_GpuQualifiers.H>
12 #include <AMReX_MultiFab.H>
13 #include <AMReX_REAL.H>
14 #include <AMReX_Vector.H>
15 
16 
17 #include <array>
18 #include <cmath>
19 
20 
21 namespace ablastr::fields
22 {
23 
34  amrex::Real
35  IntegratedPotential (amrex::Real x, amrex::Real y, amrex::Real z)
36  {
37  using namespace amrex::literals;
38 
39  amrex::Real const r = std::sqrt( x*x + y*y + z*z );
40  amrex::Real const G =
41  - 0.5_rt * z*z * std::atan( x*y/(z*r) )
42  - 0.5_rt * y*y * std::atan( x*z/(y*r) )
43  - 0.5_rt * x*x * std::atan( y*z/(x*r) )
44  + y*z*std::asinh( x/std::sqrt(y*y + z*z) )
45  + x*z*std::asinh( y/std::sqrt(x*x + z*z) )
46  + x*y*std::asinh( z/std::sqrt(x*x + y*y) );
47  return G;
48  }
49 
58  void
59  computePhiIGF (amrex::MultiFab const & rho,
60  amrex::MultiFab & phi,
61  std::array<amrex::Real, 3> const & cell_size,
62  amrex::BoxArray const & ba);
63 
64 } // namespace ablastr::fields
65 
66 #endif // ABLASTR_IGF_SOLVER_H
#define AMREX_FORCE_INLINE
#define AMREX_GPU_HOST_DEVICE
Definition: IntegratedGreenFunctionSolver.cpp:33
void computePhiIGF(amrex::MultiFab const &rho, amrex::MultiFab &phi, std::array< amrex::Real, 3 > const &cell_size, amrex::BoxArray const &ba)
Compute the electrostatic potential using the Integrated Green Function method as in http://dx....
Definition: IntegratedGreenFunctionSolver.cpp:36
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real IntegratedPotential(amrex::Real x, amrex::Real y, amrex::Real z)
Implements equation 2 in https://doi.org/10.1103/PhysRevSTAB.10.129901 with some modification to symm...
Definition: IntegratedGreenFunctionSolver.H:35
cell_size
Definition: compute_domain.py:37