1 /* ----------------------------------------------------------------------
2 SPARTA - Stochastic PArallel Rarefied-gas Time-accurate Analyzer
3 http://sparta.sandia.gov
4 Steve Plimpton, sjplimp@sandia.gov, Michael Gallis, magalli@sandia.gov
5 Sandia National Laboratories
6
7 Copyright (2014) Sandia Corporation. Under the terms of Contract
8 DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
9 certain rights in this software. This software is distributed under
10 the GNU General Public License.
11
12 See the README file in the top-level SPARTA directory.
13 ------------------------------------------------------------------------- */
14
15 #include "string.h"
16 #include "compute_boundary_kokkos.h"
17 #include "particle_kokkos.h"
18 #include "mixture.h"
19 #include "grid.h"
20 #include "surf.h"
21 #include "update.h"
22 #include "domain.h"
23 #include "modify.h"
24 #include "memory_kokkos.h"
25 #include "error.h"
26 #include "sparta_masks.h"
27 #include "kokkos.h"
28
29 using namespace SPARTA_NS;
30
31 /* ---------------------------------------------------------------------- */
32
ComputeBoundaryKokkos(SPARTA * sparta,int narg,char ** arg)33 ComputeBoundaryKokkos::ComputeBoundaryKokkos(SPARTA *sparta, int narg, char **arg) :
34 ComputeBoundary(sparta, narg, arg)
35 {
36 kokkos_flag = 1;
37
38 memory->destroy(array);
39 memory->destroy(myarray);
40 memoryKK->create_kokkos(k_array,array,size_array_rows,size_array_cols,"boundary:array");
41 memoryKK->create_kokkos(k_myarray,myarray,size_array_rows,size_array_cols,"boundary:myarray");
42 d_myarray = k_myarray.d_view;
43 d_array = k_array.d_view;
44 d_which = DAT::t_int_1d("boundary:which",nvalue);
45 }
46
ComputeBoundaryKokkos(SPARTA * sparta)47 ComputeBoundaryKokkos::ComputeBoundaryKokkos(SPARTA *sparta) :
48 ComputeBoundary(sparta)
49 {
50 array = NULL;
51 myarray = NULL;
52 which = NULL;
53 id = NULL;
54 style = NULL;
55 tlist = NULL;
56 }
57
58 /* ---------------------------------------------------------------------- */
59
~ComputeBoundaryKokkos()60 ComputeBoundaryKokkos::~ComputeBoundaryKokkos()
61 {
62 if (copy || copymode) return;
63
64 memoryKK->destroy_kokkos(k_array,array);
65 memoryKK->destroy_kokkos(k_myarray,myarray);
66 }
67
68 /* ---------------------------------------------------------------------- */
69
init()70 void ComputeBoundaryKokkos::init()
71 {
72 ComputeBoundary::init();
73
74 auto h_which = Kokkos::create_mirror_view(d_which);
75 for (int n=0; n<nvalue; n++)
76 h_which(n) = which[n];
77 Kokkos::deep_copy(d_which,h_which);
78 }
79
80 /* ---------------------------------------------------------------------- */
81
compute_array()82 void ComputeBoundaryKokkos::compute_array()
83 {
84 invoked_array = update->ntimestep;
85
86 // sum tallies across processors
87
88 if (sparta->kokkos->gpu_direct_flag) {
89 MPI_Allreduce(d_myarray.data(),d_array.data(),nrow*ntotal,
90 MPI_DOUBLE,MPI_SUM,world);
91 k_array.modify_device();
92 k_array.sync_host();
93 } else {
94 k_myarray.modify_device();
95 k_myarray.sync_host();
96 MPI_Allreduce(k_myarray.h_view.data(),k_array.h_view.data(),nrow*ntotal,
97 MPI_DOUBLE,MPI_SUM,world);
98 }
99
100 // normalize tallies
101
102 int m;
103 for (int j = 0; j < ntotal; j++) {
104 m = j % nvalue;
105 if (which[m] != NUM && which[m] != NUMWT)
106 for (int i = 0; i < size_array_rows; i++)
107 array[i][j] /= normflux[i];
108 }
109 }
110
111 /* ---------------------------------------------------------------------- */
112
clear()113 void ComputeBoundaryKokkos::clear()
114 {
115 // reset tally values to zero
116 // called by Update at beginning of timesteps boundary tallying is done
117
118 Kokkos::deep_copy(d_myarray,0.0);
119 }
120
pre_boundary_tally()121 void ComputeBoundaryKokkos::pre_boundary_tally()
122 {
123 mvv2e = update->mvv2e;
124
125 ParticleKokkos* particle_kk = (ParticleKokkos*) particle;
126 particle_kk->sync(Device,PARTICLE_MASK|SPECIES_MASK);
127 d_species = particle_kk->k_species.d_view;
128 d_s2g = particle_kk->k_species2group.d_view;
129
130 need_dup = sparta->kokkos->need_dup<DeviceType>();
131 if (need_dup)
132 dup_myarray = Kokkos::Experimental::create_scatter_view<typename Kokkos::Experimental::ScatterSum, typename Kokkos::Experimental::ScatterDuplicated>(d_myarray);
133 else
134 ndup_myarray = Kokkos::Experimental::create_scatter_view<typename Kokkos::Experimental::ScatterSum, typename Kokkos::Experimental::ScatterNonDuplicated>(d_myarray);
135 }
136
post_boundary_tally()137 void ComputeBoundaryKokkos::post_boundary_tally()
138 {
139 if (need_dup) {
140 Kokkos::Experimental::contribute(d_myarray, dup_myarray);
141 dup_myarray = decltype(dup_myarray)(); // free duplicated memory
142 }
143 }
144