phase_space_dg_geometry.h

Go to the documentation of this file.

#ifndef geometry_phase_space_dg_geometry_h
#define geometry_phase_space_dg_geometry_h
 
// WARPXM includes
#include "warpxm/warpxm_config.h"
#include "basis/wmbasisarrayset.h"
#include "lib/geometry/unstructured_geometry/wmudggeometry.h"
#include "lib/wmunstructuredpatch.h"
 
// STL includes
#include <vector>
#include <string>
#include <array>
 
namespace geometry
{
namespace phase_space
{
 
class PhaseSpaceDGGeometry
{
public:
    PhaseSpaceDGGeometry(const std::string basis_set_name,
                         const std::string phase_space_order,
                         const std::vector<int> vspace_nelems,
                         const std::vector<real> vspace_minimums,
                         const std::vector<real> vspace_dv);
 
    ~PhaseSpaceDGGeometry();
 
    void setupBases(const WmUnstructuredPatch& patch);
 
    void setupWmUDGGeometry(const WmUnstructuredPatch& patch,
                            const std::vector<std::string>& boundary_names,
                            const bool has_fluxes);
 
    // /**
    //  * Setup:
    //  * 1. phase space position coordinates for each phase space node in
    //  * physical_element_scope_range range given of physical elements
    //  * 2. matching global physical space node indexes for each phase space node in
    //  * physical_element_scope_range
    //  * @param physical_element_scope_range
    //  *
    //  */
    // void setupPhaseSpaceDGGeometry(const WxRange& physical_element_scope_range);
 
    int getFlattenedVelocitySpaceElementIndex(const int vx_element_index,
                                              const int vy_element_index,
                                              const int vz_element_index) const;
 
    void getVelocityElementCornerCoordinates(
        const int velocity_space_element_indexes[3],
        real velocity_element_corner_coordinates[3]) const;
 
    void calcPhaseSpaceInverseJacobian(
        const int physical_element_index,
        std::vector<std::vector<real>>& inverse_jacobian_phase_space) const;
 
    // /**
    //  * Calculate all phase space position coordinates for a physical element scope range
    //  * and write into _phase_space_coordinates
    //  * @param physical_element_scope_range
    //  *
    //  */
    // void setAllPhaseSpacePositionCoordinates(const WxRange& physical_element_scope_range);
 
    // /**
    //  * Calculate the phase space position coordinate at a particular phase space node
    //  * @param phase_space_element_local_node_index phase space local node index where we
    //  * are calculating the phase space position
    //  * @param physical_element_index physical element index
    //  * @param phase_space_element_base_velocity_coordinates phase space element corder
    //  * coordinates
    //  * @return position phase space position at this phase_space_element_local_node_index
    //  * {x, y, z, vx, vy, vz}
    //  *
    //  */
    // void setPhaseSpaceCoords(const int phase_space_element_local_node_index,
    //                          const int physical_element_index,
    //                          const real phase_space_element_base_velocity_coordinates[3],
    //                          real position[6]);
 
    // /**
    //  * Calculate the matching physical space node index for each phase space node index,
    //  * writing into _matching_global_physical_space_indexes
    //  * @param physical_element_scope_range
    //  *
    //  */
    // void
    // setMatchingGlobalPhysicalSpaceIndexes(const WxRange& physical_element_scope_range);
 
    WmUDGGeometry const* getWmUDGGeometry() const
    {
        return _physical_space_dg_geometry.get();
    }
 
    WmBasisArraySet const* getPhysicalSpaceBasisArraySet() const
    {
        return _physical_space_basis_set.get();
    }
 
    WmBasisArraySet const* getPhaseSpaceBasisArraySet() const
    {
        return _phase_space_basis_set.get();
    }
 
    int getVspaceNdims() const
    {
        return _vspace_ndims;
    }
 
    const std::vector<int>& getVspaceNelems() const
    {
        return _vspace_nelems;
    }
 
    const std::vector<real>& getVspaceDv() const
    {
        return _vspace_dv;
    }
 
