![]() |
Tesseroids 1.0: User Manual and API Documentation |
Functions for implementing a Gauss-Legendre Quadrature numerical integration. More...
#include <stdlib.h>
#include <math.h>
#include "constants.h"
#include "utils.h"
#include <stdio.h>
#include "logger.h"
Defines | |
#define | GLQ_ABS(x) ((x) < 0 ? -1*(x) : (x)) |
Functions | |
GLQ * | glq_new (int order, double lower, double upper) |
Make a new GLQ structure and set all the parameters needed. | |
void | glq_free (GLQ *glq) |
Free the memory allocated to make a GLQ structure. | |
int | glq_nodes (int order, double *nodes) |
Calculates the GLQ nodes using glq_next_root. | |
int | glq_set_limits (double lower, double upper, GLQ *glq) |
Put the GLQ nodes to the integration limits IN PLACE. | |
int | glq_next_root (double initial, int root_index, int order, double *roots) |
Calculate the next Legendre polynomial root given the previous root found. | |
int | glq_weights (int order, double *nodes, double *weights) |
Calculates the weighting coefficients for the GLQ integration. | |
Variables | |
const int | GLQ_MAXIT = 1000 |
Max iterations of the root-finder algorithm. | |
const double | GLQ_MAXERROR = 0.000000000000001 |
Max error allowed for the root-finder algorithm. |
Functions for implementing a Gauss-Legendre Quadrature numerical integration.
void glq_free | ( | GLQ * | glq | ) |
Free the memory allocated to make a GLQ structure.
glq | pointer to the allocated memory |
GLQ* glq_new | ( | int | order, | |
double | lower, | |||
double | upper | |||
) |
Make a new GLQ structure and set all the parameters needed.
WARNING: Don't forget to free the memory malloced by this function using glq_free()!
Prints error and warning messages using the logging.h module.
order | order of the quadrature, ie number of nodes | |
lower | lower integration limit | |
upper | upper integration limit |
int glq_next_root | ( | double | initial, | |
int | root_index, | |||
int | order, | |||
double * | roots | |||
) |
Calculate the next Legendre polynomial root given the previous root found.
Uses the root-finder algorithm of:
Barrera-Figueroa, V., Sosa-Pedroza, J. and López-Bonilla, J., 2006, "Multiple root finder algorithm for Legendre and Chebyshev polynomials via Newton's method", 2006, Annales mathematicae et Informaticae, 33, pp 3-13
initial | initial estimate of the next root. I recommend the use of ![]() ![]() | |
root_index | index of the desired root, starting from 0 | |
order | order of the Legendre polynomial, ie number of roots. | |
roots | array with the roots found so far. Will return the next root in roots[root_index], so make sure to malloc enough space. |
Compute the absolute value of x
int glq_nodes | ( | int | order, | |
double * | nodes | |||
) |
Calculates the GLQ nodes using glq_next_root.
Nodes will be in the [-1,1] interval. To convert them to the integration limits use glq_scale_nodes
order | order of the quadrature, ie how many nodes. Must be >= 2. | |
nodes | pre-allocated array to return the nodes. |
int glq_set_limits | ( | double | lower, | |
double | upper, | |||
GLQ * | glq | |||
) |
Put the GLQ nodes to the integration limits IN PLACE.
Will replace the values of glq.nodes with ones in the specified integration limits.
In case the GLQ structure was created with glq_new(), the integration limits can be reset using this function.
lower | lower integration limit | |
upper | upper integration limit | |
glq | pointer to a GLQ structure created with glq_new() and with all necessary memory allocated |
int glq_weights | ( | int | order, | |
double * | nodes, | |||
double * | weights | |||
) |
Calculates the weighting coefficients for the GLQ integration.
order | order of the quadrature, ie number of nodes and weights. | |
nodes | array containing the GLQ nodes calculated by glq_nodes. IMPORTANT: needs the nodes in [-1,1] interval! Scaled nodes will result in wrong weights! | |
weights | pre-allocated array to return the weights |