diff --git a/Tutorial1_Introduction/.ipynb_checkpoints/Introduction-to-BO-checkpoint.ipynb b/Tutorial1_Introduction/.ipynb_checkpoints/Introduction-to-BO-checkpoint.ipynb
new file mode 100755
index 0000000000000000000000000000000000000000..bde76177e47e7661b7edfb1d95f740488e852e6c
--- /dev/null
+++ b/Tutorial1_Introduction/.ipynb_checkpoints/Introduction-to-BO-checkpoint.ipynb
@@ -0,0 +1,363 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "QHJs6azjC8Il"
+ },
+ "source": [
+ "# Tutorial 1: Introduction to active learning with Bayesian Optimization"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "i_kyWo5KFGQj"
+ },
+ "source": [
+ "## The principles Bayesian Optimization"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "On binder, we can just render any image (from the git repo or generated by BOSS) using markdown like so:\n",
+ ""
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "hlWKlGryC8In",
+ "lines_to_next_cell": 2
+ },
+ "source": [
+ "## 1D example of Bayesian Optimization\n",
+ "A one dimensional problem is the simplest case to solve with BOSS. By looking at the plotted GP models step-by-step,\n",
+ "it is transparent what is going on in BOSS. Studying a 1D case is also often a part of solving a bigger optimization\n",
+ "problem, as short 1D searches of each simulation variable separately gives important insight on how to set up the large problem. In this demo the target function is"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "Brbyz6tyC8Io"
+ },
+ "source": [
+ "This is a Gaussian added to a sine wave, which makes up a fairly simple but not trivial non-periodic 1D function. The\n",
+ "domain of the variable x is [0, 7] by defition, and because of the sine, the values of f(x) are known to be roughly within [ô€€€1; 1]. Since we have an anlytic expression for the true function we can go ahead and plot it:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "metadata": {
+ "colab": {
+ "base_uri": "https://localhost:8080/"
+ },
+ "executionInfo": {
+ "elapsed": 6370,
+ "status": "ok",
+ "timestamp": 1630579817024,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "Q1wK8_INC8Ip",
+ "lines_to_next_cell": 1,
+ "outputId": "dea01038-08b2-48c0-be4c-bf888a532db2"
+ },
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import matplotlib.pyplot as plt\n",
+ "%matplotlib inline\n",
+ "from boss.bo.bo_main import BOMain\n",
+ "from boss.pp.pp_main import PPMain"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 225,
+ "status": "ok",
+ "timestamp": 1630579841548,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "YemkoQ0LC8Ip",
+ "lines_to_next_cell": 1
+ },
+ "outputs": [],
+ "source": [
+ "def f(X):\n",
+ " x = X[0, 0]\n",
+ " y = np.sin(x) + 1.5*np.exp(-(x - 4.3)**2)\n",
+ " return y"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "metadata": {
+ "colab": {
+ "base_uri": "https://localhost:8080/",
+ "height": 313
+ },
+ "executionInfo": {
+ "elapsed": 709,
+ "status": "ok",
+ "timestamp": 1630579843739,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "5VUbwTOtC8Iq",
+ "lines_to_next_cell": 2,
+ "outputId": "5e914b18-3cd1-41f1-f7f9-19975c3bb888"
+ },
+ "outputs": [],
+ "source": [
+ "bounds = np.array([[0., 7.]])\n",
+ "fig, ax = plt.subplots()\n",
+ "x = np.linspace(bounds[0, 0], bounds[0, 1], 100)\n",
+ "y_true = np.array([f(np.atleast_2d(xi)) for xi in x])\n",
+ "ax.plot(x, y_true)\n",
+ "ax.set_xlabel('x')\n",
+ "ax.set_ylabel('f(x)')\n",
+ "ax.set_title('True function')\n",
+ "plt.show()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "dFzs0bD08gep"
+ },
+ "source": [
+ "# Running instructions"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "8EA4GyoCC8Ir",
+ "lines_to_next_cell": 2
+ },
+ "source": [
+ "The first step to running BOSS typically consists of defining an objective function and\n",
+ "the optimization bounds, where the latter should be specified as a hypercube.\n",
+ "For the Forrester problem, we define the function and bounds as follows"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "5ndWs_1kC8Ir",
+ "lines_to_next_cell": 2
+ },
+ "source": [
+ "Note that BOSS expects the objective function to take a single 2D numpy array \n",
+ "as argument and return a scalar value (this behaviour can be modified).\n",
+ "Next, we import BOMain, which will be used to launch and configure the optimization.\n",
+ "When creating this object we can supply any number of BOSS *keywords*,\n",
+ "these are used to provide essential input information and modify BOSS's behavior.\n",
+ "In the following, only a minimal set of keywords are provided for brevity."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 208,
+ "status": "ok",
+ "timestamp": 1630579846701,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "FbkFz5egC8It"
+ },
+ "outputs": [],
+ "source": [
+ "bo = BOMain(\n",
+ " f, \n",
+ " bounds, \n",
+ " yrange=[-1, 1],\n",
+ " kernel='rbf',\n",
+ " initpts=5,\n",
+ " iterpts=5,\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "eJ9_AjA3C8It"
+ },
+ "source": [
+ "Consider the variable boundary and estimated target function values range (keyword yrange). Are they set properly?\n",
+ "The kernel type should be chosen to be the non-periodic rbf (why?)."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 1627,
+ "status": "ok",
+ "timestamp": 1630579850682,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "Wqmj8WqrC8Iu"
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting file 'boss.out'\n",
+ "warning: overwriting file 'boss.rst'\n"
+ ]
+ }
+ ],
+ "source": [
+ "res = bo.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "zKy8LL1oC8Iu"
+ },
+ "source": [
+ "# Results and Analysis\n",
+ "Let's visualize the quality of the GP regression fit and the estimated minimum using the information provided the by the results object."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 10,
+ "metadata": {
+ "colab": {
+ "base_uri": "https://localhost:8080/",
+ "height": 295
+ },
+ "executionInfo": {
+ "elapsed": 764,
+ "status": "ok",
+ "timestamp": 1630579855186,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "MBWdGerPC8Iu",
+ "outputId": "d3697054-c4cf-4a79-84bb-1b5baf49ac3e"
+ },
+ "outputs": [],
+ "source": [
+ "import matplotlib.pyplot as plt\n",
+ "fig, ax = plt.subplots()\n",
+ "ax.plot(x, y_true, color='tab:blue', label='True') # Plot the true function.\n",
+ "y_fit = res.f(np.atleast_2d(x).T)\n",
+ "ax.plot(x, y_fit, color='tab:red', ls='--', label='Fit') # Plot the GP regression fit.\n",
+ "# Plot the estimated minimum. Note the use the xmin, and fmin members of the results object.\n",
+ "ax.plot(res.xmin, res.fmin, 'ro', label='Minimum')\n",
+ "ax.set_xlabel('x')\n",
+ "ax.set_ylabel('f(x)')\n",
+ "ax.set_title('BO minimization of Forrester function')\n",
+ "ax.legend()\n",
+ "plt.show()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 11,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 10651,
+ "status": "ok",
+ "timestamp": 1630579868087,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "jtb44sk8EfoN"
+ },
+ "outputs": [],
+ "source": [
+ "pp = PPMain(res, pp_models=True, pp_acq_funcs=True, pp_truef_npts=100) # add some more pp keywords here if needed\n",
+ "pp.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We can include images generated by BOSS the same way we do uploaded images \n",
+ "\n",
+ ""
+ ]
+ }
+ ],
+ "metadata": {
+ "colab": {
+ "collapsed_sections": [],
+ "name": "tutorial1.ipynb",
+ "provenance": []
+ },
+ "jupytext": {
+ "cell_metadata_filter": "-all",
+ "main_language": "python",
+ "notebook_metadata_filter": "-all"
+ },
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.6.3"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/Tutorial1_Introduction/Introduction-to-BO.ipynb b/Tutorial1_Introduction/Introduction-to-BO.ipynb
new file mode 100755
index 0000000000000000000000000000000000000000..2a08fcd47302b0513073ba65350ba93b7ba191a6
--- /dev/null
+++ b/Tutorial1_Introduction/Introduction-to-BO.ipynb
@@ -0,0 +1,345 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "QHJs6azjC8Il"
+ },
+ "source": [
+ "# Tutorial 1: Introduction to active learning with Bayesian Optimization"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Bayesian Optimization (BO) is frequently based on Gaussian Process regression (GPR), a probabilistic regression tool. GPR is used to fit *the statistically most likely function*, given data. A Gaussian process (GP) had two forms: the **prior** (encoding any prior belief about the function) and the **posterior** (see image below, from [Rasmussen & Williams, MIT Press (2006)](http://www.gaussianprocess.org/gpml/)). Given the GP prior and some data points (observations), the Bayes' rule is used to compute the GP posterior. **GP posterior mean** is statistically the most likely function that fits the data. **GP posterior variance** is the measure of confidence on the GP mean, with high confidence where variance vanishes.\n",
+ "\n",
+ "<br><img src=\"../figures/GPR_example.png\" width=\"650px\"/><br>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "In BO, GPR is combined with an acquisition function to sample further data \"on-the-fly\". Acquisition functions sample data points with high information content, which are added to the dataset to iteratively refine the model until convergence. Here, model training and refinement is carried out concurrently with dataset building, which makes this an *active learning* approach. The acquisition functions are typically set up to find the global minimum of the unknown function in as few data points as possible.\n",
+ "\n",
+ "In this tutorial, we use the free [Bayesian Optimization Structure Search (BOSS)](https://gitlab.com/cest-group/boss) code. See the [documentation website](https://cest-group.gitlab.io/boss/) for tutorials, manual and keywords. "
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "hlWKlGryC8In",
+ "lines_to_next_cell": 2
+ },
+ "source": [
+ "## 1D example of Bayesian Optimization\n",
+ "Let us consider an unknown, non-periodic, 1-dimensional function, and infer its global minimum with BO.\n",
+ "First, we load the required python packages."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 65,
+ "metadata": {
+ "colab": {
+ "base_uri": "https://localhost:8080/"
+ },
+ "executionInfo": {
+ "elapsed": 6370,
+ "status": "ok",
+ "timestamp": 1630579817024,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "Q1wK8_INC8Ip",
+ "lines_to_next_cell": 1,
+ "outputId": "dea01038-08b2-48c0-be4c-bf888a532db2"
+ },
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import matplotlib.pyplot as plt\n",
+ "%matplotlib inline\n",
+ "\n",
+ "from boss.bo.bo_main import BOMain\n",
+ "from boss.pp.pp_main import PPMain"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Next, let's define a simple but non-trivial 1D function by adding a Gaussian to a sine wave. We chose the domain of the variable x to be [0, 7], and because of the sine, the values of f(x) are known to be roughly within [-1; 1]. Since we have an analytic expression for the true function we can go ahead and plot it:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 66,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 225,
+ "status": "ok",
+ "timestamp": 1630579841548,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "YemkoQ0LC8Ip",
+ "lines_to_next_cell": 1
+ },
+ "outputs": [
+ {
+ "data": {
+ "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEWCAYAAAB42tAoAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAA210lEQVR4nO3deXxU9b3/8dcnk41AIGSFBEJAkrDJGgE3RAEFBNHWutelWq9t7dXb9rZa29r13ra3rXbxVr3W1lYrtVoVBUQWFTeWsGdhCWvIngBZCVnm8/sjg7+ISQiQmTMz+TwfjzyY5cw5b9Dkk/NdRVUxxhhjOhPidABjjDH+zQqFMcaYLlmhMMYY0yUrFMYYY7pkhcIYY0yXrFAYY4zpkhUKY86SiFwnIoUiUicik3x43VtF5G1fXc8YKxSmVxCR+0UkW0ROiMhfTnlvpoi4PT/w60TksIi8JCIXnOa0vwLuV9V+qrrFS7nTRERFJPTka6r6gqpe6Y3rGdMRKxSmtygGfgo829n7qtoPiAamAzuB90VkVhfnHAbk9mhKY/yQFQrTK6jqv1T1NaDqNMepqh5W1R8AzwC/OPUYEYkQkTrABWwTkb2e11VERrY77i8i8lPP45meO5Vviki5iJSIyF3tju0jIr8WkYMiUi0iH4hIH2Ct55BjnrudC0XkThH5oN1nLxKRjZ7PbRSRi9q9966I/EREPhSRWhF5W0Tiz/xf0PRmViiM6dy/gMki0rf9i6p6wnP3ATBBVc/r5vkGAQOAFOBu4AkRGeh571fAFOAiIBb4NuAGZnjej/E0cX3c/oQiEgssBX4HxAG/AZaKSFy7w24B7gISgXDgW93MawxghcKYrhQDAsT00PmagR+rarOqLgPqgEwRCQG+BDygqkWq2qqqH6nqiW6c82pgj6r+TVVbVPVF2prNFrY75s+qultVjwMvARN76O9jegkrFMZ0LgVQ4FgPna9KVVvaPW8A+gHxQCSw9yzOmQwcPOW1g7RlP6m0g2sa021WKIzp3HXAZlWt7+bxDUBUu+eDuvm5SqAR6KgJ63TLOxfT1qneXipQ1M1rG3NaVihMryAioSISSVsHtEtEItsPOW13nIhIiog8CtwDfPcMLrMVuEVEXCIyF7isOx9SVTdto7F+IyLJns9fKCIRQAVtfRUjOvn4MiBDRG7x/B1vBMYAb55BbmO6ZIXC9BbfA44DDwG3eR5/r937yZ6RTHXARuB8YKaqnsnEtgdo6xs4BtwKvHYGn/0WsMNz7SO0jbYKUdUG4GfAhyJyTESmt/+QqlYBC4Bv0jai69vAAlWtPINrG9MlsY2LjDHGdMXuKIwxxnTJCoUxxpguWaEwxhjTJSsUxhhjuvSZ4YGBLj4+XtPS0pyOYYwxAWXTpk2VqprQ0XtBVyjS0tLIzs52OoYxxgQUETl1hv8nrOnJGGNMl6xQGGOM6ZIVCmOMMV2yQmGMMaZLViiMMcZ0ydFCISLPeraFzOnkfRGR34lIgYhsF5HJvs5ojDG9ndN3FH8B5nbx/jwg3fN1L/BHH2QyxhjTjqPzKFR1rYikdXHIIuCv2rbE7ToRiRGRwapa4puEvU9jcysF5XXkl9RQfKyRqHAXfSNCGRgVxtThscT1i3A6ojHGx/x9wl0KUNju+WHPa58qFCJyL213HKSmpvosXLBodSvv7irnb+sO8sGeSlrcHS89LwLnpwxgZmYit01LJbF/pI+TGmOc4O+FoltU9WngaYCsrCzbYKOb3G7lpexCfr+mgKJjx0mMjuDuS4YzfkgMowZHkxobxYkWN/UnWig6dpwP9lSydncFf1izh6fe28ut04Zx38wRJEZbwTAmmPl7oSgChrZ7PgTbC7hHHKis56F/bWfdviNMGTaQR64ezZwxSYS5Pt1tFeYKoV9EKEn9I5mcOpB/n5XOgcp6/vBOAc99fIC/bzjIw/NGc/uFwxARh/42xhhvcroz+3SWALd7Rj9NB6qtf+Lc/e3jA1z1+Fpyi2v4xefP5+X7LmT++YM/UyQ6kxbfl199YQKrv3EZF46I49EluXz5r9kcqW/ycnJjjBMc3QpVRF4EZgLxQBnwKBAGoKpPStuvqH+gbWRUA3CXqna54l9WVpbaooAda3UrP1uaz7Mf7ufyzAR+/vnxJJ1jP4Oq8pePDvDfy3YSExXGU1+cwqTUgT2U2BjjKyKySVWzOnwv2PbMtkLRscbmVh5cvJW3cku56+I0vnf1GFwhPddUlFdcw33Pb6Ky7gRPfzGLS9Lje+zcxhjv66pQ+HvTk+kBjc2t3PnnDazIK+X7C8bw6MKxPVokAMYk9+fl+y4kNTaKL/1lI2/llPbo+Y0xzrFCEeRa3coDi7ewfv8RHrthIndfMtxr10rsH8nie6czNqU/X31hE29uL/batYwxvmOFIoipKt97bQcrcst4dMEYrp2U4vVrxkSF88I908gaFss3/rGNj/dWef2axhjvskIRxB5ftYcXNxTytcvP486LvXcncaqo8FCevn0KqXFR3Pu3bHaX1frs2saYnmeFIkitzi/jt6v3cP2UIXzrykyfXz8mKpy/3HUBfcJc3PnsBkqrG32ewRjTM6xQBKHDRxv4xkvbGDO4Pz+9dpxjE+GGDIziz3ddQPXxZr7ywiaaWtyO5DDGnBsrFEGmqcXN1/6+Bbdb+d9bJxMZ5nI0z9jkAfzy+glsOXSMX76109EsxpizY4UiyPx8+U62FR7jl9ePJy2+r9NxALh6/GDuvCiNZz7Yb8NmjQlAViiCyMd7q3j2w/3cceEw5p0/2Ok4n/Lw/FFMGDKA//znNg5W1TsdxxhzBqxQBInjTa089K/tDIuL4qF5o52O8xkRoS7+cMtkRODBf2yltZOlzI0x/scKRZD4zcpdHKxq4OefG0+fcGf7JTozNDaKn1w7ji2HjvF/7+9zOo4xppusUASBrYXH+NMH+7llWioXnhfndJwuXTMhmbljB/Gbt3ezx+ZXGBMQrFAEuOZWN995eTtJ/SN5eN4op+Oclojw0+vG0S8ylG/+cxstrTZk1hh/Z4UiwL2w7iC7ymr50TVjiY4MczpOt8T3i+Ani8ax/XA1T621Jihj/J0VigB2tL6Jx1bt4ZKR8cwZk+R0nDNy9fjBzD9/EL9bvYdDVQ1OxzHGdMEKRQB7fNVuahub+f6CMQG5DekPFowlNER4dEkOwbYvijHBxApFgNpdVsvz6w9x67RhZA6KdjrOWRk0IJL/mJPBO7sqWJFrE/GM8VdWKAKQqvKTN/PoG+7iP+ZkOB3nnNx5URqjBkXzozfyqD/R4nQcY0wHrFAEoI/2VvH+nkoemJ1BbN9wp+Ock1BXCD+7bhwl1Y38dvUep+MYYzpghSLAqCq/ensXyQMiuW16qtNxesSUYbHckDWEP3+4n/2VtryHMf7G0UIhInNFZJeIFIjIQx28nyoi74jIFhHZLiLzncjpT97ZVc6WQ8f4+qx0IkL9cwb22fjWVZmEu0L472X5TkcxxpzCsUIhIi7gCWAeMAa4WUTGnHLY94CXVHUScBPwv75N6V/cbuVXK3aTGhvF9VOGOB2nRyVGR/LVy0fydl4ZH+2tdDqOMaYdJ+8opgIFqrpPVZuAxcCiU45RoL/n8QCg2If5/M5buaXkldTw4Ox0wlzB12p49yXDSYnpw0/ezLdFA43xI07+tEkBCts9P+x5rb0fAreJyGFgGfB130TzP61u5Tcrd3NeQl8WTTz1nyk4RIa5+M68UeSX1PDypsLTf8AY4xP+/mvpzcBfVHUIMB/4m4h8JrOI3Csi2SKSXVFR4fOQvvBWTikF5XU8ODsDV0jgTa7rroXjBzM5NYZfvb2bhiYbLmuMP3CyUBQBQ9s9H+J5rb27gZcAVPVjIBKIP/VEqvq0qmapalZCQoKX4jpHVfnjewUMj+/LfD/bkKiniQjfnT+aitoT/PnDA07HMcbgbKHYCKSLyHARCaets3rJKcccAmYBiMho2gpFcN4ydOGDgkpyimr4txkjgvpu4qSstFhmj07iyXf3crS+yek4xvR6jhUKVW0B7gdWAPm0jW7KFZEfi8g1nsO+CXxZRLYBLwJ3ai9cFOh/39lLUv8IrpscnH0THfn23Ezqm1p44p0Cp6MY0+uFOnlxVV1GWyd1+9d+0O5xHnCxr3P5ky2HjvLxvioemT86qOZNnE5GUjTXTxnCXz8+yJ0XpzFkYJTTkYzptfy9M7vXe/K9vQzoE8bN04JjFvaZeHB2Bgg8ttKW9jDGSVYo/NjeijpW5JZxx4XD6Bfh6M2fI5Jj+nDnRWm8uuUwBeW2baoxTrFC4cf+/OF+wl0h3H5RmtNRHHPfZecRFR5qdxXGOMgKhZ+qbmjmlU1FXDMxmfh+EU7HcUxs33C+dMlwlu4oIaeo2uk4xvRKVij81D+yD3G8uZW7Lk5zOorj7rl0OAP6hPHYyt1ORzGmV7JC4YdaWt0899FBpg6PZWzyAKfjOK5/ZBj/dtkIVu8sZ/Oho07HMabXsULhh1bll1F07DhfsruJT9x5URrx/cL51YpdTkcxptexQuGHnv3gACkxfZgzZpDTUfxGVHgoX5k5ko/2VrFuX5XTcYzpVaxQ+Jnc4mo2HDjCnRel9YrlOs7ErdNSSYyO4Dcrd9MLJ+gb4xgrFH7mhfWHiAgN4Yasoac/uJeJDHPxtctHsmH/ET7ea3cVxviKFQo/Uneihde3FLFwQjIDosKcjuOXbrxgKIMHRNpdhTE+ZIXCj7y2pYj6plZu7YXLdXRXZJiLr14+kuyDR3l/j22ZaowvWKHwE6rKC+sPMXpwfyYOjXE6jl+7IWsIKTF9+LXdVRjjE1Yo/MTWwmPkl9Rw67RURKwTuysRoS7uv2Ik2wqP8e6uXrc9iTE+Z4XCT/x9/SH6hru4dlLv2XPiXFw/ZQhDY/vw2Cq7qzAGYPGGQ/zt4wNe+X6wQuEHqo8388b2YhZNSumVq8SejTBXCF+/Ip3th6tZnV/udBxjHNXqVn63eg9v55V5pUXCCoUfeH1rEY3Nbm6Zap3YZ+Jzk1IYFhdlI6BMr/f+ngqKqxu58QLvDKu3QuEHXsouZGxyf8al2LpOZyLUFcK/X5FOXkkNK3LLnI5jjGNeyi5kYFQYc8YkeeX8VigclldcQ05RjU2wO0uLJiYzIr4vj6/ajdttdxWm96mqO8HKvDKunZTite2SrVA47J+bCgl3hbBoYrLTUQJSqCuEB2ans7O0luU5pU7HcZSqcrCqnje2FbNsRwml1Y1ORzI+8OqWIppb1WvNTgDWc+qgphY3r20pYs7YJGKiwp2OE7AWjE/miXcK+M3KXcwdN6jXrZFVdOw4//PWTt7dXcGxhuZPvZcS04e54wbxzSsziAq3b/dgo6q8lF3IhKExjBrU32vXcfT/HBGZC/wWcAHPqOrPOzjmBuCHgALbVPUWn4b0otX5ZRxtaOYLU4Y4HSWguUKEb8zJ4L7nN/P61iI+N7l3/Hs2Nrfy1Hv7+ON7BQAsHJ/MpNSBjB8ygFa3sungUTYeOMKzH+5nzc5yHrtxok3mDDJbC4+xu6yO/7rufK9ex7FCISIu4AlgDnAY2CgiS1Q1r90x6cDDwMWqelREEp1J6x0vZRcyqH8kl6YnOB0l4F01dhBjk/vz+Ko9LJyQTJgruFtVD1U1cOdfNrCvop6rxw/mu/NHkxLT51PHTBgaw5cuGc7He6v41j+38fk/fsS3rszkKzPPcyi16WkvZRfSJ8zFwgmDvXodJ7+bpgIFqrpPVZuAxcCiU475MvCEqh4FUNWgGTBfWt3Ie7sr+PyUlF7XVOINIsK3rszk0JEG/pl92Ok4XpVTVM3n/vgRR+qbeP7uaTxxy+TPFIn2LjwvjuUPXsrccYP4xVs7eXHDIR+mNd7S2NzKG9tKmHf+IKIjvbuIqJOFIgUobPf8sOe19jKADBH5UETWeZqqPkNE7hWRbBHJrqgIjCUdXt1ShFvh+ik22qmnzMxMYHJqDL9fs4fG5lan43jFhwWV3PT0OiJCQ3j5vgu5JD2+W5/rHxnGb2+cyMzMBL7/Wg4f7bUFFQPditxS6k60cL0Pmlr9/f48FEgHZgI3A/8nIjGnHqSqT6tqlqpmJST4fzOOqvLqlsNMGTaQ4fF9nY4TNESEb12VSUl1I8+vO+h0nB638cAR7vrLRlJi+vDKVy5iZGL0GX0+1BXC72+exIiEvnzl+c3sq6jzUlLjC//aXETygEimj4jz+rWcLBRFQPtfp4d4XmvvMLBEVZtVdT+wm7bCEdDySmrYXVZn6zp5wUXnxTMjI4E/vFNATWPz6T8QIArK67jnuWyGxPRh8b3TGTQg8qzOEx0Zxp/uuIDQEOGev2YH7Z1XsCuvaeT9PRVcNzmFEB80XTtZKDYC6SIyXETCgZuAJacc8xptdxOISDxtTVH7fJjRK17dXESYS1hwvnc7oHqrb1+VybGGZp56b6/TUXpEeW0jd/55A2Eu4bkvTWVg33MbSj00NorHb5rIvop6nl4b8N9OvdLrW4txK1w3yTcj/BwrFKraAtwPrADygZdUNVdEfiwi13gOWwFUiUge8A7wn6oa0HtgtrqV17cVMzMz8Zy/4U3HxqUMYNHEZP70wX7KagJ70lljcyv3PJdNVV0Tz955AUNjo3rkvJemJ3D1+YN54p0CCo809Mg5jW+oKq9sPsyEoTGMTOznk2s62kehqstUNUNVz1PVn3le+4GqLvE8VlX9hqqOUdXzVXWxk3l7wocFlVTUnuBz1uzkVd+ck0mrW/nt6j1ORzlrqsr3XsthR1E1v795EuOHxPTo+R+5ejQhIvx0ad7pDzZ+I6+khp2ltVw/2Xc/Q/y9MzvovLaliP6RoVw+KqimhPid1Lgobp02jH9sLKSgPDA7bRdvLOTlTYf59yvSme2Fxd6SY/pw/xUjWZFbxnu7A2O0oGnrxA5zCQvG+27ZHysUPtTQ1MJbuaVcPX4wkWHeWbzL/H/3XzGSqDAX/7Us3+koZ2z74WM8+nouMzIS+PdZ3hu/cc+lw0mLi+JHS3JpaXV77TqmZ7S0unl9azGX+7jp2gqFD72dW0ZDU6vPOqB6u/h+EXx91kjW7CwPqN+Yqxua+crzm4nvF87jN0706oTMiFAXD88fzb7KepbuKPHadUzP+GhvFZV1J7jOx03XVih86PWtRaTE9CFr2ECno/Qad1yUxrC4KH7yZh7NAfAbs6ry7Ve2UVbTyBO3TibWB781zhmdRHpiP558b59tAOXnXt9aTHSE75uurVD4yNH6Jt7fU8mCCYN9Mu7ZtIkIdfHI/NEUlNfx9/X+v3TF8+sOsiK3jO/MHcWkVN/8QhESItw7YwT5JTUBdefV2zQ2t7Iit5S54wb5vOnaCoWPLM8ppcWtLPRhB5RpM2dMEhePjOOxVbs51tDkdJxO5RXX8JOl+czMTODuS4b79NqLJqYweEAkTwbJ3JNgtDq/nLoTLY5M1LVC4SNvbCtmREJfxiZ7b8140zER4fsLxlDb2MLPl+90Ok6H6k+0cP+Lm4npE8avvzDB53ed4aEh3H3JcNbtO8LWwmM+vbbpnte3FpEYHeGTJTtOZYXCB8pqGlm3v4qF45MRsWYnJ4wa1J97LhnO4o2FrN/nf3M2f/B6Lvsr63n8ponE9YtwJMNNU1PpHxnKk+/aXYW/qW5o5t1dFSyckOzIatNWKHxg6fYSVGHhBGt2ctIDs9MZMrAPD7+6gxMt/rPG0cubDvPK5rb5Ehed173VYL2hX0Qot1+Yxoq8Ug5V2Wxtf7I8p4SmVjfXTnRmoq4VCh9Ysq2YMYP7+2y6velYVHgoP712HPsq6vmjn/zWXFBex/dfy2H6iFivzpforlunpyK07eVu/MfrW4sZEd+XcSnONF1bofCywiMNbC08xjUT7W7CH8zMTOSaCcn87zt72VNW62iWxuZW7v/7ZqLCXfz2pkl+sYHV4AF9mJGRwMubDtPqtqGy/uBk0/U1E51rurZC4WVvbC8GYMF4WynWX3x/wRj6Rrh4YPFWx5qgVJXvvrqDXWW1/PqGCST1P7tlw73hxqyhlFQ3snaPDZX1B8t2tDVd+3LJjlNZofCypdtLmJQaw5CBPbPqpzl3CdER/PL6CeSV1PA/b+1yJMNzHx3gX5uLeHBWBjMz/Wvdr1mjk4jrG85LG635yR+8sa2Y0Q43XVuh8KIDlfXkFtdwte074XfmjEnii9OH8cwH+30+yWzD/iP8dGk+s0cn8fUrRvr02t0RHhrCdZNSWJlXRmXdCafj9GqHjzaw+dAxx1skrFB40cm1c+ZZofBLj1w9moykfnzzpW1U1PrmB2LxseN89YVNpMZG8ZsbfT9fortuvGAoLW7l1c2nbjppfGnp9rafIU5P1LVC4UXLdrQ1O6XE9HE6iulAZJiL3988mdrGZr7y/Cavbwt6rKGJO57dQGOzm6e+OIX+kWFevd65SE+KZnJqDP/ILrT1nxz05vYSJgwZQGqcs03XVii8xJqdAkPmoGgeu3Ei2QeP8s1/bsPtpZE+x5taufu5bA5WNfD07VNIT4r2ynV60o0XDKWgvI5th6udjtIrHaisZ0dRtaOd2CdZofASa3YKHPPPH8x3549i6fYSfrmi5zu3m1vd3P/3zWw+dJTf3jTR0Ul1Z2LuuMGEu0J4c1ux01F6pTc9Iyav9oMRk1YovGTZjhImDrVmp0Dx5UtHcNv0VJ58by//t3Zfj523sbmVr72wmdU7y/nJonEB9YvDgD5hzMiIZ+mOEq/daZnOvbm9hKxhA0n2g58hVii84GCVNTsFGhHhhwvHMm/cIH62LJ+fvpl3zj8cj9Y3ccv/rWNlfhk/XDiG26YP66G0vrNgfDIl1Y1sKTzqdJRepaC8lp2ltY6PdjrJ0UIhInNFZJeIFIjIQ10c93kRURHJ8mW+s/X/m50GOZzEnIlQVwh/uGUyd1zYNmz264u3nHUH976KOj7/5EfkFNfwv7dM5s6LfbtseE+ZNTqR8NAQ3thmu9/50tLtpYj4T9O1Y4VCRFzAE8A8YAxws4iM6eC4aOABYL1vE5695TtKmTDUJtkFIleI8MNrxvLwvLY+i0V/+PCMVpttanHzhzV7mPvb96mqa+KFe6b5zTf72YiODOPyzASW7SixJT18aNmOtmYnf5mx7+QdxVSgQFX3qWoTsBhY1MFxPwF+ATT6MtzZOny0gR1F1cwbZ3cTgUpE+LfLzuOZ27OoO9HCjU+v44HFWzhQWd/pZ443tbJsRwlX/+59fvX2buaMTmLlf8zggrRYHyb3jgXjkymvPUH2gSNOR+kVCsrr2FVW61dN16EOXjsFaL9GwGFgWvsDRGQyMFRVl4rIf3Z2IhG5F7gXIDU11QtRu++tnFIA5o61QhHoZo9J4uKR8fzx3QKeXLvvkxU8Z2YmMjKxH00trTS1utlRVMPq/DIamlpJienDn+7IYtboJKfj95grRiUSGRbCm9tLmObApjm9zbIdJX7V7ATOFoouiUgI8BvgztMdq6pPA08DZGVlOXp//FZOKaMGRZMW39fJGKaH9Al38Y0rM7lpaiorckt5d1cFL6w/yIkW9yfHDIwKY9HEFBaOH8zU4bGEuoJrjEjfiFBmjUpieU4Jjy4cE3R/P3/jb81O4GyhKAKGtns+xPPaSdHAOOBdz9K6g4AlInKNqmb7LOUZKK9pZNOhozw4K8PpKKaHJcf04a6Lh3PXxcNpbG7laEMTEaEuIkJD6BPm8tulOHrKgvGDWbqjhA0HjgTMPJBAtLeijp2ltTy68DPdtY5y8leDjUC6iAwXkXDgJmDJyTdVtVpV41U1TVXTgHWA3xYJgBV5ZajCXOufCGqRYS4GD+hDbN9w+kaEBn2RALgsM4GI0BBW5pU5HSWoLfOs7TRvnP80O4GDhUJVW4D7gRVAPvCSquaKyI9F5Bqncp2LFTmljIjvS0aS7WRngktUeCiXjIxnZV6Zrf3kRUs9zU6DBvhPsxM4PI9CVZepaoaqnqeqP/O89gNVXdLBsTP9+W7iaH0TH++r4qpxgxzbhcoYb5o9JonDR4+zy+GdAYPVPk+z03w/6sQ+qVt9FCKSCFwMJAPHgRwgW1XdXX6wF1mVX0arW21YrAlas0YnIgIrc8sYNciZvZuD2XLPiEl/nKjb5R2FiFwuIiuApbRNjBtM2+S47wE7RORHImL/x9A22iklpg/npwxwOooxXpEYHcnEoTGszLd+Cm9YntO2LcHgAc6v7XSq091RzAe+rKqHTn1DREKBBcAc4BUvZAsY9SdaeL+gklunpVqzkwlqs0cn8T8rdlFa3eh37eiB7FBVAzlFNTwyf7TTUTrU5R2Fqv5nR0XC816Lqr6mqr26SACs3V1BU4ubK8f43y2jMT3pyjFtEwlX2V1Fj1qe0zbayV9HTHarM1tE/iYiA9o9TxOR1d6LFVhW5JYyMCqMC9IGOh3FGK8amdiPtLgoGybbw5bnlHJ+ygCGxvrn+nDdHfX0AbBeROaLyJeBt4HHvZYqgDS3ulm9s5xZo5NsxqoJeiLCnDFJfLy3iroTLU7HCQpFx46ztfCYX3Zin9Stn2yq+hRwD/A68GNghqq+4c1ggWL9viPUNrZ8cktuTLCbPTqJplY3a3dXOB0lKJxcH87fJtm1192mpy8CzwK3A38BlonIBC/mChgrckvpE+ZiRkaC01GM8YkpwwbSPzKUd3eVOx0lKLyVU8KoQdEM9+P14brbVvJ54BJVfVFVHwbuo61g9Gput7Iyr4wZGfFEhrmcjmOMT4S6Qrg0I4F3dlXYLO1zVF7TSPbBo345ya697jY9Xauq5e2eb+CUJcF7ox1F1ZTWNNpoJ9PrXJ6ZSEXtCXKLa5yOEtBW5JYGxPpwp5tw9z0R6XDnFVVtEpErRGSBd6L5vxW5pbhChFmjE52OYoxPXeZparXmp3OzPKeUEQl9SU/07/XhTjfhbgfwhog0ApuBCiASSAcmAquA//JmQH/2dl4ZU9NiiYkKdzqKMT6VEB3B+CEDeHdXBfdfke50nIB0pL6J9fuPcN9lI/x+ou7pmp6uV9WLaVvhNRdwATXA88BUVf0PVe2VQx/2V9ZTUF7HlWNttJPpnWZmJLD50FGONTQ5HSUgrco7uT6cf/dPwOkLxRQRSQZupW2viKeAv9K2l4T/LUjiQ6s8E45mB9GWl8aciZmjEnErrN1T6XSUgLQ8p4QhA/swNtn/l8s7XaF4ElgNjAKy231t8vzZa63MK2PUoGi/nUlpjLdNGBLDwKgw3t1p/RRnqqaxmQ8KKpk7NjC2JTjdWk+/U9XRwLOqOqLd13BVHeGjjH7nSH0T2QeP2CQ706u5QoQZGQm8t7sCt9uGyZ6Jd3aW09yqfj0bu73uDo/9ireDBJI1O8txK8yxYbGml7s8M5Gq+iZ2FFU7HSWgLN9RSmJ0BJOGBsb6cLY40VlYmVfKoP6RjEvx/7ZFY7xpRkYCIvDurl45puWsNDS18O7ucq4aOyhg9lu3QnGGGptbWbu7ktljEgOibdEYb4rtG874lAGs3WOForve21VBY7Pb7yfZtWeF4gx9WFDJ8eZWa3YyxmNGRgJbDh2luqHZ6SgB4S3PtgTThnc4l9kvOVooRGSuiOwSkQIReaiD978hInkisl1EVovIMCdytrcqv4x+EaFMHxE4/5GN8aYZGQm4FT7ca8NkT+dESytr8suZMyawtiVwLKmIuIAnaNuLewxws4iMOeWwLUCWqo4HXgZ+6duUn+Z2K6vyy7ksI4GIUFsE0BiAiUNjiI4ItWXHu+HDgkpqT7QExCS79pwsaVOBAlXdp6pNwGJgUfsDVPUdVW3wPF0HDPFxxk/ZXlRNRe0J5tiwWGM+EeYK4aKRcazdbavJns5bOaVER4Ry0cg4p6OcEScLRQpQ2O75Yc9rnbkbWN7RGyJyr4hki0h2RYX3fqtZlVeGK0SYmWl7TxjT3mUZiRRXN7K3os7pKH6rpdXNyrwyZo1ODLgWiYBoJBOR24As4H86el9Vn1bVLFXNSkjw3g/xVfllZA0baIsAGnOKGRnxALy32/opOrN+/xGONjQH1Gink5wsFEXA0HbPh3he+xQRmQ08Alyjqid8lO0zCo80sLO01pqdjOnAkIFRjEjoa/0UXXgrp203zMsyAm9bAicLxUYgXUSGi0g4cBNtCw9+QkQm0bYQ4TXtN05ywur8tkUAZ9kigMZ0aEZ6Auv3V9HY3Op0FL/jdisrckuZmZlAn/DAanYCBwuFqrYA99O2hHk+8JKq5orIj0XkGs9h/wP0A/4pIltFZEknp/O61TvLOS+hr1/va2uMky7LSKCx2c2G/UecjuJ3Nh86SnntiYBsdoLTb1zkVaq6DFh2yms/aPd4ts9DdaCmsZl1+6r40sXDnY5ijN+aPiKOcFcIa3dXMCPDBny0tzynlHBXCFeMCrxmJwiQzmynrd1dQXOrMtv6J4zpVJ9wFxcMH8j7tj/Fp6gqb+WUcml6PNGRYU7HOStWKLphdX45A6PCmJwaGCs9GuOUGekJ7CqrpbS60ekofmP74WqKjh0P2GYnsEJxWi2tbtbsLOfyzERcAbLSozFOuTS9rcnpfVsk8BPLc0oJDZGAHjFpheI0Nh08SvXxZmt2MqYbRg+OJiE6wrZH9WhrdirhwvPiAnr+lRWK01i9s5wwl3BperzTUYzxeyJt3ysf7LFd7wDyS2o5UNUQcGs7ncoKxWmsyi9j+oi4gO2EMsbXZqQncLShmZxi2/XurZwSQgSuHBvYLRJWKLqwv7KefRX1ATukzRgnXOK5+7ZZ2m39E1OHxxLfL8LpKOfECkUXTs7Gnm2zsY3ptvh+EYxN7t/r+ykKymvZU14X8M1OYIWiS6vzy8lI6sfQ2CinoxgTUGZkJLD54FFqG3vvrnfLdpQCcNXYwB0We5IVik5UH29m44EjtraTMWfh0vR4WtzKun29dzmPZTtKyBo2kEEDIp2Ocs6sUHRi7e4KWtzK7NHWP2HMmcoaFktUuIv3dju6lqdj9lbUsbO0lvnnB36zE1ih6NTq/DJi+4YzcajNxjbmTIWHhnDReXG810t3vVu+owSAeecHfrMTWKHoUEurm3d3VzAzM8FmYxtzli7LSKDwyHEOVDWc/uAgs3RHKVOGDWTwgD5OR+kRVig6sPnQMY41NDNrlPVPGHO2Tm7Q896u3tX8tL+ynvySGuYF8NpOp7JC0YHV+WWEueST7R2NMWcuNS6KtLioXjdMdpmn2SlY+ifACkWHVu8sZ9pwm41tzLm6LCOBj/f2rl3vlm4vYVJqDMkxwdHsBFYoPuNgVT0F5XU2G9uYHnBZZgLHm1vJPnDU6Sg+caCynrySGq4OorsJsELxGavz29pTZ9mwWGPO2cld73rLMNmln4x2skIR1NbsLGdkYj+Gxdne2Macq6jwUKYOj+W9XrLu05vbS5icGkNKEDU7gRWKT6ltbGb9/ipmWbOTMT3msowEdpfVUXzsuNNRvKqgvI78khoWjE92OkqPc7RQiMhcEdklIgUi8lAH70eIyD88768XkTRv5nl/TyXNrWrLdhjTg2Zk9I5d797cXowIXD0+uJqdwMFCISIu4AlgHjAGuFlExpxy2N3AUVUdCTwG/MKbmVbllzGgTxiTU2O8eRljepWMpH4MHhDJOzuDt1CoKm9sK2ba8FiS+gf+2k6ncvKOYipQoKr7VLUJWAwsOuWYRcBznscvA7NExCtTpVvdyru7Krg8M4FQl7XIGdNTRISZmYl8UFBJU4vb6ThekV9Sy96K+qBsdgJnC0UKUNju+WHPax0eo6otQDUQd+qJROReEckWkeyKirP7raWk+jh9wlxcYc1OxvS4K0YlUneihY0HgnM12Te3F+MKkaCajd1eUPzqrKpPq2qWqmYlJCSc1TmGDIzig+9cHnTjn43xBxePjCM8NIQ1O4NvmKyq8sb2Yi46L464AN/JrjNOFooiYGi750M8r3V4jIiEAgOAKm8FEhFbBNAYL4gKD2X6iLigLBTbD1dTeOQ4CycEZ7MTOFsoNgLpIjJcRMKBm4AlpxyzBLjD8/h6YI32xjWLjQkCs0Ylsr+ynv2V9U5H6VFLthUT5hKuGhOczU7gYKHw9DncD6wA8oGXVDVXRH4sItd4DvsTECciBcA3gM8MoTXGBIaTy+IE011Fq1tZsq2YyzMTGRAVvGvDhTp5cVVdBiw75bUftHvcCHzB17mMMT1vaGwUIxP78c7Ocu6+ZLjTcXrER3srqag9wbWTTh2HE1yCojPbGBMYrhiVyPr9VdSdaHE6So94bUsx0RGhQb+IqBUKY4zPXJ6ZSHOr8kEQzNI+3tTKWzklzDt/EJFhLqfjeJUVCmOMz2SlDSQ6MvSTVZoD2ar8MuqbWoO+2QmsUBhjfCjMFcLMzETW7Cyn1R3YAxhf21LEoP6RTB/+mTnAQccKhTHGp64am0RVfRPZATxL+0h9E+/trmDRxGRCesHcKysUxhifmpmZSHhoCCtyy5yOctaW7iihxa0smhj8zU5ghcIY42P9IkK5ZGQ8K3JLCdT5s69sOkxmUjSjB0c7HcUnrFAYY3zuqrFJFB07Tm5xjdNRztieslq2Fh7jC1lD8NJi1n7HCoUxxudmj04iRODt3FKno5yxf246TGiIcF0vGO10khUKY4zPxfWLICstlrfzAqufornVzb82H+aKUYlBu1JsR6xQGGMccdXYQewsreVgVeAsEvjurgoq65q4IWvo6Q8OIlYojDGOuHJM2yZhKwKo+emf2YXE94tgZubZ7XsTqKxQGGMcMTQ2irHJ/VmeExiForLuBGt2lvO5ySm9brvk3vW3Ncb4lQXjk9ly6BiHqhqcjnJar20posWtfGHKEKej+JwVCmOMYxZOaNt6+I3txQ4n6Zqq8uKGQ0xKjSE9qXfMnWjPCoUxxjFDBkaRNWwgS7b6d6H4eF8VeyvquW3aMKejOMIKhTHGUYsmJrOrrJadpf47+e75dQeJiQrj6vGDnY7iCCsUxhhHzT9/MK4Q4XU/vasoq2nk7dwybsgaGvT7TnTGCoUxxlFx/SK4ND2eJVuLcfvh0uOLNxTS4lZumZrqdBTHWKEwxjjumgnJFB07zuZDR52O8iktrW5e3HCIGRkJpMX3dTqOYxwpFCISKyIrRWSP58+BHRwzUUQ+FpFcEdkuIjc6kdUY431Xjh1ERGiI3zU/rcovp7Smkdum9d67CXDujuIhYLWqpgOrPc9P1QDcrqpjgbnA4yIS47uIxhhf6RcRypwxSbyxvZjG5lan43zirx8fIHlAJFeMSnQ6iqOcKhSLgOc8j58Drj31AFXdrap7PI+LgXKgd82bN6YXuXlqKscamlmeU+J0FAByiqr5aG8Vt1+U1utmYp/Kqb99kqqe/L+hFEjq6mARmQqEA3s7ef9eEckWkeyKioqeTWqM8YkLR8SRFhfFC+sOOR0FgKfW7qNfRCi39PJmJ/BioRCRVSKS08HXovbHadsWV50OdRCRwcDfgLtU1d3RMar6tKpmqWpWQoLddBgTiEJChFumpZJ98Ci7SmsdzVJ4pIGl24u5dVoq/SPDHM3iD7xWKFR1tqqO6+DrdaDMUwBOFoLyjs4hIv2BpcAjqrrOW1mNMf7h+ilDCXeF8Pf1Bx3N8cz7+3CFCHddPNzRHP7CqaanJcAdnsd3AK+feoCIhAOvAn9V1Zd9mM0Y45DYvuHMO38Q/9pSRENTiyMZjtQ38Y/sQhZNTGHQgEhHMvgbpwrFz4E5IrIHmO15johkicgznmNuAGYAd4rIVs/XREfSGmN85tZpw6htbOHNbc50av/14wM0Nru5d8YIR67vj0KduKiqVgGzOng9G7jH8/h54HkfRzPGOOyCtIGkJ/bj+fUH+ULWEETEZ9euPt7Mnz88wKxRiWT0wlViO9O7x3wZY/yOiHD7RWlsP9w2PNWX/m/tPqqPN/MfczJ8el1/Z4XCGON3vjBlCInREfxu9R6fXbO8tpE/fbCfhROSGZcywGfXDQRWKIwxficyzMV9l53H+v1H2LD/iE+u+Yc1BTS3uvmm3U18hhUKY4xfunlqKvH9wvn9Gu/fVRyqauDv6w9x4wVDe/Xif52xQmGM8Ut9wl18+dIRvL+nki1eXlX21yt3EeoS/n1WulevE6isUBhj/NZt04cxMCqM368p8No1Piqo5PWtxdx9yXCS+tu8iY5YoTDG+K2+EaF8ecYI1uws573dPb+O2/GmVh761w7S4qL4+hV2N9EZKxTGGL929yXDGZHQl++9toPjTT27BPljq3Zz6EgD//258b12m9PusEJhjPFrEaEu/uu68yk8cpzHV+/usfNuKzzGM+/v4+apqVx4XlyPnTcYWaEwxvi96SPiuCFrCM+8v5+84ppzPl9DUwvfeWU7CdERPDx/VA8kDG5WKIwxAeG780cT0yeMh1/dQVNLhzsOdIvbrTy4eCu7ymr5+efH2zLi3WCFwhgTEGKiwvnxonFsKzzGd17Zjtvd6TY2XfrFip28nVfG964ew+WZvXuL0+5yZFFAY4w5G1ePH8z+ygx+9fZuEqMjeHj+6DP6/OINh3jqvX3cNj2VL12c5p2QQcgKhTEmoHzt8pGU157gqbX7SIiO4J5LT78cuNutPLl2L79asYtL0+P54cKxPl2VNtBZoTDGBBQR4dGFY6msO8FPl+aTV1LDI/NHE9cvosPjqxua+cZLW1m9s5yrxw/ml58fT6jLWt3PhBUKY0zAcYUIj904kRHxBTy1di9rdpbzn1dlcunIBJJjInGFCHklNbyVU8rLmw5TWXeCH10zltsvHGZ3EmdBVM+uQ8hfZWVlaXZ2ttMxjDE+sqeslkdezWHDgbZVZkNDhAF9wqiqbyJEYOrwWL4zdxSTUgc6nNS/icgmVc3q6D27ozDGBLT0pGgW3zudLYVH2Vtez/6qespqGrkgLZY5Y5KI76RJynSfFQpjTMALCRGmDItlyrBYp6MEJevRMcYY0yVHCoWIxIrIShHZ4/mz08ZDEekvIodF5A++zGiMMaaNU3cUDwGrVTUdWO153pmfAGt9ksoYY8xnOFUoFgHPeR4/B1zb0UEiMgVIAt72TSxjjDGncqpQJKlqiedxKW3F4FNEJAT4NfCt051MRO4VkWwRya6o6PnNTYwxpjfz2qgnEVkFDOrgrUfaP1FVFZGOJnN8FVimqodPN0FGVZ8Gnoa2eRRnl9gYY0xHvFYoVHV2Z++JSJmIDFbVEhEZDJR3cNiFwKUi8lWgHxAuInWq2lV/hjHGmB7m1DyKJcAdwM89f75+6gGqeuvJxyJyJ5BlRcIYY3zPkSU8RCQOeAlIBQ4CN6jqERHJAu5T1XtOOf5O2grF/d04d4XnnGcrHqg8h8/7UiBlhcDKG0hZIbDyBlJWCKy855J1mKomdPRG0K31dK5EJLuz9U78TSBlhcDKG0hZIbDyBlJWCKy83spqM7ONMcZ0yQqFMcaYLlmh+KynnQ5wBgIpKwRW3kDKCoGVN5CyQmDl9UpW66MwxhjTJbujMMYY0yUrFMYYY7pkhcJDROaKyC4RKRARv57YJyLPiki5iOQ4neV0RGSoiLwjInkikisiDzidqSsiEikiG0Rkmyfvj5zOdDoi4hKRLSLyptNZTkdEDojIDhHZKiJ+vWexiMSIyMsislNE8kXkQqczdUZEMj3/pie/akTkwR47v/VRtH2jAbuBOcBhYCNws6rmORqsEyIyA6gD/qqq45zO0xXPEi2DVXWziEQDm4Br/fjfVoC+qlonImHAB8ADqrrO4WidEpFvAFlAf1Vd4HSerojIAdomz/r9BDYReQ54X1WfEZFwIEpVjzkc67Q8P8+KgGmqei6Tjz9hdxRtpgIFqrpPVZuAxbQthe6XVHUtcMTpHN2hqiWqutnzuBbIB1KcTdU5bVPneRrm+fLb36ZEZAhwNfCM01mCiYgMAGYAfwJQ1aZAKBIes4C9PVUkwArFSSlAYbvnh/HjH2aBSkTSgEnAeoejdMnTlLOVtsUqV6qqP+d9HPg24HY4R3cp8LaIbBKRe50O04XhQAXwZ0+z3jMi0tfpUN10E/BiT57QCoXxCRHpB7wCPKiqNU7n6YqqtqrqRGAIMFVE/LJ5T0QWAOWqusnpLGfgElWdDMwDvuZpRvVHocBk4I+qOgmop+udOP2Cp4nsGuCfPXleKxRtioCh7Z4P8bxmeoCnrf8V4AVV/ZfTebrL09TwDjDX4SiduRi4xtPuvxi4QkSedzZS11S1yPNnOfAqbc2+/ugwcLjd3eTLtBUOfzcP2KyqZT15UisUbTYC6SIy3FORb6JtKXRzjjydw38C8lX1N07nOR0RSRCRGM/jPrQNcNjpaKhOqOrDqjpEVdNo+392jare5nCsTolIX8+ABjzNOFcCfjlyT1VLgUIRyfS8NAvwywEYp7iZHm52Auf2o/ArqtoiIvcDKwAX8Kyq5jocq1Mi8iIwE4gXkcPAo6r6J2dTdepi4IvADk+7P8B3VXWZc5G6NBh4zjNyJAR4SVX9fthpgEgCXvXsWBkK/F1V33I2Upe+Drzg+eVxH3CXw3m65Cm+c4B/6/Fz2/BYY4wxXbGmJ2OMMV2yQmGMMaZLViiMMcZ0yQqFMcaYLlmhMMYY0yUrFMYYY7pkhcIYY0yXrFAY42UicoGIbPfsddHXs8+FX64fZUxHbMKdMT4gIj8FIoE+tK0h9N8ORzKm26xQGOMDnmUgNgKNwEWq2upwJGO6zZqejPGNOKAfEE3bnYUxAcPuKIzxARFZQttS4MNp2xr2focjGdNttnqsMV4mIrcDzar6d8+qtB+JyBWqusbpbMZ0h91RGGOM6ZL1URhjjOmSFQpjjDFdskJhjDGmS1YojDHGdMkKhTHGmC5ZoTDGGNMlKxTGGGO69P8A2+U8PKXV20cAAAAASUVORK5CYII=\n",
+ "text/plain": [
+ "<Figure size 432x288 with 1 Axes>"
+ ]
+ },
+ "metadata": {
+ "needs_background": "light"
+ },
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "# Here we define the functional form\n",
+ "def f(X):\n",
+ " x = X[0, 0]\n",
+ " y = np.sin(x) + 1.5*np.exp(-(x - 4.3)**2)\n",
+ " return y\n",
+ "\n",
+ "# Here we define the domain and examine the function\n",
+ "bounds = np.array([[0., 7.]])\n",
+ "fig, ax = plt.subplots()\n",
+ "x = np.linspace(bounds[0, 0], bounds[0, 1], 100)\n",
+ "y_true = np.array([f(np.atleast_2d(xi)) for xi in x])\n",
+ "ax.plot(x, y_true)\n",
+ "ax.set_xlabel('x')\n",
+ "ax.set_ylabel('f(x)')\n",
+ "ax.set_title('1D function')\n",
+ "plt.show()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "8EA4GyoCC8Ir",
+ "lines_to_next_cell": 2
+ },
+ "source": [
+ "## Running BO\n",
+ "The first step in applying BO to a problem is to consider problem dimensionality and the optimization domain, which involves setting the upper and lower bound for each dimension. In BOSS software, this can be set by keywords in an input file, as listed below.\n",
+ "\n",
+ "From our problem above, we pass the sampling function and the bounds as an argument. A single entry in the `kernel` array indicates a 1D problem, with an expected `yrange` of approximately [-1; 1]. We begin BO with 2 initial points and acquire 5 more.\n",
+ "\n",
+ "Once the BO run is defined, we can execute it using the run function, which produces the the `boss.out` (main output file), `boss.rst` (restart file). "
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 62,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 208,
+ "status": "ok",
+ "timestamp": 1630579846701,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "FbkFz5egC8It"
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting file 'boss.out'\n",
+ "warning: overwriting file 'boss.rst'\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Here we define the BO run\n",
+ "bo = BOMain(\n",
+ " f, \n",
+ " bounds, \n",
+ " yrange=[-1, 1],\n",
+ " kernel='rbf',\n",
+ " initpts=2,\n",
+ " iterpts=5,\n",
+ ")\n",
+ "\n",
+ "# Here we perform BO\n",
+ "result = bo.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "eJ9_AjA3C8It"
+ },
+ "source": [
+ "To check the quality of the BO, we set up an run the BO postprocessing (PP) routines. The result is the `postprocessing` folder with graphs and raw data. "
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 63,
+ "metadata": {
+ "executionInfo": {
+ "elapsed": 1627,
+ "status": "ok",
+ "timestamp": 1630579850682,
+ "user": {
+ "displayName": "Joakim Löfgren",
+ "photoUrl": "https://lh3.googleusercontent.com/a/default-user=s64",
+ "userId": "04244087693112310568"
+ },
+ "user_tz": -180
+ },
+ "id": "Wqmj8WqrC8Iu"
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting directory 'postprocessing'\n"
+ ]
+ }
+ ],
+ "source": [
+ "pp = PPMain(result, pp_models=True, pp_acq_funcs=True, pp_truef_npts=100) # add more pp keywords here if needed\n",
+ "pp.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "zKy8LL1oC8Iu"
+ },
+ "source": [
+ "## Results and Analysis\n",
+ "Let's visualize the quality of the GP regression fit by checking the results in the `postprocessing` folder. This typically contains:\n",
+ "* a record of data acquisitions (locations and values), with the global minimum predictions\n",
+ "* a set of figures showing the GPR fit at every iteration\n",
+ "* a set of figures showing the GPR acquisition function at every iteration\n",
+ "* global minimum convergence monitoring\n",
+ "* a record of hyperparameter evolution throughout the run.\n",
+ "\n",
+ "You can tailor the amount of output data with postprocessing `keywords`."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### 1. Acquisition data check\n",
+ "For the acquired points, let's review the **acquisition locations** (in x, lower graph and y, upper graph), alongside the (x,y) tracking of the **global minimum prediction** (full line).\n",
+ "The key quality check questions to ask are:\n",
+ "\n",
+ "* is the sampling distributed across the entire x domain?\n",
+ "* are there any errors in y function evaluation (are the y values distributed across the expected domain)?\n",
+ "* is the predicted global minimum location (full line) converging?\n",
+ "\n",
+ "<img src=\"postprocessing/acquisition_locations.png\" width=\"400px\">"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### 2. Model refinement check\n",
+ "Let us review the **progression of the model fitting** with 1, 3 and 5 BO acquisitions (see below) against the true function.<br> The key quality check questions to ask are:\n",
+ "\n",
+ "* Is the model converging? \n",
+ "* Has the global minimum been found?\n",
+ "* Is the model uncertainty reducing?\n",
+ "* Has the entire model been converged to the true function?\n",
+ "\n",
+ "<table><tr><td><img src='postprocessing/graphs_models/it0000_npts0002.png' width=\"400px\"></td><td><img src='postprocessing/graphs_models/it0003_npts0005.png' width=\"400px\"></td><td><img src='postprocessing/graphs_models/it0005_npts0007.png' width=\"400px\"></td></tr></table>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### 3. Model convergence check\n",
+ "Let's consider at the **rate of convergence** of the global minimum in both x and y (left image), as well as the **GPR hyperparameter convergence** with BO iterations (right image). The key quality check questions to ask are:\n",
+ "\n",
+ "* is the **rate of convergence** reducing on the logarithmic scale?\n",
+ "* has the global minimum been evaluated with a sufficient level of accuracy, given data range?\n",
+ "* are the GPR hyperparameters - GP lengthscale and variance - converging?\n",
+ "\n",
+ "<table><tr><td><img src='postprocessing/convergence_measures.png'></td><td><img src='postprocessing/hyperparameters.png'></td></tr></table>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Exercise:\n",
+ "Re-define the BOSS run with 20 BO iterations and redo all quality check points. <br> Has the convergence improved?\n",
+ "Below you can view the model after 20 BO iterations - has the model converged to the true function? If you start the BO run with more initial points, does the model take fewer points to converge?"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "<img src=\"postprocessing/graphs_models/it0020_npts0022.png\" width=\"400px\">"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "colab": {
+ "collapsed_sections": [],
+ "name": "tutorial1.ipynb",
+ "provenance": []
+ },
+ "jupytext": {
+ "cell_metadata_filter": "-all",
+ "main_language": "python",
+ "notebook_metadata_filter": "-all"
+ },
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.6.3"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/Tutorial2_Structure-Search/.ipynb_checkpoints/BO-Structure-Search-checkpoint.ipynb b/Tutorial2_Structure-Search/.ipynb_checkpoints/BO-Structure-Search-checkpoint.ipynb
new file mode 100755
index 0000000000000000000000000000000000000000..ee0268d6ba51206d00f7640113051169309b6f91
--- /dev/null
+++ b/Tutorial2_Structure-Search/.ipynb_checkpoints/BO-Structure-Search-checkpoint.ipynb
@@ -0,0 +1,251 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Tutorial 2: Configurational Structure Search with BO"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "This tutorial demonstrates how to use BO to optimize the atomistic configuration of materials. This can be done by mapping the configurational degrees of freedom to the total energy landscape, and inferring the global minimum of energy. The energy landscapes, refined by iterative data acquisitons, can be employed to visualise low energy configurations (they are chemically interpretable). \n",
+ "\n",
+ "<img src=\"../figures/ss_scenarios.png\" width=\"800px\"/>\n",
+ "\n",
+ "At the outset, the user should consider which **key degrees of freedom** control the energetics of the research problem. In the case of adsorption of organic molecules to the surface, these could be the molecule position or orientation above the surface. The molecular conformation can often be kept fixed. In the case of flexible molecules, we need to conduct a conformer search and minimize the energy with respect to dihedral angles. The degrees of freedom and their bounds define the N-dimensional BO structure search.\n",
+ "\n",
+ "BO can proceed on any data source (experiment or computation), any computer code (molecular mechanics or ab initio approaches), on any number of cores (it has an API and command line interface) and can fit a GP surrogate model of any or many materials properties extracted from calculations (energy, band gaps, HOMO/LUMO). Flexible use of external functions is facilitated by the BOSS **user function**, where the user transforms the BO sampling across the defined search domain into an atomistic configuration, the calculation is performed, and the desired value is parsed and returned.\n",
+ "\n",
+ "<img src=\"../figures/user_fn.png\" width=\"700px\"/>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Molecular conformer structure search\n",
+ "\n",
+ "Let us consider the problem of alanine conformational structure search, as shown in the image above. We identify 4 degrees of freedom, described by dihedral angles `[d1, d2, d3, d4]`. We will use BO to refine a GPR model mapping the state vector of dihedral angles to conformer energy `E`, as computed by an external total energy code. **Our objective is to find the state vector (configuration) that minimizes the total energy.**\n",
+ "\n",
+ "First, we load the required python packages."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 71,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import GPy\n",
+ "\n",
+ "from boss.bo.bo_main import BOMain\n",
+ "from boss.pp.pp_main import PPMain"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "In the conformer structure research, we used the AMBER code and `gaff` force field to compute the total energy for each sampled configuration. To simplify this tutorial, we construted an **AMBER emulator function** to simulate AMBER energies without the code itself. We did this by building a GPR generative model based on AMBER data, then encoding it in an AMBER emulator function."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 61,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Here we load AMBER data from a file and build an GPR model. \n",
+ "def load_model(filename):\n",
+ " \"\"\"Recreates a GPR model from saved parameters and data. \"\"\"\n",
+ "\n",
+ " # load saved data\n",
+ " data = np.load(filename)\n",
+ " dim = data[\"X\"].shape[1]\n",
+ "\n",
+ " # create kernel and mean functions\n",
+ " kernel = GPy.kern.StdPeriodic(input_dim=dim, ARD1=True, ARD2=True)\n",
+ " mean_func = GPy.mappings.Constant(dim, 1)\n",
+ "\n",
+ " # create model\n",
+ " model = GPy.models.GPRegression(\n",
+ " data[\"X\"], data[\"Y\"], kernel=kernel, mean_function=mean_func\n",
+ " )\n",
+ "\n",
+ " # set model params\n",
+ " model[:] = data[\"params\"]\n",
+ " model.fix()\n",
+ " model.parameters_changed()\n",
+ "\n",
+ " return model\n",
+ "\n",
+ "# Here, we define the AMBER emulator for a 2D example\n",
+ "AMBER_emulator2D = load_model(\"../data/model_2D_E0.npz\")\n",
+ "\n",
+ "# Here, we define the utility function that retrieves data from the AMBER emulator.\n",
+ "def f(X):\n",
+ " return AMBER_emulator2D.predict(np.atleast_2d(X))[0]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Running BO"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Once more, we need to set up the BO run by defining the:\n",
+ "\n",
+ "* `bounds` of the search domain: the number of [lower, upper] bound pairs in the vector will tell BOSS that this is a 2D search. \n",
+ "* `yrange` of the expected function\n",
+ "* `kernel`: since all dihedral angles are 360-periodic, we can set the kernel to *standard periodic*\n",
+ "* number of initial points `initpts` and BO iterations `iterpts`\n",
+ "\n",
+ "After setting BO, we are ready to run!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 66,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting file 'boss.out'\n",
+ "warning: overwriting file 'boss.rst'\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Here, we define the BO run.\n",
+ "#bounds = [[-50.0, 310.0]] * 2 # It can be done like this, or explicitly below \n",
+ "bounds = [[-50.0, 310.0], [-50.0, 310.0]]\n",
+ "\n",
+ "bo = BOMain(\n",
+ " f,\n",
+ " bounds,\n",
+ " yrange=[0, 1],\n",
+ " kernel=\"stdp\",\n",
+ " initpts=5,\n",
+ " iterpts=10,\n",
+ ")\n",
+ "\n",
+ "# Here we perform BO, it should take about 2 minutes.\n",
+ "res = bo.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We immediately run post-processing to check results:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 68,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting directory 'postprocessing'\n"
+ ]
+ },
+ {
+ "ename": "ValueError",
+ "evalue": "could not broadcast input array from shape (5,1) into shape (4,1)",
+ "output_type": "error",
+ "traceback": [
+ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+ "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)",
+ "\u001b[0;32m<ipython-input-68-1e6b51289a2a>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0mpp\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mPPMain\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mres\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpp_models\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mTrue\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpp_iters\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m2\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0;31m#pp = PPMain(res, pp_models=True)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mpp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrun\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
+ "\u001b[0;32m~/Code/BOSS/BOSS_tutorial/env/lib/python3.6/site-packages/boss/pp/pp_main.py\u001b[0m in \u001b[0;36mrun\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 158\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mload_data\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 159\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdump_data\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 160\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 161\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 162\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mload_data\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+ "\u001b[0;32m~/Code/BOSS/BOSS_tutorial/env/lib/python3.6/site-packages/boss/pp/pp_main.py\u001b[0m in \u001b[0;36mplot\u001b[0;34m(self, target)\u001b[0m\n\u001b[1;32m 367\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 368\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mtarget\u001b[0m \u001b[0;32min\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m\"model\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"all\"\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 369\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_plot_model\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 370\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 371\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mtarget\u001b[0m \u001b[0;32min\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m\"truef\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"all\"\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+ "\u001b[0;32m~/Code/BOSS/BOSS_tutorial/env/lib/python3.6/site-packages/boss/pp/pp_main.py\u001b[0m in \u001b[0;36m_plot_model\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 493\u001b[0m \u001b[0mbo\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 494\u001b[0m \u001b[0mbo\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_all_params\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 495\u001b[0;31m \u001b[0mcurr_xhat\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 496\u001b[0m )\n\u001b[1;32m 497\u001b[0m mdata = np.loadtxt(\n",
+ "\u001b[0;32m~/Code/BOSS/BOSS_tutorial/env/lib/python3.6/site-packages/boss/io/dump.py\u001b[0m in \u001b[0;36mdump_model\u001b[0;34m(settings, dest_file, bo, mod_params, xhat)\u001b[0m\n\u001b[1;32m 46\u001b[0m \u001b[0mp\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0minsert\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfloat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mbo\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_mu\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 47\u001b[0m \u001b[0mp\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0minsert\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfloat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mbo\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_nu\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 48\u001b[0;31m \u001b[0mmodel_data\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0minsert\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel_data\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel_data\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mp\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0maxis\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 49\u001b[0m titleLine = \"# Model output (x mu nu)\" + \", grid of %ix%i=%i pts\" % (\n\u001b[1;32m 50\u001b[0m \u001b[0mnpts\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+ "\u001b[0;32m<__array_function__ internals>\u001b[0m in \u001b[0;36minsert\u001b[0;34m(*args, **kwargs)\u001b[0m\n",
+ "\u001b[0;32m~/Code/BOSS/BOSS_tutorial/env/lib/python3.6/site-packages/numpy/lib/function_base.py\u001b[0m in \u001b[0;36minsert\u001b[0;34m(arr, obj, values, axis)\u001b[0m\n\u001b[1;32m 4576\u001b[0m \u001b[0mnew\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mtuple\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mslobj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0marr\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mtuple\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mslobj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 4577\u001b[0m \u001b[0mslobj\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0maxis\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mslice\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mindex\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mindex\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0mnumnew\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 4578\u001b[0;31m \u001b[0mnew\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mtuple\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mslobj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mvalues\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 4579\u001b[0m \u001b[0mslobj\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0maxis\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mslice\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mindex\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0mnumnew\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 4580\u001b[0m \u001b[0mslobj2\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0mslice\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m*\u001b[0m \u001b[0mndim\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+ "\u001b[0;31mValueError\u001b[0m: could not broadcast input array from shape (5,1) into shape (4,1)"
+ ]
+ }
+ ],
+ "source": [
+ "#pp = PPMain(res, pp_models=True, pp_iters=[0, 2])\n",
+ "pp = PPMain(res, pp_models=True)\n",
+ "pp.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "zKy8LL1oC8Iu"
+ },
+ "source": [
+ "## Results and Analysis\n",
+ "We start by checking **data acqusition** plots. Since we have 2 dimensions, the lower plot now features the sampling locations in both dihedrals. Amber energy values are all positive, and a global minimum was clearly idenfitifed. On the right, the **hyperparameter plot** indicates that the model has also converged.\n",
+ "\n",
+ "<table><tr><td><img src='postprocessing/acquisition_locations.png' width=\"400px\"></td><td><img src='postprocessing/hyperparameters.png' width=\"400px\"></td></tr></table>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Let us review the potential energy surface (PES) landscapes simulated by GPR after 20 and 40 BO points.\n",
+ "<table><tr><td><img src='postprocessing/graphs_models/it0005_npts0010.png' width=\"400px\"></td><td><img src='postprocessing/graphs_models/it0010_npts0015.png' width=\"400px\"></td></tr></table>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "<img src=\"postprocessing/graphs_models/it0020_npts0025.png\" width=\"400px\">"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "jupytext": {
+ "cell_metadata_filter": "-all",
+ "main_language": "python",
+ "notebook_metadata_filter": "-all"
+ },
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.6.3"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/Tutorial2_Structure-Search/BO-Structure-Search.ipynb b/Tutorial2_Structure-Search/BO-Structure-Search.ipynb
new file mode 100755
index 0000000000000000000000000000000000000000..984c9bcadae24cb0bc90c26e6caf6bfa9f3deba2
--- /dev/null
+++ b/Tutorial2_Structure-Search/BO-Structure-Search.ipynb
@@ -0,0 +1,243 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Tutorial 2: Configurational Structure Search with BO"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "This tutorial demonstrates how to use BO to optimize the atomistic configuration of materials. This can be done by mapping the configurational degrees of freedom to the total energy landscape, and inferring the global minimum of energy. The energy landscapes, refined by iterative data acquisitons, can be employed to visualise low energy configurations (they are chemically interpretable). \n",
+ "\n",
+ "<img src=\"../figures/ss_scenarios.png\" width=\"800px\"/>\n",
+ "\n",
+ "At the outset, the user should consider which **key degrees of freedom** control the energetics of the research problem. In the case of adsorption of organic molecules to the surface, these could be the molecule position or orientation above the surface. The molecular conformation can often be kept fixed. In the case of flexible molecules, we need to conduct a conformer search and minimize the energy with respect to dihedral angles. The degrees of freedom and their bounds define the N-dimensional BO structure search.\n",
+ "\n",
+ "BO can proceed on any data source (experiment or computation), any computer code (molecular mechanics or ab initio approaches), on any number of cores (it has an API and command line interface) and can fit a GP surrogate model of any or many materials properties extracted from calculations (energy, band gaps, HOMO/LUMO). Flexible use of external functions is facilitated by the BOSS **user function**, where the user transforms the BO sampling across the defined search domain into an atomistic configuration, the calculation is performed, and the desired value is parsed and returned.\n",
+ "\n",
+ "<img src=\"../figures/user_fn.png\" width=\"700px\"/>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Molecular conformer structure search\n",
+ "\n",
+ "Let us consider the problem of alanine conformational structure search, as shown in the image above. We identify 4 degrees of freedom, described by dihedral angles `[d1, d2, d3, d4]`. We will use BO to refine a GPR model mapping the state vector of dihedral angles to conformer energy `E`, as computed by an external total energy code. **Our objective is to find the state vector (configuration) that minimizes the total energy.** \n",
+ "\n",
+ "We will first consider a 2D conformer search (other variables kept fixed) and analyse the results, then proceed with a 4D conformer search. <br>\n",
+ "To begin, let's load the required python packages."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 28,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import GPy\n",
+ "\n",
+ "from boss.bo.bo_main import BOMain\n",
+ "from boss.pp.pp_main import PPMain"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "In this example, we used the AMBER code and `gaff` force field to compute the total energy for each sampled configuration. To simplify this tutorial, we construted an **AMBER emulator function** to simulate AMBER energies without the code itself. We did this by building a GPR generative model based on AMBER data, then encoding it in an AMBER emulator function."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 31,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Here we load AMBER data from a file and build an GPR model. \n",
+ "def load_model(filename):\n",
+ " \"\"\"Recreates a GPR model from saved parameters and data. \"\"\"\n",
+ "\n",
+ " # load saved data\n",
+ " data = np.load(filename)\n",
+ " dim = data[\"X\"].shape[1]\n",
+ "\n",
+ " # create kernel and mean functions\n",
+ " kernel = GPy.kern.StdPeriodic(input_dim=dim, ARD1=True, ARD2=True)\n",
+ " mean_func = GPy.mappings.Constant(dim, 1)\n",
+ "\n",
+ " # create model\n",
+ " model = GPy.models.GPRegression(\n",
+ " data[\"X\"], data[\"Y\"], kernel=kernel, mean_function=mean_func\n",
+ " )\n",
+ "\n",
+ " # set model params\n",
+ " model[:] = data[\"params\"]\n",
+ " model.fix()\n",
+ " model.parameters_changed()\n",
+ "\n",
+ " return model\n",
+ "\n",
+ "# Here, we define the AMBER emulator for a 2D structure search in d1-d4\n",
+ "AMBER_emulator = load_model(\"../data/model_2D_E0.npz\")\n",
+ "# Here, we define the AMBER emulator for a 4D structure search in d1-d2-d3-d4 \n",
+ "#AMBER_emulator = load_model(\"../data/model_4D_E0.npz\")\n",
+ "\n",
+ "# Here, we define the utility function that retrieves data from the AMBER emulator.\n",
+ "def f(X):\n",
+ " return AMBER_emulator.predict(np.atleast_2d(X))[0]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Running BO"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Once more, we need to set up the BO run by defining the:\n",
+ "\n",
+ "* `bounds` of the search domain: the number of [lower, upper] bound pairs in the vector will tell BOSS that this is a 2D search. \n",
+ "* `yrange` of the expected function\n",
+ "* `kernel`: since all dihedral angles are 360-periodic, we can set the kernel to *standard periodic*\n",
+ "* number of initial points `initpts` and BO iterations `iterpts`\n",
+ "\n",
+ "After setting BO, we are ready to run!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 32,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting file 'boss.out'\n",
+ "warning: overwriting file 'boss.rst'\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Here, we define the BO run.\n",
+ "bounds = [[-50.0, 310.0]] * 2 # It can be done like this, or explicitly below \n",
+ "#bounds = [[-50.0, 310.0], [-50.0, 310.0]]\n",
+ "\n",
+ "bo = BOMain(\n",
+ " f,\n",
+ " bounds,\n",
+ " yrange=[0, 1],\n",
+ " kernel=\"stdp\",\n",
+ " initpts=5,\n",
+ " iterpts=40,\n",
+ ")\n",
+ "\n",
+ "# Here we perform BO, it should take about 2.5 minutes.\n",
+ "res = bo.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We immediately run post-processing to check results. We set additional options to:\n",
+ "\n",
+ "* `pp_iters`: postprocess only BO steps 0, 10, 20 and 40 for convenience\n",
+ "* `pp_model_slice`: explicitly request a model slice in [dim1, dim2, Npts], which is useful to visualise N-dimensional models"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 33,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "warning: overwriting directory 'postprocessing'\n"
+ ]
+ }
+ ],
+ "source": [
+ "pp = PPMain(res, pp_models=True, pp_iters=[0, 10, 20, 40], pp_model_slice=[1, 2, 50])\n",
+ "pp.run()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "zKy8LL1oC8Iu"
+ },
+ "source": [
+ "## Results and Analysis\n",
+ "We start by checking **data acqusition** plots. Since we have 2 dimensions, the lower plot now features the sampling locations in both dihedrals. Amber energy values are all positive, and a global minimum was clearly idenfitifed. On the right, the **hyperparameter plot** indicates that the model has also converged.\n",
+ "\n",
+ "<table><tr><td><img src='postprocessing/acquisition_locations.png' width=\"400px\"></td><td><img src='postprocessing/hyperparameters.png' width=\"400px\"></td></tr></table>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Let us review the potential energy surface (PES) landscapes simulated by GPR after 10, 20 and 40 BO points. We can see that **the global mininum was found with few datapoints**, but more sampling is needed to converge the energy landscapes.\n",
+ "<table><tr><td><img src='postprocessing/graphs_models/it0010_npts0015.png' width=\"400px\"></td><td><img src='postprocessing/graphs_models/it0020_npts0025.png' width=\"400px\"></td><td><img src='postprocessing/graphs_models/it0040_npts0045.png' width=\"400px\"></td></tr></table>"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Exercise:\n",
+ "Re-do the conformer search with **all 4 degrees of freedom.** Note that you have to load the 4D dataset for the AMBER emulator function, and adjust the bounds for 4 dimensions instead of 2. You can execute the BO run with the same settings. In 4D, BO should take about 5min.\n",
+ "\n",
+ "Re-plot all the graphs above. In a 4D structure search, domain volume to cover by sampling is much bigger than in 2D, so the same amount of sampling results in poorer models. \n",
+ "\n",
+ "* How does the convergence look for the global minimum prediction and the model hyperparameters? \n",
+ "* What about the energy landscapes?\n",
+ "\n",
+ "When post-processing, you can now choose to slice the 4D surrogate model in the plane of any dihedral pair with the `pp_model_slice` keyword. The cross-section is always extracted in the plane of the global minimum solution (you can choose other planes via keywords). Once the BO run is complete, you can rerun the post-processing in different ways, but note that the `postprocessing` directory will be overwritten unless you rename it.\n",
+ "\n",
+ "* Dihedral variable d7 has a 120deg symmetry which can be accounted for with a shorter domain boundary - is this visible in the energy landscapes? \n",
+ "* The 2D model was produced for d1-d4 pair, which dominates conformer energetics. If you slice the 4D model in the d1-d4 plane, does it look similar? "
+ ]
+ }
+ ],
+ "metadata": {
+ "jupytext": {
+ "cell_metadata_filter": "-all",
+ "main_language": "python",
+ "notebook_metadata_filter": "-all"
+ },
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.6.3"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/data/model_2D_E0.npz b/data/model_2D_E0.npz
new file mode 100644
index 0000000000000000000000000000000000000000..5d0df6c570b8bd00c9e47347be9f3205f2dc394d
Binary files /dev/null and b/data/model_2D_E0.npz differ
diff --git a/data/model_2D_E1.npz b/data/model_2D_E1.npz
new file mode 100644
index 0000000000000000000000000000000000000000..cb2642c886b9ed2a07340f91e350192dc0283dc4
Binary files /dev/null and b/data/model_2D_E1.npz differ
diff --git a/data/model_4D_E0.npz b/data/model_4D_E0.npz
new file mode 100644
index 0000000000000000000000000000000000000000..374dd6aa49f1ac51ba70fae87928ce4d4872f0f8
Binary files /dev/null and b/data/model_4D_E0.npz differ
diff --git a/data/model_4D_E1.npz b/data/model_4D_E1.npz
new file mode 100644
index 0000000000000000000000000000000000000000..e33b4fe979c027097bc66d298b0e9cc98ca8a7c9
Binary files /dev/null and b/data/model_4D_E1.npz differ
diff --git a/figures/GPR_example.png b/figures/GPR_example.png
new file mode 100644
index 0000000000000000000000000000000000000000..b463cec1f90bc2f560cede0dd644930e1ff3751b
Binary files /dev/null and b/figures/GPR_example.png differ
diff --git a/figures/ss_scenarios.png b/figures/ss_scenarios.png
new file mode 100755
index 0000000000000000000000000000000000000000..cd4c4e9d5f4051d04ae7f4a2bbbcddc0b56c7d9a
Binary files /dev/null and b/figures/ss_scenarios.png differ
diff --git a/figures/user_fn.png b/figures/user_fn.png
new file mode 100644
index 0000000000000000000000000000000000000000..67d4d03e7a679cbafbe1ee24787e41a50ea0747c
Binary files /dev/null and b/figures/user_fn.png differ