1 Introduction

1.1 What is ALGLIB

ALGLIB is a cross-platform numerical analysis and data mining library. It supports several programming languages (C++, C#, Delphi, VB.NET, Python) and several operating systems (Windows, *nix family).

ALGLIB features include:

Data analysis (classification/regression, including neural networks)
Optimization and nonlinear solvers
Interpolation and linear/nonlinear least-squares fitting
Linear algebra (direct algorithms, EVD/SVD), direct and iterative linear solvers, Fast Fourier Transform and many other algorithms (numerical integration, ODEs, statistics, special functions)

ALGLIB Project (the company behind ALGLIB) delivers to you several editions of ALGLIB:

ALGLIB Free Edition - full functionality but limited performance and license
ALGLIB Commercial Edition - high-performance version of ALGLIB with business-friendly license

Free Edition is a serial version without multithreading support or extensive low-level optimizations (generic C or C# code). Commercial Edition is a heavily optimized version of ALGLIB. It supports multithreading, it was extensively optimized, and (on Intel platforms) - our commercial users may enjoy precompiled version of ALGLIB which internally calls Intel MKL to accelerate low-level tasks. We obtained license from Intel corp., which allows us to integrate Intel MKL into ALGLIB, so you don't have to buy separate license from Intel.

1.2 ALGLIB license

ALGLIB Free Edition is distributed under license which favors non-commmercial usage, but is not well suited for commercial applications:

ALGLIB for C++ and ALGLIB for C# are distributed under GPL 2+, GPL license version 2 or at your option any later version. A copy of the GNU General Public License is available at http://www.fsf.org/licensing/licenses
ALGLIB for Delphi and ALGLIB for CPython are distributed under ALGLIB Personal and Academic Use License Agreement.

ALGLIB Commercial Edition is distributed under license which is friendly to commericial users. A copy of the commercial license can be found at http://www.alglib.net/commercial.php.

1.3 Documentation license

This reference manual is licensed under BSD-like documentation license:

Redistribution and use of this document (ALGLIB Reference Manual) with or without modification, are permitted provided that such redistributions will retain the above copyright notice, this condition and the following disclaimer as the first (or last) lines of this file.

THIS DOCUMENTATION IS PROVIDED BY THE ALGLIB PROJECT "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE ALGLIB PROJECT BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS DOCUMENTATION, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

1.4 Reference Manual and User Guide

ALGLIB Project provides two sources of information: ALGLIB Reference Manual (this document) and ALGLIB User Guide.

ALGLIB Reference Manual contains full description of all publicly accessible ALGLIB units accompanied with examples. Reference Manual is focused on the source code: it documents units, functions, structures and so on. If you want to know what unit YYY can do or what subroutines unit ZZZ contains Reference Manual is a place to go. Free software needs free documentation - that's why ALGLIB Reference Manual is licensed under BSD-like documentation license.

Additionally to the Reference Manual we provide you User Guide. User Guide is focused on more general questions: how fast ALGLIB is? how reliable it is? what are the strong and weak sides of the algorithms used? We aim to make ALGLIB User Guide an important source of information both about ALGLIB and numerical analysis algorithms in general. We want it to be a book about algorithms, not just software documentation. And we want it to be unique - that's why ALGLIB User Guide is distributed under less-permissive personal-use-only license.

1.5 Acknowledgements

ALGLIB was not possible without contribution of the next open source projects:

We also want to thank developers of the Intel's local development center (Nizhny Novgorod branch) for their help during MKL integration.

2 ALGLIB structure

2.1 Packages

ALGLIB is a C++ interface to the computational core written in C. Both C library and C++ wrapper are automatically generated by code generation tools developed within ALGLIB project. Pre-3.0 versions of ALGLIB included more than 100 units, but it was difficult to work with such large number of files. Since ALGLIB 3.0 all units are merged into 11 packages and two support units:

alglibmisc.cpp - contains different algorithms which are hard to classify
dataanalysis.cpp - contains data mining algorithms
diffequations.cpp - contains differential equation solvers
fasttransforms.cpp - contains FFT and other related algorithms
integration.cpp - contains numerical integration algorithms
interpolation.cpp - contains interpolation algorithms
linalg.cpp - contains linear algebra algorithms
optimization.cpp - contains optimization algorithms
solvers.cpp - contains linear and nonlinear solvers
specialfunctions.cpp - contains special functions
statistics.cpp - statistics
alglibinternal.cpp - contains internal functions which are used by other packages, but not exposed to the external world
ap.cpp - contains publicly accessible vector/matrix classes, most important and general functions and other "basic" functionality

One package may rely on other ones, but we have tried to reduce number of dependencies. Every package relies on ap.cpp and many packages rely on alglibinternal.cpp. But many packages require only these two to work, and many other packages need significantly less than 13 packages. For example, statistics.cpp requires two packages mentioned above and only one additional package - specialfunctions.cpp.

2.2 Subpackages

There is one more concept to learn - subpackages. Every package was created from several source files. For example (as of ALGLIB 3.0.0), linalg.cpp was created by merging together 14 .cpp files (C++ interface) and 14 .c files (computational core). These files provide different functionality: one of them calculates triangular factorizations, another generates random matrices, and so on. We've merged source code, but what to do with their documentation?

Of course, we can merge their documentation (as we've merged units) in one big list of functions and data structures, but such list will be hard to read. Instead, we have decided to merge source code, but leave documentation separate.

If you look at the list of ALGLIB packages, you will see that each package includes several subpackages. For example, linalg.cpp includes trfac, svd, evd and other subpackages. These subpackages do no exist as separate files, namespaces or other entities. They are just subsets of one large unit which provide significantly different functionality. They have separate documentation sections, but if you want to use svd subpackage, you have to include linalg.h, not svd.h.

2.3 Open Source and Commercial versions

ALGLIB comes in two versions - open source (GPL-licensed) and commercial (closed source) one. Both versions have same functionality, i.e. may solve same set of problems. However, commercial version differs from open source one in following aspects:

License. Commercial ALGLIB is licensed under non-copyleft license which is friendly to commercial users.
Performance. Many algorithms in commercial ALGLIB are multi-threaded and SIMD-optimized (when used on Intel systems). Open source ALGLIB is single-threaded and can not efficiently use modern multicore CPU's.
You have to study comments on specific functions if you want to know whether they have multithreaded versions or not.

This documentation applies to both versions of ALGLIB. Detailed description of commercial version can be found below.

3 Compatibility

3.1 CPU

ALGLIB is compatible with any CPU which:

supports double precision arithmetics
complies with IEEE 754 floating point standard (especially in its handling of IEEE special values)
either big-endian or little-endian (but not mixed-endian)

Most mainstream CPU's (in particular, x86, x86_64, ARM and SPARC) satisfy these requirements.

As for Intel architecture, ALGLIB works with both FPU-based and SIMD-based implementations of floating point math. 80-bit internal representation used by Intel FPU is not a problem for ALGLIB.

3.2 OS

ALGLIB for C++ (both open source and commercial versions) can be compiled in OS-agnostic mode (no OS-specific preprocessor definitions), when it is compatible with any OS which supports C++98 standard library. In particular, it will work under any POSIX-compatible OS and under Windows.

If you want to use multithreaded capabilities of commercial version of ALGLIB, you should compile it in OS-specific mode by #defining either AE_OS=AE_WINDOWS or AE_OS=AE_POSIX at compile time, depending on OS being used. Former corresponds to any modern OS (32/64-bit Windows XP and later) from Windows family, while latter means almost any POSIX-compatible OS. It applies only to commercial version of ALGLIB. Open source version is always OS-agnostic, even in the presence of OS-specific definitions.

3.3 Compiler

ALGLIB is compatible with any C++ compiler which:

supports 32-bit and 64-bit signed integer datatypes
emits code which handles comparisons with IEEE special values without raising exception. We don't require that x/0 will return INF. But at least we must be able to compare double precision value with infinity or NAN without raising exception.

All modern compilers satisfy these requirements.

However, some very old compilers (ten years old version of Borland C++ Builder, for example) may emit code which does not correctly work with IEEE special values. If you use one of these old compilers, we recommend you to run ALGLIB test suite to ensure that library works.

3.4 Optimization settings

ALGLIB is compatible with any kind of optimizing compiler as long as:

volatile modifier is correctly handled (i.e. compiler does not optimize volatile reads/writes)
optimized code correctly handles IEEE special values

Generally, all kinds of optimization that were marked by compiler vendor as "safe" are possible. For example, ALGLIB can be compiled:

under MSVC: with /O1, /O2, /Og, /Os, /Ox, /Ot, /Oy, /fp:precise, /fp:except, /fp:strict
under GCC: with -O1, -O2, -O3, -Os

From the other side, following "unsafe" optimizations will break ALGLIB:

under MSVC: /fp:fast
under GCC: -Ofast, -ffast-math

4 Compiling ALGLIB

4.1 Adding to your project

Adding ALGLIB to your project is easy - just pick packages you need and... add them to your project! Under most used compilers (GCC, MSVC) it will work without any additional settings. In other cases you will need to define several preprocessor definitions (this topic will be detailed below), but everything will still be simple.

By "adding to your project" we mean that you should a) compile .cpp files with the rest of your project, and b) include .h files you need. Do not include .cpp files - these files must be compiled separately, not as part of some larger source file. The only files you should include are .h files, stored in the /src folder of the ALGLIB distribution.

As you see, ALGLIB has no project files or makefiles. Why? There are several reasons:

first, because many ALGLIB users don't need separate static library (which will be created by invoking makefile) - they prefer to integrate source code in their projects. We have provided script-based build system before, but majority of our users prefer to build ALGLIB themselves.
second, because we want ALGLIB to be usable in any programming environment, whether it is Visual Studio, GNU build system or something else. The best solution is to write package which doesn't depend on any particular programming environment.

In any case, compiling ALGLIB is so simple that even without project file you can do it in several minutes.

4.2 Configuring for your compiler

If you use modern versions of MSVC or GCC, you don't need to configure ALGLIB at all. But if you use outdated versions of these compilers (or something else), then you may need to tune definitions of several data types:

alglib_impl::ae_int32_t - signed integer which is 32 bits wide
alglib_impl::ae_int64_t - signed integer which is 64 bits wide
alglib_impl::ae_uint64_t - unsigned integer which is 64 bits wide
alglib_impl::ae_int_t - signed integer which has same width as pointer

ALGLIB tries to autodetect your compiler and to define these types in compiler-specific manner:

ae_int32_t is defined as int, because this type is 32 bits wide in all modern compilers.
ae_int64_t is defined as _int64 (MSVC) or as signed long long (GCC, Sun Studio).
ae_uint64_t is defined as unsigned _int64 (MSVC) or as unsigned long long (GCC, Sun Studio).
ae_int_t is defined as ptrdiff_t.

In most cases, it is enough. But if anything goes wrong, you have several options:

if your compiler provides stdint.h, you can define AE_HAVE_STDINT conditional symbol
alternatively, you can manually define AE_INT32_T and/or AE_INT64_T and/or AE_UINT64_T and/or AE_INT_T symbols. Just assign datatype name to them, and ALGLIB will automatically use your definition. You may define all or just one/two types (those which are not detected automatically).

4.3 Improving performance (CPU-specific and OS-specific optimizations)

You can improve performance of ALGLIB in a several ways:

by compiling with advanced optimization turned on
by telling ALGLIB about CPU it will run on

ALGLIB has two-layered structure: some set of basic performance-critical primitives is implemented using optimized code, and the rest of the library is built on top of these primitives. By default, ALGLIB uses generic C code to implement these primitives (matrix multiplication, decompositions, etc.). This code works everywhere from Intel to SPARC. However, you can tell ALGLIB that it will work under particular architecture by defining appropriate macro at the global level:

defining AE_CPU=AE_INTEL - to tell ALGLIB that it will work under Intel

When AE_CPU macro is defined and equals to the AE_INTEL, it enables SIMD support. ALGLIB will use cpuid instruction to determine SIMD presence at run-time and use SIMD-capable code. ALGLIB uses SIMD intrinsics which are portable across different compilers and efficient enough for most practical purposes.

If you want to use multithreaded capabilities of commercial version of ALGLIB, you should compile it in OS-specific mode by #defining either AE_OS=AE_WINDOWS, AE_OS=AE_POSIX or AE_OS=AE_LINUX (POSIX with Linux-specific extensions) at compile time, depending on OS being used. Former corresponds to any modern OS (32/64-bit Windows XP and later) from Windows family, while latter means almost any POSIX-compatible OS (or any OS from the Linux family). It applies only to commercial version of ALGLIB. Open source version is always OS-agnostic, even in the presence of OS-specific definitions.

4.4 Examples (free and commercial editions)

4.4.1 Introduction

In this section we'll consider different compilation scenarios for free and commercial versions of ALGLIB - from simple platform-agnostic compilation to compiling/linking with MKL extensions.

We assume that you unpacked ALGLIB distribution in the current directory and saved here demo.cpp file, whose code is given below. Thus, in the current directory you should have exactly one file (demo.cpp) and exactly one subdirectory (cpp folder with ALGLIB distribution).

4.4.2 Compiling under Windows

File listing below contains the very basic program which uses ALGLIB to perform matrix-matrix multiplication. After that program evaluates performance of GEMM (function being called) and prints result to console. We'll show how performance of this program continually increases as we add more and more sophisticated compiler options.

demo.cpp (WINDOWS EXAMPLE)

#include <stdio.h>
#include <windows.h>
#include "LinAlg.h"

double counter()
{
    return 0.001*GetTickCount();
}

int main()
{
    alglib::real_2d_array a, b, c;
    int n = 2000;
    int i, j;
    double timeneeded, flops;
    
    // Initialize arrays
    a.setlength(n, n);
    b.setlength(n, n);
    c.setlength(n, n);
    for(i=0; i<n; i++)
        for(j=0; j<n; j++)
        {
            a[i][j] = alglib::randomreal()-0.5;
            b[i][j] = alglib::randomreal()-0.5;
            c[i][j] = 0.0;
        }
    
    // Set global threading settings (applied to all ALGLIB functions);
    // default is to perform serial computations, unless parallel execution
    // is activated. Parallel execution tries to utilize all cores; this
    // behavior can be changed with alglib::setnworkers() call.
    alglib::setglobalthreading(alglib::parallel);
    
    // Perform matrix-matrix product.
    flops = 2*pow((double)n, (double)3);
    timeneeded = counter();
    alglib::rmatrixgemm(
        n, n, n,
        1.0,
        a, 0, 0, 0,
        b, 0, 0, 1,
        0.0,
        c, 0, 0);
    timeneeded = counter()-timeneeded;
    
    // Evaluate performance
    printf("Performance is %.1f GFLOPS\n", (double)(1.0E-9*flops/timeneeded));
    
    return 0;
}

Examples below cover Windows compilation from command line with MSVC. It is very straightforward to adapt them to compilation from MSVC IDE - or to another compilers. We assume that you already called %VCINSTALLDIR%\bin\amd64\vcvars64.bat batch file which loads 64-bit build environment (or its 32-bit counterpart). We also assume that current directory is clean before example is executed (i.e. it has ONLY demo.cpp file and cpp folder). We used 3.2 GHz 4-core CPU for this test.

First example covers platform-agnostic compilation without optimization settings - the most simple way to compile ALGLIB. This step is same in both open source and commercial editions. However, in platform-agnostic mode ALGLIB is unable to use all performance related features present in commercial edition.

We starts from copying all cpp and h files to current directory, then we will compile them along with demo.cpp. In this and following examples we will omit compiler output for the sake of simplicity.

OS-agnostic mode, no compiler optimizations

> copy cpp\src\*.* .
> cl /I. /EHsc /Fedemo.exe *.cpp
> demo.exe
Performance is 0.7 GFLOPS

Well, 0.7 GFLOPS is not very impressing for a 3.2GHz CPU... Let's add /Ox to compiler parameters.

OS-agnostic mode, /Ox optimization

> cl /I. /EHsc /Fedemo.exe /Ox *.cpp
> demo.exe
Performance is 0.9 GFLOPS

Still not impressed. Let's turn on optimizations for x86 architecture: define AE_CPU=AE_INTEL. This option provides some speed-up in both free and commercial editions of ALGLIB.

OS-agnostic mode, ALGLIB knows it is x86/x64

> cl /I. /EHsc /Fedemo.exe /Ox /DAE_CPU=AE_INTEL *.cpp
> demo.exe
Performance is 4.5 GFLOPS

It is good, but we have 4 cores - and only one of them was used. Defining AE_OS=AE_WINDOWS allows ALGLIB to use Windows threads to parallelize execution of some functions. Starting from this moment, our example applies only to Commercial Edition.

ALGLIB knows it is Windows on x86/x64 CPU (COMMERCIAL EDITION)

> cl /I. /EHsc /Fedemo.exe /Ox /DAE_CPU=AE_INTEL /DAE_OS=AE_WINDOWS *.cpp
> demo.exe
Performance is 16.0 GFLOPS

Not bad. And now we are ready to the final test - linking with MKL extensions.

Linking with MKL extensions differs a bit from standard way of linking with ALGLIB. ALGLIB itself is compiled with one more preprocessor definition: we define AE_MKL symbol. We also link ALGLIB with appropriate (32-bit or 64-bit) alglib???_??mkl.lib static library, which is an import library for special lightweight MKL distribution, shipped with ALGLIB. We also should copy to current directory appropriate alglib???_??mkl.dll binary file which contains Intel MKL.

Linking with MKL extensions (COMMERCIAL EDITION)

> copy cpp\addons-mkl\alglib*64mkl.lib .
> copy cpp\addons-mkl\alglib*64mkl.dll .
> cl /I. /EHsc /Fedemo.exe /Ox /DAE_CPU=AE_INTEL /DAE_OS=AE_WINDOWS /DAE_MKL *.cpp alglib*64mkl.lib
> demo.exe
Performance is 33.1 GFLOPS

From 0.7 GFLOPS to 33.1 GFLOPS - you may see that commercial version of ALGLIB is really worth it!

4.4.3 Compiling under Linux

demo.cpp (LINUX EXAMPLE)

#include <stdio.h>
#include <sys/time.h>
#include "LinAlg.h"

double counter()
{
    struct timeval now;
    alglib_impl::ae_int64_t r, v;
    gettimeofday(&now, NULL);
    v = now.tv_sec;
    r = v*1000;
    v = now.tv_usec/1000;
    r = r+v;
    return 0.001*r;
}

int main()
{
    alglib::real_2d_array a, b, c;
    int n = 2000;
    int i, j;
    double timeneeded, flops;
    
    // Initialize arrays
    a.setlength(n, n);
    b.setlength(n, n);
    c.setlength(n, n);
    for(i=0; i<n; i++)
        for(j=0; j<n; j++)
        {
            a[i][j] = alglib::randomreal()-0.5;
            b[i][j] = alglib::randomreal()-0.5;
            c[i][j] = 0.0;
        }
    
    // Set global threading settings (applied to all ALGLIB functions);
    // default is to perform serial computations, unless parallel execution
    // is activated. Parallel execution tries to utilize all cores; this
    // behavior can be changed with alglib::setnworkers() call.
    alglib::setglobalthreading(alglib::parallel);
    
    // Perform matrix-matrix product.
    flops = 2*pow((double)n, (double)3);
    timeneeded = counter();
    alglib::rmatrixgemm(
        n, n, n,
        1.0,
        a, 0, 0, 0,
        b, 0, 0, 1,
        0.0,
        c, 0, 0);
    timeneeded = counter()-timeneeded;
    
    // Evaluate performance
    printf("Performance is %.1f GFLOPS\n", (double)(1.0E-9*flops/timeneeded));
    
    return 0;
}

Examples below cover x64 Linux compilation from command line with GCC. We assume that current directory is clean before example is executed (i.e. it has ONLY demo.cpp file and cpp folder). We used 2.3 GHz 2-core Skylake CPU with 2x Hyperthreading enabled for this test.

OS-agnostic mode, no compiler optimizations

> cp cpp/src/* .
> g++ -I. -o demo.out *.cpp
> ./demo.out
Performance is 0.9 GFLOPS

Let's add -O3 to compiler parameters.

OS-agnostic mode, -O3 optimization

> g++ -I. -o demo.out -O3 *.cpp
> ./demo.out
Performance is 2.8 GFLOPS

Better, but not impressed. Let's turn on optimizations for x86 architecture: define AE_CPU=AE_INTEL. This option provides some speed-up in both free and commercial editions of ALGLIB.

OS-agnostic mode, ALGLIB knows it is x86/x64

> g++ -I. -o demo.out -O3 -DAE_CPU=AE_INTEL *.cpp
> ./demo.out
Performance is 5.0 GFLOPS

It is good, but we have 4 cores (in fact, 2 cores - it is 2-way hyperthreaded system) and only one of them was used. Defining AE_OS=AE_POSIX allows ALGLIB to use POSIX threads to parallelize execution of some functions. You should also specify -pthread flag to link with pthreads standard library. Starting from this moment, our example applies only to Commercial Edition.

ALGLIB knows it is POSIX OS on x86/x64 CPU (COMMERCIAL EDITION)

> g++ -I. -o demo.out -O3 -DAE_CPU=AE_INTEL -DAE_OS=AE_POSIX -pthread *.cpp
> ./demo.out
Performance is 9.0 GFLOPS

Not bad. You may notice that performance growth was ~2x, not 4x. The reason is that we tested ALGLIB on hyperthreaded system: although we have 4 logical cores, they share computational resources of just 2 physical cores. And now we are ready to the final test - linking with MKL extensions.

We should note that on typical Linux system shared libraries are not loaded from current directory by default. Either you install them into one of the system directories, or use some way to tell linker/loader that you want to load shared library from some specific directory. For our example we choose to update LD_LIBRARY_PATH environment variable.

Linking with MKL extensions (COMMERCIAL EDITION, relevant for ALGLIB 3.13)

> cp cpp/addons-mkl/libalglib*64mkl.so .
> ls *.so
libalglib313_64mkl.so
> g++ -I. -o demo.out -O3 -DAE_CPU=AE_INTEL -DAE_OS=AE_POSIX -pthread -DAE_MKL -L. *.cpp -lalglib313_64mkl
> LD_LIBRARY_PATH=.
> export LD_LIBRARY_PATH
> ./demo.out
Performance is 33.8 GFLOPS

Final result: from 0.9 GFLOPS to 33.8 GFLOPS!

5 Using ALGLIB

5.1 Thread-safety

Both open source and commercial versions of ALGLIB are 100% thread-safe as long as different user threads work with different instances of objects/arrays. Thread-safety is guaranteed by having no global shared variables.

However, any kind of sharing ALGLIB objects/arrays between different threads is potentially hazardous. Even when this object is seemingly used in read-only mode!

Say, you use ALGLIB neural network NET to process two input vectors X0 and X1, and get two output vectors Y0 and Y1. You may decide that neural network is used in read-only mode which does not change state of NET, because output is written to distinct arrays Y. Thus, you may want to process these vectors from parallel threads.

But it is not read-only operation, even if it looks like that! Neural network object NET allocates internal temporary buffers, which are modified by neural processing functions. Thus, sharing one instance of neural network between two threads is thread-unsafe!

5.2 Global definitions

ALGLIB defines several conditional symbols (all start with "AE_" which means "ALGLIB environment") and two namespaces: alglib_impl (contains computational core) and alglib (contains C++ interface).

Although this manual mentions both alglib_impl and alglib namespaces, only alglib namespace should be used by you. It contains user-friendly C++ interface with automatic memory management, exception handling and all other nice features. alglib_impl is less user-friendly, is less documented, and it is too easy to crash your system or cause memory leak if you use it directly.

5.3 Datatypes

ALGLIB (ap.h header) defines several "basic" datatypes (types which are used by all packages) and many package-specific datatypes. "Basic" datatypes are:

alglib::ae_int_t - signed integer type used by library
alglib::complex - double precision complex datatype, safer replacement for std::complex
alglib::ap_error - exception which is thrown by library
boolean_1d_array - 1-dimensional boolean array
integer_1d_array - 1-dimensional integer array
real_1d_array - 1-dimensional real (double precision) array
complex_1d_array - 1-dimensional complex array
boolean_2d_array - 2-dimensional boolean array
integer_2d_array - 2-dimensional integer array
real_2d_array - 2-dimensional real (double precision) array
complex_2d_array - 2-dimensional complex array

Package-specific datatypes are classes which can be divided into two distinct groups:

"struct-like" classes with public fields which you can access directly.
"object-like" classes which have no public fields. You should use ALGLIB functions to work with them.

5.4 Constants

The most important constants (defined in the ap.h header) from ALGLIB namespace are:

alglib::machineepsilon - small number which is close to the double precision &eps;, but is slightly larger
alglib::maxrealnumber - very large number which is close to the maximum real number, but is slightly smaller
alglib::minrealnumber - very small number which is close to the minimum nonzero real number, but is slightly larger
alglib::fp_nan - NAN (non-signalling under most platforms except for PA-RISC, where it is signalling; but when PA-RISC CPU is in its default state, it is silently converted to the quiet NAN)
alglib::fp_posinf - positive infinity
alglib::fp_neginf - negative infinity

5.5 Functions

The most important "basic" functions from ALGLIB namespace (ap.h header) are:

alglib::randomreal() - returns random real number from [0,1)
alglib::randominteger(mx) - returns random integer number from [0,nx); mx must be less than RAND_MAX
alglib::fp_eq(v1,v2) - makes IEEE-compliant comparison of two double precision numbers. If numbers are represented with greater precision than specified by IEEE 754 (as with Intel 80-bit FPU), this functions converts them to 64 bits before comparing.
alglib::fp_neq(v1,v2) - makes IEEE-compliant comparison of two double precision numbers.
alglib::fp_less(v1,v2) - makes IEEE-compliant comparison of two double precision numbers.
alglib::fp_less_eq(v1,v2) - makes IEEE-compliant comparison of two double precision numbers.
alglib::fp_greater(v1,v2) - makes IEEE-compliant comparison of two double precision numbers.
alglib::fp_greater_eq(v1,v2) - makes IEEE-compliant comparison of two double precision numbers.
alglib::fp_isnan - checks whether number is NAN
alglib::fp_isposinf - checks whether number is +INF
alglib::fp_isneginf - checks whether number is -INF
alglib::fp_isinf - checks whether number is +INF or -INF
alglib::fp_isfinite - checks whether number is finite value (possibly subnormalized)

5.6 Working with vectors and matrices

ALGLIB (ap.h header) supports matrixes and vectors (one-dimensional and two-dimensional arrays) of variable size, with numeration starting from zero.

Everything starts from array creation. You should distinguish the creation of array class instance and the memory allocation for the array. When creating the class instance, you can use constructor without any parameters, that creates an empty array without any elements. An attempt to address them may cause the program failure.

You can use copy and assignment constructors that copy one array into another. If, during the copy operation, the source array has no memory allocated for the array elements, destination array will contain no elements either. If the source array has memory allocated for its elements, destination array will allocate the same amount of memory and copy the elements there. That is, the copy operation yields into two independent arrays with indentical contents.

You can also create array from formatted string like "[]", "[true,FALSE,tRUe]", "[[]]]" or "[[1,2],[3.2,4],[5.2]]" (note: '.' is used as decimal point independently from locale settings).

alglib::boolean_1d_array b1;
b1 = "[true]";

alglib::real_2d_array r2("[[2,3],[3,4]]");
alglib::real_2d_array r2_1("[[]]");
alglib::real_2d_array r2_2(r2);
r2_1 = r2;

alglib::complex_1d_array c2;
c2 = "[]";
c2 = "[0]";
c2 = "[1,2i]";
c2 = "[+1-2i,-1+5i]";
c2 = "[ 4i-2,  8i+2]";
c2 = "[+4i-2, +8i+2]";
c2 = "[-4i-2, -8i+2]";

After an empty array has been created, you can allocate memory for its elements, using the setlength() method. The content of the created array elements is not defined. If the setlength method is called for the array with already allocated memory, then, after changing its parameters, the newly allocated elements also become undefined and the old content is destroyed.

alglib::boolean_1d_array b1;
b1.setlength(2);

alglib::integer_2d_array r2;
r2.setlength(4,3);

Another way to initialize array is to call setcontent() method. This method accepts pointer to data which are copied into newly allocated array. Vectors are stored in contiguous order, matrices are stored row by row.

alglib::real_1d_array r1;
double _r1[] = {2, 3};
r1.setcontent(2,_r1);

alglib::real_2d_array r2;
double _r2[] = {11, 12, 13, 21, 22, 23};
r2.setcontent(2,3,_r2);

You can also attach real vector/matrix object to already allocated double precision array (attaching to boolean/integer/complex arrays is not supported). In this case, no actual data is copied, and attached vector/matrix object becomes a read/write proxy for external array.

alglib::real_1d_array r1;
double a1[] = {2, 3};
r1.attach_to_ptr(2,a1);

alglib::real_2d_array r2;
double a2[] = {11, 12, 13, 21, 22, 23};
r2.attach_to_ptr(2,3,_r2);

To access the array elements, an overloaded operator() or operator[] can used. That is, the code addressing the element of array a with indexes [i,j] can look like a(i,j) or a[i][j].

alglib::integer_1d_array a("[1,2,3]");
alglib::integer_1d_array b("[3,9,27]");
a[0] = b(0);

alglib::integer_2d_array c("[[1,2,3],[9,9,9]]");
alglib::integer_2d_array d("[[3,9,27],[8,8,8]]");
d[1][1] = c(0,0);

You can access contents of 1-dimensional array by calling getcontent() method which returns pointer to the array memory. For historical reasons 2-dimensional arrays do not provide getcontent() method, but you can use create reference to any element of array. 2-dimensional arrays store data in row-major order with aligned rows (i.e. generally distance between rows is not equal to number of columns). You can get stride (distance between consequtive elements in different rows) with getstride() call.

alglib::integer_1d_array a("[1,2]");
alglib::real_2d_array b("[[0,1],[10,11]]");

alglib::ae_int_t *a_row = a.getcontent();

// all three pointers point to the same location
double *b_row0 = &b[0][0];
double *b_row0_2 = &b(0,0);
double *b_row0_3 = b[0];

// advancing to the next row of 2-dimensional array
double *b_row1 = b_row0 + b.getstride();

Finally, you can get array size with length(), rows() or cols() methods:

alglib::integer_1d_array a("[1,2]");
alglib::real_2d_array b("[[0,1],[10,11]]");

printf("%ld\n", (long)a.length());
printf("%ld\n", (long)b.rows());
printf("%ld\n", (long)b.cols());

5.7 Using functions: 'expert' and 'friendly' interfaces

Most ALGLIB functions provide two interfaces: 'expert' and 'friendly'. What is the difference between two? When you use 'friendly' interface, ALGLIB:

tries to automatically determine size of input arguments
throws an exception when arguments have inconsistent size (for example, square matrix is expected, but non-square is passed; another example - two parameters must have same size, but have different size)
if semantics of input parameter assumes that it is symmetric/Hermitian matrix, checks that lower triangle is equal to upper triangle (conjugate of upper triangle) and throws an exception otherwise
if semantics of output parameter assumes that it is symmetric/Hermitian matrix, returns full matrix (both upper and lower triangles)

When you use 'expert' interface, ALGLIB:

requires caller to explicitly specify size of input arguments. If vector/matrix is larger than size being specified (say, N), only N leading elements are used
if semantics of input parameter assumes that it is symmetric/Hermitian matrix, uses only upper or lower triangle of input matrix and requires caller to specify what triangle to use
if semantics of output parameter assumes that it is symmetric/Hermitian matrix, returns only upper or lower triangle (when you look at specific function, it is clear what triangle is returned)

Here are several examples of 'friendly' and 'expert' interfaces:

#include "interpolation.h"

...

alglib::real_1d_array    x("[0,1,2,3]");
alglib::real_1d_array    y("[1,5,3,9]");
alglib::real_1d_array   y2("[1,5,3,9,0]");
alglib::spline1dinterpolant s;

alglib::spline1dbuildlinear(x, y, 4, s);  // 'expert' interface is used
alglib::spline1dbuildlinear(x, y, s);     // 'friendly' interface - input size is
                                          // automatically determined

alglib::spline1dbuildlinear(x, y2, 4, s); // y2.length() is 5, but it will work

alglib::spline1dbuildlinear(x, y2, s);    // it won't work because sizes of x and y2
                                          // are inconsistent

'Friendly' interface - matrix semantics:

#include "linalg.h"

...

alglib::real_2d_array a;
alglib::matinvreport  rep;
alglib::ae_int_t      info;

// 
// 'Friendly' interface: spdmatrixinverse() accepts and returns symmetric matrix
// 

// symmetric positive definite matrix
a = "[[2,1],[1,2]]";

// after this line A will contain [[0.66,-0.33],[-0.33,0.66]]
// which is symmetric too
alglib::spdmatrixinverse(a, info, rep); 

// you may try to pass nonsymmetric matrix
a = "[[2,1],[0,2]]";

// but exception will be thrown in such case
alglib::spdmatrixinverse(a, info, rep);

Same function but with 'expert' interface:

#include "linalg.h"

...

alglib::real_2d_array a;
alglib::matinvreport  rep;
alglib::ae_int_t      info;

// 
// 'Expert' interface, spdmatrixinverse()
// 

// only upper triangle is used; a[1][0] is initialized by NAN,
// but it can be arbitrary number
a = "[[2,1],[NAN,2]]";

// after this line A will contain [[0.66,-0.33],[NAN,0.66]]
// only upper triangle is modified
alglib::spdmatrixinverse(a, 2 /* N */, true /* upper triangle is used */, info, rep);

