+ * The {@code function} package contains function objects that wrap the + * methods contained in {@link java.lang.Math}, as well as common + * mathematical functions such as the gaussian and sinc functions. + *

+ * Parent package for common numerical analysis procedures, including root finding, + * function interpolation and integration. Note that the optimization (i.e. minimization + * and maximization) is a huge separate top package, despite it also operate on functions + * as defined by this top-level package. + *

+ *

+ * Functions interfaces are intended to be implemented by user code to represent their + * domain problems. The algorithms provided by the library will then operate on these + * function to find their roots, or integrate them, or ... Functions can be multivariate + * or univariate, real vectorial or matrix valued, and they can be differentiable or not. + *

Another floating point class. This one is built using radix 10000 + *which is 104, so its almost decimal.

+ * + *

The design goals here are: + *

+ *
1. Decimal math, or close to it
2. + *
3. Settable precision (but no mix between numbers using different settings)
4. + *
5. Portability. Code should be keep as portable as possible.
6. + *
7. Performance
8. + *
9. Accuracy - Results should always be +/- 1 ULP for basic + * algebraic operation
10. + *
11. Comply with IEEE 854-1987 as much as possible. + * (See IEEE 854-1987 notes below)
12. + *

+ * + *

+ *
1. Memory foot print. I'm using more memory than necessary to + * represent numbers to get better performance.
2. + *
3. Digits are bigger, so rounding is a greater loss. So, if you + * really need 12 decimal digits, better use 4 base 10000 digits + * there can be one partially filled.
4. + *

+ * + *

Numbers are represented in the following form: + *

```+ *n  =  sign × mant × (radix)exp;
+ *```
+ *where sign is ±1, mantissa represents a fractional number between + *zero and one. mant[0] is the least significant digit. + *exp is in the range of -32767 to 32768

+ * + *

IEEE 854-1987 Notes and differences

+ * + *

IEEE 854 requires the radix to be either 2 or 10. The radix here is + *10000, so that requirement is not met, but it is possible that a + *subclassed can be made to make it behave as a radix 10 + *number. It is my opinion that if it looks and behaves as a radix + *10 number then it is one and that requirement would be met.

+ * + *

The radix of 10000 was chosen because it should be faster to operate + *on 4 decimal digits at once instead of one at a time. Radix 10 behavior + *can be realized by add an additional rounding step to ensure that + *the number of decimal digits represented is constant.

+ * + *

The IEEE standard specifically leaves out internal data encoding, + *so it is reasonable to conclude that such a subclass of this radix + *10000 system is merely an encoding of a radix 10 system.

+ * + *

IEEE 854 also specifies the existence of "sub-normal" numbers. This + *class does not contain any such entities. The most significant radix + *10000 digit is always non-zero. Instead, we support "gradual underflow" + *by raising the underflow flag for numbers less with exponent less than + *expMin, but don't flush to zero until the exponent reaches MIN_EXP-digits. + *Thus the smallest number we can represent would be: + *1E(-(MIN_EXP-digits-1)∗4), eg, for digits=5, MIN_EXP=-32767, that would + *be 1e-131092.

+ * + *

IEEE 854 defines that the implied radix point lies just to the right + *of the most significant digit and to the left of the remaining digits. + *This implementation puts the implied radix point to the left of all + *digits including the most significant one. The most significant digit + *here is the one just to the right of the radix point. This is a fine + *detail and is really only a matter of definition. Any side effects of + *this can be rendered invisible by a subclass.

+ *This package provides Genetic Algorithms components and implementations. + *

+ *This package provides basic 1D geometry components. + *

+ *This package provides basic 3D geometry components. + *

+ *This package provides basic 2D geometry components. + *

+ *This package is the top level package for geometry. It provides only a few interfaces + *related to vectorial/affine spaces that are implemented in sub-packages. + *

+ *{@link org.apache.commons.math.geometry.partitioning.BSPTree BSP trees} + *are an efficient way to represent parts of space and in particular + *polytopes (line segments in 1D, polygons in 2D and polyhedrons in 3D) + *and to operate on them. The main principle is to recursively subdivide + *the space using simple hyperplanes (points in 1D, lines in 2D, planes + *in 3D). + *

+ * + *

