BinnedFieldSolver.h

Go to the documentation of this file.

#ifndef OPALX_BINNED_FIELD_SOLVER_H
#define OPALX_BINNED_FIELD_SOLVER_H
 
#include <algorithm>
#include <cmath>
#include <functional>
#include <iomanip>
#include <memory>
#include <numeric>
#include <string>
#include <vector>
 
#include "PartBunch/FieldSolver.hpp"
#include "PartBunch/ImageChargeScatterController.h"
#include "PartBunch/PartBunch.h"
#include "Utilities/OpalException.h"
 
template <typename T, unsigned Dim>
class BinnedFieldSolver : public FieldSolver<T, Dim> {
    static_assert(Dim == 3, "BinnedFieldSolver currently supports Dim == 3 only.");
 
public:
    using PartBunch_t    = PartBunch<T, Dim>;
    using ParticleCtr_t  = typename PartBunch_t::ParticleContainer_t;
    using AdaptBins_t    = typename PartBunch_t::AdaptBins_t;
    using BCHandler_t    = BCHandler<Dim>;
    using bin_index_type = typename AdaptBins_t::bin_index_type;
    using size_type      = typename AdaptBins_t::size_type;
 
    using particle_position_type = typename PartBunch_t::Base::particle_position_type;
 
    enum class ScatterAttribute { ChargeQ };
 
    enum class ImageScatterMode { PrimaryAndImage, PrimaryOnly, ImageOnly };
 
    enum class GatherAttribute { ElectricFieldE };
 
    BinnedFieldSolver(
            std::string solver, Field_t<Dim>* rho, VField_t<T, Dim>* E, Field_t<Dim>* phi,
            std::shared_ptr<BCHandler_t> bcHandler, int tablePrintFrequency, bool adaptiveBinning);
 
    void computeSelfFields(PartBunch_t& bunch);
 
    void setScatterAttribute(const ScatterAttribute attr);
 
    void setGatherAttribute(const GatherAttribute attr);
 
    void setImageChargeConfiguration(bool enabled, double zPlane);
    bool isImageChargeEnabled() const { return imageScatterController_m.isEnabled(); }
    double getImageChargePlaneZ() const { return imageScatterController_m.getZPlane(); }
 
    void setShiftedGreensConfiguration(bool enabled, double zPlane);
    bool isShiftedGreensEnabled() const { return shiftedGreensEnabled_m; }
    double getShiftedGreensPlaneZ() const { return shiftedGreensPlaneZ_m; }
 
    void setZerofaceMaxSteps(int maxSteps);
    int getZerofaceMaxSteps() const { return zerofaceMaxSteps_m; }
 
    bool isImageChargeActiveForStep(size_t step) const;
 
    bool isShiftedGreensActiveForStep(size_t step) const;
 
    void setZeroFacePlaneDumpFrequency(int frequency);
    int getZeroFacePlaneDumpFrequency() const { return zeroFacePlaneDumpFrequency_m; }
 
    struct BinKinematics {
        Vector_t<double, Dim> pmean = Vector_t<double, Dim>(0.0);
        double gammaBin             = 1.0;
    };
 
private:
    ScatterAttribute scatterAttribute_m;
    GatherAttribute gatherAttribute_m;
    int tablePrintFrequency_m        = 0;
    bool adaptiveBinning_m           = true;
    int zeroFacePlaneDumpFrequency_m = 0;
    int zerofaceMaxSteps_m           = 0;
    ImageChargeScatterController<T, Dim> imageScatterController_m;
    bool warnedPlaneDumpParallelUnsupported_m = false;
 
    // Shifted Green's function Dirichlet correction (alternative to image charges).
    // Mutually exclusive with the image-charge path (enforced at config time).
    bool shiftedGreensEnabled_m  = false;
    double shiftedGreensPlaneZ_m = 0.0;
 
    // Scratch field holding the axis-flipped version of E' for the shifted-GF
    // correction pass. Populated by buildFlippedZSlab (which delegates to
    // opalx::detail::mirrorField) before accumulateFieldToTemp's flipped branch.
    // Same layout/mesh/ghost count as *(this->getE()). Allocated lazily on first
    // use; reused across bins and timesteps.
    std::shared_ptr<VField_t<T, Dim>> flippedZSlabField_m;
 
    struct BinStatsRow {
        long long binNumber;            
        unsigned long long nParticles;  
        double gammaBin;                
    };
 
    void printBinStatsTable(
            const std::string& binningCmdName, const std::vector<BinStatsRow>& rows);
 
public:
    static void setScalarField(Field_t<Dim>& field, double value);
 
    static void scaleAndShiftScalarField(Field_t<Dim>& field, double scale, double shift);
 
    static void setVectorField(VField_t<T, Dim>& field, const Vector_t<T, Dim>& value);
 
private:
    void dumpDirichletPlaneDiagnosticsIfRequested(PartBunch_t& bunch, const std::string& solveTag);
 
    void computeBinnedSelfFields(PartBunch_t& bunch);
 
    void computeLegacySelfFields(PartBunch_t& bunch);
 
    void rebinAndPrepare(PartBunch_t& bunch, std::shared_ptr<AdaptBins_t> bins);
 
public:
    BinKinematics computeGammaBinGlobal(
            PartBunch_t& bunch, std::shared_ptr<AdaptBins_t> bins, const bin_index_type binIndex,
            const size_type nPartGlobal) const;
 
    void prepareRhoForBin(
            PartBunch_t& bunch, std::shared_ptr<AdaptBins_t> bins, const bin_index_type binIndex,
            const size_type nPartGlobal, const double gammaBin,
            ImageScatterMode mode = ImageScatterMode::PrimaryAndImage);
 
    void accumulateFieldToTemp(
            const double gammaBin, const Vector_t<double, Dim>& pmean,
            std::shared_ptr<VField_t<T, Dim>> EtmpSP, std::shared_ptr<VField_t<T, Dim>> BtmpSP,
            double bFieldSign = 1.0, int flipAxis = -1);
 
private:
    void buildFlippedZSlab(const VField_t<T, Dim>& src);
 
    void gatherFromTempToParticles(
            PartBunch_t& bunch, std::shared_ptr<VField_t<T, Dim>> EtmpSP,
            std::shared_ptr<VField_t<T, Dim>> BtmpSP);
};
 
// Reduce compile-time churn: instantiate the only supported concrete solver in one TU.
extern template class BinnedFieldSolver<double, 3>;
 
#endif  // OPALX_BINNED_FIELD_SOLVER_H