5.8 Handling errors

ALGLIB uses two error handling strategies:

returning error code
throwing exception

What is actually done depends on function being used and error being reported:

if function returns some error code and has corresponding value for this kind of error, ALLGIB returns error code
if function does not return error code (or returns error code, but there is no code for error being reported), ALGLIB throws alglib::ap_error exception. Exception object has msg parameter which contains short description of error.

To make things clear we consider several examples of error handling.

Example 1. mincgreate function creates nonlinear CG optimizer. It accepts problem size N and initial point X. Several things can go wrong - you may pass array which is too short, filled by NAN's, or otherwise pass incorrect data. However, this function returns no error code - so it throws an exception in case something goes wrong. There is no other way to tell caller that something went wrong.

Example 2. rmatrixinverse function calculates inverse matrix. It returns error code, which is set to +1 when problem is solved and is set to -3 if singular matrix was passed to the function. However, there is no error code for matrix which is non-square or contains infinities. Well, we could have created corresponding error codes - but we didn't. So if you pass singular matrix to rmatrixinverse, you will get completion code -3. But if you pass matrix which contains INF in one of its elements, alglib::ap_error will be thrown.

First error handling strategy (error codes) is used to report "frequent" errors which can occur during normal execution of user program. Second error handling strategy (exceptions) is used to report "rare" errors which are result of serious flaws in your program (or ALGLIB) - infinities/NAN's in the inputs, inconsistent inputs, etc.

