### DREAM Example 5 Setup and Model Definition

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianDREAMExamples5.html

Initializes the DREAM Example 5, sets parameters (5 dimensions, 50 chains, 1000 iterations), defines a custom model lambda and batch model lambda to evaluate the sin-wave decomposition.

```cpp
#endif
// using the default random engine, but must reset the random number generator
srand((int) time(nullptr));
// EXAMPLE 5:
cout << "\n" << "---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(5);
cout << "EXAMPLE 5: infer the frequency and magnitude of two signals from noisy data\n"
<< " the model has 5 parameters: f(x_1 ... x_5) = sum x_k sin(k * pi * t)\n"
<< " data = 2.0 * sin(2 * pi * t) + sin(4 * pi * t) + noise\n"
<< " t in [0, 1], t is discretized using 64 equidistant nodes\n"
<< " we use two different likelihood functions,"
<< "corresponding to l-2 and l-1 norms\n"
<< " we are looking for the mode of the posterior,"
<<"i.e., the optimal fit to the data\n\n";
constexpr double p = to be defined; // half-period of the std::sin()
int num_dimensions = 5;
int num_chains = 50;
int num_burnup_iterations = 1000;
int num_sample_iterations = 1000;
int num_discrete_nodes = 64;

```

```cpp
// create a lambda function that represents the model
// normally this would be a call to an external code
auto model = [&](std::vector<double> const &x, std::vector<double> &y)->
void{
double dt = 1.0 / ((double) y.size());
double t = 0.5 * dt;
for(auto &output : y){
output = 0.0;
int signal_number = 1;
for(auto const &weight : x)
output += weight * std::sin(double(signal_number++) * t * pi);

t += dt;
}
};

```

--------------------------------

### Example virtual environment install with CUDA

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2Installation.html

Complete workflow to create a virtual environment, enable CUDA, install Tasmanian with pip, and verify the installation.

```bash
python3 -m venv tasmanian_virtual_env # create a virtual environment
source ./tasmanian_virtual_env/bin/activate # activate the virtual environment
export Tasmanian_ENABLE_CUDA=/usr/local/cuda/bin/nvcc # specify the CUDA compiler
python -m pip install Tasmanian --no-cache-dir # will install Tasmanian with CUDA
python -m Tasmanian # print the install log
tasgrid -v # print the available CUDA devices
```

--------------------------------

### Build and install on MS Windows with CMake

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2Installation.html

After configuring with the CMake GUI, these commands build the Release configuration, run tests, and install into the configured folder. Use only one configuration (Release or Debug) at a time.

```cmd
cd <cmake-build-folder>
```

```cmd
cmake --build . --config Release
```

```cmd
ctest -C Release
```

```cmd
cmake --build . --config Release --target install
```

--------------------------------

### MPIGridBcast example: broadcast grid from root to all processes

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianAddonsMPIGridSend.html

Broadcast a grid from the root process to all processes in an MPI communicator. The root process reads the grid from a file, while other processes start with an empty grid. After the call, every process has the same grid as if they all read it from the file.

```cpp
int me;
MPI_Comm_rank(MPI_COMM_WORLD, &me);
auto grid = (me == 0) ? TasGrid::readGrid("foo") : TasGrid::TasmanianSparseGrid();
MPIGridBcast(grid, 0, MPI_COMM_WORLD);
// at this line, every process has the same grid as if they all read it from "foo"
```

--------------------------------

### sparse_grids_example_04: Build and evaluate a surrogate model

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianSGExamples4.html

Constructs a global grid with the Clenshaw-Curtis rule and total-degree polynomial space, evaluates the model at needed points, loads the values, and evaluates the interpolant at a point of interest. The example loops over two polynomial orders (6 and 12) and prints the number of grid points, the interpolant value, and the error against the reference solution. Note that the code snippet is incomplete as it starts with '#endif' and omits the function header and includes.