    int getVspaceNVxElems() const
    {
        return _vspace_n_vx_elems;
    }
 
    int getVspaceNVyElems() const
    {
        return _vspace_n_vy_elems;
    }
 
    int getVspaceNVzElems() const
    {
        return _vspace_n_vz_elems;
    }
 
    const std::vector<real>& getVspaceMinimums() const
    {
        return _vspace_minimums;
    }
 
    const std::vector<real>& getVspaceMaximums() const
    {
        return _vspace_maximums;
    }
 
    int getMatchingGlobalPhysicalSpaceNodeIndex(
        const int physical_element_index,
        const int phase_space_element_local_node_index) const
    {
        // Phase space local nodes should line up with physical space such that the
        // modulus operator should just map the phase space node to the physical space
        // node, and then add the total nodes up to the physical element to get the global
        // node index of the physical space node.
        const int N_p = _physical_space_basis_set->getNumTotalNodes();
        return physical_element_index * N_p + phase_space_element_local_node_index % N_p;
    }
 
    // const std::vector<real>&
    // getPhaseSpaceCoordinatesOld(const int physical_element_index,
    //                          const int flattened_velocity_space_element_index,
    //                          const int phase_space_element_local_node_index) const
    // {
    //     return _phase_space_coordinates[physical_element_index]
    //                                    [flattened_velocity_space_element_index]
    //                                    [phase_space_element_local_node_index];
    // }
 
    void getPhaseSpaceCoordinates(const int physical_element_index,
                                  const int velocity_space_element_indexes[3],
                                  const int phase_space_element_local_node_index,
                                  std::array<real, 6>& position) const;
 
    int getGlobalPhaseSpaceElementIndex(const int physical_element_index,
                                        const int vx_element_index,
                                        const int vy_element_index,
                                        const int vz_element_index);
 
    int getGlobalPhaseSpaceNodeIndex(const int physical_element_index,
                                     const int vx_element_index,
                                     const int vy_element_index,
                                     const int vz_element_index,
                                     const int local_phase_space_element_node_index);
 
    const std::array<int, 18>& getVfaceLookupTable()
    {
        return _vface_lookup;
    }
 
protected:
private:
    // To hold phase space information, we wan to hold both:
    // -> physical space information
    // -> velocity space information
 
    // physical space information:
 
    // Spatial basis set name
    std::string _basis_set_name;
 
    // Spatial basis set
    std::unique_ptr<WmBasisArraySet> _physical_space_basis_set;
 
    // WmUDGGeometry
    std::unique_ptr<WmUDGGeometry> _physical_space_dg_geometry;
 
    // velocity space information
 
    std::vector<int> _vspace_nelems;
    int _vspace_ndims;
 
    std::vector<real> _vspace_minimums;
    std::vector<real> _vspace_maximums;
 
    std::vector<real> _vspace_dv;
 
    // number of velocity elements (will be minimum 1)
    int _vspace_n_vx_elems;
    int _vspace_n_vy_elems;
    int _vspace_n_vz_elems;
 
    // Then we can hold the phase space
    // Phase space basis set
    std::string _phase_space_order;
    std::unique_ptr<WmBasisArraySet> _phase_space_basis_set;
 
    // lookup table for neighboring velocity face phase space elements
    // there should be up to 6 velocity faces and vx, vy, or vz neighbors are +-1,
    // depending on the face
    // faces 0, 1 are -+ vx respectively
    // faces 2, 3 are -+ vy respectively
    // faces 4, 5 are -+ vz respectively
    // The lookup table should thus be a 6x3 array looking like:
    // [[-1, 0, 0],
    // [+1, 0, 0],
    // [0, -1, 0],
    // [0, +1, 0],
    // [0, 0, -1],
    // [0, 0, +1]]
    // but is flattened c-style row-wise
    std::array<int, 18> _vface_lookup;
};
 
} // namespace phase_space
} // namespace geometry
 
#endif