5.9 Working with Level 1 BLAS functions

ALGLIB (ap.h header) includes following Level 1 BLAS functions:

alglib::vdotproduct() family, which allows to calculate dot product of two real or complex vectors
alglib::vmove() family, which allows to copy real/complex vector to another location with optimal multiplication by real/complex value
alglib::vmoveneg() family, which allows to copy real/complex vector to another location with multiplication by -1
alglib::vadd() and alglib::vsub() families, which allows to add or subtract two real/complex vectors with optimal multiplication by real/complex value
alglib::vmul() family, which implements in-place multiplication of real/complex vector by real/complex value

Each Level 1 BLAS function accepts input stride and output stride, which are expected to be positive. Input and output vectors should not overlap. Functions operating with complex vectors accept additional parameter conj_src, which specifies whether input vector is conjugated or not.

For each real/complex function there exists "simple" companion which accepts no stride or conjugation modifier. "Simple" function assumes that input/output stride is +1, and no input conjugation is required.

alglib::real_1d_array    rvec("[0,1,2,3]");
alglib::real_2d_array    rmat("[[1,2],[3,4]]");
alglib::complex_1d_array cvec("[0+1i,1+2i,2-1i,3-2i]");
alglib::complex_2d_array cmat("[[3i,1],[9,2i]]");