```cpp
#endif
cout << "\n---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(4);
cout << "Example 4: interpolate f(x,y) = exp(-x^2) * cos(y),"
<< "using clenshaw-curtis iptotal rule\n";
const int num_inputs = 2; // just x and y as inputs
const int num_outputs = 1; // using single-value output
const int num_samples = 1; // this is a huge number in an actual application
std::vector<double> pnts = { 0.3, 0.7 }; // points of interest where the model is not known
double reference_solution = std::exp(-0.3 * 0.3) * std::cos(0.7);
std::vector<int> polynomial_order = {6, 12}; // use two grids with different polynomial order
for(auto prec : polynomial_order){
// the first 3 inputs to make grid are always: inputs, outputs, depth
// TasGrid::type_iptotal means interpret the precision as "total degree polynomial space"
auto grid = TasGrid::makeGlobalGrid(num_inputs, num_outputs, prec,
TasGrid::type_iptotal,
TasGrid::rule_clenshawcurtis);
// having constructed the initial grid
// to complete the interpolant, Tasmanian needs the model values at the needed points
auto points = grid.getNeededPoints();
// allocate a vector for the associated model values
std::vector<double> model_values(grid.getNumNeeded() * num_outputs);
for(int i=0; i<grid.getNumNeeded(); i++){
double x = points[i*num_inputs];
double y = points[i*num_inputs+1];
model_values[i] = std::exp(-x*x) * std::cos(y);
}
// feed the data to the grid
grid.loadNeededValues(model_values);
// interpolation can be performed after loadNeededPoints() is called
// note that we don't need to set the size, Tasmanian will automatically set the right size
std::vector<double> result(num_samples * num_outputs);
// evaluate() deals with a single point
// if num_samples is large, evaluateBatch() is much more efficient
grid.evaluate(pnts, result);
cout << "\n using polynomials of total degree up to: " << prec << "\n"
<< " the grid has: " << grid.getNumPoints() << " poitns\n"
<< " interpolant at (0.3,0.7): " << result[0] << "\n"
<< " error: "
<< std::abs(result[0] - reference_solution) << "\n";
}
```

--------------------------------

### Example 11: Construct sparse grid from unstructured data with loadUnstructuredDataL2

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianSGExamples11.html

Constructs a global sparse grid from 2000 random data points using a least-squares fit (loadUnstructuredDataL2) with a tolerance of 1.E-4. Requires BLAS, CUDA, or MAGMA acceleration; the example enables accel_gpu_magma and prints the approximation error on 1000 random test points. The code is part of Example 11 and uses a simple exponential model.

```cpp
#endif
cout << "\n---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(4);
cout << "Example 11: construction using unstructured data\n\n";
int const num_inputs = 2;
// using random points to test the error
int const num_test_points = 1000;
std::vector<double> test_points(num_test_points * num_inputs);
std::minstd_rand park_miller(42);
std::uniform_real_distribution<double> domain(-1.0, 1.0);
for(auto &t : test_points) t = domain(park_miller);
// computes the error between the gird surrogate model and the actual model
// using the test points, finds the largest absolute error
auto get_error = [&](TasGrid::TasmanianSparseGrid const &grid,
std::function<void(double const x[], double y[], size_t)> model)->
double{
std::vector<double> grid_result;
grid.evaluateBatch(test_points, grid_result);
double err = 0.0;
for(int i=0; i<num_test_points; i++){
double model_result; // using only one output
model(&test_points[i*num_inputs], &model_result, 0);
err = std::max(err, std::abs(grid_result[i] - model_result));
}
return err;
};
// using a simple model
auto model = [](double const x[], double y[], size_t)->
void{ y[0] = std::exp(-x[0]*x[0] -x[1]-x[1]); };
auto grid = TasGrid::makeGlobalGrid(num_inputs, 1, 4, TasGrid::type_level, TasGrid::rule_clenshawcurtis);
// generate random data for the inputs, and compute the corresponding outputs
int const num_data_points = 2000;
std::vector<double> data_input(num_inputs * num_data_points);
std::vector<double> data_output(num_data_points);
for(auto &d : data_input) d = domain(park_miller);
for(int i=0; i<num_data_points; i++) model(&data_input[i * num_inputs], &data_output[i], 0);
// check if capability is available
if (not grid.isAccelerationAvailable(TasGrid::accel_cpu_blas) and
not grid.isAccelerationAvailable(TasGrid::accel_gpu_cuda) and
not grid.isAccelerationAvailable(TasGrid::accel_gpu_magma)){
cout << "Skipping example 11, BLAS, CUDA, or MAGMA acceleration required.\n";
return;
}
// accel_cpu_blas: works on the CPU and can utilize all available RAM
// accel_gpu_cuda: works on the GPU but it is restricted to the case
// where the data can fit in GPU memory
// accel_gpu_magma: works out-of-core, the data is stored in CPU RAM
// while computations are still done on the GPU
grid.enableAcceleration(TasGrid::accel_gpu_magma);
// constructs the grid, depending on the amount of data data,
// the side of the grid and the enabled acceleration
// this can take significant amount of time
TasGrid::loadUnstructuredDataL2(data_input, data_output, 1.E-4, grid);
cout << "Using construction from unstructured (random) data\n";
cout << " approximatino error = " << get_error(grid, model) << "\n\n";
#ifndef __TASMANIAN_DOXYGEN_SKIP
TasGrid::TasmanianSparseGrid
```

--------------------------------

### Local Polynomial Refinement Example (C++)

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianSGExamples9.html

Full example code for Tasmanian Sparse Grids Example 9. It sets up a sharp model, creates local polynomial grids with rule_localp, and compares classic vs FDS refinement using batch refinement and construction. Requires Tasmanian library and standard headers.

