XGC1
grid_files.hpp
Go to the documentation of this file.
1 #ifndef GRID_FILES_HPP
2 #define GRID_FILES_HPP
3 
4 #ifdef USE_MPI
5 #include "my_mpi.hpp"
6 #endif
7 #include "file_reader.hpp"
8 #include "magnetic_field.hpp"
9 #include "grid_structs.hpp"
10 #include "space_settings.hpp"
11 #include "NamelistReader.hpp"
12 
14  public:
15  virtual ~PlaneFilesInterface() = default;
16 
17  // Functions to query data from PlaneFiles
18  [[nodiscard]] virtual int get_nnodes() const = 0;
19  [[nodiscard]] virtual View<RZPair*, HostType> get_x() const = 0;
20  [[nodiscard]] virtual View<int*, HostType> get_rgn() const = 0;
21 
22  [[nodiscard]] virtual int get_ntriangles() const = 0;
23  [[nodiscard]] virtual View<Vertex*, HostType> get_nodes() const = 0;
24 
25  [[nodiscard]] virtual View<int*, HostType> get_i_x() const = 0;
26 
27  [[nodiscard]] virtual int get_nsurfaces() const = 0;
28  [[nodiscard]] virtual int get_nsurfaces1() const = 0;
29  [[nodiscard]] virtual int get_nsurfaces2() const = 0;
30  [[nodiscard]] virtual int get_nsurfaces3a() const = 0;
31  [[nodiscard]] virtual int get_nsurfaces3b() const = 0;
32  [[nodiscard]] virtual View<int*, HostType> get_i_surf_separatrices() const = 0;
33  [[nodiscard]] virtual View<int*, HostType> get_nnodes_on_surface() const = 0;
34  [[nodiscard]] virtual int get_sum_nnodes_on_surfaces() const = 0;
35  [[nodiscard]] virtual int get_num_non_aligned ()const = 0;
36  [[nodiscard]] virtual View<int*, HostType> get_non_aligned_vert() const = 0;
37  [[nodiscard]] virtual View<int*, HostType> get_non_aligned_nsurf() const = 0;
38  [[nodiscard]] virtual View<int**, CLayout, HostType> get_non_aligned_surf_idx() const = 0;
39  [[nodiscard]] virtual View<int**, CLayout, HostType> get_surf_idx() const = 0;
40 
41  [[nodiscard]] virtual View<int*, HostType> pack_scalars(int max_nodes_on_surface, int n_xpts, int n_separatrices) const = 0;
42  virtual void unpack_scalars(const View<int*, HostType>& scalars, int& max_nodes_on_surface, int& n_xpts, int& n_separatrices) = 0;
43 
44  [[nodiscard]] virtual int n_specified_xpts() const = 0;
45  [[nodiscard]] virtual View<Equil::XPoint*, HostType> get_xpts() const = 0;
46  [[nodiscard]] virtual RZPair get_axis() const = 0;
47  [[nodiscard]] virtual RZPair get_reference_point() const = 0;
48  [[nodiscard]] virtual View<double*, HostType> get_psiN() const = 0;
49 };
50 
52 
53  protected:
54  // These are the regions designated in the .node file
55  enum FileRegion{
56  Inside=0,
57  OnWall
58  };
59 
60  // Contents of .node
61  int nnodes;
62  View<RZPair*, HostType> x;
63  View<int*, HostType> rgn;
64 
65  // Contents of .ele
67  View<Vertex*, HostType> nodes;
68 
69  // Contents of .flx.aif
70  View<int*, HostType> i_x;
71 
72  int nsurfaces;
74  View<int*, HostType> i_surf_separatrices;
75  View<int*, HostType> nnodes_on_surface;
78  View<int*, HostType> non_aligned_vert;
79  View<int*, HostType> non_aligned_nsurf;
80  View<int**, CLayout, HostType> non_aligned_surf_idx;
81  View<int**, CLayout, HostType> surf_idx;
82 
83  View<double*, HostType> psiN;
84 };
85 
86 struct PlaneFiles: public PlaneFilesBase{
87 
88  void convert_to_x_rgn_format(const View<double**, CLayout, HostType>& tmp_nodes, View<RZPair*, HostType>& x, View<int*, HostType>& rgn, bool is_stellarator);
89 
90  void convert_to_node_format(const View<int**, CLayout, HostType>& tmp_ele, View<Vertex*, HostType>& nodes);
91 
92  // Default constructor, to use to deallocate the file arrays after use
94 
95 // Should remove this USE_MPI at some point
96 #ifdef USE_MPI
97 
98  // Copy data from a PlaneFiles object on a neighboring rank
99  PlaneFiles(const MyMPI& my_mpi, const PlaneFilesInterface& plane_files);
100 
101  // Read the files
102  PlaneFiles(const std::string& node_file, const std::string& element_file, const std::string& flx_aif_file, bool is_stellarator, const MPI_Comm& comm);
103 #endif
104 
105  // Options: hex vs circular
108  Circular=1
109  };
110 
111  void construct_hex_grid(int nshells, double raxis, double zaxis, double rscale, double zscale);
112 
113  double node_to_node_dist(RZPair a, RZPair b);
114 
115  void construct_circular_grid(int nshells, double raxis, double zaxis, double rscale, double zscale);
116 
117  static int nshells_from_nnodes(int nnodes);
118 
119  PlaneFiles(int nshells, double raxis, double zaxis, double rscale, double zscale, GridCreateOpt grid_opt=Hexagonal);
120 
122 
123  // Functions
124  int get_nnodes() const override;
125  View<RZPair*, HostType> get_x() const override;
126  View<int*, HostType> get_rgn() const override;
127 
128  int get_ntriangles() const override;
129  View<Vertex*, HostType> get_nodes() const override;
130 
131  View<int*, HostType> get_i_x() const override;
132 
133  int get_nsurfaces() const override;
134  int get_nsurfaces1() const override;
135  int get_nsurfaces2() const override;
136  int get_nsurfaces3a() const override;
137  int get_nsurfaces3b() const override;
138  View<int*, HostType> get_i_surf_separatrices() const override;
139  View<int*, HostType> get_nnodes_on_surface() const override;
140  int get_sum_nnodes_on_surfaces() const override;
141  int get_num_non_aligned() const override;
142  View<int*, HostType> get_non_aligned_vert() const override;
143  View<int*, HostType> get_non_aligned_nsurf() const override;
144  View<int**, CLayout, HostType> get_non_aligned_surf_idx() const override;
145  View<int**, CLayout, HostType> get_surf_idx() const override;
146 
147  View<int*, HostType> pack_scalars(int max_nodes_on_surface, int n_xpts, int n_separatrices) const override;
148  void unpack_scalars(const View<int*, HostType>& scalars, int& max_nodes_on_surface, int& n_xpts, int& n_separatrices) override;
149 
150  int n_specified_xpts() const override;
151  View<Equil::XPoint*, HostType> get_xpts() const override;
152  RZPair get_axis() const override;
153  RZPair get_reference_point() const override;
154  View<double*, HostType> get_psiN() const override;
155 };
156 
157 struct GridFiles{
161 
162  std::vector<std::unique_ptr<PlaneFilesInterface>> plane_files_vec;
163 
164  // For test grids
165  GridFiles(int nshells, double raxis, double zaxis, double rscale, double zscale, GridCreateOpt grid_opt=Hexagonal);
166 
168 
169  std::string ifile_int2str(int ifile) const;
170 
171 #ifdef USE_MPI
172  GridFiles(NLReader::NamelistReader& nlr, bool is_stellarator, const MyMPI& my_mpi);
173 #endif
174 
175  const PlaneFilesInterface& lplane() const;
176 
177  const PlaneFilesInterface& midplane() const;
178 
179  const PlaneFilesInterface& rplane() const;
180 
181  int n_planes() const;
182 
183  void clear();
184 
185  static int nshells_from_nnodes(int nnodes);
186 };
187 
188 #endif
Definition: magnetic_field.hpp:12
Definition: NamelistReader.hpp:193
Definition: grid_files.hpp:51
int nsurfaces3b
Definition: grid_files.hpp:73
View< int *, HostType > nnodes_on_surface
Definition: grid_files.hpp:75
int nsurfaces1
Definition: grid_files.hpp:73
View< int *, HostType > non_aligned_nsurf
Definition: grid_files.hpp:79
View< RZPair *, HostType > x
Definition: grid_files.hpp:62
int num_non_aligned
Definition: grid_files.hpp:77
int ntriangles
Definition: grid_files.hpp:66
View< Vertex *, HostType > nodes
Definition: grid_files.hpp:67
View< int *, HostType > i_x
Definition: grid_files.hpp:70
int nnodes
Definition: grid_files.hpp:61
View< int *, HostType > rgn
Definition: grid_files.hpp:63
int nsurfaces2
Definition: grid_files.hpp:73
View< int **, CLayout, HostType > non_aligned_surf_idx
Definition: grid_files.hpp:80
View< double *, HostType > psiN
An array of normalized psi coordinates temporarily needed for the stellarator version.
Definition: grid_files.hpp:83
View< int **, CLayout, HostType > surf_idx
Definition: grid_files.hpp:81
int nsurfaces3a
Definition: grid_files.hpp:73
int nsurfaces
Definition: grid_files.hpp:72
int sum_nnodes_on_surfaces
Definition: grid_files.hpp:76
View< int *, HostType > i_surf_separatrices
Definition: grid_files.hpp:74
View< int *, HostType > non_aligned_vert
Definition: grid_files.hpp:78
FileRegion
Definition: grid_files.hpp:55
@ Inside
Definition: grid_files.hpp:56
@ OnWall
Definition: grid_files.hpp:57
Definition: grid_files.hpp:13
virtual View< int *, HostType > get_rgn() const =0
virtual int get_nsurfaces3a() const =0
virtual View< int **, CLayout, HostType > get_non_aligned_surf_idx() const =0
virtual int get_nsurfaces1() const =0
virtual int n_specified_xpts() const =0
virtual int get_nsurfaces3b() const =0
virtual View< RZPair *, HostType > get_x() const =0
virtual int get_nsurfaces() const =0
virtual int get_nnodes() const =0
virtual View< int *, HostType > get_non_aligned_nsurf() const =0
virtual RZPair get_reference_point() const =0
virtual View< Vertex *, HostType > get_nodes() const =0
virtual int get_ntriangles() const =0
virtual int get_num_non_aligned() const =0
virtual View< int *, HostType > get_nnodes_on_surface() const =0
virtual View< Equil::XPoint *, HostType > get_xpts() const =0
virtual View< int *, HostType > get_i_surf_separatrices() const =0
virtual int get_sum_nnodes_on_surfaces() const =0
virtual RZPair get_axis() const =0
virtual ~PlaneFilesInterface()=default
virtual View< double *, HostType > get_psiN() const =0
virtual View< int **, CLayout, HostType > get_surf_idx() const =0
virtual View< int *, HostType > get_i_x() const =0
virtual void unpack_scalars(const View< int *, HostType > &scalars, int &max_nodes_on_surface, int &n_xpts, int &n_separatrices)=0
virtual int get_nsurfaces2() const =0
virtual View< int *, HostType > get_non_aligned_vert() const =0
virtual View< int *, HostType > pack_scalars(int max_nodes_on_surface, int n_xpts, int n_separatrices) const =0
Definition: magnetic_field.F90:1
Definition: my_mpi.F90:1
Definition: grid_files.hpp:157
static int nshells_from_nnodes(int nnodes)
Definition: grid_files.cpp:918
int n_planes() const
Definition: grid_files.cpp:909
std::vector< std::unique_ptr< PlaneFilesInterface > > plane_files_vec
Definition: grid_files.hpp:162
static constexpr GridCreateOpt Hexagonal
Definition: grid_files.hpp:159
const PlaneFilesInterface & midplane() const
Definition: grid_files.cpp:901
const PlaneFilesInterface & rplane() const
Definition: grid_files.cpp:905
const PlaneFilesInterface & lplane() const
Definition: grid_files.cpp:897
GridFiles(int nshells, double raxis, double zaxis, double rscale, double zscale, GridCreateOpt grid_opt=Hexagonal)
Definition: grid_files.cpp:776
static constexpr GridCreateOpt Circular
Definition: grid_files.hpp:160
std::string ifile_int2str(int ifile) const
Definition: grid_files.cpp:784
void clear()
Definition: grid_files.cpp:913
Definition: my_mpi.hpp:19
Definition: grid_files.hpp:86
int get_nsurfaces() const override
Definition: grid_files.cpp:694
void unpack_scalars(const View< int *, HostType > &scalars, int &max_nodes_on_surface, int &n_xpts, int &n_separatrices) override
Definition: grid_files.cpp:50
PlaneFiles()
Definition: grid_files.hpp:93
static int nshells_from_nnodes(int nnodes)
Definition: grid_files.cpp:589
View< Vertex *, HostType > get_nodes() const override
Definition: grid_files.cpp:686
View< int *, HostType > get_non_aligned_nsurf() const override
Definition: grid_files.cpp:734
int get_nsurfaces3b() const override
Definition: grid_files.cpp:710
View< int *, HostType > get_non_aligned_vert() const override
Definition: grid_files.cpp:730
void convert_to_x_rgn_format(const View< double **, CLayout, HostType > &tmp_nodes, View< RZPair *, HostType > &x, View< int *, HostType > &rgn, bool is_stellarator)
Definition: grid_files.cpp:15
int get_num_non_aligned() const override
Definition: grid_files.cpp:726
View< int **, CLayout, HostType > get_non_aligned_surf_idx() const override
Definition: grid_files.cpp:738
View< int **, CLayout, HostType > get_surf_idx() const override
Definition: grid_files.cpp:742
RZPair get_axis() const override
Definition: grid_files.cpp:760
int get_ntriangles() const override
Definition: grid_files.cpp:682
int get_nsurfaces3a() const override
Definition: grid_files.cpp:706
View< int *, HostType > get_i_x() const override
Definition: grid_files.cpp:690
RZPair get_reference_point() const override
Definition: grid_files.cpp:765
void construct_hex_grid(int nshells, double raxis, double zaxis, double rscale, double zscale)
Definition: grid_files.cpp:372
View< int *, HostType > get_i_surf_separatrices() const override
Definition: grid_files.cpp:714
View< int *, HostType > pack_scalars(int max_nodes_on_surface, int n_xpts, int n_separatrices) const override
Definition: grid_files.cpp:33
View< double *, HostType > get_psiN() const override
Definition: grid_files.cpp:770
View< RZPair *, HostType > get_x() const override
Definition: grid_files.cpp:674
double node_to_node_dist(RZPair a, RZPair b)
Definition: grid_files.cpp:475
int get_nsurfaces2() const override
Definition: grid_files.cpp:702
void construct_circular_grid(int nshells, double raxis, double zaxis, double rscale, double zscale)
Definition: grid_files.cpp:488
void convert_to_node_format(const View< int **, CLayout, HostType > &tmp_ele, View< Vertex *, HostType > &nodes)
Definition: grid_files.cpp:24
View< int *, HostType > get_rgn() const override
Definition: grid_files.cpp:678
GridCreateOpt
Definition: grid_files.hpp:106
@ Hexagonal
Definition: grid_files.hpp:107
@ Circular
Definition: grid_files.hpp:108
int n_specified_xpts() const override
Definition: grid_files.cpp:746
int get_nsurfaces1() const override
Definition: grid_files.cpp:698
int get_nnodes() const override
Definition: grid_files.cpp:670
int get_sum_nnodes_on_surfaces() const override
Definition: grid_files.cpp:722
View< Equil::XPoint *, HostType > get_xpts() const override
Definition: grid_files.cpp:750
View< int *, HostType > get_nnodes_on_surface() const override
Definition: grid_files.cpp:718
Definition: grid_structs.hpp:28