ASCOT5
|
Mathematical utility functions. More...
Go to the source code of this file.
Functions | |
int | rcomp (const void *a, const void *b) |
Helper comparison routine for "math_rsearch". | |
double | math_simpson_helper (double(*f)(double), double a, double b, double eps, double S, double fa, double fb, double fc, int bottom) |
Helper routine for "math_simpson". | |
real | fmod (real x, real y) |
Compute the modulus of two real numbers. | |
void | math_jac_rpz2xyz (real *rpz, real *xyz, real r, real phi) |
Convert a Jacobian from cylindrical to cartesian coordinates. | |
void | math_jac_xyz2rpz (real *xyz, real *rpz, real r, real phi) |
Convert a Jacobian from cartesian to cylindrical coordinates. | |
void | math_matmul (real *matA, real *matB, int d1, int d2, int d3, real *matC) |
Matrix multiplication. | |
real | math_normal_rand (void) |
Generate normally distributed random numbers. | |
int | math_ipow (int a, int p) |
Calculate a^p where both a and p are integers (p >= 0) | |
double | math_simpson (double(*f)(double), double a, double b, double eps) |
Adaptive Simpsons rule for integral. | |
void | math_linspace (real *vec, real a, real b, int n) |
Generate linearly spaced vector. | |
void | math_uniquecount (int *in, int *unique, int *count, int n) |
Find unique numbers and their frequency in given array. | |
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. | |
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. | |
real * | math_rsearch (const real key, const real *base, int num) |
Search for array element preceding a key value. | |
int | math_point_in_polygon (real r, real z, real *rv, real *zv, int n) |
Check if coordinates are within polygon. | |
Mathematical utility functions.
Definition in file math.c.
int rcomp | ( | const void * | a, |
const void * | b ) |
Helper comparison routine for "math_rsearch".
This function checks if a key value is between two consecutive values in a real array. The array to be searched has to be in ascending sorted order.
a | first comparison value (key) |
b | first value in array |
double math_simpson_helper | ( | double(*)(double) | f, |
double | a, | ||
double | b, | ||
double | eps, | ||
double | S, | ||
double | fa, | ||
double | fb, | ||
double | fc, | ||
int | bottom ) |
Convert a Jacobian from cylindrical to cartesian coordinates.
This function converts a Jacobian located at angle phi and radius r from cylindrical to cartesian coordinates. The input Jacobian is an array with
[Br dBr/dr dBr/dphi dBr/dz Bphi dBphi/dr dBphi/dphi dBphi/dz Bz dBz/dr dBz/dphi dBz/dz]
an the output is
[Bx dBx/dx dBx/dy dBx/dz By dBy/dx dBy/dy dBy/dz Bz dBz/dx dBz/dy dBz/dz].
rpz | input rpz coordinates in a 3-length array |
xyz | output xyz coordinates in a 3-length array |
r | r coordinate of the gradient |
phi | phi coordinate of the gradient |
Convert a Jacobian from cartesian to cylindrical coordinates.
This function converts a Jacobian located at angle phi and radius r from cartesian to cylindrical coordinates. The input Jacobian is an array with
[Bx dBx/dx dBx/dy dBx/dz By dBy/dx dBy/dy dBy/dz Bz dBz/dx dBz/dy dBz/dz]
an the output is
[Br dBr/dr dBr/dphi dBr/dz Bphi dBphi/dr dBphi/dphi dBphi/dz Bz dBz/dr dBz/dphi dBz/dz].
xyz | output xyz coordinates in a 3-length array |
rpz | input rpz coordinates in a 3-length array |
r | r coordinate of the gradient |
phi | phi coordinate of the gradient |
Matrix multiplication.
This function performs matrix multiplication. The matrices are given as arrays, the indexing being row-wise.
matA | input array representing a d1 x d2 matrix |
matB | input array representing a d2 x d3 matrix |
d1 | input scalar dimension (columns in matA and matC) |
d2 | input scalar dimension (rows in matA and columns in matB) |
d3 | input scalar dimension (rows in matB and matC) |
matC | output array representing a d1 x d3 matrix |
real math_normal_rand | ( | void | ) |
Generate normally distributed random numbers.
This function generates normally distributed random numbers with expected value of 0 and variance of 1 using the polar method [1]. The rand48 random number generator should be initialized with srand48() before calling this function.
[1] G. Marsaglia, T. A. Bray. A convenient method for generating normal variables. Siam Review 6.3 (1964): 260-264. http://epubs.siam.org/doi/abs/10.1137/1006063:w
Currently only one of the generated numbers is returned; a second independent variate is Y = v2 * sqrt(-2*log(s) / s)
Try other random number generators such as those in Intel MKL
int math_ipow | ( | int | a, |
int | p ) |
double math_simpson | ( | double(*)(double) | f, |
double | a, | ||
double | b, | ||
double | eps ) |
Adaptive Simpsons rule for integral.
This function uses recursive splitting of the interval until either maximum number of intervals is reached or given relative error tolerance is obtained
f | pointer to the integrand function (of one variable) |
a | lower limit for the interval |
b | upper limit for the interval |
eps | absolute tolerance |
void math_uniquecount | ( | int * | in, |
int * | unique, | ||
int * | count, | ||
int | n ) |
Find unique numbers and their frequency in given array.
All argument arrays are assumed to have the same length as the array holding the input numbers. If not all numbers are unique, the remaining entries in unique and count arguments are zero.
in | input array of length n |
unique | array where unique numbers are stored, must be of length n |
count | count[i] is how many times unique[i] appears in array in, must be of length n |
n | length of the argument arrays |
Find if a point is on a given plane.
Let t1=(x1, y1, z1), t2=(x2, y2, z2), and t3=(x3, y3, z3) be non-collinear points on the plane. Let q =(x, y, z) be an arbitrary point. It lies in the plane if:
det( |x -x1 y -y1 z -z1| |x2-x1 y2-y1 z2-z1| |x3-x1 y3-y1 z3-z1| ) == 0
det( |x -x1 y -y1 z -z1| |x -x2 y -y2 z -z2| |x -x3 y -y3 z -z3| ) == 0
q | xyz coordinates of a query point |
t1 | xyz coordinates of the first point defining the plane |
t2 | xyz coordinates of the second point defining the plane |
t3 | xyz coordinates of the third point defining the plane |
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.
Let a,b,c be the vertices of a triangle and p an arbitrary point. The s,t are the barycentric coordinates of p if: p = (1-s-t)*a + s*b + t*c
AP | point position vector [x,y,z] |
AB | triangle edge AB vector |
AC | triangle edge AC vector |
n | normal vector [nx, ny, nz] of the triangle |
s | evaluated s coordinate |
t | evaluated t coordinate |
Search for array element preceding a key value.
This function takes an array of real values and a key value, and finds the array element which precedes the key value. The array to be searched has to be in ascending sorted order (according to comparison function rcomp).
key | real value that serves as key for the search |
base | pointer to the first object of the array to be searched |
num | number of elements in array |
Check if coordinates are within polygon.
This function checks if the given coordinates are within a 2D polygon using a modified axis crossing method [1]. Origin is moved to the coordinates and the number of wall segments crossing the positive r-axis are calculated. If this is odd, the point is inside the polygon.
[1] D.G. Alciatore, R. Miranda. A Winding Number and Point-in-Polygon Algorithm. Technical report, Colorado State University, 1995. http://www.engr.colostate.edu/~dga/dga/papers/point_in_polygon.pdf
r | r coordinate |
z | z coordinate |
rv | array of r points for the polygon |
zv | array of z points for the polygon |
n | number of points in the polygon |