```cpp
cout << "\n---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(4);
cout << "Example 9: comparison between local polynomial refinement strategies\n\n";

int const num_inputs = 2;
// using random points to test the error
int const num_test_points = 1000;
std::vector<double> test_points(num_test_points * num_inputs);
std::minstd_rand park_miller(42);
std::uniform_real_distribution<double> domain(-1.0, 1.0);
for(auto &t : test_points) t = domain(park_miller);

// computes the error between the gird surrogate model and the actual model
// using the test points, finds the largest absolute error
auto get_error = [&](TasGrid::TasmanianSparseGrid const &grid,
                     std::function<void(double const x[], double y[], size_t)> model)->
double{
    std::vector<double> grid_result;
    grid.evaluateBatch(test_points, grid_result);
    double err = 0.0;
    for(int i=0; i<num_test_points; i++){
        double model_result; // using only one output
        model(&test_points[i*num_inputs], &model_result, 0);
        err = std::max(err, std::abs(grid_result[i] - model_result));
    }
    return err;
};

// using a model with sharp transition, the modified rules are not
// suitable for this model, using rule_localp
auto sharp_model = [](double const x[], double y[], size_t)->
void{ y[0] = exp(-x[0]) / (1.0 + 100.0 * exp(-10.0 * x[1])); };

// using maximal order
int const order = -1;
auto grid_classic = TasGrid::makeLocalPolynomialGrid(num_inputs, 1, 2, order,
                                                     TasGrid::rule_localp);
auto grid_fds = TasGrid::copyGrid(grid_classic);

// setup the top rows of the output table
cout << "Using batch refinement:\n"
     << setw(22) << "classic" << setw(22) << "fds\n"
     << setw(8) << "points" << setw(14) << "error"
     << setw(8) << "points" << setw(14) << "error\n";

double const tolerance = 1.E-5;
while((grid_classic.getNumNeeded() > 0) || (grid_fds.getNumNeeded() > 0)){
    TasGrid::loadNeededValues(sharp_model, grid_classic, 4);
    TasGrid::loadNeededValues(sharp_model, grid_fds, 4);

    // output the grids and results at the current stage
    cout << setw(8) << grid_classic.getNumLoaded()
         << setw(14) << get_error(grid_classic, sharp_model)
         << setw(8) << grid_fds.getNumLoaded()
         << setw(14) << get_error(grid_fds, sharp_model) << "\n";

    // setting refinement for each grid
    grid_classic.setSurplusRefinement(tolerance, TasGrid::refine_classic);
    grid_fds.setSurplusRefinement(tolerance, TasGrid::refine_fds);
}

// Repeat the example using construction
auto vector_model = [=](std::vector<double> const &x,
                        std::vector<double> &y,
                        size_t tid)->
void{
    y.resize(1);
    sharp_model(x.data(), y.data(), tid);
};

// reset the grids
grid_classic = TasGrid::makeLocalPolynomialGrid(num_inputs, 1, 2, order,
                                                TasGrid::rule_localp);
grid_fds = TasGrid::copyGrid(grid_classic);

// setup the top rows of the output table
cout << "\nUsing construction:\n"
     << setw(22) << "classic" << setw(22) << "fds\n"
     << setw(8) << "points" << setw(14) << "error"
     << setw(8) << "points" << setw(14) << "error\n";

for(size_t budget = 50; budget < 800; budget += 100){
    TasGrid::constructSurrogate(vector_model, budget, 1, 1, grid_classic,
                                tolerance, TasGrid::refine_classic);
    TasGrid::constructSurrogate(vector_model, budget, 1, 1, grid_fds,
                                tolerance, TasGrid::refine_fds);

    // output the grids and results at the current stage
    cout << setw(8) << grid_classic.getNumLoaded()
         << setw(14) << get_error(grid_classic, sharp_model)
         << setw(8) << grid_fds.getNumLoaded()
         << setw(14) << get_error(grid_fds, sharp_model) << "\n";
}

// The FDS strategy does not win at every step of the way, but at the final
// stage the FDS approach results in fewer nodes while keeping the error
// at near the same magnitude
```

--------------------------------

### ASCII matrix file format example (1x2)

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2InterfaceCLI.html

Example of an ASCII matrix file with one row and two columns. A matrix file may contain only one row or column.

```text
1 2
13.0 14.0
```

--------------------------------

### MPIGridSend and MPIGridRecv example

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianAddonsMPIGridSend.html

Use this MPI example to send a TasmanianSparseGrid from process 0 to process 1. Process 0 creates a global grid with makeGlobalGrid and sends it with MPIGridSend; process 1 receives it with MPIGridRecv. Requires Tasmanian_ENABLE_MPI=ON and the send must be mirrored by a receive on the destination.

