Patterns in static

Apophenia

Macros | Functions
Things to make life easier with the GSL

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)
 

Detailed Description

Macro Definition Documentation

#define apop_sum (   in)    apop_vector_sum(in)

An alias for apop_vector_sum. Returns the sum of the data in the given vector.

Function Documentation

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:

1 gsl_matrix *a_copy = apop_matrix_copy(original);
Parameters
inthe input data
Returns
a structure that this function will allocate and fill. If 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.

Parameters
datathe 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.

  • If you want sample variance, multiply the result returned by (data->size1*data->size2)/(data->size1*data->size2-1.0).
Parameters
datathe matrix to be averaged.
meanwhere to put the mean to be calculated.
varwhere to put the variance to be calculated.
  • If 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.

Parameters
m1the upper/rightmost matrix (default=NULL, in which case this basically copies m2)
m2the second matrix (default = NULL, in which case m1 is returned)
posnif 'r', stack rows on top of other rows, else, e.g. 'c' stack columns next to columns. (default ='r')
inplaceIf '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')
Returns
the stacked data, either in a new matrix or a pointer to 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.

1 gsl_matrix *apop_stack_two_by_two(gsl_matrix *ul, gsl_matrix *ur, gsl_matrix *dl, gsl_matrix *dr){
2  gsl_matrix *output, *t;
3  output = apop_matrix_stack(ul, ur, 'c');
4  t = apop_matrix_stack(dl, dr, 'c');
5  apop_matrix_stack(output, t, 'r', .inplace='y');
6  gsl_matrix_free(t);
7  return output;
8 }
long double apop_matrix_sum ( const gsl_matrix *  m)

Returns the sum of the elements of a matrix. Occasionally convenient.

Parameters
mthe matrix to be summed.
gsl_rng* apop_rng_alloc ( int  seed)

Initialize a gsl_rng.

Uses the Tausworth routine.

Parameters
seedThe seed. No need to get funny with it: 0, 1, and 2 will produce wholly different streams.
Returns
The RNG ready for your use.
  • If you are confident that your code is debugged and would like a new stream of values every time your program runs (provided your runs are more than a second apart), seed with the time:
#include <apop.h>
#include <time.h>
int main(){
apop_opts.rng_seed = time(NULL);
.count=10,
)
);
}
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.

Parameters
inA gsl_vector
maxAn upper and lower bound to the elements of the vector. (default: GSL_POSINF)
Returns
1 if everything is bounded: not Inf, -Inf, or NaN, and $-\max < x < \max$; zero otherwise.
  • A NULL vector has no unbounded elements, so NULL input returns 1. You get a warning if apop_opts.verbosity >=1.
  • This function uses the Designated initializers syntax for inputs.
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:

1 gsl_vector *a_copy = apop_vector_copy(original);
Parameters
inthe input data
Returns
a structure that this function will allocate and fill. If 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:

  • 'e' or 'E' (the default): scalar distance (standard Euclidean metric) between two vectors. Simply $\sqrt{\sum_i{(a_i - b_i)^2}},$ where $i$ iterates over dimensions.
  • 'm' or 'M' Returns the Manhattan metric distance between two vectors: $\sum_i{|a_i - b_i|},$ where $i$ iterates over dimensions.
  • 'd' or 'D' The discrete norm: if $a = b$, return zero, else return one.
  • 's' or 'S' The sup norm: find the dimension where $|a_i - b_i|$ is largest, return the distance along that one dimension.
  • 'l' or 'L' The $L_p$ norm, $\left(\sum_i{(a_i - b_i)^2}\right)^{1/p},$. The value of $p$ is set by the fourth (optional) argument.
Parameters
inaFirst vector (No default, must not be NULL)
inbSecond vector (Default = zero)
metricThe type of metric, as above.
normIf you are using an $L_p$ norm, this is $p$. Must be strictly greater than zero. (default = 2)

Notice that the defaults are such that

1 apop_vector_distance(v);
2 apop_vector_distance(v, .metric = 's');

gives you the standard Euclidean length of v and its longest element.

#include <apop.h>
/* Test distance calculations using a 3-4-5 triangle */
int main(){
gsl_vector *v1 = gsl_vector_alloc(2);
gsl_vector *v2 = gsl_vector_alloc(2);
gsl_vector_set(v1, 0,2);
gsl_vector_set(v1, 1,2);
gsl_vector_set(v2, 0,5);
gsl_vector_set(v2, 1,6);
assert(apop_vector_distance(v1, v1, 'd') == 0); //discrete: if vectors are equal d==0;
assert(apop_vector_distance(v1, v2, 'd') == 1); // if vectors differ d ==1
assert(apop_vector_distance(v1, NULL, 'm') == 4.); //length of v1, Manhattan metric
assert(apop_vector_distance(v1,v2) == 5.); //the hypotenuse of the 3-4-5 triangle
assert(apop_vector_distance(v2, NULL, 's') == 6); //length of v2, sup norm
assert(apop_vector_distance(v1,v2, 'm') == 7.); //distance v1 to v2 in Manhattan
assert(apop_vector_distance(v1,v2, 'L', 2) == 5.); //L_2 norm == Euclidean
}
void apop_vector_exp ( gsl_vector *  v)

Replace every vector element $v_i$ with exp( $v_i$).

  • If the input vector is NULL, do nothing.
void apop_vector_log ( gsl_vector *  v)

Take the natural log of every element in a vector.

  • If the input vector is NULL, do nothing.
void apop_vector_log10 ( gsl_vector *  v)

Take the log (base ten) of every element in a vector.

  • If the input vector is 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.

Parameters
v1the upper vector (default=NULL, in which case this basically copies v2)
v2the second vector (default=NULL, in which case nothing is added)
inplaceIf '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')
Returns
the stacked data, either in a new vector or a pointer to v1.
long double apop_vector_sum ( const gsl_vector *  in)

Returns the sum of the data in the given vector.

Autogenerated by doxygen on Sun Oct 26 2014 (Debian 0.999b+ds3-2).