Highly Efficient FFT for Exascale: HeFFTe v2.4
Loading...
Searching...
No Matches
heffte_utils.h
1/*
2 -- heFFTe --
3 Univ. of Tennessee, Knoxville
4 @date
5*/
6
7#ifndef HEFFTE_UTILS_H
8#define HEFFTE_UTILS_H
9
10#include <algorithm>
11#include <vector>
12#include <complex>
13#include <memory>
14#include <numeric>
15#include <algorithm>
16#include <functional>
17#include <cassert>
18#include <utility>
19#include <iostream>
20#include <ostream>
21#include <iomanip>
22#include <string>
23#include <deque>
24#include <fstream>
25#include <mpi.h>
26#include <limits>
27
28#include "heffte_config.h"
29
30namespace heffte {
31
36template<class T, class U = T>
38{
39 T old_value = std::move(obj);
40 obj = std::forward<U>(new_value);
41 return old_value;
42}
43
57namespace mpi {
58
79inline int comm_rank(MPI_Comm const comm){
80 int me;
81 MPI_Comm_rank(comm, &me);
82 return me;
83}
98inline bool world_rank(int rank){ return (comm_rank(MPI_COMM_WORLD) == rank); }
105inline int world_rank(){ return comm_rank(MPI_COMM_WORLD); }
106
117template<typename vector_like>
118void dump(int me, vector_like const &x, std::string const &message){
119 if (me < 0 or world_rank(me)){
120 std::cout << message << "\n";
121 for(auto i : x) std::cout << i << " ";
122 std::cout << std::endl;
123 }
124}
125
136inline int comm_size(MPI_Comm const comm){
137 int nprocs;
138 MPI_Comm_size(comm, &nprocs);
139 return nprocs;
140}
152inline MPI_Comm new_comm_from_group(std::vector<int> const &ranks, MPI_Comm const comm){
155 MPI_Group_incl(orig_group, (int) ranks.size(), ranks.data(), &new_group);
160 return result;
161}
162
175inline void comm_free(MPI_Comm const comm){
176 if (MPI_Comm_free(const_cast<MPI_Comm*>(&comm)) != MPI_SUCCESS)
177 throw std::runtime_error("Could not free a communicator.");
178}
179
190template<typename scalar> inline MPI_Datatype type_from(){
191 // note that "!std::is_same<scalar, scalar>::value" is always false,
192 // but will not be checked until the template is instantiated
193 static_assert(!std::is_same<scalar, scalar>::value, "The C++ type has unknown MPI equivalent.");
194 return MPI_BYTE; // some compilers complain about lack of return statement.
195}
200template<> inline MPI_Datatype type_from<int>(){ return MPI_INT; }
205template<> inline MPI_Datatype type_from<float>(){ return MPI_FLOAT; }
210template<> inline MPI_Datatype type_from<double>(){ return MPI_DOUBLE; }
221
222}
223
252template<typename scalar_type> struct is_ccomplex : std::false_type{};
270template<typename scalar_type> struct is_zcomplex : std::false_type{};
271
276template<> struct is_ccomplex<std::complex<float>> : std::true_type{};
281template<> struct is_zcomplex<std::complex<double>> : std::true_type{};
282
315template<typename, typename = void> struct define_standard_type{};
316
321template<> struct define_standard_type<float, void>{
323 using type = float;
324};
329template<> struct define_standard_type<double, void>{
331 using type = double;
332};
333
338template<typename scalar_type> struct define_standard_type<scalar_type, typename std::enable_if<is_ccomplex<scalar_type>::value>::type>{
340 using type = std::complex<float>;
341};
346template<typename scalar_type> struct define_standard_type<scalar_type, typename std::enable_if<is_zcomplex<scalar_type>::value>::type>{
348 using type = std::complex<double>;
349};
350
354template<typename scalar_type>
362template<typename scalar_type>
366
380template<typename some_class>
381int get_last_active(std::array<std::unique_ptr<some_class>, 4> const &shaper){
382 int last = -1;
383 for(int i=0; i<4; i++) if (shaper[i]) last = i;
384 return last;
385}
386
391template<typename some_class>
392int count_active(std::array<std::unique_ptr<some_class>, 4> const &shaper){
393 int num = 0;
394 for(int i=0; i<4; i++) if (shaper[i]) num++;
395 return num;
396}
397
402template<typename some_class>
403size_t get_max_box_size(std::array<some_class, 3> const &executors){
404 size_t max_size = (executors[0]) ? executors[0]->box_size() : 0;
405 max_size = std::max(max_size, (executors[1]) ? executors[1]->box_size() : static_cast<size_t>(0));
406 return std::max(max_size, (executors[2]) ? executors[2]->box_size() : static_cast<size_t>(0));
407}
408
413template<typename some_class>
414size_t get_max_box_size_r2c(std::array<some_class, 3> const &executors){
415 size_t max_size = (executors[0]) ? executors[0]->complex_size() : 0;
416 max_size = std::max(max_size, (executors[1]) ? executors[1]->box_size() : static_cast<size_t>(0));
417 return std::max(max_size, (executors[2]) ? executors[2]->box_size() : static_cast<size_t>(0));
418}
423template<typename some_class>
424size_t get_max_work_size(std::array<some_class, 3> const &executors){
425 size_t max_size = (executors[0]) ? executors[0]->workspace_size() : 0;
426 max_size = std::max(max_size, (executors[1]) ? executors[1]->workspace_size() : static_cast<size_t>(0));
427 return std::max(max_size, (executors[2]) ? executors[2]->workspace_size() : static_cast<size_t>(0));
428}
429
434template<typename some_class_r2c, typename some_class>
435size_t get_max_work_size(some_class_r2c const &executors_r2c, std::array<some_class, 2> const &executors){
436 return std::max(executors_r2c->workspace_size(), std::max(executors[0]->workspace_size(), executors[1]->workspace_size()));
437}
438
439}
440
441#endif /* HEFFTE_UTILS_H */
T c11_exchange(T &obj, U &&new_value)
Replace with the C++ 2014 std::exchange later.
Definition heffte_utils.h:37
int get_last_active(std::array< std::unique_ptr< some_class >, 4 > const &shaper)
Return the index of the last active (non-null) unique_ptr.
Definition heffte_utils.h:381
size_t get_max_box_size_r2c(std::array< some_class, 3 > const &executors)
Returns the max of the box_size() for each of the executors.
Definition heffte_utils.h:414
size_t get_max_box_size(std::array< some_class, 3 > const &executors)
Returns the max of the box_size() for each of the executors.
Definition heffte_utils.h:403
int count_active(std::array< std::unique_ptr< some_class >, 4 > const &shaper)
Return the number of active (non-null) unique_ptr.
Definition heffte_utils.h:392
size_t get_max_work_size(std::array< some_class, 3 > const &executors)
Returns the max of the workspace_size() for each of the executors.
Definition heffte_utils.h:424
MPI_Datatype type_from< float >()
Specialization to hand the float type.
Definition heffte_utils.h:205
MPI_Datatype type_from< int >()
Specialization to hand the int type.
Definition heffte_utils.h:200
int comm_rank(MPI_Comm const comm)
Returns the rank of this process within the specified comm.
Definition heffte_utils.h:79
void comm_free(MPI_Comm const comm)
Calls free on the MPI comm.
Definition heffte_utils.h:175
int comm_size(MPI_Comm const comm)
Returns the size of the specified communicator.
Definition heffte_utils.h:136
MPI_Datatype type_from< double >()
Specialization to hand the double type.
Definition heffte_utils.h:210
void dump(int me, vector_like const &x, std::string const &message)
Write the message and the data from the vector-like x, performed only on rank me (if positive),...
Definition heffte_utils.h:118
MPI_Datatype type_from()
Returns the MPI equivalent of the scalar C++ type.
Definition heffte_utils.h:190
MPI_Comm new_comm_from_group(std::vector< int > const &ranks, MPI_Comm const comm)
Creates a new sub-communicator from the provided processes in comm.
Definition heffte_utils.h:152
int world_rank()
Returns the rank of this process within the MPI_COMM_WORLD (useful for debugging).
Definition heffte_utils.h:105
Namespace containing all HeFFTe methods and classes.
Definition heffte_backend_cuda.h:38
define_standard_type< scalar_type >::type * convert_to_standard(scalar_type input[])
Converts an array of some type to an array of the C++ equivalent type.
Definition heffte_utils.h:355
std::complex< double > type
If heffte::is_ccomplex is true_type, then the type is equivalent to std::complex<double>.
Definition heffte_utils.h:348
std::complex< float > type
If heffte::is_ccomplex is true_type, then the type is equivalent to std::complex<float>.
Definition heffte_utils.h:340
Struct to specialize that returns the C++ equivalent of each type.
Definition heffte_utils.h:315
Struct to specialize to allow HeFFTe to recognize custom single precision complex types.
Definition heffte_utils.h:252
Struct to specialize to allow HeFFTe to recognize custom double precision complex types.
Definition heffte_utils.h:270
Wrapper around cufftHandle plans, set for float or double complex.
Definition heffte_backend_cuda.h:346