```cpp
int me, tag = 0;
MPI_Comm_rank(MPI_COMM_WORLD, &me);
MPI_Status status;
TasGrid::TasmanianSparseGrid grid;
if (me == 0) grid.makeGlobalGrid(3, 4, 5, TasGrid::type_level, TasGrid::rule_clenshawcurtis);
if (me == 0) MPIGridSend(grid, 1, tag, MPI_COMM_WORLD);
else if (me == 1) MPIGridRecv(grid, 0, tag, MPI_COMM_WORLD, &status);
// at this line, process 1 has a grid equivalent to that of process 0
// processes with rank 2 and above do nothing, i.e., they have an empty grid
```

--------------------------------

### ASCII matrix file format example (3x4)

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2InterfaceCLI.html

Example of an ASCII matrix file: the first line contains the number of rows and columns, followed by the matrix entries row by row. This file represents a 3x4 matrix.

```text
3 4
1.0 2.0 3.0 4.0
5.0 6.0 7.0 8.0
9.0 10.0 11.0 12.0
```

--------------------------------

### sparse_grids_example_08()

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianSGExamples8.html

Example function that demonstrates creating a local polynomial sparse grid. It calls the makeLocalPolynomialGrid factory method.

```APIDOC
## sparse_grids_example_08()

### Description
Example function that demonstrates creating a local polynomial sparse grid using the TasmanianSparseGrid factory method makeLocalPolynomialGrid.

### Method
Function call

### Signature
```cpp
void sparse_grids_example_08();
```

### Related Factory Method
TasmanianSparseGrid makeLocalPolynomialGrid(int dimensions, int outputs, int depth, int order=1, TypeOneDRule rule=rule_localp, std::vector< int > const &level_limits=std::vector< int >())

### Parameters
- **dimensions** (int) - Required - Number of dimensions for the grid.
- **outputs** (int) - Required - Number of outputs for the grid.
- **depth** (int) - Required - Depth of the grid.
- **order** (int) - Optional - Order of the local polynomial rule. Default is 1.
- **rule** (TypeOneDRule) - Optional - One-dimensional rule. Default is rule_localp.
- **level_limits** (std::vector<int> const&) - Optional - Level limits for the grid. Default is empty vector.

### Description of Factory Method
Factory method, creates a new grid and calls TasmanianSparseGrid::makeLocalPolynomialGrid().

### Definition
TasmanianSparseGrid.hpp:2316
```

--------------------------------

### sparse_grids_example_10()

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianSGExamples10.html

Example function that demonstrates creating a local polynomial sparse grid. It likely calls the makeLocalPolynomialGrid factory method.

```APIDOC
## sparse_grids_example_10()

### Description
Example function that demonstrates creating a local polynomial sparse grid. It likely calls the makeLocalPolynomialGrid factory method.

### Method
Function call

### Signature
```cpp
void sparse_grids_example_10();
```

### Related Factory Method
TasmanianSparseGrid makeLocalPolynomialGrid(int dimensions, int outputs, int depth, int order=1, TypeOneDRule rule=rule_localp, std::vector< int > const &level_limits=std::vector< int >())

### Description
Factory method, creates a new grid and calls TasmanianSparseGrid::makeLocalPolynomialGrid().

### Parameters
- **dimensions** (int) - Required - Number of dimensions.
- **outputs** (int) - Required - Number of outputs.
- **depth** (int) - Required - Depth of the grid.
- **order** (int) - Optional - Order of the polynomial, default is 1.
- **rule** (TypeOneDRule) - Optional - One-dimensional rule, default is rule_localp.
- **level_limits** (std::vector<int> const&) - Optional - Level limits, default is empty vector.

### Returns
- **TasmanianSparseGrid** - The created sparse grid.

### Definition
TasmanianSparseGrid.hpp:2316
```

--------------------------------

### Create a sparse grid instance

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2InterfacePython.html

Example for creating an instance of a sparse grid. The grid is created with Tasmanian.SparseGrid() and then configured with makeGlobalGrid(...).

```python
import Tasmanian
grid = Tasmanian.SparseGrid()
grid.makeGlobalGrid(...)
```

--------------------------------

### DREAM Example 3: signal decomposition with bi-modal posterior (C++)

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianDREAMExamples3.html

C++ example for DREAM Example 3: signal decomposition with a bi-modal posterior. It sets up a sparse grid surrogate, defines a likelihood, runs Tasmanian DREAM sampling, and splits the bi-modal history at frequency 6.5 to separate low and high frequency components. Requires the Tasmanian DREAM and Sparse Grid libraries.

```cpp
srand((int) time(nullptr));
```