alglib::vmove(&rvec[0],  1, &rmat[0][0], rmat.getstride(), 2); // now rvec is [1,3,2,3]

alglib::vmove(&cvec[0],  1, &cmat[0][0], rmat.getstride(), "No conj", 2); // now cvec is [3i, 9, 2-1i, 3-2i]
alglib::vmove(&cvec[2],  1, &cmat[0][0], 1,                "Conj", 2);    // now cvec is [3i, 9, -3i,  1]

Here is full list of Level 1 BLAS functions implemented in ALGLIB:

double vdotproduct(
    const double *v0,
     ae_int_t stride0,
     const double *v1,
     ae_int_t stride1,
     ae_int_t n);
double vdotproduct(
    const double *v1,
     const double *v2,
     ae_int_t N);

alglib::complex vdotproduct(
    const alglib::complex *v0,
     ae_int_t stride0,
     const char *conj0,
     const alglib::complex *v1,
     ae_int_t stride1,
     const char *conj1,
     ae_int_t n);
alglib::complex vdotproduct(
    const alglib::complex *v1,
     const alglib::complex *v2,
     ae_int_t N);

void vmove(
    double *vdst,
      ae_int_t stride_dst,
     const double* vsrc,
      ae_int_t stride_src,
     ae_int_t n);
