1 #ifndef FIELD_FOLLOWING_COORDINATES_HPP
2 #define FIELD_FOLLOWING_COORDINATES_HPP
10 View<int**,CLayout,DeviceType>
tr;
11 View<double***,CLayout,DeviceType>
p;
14 View<double**,CLayout,DeviceType>
result;
22 p(
NoInit(
"p"),projection.
p.layout()),
26 Kokkos::deep_copy(
tr, projection.
tr);
27 Kokkos::deep_copy(
p, projection.
p);
37 if(input.extent(0) != 2 || input.extent(1) != grid.
midplane.nnode ||
38 output.extent(0)!= 2 || output.extent(1)!= std::max(grid.
rplane.nnode, grid.
lplane.nnode) )
exit_XGC(
"\nError in cnvt_grid_ff2real: Unexpected View sizes. \n");
41 Kokkos::deep_copy(output, 0.0);
44 View<double**,CLayout,DeviceType> vol(
"vol",input.layout());
47 View<double**,CLayout,DeviceType> node_vol_ff(
NoInit(
"node_vol_ff"),input.layout());
56 constexpr
int LEFTWARD=0;
58 for (
int j=0;j<3;j++){
60 int nd=grid.
lplane.get_node_index(ff.
tr(LEFTWARD,i), j);
61 double wt=ff.
p(LEFTWARD,i,j)*node_vol_ff(LEFTWARD,nd);
64 Kokkos::atomic_add(&(output(LEFTWARD,nd)), wt*input(LEFTWARD,i));
65 Kokkos::atomic_add(&(vol(LEFTWARD,nd)), wt);
68 constexpr
int RIGHTWARD=1;
70 for (
int j=0;j<3;j++){
72 int nd=grid.
rplane.get_node_index(ff.
tr(RIGHTWARD,i), j);
73 double wt=ff.
p(RIGHTWARD,i,j)*node_vol_ff(RIGHTWARD,nd);
76 Kokkos::atomic_add(&(output(RIGHTWARD,nd)), wt*input(RIGHTWARD,i));
77 Kokkos::atomic_add(&(vol(RIGHTWARD,nd)), wt);
82 Kokkos::parallel_for(
"cnvt_grid_ff2real_norm", Kokkos::RangePolicy<ExSpace>(0,output.extent(1)), KOKKOS_LAMBDA(
const int i ){
83 for (
int dir=0;dir<2;dir++){
84 if(vol(dir, i) == 0.0){
87 output(dir, i) /= vol(dir, i);
98 Kokkos::parallel_for(
"cnvt_grid_ff2real_fillempty", Kokkos::RangePolicy<ExSpace>(0,output.extent(1)), KOKKOS_LAMBDA(
const int i ){
99 if(vol(LEFTWARD,i)==0.0){
100 output(LEFTWARD,i)=0.0;
101 for (
int j=0;j<3;j++){
102 int nd=grid.
lplane.get_node_index(ff.
tr(RIGHTWARD,i), j);
103 double wt=ff.
p(RIGHTWARD,i,j);
106 output(LEFTWARD,i) += wt*input(LEFTWARD,nd);
110 Kokkos::parallel_for(
"cnvt_grid_ff2real_fillempty", Kokkos::RangePolicy<ExSpace>(0,output.extent(1)), KOKKOS_LAMBDA(
const int i ){
111 if(vol(RIGHTWARD,i)==0.0){
112 output(RIGHTWARD,i)=0.0;
113 for (
int j=0;j<3;j++){
114 int nd=grid.
rplane.get_node_index(ff.
tr(LEFTWARD,i), j);
115 double wt=ff.
p(LEFTWARD,i,j);
118 output(RIGHTWARD,i) += wt*input(RIGHTWARD,nd);
133 void cnvt_grid_ff2real(
const Grid<DeviceType>& grid,
const View<double**,CLayout,HostType>& input,
const View<double**,CLayout,HostType>& output,
const View<double**,CLayout,DeviceType>& tmp)
const{
135 Kokkos::deep_copy(tmp, input);
141 Kokkos::deep_copy(output,
result);
153 Kokkos::deep_copy(output, input);
159 Kokkos::deep_copy(output,
result);
171 if(input.extent(0) != 2 || input.extent(1) != std::max(grid.
rplane.nnode, grid.
lplane.nnode) ||
172 output.extent(0)!= 2 || output.extent(1)!= grid.
midplane.nnode )
exit_XGC(
"\nError in cnvt_grid_real2ff: Unexpected View sizes. \n");
175 Kokkos::deep_copy(output, 0.0);
181 constexpr
int LEFTWARD=0;
183 for (
int j=0;j<3;j++){
184 int nd=grid.
lplane.get_node_index(ff.
tr(LEFTWARD,i), j);
186 double wt=ff.
p(LEFTWARD,i,j);
189 output(LEFTWARD,i) += wt*input(LEFTWARD,nd);
192 constexpr
int RIGHTWARD=1;
194 for (
int j=0;j<3;j++){
195 int nd=grid.
rplane.get_node_index(ff.
tr(RIGHTWARD,i), j);
197 double wt=ff.
p(RIGHTWARD,i,j);
200 output(RIGHTWARD,i) += wt*input(RIGHTWARD,nd);
215 Kokkos::deep_copy(input,
result);
224 for(
int i = 0; i<input.extent(0); i++){
238 for(
int i = 0; i<input.extent(0); i++){
252 Kokkos::View<double**,Kokkos::LayoutRight,DeviceType> dest_view(
NoInit(
"dest_view"), nphi, field_in.nnode());
255 field_in.copy_to_double_view(
my_subview(dest_view, 0));
256 field_in.copy_to_double_view(
my_subview(dest_view, 1));
262 field_out.transpose_copy_from_double_view(
result);
Definition: field_following_coordinates.hpp:9
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const View< double **, CLayout, DeviceType > &input, const View< double **, Kokkos::LayoutRight, DeviceType > &output) const
Definition: field_following_coordinates.hpp:169
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const ScalarGridField &field_in, const GridField< DeviceType, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::DriftKin > &field_out) const
Definition: field_following_coordinates.hpp:247
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const View< double ***, CLayout, DeviceType > &input) const
Definition: field_following_coordinates.hpp:237
View< double **, CLayout, DeviceType > result
Definition: field_following_coordinates.hpp:14
View< int **, CLayout, DeviceType > tr
Definition: field_following_coordinates.hpp:10
void cnvt_grid_ff2real(const Grid< DeviceType > &grid, const View< double **, CLayout, HostType > &input, const View< double **, CLayout, HostType > &output, const View< double **, CLayout, DeviceType > &tmp) const
Definition: field_following_coordinates.hpp:133
void cnvt_grid_ff2real(const Grid< DeviceType > &grid, const View< double **, CLayout, DeviceType > &input, const View< double **, Kokkos::LayoutRight, DeviceType > &output) const
Definition: field_following_coordinates.hpp:35
FieldFollowingCoordinates(const Projection< HostType > &projection)
Definition: field_following_coordinates.hpp:20
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const View< double **, CLayout, DeviceType > &input) const
Definition: field_following_coordinates.hpp:211
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const View< double ***, CLayout, DeviceType > &input, const View< double ***, Kokkos::LayoutRight, DeviceType > &output) const
Definition: field_following_coordinates.hpp:223
View< double ***, CLayout, DeviceType > p
Definition: field_following_coordinates.hpp:11
void cnvt_grid_ff2real(const Grid< DeviceType > &grid, const View< double **, CLayout, HostType > &input, const View< double **, CLayout, DeviceType > &output) const
Definition: field_following_coordinates.hpp:151
FieldFollowingCoordinates()
Definition: field_following_coordinates.hpp:18
VolumesAndAreas volumes_and_areas
Definition: grid.hpp:191
Plane< Device > lplane
Definition: grid.hpp:170
Plane< Device > rplane
Definition: grid.hpp:171
Plane< Device > midplane
Definition: grid.hpp:169
void exit_XGC(std::string msg)
Definition: globals.hpp:37
Kokkos::View< T *, Kokkos::LayoutRight, Device > my_subview(const Kokkos::View< T ****, Kokkos::LayoutRight, Device > &view, int i, int j, int k)
Definition: my_subview.hpp:8
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
Kokkos::ViewAllocateWithoutInitializing NoInit
Definition: space_settings.hpp:69
Definition: grid_field.hpp:22
View< double ***, CLayout, Device > p
Definition: ff_projection.hpp:24
View< int **, CLayout, Device > tr
The triangles that are mapped to by the projection. One-indexed.
Definition: ff_projection.hpp:22
View< double **, CLayout, HostType > node_vol_ff_l_h
Definition: volumes_and_areas.hpp:13
View< double **, CLayout, HostType > node_vol_ff_r_h
Definition: volumes_and_areas.hpp:14