```cpp
cout << "\n" << "---------------------------------------------------------------------------------------------------\n";
cout << "EXAMPLE 3: set the inference problem: identify x_0 and x_1 model parameters\n"
<< " from data (noise free example)\n"
<< " model: f(x) = sin(x_0*M_PI*t + x_1),\n"
<< " data: d = sin(5*M_PI*t + 0.3*M_PI) + sin(10*M_PI*t + 0.1*M_PI)\n"
<< " compared to Example 2, the data is a superposition of two signals\n"
<< " and the posterior is multi-modal\n"
<< " -- problem setup --\n"
<< " t in [0,1], t is discretized with 32 equidistant nodes\n"
<< " the likelihood is exp(- 16 * (f(x) - d)^2)\n"
<< " using a sparse grid to interpolate the model\n"
<< " NOTE: 16 = 32/2 corresponds to the discretization error in t\n\n";
```

```cpp
constexpr double pi = 3.14159265358979323846;
```

```cpp
int num_chains = 500;
int num_burnup_iterations = 1000;
int num_sample_iterations = 300;
int num_discrete_nodes = 32;
```

```cpp
auto model = [&](double x0, double x1, std::vector<double> &data)->
void{
double dt = 1.0 / ((double) data.size());
double t = 0.5 * dt;
for(auto &d : data){
d = std::sin(x0 * pi * t + x1);
t += dt;
}
};
```

```cpp
std::vector<double> signal1(num_discrete_nodes),
signal2(num_discrete_nodes),
data(num_discrete_nodes);
model( 5.0, 0.3 * pi, signal1);
model(10.0, 0.1 * pi, signal2);
std::transform(signal1.begin(), signal1.end(), signal2.begin(),
data.begin(), std::plus<double>());
```

```cpp
auto grid = TasGrid::makeSequenceGrid(2, num_discrete_nodes, 30, // 30-th order polynomial
TasGrid::type_iptotal, TasGrid::rule_leja);
```

```cpp
std::vector<double> domain_a = { 1.0, -0.1};
std::vector<double> domain_b = {12.0, 1.7};
grid.setDomainTransform(domain_a, domain_b);
```

```cpp
TasGrid::loadNeededPoints<TasGrid::mode_sequential>(
[&](std::vector<double> const &x, std::vector<double> &y, size_t)->
void{
model(x[0], x[1], y);
},
grid, 1);
```

```cpp
// when working with grids with large N, acceleration becomes very useful
// if BLAS is enabled on compile time, the grid will use BLAS by default
// GPU acceleration can be enabled here using:
// grid.enableAcceleration(TasGrid::accel_gpu_cuda);
// grid.setGPUID(0);
```

```cpp
TasDREAM::LikelihoodGaussIsotropic likely(1.0 / ((double) num_discrete_nodes), data);
TasDREAM::TasmanianDREAM state(num_chains, grid); // get the dimensions from the sparse grid
state.setState(TasDREAM::genUniformSamples(domain_a, domain_b, num_chains));
```

```cpp
TasDREAM::SampleDREAM<TasDREAM::logform>
(num_burnup_iterations, num_sample_iterations,
TasDREAM::posterior<TasDREAM::logform>
(grid, // provide the surrogate model
likely, // provide the likelihood
TasDREAM::uniform_prior), // assume non-informative prior
grid.getDomainInside(),
state,
TasDREAM::dist_gaussian, 0.01, // independent update of magnitude 0.01
TasDREAM::const_percent<90> // use 90% differential update
);
```

```cpp
const std::vector<double> &history = state.getHistory();
```

```cpp
double frequency_low = 0.0, correction_low = 0.0;
double frequency_high = 0.0, correction_high = 0.0;
int num_low = 0, num_high = 0;
for(size_t i=0; i<history.size(); i+=2){ // 2 parameters, the samples in a stride of 2
if (history[i] < 6.5){
frequency_low += history[i];
correction_low += history[i+1];
num_low++;
}else{
frequency_high += history[i];
correction_high += history[i+1];
num_high++;
}
}
```

--------------------------------

### dream_example_02: Sparse Grid Likelihood Calibration

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianDREAMExamples2.html

This C++ code implements the DREAM Example 2, which performs Bayesian parameter calibration using a sparse grid approximation to the likelihood. It sets up a model and likelihood, constructs a sparse grid, and runs the DREAM sampling algorithm. The example uses a 10th-order polynomial global grid and then a 30th-order sequence grid for comparison.