void vmove(
    double *vdst,
     const double* vsrc,
     ae_int_t N);

void vmove(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex* vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n);
void vmove(
    alglib::complex *vdst,
     const alglib::complex* vsrc,
     ae_int_t N);

void vmoveneg(
    double *vdst,
      ae_int_t stride_dst,
     const double* vsrc,
      ae_int_t stride_src,
     ae_int_t n);
void vmoveneg(
    double *vdst,
     const double *vsrc,
     ae_int_t N);

void vmoveneg(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex* vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n);
void vmoveneg(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N);

void vmove(
    double *vdst,
      ae_int_t stride_dst,
     const double* vsrc,
      ae_int_t stride_src,
     ae_int_t n,
     double alpha);
void vmove(
    double *vdst,
     const double *vsrc,
     ae_int_t N,
     double alpha);

void vmove(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex* vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n,
     double alpha);
void vmove(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N,
     double alpha);

void vmove(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex* vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n,
     alglib::complex alpha);
void vmove(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N,
     alglib::complex alpha);

void vadd(
    double *vdst,
      ae_int_t stride_dst,
     const double *vsrc,
      ae_int_t stride_src,
     ae_int_t n);
void vadd(
    double *vdst,
     const double *vsrc,
     ae_int_t N);

void vadd(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex *vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n);
void vadd(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N);

void vadd(
    double *vdst,
      ae_int_t stride_dst,
     const double *vsrc,
      ae_int_t stride_src,
     ae_int_t n,
     double alpha);
void vadd(
    double *vdst,
     const double *vsrc,
     ae_int_t N,
     double alpha);

void vadd(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex *vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n,
     double alpha);
void vadd(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N,
     double alpha);

void vadd(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex *vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n,
     alglib::complex alpha);
void vadd(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N,
     alglib::complex alpha);

void vsub(
    double *vdst,
      ae_int_t stride_dst,
     const double *vsrc,
      ae_int_t stride_src,
     ae_int_t n);
void vsub(
    double *vdst,
     const double *vsrc,
     ae_int_t N);

void vsub(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex *vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n);
void vsub(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N);

void vsub(
    double *vdst,
      ae_int_t stride_dst,
     const double *vsrc,
      ae_int_t stride_src,
     ae_int_t n,
     double alpha);
void vsub(
    double *vdst,
     const double *vsrc,
     ae_int_t N,
     double alpha);

void vsub(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex *vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n,
     double alpha);
void vsub(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N,
     double alpha);

void vsub(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     const alglib::complex *vsrc,
     ae_int_t stride_src,
     const char *conj_src,
     ae_int_t n,
     alglib::complex alpha);
void vsub(
    alglib::complex *vdst,
     const alglib::complex *vsrc,
     ae_int_t N,
     alglib::complex alpha);

void vmul(
    double *vdst,
      ae_int_t stride_dst,
     ae_int_t n,
     double alpha);
void vmul(
    double *vdst,
     ae_int_t N,
     double alpha);

void vmul(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     ae_int_t n,
     double alpha);
void vmul(
    alglib::complex *vdst,
     ae_int_t N,
     double alpha);

void vmul(
    alglib::complex *vdst,
     ae_int_t stride_dst,
     ae_int_t n,
     alglib::complex alpha);
void vmul(
    alglib::complex *vdst,
     ae_int_t N,
     alglib::complex alpha);

5.10 Reading data from CSV files

ALGLIB (ap.h header) has alglib::read_csv() function which allows to read data from CSV file. Entire file is loaded into memory as double precision 2D array (alglib::real_2d_array object). This function provides following features:

support for ASCI encoding and UTF-8 without BOM (in header names)
ability to use any character (comma/tab/space) as field separator (as long as it is distinct from one used for decimal point)
ability to use both comma and full stop as decimal point (as long is decimal point is distinct from field separator).
support for Unix and Windows text files (CR vs CRLF)

See comments on alglib::read_csv() function for more information about its functionality.

6 Working with commercial version

6.1 Benefits of commercial version

Commercial version of ALGLIB for C++ features four important improvements over open source one:

