{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Model selection\n",
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "# Problem definition\n",
    "from probeye.definition.inverse_problem import InverseProblem\n",
    "from probeye.definition.forward_model import ForwardModelBase\n",
    "from probeye.definition.distribution import Normal\n",
    "from probeye.definition.sensor import Sensor\n",
    "from probeye.definition.likelihood_model import GaussianLikelihoodModel\n",
    "from probeye.definition.correlation_model import ExpModel\n",
    "\n",
    "# Inference\n",
    "from probeye.inference.dynesty.solver import DynestySolver"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problem setup"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fixed parameters\n",
    "I = 1e9  # mm^4\n",
    "L = 10_000  # mm\n",
    "Q_loc = 2500  # mm\n",
    "\n",
    "# Measurements\n",
    "x_sensors = [2500, 5000]  # mm  (always from lower to higher, bug in inv_cov_vec_1D in tripy)\n",
    "d_sensors = [35, 50]  # mm\n",
    "sigma_model = 2.5  # mm\n",
    "pearson = 0.5\n",
    "l_corr = -np.abs(x_sensors[1] - x_sensors[0]) / np.log(pearson)  # mm (assuming exponential correlation)\n",
    "\n",
    "# Prior\n",
    "E_mean = 60  # GPa\n",
    "E_std = 20  # GPa\n",
    "Q_mean = 60  # kN\n",
    "Q_std = 30  # kN"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define forward model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def beam_deflection(E, Q, a, x):  # a is load position, x is sensor position\n",
    "    if x < a:\n",
    "        b = L - a\n",
    "        return Q * b * x * (L ** 2 - b ** 2 - x ** 2) / (6 * E * I * L)\n",
    "    \n",
    "    return Q * a * (L - x) * (2 * L * x - x ** 2 - a ** 2) / (6 * E * I * L)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class BeamModel(ForwardModelBase):\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define inverse problem"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Instantiate the inverse problem\n",
    "problem = InverseProblem(\"Beam model with two sensors\", print_header=False)\n",
    "\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Solve with Dynesty"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dynesty_solver = DynestySolver(problem, show_progress=True)\n",
    "inference_data = dynesty_solver.run(estimation_method='static', nlive=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model comparison"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "log_evidence_1 = ...  # See value in previous notebook\n",
    "log_evidence_2 = ...  # See value in this notebook\n",
    "\n",
    "prior_M_1 = ...\n",
    "prior_M_2 = ...\n",
    "\n",
    "denominator = ...\n",
    "\n",
    "posterior_M_1 = ...\n",
    "posterior_M_2 = ...\n",
    "\n",
    "print(f\"Posterior probability model 1: {posterior_M_1}\")\n",
    "print(f\"Posterior probability model 2: {posterior_M_2}\")\n",
    "\n",
    "bayes_factor = ...\n",
    "print(f\"Bayes factor: {bayes_factor}\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "venv",
   "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.11.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}