```cpp
#endif
// using the default random engine, but must reset the random number generator
srand((int) time(nullptr));
// Example 2
cout << "\n" << "---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(5);
cout << "EXAMPLE 2: set the inference problem: identify model parameters x_0 and x_1\n"
<< " from data (noise free example)\n"
<< " model: f(x) = sin(x_0*M_PI*t + x_1), data: d = sin(M_PI*t + 0.3*M_PI)\n"
<< " t in [0,1], t is discretized with 32 equidistant nodes\n"
<< " the likelihood is exp(- 16 * (f(x) - d)^2)\n"
<< " using a sparse grid to interpolate the likelihood\n"
<< " NOTE: 16 = 32/2 corresponds to the discretization error in t\n\n";
constexpr double pi = 3.14159265358979323846;
int num_dimensions = 2;
int num_chains = 100;
int num_burnup_iterations = 3000;
int num_sample_iterations = 100;
// the total number of samples is num_chains * num_iterations = 100 * 100 = 10,000
int num_discrete_nodes = 32;
// create a lambda function that represents the model,
// normally this would be a call to an external code
auto model = [&](double x0, double x1, std::vector<double> &data)->
void{
double dt = 1.0 / ((double) data.size());
double t = 0.5 * dt;
for(auto &d : data){
d = std::sin(x0 * pi * t + x1);
t += dt;
}
};
// create a lambda function that represents the likelihood
// the Tasmanian LikelihoodGaussIsotropic class can be used here,
// but the objective of the example is to use sparse grid interpolant for the likelihood
// thus we assume that the likelihood formula is potentially much more complex
auto likelihood = [&](double scale,
std::vector<double> &model_values,
std::vector<double> &data)->
double{
double likelihood_value = 0.0;
for(size_t j=0; j<model_values.size(); j++){
likelihood_value += (model_values[j] - data[j]) * (model_values[j] - data[j]);
}
return - 0.5 * scale * likelihood_value;
};
// generate the data, normally this will come from an external source
std::vector<double> data(num_discrete_nodes);
model(1.0, 0.3 * pi, data); // "true" values: x0 = 1.0, x1 = 0.3 * pi
// construct a sparse grid approximation for the log of the likelihood
// the sampling will be done in logform
// for details consult the documentation of the sparse grids module
// define the domain over which the sparse grid will be constructed
std::vector<double> domain_a = {0.5, -0.1}, domain_b = {8.0, 1.7};
auto grid = TasGrid::makeGlobalGrid(num_dimensions, 1, 10, // using 10-th order polynomial
TasGrid::type_iptotal, TasGrid::rule_clenshawcurtis);
grid.setDomainTransform(domain_a, domain_b);
TasGrid::loadNeededPoints<TasGrid::mode_sequential>(
[&](std::vector<double> const &x, std::vector<double> &y, size_t)->void{
std::vector<double> model_at_point(num_discrete_nodes);
model(x[0], x[1], model_at_point);
y[0] = likelihood((double) num_discrete_nodes, model_at_point, data);
},
grid, 1);
TasDREAM::TasmanianDREAM state(num_chains, grid); // get the dimensions from the sparse grid
std::vector<double> init_chains = TasDREAM::genUniformSamples(domain_a, domain_b, num_chains);
state.setState(init_chains); // use chains distributed uniformly over the domain
// Call to Tasmanian DREAM Sampling algorithm
TasDREAM::SampleDREAM<TasDREAM::logform>
(num_burnup_iterations, num_sample_iterations,
TasDREAM::posterior(grid, TasDREAM::uniform_prior),
grid.getDomainInside(),
state,
TasDREAM::dist_uniform, 0.5, // independent update of magnitude 0.5
TasDREAM::const_percent<90> // use 90% differential update
);
std::vector<double> expectation, variance;
state.getHistoryMeanVariance(expectation, variance);
cout << "Inferred values (using 10th order polynomial sparse grid):\n";
cout << " frequency:" << setw(12) << std::fixed << expectation[0]
<< " error:" << setw(12) << std::scientific << std::abs(expectation[0] - 1.0) << "\n";
cout << " correction:" << setw(12) << std::fixed << expectation[1]
<< " error:" << setw(12) << std::scientific << std::abs(expectation[1] - 0.3 * pi)
<< "\n\n";
// Solve the same example, but switch to 30-th order polynomial with a Sequence grid
// rebuild the sparse grid interpolant
// use 30-th order polynomial and the faster Sequence grid
grid.makeSequenceGrid(2, 1, 30, TasGrid::type_iptotal, TasGrid::rule_leja);
grid.setDomainTransform(domain_a, domain_b);
TasGrid::loadNeededPoints<TasGrid::mode_sequential>(
```

--------------------------------

### Sparse Grids Example 2: Integral with Gauss-Patterson nodes (C++)

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianSGExamples2.html

This code computes the integration of f(x,y) = exp(-x^2)*cos(y) over [-5,5] x [-2,3] using a global grid with total-degree polynomial space and Gauss-Patterson rule. It sets a domain transform, obtains quadrature points and weights, and evaluates the weighted sum. The example runs for exactness values 20 and 40, printing the integral and error compared to the exact value.