License. Commercial license used by ALGLIB is friendly to closed source applications. Unlike GPL, it does not require you to open source your application. Thus, almost any commercial software developer is interested in obtaining commercial license.
Low-level optimizations. Commercial version of ALGLIB includes SIMD-optimized versions of many computationally intensive functions. It allows to increase performance on Intel/AMD platforms while still being able to use software under non-x86 CPU's.
Multithreading. Commercial version of ALGLIB can utilize multicore capabilities of modern CPU's. Large computational problems can be automatically split between different cores. ALGLIB uses its own multithreading framework which does not depend on vendor/compiler support for technologies like OpenMP/MPI/... It gives ALGLIB unprecedented portability across operating systems and compilers.
Integrated Intel MKL. Commercial version of ALGLIB includes special MKL extensions - special lightweight distribution of Intel MKL, high-performance numerical analysis library, accompanied by ALGLIB-MKL interface. We obtained license from Intel which allows us to integrate MKL into ALGLIB distribution. Linking with MKL accelerates many ALGLIB functions, however due to license restrictions you can not use MKL directly (i.e. bypass ALGLIB interface between your program and MKL).

6.2 Working with SIMD support (Intel/AMD users)

ALGLIB for C++ can utilize SIMD instructions supported by Intel and AMD processors. This feature is optional and must be explicitly turned on during compile-time. If you do not activate it, ALGLIB will use generic C code, without any processor-specific assembly/intrinsics.

Thus, if you turn on this feature, your code will run faster on x86_32 and x86_64 processors, but will be unportable to non-x86 platforms (and Intel MIC platform, which is not exactly x86!). From the other side, if you do not activate this feature, your code will be portable to almost any modern CPU (SPARC, ARM, ...).

In order to turn on x86-specific optimizations, you should define AE_CPU=AE_INTEL preprocessor definition at global level. It will tell ALGLIB to use SIMD intrinsics supported by GCC, MSVC and Intel compilers. Additionally you should tell compiler to generate SIMD-capable code. It can be done in the project settings of your IDE or in the command line:


GCC example:
> g++ -msse2 -I. -DAE_CPU=AE_INTEL *.cpp -lm

MSVC example:
> cl /I. /EHsc /DAE_CPU=AE_INTEL *.cpp

6.3 Using multithreading

6.3.1 General information

Commercial version of ALGLIB includes out-of-the-box support for multithreading. Many (not all) computationally intensive problems can be solved in multithreaded mode. You should read comments on specific ALGLIB functions to determine what can be multithreaded and what can not.

ALGLIB does not depend on vendor/compiler support for technologies like OpenMP/MPI/... Under Windows ALGLIB uses OS threads and custom synchronization framework. Under POSIX-compatible OS (Solaris, Linux, FreeBSD, NetBSD, OpenBSD, ...) ALGLIB uses POSIX Threads (standard *nix library which is shipped with any POSIX system) with its threading and synchronization primitives. It gives ALGLIB unprecedented portability across operating systems and compilers. ALGLIB does not depend on presence of any custom multithreading library or compiler support for any multithreading technology.

If you want to use multithreaded capabilities of ALGLIB, you should:

compile it in OS-specific mode (ALGLIB have to know what OS it is running on)
activate multithreading at global level with alglib::setglobalthreading function call or enable it at per-function basis by passing alglib::parallel to the specific computational function
optionally, tell ALGLIB about number of worker threads to use (default is to utilize all cores)

Let explain it in more details...

1. You should compile ALGLIB in OS-specific mode by #defining either AE_OS=AE_WINDOWS or AE_OS=AE_POSIX (or AE_OS=AE_LINUX, which means "POSIX with Linux extensions") at compile time, depending on OS being used. Former corresponds to any modern OS (32/64-bit Windows XP and later) from Windows family, while latter two mean almost any POSIX-compatible OS or any OS from the Linux family. When compiling on POSIX/Linux, do not forget to link ALGLIB with libpthread library.

2. By default, ALGLIB is configured to perform all calculations in serial manner. Parallel execution can be manually enabled at either global or per-function level. Former means that all ALGLIB functions will be able to parallelize themselves at their discretion. Similarly, you can enable global parallelism, but selectively disable it for specific function calls.

Global parallelism is enabled by alglib::setglobalthreading(alglib::parallel) and disabled by alglib::setglobalthreading(alglib::serial) call. Function-level parallelism can be enabled (or disabled, if global default is to parallelize) by adding alglib::parallel (or alglib::serial) to the end of the parameters list of specific ALGLIB function being called.

Enabling parallelism does not guarantee that ALGLIB will parallelize its computations. Small problems (say, products of 64x64 matrices) can not be efficiently parallelized. ALGLIB will automatically decide whether your problem is large enough or not for efficient parallelization.

3. ALGLIB automatically determines number of cores on application startup. On Windows it is done using GetSystemInfo() call. On POSIX systems ALGLIB performs sysconf(_SC_NPROCESSORS_ONLN) system call. This system call is supported by all modern POSIX-compatible systems: Solaris, Linux, FreeBSD, NetBSD, OpenBSD.

By default, ALGLIB uses all available cores (when told to parallelize calculations). Such behavior may be changed with setnworkers() call:

alglib::setnworkers(0) = use all cores
alglib::setnworkers(-1) = leave one core unused
alglib::setnworkers(-2) = leave two cores unused
alglib::setnworkers(+2) = use 2 cores (even if you have more)

You may want to specify maximum number of worker threads during compile time by means of preprocessor definition AE_NWORKERS=N. You can add this definition to compiler command line or change corresponding project settings in your IDE. Here N can be any positive number. ALGLIB will use exactly N worker threads, unless being told to use less by setnworkers() call.

Some old POSIX-compatible operating systems do not support sysconf(_SC_NPROCESSORS_ONLN) system call which is required in order to automatically determine number of active cores. On these systems you should specify number of cores manually at compile time. Without it ALGLIB will run in single-threaded mode.

6.3.2 SMT (CMT/hyper-threading) issues

Simultaneous multithreading (SMT) also known as Hyper-threading (Intel) and Cluster-based Multithreading (AMD) is a CPU design where several (usually two) logical cores share resources of one physical core. Say, on dual-core system with 2x HT scale factor you will see 4 logical cores. Each pair of these 4 cores, however, share same hardware resources. Thus, you may get only marginal speedup when running highly optimized software which fully utilizes CPU resources.

Say, if one thread occupies floating-point unit, another thread on the same physical core may work with integer numbers at the same time without any performance penalties. In this case you may get some speedup due to having additional cores. But if both threads keep FPU unit 100% busy, they won't get any multithreaded speedup.

So, if 2 math-intensive threads are dispatched by OS scheduler to different physical cores, you will get 2x speedup due to use of multithreading. But if these threads are dispatched to different logical cores - but same physical core - you won't get any speedup at all! One physical core will be 100% busy, and another one will be 100% idle. From the other side, if you start four threads instead of two, your system will be 100% utilized independently of thread scheduling details.

Let we stress it one more time - multithreading speedup on SMT systems is highly dependent on number of threads you are running and decisions made by OS scheduler. It is not 100% deterministic! With "true SMP" when you run 2 threads, you get 2x speedup (or 1.95, or 1.80 - it depends on algorithm, but this factor is always same). With SMT when you run 2 threads you may get your 2x speedup - or no speedup at all. Modern OS schedulers do a good job on single-socket hardware, but even in this "simple" case they give no guarantees of fair distribution of hardware resources. And things become a bit tricky when you work with multi-socket hardware. On SMT systems the only guaranteed way to 100% utilize your CPU is to create as many worker threads as there are logical cores. In this case OS scheduler has no chance to make its work in a wrong way.

6.4 Linking with Intel MKL

6.4.1 Using lightweight Intel MKL supplied by ALGLIB Project

Commercial edition of ALGLIB includes MKL extensions - special lightweight distribution of Intel MKL, highly optimized numerical library from Intel - and precompiled ALGLIB-MKL interface libraries. Linking your programs with MKL extensions allows you to run ALGLIB with maximum performance. MKL binaries are delivered for x86/x64 Windows and x64 Linux platforms.

Unlike the rest of the library, MKL extensions are distributed in binary-only form. ALGLIB itself is still distributed in source code form, but Intel MKL and ALGLIB-MKL interface are distributed as precompiled dynamic/static libraries. We can not distribute them in source because of license restrictions associated with Intel MKL. Also due to license restrictions we can not give you direct access to MKL functionality. You may use MKL to accelerate ALGLIB - without paying for MKL license - but you may not call its functions directly. It is technically possible, but strictly prohibited by both MKL's EULA and ALGLIB License Agreement. If you want to work with MKL, you should obtain separate license from Intel (as of 2018, free licenses are available).

MKL extensions are located in the /cpp/addons-mkl subdirectory of the ALGLIB distribution. This directory includes following files:

