grid.hpp

Go to the documentation of this file.

#ifndef GRID_HPP
#define GRID_HPP

#include "magnetic_field.hpp"
#include "grid_structs.hpp"
#include "grid_files.hpp"

// Grid class 
template<class Device>
class Grid {
    public:

    // Constructor from inputs sent through fortran interface
    Grid(int ntriangle, int nnode, int nplanes, double delta_phi, int nguess, int guess_list_len, double guess_min1, double guess_min2, double guess_max1,
         double guess_max2, double inv_guess_d1, double inv_guess_d2, int npsi_surf2, int iphi_offset, int npsi00, double psi00min, double psi00max,
         double dpsi00, int nwall, double wedge_angle, double eq_x_psi,
         double* psi, int* basis, Vertex* nodes, VertexMap* mapping, int* guess_list, int* guess_xtable, int* guess_count,
         int* guess_list_1d, double* psi_surf2,double minval_psi_surf2, double maxval_psi_surf2, 
         int* wall_nodes, int* node_to_wall, int* rgn, RZPair* gx, int nrho, double rhomax, double phimax);

    // Constructor from file inputs
    Grid(GridFiles& grid_files, const MagneticField<Device>& magnetic_field, int nplanes_in, int wedge_n_in, int iphi_offset_in, int guess_table_size, int npsi00_in, int nrho_in, double rhomax_in);

    // Default constructor
    Grid(){}

    KOKKOS_INLINE_FUNCTION void t_coeff(const SimdVector2D &x, const Simd<int>& itr, SimdGridVec &p) const;

    KOKKOS_INLINE_FUNCTION void t_coeff_mod(const MagneticField<Device> &magnetic_field, const SimdVector2D &xy, const Simd<double>& psiin, const Simd<int>& itr, SimdGridVec &p ) const;

    KOKKOS_INLINE_FUNCTION void search_tr2( const SimdVector2D &xy, Simd<int>& itr, SimdGridVec &pout ) const;

    KOKKOS_INLINE_FUNCTION void search_tr_check_guess(const SimdVector2D &x, const Simd<int>& old_itr,Simd<int>& itr,SimdGridVec &p ) const;

    KOKKOS_INLINE_FUNCTION void search_ptl_1d(double psi,double &wp,int &ip) const;

    KOKKOS_INLINE_FUNCTION void charge_search_index(const MagneticField<Device> &magnetic_field, const SimdVector2D &x, const Simd<double>& phi, const Simd<double>& psi, SimdVector2D &xff,Simd<int>& itr, SimdGridVec &p) const;

    KOKKOS_INLINE_FUNCTION void get_wall_index(const Simd<bool>& just_left_the_grid,const SimdVector2D &x,const Simd<int>& itr, const SimdGridVec &p, Simd<int>& widx) const;

    KOKKOS_INLINE_FUNCTION void wedge_modulo_phi(Simd<double>& phi_mod) const;

    KOKKOS_INLINE_FUNCTION void check_triangle(const Simd<int>& itr,Simd<bool>& wasnt_in_triangle) const;

    KOKKOS_INLINE_FUNCTION bool node_is_inside_psi_range(const MagneticField<Device> &magnetic_field, const int node) const;
    
    void draw_ascii_grid(MagneticField<Device> magnetic_field, double rmin, double rmax, double dr, double zmin, double zmax, double dz){
        printf("\n\n******** PLOT OF GRID ******* \n");
        Kokkos::View<int*,HostType> save_vert("save_vert", 1+((rmax-rmin)/dr));
        for (double zplot = zmax; zplot>=zmin; zplot-=dz){
            printf("\n");
            int save_num = 0;
            int rsave_ind = -1;
            for (double rplot = rmin; rplot<=rmax; rplot+=dr){
                rsave_ind++;
                SimdVector2D x;
                x.r = rplot;
                x.z = zplot;
                Simd<int> itr = -1;

                Simd<double> phi = 0.001;
                Simd<double> psi;
                //magnetic_field.psi_bicub.interpol(x.r,x.z,0,0,psi);

                // Get itr (triangle index)
                SimdVector2D xff; // Unused
                SimdGridVec p; // Unused
                charge_search_index(magnetic_field, x, phi, psi,xff,itr, p);
                if(itr[0]!=save_num || itr[0]!=save_vert[rsave_ind]){
                    if(itr[0]<10 && itr[0]>=0) printf("%d  ",itr[0]);
                    else printf("%d ",itr[0]);
                    save_num = itr[0];
                    save_vert[rsave_ind] = itr[0];
                } else {printf("   ");}
            }
        }
    }

    int ntriangle;         
    int nnode;             
    int nplanes;           
    double wedge_angle;    
    double delta_phi;      
    double inv_delta_phi;  

    Kokkos::View<double*,Kokkos::LayoutRight,Device> psi;   
    Kokkos::View<int*,   Kokkos::LayoutRight,Device> basis; 
    Kokkos::View<Vertex*,Kokkos::LayoutRight,Device> nodes;
    Kokkos::View<VertexMap*,Kokkos::LayoutRight,Device> mapping;

    int nguess;
    int guess_list_len;
    double guess_min1;
    double guess_min2;
    double guess_max1;
    double guess_max2;
    double inv_guess_d1;
    double inv_guess_d2;
    Kokkos::View<int*,Kokkos::LayoutRight,Device> guess_list;
    Kokkos::View<int**,Kokkos::LayoutRight,Device> guess_xtable;
    Kokkos::View<int**,Kokkos::LayoutRight,Device> guess_count;
    Kokkos::View<int**,Kokkos::LayoutRight,Device> guess_list_1d;

    int npsi_surf2;
    double minval_psi_surf2;
    double maxval_psi_surf2;
    Kokkos::View<double*,Kokkos::LayoutRight,Device> psi_surf2;

    int iphi_offset;
    int npsi00;
    double psi00min;
    double psi00max;
    double dpsi00;

    Kokkos::View<int*,Kokkos::LayoutRight,Device> rgn;
    Kokkos::View<RZPair*,Kokkos::LayoutRight,Device> gx;
    
    int nwall;
    Kokkos::View<int*,Kokkos::LayoutRight,Device> wall_nodes;
    Kokkos::View<int*,Kokkos::LayoutRight,Device> node_to_wall;

    double eps_flux_surface; 

    // Gyroaverage grid information
    int nrho;
    double rhomax;
    double drho;
    double phimin;
    double phimax;
};

#include "grid.tpp"
#endif