XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
field_following_coordinates.hpp
Go to the documentation of this file.
1 #ifndef FIELD_FOLLOWING_COORDINATES_HPP
2 #define FIELD_FOLLOWING_COORDINATES_HPP
3 
4 #include "grid.hpp"
5 
6 extern "C" int* get_psn_ff_hdp_tr_loc();
7 extern "C" double* get_psn_ff_hdp_p_loc();
8 extern "C" double* get_grid_node_vol_ff_loc();
9 
11 
12  Kokkos::View<int**,Kokkos::LayoutRight,DeviceType> tr;
13  Kokkos::View<double***,Kokkos::LayoutRight,DeviceType> p;
14 
15  // View for temporarily storing result data for in-place conversions
16  Kokkos::View<double**,Kokkos::LayoutRight,DeviceType> result;
17 
18  public:
19 
21 
23  : tr(Kokkos::ViewAllocateWithoutInitializing("tr"),2,grid.nnode),
24  p(Kokkos::ViewAllocateWithoutInitializing("p"),2,grid.nnode,3),
25  result(Kokkos::ViewAllocateWithoutInitializing("result"),2,grid.nnode)
26  {
27  // Copy grid conversion info to device
30  }
31 
32  /* Conversion from field-following coordinates to toroidal coordinates
33  * @param[in] grid is the Grid object needed here to map elements to nodes
34  * @param[in] input is the array one wants to convert
35  * @param[out] output is the resulting converted array
36  * */
37  void cnvt_grid_ff2real(const Grid<DeviceType>& grid, const Kokkos::View<double**,Kokkos::LayoutRight,DeviceType>& input, const Kokkos::View<double**,Kokkos::LayoutRight,DeviceType>& output){
38  // Check the input and output views are the correct size
39  if(input.extent(0) != 2 || input.extent(1) != grid.nnode ||
40  output.extent(0)!= 2 || output.extent(1)!= grid.nnode ) exit_XGC("\nError in cnvt_grid_ff2real: Unexpected View sizes. \n");
41 
42  // Initialize output to zero
43  Kokkos::deep_copy(output, 0.0);
44 
45  // Initialize vol to zero
46  Kokkos::View<double**,Kokkos::LayoutRight,DeviceType> vol("vol",input.layout());
47 
48  // Copy grid%node_vol_ff from fortran
49  Kokkos::View<double**,Kokkos::LayoutRight,DeviceType> node_vol_ff(Kokkos::ViewAllocateWithoutInitializing("node_vol_ff"),input.layout());
51 
52  // Make a (SHALLOW) copy so the lambda doesn't have to implicitly capture *this
53  FieldFollowingCoordinates ff = *this;
54 
55  // for all node point of ff
56  Kokkos::parallel_for("cnvt_grid_ff2real", Kokkos::RangePolicy<ExSpace>(0,grid.nnode), KOKKOS_LAMBDA( const int i ){
57  for (int dir=0;dir<2;dir++){
58  for (int j=0;j<3;j++){
59  // node and weight
60  int nd=grid.nodes(ff.tr(dir,i)-1).vertex[j] - 1; // 0-index
61  double wt=ff.p(dir,i,j)*node_vol_ff(dir,nd);
62 
63  // accumulate values
64  Kokkos::atomic_add(&(output(dir,nd)), wt*input(dir,i));
65  Kokkos::atomic_add(&(vol(dir,nd)), wt);
66  }
67  }
68  });
69 
70  // Normalize to vol
71  Kokkos::parallel_for("cnvt_grid_ff2real_norm", Kokkos::RangePolicy<ExSpace>(0,grid.nnode), KOKKOS_LAMBDA( const int i ){
72  for (int dir=0;dir<2;dir++){
73  if(vol(dir, i) == 0.0){
74  output(dir, i) = 0.0;
75  }else{
76  output(dir, i) /= vol(dir, i);
77  }
78  }
79  });
80 
81  // fill empty area --- replace nan values
82  Kokkos::parallel_for("cnvt_grid_ff2real_fillempty", Kokkos::RangePolicy<ExSpace>(0,grid.nnode), KOKKOS_LAMBDA( const int i ){
83  for (int dir=0; dir<2; dir++){
84  if(vol(dir,i)==0.0){
85  int cdir = 1-dir;
86  output(dir,i)=0.0;
87  for (int j=0;j<3;j++){
88  int nd=grid.nodes(ff.tr(cdir,i)-1).vertex[j] - 1; // 0-index
89  double wt=ff.p(cdir,i,j);
90 
91  // obtain interpolated values
92  output(dir,i) += wt*input(dir,nd);
93  }
94  }
95  }
96  });
97  Kokkos::fence();
98  }
99 
100  /* Conversion from real coordinates to field-following coordinates
101  * @param[in] grid is the Grid object needed here to map elements to nodes
102  * @param[in] input is the array one wants to convert
103  * @param[out] output is the resulting converted array
104  * */
105  void cnvt_grid_real2ff(const Grid<DeviceType>& grid, const Kokkos::View<double**,Kokkos::LayoutRight,DeviceType>& input, const Kokkos::View<double**,Kokkos::LayoutRight,DeviceType>& output){
106  // Check the input and output views are the correct size
107  if(input.extent(0) != 2 || input.extent(1) != grid.nnode ||
108  output.extent(0)!= 2 || output.extent(1)!= grid.nnode ) exit_XGC("\nError in cnvt_grid_real2ff: Unexpected View sizes. \n");
109 
110  // Initialize output to zero
111  Kokkos::deep_copy(output, 0.0);
112 
113 
114  // Make a (SHALLOW) copy so the lambda doesn't have to implicitly capture *this
115  FieldFollowingCoordinates ff = *this;
116 
117  Kokkos::parallel_for("cnvt_grid_real2ff", Kokkos::RangePolicy<ExSpace>(0,grid.nnode), KOKKOS_LAMBDA( const int i ){
118  for (int dir=0;dir<2;dir++){
119  for (int j=0;j<3;j++){
120  int nd=grid.nodes(ff.tr(dir,i)-1).vertex[j] - 1; // 0-index
121 
122  double wt=ff.p(dir,i,j);
123 
124  // obtain interpolated values
125  output(dir,i) += wt*input(dir,nd);
126  }
127  }
128  });
129  Kokkos::fence();
130  }
131 
132  /* In-place conversion from real coordinates to field-following coordinates, using a buffer array to store the result temporarily
133  * before copying back to the input array
134  * @param[in] grid is the Grid object needed here to map elements to nodes
135  * @param[inout] input is the array one wants to convert, and the resulting converted array
136  * */
137  void cnvt_grid_real2ff(const Grid<DeviceType>& grid, const Kokkos::View<double**,Kokkos::LayoutRight,DeviceType>& input){
138  cnvt_grid_real2ff(grid, input, result);
139 
140  // Copy back to original array
141  Kokkos::deep_copy(input, result);
142  }
143 };
144 
145 #endif
void array_deep_copy(T *array, const Kokkos::View< T *, Kokkos::LayoutRight, Device > &view)
Definition: array_deep_copy.hpp:11
Definition: field_following_coordinates.hpp:10
Kokkos::View< double **, Kokkos::LayoutRight, DeviceType > result
Definition: field_following_coordinates.hpp:16
double * get_grid_node_vol_ff_loc()
int * get_psn_ff_hdp_tr_loc()
FieldFollowingCoordinates()
Definition: field_following_coordinates.hpp:20
Kokkos::View< int **, Kokkos::LayoutRight, DeviceType > tr
Definition: field_following_coordinates.hpp:12
FieldFollowingCoordinates(const Grid< DeviceType > &grid)
Definition: field_following_coordinates.hpp:22
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const Kokkos::View< double **, Kokkos::LayoutRight, DeviceType > &input, const Kokkos::View< double **, Kokkos::LayoutRight, DeviceType > &output)
Definition: field_following_coordinates.hpp:105
Kokkos::View< double ***, Kokkos::LayoutRight, DeviceType > p
Definition: field_following_coordinates.hpp:13
void cnvt_grid_ff2real(const Grid< DeviceType > &grid, const Kokkos::View< double **, Kokkos::LayoutRight, DeviceType > &input, const Kokkos::View< double **, Kokkos::LayoutRight, DeviceType > &output)
Definition: field_following_coordinates.hpp:37
void exit_XGC(std::string msg)
Definition: globals.hpp:36
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const Kokkos::View< double **, Kokkos::LayoutRight, DeviceType > &input)
Definition: field_following_coordinates.hpp:137
double * get_psn_ff_hdp_p_loc()
void parallel_for(const std::string name, int n_ptl, Function func, Option option, HostAoSoA aosoa_h, DeviceAoSoA aosoa_d)
Definition: streamed_parallel_for.hpp:252
int nnode
Number of grid nodes.
Definition: grid.hpp:83
Kokkos::View< Vertex *, Kokkos::LayoutRight, Device > nodes
Definition: grid.hpp:92