![]() |
|
Macros | |
#define | apop_sum(in) apop_vector_sum(in) |
Functions | |
gsl_rng * | apop_rng_alloc (int seed) |
gsl_vector * | apop_vector_copy (const gsl_vector *in) |
gsl_matrix * | apop_matrix_copy (const gsl_matrix *in) |
void | apop_vector_log10 (gsl_vector *v) |
void | apop_vector_log (gsl_vector *v) |
void | apop_vector_exp (gsl_vector *v) |
gsl_vector * | apop_vector_stack (gsl_vector *v1, gsl_vector *v2, char inplace) |
gsl_matrix * | apop_matrix_stack (gsl_matrix *m1, gsl_matrix *m2, char posn, char inplace) |
int | apop_vector_bounded (const gsl_vector *in, long double max) |
long double | apop_vector_sum (const gsl_vector *in) |
double | apop_vector_distance (const gsl_vector *ina, const gsl_vector *inb, const char metric, const double norm) |
long double | apop_matrix_sum (const gsl_matrix *m) |
double | apop_matrix_mean (const gsl_matrix *data) |
void | apop_matrix_mean_and_var (const gsl_matrix *data, double *mean, double *var) |
#define apop_sum | ( | in | ) | apop_vector_sum(in) |
An alias for apop_vector_sum. Returns the sum of the data in the given vector.
gsl_matrix* apop_matrix_copy | ( | const gsl_matrix * | in | ) |
Copy one gsl_matrix
to another. That is, all data is duplicated. Unlike gsl_matrix_memcpy
, this function allocates and returns the destination, so you can use it like this:
in | the input data |
gsl_matrix_alloc
fails, returns NULL
. double apop_matrix_mean | ( | const gsl_matrix * | data | ) |
Returns the mean of all elements of a matrix.
Calculated with an eye toward avoiding overflow errors.
data | the matrix to be averaged. If NULL , return zero. |
void apop_matrix_mean_and_var | ( | const gsl_matrix * | data, |
double * | mean, | ||
double * | var | ||
) |
Returns the mean and population variance of all elements of a matrix.
(data->size1*data->size2)/(data->size1*data->size2-1.0)
.data | the matrix to be averaged. |
mean | where to put the mean to be calculated. |
var | where to put the variance to be calculated. |
NULL
, return (zero, NaN). gsl_matrix* apop_matrix_stack | ( | gsl_matrix * | m1, |
gsl_matrix * | m2, | ||
char | posn, | ||
char | inplace | ||
) |
Put the first matrix either on top of or to the right of the second matrix. The fn returns a new matrix, meaning that at the end of this function, until you gsl_matrix_free() the original matrices, you will be taking up twice as much memory. Plan accordingly.
m1 | the upper/rightmost matrix (default=NULL , in which case this basically copies m2 ) |
m2 | the second matrix (default = NULL , in which case m1 is returned) |
posn | if 'r', stack rows on top of other rows, else, e.g. 'c' stack columns next to columns. (default ='r') |
inplace | If 'y', use apop_matrix_realloc to modify m1 in place; see the caveats on that function. Otherwise, allocate a new matrix, leaving m1 unmolested. (default='n') |
m1
.For example, here is a little function to merge four matrices into a single two-part-by-two-part matrix. The original matrices are unchanged.
long double apop_matrix_sum | ( | const gsl_matrix * | m | ) |
Returns the sum of the elements of a matrix. Occasionally convenient.
m | the matrix to be summed. |
gsl_rng* apop_rng_alloc | ( | int | seed | ) |
Initialize a gsl_rng
.
Uses the Tausworth routine.
seed | The seed. No need to get funny with it: 0, 1, and 2 will produce wholly different streams. |
int apop_vector_bounded | ( | const gsl_vector * | in, |
long double | max | ||
) |
Test for a situation when a vector is diverging, so you can preempt a procedure that is about to break on infinite values.
Alternatively, set max
to INFINITY
(or GSL_INF
) to just test whether all of the matrix's elements are finite.
in | A gsl_vector |
max | An upper and lower bound to the elements of the vector. (default: GSL_POSINF) |
NULL
vector has no unbounded elements, so NULL
input returns 1. You get a warning if apop_opts.verbosity >=1
. gsl_vector* apop_vector_copy | ( | const gsl_vector * | in | ) |
Copy one gsl_vector
to another. That is, all data is duplicated. Unlike gsl_vector_memcpy
, this function allocates and returns the destination, so you can use it like this:
in | the input data |
gsl_vector_alloc
fails, returns NULL
. double apop_vector_distance | ( | const gsl_vector * | ina, |
const gsl_vector * | inb, | ||
const char | metric, | ||
const double | norm | ||
) |
Returns the distance between two vectors, where distance is defined based on the third (optional) parameter:
ina | First vector (No default, must not be NULL ) |
inb | Second vector (Default = zero) |
metric | The type of metric, as above. |
norm | If you are using an ![]() ![]() |
Notice that the defaults are such that
gives you the standard Euclidean length of v
and its longest element.
void apop_vector_exp | ( | gsl_vector * | v | ) |
Replace every vector element with exp(
).
NULL
, do nothing. void apop_vector_log | ( | gsl_vector * | v | ) |
Take the natural log of every element in a vector.
NULL
, do nothing. void apop_vector_log10 | ( | gsl_vector * | v | ) |
Take the log (base ten) of every element in a vector.
NULL
, do nothing. gsl_vector* apop_vector_stack | ( | gsl_vector * | v1, |
gsl_vector * | v2, | ||
char | inplace | ||
) |
Put the first vector on top of the second vector.
v1 | the upper vector (default=NULL , in which case this basically copies v2 ) |
v2 | the second vector (default=NULL , in which case nothing is added) |
inplace | If 'y', use apop_vector_realloc to modify v1 in place; see the caveats on that function. Otherwise, allocate a new vector, leaving v1 unmolested. (default='n') |
v1
.long double apop_vector_sum | ( | const gsl_vector * | in | ) |
Returns the sum of the data in the given vector.