+ *We start with a tree composed of a single node without any cut + *hyperplane: it represents the complete space, which is a convex + *part. If we add a cut hyperplane to this node, this represents a + *partition with the hyperplane at the node level and two half spaces at + *each side of the cut hyperplane. These half-spaces are represented by + *two child nodes without any cut hyperplanes associated, the plus child + *which represents the half space on the plus side of the cut hyperplane + *and the minus child on the other side. Continuing the subdivisions, we + *end up with a tree having internal nodes that are associated with a + *cut hyperplane and leaf nodes without any hyperplane which correspond + *to convex parts. + *

+ * + *

+ *When BSP trees are used to represent polytopes, the convex parts are + *known to be completely inside or outside the polytope as long as there + *is no facet in the part (which is obviously the case if the cut + *hyperplanes have been chosen as the underlying hyperplanes of the + *facets (this is called an autopartition) and if the subdivision + *process has been continued until all facets have been processed. It is + *important to note that the polytope is not defined by a + *single part, but by several convex ones. This is the property that + *allows BSP-trees to represent non-convex polytopes despites all parts + *are convex. The {@link + *org.apache.commons.math.geometry.partitioning.Region Region} class is + *devoted to this representation, it is build on top of the {@link + *org.apache.commons.math.geometry.partitioning.BSPTree BSPTree} class using + *boolean objects as the leaf nodes attributes to represent the + *inside/outside property of each leaf part, and also adds various + *methods dealing with boundaries (i.e. the separation between the + *inside and the outside parts). + *

+ * + *

+ *Rather than simply associating the internal nodes with an hyperplane, + *we consider sub-hyperplanes which correspond to the part of + *the hyperplane that is inside the convex part defined by all the + *parent nodes (this implies that the sub-hyperplane at root node is in + *fact a complete hyperplane, because there is no parent to bound + *it). Since the parts are convex, the sub-hyperplanes are convex, in + *3D the convex parts are convex polyhedrons, and the sub-hyperplanes + *are convex polygons that cut these polyhedrons in two + *sub-polyhedrons. Using this definition, a BSP tree completely + *partitions the space. Each point either belongs to one of the + *sub-hyperplanes in an internal node or belongs to one of the leaf + *convex parts. + *

+ * + *

+ *In order to determine where a point is, it is sufficient to check its + *position with respect to the root cut hyperplane, to select the + *corresponding child tree and to repeat the procedure recursively, + *until either the point appears to be exactly on one of the hyperplanes + *in the middle of the tree or to be in one of the leaf parts. For + *this operation, it is sufficient to consider the complete hyperplanes, + *there is no need to check the points with the boundary of the + *sub-hyperplanes, because this check has in fact already been realized + *by the recursive descent in the tree. This is very easy to do and very + *efficient, especially if the tree is well balanced (the cost is + *`O(log(n))` where `n` is the number of facets) + *or if the first tree levels close to the root discriminate large parts + *of the total space. + *

+ * + *

+ *One of the main sources for the development of this package was Bruce + *Naylor, John Amanatides and William Thibault paper Merging + *BSP Trees Yields Polyhedral Set Operations Proc. Siggraph '90, + *Computer Graphics 24(4), August 1990, pp 115-124, published by the + *Association for Computing Machinery (ACM). The same paper can also be + *found here. + *

+ *This package provides multidimensional ordering features for partitioning. + *

+ *This package provides classes to handle discrete events occurring during + *Ordinary Differential Equations integration. + *

+ * + *

+ *Discrete events detection is based on switching functions. The user provides + *a simple {@link org.apache.commons.math.ode.events.EventHandler#g g(t, y)} + *function depending on the current time and state. The integrator will monitor + *the value of the function throughout integration range and will trigger the + *event when its sign changes. The magnitude of the value is almost irrelevant, + *it should however be continuous (but not necessarily smooth) for the sake of + *root finding. The steps are shortened as needed to ensure the events occur + *at step boundaries (even if the integrator is a fixed-step integrator). + *

+ * + *

+ *When an event is triggered, several different options are available: + *

+ *
+ *
• integration can be stopped (this is called a G-stop facility),
• + *
• the state vector or the derivatives can be changed,
• + *
• or integration can simply go on.
• + *
+ * + *

+ *The first case, G-stop, is the most common one. A typical use case is when an + *ODE must be solved up to some target state is reached, with a known value of + *the state but an unknown occurrence time. As an example, if we want to monitor + *a chemical reaction up to some predefined concentration for the first substance, + *we can use the following switching function setting: + *

```+ *  public double g(double t, double[] y) {
+ *    return y[0] - targetConcentration;
+ *  }
+ *
+ *  public int eventOccurred(double t, double[] y) {
+ *    return STOP;
+ *  }
+ *```
+ *

+ * + *

+ *The second case, change state vector or derivatives is encountered when dealing + *with discontinuous dynamical models. A typical case would be the motion of a + *spacecraft when thrusters are fired for orbital maneuvers. The acceleration is + *smooth as long as no maneuver are performed, depending only on gravity, drag, + *third body attraction, radiation pressure. Firing a thruster introduces a + *discontinuity that must be handled appropriately by the integrator. In such a case, + *we would use a switching function setting similar to this: + *

```+ *  public double g(double t, double[] y) {
+ *    return (t - tManeuverStart) ∗ (t - tManeuverStop);
+ *  }
+ *
+ *  public int eventOccurred(double t, double[] y) {
+ *    return RESET_DERIVATIVES;
+ *  }
+ *```
+ *

+ * + *

+ *The third case is useful mainly for monitoring purposes, a simple example is: + *

```+ *  public double g(double t, double[] y) {
+ *    return y[0] - y[1];
+ *  }
+ *
+ *  public int eventOccurred(double t, double[] y) {
+ *    logger.log("y0(t) and y1(t) curves cross at t = " + t);
+ *    return CONTINUE;
+ *  }
+ *```
+ *

+ *This package provides classes to solve non-stiff Ordinary Differential Equations problems. + *

+ *This package provides classes to solve Ordinary Differential Equations problems. + *

+ * + *

+ *This package solves Initial Value Problems of the form + *`y'=f(t,y)` with `t0` and + *`y(t0)=y0` known. The provided + *integrators compute an estimate of `y(t)` from + *`t=t0` to `t=t1`. + *It is also possible to get thederivatives with respect to the initial state + *`dy(t)/dy(t0)` or the derivatives with + *respect to some ODE parameters `dy(t)/dp`. + *

+ * + *

+ *All integrators provide dense output. This means that besides + *computing the state vector at discrete times, they also provide a + *cheap mean to get the state between the time steps. They do so through + *classes extending the {@link + *org.apache.commons.math.ode.sampling.StepInterpolator StepInterpolator} + *abstract class, which are made available to the user at the end of + *each step. + *

+ * + *

+ *All integrators handle multiple discrete events detection based on switching + *functions. This means that the integrator can be driven by user specified + *discrete events. The steps are shortened as needed to ensure the events occur + *at step boundaries (even if the integrator is a fixed-step + *integrator). When the events are triggered, integration can be stopped + *(this is called a G-stop facility), the state vector can be changed, + *or integration can simply go on. The latter case is useful to handle + *discontinuities in the differential equations gracefully and get + *accurate dense output even close to the discontinuity. + *

+ * + *

+ *The user should describe his problem in his own classes + *(`UserProblem` in the diagram below) which should implement + *the {@link org.apache.commons.math.ode.FirstOrderDifferentialEquations + *FirstOrderDifferentialEquations} interface. Then he should pass it to + *the integrator he prefers among all the classes that implement the + *{@link org.apache.commons.math.ode.FirstOrderIntegrator + *FirstOrderIntegrator} interface. + *

+ * + *

+ *The solution of the integration problem is provided by two means. The + *first one is aimed towards simple use: the state vector at the end of + *the integration process is copied in the `y` array of the + *{@link org.apache.commons.math.ode.FirstOrderIntegrator#integrate + *FirstOrderIntegrator.integrate} method. The second one should be used + *when more in-depth information is needed throughout the integration + *process. The user can register an object implementing the {@link + *org.apache.commons.math.ode.sampling.StepHandler StepHandler} interface or a + *{@link org.apache.commons.math.ode.sampling.StepNormalizer StepNormalizer} + *object wrapping a user-specified object implementing the {@link + *org.apache.commons.math.ode.sampling.FixedStepHandler FixedStepHandler} + *interface into the integrator before calling the {@link + *org.apache.commons.math.ode.FirstOrderIntegrator#integrate + *FirstOrderIntegrator.integrate} method. The user object will be called + *appropriately during the integration process, allowing the user to + *process intermediate results. The default step handler does nothing. + *

+ * + *

+ *{@link org.apache.commons.math.ode.ContinuousOutputModel + *ContinuousOutputModel} is a special-purpose step handler that is able + *to store all steps and to provide transparent access to any + *intermediate result once the integration is over. An important feature + *of this class is that it implements the `Serializable` + *interface. This means that a complete continuous model of the + *integrated function throughout the integration range can be serialized + *and reused later (if stored into a persistent medium like a filesystem + *or a database) or elsewhere (if sent to another application). Only the + *result of the integration is stored, there is no reference to the + *integrated problem by itself. + *

+ * + *

+ *Other default implementations of the {@link + *org.apache.commons.math.ode.sampling.StepHandler StepHandler} interface are + *available for general needs ({@link + *org.apache.commons.math.ode.sampling.DummyStepHandler DummyStepHandler}, {@link + *org.apache.commons.math.ode.sampling.StepNormalizer StepNormalizer}) and custom + *implementations can be developed for specific needs. As an example, + *if an application is to be completely driven by the integration + *process, then most of the application code will be run inside a step + *handler specific to this application. + *

+ * + *

+ *Some integrators (the simple ones) use fixed steps that are set at + *creation time. The more efficient integrators use variable steps that + *are handled internally in order to control the integration error with + *respect to a specified accuracy (these integrators extend the {@link + *org.apache.commons.math.ode.nonstiff.AdaptiveStepsizeIntegrator + *AdaptiveStepsizeIntegrator} abstract class). In this case, the step + *handler which is called after each successful step shows up the + *variable stepsize. The {@link + *org.apache.commons.math.ode.sampling.StepNormalizer StepNormalizer} class can + *be used to convert the variable stepsize into a fixed stepsize that + *can be handled by classes implementing the {@link + *org.apache.commons.math.ode.sampling.FixedStepHandler FixedStepHandler} + *interface. Adaptive stepsize integrators can automatically compute the + *initial stepsize by themselves, however the user can specify it if he + *prefers to retain full control over the integration or if the + *automatic guess is wrong. + *

+ * + *

+ * + * + * + * + * + * + * + * + *
+ *

+ * + * + * + * + * + * + * + * + * + * + *
+ *

+ * + *

+ *In the table above, the {@link org.apache.commons.math.ode.nonstiff.AdamsBashforthIntegrator + *Adams-Bashforth} and {@link org.apache.commons.math.ode.nonstiff.AdamsMoultonIntegrator + *Adams-Moulton} integrators appear as variable-step ones. This is an experimental extension + *to the classical algorithms using the Nordsieck vector representation. + *

+ *This package provides classes to handle sampling steps during + *Ordinary Differential Equations integration. + *

+ * + *

+ *In addition to computing the evolution of the state vector at some grid points, all + *ODE integrators also build up interpolation models of this evolution inside the + *last computed step. If users are interested in these interpolators, they can register a + *{@link org.apache.commons.math.ode.sampling.StepHandler StepHandler} instance using the + *{@link org.apache.commons.math.ode.FirstOrderIntegrator#addStepHandler addStepHandler} + *method which is supported by all integrators. The integrator will call this instance + *at the end of each accepted step and provide it the interpolator. The user can do + *whatever he wants with this interpolator, which computes both the state and its + *time-derivative. A typical use of step handler is to provide some output to monitor + *the integration process. + *

+ * + *

+ *In a sense, this is a kind of Inversion Of Control: rather than having the master + *application driving the slave integrator by providing the target end value for + *the free variable, we get a master integrator scheduling the free variable + *evolution and calling the slave application callbacks that were registered at + *configuration time. + *

+ * + *

+ *Since some integrators may use variable step size, the generic {@link + *org.apache.commons.math.ode.sampling.StepHandler StepHandler} interface can be called + *either at regular or irregular rate. This interface allows to navigate to any location + *within the last computed step, thanks to the provided {@link + *org.apache.commons.math.ode.sampling.StepInterpolator StepInterpolator} object. + *If regular output is desired (for example in order to write an ephemeris file), then + *the simpler {@link org.apache.commons.math.ode.sampling.FixedStepHandler FixedStepHandler} + *interface can be used. Objects implementing this interface should be wrapped within a + *{@link org.apache.commons.math.ode.sampling.StepNormalizer StepNormalizer} instance + *in order to be registered to the integrator. + *

Curve fitting is a special case of a least squares problem + *were the parameters are the coefficients of a function `f` + *whose graph `y=f(x)` should pass through sample points, and + *were the objective function is the squared sum of residuals + *`f(xi)-yi` for observed points + *(xi, yi).