ASCOT5
Loading...
Searching...
No Matches
Functions
math.c File Reference

Mathematical utility functions. More...

#include <math.h>
#include <stdlib.h>
#include "ascot5.h"

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.
 
realmath_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.
 

Detailed Description

Mathematical utility functions.

Definition in file math.c.

Function Documentation

◆ rcomp()

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.

Parameters
afirst comparison value (key)
bfirst value in array
Returns
0 if key is between the consecutive values, 1 if key is greater than the first value, and -1 if key is smaller than the first value.

Definition at line 381 of file math.c.

◆ math_simpson_helper()

double math_simpson_helper ( double(*)(double) f,
double a,
double b,
double eps,
double S,
double fa,
double fb,
double fc,
int bottom )

Helper routine for "math_simpson".

Definition at line 449 of file math.c.

◆ fmod()

real fmod ( real x,
real y )

Compute the modulus of two real numbers.

Parameters
xThe dividend
yThe divisor
Returns
The modulus (remainder) of x and y

Definition at line 22 of file math.c.

◆ math_jac_rpz2xyz()

void math_jac_rpz2xyz ( real * rpz,
real * xyz,
real r,
real phi )

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].

Parameters
rpzinput rpz coordinates in a 3-length array
xyzoutput xyz coordinates in a 3-length array
rr coordinate of the gradient
phiphi coordinate of the gradient

Definition at line 44 of file math.c.

◆ math_jac_xyz2rpz()

void math_jac_xyz2rpz ( real * xyz,
real * rpz,
real r,
real phi )

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].

Parameters
xyzoutput xyz coordinates in a 3-length array
rpzinput rpz coordinates in a 3-length array
rr coordinate of the gradient
phiphi coordinate of the gradient

Definition at line 103 of file math.c.

◆ math_matmul()

void math_matmul ( real * matA,
real * matB,
int d1,
int d2,
int d3,
real * matC )

Matrix multiplication.

This function performs matrix multiplication. The matrices are given as arrays, the indexing being row-wise.

Parameters
matAinput array representing a d1 x d2 matrix
matBinput array representing a d2 x d3 matrix
d1input scalar dimension (columns in matA and matC)
d2input scalar dimension (rows in matA and columns in matB)
d3input scalar dimension (rows in matB and matC)
matCoutput array representing a d1 x d3 matrix

Definition at line 159 of file math.c.

◆ math_normal_rand()

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

Todo

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

Definition at line 188 of file math.c.

◆ math_ipow()

int math_ipow ( int a,
int p )

Calculate a^p where both a and p are integers (p >= 0)

Parameters
aargument
ppower
Returns
a^b

Definition at line 210 of file math.c.

◆ math_simpson()

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

Parameters
fpointer to the integrand function (of one variable)
alower limit for the interval
bupper limit for the interval
epsabsolute tolerance
Bug
Probably does not work with SIMD

Definition at line 232 of file math.c.

◆ math_linspace()

void math_linspace ( real * vec,
real a,
real b,
int n )

Generate linearly spaced vector.

Generates linearly space vector with n elements whose end points are a and b. If n = 1, return b.

Parameters
vecn length array where result is stored
astart point
bend point
nnumber of elements

Definition at line 251 of file math.c.

◆ math_uniquecount()

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.

Parameters
ininput array of length n
uniquearray where unique numbers are stored, must be of length n
countcount[i] is how many times unique[i] appears in array in, must be of length n
nlength of the argument arrays

Definition at line 276 of file math.c.

◆ math_point_on_plane()

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.

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

Parameters
qxyz coordinates of a query point
t1xyz coordinates of the first point defining the plane
t2xyz coordinates of the second point defining the plane
t3xyz coordinates of the third point defining the plane

Definition at line 324 of file math.c.

◆ math_barycentric_coords_triangle()

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

Parameters
APpoint position vector [x,y,z]
ABtriangle edge AB vector
ACtriangle edge AC vector
nnormal vector [nx, ny, nz] of the triangle
sevaluated s coordinate
tevaluated t coordinate

Definition at line 360 of file math.c.

◆ math_rsearch()

real * math_rsearch ( const real key,
const real * base,
int num )

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).

Parameters
keyreal value that serves as key for the search
basepointer to the first object of the array to be searched
numnumber of elements in array
Returns
pointer to array element preceding the key, or NULL if search failed

Definition at line 405 of file math.c.

◆ math_point_in_polygon()

int math_point_in_polygon ( real r,
real z,
real * rv,
real * zv,
int n )

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

Parameters
rr coordinate
zz coordinate
rvarray of r points for the polygon
zvarray of z points for the polygon
nnumber of points in the polygon

Definition at line 427 of file math.c.