{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "# Custom Models in pycalphad: Viscosity"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Viscosity Model Background\n",
    "\n",
    "We are going to take a CALPHAD-based property model from the literature and use it to predict the viscosity of Al-Cu-Zr liquids.\n",
    "\n",
    "For a binary alloy liquid under small undercooling, Gąsior suggested an entropy model of the form\n",
    "$$\\eta = (\\sum_i x_i \\eta_i ) (1 - 2\\frac{S_{ex}}{R})$$\n",
    "\n",
    "where $\\eta_i$ is the viscosity of the element $i$, $x_i$ is the mole fraction, $S_{ex}$ is the excess entropy, and $R$ is the gas constant.\n",
    "\n",
    "For more details on this model, see \n",
    "\n",
    "1. M.E. Trybula, T. Gancarz, W. Gąsior, *Density, surface tension and viscosity of liquid binary Al-Zn and ternary Al-Li-Zn alloys*, Fluid Phase Equilibria 421 (2016) 39-48, [doi:10.1016/j.fluid.2016.03.013](http://dx.doi.org/10.1016/j.fluid.2016.03.013).\n",
    "\n",
    "2. Władysław Gąsior, *Viscosity modeling of binary alloys: Comparative studies*, Calphad 44 (2014) 119-128, [doi:10.1016/j.calphad.2013.10.007](http://dx.doi.org/10.1016/j.calphad.2013.10.007).\n",
    "\n",
    "3. Chenyang Zhou, Cuiping Guo, Changrong Li, Zhenmin Du, *Thermodynamic assessment of the phase equilibria and prediction of glass-forming ability of the Al–Cu–Zr system*, Journal of Non-Crystalline Solids 461 (2017) 47-60, [doi:10.1016/j.jnoncrysol.2016.09.031](https://doi.org/10.1016/j.jnoncrysol.2016.09.031)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:46.402127Z",
     "iopub.status.busy": "2022-02-19T19:24:46.402127Z",
     "iopub.status.idle": "2022-02-19T19:24:49.201546Z",
     "shell.execute_reply": "2022-02-19T19:24:49.201546Z"
    }
   },
   "outputs": [],
   "source": [
    "from pycalphad import Database"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## TDB Parameters\n",
    "We can calculate the excess entropy of the liquid using the Al-Cu-Zr thermodynamic database from Zhou et al.\n",
    "\n",
    "We add three new parameters to describe the viscosity (in Pa-s) of the pure elements Al, Cu, and Zr:\n",
    "```\n",
    "   $ Viscosity test parameters\n",
    "   PARAMETER ETA(LIQUID,AL;0) 2.98150E+02  +0.000281*EXP(12300/(8.3145*T));   6.00000E+03   \n",
    "  N REF:0 !\n",
    "   PARAMETER ETA(LIQUID,CU;0) 2.98150E+02  +0.000657*EXP(21500/(8.3145*T));   6.00000E+03   \n",
    "  N REF:0 !\n",
    "  PARAMETER ETA(LIQUID,ZR;0) 2.98150E+02  +4.74E-3 - 4.97E-6*(T-2128) ;   6.00000E+03   \n",
    "    N REF:0 !\n",
    "```\n",
    "\n",
    "Great! However, if we try to load the database now, we will get an error. This is because `ETA` parameters are not supported by default in pycalphad, so we need to tell pycalphad's TDB parser that \"ETA\" should be on the list of supported parameter types."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.210107Z",
     "iopub.status.busy": "2022-02-19T19:24:49.201546Z",
     "iopub.status.idle": "2022-02-19T19:24:49.261598Z",
     "shell.execute_reply": "2022-02-19T19:24:49.261598Z"
    },
    "jupyter": {
     "outputs_hidden": false
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Failed while parsing:     PARAMETER ETA(LIQUID,AL;0) 2.98150E+02 +0.000281*EXP(12300/(8.3145*T)); 6.00000E+03 N REF:0 \n",
      "Tokens: None\n",
      "Expected {{\"ELEMENT\" W:(ABCD...) W:(ABCD...) Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)') Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)') Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)') LineEnd} | {\"SPECIES\" W:(ABCD...) [Suppress:(\"%\")] Group:({{W:(ABCD...) [Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)')]}}...) [{Suppress:(\"/\") W:(+-01...)}] LineEnd} | {\"TYPE_DEFINITION\" Suppress:(<SP><TAB><CR><LF>) !W:( !) SkipTo:(LineEnd)} | {\"FUNCTION\" W:(ABCD...) {{Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)') | [\",\"]...} {{SkipTo:(\";\") Suppress:(\";\") [Suppress:(\",\")]... [Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)')] Suppress:([W:(Yy)])}}... Suppress:([W:(Nn)])} [Suppress:(W:(ABCD...))] LineEnd} | {\"ASSESSED_SYSTEMS\" SkipTo:(LineEnd)} | {\"DEFINE_SYSTEM_DEFAULT\" SkipTo:(LineEnd)} | {\"DEFAULT_COMMAND\" SkipTo:(LineEnd)} | {\"DATABASE_INFO\" SkipTo:(LineEnd)} | {\"VERSION_DATE\" SkipTo:(LineEnd)} | {\"REFERENCE_FILE\" SkipTo:(LineEnd)} | {\"ADD_REFERENCES\" SkipTo:(LineEnd)} | {\"LIST_OF_REFERENCES\" SkipTo:(LineEnd)} | {\"TEMPERATURE_LIMITS\" SkipTo:(LineEnd)} | {\"PHASE\" W:(ABCD...) Suppress:(<SP><TAB><CR><LF>) !W:( !) Suppress:(<SP><TAB><CR><LF>) Suppress:(W:(0123...)) Group:({Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)')}...) Suppress:(SkipTo:(LineEnd))} | {\"CONSTITUENT\" W:(ABCD...) Suppress:(<SP><TAB><CR><LF>) Suppress:(\":\") Group:(Group:({{[Suppress:(\",\")] {W:(ABCD...) [Suppress:(\"%\")]}}}...) [: Group:({{[Suppress:(\",\")] {W:(ABCD...) [Suppress:(\"%\")]}}}...)]...) Suppress:(\":\") LineEnd} | {\"PARAMETER\" {\"BMAGN\" | \"DF\" | \"DQ\" | \"ELRS\" | \"G\" | \"GD\" | \"L\" | \"MF\" | \"MQ\" | \"NT\" | \"SIGM\" | \"TC\" | \"THCD\" | \"THETA\" | \"V0\" | \"VA\" | \"VC\" | \"VISC\" | \"VK\" | \"VS\" | \"XI\"} Suppress:(\"(\") W:(ABCD...) [{Suppress:(\"&\") W:(ABCD...)}] Suppress:(\",\") Group:(Group:({{[Suppress:(\",\")] {W:(ABCD...) [Suppress:(\"%\")]}}}...) [: Group:({{[Suppress:(\",\")] {W:(ABCD...) [Suppress:(\"%\")]}}}...)]...) [{Suppress:(\";\") W:(0123...)}] Suppress:(\")\") {{Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)') | [\",\"]...} {{SkipTo:(\";\") Suppress:(\";\") [Suppress:(\",\")]... [Re:('[-+]?([0-9]+\\\\.(?!([0-9]|[eE])))|([0-9]*\\\\.?[0-9]+([eE][-+]?[0-9]+)?)')] Suppress:([W:(Yy)])}}... Suppress:([W:(Nn)])} [Suppress:(W:(ABCD...))] LineEnd}}, found '('  (at char 17), (line:1, col:18)\n"
     ]
    }
   ],
   "source": [
    "try:\n",
    "    dbf = Database('alcuzr-viscosity.tdb')\n",
    "except Exception as e:\n",
    "    print(e)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Adding the `ETA` parameter to the TDB parser "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.261598Z",
     "iopub.status.busy": "2022-02-19T19:24:49.261598Z",
     "iopub.status.idle": "2022-02-19T19:24:49.281502Z",
     "shell.execute_reply": "2022-02-19T19:24:49.281502Z"
    },
    "jupyter": {
     "outputs_hidden": false
    }
   },
   "outputs": [],
   "source": [
    "import pycalphad.io.tdb_keywords\n",
    "pycalphad.io.tdb_keywords.TDB_PARAM_TYPES.append('ETA')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "Now the database will load:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.281502Z",
     "iopub.status.busy": "2022-02-19T19:24:49.281502Z",
     "iopub.status.idle": "2022-02-19T19:24:49.341526Z",
     "shell.execute_reply": "2022-02-19T19:24:49.341526Z"
    },
    "jupyter": {
     "outputs_hidden": false
    }
   },
   "outputs": [],
   "source": [
    "dbf = Database('alcuzr-viscosity.tdb')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Writing the Custom Viscosity Model\n",
    "\n",
    "Now that we have our `ETA` parameters in the database, we need to write a `Model` class to tell pycalphad how to compute viscosity. All custom models are subclasses of the pycalphad `Model` class.\n",
    "\n",
    "When the `ViscosityModel` is constructed, the `build_phase` method is run and we need to construct the viscosity model after doing all the other initialization using a new method `build_viscosity`. The implementation of `build_viscosity` needs to do four things:\n",
    "1. Query the Database for all the `ETA` parameters\n",
    "2. Compute their weighted sum\n",
    "3. Compute the excess entropy of the liquid\n",
    "4. Plug all the values into the Gąsior equation and return the result\n",
    "\n",
    "Since the `build_phase` method sets the attribute `viscosity` to the `ViscosityModel`, we can access the property using `viscosity` as the output in pycalphad caluclations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.351408Z",
     "iopub.status.busy": "2022-02-19T19:24:49.351408Z",
     "iopub.status.idle": "2022-02-19T19:24:49.361677Z",
     "shell.execute_reply": "2022-02-19T19:24:49.361677Z"
    }
   },
   "outputs": [],
   "source": [
    "from tinydb import where\n",
    "from pycalphad import Model, variables as v\n",
    "\n",
    "class ViscosityModel(Model):\n",
    "    def build_phase(self, dbe):\n",
    "        super(ViscosityModel, self).build_phase(dbe)\n",
    "        self.viscosity = self.build_viscosity(dbe)\n",
    "\n",
    "    def build_viscosity(self, dbe):\n",
    "        if self.phase_name != 'LIQUID':\n",
    "            raise ValueError('Viscosity is only defined for LIQUID phase')\n",
    "        phase = dbe.phases[self.phase_name]\n",
    "        param_search = dbe.search\n",
    "        # STEP 1\n",
    "        eta_param_query = (\n",
    "            (where('phase_name') == phase.name) & \\\n",
    "            (where('parameter_type') == 'ETA') & \\\n",
    "            (where('constituent_array').test(self._array_validity))\n",
    "        )\n",
    "        # STEP 2\n",
    "        eta = self.redlich_kister_sum(phase, param_search, eta_param_query)\n",
    "        # STEP 3\n",
    "        excess_energy = self.GM - self.models['ref'] - self.models['idmix']\n",
    "        #liquid_mod = Model(dbe, self.components, self.phase_name)\n",
    "        ## we only want the excess contributions to the entropy\n",
    "        #del liquid_mod.models['ref']\n",
    "        #del liquid_mod.models['idmix']\n",
    "        excess_entropy = -excess_energy.diff(v.T)\n",
    "        ks = 2\n",
    "        # STEP 4\n",
    "        result = eta * (1 - ks * excess_entropy / v.R)\n",
    "        self.eta = eta\n",
    "        return result"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Performing Calculations\n",
    "\n",
    "Now we can create an instance of `ViscosityModel` for the liquid phase using the `Database` object we created earlier. We can verify this model has a `viscosity` attribute containing a symbolic expression for the viscosity."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.372613Z",
     "iopub.status.busy": "2022-02-19T19:24:49.372613Z",
     "iopub.status.idle": "2022-02-19T19:24:49.381328Z",
     "shell.execute_reply": "2022-02-19T19:24:49.381328Z"
    },
    "jupyter": {
     "outputs_hidden": false
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(0.000657*LIQUID0CU*2.71828182845905**(2585.84400745685*T**(-1.0)) + LIQUID0ZR*(0.00474 - 4.97e-06*(-2128.0 + T)))*(1 + 0.240543628600637*(LIQUID0CU*LIQUID0ZR*(392.8485 - 51.3121*log(T)) + (LIQUID0CU - LIQUID0ZR)*LIQUID0CU*LIQUID0ZR*(75.3798 - 9.6125*log(T)) + (LIQUID0CU - LIQUID0ZR)**2*LIQUID0CU*LIQUID0ZR*(-270.5305 + 36.8512*log(T)) + (LIQUID0CU - LIQUID0ZR)**3*LIQUID0CU*LIQUID0ZR*(105.895 - 13.6488*log(T)))/(LIQUID0CU + LIQUID0ZR))\n"
     ]
    }
   ],
   "source": [
    "mod = ViscosityModel(dbf, ['CU', 'ZR'], 'LIQUID')\n",
    "print(mod.viscosity)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "Finally we calculate and plot the viscosity."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.431215Z",
     "iopub.status.busy": "2022-02-19T19:24:49.411882Z",
     "iopub.status.idle": "2022-02-19T19:24:49.657823Z",
     "shell.execute_reply": "2022-02-19T19:24:49.657823Z"
    },
    "jupyter": {
     "outputs_hidden": false
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAGDCAYAAAA1cVfYAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAokUlEQVR4nO3dfbhcVX328e/NCUkUAwoExEAMJFirEIJNJRRrUYxFQIgFBUtULBK15VHER0uaCELJA/gS0bYKRG3RoAbfTiOJsWklWijEBhMSUJQAEUxpEkEggAESfs8fex8zOcyZWTNn9rzen+s6FzOz1+z5nX2Fuc9ea+21FRGYmZkNZbdWF2BmZu3NQWFmZhU5KMzMrCIHhZmZVeSgMDOzihwUZmZWkYPC2o6kOyUd24LPfVzSIc3+XLN256CwppK0TNIlZV4/RdL/ShoREa+MiBXNri0iXhAR9+b1/IukS4v4HElnSbqpSptPSbpb0lZJd0l656Dt10j6haRnJZ1V5v0fyo/nY5K+LGlUybYJkm6U9GS+7zeUbPu4pIUlz8flbT4nScP6xa1jOSis2a4FZpb50nkHcF1EbG9BTe3oCeDNwF7Au4DPSvqTku23A38N/HTwGyX9OXABcBzwUuAQ4OKSJl8HVgP7AHOAb0kaW2Y/LwV+DCyOiA+Er87tXRHhH/807Qd4HvAo8NqS114EbAOOyJ9vAN6QP341sAp4DNgEzC9532uA/wIeAR4Azspf3wv4CrAF+BUwF9gt3zYJ+FFew2+ARSX7i3z7LOAZ4GngceB7wEeAbw/6XT4HfHaI3/MC4B5gK/Az4C3563+Y/6478n0/knjcFgMfLvP6TQO/d8lrXwP+X8nz44D/zR+/DHgKGFOy/T+B9+WPPw4sBCbmx+6SVv+b8U/rf3xGYU0VEb8DrgdKu1LeBtwVEbeXectnyb6M9yT78roefv/X7veBfwDGAlOANfl7/oEsLA4B/iz/rHfn2/4e+DeycDowbzu4xmuA64BPRNYd9WayL8/jJb0w//wRwBlkgVTOPcCf5nVcDCyUdEBE/Bx4H3BLvu8XDvH+35P0POCPgTurtc29kuyMY8DtwP6S9sm33RsRWwdtf2XJ80PIziSujogLEz/TupiDwlrhWuA0SaPz5+/MXyvnGWCSpH0j4vGIuDV//S+Bf4+Ir0fEMxHxUESskdRH9gU+OyK2RsQG4NNkXVsD+3sp8JKI2BYRFccKBkTEg2Rfnm/NXzoe+E1E3DZE+29GxP9ExLMRsQi4m+zsqB5XkX2Z/yCx/QvIzpgGDDweU2bbwPYxJc8PA/YAFtVcqXUlB4U1Xf7l/BtghqSJZF+gXxui+dlk3SV3SfpvSSflrx9E9lf7YPsCu5N1mwz4FTAuf/xRQMBP8tlVf1VD6dcCM/PHM4GvDtVQ0jslrZH0iKRHyL58963hswb288n8vW+LiNQxgseBPUueDzzeWmbbwPbSM4zFwJeBH+ZnbtbjHBTWKl8hO5OYCfwgIjaVaxQRd0fE24H9gCvIBl73IBuTmFjmLb9h51nDgPHAxnx//xsR50TES4D3Ap+XNKncR5d5rR+YLOkw4CSy7qnnyL9cFwDnAvvk3Ut3kAXUUPsut5+LgTcBb4yIx1Lek7sTOKLk+RHApoh4KN92iKQxg7bv0q0VEecDN5CFxTispzkorFW+ArwBOIehu52QNFPS2Ih4lmzQGuBZsi/pN0h6m6QRkvaRNCUidpCNY8yTNCb/0j6fbIwBSW+VdGC+n9+SfWk/W+ajN5H11f9eRGwDvkV29vOTiLh/iLL3yPe7Jf/Md5OdFZTu+0BJIyv83rPJutfekH/BD94+Mu+6E7C7pNGSBv5//gpwtqRX5GMqc4F/yX+HX5KN5VyUv+ctwGTg22XKOBe4EfgPSfsPVav1gFaPpvund3+AFWRf1qMGvb6BnbOeFgKbybpM7gRmlLT7U2Al2YyoB4B35a+/KH/flvz1C9k56+kTZGcXj5N1Xc0q2V8Ak/LHh5J9oT4C9Je0eU3e7t1Vfrd5wMNkZzjzyWZavSffNhJYMrB9iPcH2eykx0t+/m7QsYtBP8eWbD+fLJAeA/659BgDE/L3/w74xcCxzrd9HFhY8nw3suBZB+zb6n8z/mnNj/J/DGaWQNJ44C7gxVFbd5BZx3LXk1mivGvnfOAbDgnrJSNaXYBZJ8gH0DeRzaA6vsXlmDWVu57MzKwidz2ZmVlFDgozM6uo48Yo9t1335gwYUKryzAz6yi33XbbbyLiOasEp+i4oJgwYQKrVq1qdRlmZh1F0q+qtyrPXU9mZlaRg8LMzCpyUJiZWUUOCjMzq8hBYWZmFTkozMysIgeFmZlV5KAwM7OKHBRmZlZRx12ZbWbWieb2r2PhreXvnrvh8hObXE1tHBRmZg125oJbuPmeh5PbT7hgyXNea6fwcFCYmTXIpNlL2N6gW/wMhMeVp09hxpHjGrPTOnXcjYumTp0aXhTQzNpJIwNiKMM9w5B0W0RMree9PqMwM6vT5IuW8dhTO5ryWRMuWMLMaeO5dMbhTfm8Up71ZGZWo/7VG5lwwZKmhcSAhbfeX3Y8o2gOCjOzGhw1bznnLVrT0hqaHRYOCjOzRBMuWMKmrU+3ugyguWHhMQozsypqne462AjB+suGHoyePn8Fd29+oub9TrhgSVOm0XrWk5lZBcMZsK5n8LmeM4WUsBjOrCd3PZmZDaHekJg5bTwbLj+xrhlKGy4/kStPn1LTe4ruhnJQmJmVUU9I7D9mZN0BUWrGkeNq7lJ6+Zylw/rMShwUZmaDvHzO0ppDYsPlJ7JyzvSG1rHh8hMZobS223YER81b3tDPH+CgMDMrMeGCJWzbkT52u+eovkIHlNdfdiJ7jupLartp69PM7V/X8BocFGZmuVr7+q88fQprLz6+oGp2Wnvx8ew/ZmRS26FWqB0OB4WZGbWHxIbLT2zqYn0r50xP7oZq9HhFoUEhaYOkdZLWSHrOnFZlPidpvaS1kl5VZD1mZuXUExKtUOlajFLbdgT9qzc27HObcUbxuoiYMsT83TcBh+Y/s4AvNKEeM7Pfq/Wv71bfJyL18xu5zEiru55OAb4SmVuBF0o6oMU1mVmPOHPBLTUNXLc6JAak1jH5omUN+byigyKAf5N0m6RZZbaPAx4oef7r/DUzs0LN7V9X07Ic7RISA2ZOG1+1TaNWty06KF4TEa8i62L6G0mvrWcnkmZJWiVp1ZYtWxpboZn1pFpmB7VbSADJF/Ud3ICrtgsNiojYmP93M/Bd4NWDmmwEDip5fmD+2uD9XBMRUyNi6tixY4sq18x6RC2D1+0YEgNSagsY9rUVhQWFpD0kjRl4DLwRuGNQs8XAO/PZT9OARyPiwaJqMjOr5S/sdg6JASldUMO9tqLIM4r9gZsk3Q78BFgSEcskvU/S+/I2S4F7gfXAAuCvC6zHzHrcpNlLSB267oSQgKwLKuX6ihF7jzu03s8o7H4UEXEvcESZ168qeRzA3xRVg5nZgMkXLWN7Ykp0SkgMWH/ZiVW703Yb+bw9691/q6fHmpkVbm7/uuQZQJ0WEgMO3W+PwvbtoDCzrpfaR1/kl23Rlp9/bGH7dlCYWVdLneE0uk+Fftk2Q+oqs7VyUJhZ16pleY675p1QYCXNUdRKtg4KM+tKc/vXJS/P0anjEuXUehvVFA4KM+tKqeMS3RQSQCFLnzsozKzrpI5LFPHXdzto9O/loDCzrpJ63+j9x4xs6o2HmqnRv5eDwsy6Rv/qjWza+nRS25VzphdcTWsdM3Hvhu3LQWFmXSP1Zj3dNi5RznXnHM3ovsR7p1bhoDCzrtDr4xLl3DXvhIZcRFjYWk9mZs0yff6KpHbHTNy7a8clhjJwEaGuOOm2evfhMwoz62j9qzdy9+Ynktped87RBVfTnRwUZtbRPC5RPAeFmXWsyRctS2rnkBgeB4WZdaT+1RuTlg4vaqG8XuKgMLOOlNrlVNRCeb3EQWFmHSd1Kqy7nBrDQWFmHeXMBbcktZs5bXzBlfQOB4WZdZSb73m4apvRfeLSGYc3oZre4KAws46R2uXUDTchaicOCjPrCKlTYd3l1HgOCjPrCClTYUcIdzkVwEFhZm3v4MQup/WXeZZTERwUZtbWzlxwCyl3vvZU2OI4KMysraXMctp/zMgmVNK7HBRm1rZSZzl1+93qWs1BYWZtKfXCOnc5Fc9BYWZtKaXLaURj7vRpVTgozKztpHY5eZZTczgozKytHDVveVI7dzk1j4PCzNrKpq1PV21zzMS9m1CJDXBQmFnbSL2wzve+bi4HhZm1hbn965IurPNaTs3noDCztrDw1vurttl/zEiv5dQCDgoza7mXz1ma1M4X1rWGg8LMWurMBbewbUf1Tid3ObWOg8LMWir1wjp3ObWOg8LMWib1ZkS+sK61HBRm1hL9qzcm3YzIXU6t56Aws5Y4b9Gaqm32HNXnLqc24KAws6ZLXRl27cXHF1yJpSg8KCT1SVot6YYy286StEXSmvznPUXXY2atlzKA7WU62seIJnzGB4GfA3sOsX1RRJzbhDrMrA2krgzrZTraR6FnFJIOBE4Evljk55hZZ/DNiDpT0V1PVwIfBZ6t0OZUSWslfUvSQQXXY2YtlNLldOh+ezShEqtFYUEh6SRgc0TcVqHZ94AJETEZWA5cO8S+ZklaJWnVli1bCqjWzIqW2uW0/Pxjiy3EalbkGcUxwMmSNgDfAF4vaWFpg4h4KCKeyp9+EfijcjuKiGsiYmpETB07dmyBJZtZEfpXb0xq5y6n9lRYUETE7Ig4MCImAGcAP4yImaVtJB1Q8vRkskFvM+syKddMuMupfTVj1tMuJF0CrIqIxcAHJJ0MbAceBs5qdj1mVqxJs93l1OmaEhQRsQJYkT++sOT12cDsZtRgZs03t38d2xPuRuQup/bmK7PNrDApNyNyl1P7c1CYWSHc5dQ9HBRm1nD9qzcmdTl5ZdjO4KAws4ZLmeU0uk9eGbZDOCjMrKGOmrc8qd1d804ouBJrFAeFmTVM/+qNbNr6dNV27nLqLA4KM2sYdzl1JweFmTWEu5y6l4PCzIYttcvJNyPqTA4KMxu2lC4n8M2IOpWDwsyGJbXLyct0dC4HhZnVzV1OvcFBYWZ1S+lyEu5y6nQOCjOry8vnLE1qd5+7nDqeg8LMaja3fx3bdlRfzOnK06cUX4wVzkFhZjVLWT5cwIwjxxVfjBXOQWFmNZl80bKkdu5y6h4OCjNL1r96I489taNqO9+MqLs4KMwsWeqFdb4ZUXdxUJhZkgkXpN2xzhfWdR8HhZlVdeaCW5LaeZZTd3JQmFlVN9/zcNU2e47q8yynLuWgMLOKUruc1l58fMGVWKs4KMxsSKlTYd3l1N0cFGY2pJSpsCPkC+u6nYPCzMpK7XJaf5lnOXU7B4WZPYe7nKyUg8LMdpF69bVnOfUOB4WZ7SL16mvPcuodDgoz+71Js331tT2Xg8LMgOweE9ur32KC/ceMLL4YaysOCjMD0u4xAbByzvSCK7F246AwMy/4ZxU5KMx63MGJIeGpsL3LQWHWw6bPX0HCsASj++SpsD1sRLUGknYDjgBeAvwOuCMiNhddmJkV7+7NTyS1u2veCQVXYu1syKCQNBH4W+ANwN3AFmA08DJJTwJXA9dGxLPNKNTMGsvjEpaq0hnFpcAXgPdGxC5np5L2A/4SeAdwbXHlmVkRUpfoGKGCC7GOMGRQRMTbK2zbDFxZREFmVqy5/euSlugAL/hnmaqD2ZLeKmlM/vhjkr4j6VXFl2ZmRUi9XsJdTjYgZdbTxyJiq6TXAMcBXyLrkjKzDpM6LuGpsFYqJSgGzlFPBK6JiCWAr+E36zCp6zgdut8engpru0gJio2SrgZOB5ZKGpX4PgAk9UlaLemGMttGSVokab2klZImJFduZskmzV6StI6TgOXnH1t0OdZhUr7w3wb8APjziHgE2Bv4SA2f8UHg50NsOxv4bURMAj4DXFHDfs0swZkLbkkKCYD7PC5hZVQNioh4MiK+ExF3S5oVEQ9GxL+l7FzSgWRdVl8coskp7Jxe+y3gOEmekGfWQDff83BSOw9e21BqXcLjfTW2vxL4KDDURXnjgAcAImI78CiwT42fYWZDSB28njltfMGVWCerNSiS/9qXdBKwOSJuq/Ezyu1rlqRVklZt2bJluLsz6wmpITFCcOmMwwuuxjpZrUHx5hraHgOcLGkD8A3g9ZIWDmqzETgIQNIIYC/gocE7iohrImJqREwdO3ZsjSWb9Z7UkABfVGfVVV0UEEDSicArgdEDQwgRcUml90TEbGB2/v5jgf8bETMHNVsMvAu4BTgN+OHg5ULMrDap02DB4xKWJmX12KuA5wOvIxuUPg34Sb0fKOkSYFVELCa7eO+rktYDDwNn1LtfM6tthpNDwlKlnFH8SURMlrQ2Ii6W9Gng+7V8SESsAFbkjy8seX0b8NZa9mVmQ0ud4XTMxL0LrsS6ScoYxbb8v09KegnwDHBAcSWZWT1SxyVG94nrzjm64Gqsm6ScUXxP0guBTwI/BQJYUGRRZlabWgavfRMiq1XFoJA0heyq6gMi4tv5MhyjI+LRZhRnZtWl3vMaPC5h9Rmy60nShcD1wKnAEknnRMRTDgmz9jFp9pKke16DQ8LqV+mM4nRgSkQ8KWkfYBnucjJrG7V0N3nZcBuOSoPZT0XEkwAR8VCVtmbWRLWExDET9/ay4TYslc4oDpG0OH8sYGLJcyLi5EIrM7OyUu93DbD/mJGe4WTDVikoThn0/FNFFmJm1U2+aFny/a73HNXHyjnTC67IesGQQRERP2pmIWZWWS0hIWDtxccXW5D1jKrjDpJOyu9Q97CkxyRtlfRYM4ozs8xR85YnhwT4BkTWWCkX3F0J/AWwzgv2mTXf9Pkr2LT16eT2ngZrjZYSFA8AdzgkzJrv5XOWsm1H+v96DgkrQkpQfBRYKulHwFMDL0bE/MKqMjMOviD9YjpwSFhxUoJiHvA4MBoYWWw5ZgYOCWsvKUHxkog4rPBKzAyo7WI6cEhY8VKutl4q6Y2FV2JmDglrSylB8X5gmaTfeXqsWXEcEtauqnY9RcSYZhRi1sscEtbOKi0zPqHSG5U5sOEVmfWQ/tUbHRLW9iqdUXxS0m7AvwK3AVvIZj5NAl4HHAdcBPy66CLNutH0+Su4e/MTNb3HIWGtUGmtp7dKegVwJvBXZPfJfpLsjndLgXkRsW2o95vZ0Gq9kA4cEtY6FccoIuJnwJwm1WLWE2rtagKHhLVWynUUZtYAc/vXsfDW+2t+n0PCWs1BYdYEtSwRPuDQ/fZg+fnHFlOQWQ0cFGYFq6er6crTp/j2pdY2qgaFpO8AXwK+HxHPFl+SWXeo5ywC3NVk7SflyuzPA38J3C3pckl/UHBNZh1vwgVLHBLWNVKuzP534N8l7QW8PX/8ALAAWBgRzxRco1nHqGfaK2T3t/atS61dpZxRIGkf4CzgPcBq4LPAq4DlhVVm1kEGrrCuJyRmThvvkLC2ljJG8V3gD4CvAm+OiAfzTYskrSqyOLN21796I+ctWlP3+93VZJ0gZdbTgohYWvqCpFER8VRETC2oLrO2V89spgGe+mqdJCUoLiVbsqPULWRdT2Y9ZzgBAT6LsM4zZFBIejEwDniepCMB5Zv2BJ7fhNrM2spwA0LAfQ4J60CVzij+nGwA+0BgfsnrW4G/K7Ams7ZR77Ibg/kswjpZpdVjrwWulXRqRHy7iTWZtVy9F8sNdszEvbnunKMbUJFZ61TqepoZEQuBCZLOH7w9IuaXeZtZx6rn/hBDGSFYf5nPIqw7VOp62iP/7wuaUYhZKzQyHAa4m8m6TaWup6vz/17cvHLMinfwBUuo/bK46mZOG8+lMw4vYM9mrZVywd0nyKbI/g5YBkwGPpR3S5m1vSLOGkr5DMK6Xcp1FG+MiI9KeguwAfgL4MeAg8LaTqNmKaXwUuDWK1KCYqDNicA3I+JRSZXamxVu0uwlbC+i/6iK0X3irnknNP+DzVooJShukHQXWdfT+yWNBbYVW5b1qqLGD4bLS25YL0tZZvyCfJzi0YjYIekJ4JRq75M0mqyLalT+Od+KiIsGtTkL+CSwMX/pHyPii7X9CtZOWvWXflE8/mCWNpi9OzATeG3e5fQj4KqEfT8FvD4iHs/3cZOk70fErYPaLYqIc2us25psuMtXdBKHg9muUrqevgDsTnanO4B35K+9p9KbIiKAx/Onu+c/XfS3ZvcZ7pLZnco3DTKrLCUo/jgijih5/kNJt6fsXFIfcBswCfiniFhZptmpkl4L/JJs2u0DKfu24eu2bqJUHm8wq01KUOyQNDEi7gGQdAiQtAhOROwApkh6IfBdSYdFxB0lTb4HfD0inpL0XuBa4PWD9yNpFjALYPz48SkfbWX06hmDu5LMhkdZD1GFBtJxwD8D95KtlPxS4N0RcWNNHyRdCDwZEZ8aYnsf8HBE7FVpP1OnTo1Vq3xjvVS9Eg4+SzCrTNJt9d5sLmXW039IOpTsdqgAv4iIpxKKGgs8ExGPSHoeMB24YlCbA0purXoy8POaqrchvXzO0rru39xufDZg1nops57+BrguItbmz18k6eyI+HyVtx5Atkx5H7AbcH1E3CDpEmBVRCwGPiDpZGA78DDZ/S+sTmcuuIWb73m41WXsYv8xI1k5Z3qryzCzYUjpeloTEVMGvbY6Io4ssrChuOvpuYpey6gcL19h1lkK7XoC+iQpn+46MJYwsp4Ps8ZqxrpG7voxs5SgWAYsknR1/vy9+WvWIkUNUHuZbDMrJyUo/pZsaur78+fLAS+z0SKNukXnAJ8xmFk1KbOeniVbsuMqSXsDB+bXR1gTNXKg2uFgZrVImfW0gmzq6giyq6w3S/qviPhQwbVZrhHrLLlbyczqldL1tFdEPCbpPcBXIuIiSWuLLsyGP5tJwH0+ezCzYUq6cZGkA4C3AXMKrsdywzmLGCFYf5kDwswaIyUoLgF+ANwUEf+dr/V0d7Fl9bbhhITHH8ys0VIGs78JfLPk+b3AqUUW1auGc12EA8LMijJkUEj6aER8QtI/UOY+EhHxgUIr6zFHzVvOpq1P1/w+L4ZnZkWrdEbxs/y/Xi+jYPVeG+GzCDNrhkpB8SZJv42Ia5tWTQ+q5+ZBPosws2aqFBS/BD6Vz3i6nuwGQ6ubU1ZvqGfQ2mcRZtZsuw21ISI+GxFHA38GPAR8WdJdki6S9LKmVdilag0J4ZAws9YYMigGRMSvIuKKfFnxtwMz8A2GhqXWkDhm4t6+cM7MWiZlCY8RwJuAM4DjgBXAxwutqovVGhI+izCzVqs0PXY62RnECcBPgG8AsyKiuXfI6SIOCTPrRJXOKGYDXwM+HBG/bVI9XWvSbIeEmXWmIYMiIl7fzEK62eSLltU0BdYhYWbtpOpgtg3PmQtuqeliOoeEmbUbB0XBarnZkEPCzNqRg6JAtQxeOyTMrF05KArikDCzbuGgKIBDwsy6iYOiwabPX5Hc1iFhZp3AQdFA/as3Jt/j2iFhZp3CQdFA5y1ak9TuytOnFFqHmVkjOSgaZPJFy5Laje4TM44cV3A1ZmaN46BogLn965IuqhNw17wTii/IzKyBHBQNsPDW+5PaealwM+tEDophSp0K68FrM+tUDophOGre8qR2x0zcu+BKzMyK46AYhk1bn05qd905RxdciZlZcRwUdXKXk5n1CgdFHV4+Z2lSO4eEmXUDB0Udtu2ofheiPUf1NaESM7PiOShqlHo2sfbi4wuuxMysORwUNehfvTHpbMJdTmbWTRwUNUhZy8lTYc2s2zgoEqWu5eSpsGbWbRwUiVLWcvKqsGbWjRwUCQ5OuGbCq8KaWbcqLCgkjZb0E0m3S7pT0sVl2oyStEjSekkrJU0oqp7hqD587VVhzax7FXlG8RTw+og4ApgCHC9p2qA2ZwO/jYhJwGeAKwqspy4pV2DPnDa+CZWYmbVGYUERmcfzp7vnP4P/OD8FuDZ//C3gOEkqqqZapQ5gXzrj8IIrMTNrnULHKCT1SVoDbAaWR8TKQU3GAQ8ARMR24FFgnyJrqkXKALbPJsys2xUaFBGxIyKmAAcCr5Z0WD37kTRL0ipJq7Zs2dLQGocyff6Kqm32HNXnswkz63pNmfUUEY8ANwKD17XYCBwEIGkEsBfwUJn3XxMRUyNi6tixYwuuNnP35ieqtvEyHWbWC4qc9TRW0gvzx88DpgN3DWq2GHhX/vg04IcRkTLJqFApYxO+AtvMesWIAvd9AHCtpD6yQLo+Im6QdAmwKiIWA18CvippPfAwcEaB9SRLGZvwFdhm1isKC4qIWAscWeb1C0sebwPeWlQN9fDZhJnZrnxl9iA+mzAz25WDosRR85ZXbePpsGbWaxwUJTZtfbpqG0+HNbNe46DIpdy5zjckMrNe5KDIVbtz3f5jRjapEjOz9uKgIO1sYuWc6U2oxMys/TgoqH42seeoviZVYmbWfno+KM5ccEvVNl6qw8x6Wc8Hxc33PFxxu8cmzKzX9XRQpJxNeGzCzHpdTwdFtbOJtrmDkplZC/VsUKRchf2Z06cUX4iZWZvr2aCodhV2n2DGkeOaVI2ZWfvqyaBIGZv49NumFF+ImVkH6MmgqDY2AT6bMDMb0JNBUY3vN2FmtlPPBUXKjYl8vwkzs516Liiq3ZjIZxNmZrvqqaBImRLrswkzs131VFBUmxJ76H57NKkSM7PO0TNBMX3+iqptlp9/bOF1mJl1mp4Jirs3P1Fxu5cSNzMrr2eCohovJW5mVl5PBEXKldhmZlZeTwRFtSuxZ04b36RKzMw6T9cHRcrZxKUzDm9CJWZmnanrg6La2cSoEV1/CMzMhqXnvyWvOHVyq0swM2trXR0UKd1OXiXWzKyyrg6Kat1OXtfJzKy6rg2Kuf3rqrbxuk5mZtV1bVAsvPX+itu9rpOZWZquDYpqvK6TmVmargyKlG4nMzNL05VBUa3b6crTpzSnEDOzLtCVQVGNp8SamaXruqDwAoBmZo3VdUFR7doJdzuZmdWmq4Kif/XGqm3c7WRmVpuuCorzF62puN3XTpiZ1a6rguLZKtt97YSZWe26Jiimz19RcbuaU4aZWdcpLCgkHSTpRkk/k3SnpA+WaXOspEclrcl/Lqz38+7e/ETF7Z/xILaZWV1GFLjv7cCHI+KnksYAt0laHhE/G9TuPyPipALrADyIbWZWr8LOKCLiwYj4af54K/BzoJBv66PmLS9it2ZmRpPGKCRNAI4EVpbZfLSk2yV9X9Ir69n/pq1PV9zu2U5mZvUrsusJAEkvAL4NnBcRjw3a/FPgpRHxuKQTgH7g0DL7mAXMAhg/fnzNNXi2k5lZ/Qo9o5C0O1lIXBcR3xm8PSIei4jH88dLgd0l7Vum3TURMTUipo4dO3aXbdWW7PDZhJnZ8BQ560nAl4CfR8T8Idq8OG+HpFfn9TxUy+dUW7LDZxNmZsNTZNfTMcA7gHWS1uSv/R0wHiAirgJOA94vaTvwO+CMiIgCazIzsxoVFhQRcRNVrnOLiH8E/rHez3C3k5lZ8Tr6ymx3O5mZFa9jgyJlpVgzMxu+jg2K86qsFHvMxL2bU4iZWZfr2KCo5rpzjm51CWZmXaFrg8LMzBqjI4Oi2mynEV5T3MysYToyKKrNdlp/2YlNqsTMrPt1ZFCYmVnzdFxQ/PzBwesK7soX2ZmZNVbHBcX2Zyuv8OGL7MzMGqvjgsLMzJqrq4LCF9mZmTVeVwWFL7IzM2u8rgoKMzNrvK4JCnc7mZkVo2uCwt1OZmbF6JqgMDOzYnRFUIzs8+JOZmZF6Yqg+MRpR7S6BDOzrtUVQTHjyHGtLsHMrGt1RVCYmVlxOj4oPC3WzKxYHR0U+48Z6WmxZmYFG9HqAurRJ/H2ow7i0hmHt7oUM7Ou13FBcfi4vVh12QmtLsPMrGd0dNeTmZkVz0FhZmYVOSjMzKwiB4WZmVXkoDAzs4ocFGZmVpGDwszMKnJQmJlZRQ4KMzOryEFhZmYVKSJaXUNNJG0FftHqOtrEvsBvWl1Em/Cx2MnHYicfi53+ICLG1PPGjlvrCfhFRExtdRHtQNIqH4uMj8VOPhY7+VjsJGlVve9115OZmVXkoDAzs4o6MSiuaXUBbcTHYicfi518LHbysdip7mPRcYPZZmbWXJ14RmFmZk3UtkEh6XhJv5C0XtIFZbaPkrQo375S0oQWlNkUCcfifEk/k7RW0n9Iemkr6myGaseipN2pkkJS1854STkWkt6W/9u4U9LXml1jsyT8PzJe0o2SVuf/n3TlbTIlfVnSZkl3DLFdkj6XH6e1kl6VtOOIaLsfoA+4BzgEGAncDrxiUJu/Bq7KH58BLGp13S08Fq8Dnp8/fn8vH4u83Rjgx8CtwNRW193CfxeHAquBF+XP92t13S08FtcA788fvwLY0Oq6CzoWrwVeBdwxxPYTgO8DAqYBK1P2265nFK8G1kfEvRHxNPAN4JRBbU4Brs0ffws4TpKaWGOzVD0WEXFjRDyZP70VOLDJNTZLyr8LgL8HrgC2NbO4Jks5FucA/xQRvwWIiM1NrrFZUo5FAHvmj/cC/qeJ9TVNRPwYeLhCk1OAr0TmVuCFkg6ott92DYpxwAMlz3+dv1a2TURsBx4F9mlKdc2VcixKnU32F0M3qnos8lPpgyJiSTMLa4GUfxcvA14m6WZJt0o6vmnVNVfKsfg4MFPSr4GlwP9pTmltp9bvE6Azr8y2IUiaCUwF/qzVtbSCpN2A+cBZLS6lXYwg6346luws88eSDo+IR1pZVIu8HfiXiPi0pKOBr0o6LCKebXVhnaBdzyg2AgeVPD8wf61sG0kjyE4nH2pKdc2VciyQ9AZgDnByRDzVpNqardqxGAMcBqyQtIGsD3Zxlw5op/y7+DWwOCKeiYj7gF+SBUe3STkWZwPXA0TELcBosnWgek3S98lg7RoU/w0cKulgSSPJBqsXD2qzGHhX/vg04IeRj9Z0marHQtKRwNVkIdGt/dBQ5VhExKMRsW9ETIiICWTjNSdHRN1r3LSxlP9H+snOJpC0L1lX1L1NrLFZUo7F/cBxAJL+kCwotjS1yvawGHhnPvtpGvBoRDxY7U1t2fUUEdslnQv8gGxGw5cj4k5JlwCrImIx8CWy08f1ZIM3Z7Su4uIkHotPAi8AvpmP598fESe3rOiCJB6LnpB4LH4AvFHSz4AdwEciouvOuhOPxYeBBZI+RDawfVY3/mEp6etkfxzsm4/HXATsDhARV5GNz5wArAeeBN6dtN8uPFZmZtZA7dr1ZGZmbcJBYWZmFTkozMysIgeFmZlV5KAwM7OKHBRmZlaRg8KsDEkHSbpP0t758xflzydIOkDSDfnr8yStKfn5paQdkl4g6SxJW/LX78rn8A/s/1xJf9Wq38+sFr6OwmwIkj4KTIqIWZKuJlua+jJJnwRuioh/LfOe64D7ImKupLPIljk/V9I+wC+AIyPiAUnPB26OiCOb+CuZ1cVnFGZD+wwwTdJ5wGuAT+WvnwosG9w4X5RxEtlKpbvIr4heDxyQP38S2CDp1UUUbtZIbbmEh1k7iIhnJH2ELBTemD8/GPjt4IUXld1h8XLg2HzZewZtH0+2vtDakpdXAX8K/KSgX8GsIXxGYVbZm4AHyValheyMYJfF5CT1AQuBj0XE+kHvP13SWrKzic9HROnNlDYDLymkarMGclCYDUHSFGA62XLlH8rvBPY7sjODUnOBByPin8vsZlFETAb+BLhc0otLto3O92fW1hwUZmXkt9X9AnBeRNxPtkLvp8ju6TChpN00shslzaq0v3yp868CHyx5+WXAHY2s26wIDgqz8s4hW659ef7888Afkt1B8B5Jk/LXLwaeD9w4aJrsxDL7vAJ4t6Qx+fNjgOVl2pm1FU+PNauRpLcAfxQRc4exjyOB8yPiHY2rzKwYnvVkVqOI+G5+XcRw7At8rBH1mBXNZxRmZlaRxyjMzKwiB4WZmVXkoDAzs4ocFGZmVpGDwszMKvr/6+FRmzXpsIUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "from pycalphad import calculate\n",
    "\n",
    "mod = ViscosityModel(dbf, ['CU', 'ZR'], 'LIQUID')\n",
    "\n",
    "temp = 2100\n",
    "# NOTICE: we need to tell pycalphad about our model for this phase\n",
    "models = {'LIQUID': mod}\n",
    "res = calculate(dbf, ['CU', 'ZR'], 'LIQUID', P=101325, T=temp, model=models, output='viscosity') \n",
    "\n",
    "fig = plt.figure(figsize=(6,6))\n",
    "ax = fig.gca()\n",
    "ax.scatter(res.X.sel(component='ZR'), 1000 * res.viscosity.values)\n",
    "ax.set_xlabel('X(ZR)')\n",
    "ax.set_ylabel('Viscosity (mPa-s)')\n",
    "ax.set_xlim((0,1))\n",
    "ax.set_title('Viscosity at {}K'.format(temp));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "We repeat the calculation for Al-Cu."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "execution": {
     "iopub.execute_input": "2022-02-19T19:24:49.701334Z",
     "iopub.status.busy": "2022-02-19T19:24:49.701334Z",
     "iopub.status.idle": "2022-02-19T19:24:49.911692Z",
     "shell.execute_reply": "2022-02-19T19:24:49.911692Z"
    },
    "jupyter": {
     "outputs_hidden": false
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAGDCAYAAAAiU8cRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAmLklEQVR4nO3de7xcVX338c+XE5IoHlQgSRGMEYJS5WpTgeJjUYzFxEu8giZarUJrrxQfNTxJn4glJbU2xbbeiPoSBTXeekoNxkYl8ogBG3qQAEYhENEUk3BNuAUIv+ePvY85hHPOrJkze8/M3t/365UXZ2bW7PmdTbJ+s39r7bUUEZiZWf3s0+kAzMysM5wAzMxqygnAzKymnADMzGrKCcDMrKacAMzMasoJwEoj6UZJp3Tgc++XdFjZn2vW7ZwArC0krZb04RGef52kX0uaEBEvjIi1ZccWEU+LiFvzeD4v6fwiPkfSOyX9sEGbt0j6kaQHJa3d67WDJF0l6S5J90paJ+nkvdr8dX4+d0j6nKRJw16bIemK/NgbJb1i2GsfknTJsMeH5G3+WZLG/ctbT3ICsHa5GFgwQmfyduDSiHisAzF1o7uBC4FlI7x2P/BHwBTgmcDfA/8haQKApD8AFgKnAs8BDgPOG/b+LwODwIHAIuDrkqbs/SGSngNcCVwWEX8Zvhu0viLCf/xn3H+ApwD3AS8d9twzgYeBY/PHm4FX5D+/GFgP7AC2AsuHve8lwI+Ae4FfAu/Mn3868AVgO/ALYDGwT/7aTOAHeQx3AiuHHS/y188CHgUeIets/wN4P/CNvX6XfwY+NsrvuRDYBOwEbgJenz//2/nvujs/9r0Nztd7gLVjvL4P8Jo89qn5c18C/m5Ym1OBX+c/Pw/YBfQPe/3/AX+S//wh4BLg8PzcfbjTf2f8p/N/JiTmCbMxRcRDkr4KvIPs2yXAW4CNEfGTEd7yMbJO9ouSngYcBb/5dvptss7668D+wLPz9/wLWRI4jOxb7n8CdwCfBf42f/wyYCIwa4QYL5L0e8CvImJx/nkHAx+S9IyIuDf/tn0G8KpRftVNwP8Cfg28GbhE0syI+KmkPwHeExEvaXzGRifpeuBIYF/gMxGxLX/phcC/D2v6E2CapAPz126NiJ17vf7CYY8PI/t/8/GI+LvxxGjV4BKQtdPFwJskTc4fvyN/biSPAjMlHRQR90fE1fnzbwO+GxFfjohHI+KuiLhOUh9Zx3xuROyMiM3AP5KVmIaO9xzgWRHxcESMWYsfEhF3kHWKb86fOg24MyKuHaX91yLifyLi8YhYCdxMdjXTNhFxDFniexsw/Pd4GtkVzpChn/tHeG3o9f5hj48C9gNWtjNe611OANY2ead7JzBP0uFkHeOXRmn+brKyxUZJ/yXp1fnzzyb7lr23g8i+Ef9i2HO/AA7Jf/4AIODH+WyjP2oi9IuBBfnPC4AvjtZQ0jskXZcP0t5L1qke1MRnJcmT2JeBhZKOzZ++nywxDBn6eecIrw29PvyK4DLgc8D38ystqzknAGu3L5B9818AfCcito7UKCJujoi3AlPJBju/Lmk/spr/4SO85U72fMsfMh3Ykh/v1xFxZkQ8C/hj4BOSZo700SM8NwAcI+ko4NXApSPFnHeaK4A/Bw6MiGcAN5AlntGOPV77kpVuAG4Ejh322rHA1oi4K3/tMEn9e71+4/CDRcQ5wLfIksAhWK05AVi7fQF4BXAmo5d/kLRA0pSIeJxssBfgcbLO9xX5dMkJkg6UdFxE7Aa+CiyV1J93xueQDWwi6c2SDs2Pcw9ZZ/z4CB+9lT0dKpB92yYbb/gS8OOIuH2UsPfLj7s9/8x3kY9dDDv2oZImjvF79+UlsgnAPpImS9o3f+1ESS+RNFHSUyR9EJgGXJO//QvAuyW9QNIzyAbBP5//Dj8HrgOW5Md8PXAM8I0Rwvhz4Arge5KmjRarVZ8TgLVVXpv/EVlnedkYTU8DbpR0P9mA8BkR8VDe+c4B3kc2ZfI69nzr/QvgAeBWstr4l8hKGgC/C1yTH+8y4K8in/u/l88CL8hLOAPDnr8YOJoxyj8RcRPZuMM6ss7+aOCqYU2+T/aN+9eS7hzlMG8HHgI+STaY/BDZVQXAJODjwF1kVzZzgLkR8T/5568GPkLWed9OVgJbMuzYZ5ANft9DNs30TRGxfYTfI8gG2X8MfFdS20tY1huU/V0wqzdJ04GNwG9FxI5Ox2NWBl8BWO1J2oesnPQVd/5WJ74PwGotH3jeSlZOOa3D4ZiVyiUgM7OacgnIzKymnADMzGqqq8YADjrooJgxY0anwzAz6xnXXnvtnRHxpFVfU3RVApgxYwbr16/vdBhmZj1D0i8atxqZS0BmZjXlBGBmVlNOAGZmNeUEYGZWU04AZmY15QRgZlZTTgBmZjXlBGBmVlNOAGZmNdVVdwKbWXU9d+GqMTdN3rxsbmmxWMYJwMzaZsbCVYW8d1r/RK5ZNLvlY9vInADMrGUzz13FYyVsKbJ15yNPSBC+WmgPJwAza8rs5Wu5edsDHY3ByaA9nADMLEmjGn6nDCWDI6bux5pzTulsMD3GCcDMxlRWmWe8bt72ADMWrmKC4JYLfFWQwgnAzEY0ngHdTnosstgF3Oby0JicAMzsCbq11NOsIEsE+0/q4/rzTut0OF3JCcDMAJi/Yh1Xbbq7rcdsdYD2yEWX8/Du9qShHbt2M2PhKo8RjEAR3ZPrZ82aFd4S0qx87Sr3nHz4AVx65kltOdZw7SxHVW3WkKRrI2JWS+91AjCrr3Z86y+7xNKOZDC5T2xcOqcN0XSeE4CZNW28HemFpx/HvOMPaVM0zWvHWEWnf4d2cAIws2QDg1s4e+V1Lb+/20oo470xrdeXmRhPAvBqoGY1Mn/FupY7/wUnTu+6zh9gzTmnsHnZXE4+/ICW3r/3MhN14isAs5potZPrtdkz45lBVNQgdpFcAjKzMbXa+XfjN/5U47mDuZd+b5eAzGxUrXT+3VruacYtF8zlwtOPa+m9dSkJ+QrArKJaHezt9Y5/JK2WhXqhJOQrADN7gsUDG5ru/KvwrX80G5fOael3u2rT3ZW+GnACMKuYE5au4ZKrb2/qPZuXzeX8eUcXFFH32LxsLpP71PT7qpoEnADMKuTIRZezdecjye0nqJoln7G0ejUwY+EqFg9sKCCiznECMKuImeeuaqrOfcTU/Wq9bv7mZXNp9lrgkqtvZ+a51bkacAIwq4BmpzwuOHF6T83tL8pty+YyrX9iU+8Z2m+gCjwLyKzHNTvDpW4ln1StdOrdcC49C8isptz5t0+r4wLzV6wrIJpyOAGY9ahjlqx2599mm5fNZf9JfU2956pNd3PMktUFRVQsJwCzHnTMktXs2LU7ub07/3TXn3da0+draNexXuMEYNZjZi9f686/BK2WhHqJE4BZDxkY3NLU2vfu/Men6knACcCshzSzvIM7//Zo5e7hGQtXMTC4paCI2scJwKxHNPPN0p1/e21cOqfpDWfOXnkdJyxdU1BE7eEEYNYD3Pl33qVnntT0ue323cacAMy6nDv/7lKlcQEnALMu5s6/O1UlCRSaACRtlrRB0nWSvMaDWROaucPUnX/5WrlprNsGh8u4AnhZRBzX6loVZnU0MLiFqzbdndS21W0PbfyuP++0ps//2Suv48hFlxcTUJNcAjLrQqnTPaf1T2Te8YcUG4yNad7xhzR9Bfbw7uiKklDRCSCA/5R0raSzCv4ss0popmO4ZtHsAiOxZrQ6LtDJklDRCeAlEfEi4FXAn0l66d4NJJ0lab2k9du3by84HLPu5kHf3tbK/5OzV17XsU1mCk0AEbEl/+824N+AF4/Q5qKImBURs6ZMmVJkOGZdrZmbhtz5d6/Ny+ay4MTpTb2nU5vMFJYAJO0nqX/oZ+CVwA1FfZ5ZLxsY3JK8l687/+53/ryjWy4JzV6+tv0BjaLIK4BpwA8l/QT4MbAqInpz0WyzgqUO+rrz7y2t/P+6edsDpV0NFJYAIuLWiDg2//PCiFha1GeZ9bLUf+zNlhWsO7RyvwBkfy+Kni7qaaBmHZTa+U/rn8j5844uOBoryvXnndZSAh+aLrp4YEMBUXlTeLOOmb18bfLa/i79VMd4yjsnH34Al5550hOe86bwZj2mmY1d3PlXy+Zlc5nWP7Gl91616e62DhQ7AZh1gAd96+2aRbPHNaYzNFA8Y+Eq+p4+reUDOQGYlSy1BNDsBiTWW4amija319iT9T2lv+UbqJwAzEqUeuk+QTyp1mvVdNuyuR1b0M8JwKwkzdT9b7nApZ86GVpQ7oip+5X6uU4AZiVx3d8aWXPOKaX+/3cCMCtB6g09rvsbZF8CykgETgBmBZu/Yh0P7258v41w3d+eaCgRtDpttJEJhRzVzH4jdWev21z6sVEM3/ehnesE+U5gswKl/mN13d9a5TuBzbpQ6votrvtbpzgBmBXkkqtvb9jG8/2tk5wAzAqQWvrxfH/rJCcAszZL3d/VdX/rNCcAszYaGNzCYwnzKsq+49NsJE4AZm2UerfvmnNOKTQOsxROAGZtkrrQm0s/1i2cAMzaIHWhN0/5tG7iBGDWBqmlH0/5tG7iBGA2Tp71Y73KCcBsHBYPbEia9dOpDT/MxuIEYDYOKXf7Tu4T844/pIRozJrjBGDWohOWrklqt3HpnIIjMWuNE4BZCwYGt7B15yMN2y04cXoJ0Zi1xgnArAUps34m94nz5x1dfDBmLXICMGuSSz9WFU4AZk1KKf141o/1AicAsyakLPO8/6Q+z/qxnuAEYJYotfRz/XmnFRyJWXs4AZgl8qwfqxonALMEz03c4cuzfqyXOAGYNTAwuIWE1R681o/1HCcAswZS5vzvP6mv+EDM2swJwGwMqZu7e+DXepETgNkoBga3JLVz6cd6lROA2ShSSj/e3N16mROA2QiOWbI6qZ03d7de5gRgNoIdu3Y3bOPlHqzXOQGY7cXLPVhdOAGYDTN7+dqkdp71Y1XgBGA2zM3bHmjYZlr/xBIiMSueE4BZLnXg95pFswuOxKwcTgBmZHP+UwZ+PeffqsQJwAzP+bd6cgKw2ksd+PWcf6saJwCrvZSBX6/zb1XkBGC1ljLwO0Fe59+qyQnAamvxwIakgd9bLvDAr1WTE4DV1iVX396wjUs/VmVOAFZLqQO/Lv1YlTkBWC2lDPx6sTerOicAq50jF13esM20/ole7M0qzwnAamXxwAYe3t14i3cv92B1UHgCkNQnaVDSt4r+LLNGPPBrtkcZVwB/Bfy0hM8xG1PqYm8e+LW6KDQBSDoUmAt8psjPMWvEi72ZPVnRVwAXAh8AHh+tgaSzJK2XtH779u0Fh2N15cXezJ6ssAQg6dXAtoi4dqx2EXFRRMyKiFlTpkwpKhyrscUDG5LaebE3q5sirwBOBl4raTPwFeDlki4p8PPMRuSBX7ORFZYAIuLciDg0ImYAZwDfj4gFRX2e2Uiem7DBuxd7s7ryfQBWWQODW2g849+LvVl9TSjjQyJiLbC2jM8yG5Iy8Hvy4QcUH4hZl/IVgFVS6mJvl555UrGBmHUxJwCrJC/2ZtaYE4BVTsq3/wnCi71Z7TkBWOWkfPv3wK9ZwiCwpH2AY4FnAQ8BN0TEtqIDM2vFzHMbT/v0Hb9mmVETgKTDgQ8CrwBuBrYDk4HnSXoQ+DRwcUSMusyDWZnmr1jHYwnzPn3Hr1lmrCuA84FPAn8cEU/4ZyVpKvA24O3AxcWFZ5buqk13N2zjO37N9hg1AUTEW8d4bRvZQm9mXSFlly/wHb9mwzUcBJb0Zkn9+c9/I+mbkl5UfGhm6VJ2+fJSz2ZPlDIL6G8iYqeklwCnAp8lKw2ZdYWU9X6m9U8sIRKz3pKSAIZ20ZgLXBQRqwD/a7KucMLSNUnr/XiPX7MnS0kAWyR9GjgduFzSpMT3mRVu685HGrbxt3+zkaV05G8BvgP8QUTcCxwAvL/IoMxSnLB0TVI7f/s3G1nDG8Ei4kHgm5Bt3xgRFwF3FB2YWSMp3/498Gs2umZLOX9SSBRmTTpmyeqGbVz6MRtbswlAhURh1oSBwS3s2LW7YTuXfszG1mwCeE0hUZg1wRu9mLVH0o5gkuYCLwQmS9lFQER8uMC4zEY0f8W6pHbe6MWssZQ7gT9FNgX0L8hKQG8GnlNwXGYjSlnvxwO/ZmlSSkC/FxHvAO6JiPOAk4DnFRuW2ZOlLPW8/6S+EiIxq4aUBPBw/t8HJT0LeBQ4uLiQzJ5s8cCGpKWerz/vtOKDMauIlDGA/5D0DOAfgP8GAlhRZFBme7vk6tsbtvFSz2bNGTMBSDoO+ClwcER8Q9K3gMkRcV8ZwZlB2pz//Sf1ealnsyaNWgKS9H+BrwJvBFZJOjMidrnzt7KlzPl36ceseWNdAZwOHBcRD0o6EFiNSz9WspRv/7470aw1Yw0C78rXASIi7mrQ1qztUu/4vc3TPs1aMtYVwGGSLst/FnD4sMdExGsLjcxqL+WOX6/3Y9a6sRLA6/Z6/NEiAzEbbvHAhqR2Xu/HrHVjbQr/gzIDMRsuZdrnhacfV3wgZhWWshTEqyUNSrpb0g5JOyXtKCM4q6eUgV+AeccfUnAkZtWWciPYhcAbgA0RkbL9qtm4pAz8er0fs/FLmdnzS+AGd/5WhhkLG6/346Wezdoj5QrgA2Sbwf8A2DX0ZEQsLywqq6XZy9cmtfNSz2btkZIAlgL3A5MBz7mzwty87YGGbTzwa9Y+KQngWRFxVOGRWK2lLPV8xNT9PPBr1kYpYwCXS3pl4ZFYbaUu9bzmnFMKj8WsTlISwHuB1ZIe8jRQK0LKnP8jpu5XQiRm9dKwBBQR/WUEYvWUesevv/2btd9Yy0HPGOuNyhza9oisVnzHr1nnjHUF8A+S9gH+HbgW2E42E2gm8DLgVGAJ8Kuig7RqShn4Fb7j16woY60F9GZJLwDmA39Etg/wg2Q7hF0OLI2Ih0d7v9lYBga3JA38eqlns+KMOQYQETcBi0qKxWokZaln3/FrVixv8mKlm79iXVI73/FrViwnACvdVZvubthmwYnTS4jErN6cAKxUz01Y7G1ynzh/3tElRGNWbyn7AXxT0tx8RpBZyxYPbCBlSdmNS+cUHouZpV0BfAJ4G3CzpGWSnl9wTFZRKXP+PfBrVp6GCSAivhsR84EXAZuB70r6kaR3Sdq36ACtGjzwa9Z9kso6kg4E3gm8BxgEPkaWENYUFplVigd+zbpPw7WAJP0b8Hzgi8BrIuKO/KWVktYXGZxVQ8ouXx74NStfyn4AKyLi8uFPSJoUEbsiYlZBcVlFpO7y5YFfs/KllIDOH+G5tIKu1V7KLl8e+DXrjFGvACT9FnAI8BRJx5OtywWwP/DUEmKzHnfMktVJ7Tzwa9YZY5WA/oBs4PdQYPgG8DuB/1NgTFYBA4Nb2LFrd8N2XurZrHPGWg30YuBiSW+MiG80e2BJk4ErgUn553w9Ipa0HKn1lJTF3qb1T/RSz2YdNFYJaEFEXALMkHTO3q9HxPIR3jbcLuDlEXF/fr/ADyV9OyKuHl/I1u1S5/xfs2h2wZGY2VjGKgENbcL6tFYOHBEB3J8/3Df/k7ISgPU4z/k36w1jlYA+nf/3vFYPLqmPbDexmcDHI+KaVo9lvSFlzj/gOf9mXSBlMbiPSNpf0r6Svidpu6QFKQePiN0RcRzZQPKLJR01wvHPkrRe0vrt27c3/QtY9xgY3JLUbrN3+TLrCin3AbwyInYAryZbC2gm8P5mPiQi7gWuAE4b4bWLImJWRMyaMmVKM4e1LpMy8Lv/pL7iAzGzJCkJYKhMNBf4WkTcl3JgSVMkPSP/+SnAbGBjK0Fa90ud83/9eU/6DmBmHZKyFMS3JG0EHgLeK2kKkLIZ/MFk00j7yBLNVyPiW62Hat0qdc6/Sz9m3aVhAoiIhZI+AtwXEbslPQC8LuF91wPHtyFG63IppZ8jpu7XsI2ZlStlNdB9gQXASyUB/AD4VMFxWY84YWnaiuBrzjml2EDMrGkpJaBPks3h/0T++O35c+8pKijrHVt3PtKwjRd7M+tOKQngdyPi2GGPvy/pJ0UFZL0jZYN38GJvZt0qZRbQbkmHDz2QdBjQeMTPKi11g3cP/Jp1r5QrgPcDV0i6lWxJ6OcA7yo0Kut6KRu8e+DXrLulzAL6nqQjyLaFBPhZROwqNizrZjPPTSv9eODXrLulLAXxZ8BTIuL6fGrnUyX9afGhWTdaPLCBxxJqPy79mHW/lDGAM/OlHACIiHuAMwuLyLqaSz9m1ZGSAPqU3wAAv1nhc2JxIVm3Sp3149KPWW9IGQReDayU9On88R/nz1mNDAxuSZr14y0ezXpHSgL4IHAW8N788RrgM4VFZF0pZbmHyX3yFo9mPSRlFtDjZEs/fErSAcChEeH7AGpk9vK1Se02Lp1TbCBm1lYps4DW5hvCHEC2u9cKSf9UfGjWLW7e9kDDNt7i0az3pAwCPz3fEOYNwBci4gTg1GLDsm6RssXjBHmLR7NelLQhjKSDgbcAXs+/RlJX+rzlAs/5N+tFKQngw8B3gFsi4r/ytYBuLjYs6wYpK316zr9Z70oZBP4a8LVhj28F3lhkUNZ5KaUf8Jx/s142agKQ9IGI+Iikf4EnTwGPiL8sNDLrmNRZP17uway3jXUFcFP+3/VlBGLdI2XWjzd5Met9YyWAV0m6JyIuLi0a67jU0o83eTHrfWMNAv8c+KikzZI+IskbvFdcaufv0o9ZNYyaACLiYxFxEvD7wF3A5yRtlLRE0vNKi9C6iks/ZtXRcBpoRPwiIv4+Io4H3grMA35adGBWLu/va1Y/KUtBTJD0GkmXAt8GfkZ2V7BVxOzla72/r1kNjTUNdDbZN/45wI+BrwBnRUTjKSLWMwYGt3jWj1lNjTUL6FzgS8D78l3ArIJSlnkGl37MqmjUBBARLy8zECvfkYsuT2rn0o9ZNaWsBWQV9fDuxpV/l37MqssJoKZ8w5eZOQHUkEs/ZgZOALWzeGBDUulnWv/EEqIxs05yAqiZS66+PandNYtmFxyJmXWaE0CNeK0fMxvOCaAm5q9Yl9TOs37M6sMJoCau2nR3wzaT++RZP2Y14gRQA6mln41L5xQciZl1EyeAijth6ZqkdheeflyxgZhZ13ECqLCBwS1s3flIw3aT+8S84w8pISIz6yZOABWWutCbSz9m9eQEUFGe8mlmjTgBVNAxS1YntZvcp4IjMbNu5gRQMYsHNrBj1+6kti79mNWbE0DFpC714NKPmTkBVEjqxu7u/M0MnAAqY/6KdUkbu09w2d/Mck4AFTAwuCVpqQeAWy7wt38zyzgBVEDqfH+XfsxsOCeAHpc6399LPZjZ3pwAetjs5WuT2k3rn+ilHszsSZwAetTigQ3cvO2BpLbe3cvMRuIE0KM839/MxssJoAel1v29u5eZjcUJoMccuejy5Lbe3cvMxuIE0EMWD2zg4d0pt3u59GNmjTkB9BDX/c2snQpLAJKeLekKSTdJulHSXxX1WXXgur+ZtduEAo/9GPC+iPhvSf3AtZLWRMRNBX5mJaV2/hPkur+ZpSvsCiAi7oiI/85/3gn8FPDdSE1KHfQVXufHzJpTyhiApBnA8cA1I7x2lqT1ktZv3769jHB6RjODvre57m9mTSo8AUh6GvAN4OyI2LH36xFxUUTMiohZU6ZMKTqcnuJBXzMrUqEJQNK+ZJ3/pRHxzSI/q2pS6/5HTN2v4EjMrKqKnAUk4LPATyNieVGfU0WpnT/AmnNOKS4QM6u0Iq8ATgbeDrxc0nX5H+9C3kDqto7g0o+ZjU9h00Aj4odkk1Ms0eKBDUnbOoI7fzMbP98J3EU86GtmZXIC6BKpdf8FJ04vOBIzqwsngC6Q2vlP7hPnzzu64GjMrC6cADqsmRk/G5d6DN3M2scJoIM848fMOskJoEPmr1jnGT9m1lFOAB2weGADV226O6nthacfV2wwZlZbTgAlGxjckjzd84ip+zHveC+gambFcAIo2dkrr0tqN7lPXubBzArlBFCiZjZ28YwfMyuaE0BJmpnu6Y1dzKwMTgAlaKbz94wfMyuLE0DBZp7rzt/MupMTQIGOWbKaxxIn+3u6p5mVzQmgICcsXcOOXbuT2nq6p5l1QmH7AdTZzHNXJX/z93RPM+sUXwG02ZGLLk/u/MHTPc2sc5wA2mjxwAYe3p3e+3vQ18w6yQmgTRYPbEhe4gHc+ZtZ5zkBtEEz6/uAO38z6w5OAG2Qur4PuPM3s+7hBDBOvsvXzHqVE8A4uPM3s17mBNAid/5m1uucAFrgzt/MqsAJoEnu/M2sKpwAmuDO38yqxAkgkTt/M6saJ4AE7vzNrIqcABpw529mVeUEMIZmOv8FJ04vMBIzs/ZzAhjBwOCWpjr/kw8/gPPnHV1gRGZm7ecNYfZyzJLVyTt5QfbN352/mfUiJ4Bhjlx0eVPr+V94+nHeytHMepYTQK7Zzt8DvmbW65wAaG6wF9z5m1k11H4Q2J2/mdVVrROAO38zq7PaJgB3/mZWd7VMAO78zcxqmADc+ZuZZWqVANz5m5ntUYsE0OzSDuDO38yqr/L3Acxevpabtz3Q1Hvc+ZtZHVQ6ATS7rg+48zez+qhsAmi25DNBcMsF7vzNrD4qOQbQbOc/rX+iO38zq53KXQE02/l7RU8zq6vKJICBwS2cvfK6pt7jer+Z1VklEsD8Feu4atPdTb3Hnb+Z1V3PJ4CZ567isfRl/AF3/mZm0OMJoNl6P7jzNzMb0rOzgJrt/IU7fzOz4XruCqCVev+0/olcs2h2QRGZmfWmnkoAze7bC7DgxOmcP+/ogiIyM+tdhSUASZ8DXg1si4ijxns81/vNzNqryDGAzwOnteNA7vzNzNqvsAQQEVcCzRXrR+DBXjOzYnR8DEDSWcBZANOnT//N84sHNnDJ1bc3dawjpu7HmnNOaWd4ZmaV1fEEEBEXARcBzJo1K6C1m7s82Gtm1pyOJ4C9ud5vZlaOrroRbMOW+5p+jzt/M7PWFJYAJH0ZWAc8X9KvJL27ncefIHf+ZmbjUVgJKCLeWtSxXe83Mxu/rhsDaMTf+s3M2qOrxgAacedvZtY+PZEAJvfJnb+ZWZt1fQI4+fAD2Lh0TqfDMDOrnK4eA/C3fjOz4nTlFYDX8zEzK15XXQEcfcjTWe+O38ysFF15BWBmZsVzAjAzqyknADOzmnICMDOrKScAM7OacgIwM6spJwAzs5pyAjAzqyknADOzmnICMDOrKUVEp2P4DUk7gZ91Oo4ucRBwZ6eD6AI+D3v4XOzhc7HH8yOiv5U3dtVaQMDPImJWp4PoBpLW+1z4PAznc7GHz8Uekta3+l6XgMzMasoJwMysprotAVzU6QC6iM9FxudhD5+LPXwu9mj5XHTVILCZmZWn264AzMysJKUnAEmnSfqZpFskLRzh9UmSVuavXyNpRtkxliXhXJwj6SZJ10v6nqTndCLOMjQ6F8PavVFSSKrsDJCUcyHpLfnfjRslfansGMuS8G9kuqQrJA3m/07mdCLOMkj6nKRtkm4Y5XVJ+uf8XF0v6UUNDxoRpf0B+oBNwGHAROAnwAv2avOnwKfyn88AVpYZY5edi5cBT81/fm+dz0Xerh+4ErgamNXpuDv49+IIYBB4Zv54aqfj7uC5uAh4b/7zC4DNnY67wPPxUuBFwA2jvD4H+DbZtuonAtc0OmbZVwAvBm6JiFsj4hHgK8Dr9mrzOuDi/OevA6dKUokxlqXhuYiIKyLiwfzh1cChJcdYlpS/FwB/C/w98HCZwZUs5VycCXw8Iu4BiIhtJcdYlpRzEcD++c9PB/6nxPhKFRFXAneP0eR1wBciczXwDEkHj3XMshPAIcAvhz3+Vf7ciG0i4jHgPuDAUqIrV8q5GO7dZNm9ihqei/xy9tkRsarMwDog5e/F84DnSbpK0tWSTistunKlnIsPAQsk/Qq4HPiLckLrSs32KV13J7CNQNICYBbw+52OpRMk7QMsB97Z4VC6xQSyMtApZFeFV0o6OiLu7WRQHfJW4PMR8Y+STgK+KOmoiHi804H1grKvALYAzx72+ND8uRHbSJpAdll3VynRlSvlXCDpFcAi4LURsauk2MrW6Fz0A0cBayVtJqtvXlbRgeCUvxe/Ai6LiEcj4jbg52QJoWpSzsW7ga8CRMQ6YDLZOkF1lNSnDFd2Avgv4AhJz5U0kWyQ97K92lwG/GH+85uA70c+wlExDc+FpOOBT5N1/lWt80KDcxER90XEQRExIyJmkI2HvDYiWl4DpYul/BsZIPv2j6SDyEpCt5YYY1lSzsXtwKkAkn6bLAFsLzXK7nEZ8I58NtCJwH0RccdYbyi1BBQRj0n6c+A7ZCP8n4uIGyV9GFgfEZcBnyW7jLuFbMDjjDJjLEviufgH4GnA1/Jx8Nsj4rUdC7ogieeiFhLPxXeAV0q6CdgNvD8iKneVnHgu3geskPTXZAPC76zoF0YkfZks8R+Uj3ksAfYFiIhPkY2BzAFuAR4E3tXwmBU9V2Zm1oDvBDYzqyknADOzmnICMDOrKScAM7OacgIwM6spJwAzs5pyArBakfRsSbdJOiB//Mz88QxJB0v61rC2L5Z0Zb4c8aCkz0h6qqQPSfrfex13s6SDJE3M3+NlVqzrOQFYrUTEL4FPAsvyp5YBF0XEZuAcYAWApGnA14APRsTzI+J4YDXZshRjHf8R4HvA6YX8AmZt5ARgdfRPwImSzgZeAnw0f/6NZJ08wJ8BF+frywAQEV+PiK0Jxx8A5rctWrOC+DLVaiciHpX0frLO/pX54+cC9wxbcO8o9uxL0awbgN9tQ6hmhfIVgNXVq4A7yDp6gINJX0RstPVTAiAidgOPSBqzXGTWaU4AVjuSjgNmky0r/df5rkkPka0kOeRG4HdGOcRdwDP3eq4fuHfY40lUe+cyqwAnAKuVfHvRTwJnR8TtZCuufpRsTf0Zw5r+K/CHkk4Y9t435IPDVwKvHfqGL+kNwE/yb/5IOhC4MyIeLeFXMmuZxwCsbs4kW1Z7Tf74E2TL5s4CNkmaGRG3RMRWSWcAH5U0FXicrONfnb/2r8APJQWwDXjPsM94GVD1rSutArwctFlO0uuB34mIxeM8zjeBhRHx8/ZEZlYMXwGY5SLi3/LyTcvynasG3PlbL/AVgJlZTXkQ2MysppwAzMxqygnAzKymnADMzGrKCcDMrKb+P/NiE98ar1L8AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "from pycalphad import calculate\n",
    "\n",
    "temp = 1300\n",
    "models = {'LIQUID': ViscosityModel}  # we can also use Model class\n",
    "res = calculate(dbf, ['CU', 'AL'], 'LIQUID', P=101325, T=temp, model=models, output='viscosity')\n",
    "\n",
    "fig = plt.figure(figsize=(6,6))\n",
    "ax = fig.gca()\n",
    "ax.scatter(res.X.sel(component='CU'), 1000 * res.viscosity.values)\n",
    "ax.set_xlabel('X(CU)')\n",
    "ax.set_ylabel('Viscosity (mPa-s)')\n",
    "ax.set_xlim((0,1))\n",
    "ax.set_title('Viscosity at {}K'.format(temp));"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}