alglib???_32mkl.lib - x86 Windows import library (MSVC) for Intel MKL extensions
alglib???_64mkl.lib - x64 Windows import library (MSVC) for Intel MKL extensions
alglib???_32mkl.dll - x86 Windows binary with Intel MKL inside
alglib???_64mkl.dll - x64 Windows binary with Intel MKL inside
libalglib???_64mkl.so - x64 Linux binary with Intel MKL inside

Here ??? stands for specific ALGLIB version: 313 for ALGLIB 3.13, and so on. Files above are just MKL extensions - ALGLIB itself is not included in these binaries, and you still have to compile primary ALGLIB distribution.

In order to activate MKL extensions you should:

#define globally AE_MKL to activate ALGLIB-MKL connection.
additionally, #define globally AE_OS=AE_WINDOWS or AE_OS=AE_POSIX to activate multithreading capabilities
additionally, #define globally AE_CPU=AE_INTEL to use SIMD instructions provided by x86/x64 CPU's in the rest of ALGLIB
Depending on your OS, perform one of the following:
- Windows: choose 32-bit or 64-bit import library alglib???_32/64mkl.lib and link it with your application
- Linux: choose 64-bit SO file alglib???_64mkl.so and link it with your application
place DLL/SO binary into directory where it can be found during application startup (Windows: application dir; Linux: one of the system directories)

Several examples of ALGLIB+MKL usage are given in the 'compiling ALGLIB: examples' section.

6.4.2 Using your own installation of Intel MKL

If you bought separate license for Intel MKL, and want to use your own installation of MKL - and not our lightweight distribution - then you should compile ALGLIB as it was told in the previous section, with all necessary preprocessor definitions (AE_OS=AE_WINDOWS or AE_OS=AE_POSIX, AE_CPU=AE_INTEL and AE_MKL defined). But instead of linking with MKL Extensions binary, you should add to your project alglib2mkl.c file from addons-mkl directory and compile it (as C file) along with the rest of ALGLIB.

This C file implements interface between MKL and ALGLIB. Having this file in your project and defining AE_MKL preprocessor definition results in ALGLIB using MKL functions.

However, this C file is just interface! It is your responsibility to make sure that C/C++ compiler can find MKL headers, and appropriate MKL static/dynamic libraries are linked to your application.

7 Advanced topics

7.1 Exception-free mode

ALGLIB for C++ can be compiled in exception-free mode, with exceptions (throw/try/catch constructs) being disabled at compiler level. Such feature is sometimes used by developers of embedded software.

ALGLIB uses two-level model of errors: "expected" errors (like degeneracy of linear system or inconsistency of linear constraints) are reported with dedicated completion codes, and "critical" errors (like malloc failures, unexpected NANs/INFs in the input data and so on) are reported with exceptions. The idea is that it is hard to put (and handle) completion codes in every ALGLIB function, so we use exceptions to signal errors which should never happen under normal circumstances.

Internally ALGLIB for C++ is implemented as C++ wrapper around computational core written in pure C. Thus, internals of ALGLIB core use C-specific methods of error handling - completion codes and setjmp/longjmp functions. These error handling strategies are combined with sophisticated machinery of C memory management which makes sure that not even a byte of dynamic memory is lost when we make longjmp to the error handler. So, the only point where C++ exceptions are actually used is a boundary between C core and C++ interface.

If you choose to use exceptions (default mode), ALGLIB will throw an exception with short textual description of the situation. And if you choose to work without exceptions, ALGLIB will set global error flag and silently return from the current function/constructor/... instead of throwing an exception. Due to portability issues this error flag is made to be a non-TLS variable, i.e. it is shared between different threads. So, you can use exception-free error handling only in single-threaded programs - although multithreaded programs won't break, there is no way to determine which thread caused an "exception without exceptions".

Exception-free method of reporting critical errors can be activated by #defining two preprocessor symbols at global level:

AE_NO_EXCEPTIONS - to switch from exception-based to exception-free code
AE_THREADING=AE_SERIAL_UNSAFE - to confirm that you are aware of limitations associated with exception-free mode (it does not support multithreading)

We must also note that exception-free mode is incompatible with OS-aware compiling: you can not have AE_OS=??? defined together with AE_NO_EXCEPTIONS.

After you #define all the necessary preprocessor symbols, two functions will appear in alglib namespace:

bool alglib::get_error_flag(const char **p_msg = NULL), which returns current error status (true is returned on error), with optional char** parameter used to get human-readable error message.
void alglib::clear_error_flag(), which clears error flag (ALGLIB functions set flag on failure, but do not clear it on successful calls)

You must check error flag after EVERY operation with ALGLIB objects and functions. In addition to calling computational ALGLIB functions, following kinds of operations may result in "exception":

calling default constructor for ALGLIB object (simply instantiating object may result in "exception"). Due to large size ALGLIB objects are allocated dynamically (they look like value types, but internally everything is stored in the heap memory). Thus every constructor, even default one, makes at least one malloc() call which may fail.
calling copy/assignment constructors for ALGLIB objects (same reason - malloc)
resizing arrays with setlength()/setcontents() and attaching to external memory with attach_to_ptr()

7.2 Partial compilation

Due to ALGLIB modular structure it is possible to selectively enable/disable some of its subpackages along with their dependencies. Deactivation of ALGLIB source code is performed at preprocessor level - compiler does not even see disabled code. Partial compilation can be used for two purposes:

to reduce code size in cases when linker is unable to remove unused code (happens with some compilers when ALGLIB is compiled as part of the shared library)
to reduce build time (disabled code is silently removed by preprocessor, thus no time is spent in the compiler)

You can activate partial compilation by #defining at global level following symbols:

AE_PARTIAL_BUILD - to disable everything not enabled by default
AE_COMPILE_SUBPACKAGE - to selectively enable subpackage SUBPACKAGE, with its name given in upper case (case sensitive). You may combine several definitions like this in order to enable several subpackages. Subpackage names can be found in the list of ALGLIB packages and subpackages.

7.3 Testing ALGLIB

There are three test suites in ALGLIB: computational tests, interface tests, extended tests. Computational tests are located in /tests/test_c.cpp. They are focused on numerical properties of algorithms, stress testing and "deep" tests (large automatically generated problems). They require significant amount of time to finish (tens of minutes).

Interface tests are located in /tests/test_i.cpp. These tests are focused on ability to correctly pass data between computational core and caller, ability to detect simple problems in inputs, and on ability to at least compile ALGLIB with your compiler. They are very fast (about a minute to finish including compilation time).

Extended tests are located in /tests/test_x.cpp. These tests are focused on testing some special properties (say, testing that cloning object indeed results in 100% independent copy being created) and performance of several chosen algorithms.

Running test suite is easy - just

compile one of these files (test_c.cpp, test_i.cpp or test_x.cpp) along with the rest of the library
launch executable you will get. It may take from several seconds (interface tests) to several minutes (computational tests) to get final results

If you want to be sure that ALGLIB will work with some sophisticated optimization settings, set corresponding flags during compile time. If your compiler/system are not in the list of supported ones, we recommend you to run both test suites. But if you are running out of time, run at least test_i.cpp.

8.1 `AlglibMisc` package
hqrnd		High quality random numbers generator
nearestneighbor		Nearest neighbor search: approximate and exact
xdebug		Debug functions to test ALGLIB interface generator

8.2 `DataAnalysis` package
bdss		Basic dataset functions
clustering		Clustering functions (hierarchical, k-means, k-means++)
datacomp		Backward compatibility functions
dforest		Decision forest classifier (regression model)
filters		Different filters used in data analysis
knn		K Nearest Neighbors classification/regression
lda		Linear discriminant analysis
linreg		Linear models
logit		Logit models
mcpd		Markov Chains for Population/proportional Data
mlpbase		Basic functions for neural networks
mlpe		Basic functions for neural ensemble models
mlptrain		Neural network training
pca		Principal component analysis
ssa		Singular Spectrum Analysis

8.3 `DiffEquations` package
odesolver		Ordinary differential equation solver

8.4 `FastTransforms` package
conv		Fast real/complex convolution
corr		Fast real/complex cross-correlation
fft		Real/complex FFT
fht		Real Fast Hartley Transform

8.5 `Integration` package
autogk		Adaptive 1-dimensional integration
gkq		Gauss-Kronrod quadrature generator
gq		Gaussian quadrature generator

