{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Fitting multidiffusion profiles for three components\n\nA typical fitting procedure is shown in this example. Synthetic diffusion\nprofiles are generated for a system with three components, and three different\ndifferent exchange experiments. Uphill diffusion is observed for one of the\nexchange experiments. For a moderate noise, we see that the fitting procedure\nresults in an accurate estimate of the diffusion matrix.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import numpy as np\nimport matplotlib.pyplot as plt\nfrom multidiff import compute_diffusion_matrix, create_diffusion_profiles"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "First we generate synthetic data, with three components. The diffusion\nmatrix has two eigenvalues with quite different magnitudes. Exchange vectors\ncorrespond to the exchange of the three possible pairs of components\n(meaning that each endmember is enriched in one component and poorer in\nanother, compared to the other endmember).\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "n_comps = 2\n\ndiags = np.array([1, 5])\nP = np.array([[1, 1], [-1, 0]])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "It is possible to have different measurement points for the different\nexperiments. Also note that measurement points don't have to be symmetric\naround 0.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "xpoints_exp1 = np.linspace(-10, 10, 100)\nxpoints_exp2 = np.linspace(-12, 10, 100)\nxpoints_exp3 = np.linspace(-8, 10, 120)\nx_points = [xpoints_exp1, xpoints_exp2, xpoints_exp3]"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Let us now define the exchange vectors. Each column corresponds to an\nexperiment: e.g. the first column represents the exchange of the first two\ncomponents.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "exchange_vectors = np.array([[0, 1, 1], \n               [1, -1, 0],\n               [-1, 0, -1]])\n\n\nconcentration_profiles = create_diffusion_profiles((diags, P), x_points,\n                                                    exchange_vectors,\n                                                    noise_level=0.0)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The algorithm needs an initial guess for the diffusion matrix. Here we give\nan initialization that is quite far from the looked-for diffusion matrix.\nNevertheless, the result of the fit is quite good.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "diags_init = np.array([1, 1])\nP_init = np.eye(2)\ndiags_res, eigvecs, _, _, _ = compute_diffusion_matrix((diags_init, P_init), \n                               x_points,\n                               concentration_profiles, plot=True,\n                               labels=['1', '2', '3'])\n\nprint(\"True eigenvalues: \", diags)\n\nprint(\"Fitted eigenvalues: \", diags_res)\n\nprint(\"True eigenvectors: \", P)\n\nprint(\"Fitted eigenvalues: \", eigvecs)\n\nplt.show()"
      ]
    }
  ],
  "metadata": {
    "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.8.13"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}