```cpp
#endif
cout << "\n---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(17);
cout << "Example 2: integrate f(x,y) = exp(-x^2) * cos(y) over [-5,5] x [-2,3]\n"
<< " using Gauss-Patterson nodes and total degree polynomial space\n";
int dimension = 2;
int exactness = 20;
// the type_qptotal will guarantee exact integral for all polynomials with degree 20 or less
auto grid = TasGrid::makeGlobalGrid(dimension, 0, exactness,
TasGrid::type_qptotal, TasGrid::rule_gausspatterson);
grid.setDomainTransform({-5.0, -2.0}, {5.0, 3.0}); // set the non-canonical domain
auto points = grid.getPoints();
auto weights = grid.getQuadratureWeights();
double I = 0.0;
for(int i=0; i<grid.getNumPoints(); i++){
double x = points[i*dimension];
double y = points[i*dimension+1];
I += weights[i] * std::exp(-x*x) * std::cos(y);
}
double exact = 1.861816427518322e+00;
double E = std::abs(exact - I);
cout << " at polynomial exactness: " << exactness
<< "\n the grid has: " << weights.size()
<< "\n integral: " << I
<< "\n error: " << E << "\n\n";
exactness = 40;
grid = TGrid::makeGlobalGrid(dimension, 0, exactness,
TasGrid::type_qptotal, TasGrid::rule_gausspatterson);
grid.setDomainTransform({-5.0, -2.0}, {5.0, 3.0}); // set the non-canonical domain
points = grid.getPoints();
weights = grid.getQuadratureWeights();
I = 0.0;
for(int i=0; i<grid.getNumPoints(); i++){
double x = points[i*dimension];
double y = points[i*dimension+1];
I += weights[i] * std::exp(-x*x) * std::cos(y);
}
E = std::abs(exact - I);
cout << " at polynomial exactness: " << exactness
<< "\n the grid has: " << weights.size()
<< "\n integral: " << I
<< "\n error: " << E << endl;
```

--------------------------------

### CMake find_package guidance options

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2Installation.html

Additional CMake commands to guide the find_package() modules. Use these to specify the Python interpreter, CUDA compiler, HIP search path, Fortran compiler, and oneMKL installation path.

```cmake
-D Python_EXECUTABLE:PATH (specify the Python interpreter)
-D CMAKE_CUDA_COMPILER:PATH (specify the CUDA nvcc compiler)
-D Tasmanian_ENABLE_ROCM:PATH (specify the search path for HIP)
-D CMAKE_Fortran_COMPILER:PATH (specify the Fortran compiler)
-D MKLROOT:PATH (specify the path to the oneMKL installation)
```

--------------------------------

### constructSurrogate simple exponential example

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianAddonsConstruct.html

Constructs a surrogate to a simple exponential function with a budget of 100 points, running on 4 threads with one sample per thread. The grid must be initialized with the appropriate type, dimensions, outputs, and enough initial points; the procedure terminates after 100 samples or when the tolerance of 1.E-6 is reached.

```cpp
auto grid = TasGrid::makeLocalPolynomialGrid(2, ..); // sync the dimensions with number of model inputs
TasGrid::constructSurrogate([&](std::vector<double> const &x, std::vector<double> &y)
->void{ y[0] = std::exp(x[0] + x[1]); },
100, 4, 1, grid, 1.E-6, TasGrid::refine_classic);
```

--------------------------------

### Install Tasmanian with pip

Source: https://ornl.github.io/TASMANIAN/index.html

Install the Tasmanian Python package via pip. Requires the latest pip and installs to the user site-packages.

```bash
python3 -m pip install Tasmanian --user
```

--------------------------------

### Print version and info with tasgrid -version, -v, or -info

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2InterfaceCLI.html

Prints the version of the library, the available acceleration options, and (if CUDA is enabled) the visible CUDA devices. The three commands are equivalent.

```bash
./tasgrid -version
```

```bash
./tasgrid -v
```

```bash
./tasgrid -info
```

--------------------------------

### Build and install Tasmanian with CMake

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2Installation.html

Out-of-source CMake build and installation workflow for Tasmanian. Requires CMake 3.19 (3.21 for ROCm). The commands create a Build directory, configure with options, compile, run tests, install, and verify the installation.

```bash
mkdir Build
cd Build
cmake <options> <path-to-Tasmanian-source>
make
make test
make install
make test_install
```

--------------------------------

### Python Pip install commands

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2Installation.html

Install Tasmanian via pip for user or virtual environment installations. Note that --no-cache-dir is required to change build options.

```bash
python3 -m pip install Tasmanian --user (user installation)
python3 -m pip install Tasmanian (virtual env installation)
```

--------------------------------

### Run Gradient Descent with analytic gradient, surrogate gradient, and box projection

Source: https://ornl.github.io/TASMANIAN/stable/group__TasmanianOPTExamples2.html

The example begins with an objective and gradient lambda for f(x,y)=2(x-1)^2+(y-2)^2/2, executes TasOptimization::GradientDescent, then repeats using a sequence grid surrogate's differentiate, and finally performs box-constrained descent with a custom projection and step-size coefficients.

