8#include "offload_acc_omp.h"
13#define math_degrad 0.0174532925199432957692
15#define math_raddeg 57.295779513082320876798
18#define math_maxSimpsonDepth 20
21#define math_copy(a,b) do { a[0]=b[0];a[1]=b[1];a[2]=b[2]; } while(0)
24#define math_matcopy(a,b) do { a[0]=b[0];a[1]=b[1];a[2]=b[2];a[3]=b[3]; \
25 a[4]=b[4];a[5]=b[5];a[6]=b[6];a[7]=b[7];a[8]=b[8]; } while(0)
28#define math_dot(a,b) (a[0]*b[0]+a[1]*b[1]+a[2]*b[2])
31#define math_cross(a,b,c) do { c[0]=a[1]*b[2]-a[2]*b[1]; \
32 c[1]=a[2]*b[0]-a[0]*b[2];c[2]=a[0]*b[1]-a[1]*b[0]; } while(0)
39#define math_scalar_triple_product(a,b,c) \
40 ( - (a[2] * b[1] * c[0] ) \
41 + (a[1] * b[2] * c[0] ) \
42 + (a[2] * b[0] * c[1] ) \
43 - (a[0] * b[2] * c[1] ) \
44 - (a[1] * b[0] * c[2] ) \
45 + (a[0] * b[1] * c[2] ) )
48#define math_sumew(a,b) do { a[0]=a[0]+b[0];a[1]=a[1]+b[1];a[2]=a[2]+b[2] } while(0)
51#define math_prod(a,b) do { a[0]=a[0]*b;a[1]=a[1]*b;a[2]=a[2]*b; } while(0)
54#define math_matsumew(a,b) do { a[0]=a[0]+b[0];a[1]=a[1]+b[1];a[2]=a[2]+b[2]; \
55 a[3]=a[3]+b[3];a[4]=a[4]+b[4];a[5]=a[5]+b[5];a[6]=a[6]+b[6]; \
56 a[7]=a[7]+b[7];a[8]=a[8]+b[8]; } while(0)
59#define math_matprod(a,b) do { a[0]=a[0]*b;a[1]=a[1]*b;a[2]=a[2]*b; \
60 a[3]=a[3]*b;a[4]=a[4]*b;a[5]=a[5]*b;a[6]=a[6]*b;a[7]=a[7]*b; \
61 a[8]=a[8]*b } while(0)
64#define math_norm(a) (sqrt(a[0]*a[0]+a[1]*a[1]+a[2]*a[2]))
67#define math_normc(a1, a2, a3) (sqrt(a1*a1+a2*a2+a3*a3))
70#define math_unit(a, b) do {real _n=sqrt(a[0]*a[0]+a[1]*a[1]+a[2]*a[2]); \
71 b[0]=a[0]/_n;b[1]=a[1]/_n;b[2]=a[2]/_n; } while(0)
74#define math_xyz2rpz(xyz, rpz) do { rpz[0]=sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]); \
75 rpz[1]=atan2(xyz[1],xyz[0]);rpz[2]=xyz[2]; } while(0)
78#define math_rpz2xyz(rpz, xyz) do { xyz[0]=rpz[0]*cos(rpz[1]); \
79 xyz[1]=rpz[0]*sin(rpz[1]);xyz[2]=rpz[2]; } while(0)
83#define math_vec_rpz2xyz(vrpz, vxyz, phi) do { \
84 vxyz[0]=vrpz[0]*cos(phi)-vrpz[1]*sin(phi); \
85 vxyz[1]=vrpz[0]*sin(phi)+vrpz[1]*cos(phi); \
86 vxyz[2]=vrpz[2]; } while(0)
90#define math_vec_xyz2rpz(vxyz, vrpz, phi) do { \
91 vrpz[0]=vxyz[0]*cos(phi)+vxyz[1]*sin(phi); \
92 vrpz[1]=-vxyz[0]*sin(phi)+vxyz[1]*cos(phi); \
93 vrpz[2]=vxyz[2]; } while(0)
97#define math_determinant3x3( \
101 (x1) * ( (y2)*(z3) - (y3)*(z2) ) + \
102 (x2) * ( (y3)*(z1) - (y1)*(z3) ) + \
103 (x3) * ( (y1)*(z2) - (y2)*(z1) )
107#define math_deg2rad(a) (a * math_degrad)
110#define math_rad2deg(a) (a * math_raddeg)
117GPU_DECLARE_TARGET_SIMD
120GPU_DECLARE_TARGET_SIMD
127double math_simpson(
double (*f)(
double),
double a,
double b,
double epsilon);
137DECLARE_TARGET_SIMD_UNIFORM(rv,zv,n)
Main header file for ASCOT5.
void math_uniquecount(int *in, int *unique, int *count, int n)
Find unique numbers and their frequency in given array.
DECLARE_TARGET_END DECLARE_TARGET_SIMD real math_normal_rand()
Generate normally distributed random numbers.
DECLARE_TARGET_SIMD void math_barycentric_coords_triangle(real AP[3], real AB[3], real AC[3], real n[3], real *s, real *t)
Find barycentric coordinates for a given point.
double math_simpson(double(*f)(double), double a, double b, double epsilon)
Adaptive Simpsons rule for integral.
int math_point_in_polygon(real r, real z, real *rv, real *zv, int n)
Check if coordinates are within polygon.
DECLARE_TARGET_SIMD real * math_rsearch(const real key, const real *base, int num)
Search for array element preceding a key value.
GPU_DECLARE_TARGET_SIMD void math_jac_xyz2rpz(real *xyz, real *rpz, real r, real phi)
Convert a Jacobian from cartesian to cylindrical coordinates.
DECLARE_TARGET_END DECLARE_TARGET_SIMD void math_jac_rpz2xyz(real *rpz, real *xyz, real r, real phi)
Convert a Jacobian from cylindrical to cartesian coordinates.
DECLARE_TARGET_SIMD int math_ipow(int a, int p)
Calculate a^p where both a and p are integers (p >= 0)
DECLARE_TARGET_SIMD int math_point_on_plane(real q[3], real t1[3], real t2[3], real t3[3])
Find if a point is on a given plane.
DECLARE_TARGET real fmod(real x, real y)
Compute the modulus of two real numbers.
DECLARE_TARGET_END GPU_DECLARE_TARGET_SIMD void math_matmul(real *matA, real *matB, int d1, int d2, int d3, real *matC)
Matrix multiplication.
void math_linspace(real *vec, real a, real b, int n)
Generate linearly spaced vector.