FFTBoxPoissonSolver.h 4.5 KB
Newer Older
gsell's avatar
gsell committed
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
// -*- C++ -*-
/***************************************************************************
 *
 *
 * FFTBoxPoissonSolver.hh
 *
 * Open BC in x,y and z.
 *
 *
 *
 *
 *
 *
 ***************************************************************************/

////////////////////////////////////////////////////////////////////////////
// This class contains methods for solving Poisson's equation for the
// space charge portion of the calculation.
////////////////////////////////////////////////////////////////////////////

#ifndef FFT_BOXPOISSON_SOLVER_H_
#define FFT_BOXPOISSON_SOLVER_H_

//////////////////////////////////////////////////////////////
#include "PoissonSolver.h"
26
class PartBunch;
gsell's avatar
gsell committed
27 28 29 30 31
//////////////////////////////////////////////////////////////

class FFTBoxPoissonSolver : public PoissonSolver {
public:
    // constructor and destructor
32
    FFTBoxPoissonSolver(PartBunch &bunch, std::string greensFuntion);
gsell's avatar
gsell committed
33

34
    FFTBoxPoissonSolver(Mesh_t *mesh, FieldLayout_t *fl, std::string greensFunction, double boxSize);
gsell's avatar
gsell committed
35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

    ~FFTBoxPoissonSolver();

    // given a charge-density field rho and a set of mesh spacings hr,
    // compute the scalar potential with image charges at  -z
    void computePotential(Field_t &rho, Vector_t hr, double zshift);

    // given a charge-density field rho and a set of mesh spacings hr,
    // compute the scalar potential in open space
    void computePotential(Field_t &rho, Vector_t hr);

    // compute the green's function for a Poisson problem and put it in in grntm_m
    // uses grnIField_m to eliminate excess calculation in greenFunction()
    // given mesh information in nr and hr
    void greensFunction();

    /// compute the integrated Green function as described in <A HREF="http://prst-ab.aps.org/abstract/PRSTAB/v9/i4/e044204">Three-dimensional quasistatic model for high brightness beam dynamics simulation</A> by Qiang et al.
    void integratedGreensFunction();

    /// compute the shifted integrated Green function as described in <A HREF="http://prst-ab.aps.org/abstract/PRSTAB/v9/i4/e044204">Three-dimensional quasistatic model for high brightness beam dynamics simulation</A> by Qiang et al.
    void shiftedIntGreensFunction(double zshift);

    double getXRangeMin() {return -a_m;}
    double getXRangeMax() {return  a_m;}
    double getYRangeMin() {return -a_m;}
    double getYRangeMax() {return  a_m;}
61 62 63 64
    double getZRangeMin() {return -a_m; }
    double getZRangeMax() {return  a_m; }


gsell's avatar
gsell committed
65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136

    Inform &print(Inform &os) const;

private:

    // rho2_m is the charge-density field with mesh doubled in each dimension
    Field_t rho2_m;

    // grntr_m is the Fourier transformed Green's function
    Field_t grntr_m;

    // Fields used to eliminate excess calculation in greensFunction()
    // mesh2_m and layout2_m are used
    IField_t grnIField_m[3];

    SINE_t *sine_m;

    // mesh and layout objects for rho_m
    Mesh_t *mesh_m;
    FieldLayout_t *layout_m;

    // mesh and layout objects for rho2_m
    Mesh_t *mesh2_m;
    FieldLayout_t *layout2_m;

    // tmp
    Field_t tmpgreen;

    // domains for the various fields
    NDIndex<3> domain_m;             // original domain, gridsize
    // mesh and gridsize defined outside of FFT class, given as
    // parameter to the constructor (mesh and layout object).
    NDIndex<3> domain2_m;            // doubled gridsize (2*Nx,2*Ny,2*Nz)
    // (2*Nx,Ny,2*Nz)
    // mesh spacing and size values
    Vector_t hr_m;
    Vektor<int, 3> nr_m;

    std::string greensFunction_m;

    double a_m; // the box size

    IpplTimings::TimerRef GreensFunctionTimer_m;

    IpplTimings::TimerRef IntGreensFunctionTimer1_m;
    IpplTimings::TimerRef IntGreensFunctionTimer2_m;
    IpplTimings::TimerRef IntGreensFunctionTimer3_m;
    IpplTimings::TimerRef IntGreensFunctionTimer4_m;

    IpplTimings::TimerRef ShIntGreensFunctionTimer1_m;
    IpplTimings::TimerRef ShIntGreensFunctionTimer2_m;
    IpplTimings::TimerRef ShIntGreensFunctionTimer3_m;
    IpplTimings::TimerRef ShIntGreensFunctionTimer4_m;

    IpplTimings::TimerRef GreensFunctionTimer1_m;
    IpplTimings::TimerRef GreensFunctionTimer2_m;
    IpplTimings::TimerRef GreensFunctionTimer3_m;
    IpplTimings::TimerRef GreensFunctionTimer4_m;
};

inline Inform &operator<<(Inform &os, const FFTBoxPoissonSolver &fs) {
    return fs.print(os);
}



#endif

/***************************************************************************
 * $RCSfile: FFTBoxPoissonSolver.hh,v $   $Author: adelmann $
 * $Revision: 1.1.1.1 $   $Date: 2001/08/08 11:21:48 $
 ***************************************************************************/