```cpp
// Example 2:
cout << "\n" << "---------------------------------------------------------------------------------------------------\n";
cout << std::scientific; cout.precision(5);
cout << "EXAMPLE 2: use Gradient Descent algorithms on a simple quadratic\n"
<< " f(x,y) = 2.0 * (x - 1.0) * (x - 1.0) + (y - 2.0) * (y - 2.0) / 2.0\n"
<< " See the comments in example_optimization_02.cpp\n\n";
cout << "Exact solution is at (1.0, 2.0)\n\n";
auto objective = [](std::vector<double> const &x)
->double{
return 2.0 * (x[0] - 1.0) * (x[0] - 1.0) + (x[1] - 2.0) * (x[1] - 2.0) / 2.0;
};
// the gradient descent methods require the gradient of the objective
auto gradient = [](std::vector<double> const &x, std::vector<double> &grad)
->void{
grad[0] = 4.0 * (x[0] - 1.0);
grad[1] = x[1] - 2.0;
};
int num_dimensions = 2;
double initial_stepsize = 0.0;
// the {0.0, 0.0} is a vector with initial state, e.g., initial guess
// coupling particle swarm and gradient descent can happen with
// auto gradient_state =
// TasOptimization::GradientDescentState(particle_state.getBestPosition(), 0.0);
TasOptimization::GradientDescentState state({0.0, 0.0}, initial_stepsize);
int max_iterations = 200;
double tolerance = 1.E-3;
TasOptimization::OptimizationStatus status =
TasOptimization::GradientDescent(gradient, 1.0/8.0, max_iterations, tolerance, state);
std::vector<double> opt_x = state.getX();
cout << "Using the objective function\n"
<< "aiming at tolerance 1.e-3\n"
<< "performed iterations: " << status.performed_iterations << "\n"
<< "best value for x = " << opt_x[0] << "\n"
<< "best value for y = " << opt_x[1] << "\n\n";
// same problem, but using a surrogate model
// construct the surrogate
auto grid = TasGrid::makeSequenceGrid(num_dimensions, 1, 2,
TasGrid::type_iptotal, TasGrid::rule_leja);
std::vector<double> points = grid.getNeededPoints();
std::vector<double> values(grid.getNumNeeded());
for(int i=0; i<grid.getNumNeeded(); i++){
values[i] = objective(std::vector<double>(points.begin() + 2 * i,
points.begin() + 2 * i + 2));
}
grid.loadNeededValues(values);
// reset the state to zero
state = TasOptimization::GradientDescentState({0.0, 0.0}, initial_stepsize);
status = TasOptimization::GradientDescent(
[&](std::vector<double> const &x, std::vector<double> &grad)
->void{
grid.differentiate(x, grad);
},
1.0/8.0, max_iterations, tolerance, state);
opt_x = state.getX();
// actually, since the objective is quadratic, the surrogate is an exact match
// the example demonstrates syntax more than anything else
cout << "Using the surrogate function\n"
<< "aiming at tolerance 1.e-3\n"
<< "performed iterations: " << status.performed_iterations << "\n"
<< "best value for x = " << opt_x[0] << "\n"
<< "best value for y = " << opt_x[1] << "\n\n";
// finding optimum inside a box
// define the projection operator of a vector into a box
// restriction function is probably better name here, since the projection
// is not orthogonal
// Note: that such restriction is necessary when working with a sparse grid
// since the surrogates are usually restricted to a hypercube box.
auto projection = [](const std::vector<double> &x, std::vector<double> &p) {
p[0] = std::min(std::max(x[0], 0.0), 0.5);
p[1] = std::min(std::max(x[1], 0.0), 1.5);
};
// reset the state and take the initial_stepsize to 1.0
state = TasOptimization::GradientDescentState({0.0, 0.0}, 1.0);
double increase_coeff = 1.25;
double decrease_coeff = 1.25;
status = TasOptimization::GradientDescent(objective, gradient, projection,
increase_coeff, decrease_coeff,
max_iterations, tolerance, state);
opt_x = state.getX();
// actually, since the objective is quadratic, the surrogate is an exact match
// the problem demonstrates syntax more than anything else
cout << "Using the projection inside of a box\n"
<< "the solution is at the edge of the box at (0.5, 1.5)\n"
<< "performed iterations: " << status.performed_iterations << "\n"
<< "best value for x = " << opt_x[0] << "\n"
<< "best value for y = " << opt_x[1] << "\n\n";
```

--------------------------------

### Install Tasmanian with Spack

Source: https://ornl.github.io/TASMANIAN/stable/md_Doxygen_2Installation.html

Installs Tasmanian and selected features (OpenMP, BLAS, CUDA, MAGMA, Python, Fortran) using Spack's package manager.

```bash
spack install tasmanian+openmp+blas+cuda+magma+python+fortran
```