globals.hpp

Go to the documentation of this file.

#ifndef GLOBALS_HPP
#define GLOBALS_HPP
#include <limits.h>
#include <limits>
#include <string>
#include <cassert>
#include "space_settings.hpp"
#include "array_deep_copy.hpp"
#include "access_add.hpp"
#include "simd.hpp"
#ifdef USE_MPI
#include "my_mpi.hpp"
#endif
#include "constants.hpp"
 
/* Returns max number of omp threads, or 1 if no omp
 */
inline int get_num_cpu_threads(){
#ifdef USE_OMP
    return omp_get_max_threads();
#else
    return 1;
#endif
}
 
/* Return true if MPI rank is zero
 * */
inline bool is_rank_zero(){
#ifdef USE_MPI
    return SML_COMM_RANK==0;
#else
    return true;
#endif
}
 
/* Safely abort and exit the code
 * */
inline void exit_XGC(std::string msg){
    printf("%s",msg.c_str());
    fflush(stdout);
#ifdef USE_MPI
    MPI_Abort(SML_COMM_WORLD, 1);
#else
    exit(1);
#endif
}
 
/* Asserts condition
 * (Callable from Device)
 * (C++ has no standardized assert with a message)
 * */
KOKKOS_INLINE_FUNCTION void assert_XGC(bool cond, const char* msg){
    if(!cond){
        DEVICE_PRINTF("%s",msg);
        assert(false);
    }
}
 
/* Check that two integers multiplied together doesn't cause an overflow
 * */
inline bool causes_multiplication_overflow(int a,int b){
    int c = a*b;
    if(b==0) return false; // No overflow
    return !(c/b == a);
}
 
inline bool causes_multiplication_overflow(long long int a,long long int b){
    if (a==0 || b==0) return false;
#if defined(__SIZEOF_INT128__)
    __int128 prod = static_cast<__int128>(a) * static_cast<__int128>(b);
    return prod > std::numeric_limits<long long int>::max() || prod < std::numeric_limits<long long int>::lowest();
#else
    if (a > 0){
        if (b > 0){
            return a > std::numeric_limits<long long int>::max()/b;
        } else {
            return b < std::numeric_limits<long long int>::lowest()/a;
        }
    } else {
        if (b > 0){
            return a < std::numeric_limits<long long int>::lowest()/b;
        } else {
            return a != 0 && b < std::numeric_limits<long long int>::max()/a;
        }
    }
#endif
}
 
/* Check that two integers added together doesn't cause an overflow
 * */
inline bool causes_addition_overflow(int a,int b){
    if ( ((b > 0) && (a > INT_MAX - b))    // a + b would overflow
       ||((b < 0) && (a < INT_MIN - b)) ){ // a + b would underflow
        return true;
    } else {
        return false;
    }
}
 
enum class Order{
    Zero,
    One,
    Two
};
 
enum SpeciesType {
    ELECTRON = 0,
    MAIN_ION = 1
};
 
enum KinType{
    DriftKin=0,
    GyroKin=1
};
 
/* PhiInterpType specifies whether the field uses two values in order to interpolate in the phi direction
 * */
enum class PhiInterpType{
    Planes,
    None
};
 
// Eventually, nearly the entire code will be templated on PhiInterpType, which is XGC1 vs XGCa
// As an intermediate step, it is convenient to specify the template in one global location
#ifdef XGC1
constexpr PhiInterpType PIT_GLOBAL = PhiInterpType::Planes;
#else
constexpr PhiInterpType PIT_GLOBAL = PhiInterpType::None;
#endif
 
/* What type of marker weight algorithm to use, reduced delta-f or total-f
 */
enum class MarkerType{
    ReducedDeltaF,
    FullF,
    TotalF,
    None
};
enum class FAnalyticShape{
    Maxwellian,
    SlowingDown,
    None
};
enum class WeightEvoEq{
    Direct,
    PDE,
    None
};
 
enum class MagneticFieldMode{
    Electromagnetic,
    Electrostatic
};
#ifdef EXPLICIT_EM
constexpr MagneticFieldMode MFM_GLOBAL = MagneticFieldMode::Electromagnetic;
#else
constexpr MagneticFieldMode MFM_GLOBAL = MagneticFieldMode::Electrostatic;
#endif
 
enum class BFieldSymmetry{
    Tokamak,
    Stellarator
};
#ifdef STELLARATOR
constexpr BFieldSymmetry BFS_GLOBAL = BFieldSymmetry::Stellarator;
#else
constexpr BFieldSymmetry BFS_GLOBAL = BFieldSymmetry::Tokamak;
#endif
 