8.6 `Interpolation` package
fitsphere		Fitting circle/sphere to data (least squares, minimum circumscribed, maximum inscribed, minimum zone)
idw		Inverse distance weighting: interpolation/fitting with improved Shepard-like algorithm
intcomp		Backward compatibility functions
lsfit		Fitting with least squates target function (linear and nonlinear least-squares)
parametric		Parametric curves
polint		Polynomial interpolation/fitting
ratint		Rational interpolation/fitting
rbf		Scattered N-dimensional interpolation with RBF models
spline1d		1D spline interpolation/fitting
spline2d		2D spline interpolation
spline3d		3D spline interpolation

8.7 `LinAlg` package
ablas		Level 2 and Level 3 BLAS operations
bdsvd		Bidiagonal SVD
evd		Direct and iterative eigensolvers
inverseupdate		Sherman-Morrison update of the inverse matrix
matdet		Determinant calculation
matgen		Random matrix generation
matinv		Matrix inverse
normestimator		Estimates norm of the sparse matrix (from below)
ortfac		Real/complex QR/LQ, bi(tri)diagonal, Hessenberg decompositions
rcond		Condition number estimate
schur		Schur decomposition
sparse		Sparse matrices
spdgevd		Generalized symmetric eigensolver
svd		Singular value decomposition
trfac		LU and Cholesky decompositions (dense and sparse)

8.8 `Optimization` package
minbc		Box constrained optimizer with fast activation of multiple constraints per step
minbleic		Bound constrained optimizer with additional linear equality/inequality constraints
mincg		Conjugate gradient optimizer
mincomp		Backward compatibility functions
minlbfgs		Limited memory BFGS optimizer
minlm		Improved Levenberg-Marquardt optimizer
minlp		Linear programming suite
minnlc		Nonlinearly constrained optimizer
minns		Nonsmooth constrained optimizer
minqp		Quadratic programming with bound and linear equality/inequality constraints
optguardapi		OptGuard integrity checking for nonlinear models

8.9 `Solvers` package
directdensesolvers		Direct dense linear solvers
directsparsesolvers		Direct sparse linear solvers
lincg		Sparse linear CG solver
linlsqr		Sparse linear LSQR solver
nleq		Solvers for nonlinear equations
polynomialsolver		Polynomial solver

8.10 `SpecialFunctions` package
airyf		Airy functions
bessel		Bessel functions
betaf		Beta function
binomialdistr		Binomial distribution
chebyshev		Chebyshev polynomials
chisquaredistr		Chi-Square distribution
dawson		Dawson integral
elliptic		Elliptic integrals
expintegrals		Exponential integrals
fdistr		F-distribution
fresnel		Fresnel integrals
gammafunc		Gamma function
hermite		Hermite polynomials
ibetaf		Incomplete beta function
igammaf		Incomplete gamma function
jacobianelliptic		Jacobian elliptic functions
laguerre		Laguerre polynomials
legendre		Legendre polynomials
normaldistr		Univarite and bivariate normal distribution PDF and CDF
poissondistr		Poisson distribution
psif		Psi function
studenttdistr		Student's t-distribution
trigintegrals		Trigonometric integrals

8.11 `Statistics` package
basestat		Mean, variance, covariance, correlation, etc.
correlationtests		Hypothesis testing: correlation tests
jarquebera		Hypothesis testing: Jarque-Bera test
mannwhitneyu		Hypothesis testing: Mann-Whitney-U test
stest		Hypothesis testing: sign test
studentttests		Hypothesis testing: Student's t-test
variancetests		Hypothesis testing: F-test and one-sample variance test
wsr		Hypothesis testing: Wilcoxon signed rank test

ablas_d_gemm		Matrix multiplication (single-threaded)
ablas_d_syrk		Symmetric rank-K update (single-threaded)

basestat_d_base		Basic functionality (moments, adev, median, percentile)
basestat_d_c2		Correlation (covariance) between two random variables
basestat_d_cm		Correlation (covariance) between components of random vector
basestat_d_cm2		Correlation (covariance) between two random vectors

1 Introduction

1.1 What is ALGLIB

1.2 ALGLIB license

1.3 Documentation license

1.4 Reference Manual and User Guide

1.5 Acknowledgements

2 ALGLIB structure

2.1 Packages

2.2 Subpackages

2.3 Open Source and Commercial versions

3 Compatibility

3.1 CPU

3.2 OS

3.3 Compiler

3.4 Optimization settings

4 Compiling ALGLIB

4.1 Adding to your project

4.2 Configuring for your compiler

4.3 Improving performance (CPU-specific and OS-specific optimizations)

4.4 Examples (free and commercial editions)

4.4.1 Introduction

4.4.2 Compiling under Windows

4.4.3 Compiling under Linux

5 Using ALGLIB

5.1 Thread-safety

5.2 Global definitions

5.3 Datatypes

5.4 Constants

5.5 Functions

5.6 Working with vectors and matrices

5.7 Using functions: 'expert' and 'friendly' interfaces

5.8 Handling errors

5.9 Working with Level 1 BLAS functions

5.10 Reading data from CSV files

6 Working with commercial version

6.1 Benefits of commercial version

6.2 Working with SIMD support (Intel/AMD users)

6.3 Using multithreading

6.3.1 General information

6.3.2 SMT (CMT/hyper-threading) issues

6.4 Linking with Intel MKL

6.4.1 Using lightweight Intel MKL supplied by ALGLIB Project

6.4.2 Using your own installation of Intel MKL

7 Advanced topics

7.1 Exception-free mode

7.2 Partial compilation

7.3 Testing ALGLIB

8 ALGLIB packages and subpackages

8.1 AlglibMisc package

8.2 DataAnalysis package

8.3 DiffEquations package

8.4 FastTransforms package

8.5 Integration package

8.6 Interpolation package

8.7 LinAlg package

8.8 Optimization package

8.9 Solvers package

8.10 SpecialFunctions package

8.11 Statistics package

ablas subpackage

Functions

Examples

cmatrixcopy function

cmatrixgemm function

cmatrixherk function

cmatrixlefttrsm function

cmatrixmv function

cmatrixrank1 function

cmatrixrighttrsm function

cmatrixsyrk function

cmatrixtranspose function

rmatrixcopy function

rmatrixenforcesymmetricity function

rmatrixgemm function

rmatrixgemv function

rmatrixgencopy function

rmatrixger function

rmatrixlefttrsm function

rmatrixmv function

rmatrixrank1 function

8.1 `AlglibMisc` package

8.2 `DataAnalysis` package

8.3 `DiffEquations` package

8.4 `FastTransforms` package

8.5 `Integration` package

8.6 `Interpolation` package

8.7 `LinAlg` package

8.8 `Optimization` package

8.9 `Solvers` package

8.10 `SpecialFunctions` package

8.11 `Statistics` package

`ablas` subpackage

`cmatrixcopy` function

`cmatrixgemm` function

`cmatrixherk` function

`cmatrixlefttrsm` function

`cmatrixmv` function

`cmatrixrank1` function

`cmatrixrighttrsm` function

`cmatrixsyrk` function

`cmatrixtranspose` function

`rmatrixcopy` function

`rmatrixenforcesymmetricity` function

`rmatrixgemm` function

`rmatrixgemv` function

`rmatrixgencopy` function

`rmatrixger` function

`rmatrixlefttrsm` function

`rmatrixmv` function

`rmatrixrank1` function

`rmatrixrighttrsm` function

`rmatrixsymv` function

`rmatrixsyrk` function

`rmatrixsyvmv` function

`rmatrixtranspose` function

`rmatrixtrsv` function

`rvectorcopy` function

`airyf` subpackage

`airy` function

`autogk` subpackage

`autogkreport` class

`autogkstate` class

`autogkintegrate` function

`autogkresults` function

`autogksingular` function

`autogksmooth` function

`autogksmoothw` function

`basestat` subpackage

`cov2` function

`covm` function

`covm2` function

`pearsoncorr2` function

`pearsoncorrelation` function

`pearsoncorrm` function

`pearsoncorrm2` function

`rankdata` function

`rankdatacentered` function

`sampleadev` function

`samplekurtosis` function

`samplemean` function

`samplemedian` function

`samplemoments` function

`samplepercentile` function

`sampleskewness` function

`samplevariance` function

`spearmancorr2` function

`spearmancorrm` function

`spearmancorrm2` function

`spearmanrankcorrelation` function

`bdss` subpackage

`dsoptimalsplit2` function

`dsoptimalsplit2fast` function

`bdsvd` subpackage

`rmatrixbdsvd` function

`bessel` subpackage

`besseli0` function

`besseli1` function

`besselj0` function

`besselj1` function

`besseljn` function