/* Whether to use cylindrical limit geometry with periodic boundary function
 */
enum class GeometryType{
    Toroidal,
    CylindricalLimit
};
#ifdef CYLINDRICAL
constexpr GeometryType GEOMETRY = GeometryType::CylindricalLimit;
#else
constexpr GeometryType GEOMETRY = GeometryType::Toroidal;
#endif
 
template<GeometryType GT>
KOKKOS_INLINE_FUNCTION double geometry_switch(double a, double b);
 
template<>
KOKKOS_INLINE_FUNCTION double geometry_switch<GeometryType::Toroidal>(double a, double b){
    return a;
}
 
template<>
KOKKOS_INLINE_FUNCTION double geometry_switch<GeometryType::CylindricalLimit>(double a, double b){
    return b;
}
 
 
template<GeometryType GT>
KOKKOS_INLINE_FUNCTION constexpr bool use_toroidal_terms();
 
template<>
KOKKOS_INLINE_FUNCTION constexpr bool use_toroidal_terms<GeometryType::CylindricalLimit>(){
    return 0;
}
 
template<>
KOKKOS_INLINE_FUNCTION constexpr bool use_toroidal_terms<GeometryType::Toroidal>(){
    return 1;
}
 
enum class PullbackMethod{
    Electrostatic,
    IdealMHD
};
 
/***/
 
enum ParticlePhase {
    PIR = 0, 
    PIZ,     
    PIP,     
    PIRHO,   
    PIW1,    
    PIW2,    
    PTL_NPHASE
};
 
enum ParticleConsts {
    PIM = 0, 
    PIW0,    
    PIF0,    
#ifdef PTL_G2
    PIG2,    
#endif
    PTL_NCONST
};
 
// Divide two integers, then round up
KOKKOS_INLINE_FUNCTION int divide_and_round_up(int a, int b){
    return (a+b-1)/b;
}
 
KOKKOS_INLINE_FUNCTION long long int divide_and_round_up(long long int a, long long int b){
    return (a+b-1)/b;
}
 
// Positive modulo. i.e. positive_modulo(-1,3) will return 2, whereas (-1)%3 returns -1.
KOKKOS_INLINE_FUNCTION unsigned positive_modulo( int value, unsigned m) {
    int mod = value % (int)m;
    if (mod < 0) {
        mod += m;
    }
    return mod;
}
 
/* Returns the offset of an evenly distributed set of objects into even subsets. If the set can be distributed evenly,
 * then this operation is simply the local subset times the result of the division. If not, then the first k
 * subsets will have one extra object, where k is the remainder.
 *
 * n_obj      is the number of objects that need to be divided
 * n_subsets  is the number of ranks the objects are being divided into
 * i_subset   is the index of this particular subset
 *
 */
inline long long int offsets_of_even_distribution(long long int n_obj, long long int n_subsets, long long int i_subset){
    long long int obj_per_subset = n_obj/n_subsets;
    long long int remainder = n_obj%n_subsets;
    if (i_subset<remainder){
        obj_per_subset += 1;
        return i_subset*obj_per_subset;
    }else{
        return i_subset*obj_per_subset + remainder;
    }
}
 
/* Returns the count of an evenly distributed set of objects into even subsets. If the set can be distributed evenly,
 * then this operation is simply the result of the division. If not, then the first k
 * subsets will have one extra object, where k is the remainder.
 *
 * n_obj      is the number of objects that need to be divided
 * n_subsets  is the number of ranks the objects are being divided into
 * i_subset   is the index of this particular subset
 *
 */
inline long long int counts_of_even_distribution(long long int n_obj, long long int n_subsets, long long int i_subset){
    long long int obj_per_subset = n_obj/n_subsets;
    long long int remainder = n_obj%n_subsets;
    if (i_subset<remainder){
        obj_per_subset += 1;
    }
    return obj_per_subset;
}
 
// Converts integer to string and adds leading zeros
inline std::string formatted_int2str(int input, int n_digits){
    std::string string_no_leading = std::to_string(input); // no leading zeros
 
    int initial_length = string_no_leading.length();
    int zeros_to_add = n_digits - std::min(n_digits, initial_length);
 
    // Add leading zeros
    return std::string(zeros_to_add, '0') + string_no_leading;
}
 
// is_same_type is equivalent to std::is_same, but can be used in cuda kernels
template<class T, class U> struct is_same_type{static constexpr bool val=false;};
template<class T         > struct is_same_type<T, T>{static constexpr bool val=true;};
 
// Defined and initialized in my_mpi.cpp, reset in sml.tpp
extern bool global_debug_flag;
extern bool global_perf_barriers_flag;
 
#endif