{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from pystencils.session import *\n", "sp.init_printing()\n", "frac = sp.Rational" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Tutorial 05: Phase-field simulation of dentritic solidification\n", "\n", "This is the second tutorial on phase field methods with pystencils. Make sure to read the previous tutorial first. \n", "\n", "In this tutorial we again implement a model described in **Programming Phase-Field Modelling** by S. Bulent Biner.\n", "This time we implement the model from chapter 4.7 that describes dentritic growth. So get ready for some beautiful snowflake pictures.\n", "\n", "We start again by adding all required arrays fields. This time we explicitly store the change of the phase variable φ in time, since the dynamics is calculated using an Allen-Cahn formulation where a term $\\partial_t \\phi$ occurs." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "dh = ps.create_data_handling(domain_size=(300, 300), periodicity=True, \n", " default_target='cpu')\n", "φ_field = dh.add_array('phi', latex_name='φ')\n", "φ_delta_field = dh.add_array('phidelta', latex_name='φ_D')\n", "t_field = dh.add_array('T')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This model has a lot of parameters that are created here in a symbolic fashion. " ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$$\\frac{{{φ}_{C}}^{4}}{4} - {{φ}_{C}}^{3} \\left(- \\frac{m}{3} + \\frac{1}{2}\\right) + {{φ}_{C}}^{2} \\left(- \\frac{m}{2} + \\frac{1}{4}\\right) + \\frac{ε^{2}}{2} \\left({\\partial_{0} {{φ}_{C}}}^{2} + {\\partial_{1} {{φ}_{C}}}^{2}\\right)$$" ], "text/plain": [ " 4 2 ⎛ 2 2⎞\n", "φ_C 3 ⎛ m 1⎞ 2 ⎛ m 1⎞ ε ⋅⎝D(phi_C) + D(phi_C) ⎠\n", "──── - φ_C ⋅⎜- ─ + ─⎟ + φ_C ⋅⎜- ─ + ─⎟ + ──────────────────────────\n", " 4 ⎝ 3 2⎠ ⎝ 2 4⎠ 2 " ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ε, m, δ, j, θzero, α, γ, Teq, κ, τ = sp.symbols(\"ε m δ j θ_0 α γ T_eq κ τ\")\n", "εb = sp.Symbol(\"\\\\bar{\\\\epsilon}\")\n", "\n", "φ = φ_field.center\n", "T = t_field.center\n", "\n", "def f(φ, m):\n", " return φ**4 / 4 - (frac(1, 2) - m/3) * φ**3 + (frac(1,4)-m/2)*φ**2\n", "\n", "free_energy_density = ε**2 / 2 * (ps.fd.Diff(φ,0)**2 + ps.fd.Diff(φ,1)**2 ) + f(φ, m)\n", "free_energy_density" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The free energy is again composed of a bulk and interface part." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAb4AAAEWCAYAAAAZwvJqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xl8lOW9///XJ3sC2UkCCWvYERAwsmgV64pVQa22aLVotaCn1nNOz2lr23Nqj63ftqc/62mtraK1Ym3r2ioqapEqagUk7PsWIiRhyUISCNlz/f7IYNMYJCGTuWcy7+fjMY+ZuZe5P7mzvHPf93VflznnEBERCRcRXhcgIiISSAo+EREJKwo+EREJKwo+EREJKwo+EREJKwo+EREJKwo+kdNkZs7MRvheP2lmP+rkellm9q6ZHTWzB3q2ShFpL8rrAkS8YmaFQBbQDDQCHwB3OOf29/Cm5wNlQJLTjbQiAacjPgl3Vznn+gIDgEPAQwHY5hBg68lCz8xC8h/SUK1bwo+CTwRwztUBLwDjTkwzs3fM7PY2728xs/dP9Vlmlmhmb5vZL83M2s17EpgHfMvMjpnZxWb2AzN7wcyeNrNq4BYzizCze8xsj5mVm9lzZpbW5nOmm9kHZlZpZhvM7IJPqSfbzF40s1Iz22tmd7eZ9wPfZz/lO/W6xczyurBu+7rjzWyRmR0xs21m9i0zK/It/00ze7FdbQ+Z2f+dap+K+JOCTwQwswTgi8DKbn5OOrAM+Ltz7u72R3XOuVuAPwD/65zr65x7yzdrDq3Bm+KbfzdwNTATyAaOAA/7tpEDvAb8CEgD/hN40cwyOqgnAngF2ADkABcB/2Zml7VZbDbwjG/bi4FfdWHd9nXfCwwFcoFLgJvaLPs0MMvMUnyfH0XrPv/9SXanSI9Q8Em4e8nMKoFqWv9Q/6wbn5UNLAeed879VxfXXeGce8k51+KcqwUWAN9zzhU55+qBHwDX+cLiJmCJc26Jb/mlQD7wuQ4+92wgwzl3n3OuwTlXADwGzG2zzPu+z2qmNYTO7MK67ev+AvD/nHNHnHNFwC9PLOicOwC8C1zvmzQLKHPOrenivhLpFp2Tl3B3tXPuLTOLpPXoZbmZjXPOHTyNz7oCOAY8chrrtm9QMwT4i5m1tJnWTGtjnCHA9WZ2VZt50cDbHXzuECDbF+4nRALvtXnf9ms9DsT5ArYz67avO7vdtPbzFwF30hqgN6GjPfGAjvhEAOdcs3Puz7SGy2d8k2uAhDaL9T/FxzwGvAEsMbM+XS2h3fv9wOXOuZQ2jzjnXLFv3u/bzevjnPtJB5+7H9jbbtlE51xHR4ens277ug8AA9u8H9Ru/kvARDMbD1xJ6+lRkYBS8IkA1moOkAps801eD1xrZgm++/Vu68RH3QXsAF41s/hulPQIcL+ZDfHVl+GrD1qvlV1lZpeZWaSZxZnZBWY2sIPP+RCoNrNv+xqeRJrZeDM7uxM1nM66zwHfMbNU37XIu9rObNOI6I/Ah865fZ2oQ8SvFHwS7l4xs2O0XuO7H5jnnNvim/cg0EDrbQ6L6MTRia8xy3xaj5ZeNrO406zrF7Q2NPmrmR2ltdHNNN829tN6Wva7QKlvW9+kg99n33W7q4BJwF5a7x98HEjuxNdyOuveBxT5ln+L1pCrb7fMImACOs0pHjHdPysiPcXM7gTmOudmtpk2GNgO9HfOVXtWnIQtHfGJiN+Y2QAzO9d3H+Jo4D+Av7SZHwF8A3hGoSdeUatOEfGnGOBRYBhQSev9gb8G8DX4OQR8ROutDCKe0KlOEREJKzrVKSIiYSUkT3X269fPDR061OsyREQkSKxZs6bMOfeJbvs6EpLBN3ToUPLz870uQ0REgoSZfdTZZXWqU0REwoqCT0REwoqCT0REwoqCT0REwoqCT0REwoqCT0REwoqCT0REwkpYBt/Bqjoe+OsOPiqv8boUEREJsLAMvoamFh76226W7yz1uhQREQmwsAy+QWnxZCfHsbKg3OtSREQkwMIy+MyM6bnprCyoQKNTiIiEl7AMPoDpw9OpqGlg56FjXpciIiIBFLbBNyM3HUCnO0VEwkzYBt/A1HhyUuIVfCIiYSZsg+8f1/nKaWnRdT4RkXARtsEHMD03jSPHG9l5+KjXpYiISID4JfjMbJaZ7TCz3WZ2Twfzv2FmW81so5ktM7MhbebNM7Ndvsc8f9TTWdNPXOfbo9OdIiLhotvBZ2aRwMPA5cA44AYzG9dusXVAnnNuIvAC8L++ddOAe4FpwFTgXjNL7W5NnTUoLYGBqfGs0HU+EZGw4Y8jvqnAbudcgXOuAXgGmNN2Aefc28654763K4GBvteXAUudcxXOuSPAUmCWH2rqtOm56azaW6HrfCIiYcIfwZcD7G/zvsg37WRuA17v6rpmNt/M8s0sv7TUf12NTc9Np/J4IzsO6TqfiEg48EfwWQfTOjx8MrObgDzgZ11d1zm30DmX55zLy8jIOK1COzI9Nw2AFbrOJyISFvwRfEXAoDbvBwIl7Rcys4uB7wGznXP1XVm3Jw1MTWBQmu7nExEJF/4IvtXASDMbZmYxwFxgcdsFzGwy8CitoXe4zaw3gUvNLNXXqOVS37SAmj5M1/lERMJFt4PPOdcE3EVrYG0DnnPObTGz+8xstm+xnwF9gefNbL2ZLfatWwH8kNbwXA3c55sWUDOGp1NV28j2g7rOJyLS20X540Occ0uAJe2mfb/N64s/Zd0ngCf8Ucfpmua7n29FQTnjspO8LEVERHpYWPfcckJOSjyD0xJ0nU9EJAwo+Hym56bxoa7ziYj0ego+nxPX+bYeqPa6FBER6UEKPp9pwzQ+n4hIOFDw+WSnxDMkPYGVBQFvVCoiIgGk4GtjRm46q/aW06zrfCIivZaCr43puekcrWtim67ziYj0Wgq+Nj4en0/X+UREei0FXxv9k+MY1q+POqwWEenFFHztnLifT9f5RER6JwVfO9Nz0zla38TWEl3nExHpjRR87Uz/uN/OMo8rERGRnqDgaycrKY4RmX1ZvtN/o7yLiEjwUPB14KKxmawqqKC6rtHrUkRExM8UfB24eGwWTS2Od3XUJyLS6yj4OjBlcCqpCdG8tfWQ16WIiIif+SX4zGyWme0ws91mdk8H8883s7Vm1mRm17Wb1+wblf3jkdm9FhlhfHZ0Jm/vKKWpucXrckRExI+6HXxmFgk8DFwOjANuMLNx7RbbB9wC/LGDj6h1zk3yPWZ3tx5/uXhcFlW1jaz56IjXpYiIiB/544hvKrDbOVfgnGsAngHmtF3AOVfonNsIhMzh03kj+xEdaSzbftjrUkRExI/8EXw5wP4274t80zorzszyzWylmV3th3r8IjEumum56by1Tdf5RER6E38En3UwrSv9fQ12zuUBNwL/Z2bDO9yI2XxfQOaXlgamteVFYzIpKK2hoPRYQLYnIiI9zx/BVwQMavN+IFDS2ZWdcyW+5wLgHWDySZZb6JzLc87lZWRknH61XXDR2CwAlm3T6U4Rkd7CH8G3GhhpZsPMLAaYC3SqdaaZpZpZrO91P+BcYKsfavKLQWkJjOmfqNOdIiK9SLeDzznXBNwFvAlsA55zzm0xs/vMbDaAmZ1tZkXA9cCjZrbFt/pYIN/MNgBvAz9xzgVN8EFrLy75Hx2h6rh6cRER8SfnHI+9W8DOQ0cDut0of3yIc24JsKTdtO+3eb2a1lOg7df7AJjgjxp6ykVjs3j47T28s/MwcyZ1pc2OiIh8mqIjtdy/ZBsJsZGMykoM2HbVc8spTBqYQr++Mbyl63wiIn61qbgKgAk5yQHdroLvFCIijAvHZPLOjsM0qhcXERG/2VxcRVSEMbp/4I72QMHXKReNzeJoXROr91Z4XYqISK+xqbiKUVmJxEZFBnS7Cr5OOG9kP2KiInS6U0TET5xzbCmpDvhpTlDwdUpCTBTnDE9n2fZDONeVe/NFRKQjJVV1VNQ0MH6ggi9oXTw2i4/Kj7NHvbiIiHTbpqLWhi3js5MCvm0FXyddNDYTgKVbdbpTRKS7tpRUERlhjB2g4AtaA5LjOSM7iWXqxUVEpNs2FVcxMrMvcdGBbdgCCr4uuWhsFmv3HaGipsHrUkREQpZzjs3FVYz3oGELKPi65JKxWbQ4eFtj9ImInLZD1fWUHWvwpEUnKPi6ZHxOEllJseq0WkSkG0702DI+J/DX90DB1yVmxqXj+vP2jsMcrVOn1SIip2NTcRURBuMG6IgvJFwzJYe6xhZe33zQ61JERELSluIqRmT2JT4m8A1bQMHXZZMHpTCsXx/+vLbI61JERELSpuIqxmd7c7QHCr4uMzOunZzDyoIKio4c97ocEZGQcri6jsNH6z1r0QkKvtNy9eTWcfleWlfscSUiIqFlc4lvKCIPuio7QcF3GgalJTBtWBovri1W350iIl2wqagaMxjnQY8tJ/gl+MxslpntMLPdZnZPB/PPN7O1ZtZkZte1mzfPzHb5HvP8UU8gfH7KQPaW1bBuf6XXpYiIhIxNxVXk9utDn9goz2rodvCZWSTwMHA5MA64wczGtVtsH3AL8Md266YB9wLTgKnAvWaW2t2aAuHyCf2Ji45QIxcRkS7YUlLl2Y3rJ/jjiG8qsNs5V+CcawCeAea0XcA5V+ic2wi0H8L8MmCpc67COXcEWArM8kNNPS4xLprLzujPKxsOUN/U7HU5IiJBr+xYPQeq6jxt2AL+Cb4cYH+b90W+aX5d18zmm1m+meWXlpaeVqH+du2UgVTVNqoLMxGRTvhHjy2hH3zWwbTOtvjo9LrOuYXOuTznXF5GRkani+tJ5w5PJzMxlhfXqnWniMipbPEF3xkejMHXlj+CrwgY1Ob9QKAkAOt6Lioygqsn5/D29sOUH6v3uhwRkaC2qbiKYf36kBgX7Wkd/gi+1cBIMxtmZjHAXGBxJ9d9E7jUzFJ9jVou9U0LGddOyaGpxfHKhpDJaxERT2wurvb8NCf4Ificc03AXbQG1jbgOefcFjO7z8xmA5jZ2WZWBFwPPGpmW3zrVgA/pDU8VwP3+aaFjDH9kxg3IIk/62Z2EZGTqqhpoLiylvEen+YE8MuNFM65JcCSdtO+3+b1alpPY3a07hPAE/6owyvXTsnhR69tY9eho4zMSvS6HBGRoLPZd33P61sZQD23+MWcSTlERpiO+kRETuJEi84zFHy9Q0ZiLDNHZfDSumKaW9SFmYhIe1tKqhiclkByvLcNW0DB5zfXTsnhQFUdK/aUe12KiEjQ2VTsfY8tJyj4/OTisVkkxkWpCzMRkXaqjjeyv6I2KFp0goLPb+KiI7ly4gBe33yQ6rpGr8sREQkaJ4YiGp/jfYtOUPD51Y1Th1Db2MwzH+7zuhQRkaDxcVdlHo663paCz48mDExmRm46T7xfSENT+/64RUTC0+biKgamxpPaJ8brUgAFn9/Nn5nLweo6FqsnFxERoDX4guVoDxR8fnfBqAxGZyXy2LsFGp1dRMJedV0jheXHmTBQwddrmRlfPT+XHYeO8s7O4Bg+SUTEK5uDZCiithR8PWD2mdn0T4pj4fICr0sREfHUun2VAExU8PVuMVER3HruUFYUlLOpqMrrckREPLO6sIKRmX2DpmELKPh6zA3TBtM3NopH393jdSkiIp5obnGsKTzC2cPSvC7lnyj4ekhSXDQ3ThvMkk0H2F9x3OtyREQCbvvBao7WNzF1qIIvbNx67lAizPjt+3u9LkVEJODyC48AkDc01eNK/pmCrwcNSI5n9qRsnl29nyM1DV6XIyISUB8WVpCdHMfA1ASvS/knfgk+M5tlZjvMbLeZ3dPB/Fgze9Y3f5WZDfVNH2pmtWa23vd4xB/1BJP55+dS29jM0ys/8roUEZGAcc6RX1hBXpCd5gQ/BJ+ZRQIPA5cD44AbzGxcu8VuA44450YADwI/bTNvj3Nuku9xR3frCTZj+icxc1QGi1YUUtfY7HU5IiIBsb+ilkPV9UHXsAX8c8Q3FdjtnCtwzjUAzwBz2i0zB1jke/0CcJGZmR+2HRIWnJ9L2bEG/rxWI7SLSHhYXVgBwNlBdn0P/BN8OcD+Nu+LfNM6XMY51wRUAem+ecPMbJ2ZLTez8/xQT9CZMTyd8TlJPP5egUZoF5GwsLqwguT4aEZlJnpdyif4I/g6OnJr/9f9ZMscAAY75yYD3wD+aGYdDthkZvPNLN/M8ktLQ6srMDNj/vnDKSir4bVNB7wuR0Skx60urCBvSCoREcF3cs8fwVcEDGrzfiDQfmiCj5cxsyggGahwztU758oBnHNrgD3AqI424pxb6JzLc87lZWRk+KHswLpiwgDGDUjiJ0u2Uduga30i0nuVH6tnT2lNUDZsAf8E32pgpJkNM7MYYC6wuN0yi4F5vtfXAX9zzjkzy/A1jsHMcoGRQK/s4DIywvjB7DMoqarjkeXqzUVEeq/Vvvv3pg4Lvut74Ifg812zuwt4E9gGPOec22Jm95nZbN9ivwXSzWw3rac0T9zycD6w0cw20Nro5Q7nXEV3awpWU4elceXEATyyfA9FR9Sbi4j0TvmFFcRERQTViAxtWSiOGZeXl+fy8/O9LuO0lFTWcuED73DhmEx+/aWzvC5HRMTv5vzqfWKjI3luwYyAbdPM1jjn8jqzrHpuCbDslHjunDmCJZsO8sGeMq/LERHxq+MNTWwuqQ7K2xhOUPB5YMHMXHJS4rnvla00Nbd4XY6IiN+s21dJc4vj7CBt2AIKPk/ERUfyX1eMZfvBo/zpw31elyMi4jerCyuIMDhriI74pJ1Z4/szIzedB5buVAfWItJrrC6sYEz/JBLjor0u5aQUfB4xM+6dPY7q2kZ+vnSn1+WIiHRbY3ML6/ZVMjUI++dsS8HnoTH9k7hp+hD+sOojth2o9rocEZFu2VpSzfGG5qAbf689BZ/HvnHJKJLio/mfV7YQireWiIic8I+OqXXEJ58iJSGG/7hkFCsLKnh1o/rxFJHQtbqwgsFpCWQlxXldyqdS8AWBG6YOZkJOMt/7yyb2V6hHFxEJPa0Dzx4J+qM9UPAFhajICB6+cQoAd/5hjQasFZGQU1BWQ3lNQ9D2z9mWgi9IDE5P4IEvTGJzcTX3vbrV63JERLpk9d7W63vBOiJDWwq+IHLJuCzumDmcP67ax5/XFnldjohIp60uPEJ6nxhy+/XxupRTUvAFmf+8dBTThqXxvb9sZsfBo16XIyLSKasLK8gbmopZ8A08256CL8hERUbw0I2T6RsXxZ1/WMOx+iavSxIR+VSHquvYV3E8JBq2gIIvKGUmxvHQDZMpLKvh2y9u1P19IhLUQuX+vRMUfEFqem4637xsDK9tPMCiDwq9LkdE5KRWFVSQEBPJuOwkr0vpFAVfEFtwfi4Xj83i/iXbWLvviNfliIh8gnOO5TtLOWd4OtGRoREpfqnSzGaZ2Q4z221m93QwP9bMnvXNX2VmQ9vM+45v+g4zu8wf9fQWERHGA9efSf/kOOY/lc/WEvXnKSLBpaCshn0Vx7lgdKbXpXRat4PPzCKBh4HLgXHADWY2rt1itwFHnHMjgAeBn/rWHQfMBc4AZgG/9n2e+CQnRPPkrVOJjoxg7sIVOvITkaDy9vbDAFwwOsPjSjovyg+fMRXY7ZwrADCzZ4A5QNu7sOcAP/C9fgH4lbW2eZ0DPOOcqwf2mtlu3+et8ENdvcbwjL48t2AGN/12FTc9vorH5+VxzvB+XpfVqzW3OEoqaymurKW6tpHquibfcyNHfa9rGpqIioggNiqCuOhIYqMiiI2OIDaq9XVmUiyD0/owOC2Bfn1jQqKZt0hXLd9ZysjMvgxMTfC6lE7zR/DlAPvbvC8Cpp1sGedck5lVAem+6SvbrZvT0UbMbD4wH2Dw4MF+KDu0DEpL4Hlf+N3yu9U8ctMULhyT5XVZIe9gVR07Dx2lsLyGwrLjrc/lNeyvOE5jc8etafvERJIUH02f2CiaWxx1jc3UN7VQ73tuavnkegkxkQxOS2BQWgKD0xIYnZXI1GFpDElPUCBKyKqpb2JVQQW3nDvU61K6xB/B19Fvbfvf/JMt05l1Wyc6txBYCJCXlxeW7fszk+J4Zv4M5j3xIfOfWsP/zZ3ElROzvS4rZFQdb2RjcSUb9leyoaiKDfsrOXy0/uP58dGRDO3Xh9FZiVx2Rn+GpicwMDWB5PhokuKiSYqPom9sFFGnuIDf1NxCXVMLB6vq2F9xnI/Ka9hXUcs+3+v3dpVS19gCQGZiLFOHpTFtWBpTh6UzMrMvEREKQgkNH+wpp6G5hQtGhc5pTvBP8BUBg9q8HwiUnGSZIjOLApKBik6uK22k9Ynhj1+dxm1P5nP3n9ZxvL6ZL5w96NQrhqH9FcdZsaecFQXlrNt3hMLyf4x8kZvRh3NH9GPiwGTGDkhiWL8+ZCbG+uXoKyoygr6REYzI7MuIzL6fmO+cY0/pMVbtreDDvRWsajMkVUpCNDNy07lyYjYXjc0kLlqXvCV4vb3jMH1iIkOif862/BF8q4GRZjYMKKa1scqN7ZZZDMyj9drddcDfnHPOzBYDfzSznwPZwEjgQz/U1KslxkWz6CtTWfD0Gr714kaq6xq57TPDwv6U2aHqOlbsKeeDPWV8sKecoiO1APTrG8NZQ1K5Pm8QkwalMD4nmeT4aM/qNDNGZCYyIjORL00bgnOOoiO1viAs550dpby++SB9Y6O47Iz+zJmUzTnD0095pCkSSM45lu8o5TMj+xETFVo/m90OPt81u7uAN4FI4Ann3BYzuw/Id84tBn4L/N7XeKWC1nDEt9xztDaEaQK+5pzTmDydEB8TyWNfPot//dN6fvTaNv6+u4wfXTOBnJR4r0sLmLrGZlYWlLN8ZynLd5ZSUFoDQHJ8NNNz0/jqebnMGN56+jCY/ykwMwb5rv9dd9ZAmlscqwrKeWl9Ma9vPsiLa4vo1zeWKycOYM6kbCYNSgnqr0fCw67DxyiurOXrF47wupQus1DsDisvL8/l5+d7XUZQaG5xLPqgkJ+9uYMIg2/NGsPN04f0yutEzjkKy4/zzo7DvLOjlJUF5dQ3tRAbFcG03HTOG9GPGcPTGTsgiche8vXXNTbzzo5SXl5fzLLth2loauHMQSncOTOXS8b17zVfp4SeR5fv4cevb2fFdy5kQLL3/3Cb2RrnXF6nllXw9Q77K47z3b9s4r1dZUwZnMJPPz+RkVmJXpfVbbUNzawoKOOdHa1HdR/5rtPl9uvDzNEZzByVwfTc9LC4FlZd18jL60t4/L0CPio/Tm6/Pnz1/FyumZwTFl+/BJe5C1dQebyRN/7tfK9LARR8Ycs5x1/WFXPfq1upqW/ia58dwb9cMCKkzr875ygoq+GdHaW8s+Mwq/ZW0NDUQnx0JOcMT+eC0RnMHJXJ4PTQuWfI35pbHG9sPsgjy/ewqbiKjMRYbj13KF+aNsTTa5cSPo7WNTL5vqV89fxcvj1rjNflAAq+sFd2rJ77XtnK4g0ljMzsyx0zh/O5CQOIjwnOo4IjNQ2sKGhtlPLuzjL2VbQe1Q3P6MMFozO5YHQGZw9N01FNO845Vuwp55F3C3h3Zyl9Y6O47TPDWDAzl4QYf7RbE+nYG5sPcMfTa3l2/nSm5aZ7XQ6g4BOfv20/xA9f3cbeshoSY6OYPSmbL549iAk5yZ42jqipb+LDwgo+2F3G33eXs+1gNc613hg+Pbf1qO6C0ZkMSgvfo7qu2lpSzcNv7+a1TQfonxTHty8fzZwzc3rltV7x3rdf2MiSzQdY+9+XBE3H1Ao++ZhzjlV7K3hu9X6WbD5AXWMLY/on8sWzB3HN5BxSEmJ6dPstLa2nLjcWVbKxqIoNRZVsKqqiqcURExnB5MEpnDuiH+eOSGfiwJSg+SUKVfmFFfzw1a1sKKrizIHJfP+qcZw1JLTusZLg5pxj+o+XkTckjYe/NMXrcj6m4JMOVdc1snh9Cc/l72djURUxkRFMGpTC6P6JjBmQyJj+iYzKSiQxruvXiZpbHGXH6imurGV/xXG2lFSzYX8lW0qqPx5FPj46kjOyk8gbmsa5I9LJG5IWtKdfQ1lLi+Ol9cX89I3tHKqu56ozs/n2rNEh1ZeiBK+tJdV87pfv8b/XTeQLecHTeYaCT05p24FqXlxTxPr9lew4eJSjvnACGJgaz5j+iQxIjicywogwIzICIiMiWp/NaGh2HKyqpaSyjpKqWg5W1f1TH5UxURGMHZDEmQOTmZCTzMSBKQzP6KObsAPoeEMTjywv4NHlewC484LhIdfYSYLPw2/v5mdv7uDD715EZlKc1+V8TMEnXeKco7iylh0Hj7L9xONANWXH6mlucbS41iO65hZHs2t9jo40+ifHkZ0cT3ZKPAOS48hOiSc7pfU5t19f/YENEsWVtfx4yTZe3XiAsQOSeOD6M0NmpGwJPtc/8gG1jc28+vXzvC7ln3Ql+NT0SzAzBqa2dsh80djOjfjgnFPvISEiJyWeX904hdlnHuS7f9nM7F+9z9cvHMm/fHa4rqlKl1Qdb2TtvkrunDnc61K6RT/1cloUeqHn0jP6s/Tfz+eKiQN48K2dXP3w39l+sNrrsiSEvLe7lOYWx2fHhNZoDO0p+ETCSGqfGH4xdzKP3HQWh6rruOqh93lo2S4am1u8Lk1CwNvbS0lJiGbSoFSvS+kWBZ9IGJo1vj9//feZXHZGfx5YupPP/+YD9rUZtkmkvZYWx/KdhzlvZEbI9xGr4BMJU2l9YvjVjVN4+MYpFJbVcMVD7/HmloNelyVBaktJNWXHGvjs6NA+zQkKPpGwd8XEAbx293kMTe/Dgt+v4f7XturUp3zC2zsOYwbnh9ho6x1R8IkIg9ISeOHOGdw8fQiPvbeXGxau5EBVrddlSRBZtv0wE3OS6dc31utSuk3BJyIAxEZF8sOrx/OLuZPYeqCaK375Pu/tKvW6LAkCH5XXsGF/JbPGD/C6FL/oVvCZWZqZLTWzXb7nDpv6mNk83zK7zGxem+nvmNkOM1vve2R2px4R6b45k3JYfNdn6Nc3hi8/8SEPLt1Jc0vodXQh/vPy+hIA5kzK9rgS/+juEd89wDLn3Ehgme/9PzGzNOBeYBowFbi3XUB+yTk3yffhE9MAAAAUJUlEQVQ43M16RMQPRmT25aWvncs1k3P4xbJdLPj9GmradGsn4cM5x0vripmem0Z2ivcjrftDd4NvDrDI93oRcHUHy1wGLHXOVTjnjgBLgVnd3K6I9LCEmCgeuP5M/mf2Gfxt+yGue2QFxZW67hduNhZVUVBWw9WTcrwuxW+6G3xZzrkDAL7njk5V5gD727wv8k074Xe+05z/beoORCSomBnzzhnK726dSlHFceb86u+s3XfE67IkgF5aX0xMZASXT+gd1/egE8FnZm+Z2eYOHnM6uY2OwuzEBYMvOecmAOf5Hjd/Sh3zzSzfzPJLS3XBXSSQZo7K4M//cg4JMZHMXbiSl9cXe12SBEBTcwuvbCjhorGZJMd3fbiyYHXK4HPOXeycG9/B42XgkJkNAPA9d3SNrghoO2jTQKDE99nFvuejwB9pvQZ4sjoWOufynHN5GRmhfx+JSKgZmZXIS187l0kDU/jXZ9bz86U7aVGjl17t/d1llB1rYE4vOs0J3T/VuRg40UpzHvByB8u8CVxqZqm+Ri2XAm+aWZSZ9QMws2jgSmBzN+sRkR6U1ieG398+levOGsgvl+3i68+so7ah2euypIe8vL6EpLiokO+Uur3uBt9PgEvMbBdwie89ZpZnZo8DOOcqgB8Cq32P+3zTYmkNwI3AeqAYeKyb9YhID4uNiuRn103kO5ePYcmmA3zp8ZUcqWnwuizxs5r6Jt7YfJArJmYTGxXpdTl+1a3x+Jxz5cBFHUzPB25v8/4J4Il2y9QAZ3Vn+yLiDTNjwczhDE5L4F+fXc91j3zAU7dNI6eXNHcXWLr1ELWNzVzdS+7da0s9t4jIabt8wgCe+spUDh+t59pfa3y/3uQv64rJSYnn7KFpXpfidwo+EemW6bnpPH/HDACuf2QFqwrKPa5Iuqv0aD3v7y5jzqRsIkJ8CKKOKPhEpNvG9E/ixTvPISMxlpuf+JA3Nh/wuiTphlc3ltDc4rhmcu9qzXmCgk9E/GJgagIv3nEOZ2Qncecf1vL7lR95XZKcppfWFTNuQBIjsxK9LqVHKPhExG9S+8Twx9unc+HoTP77pc38fOlOnNO9fqGkoPQYG4qqeu3RHij4RMTP4mMiefTms7jed6/f/7yyVTe6h5CX1pdgBrN7YWvOE7p1O4OISEeiIiP46ecn0jcuit/9vZCa+iZ+8vmJRPbChhK9yYmRGM4Znk5WUpzX5fQYBZ+I9IiICOP7V44jMS6aXy7bxfGGZh784iRionSiKVit21/JvorjfP3CEV6X0qMUfCLSY8yMb1wyisTYKO5fso2ahiZ+86WziI/pXT2B9BYvrSsmNiqCWeP7e11Kj9K/XiLS4756fi7/75oJLN9ZyrzffcjRukavS5J2GppaeHXjAS4el0ViXO8ZiaEjCj4RCYgbpw3m/744ibUfHeFLj69S/55B5uX1xVTUNPCFvEGnXjjEKfhEJGDmTMrh0ZvPYvvBo3xx4QoOH63zuiQBWlocj71XwJj+iZw/sp/X5fQ4BZ+IBNRFY7N48pazKTpSy9xHV3KgqtbrksLeOzsPs/PQMRbMzMWs97e8VfCJSMCdM6Lfx51bf+HRFeyvOO51SWHt0eUFZCfHceXE3nvvXlsKPhHxRN7QNP5w+zSqa5v4wqMrKCg95nVJYWn9/kpW7a3gK58ZRnRkeERCeHyVIhKUzhyUwp++Op2Gpha+uHAlOw8d9bqksLPw3T0kxkUxd+pgr0sJGAWfiHhqXHYSz8yfjgFzF65kS0mV1yWFjcKyGl7ffJCbpw+hb2z43NbdreAzszQzW2pmu3zPqSdZ7g0zqzSzV9tNH2Zmq3zrP2tmMd2pR0RC08isRJ5bMIO4qAhuWLiS9fsrvS4pLDz+fgHRERHccs5Qr0sJqO4e8d0DLHPOjQSW+d535GfAzR1M/ynwoG/9I8Bt3axHRELU0H59eHbBDFISYrjp8VV8uLfC65J6tbJj9TyfX8S1U3LI7MX9cnaku8E3B1jke70IuLqjhZxzy4B/OnlvrW1mLwReONX6IhIeBqUl8NyCGWQmxTLviQ95f1eZ1yX1Wk+t+Ij6phZuPy/X61ICrrvBl+WcOwDge87swrrpQKVzrsn3vgg46QBQZjbfzPLNLL+0tPS0CxaR4NY/OY5n589gSHoCX1m0mr9tP+R1Sb3O8YYmnlpRyCXjshiR2dfrcgLulMFnZm+Z2eYOHnO6ue2O7pI86aBdzrmFzrk851xeRkZGNzctIsEsIzGWP311OmP6JzL/qTUs2XTA65J6lefzi6g83siC88PvaA86EXzOuYudc+M7eLwMHDKzAQC+58Nd2HYZkGJmJ5oSDQRKuvoFiEjvlNonhqdvn8aZg1K4649r+cu6Iq9L6hWamlt47L0CzhqSSt7QNK/L8UR3T3UuBub5Xs8DXu7sis45B7wNXHc664tI75cUF81TX5nK9Nx0vvHcBv704T6vSwp5r28+SNGRWuaH6dEedD/4fgJcYma7gEt87zGzPDN7/MRCZvYe8DxwkZkVmdllvlnfBr5hZrtpveb3227WIyK9TJ/YKJ645WxmjsrgO3/exBPv7/W6pJDlnOPRd/eQ268Pl4zN8rocz3TrjkXnXDlwUQfT84Hb27w/7yTrFwBTu1ODiPR+cdGRPHrzWdz9p3Xc9+pWjjc08bXPjgiLDpX9afnOUjYXV/PjaycQERG++049t4hISIiNiuThG6dwzeQc/r+/7uT+17bResVEOqOusZkfLN7CsH59uHbKSRvQh4Xw6aNGREJeVGQED1x/Jsnx0Tz+/l6qahv58bUTiAqTzpW745HleygsP87vb5tKbFSk1+V4SsEnIiElIsK496pxJMdH84tlu6iua+QXcycTFx3ef8w/TWFZDb9+Zw9XThzAeSN1O5j+TRKRkGNm/Pslo7j3qnG8ueUQX3lyNcfqm069YhhyzvHfL28mJjKC/75ynNflBAUFn4iErFvPHcbPv3Amq/ZW8KXHVnKkpsHrkoLOkk0HeW9XGf9x6SiywqxPzpNR8IlISLt2ykAeveksth08yvWPruBAVa3XJQWNo3WN3PfqFs7ITuLm6UO8LidoKPhEJORdPC6LRbdO5WBVHZ//9QdsP1jtdUlB4cGluzh8tJ77r1EDoLa0J0SkV5gxPJ1nF0yn2Tmu+80K3t0Z3p3Zbymp4skP9nLj1MFMGpTidTlBRcEnIr3GGdnJvPS1cxmYGs+tT64O2y7OWloc//XSZtL6xPCty8Z4XU7QUfCJSK8yIDme5++YwWdG9OM7f97ET9/YTktLeN3o/szq/azbV8l3PzeW5IRor8sJOgo+Eel1EuOi+e28PG6cNpjfvLOHu59ZR11js9dlBUTZsXp++sZ2pg1L45rJ4d1Dy8noBnYR6ZWiIiO4/+rxDElL4Mevb+dAVR2PfTmPtD4xXpfWY5qaW/jWCxupqW/iR1ePV1+mJ6EjPhHptcyMBTOH8/CNU9hUXMXVD/+dLSVVXpfVI5xz/OCVLfxt+2HuvWocI7MSvS4paCn4RKTXu2LiAJ6ZP536pmau/fUHPLd6v9cl+d3Cdwt4euU+Fpyfy80zhnpdTlBT8IlIWJgyOJXX7j6PvKGpfOvFjXzz+Q3UNvSO636vbizhx69v54qJA/j2LLXiPBUFn4iEjX59Y3nqK9O4+8IRPL+miGt+/Xf2ltV4XVa35BdW8I3nNpA3JJUHrj8zrMfZ66xuBZ+ZpZnZUjPb5XtOPclyb5hZpZm92m76k2a218zW+x6TulOPiMipREYY37h0NL+79WwOVtcx+6H3eWPzAa/LOi0Fpce4/al8BqbE89iX8zRCRSd194jvHmCZc24ksMz3viM/A24+ybxvOucm+R7ru1mPiEinfHZ0Jq/dfR65mX254+m1/OjVrdQ3hc6pz7Jj9dzyu9VEmvG7W88mtRe3VvW37gbfHGCR7/Ui4OqOFnLOLQOOdnNbIiJ+lZMSz/MLZnDLOUN5/P29fO4X77GqoNzrsk6ptqGZ2xflc/hoHY/Py2NIeh+vSwop3Q2+LOfcAQDfc+ZpfMb9ZrbRzB40s9hu1iMi0iUxURH8YPYZLPrKVBqaW/jiwpV8+4WNVB4PziGOKmoaWPD0GjYUVfKLuZOZPLjDK0zyKU4ZfGb2lplt7uAxxw/b/w4wBjgbSAO+/Sl1zDezfDPLLy0N785nRcT/Zo7K4K//NpMFM3N5YW0RF/98OS+vL8a54Onu7IPdZVz+i3dZuaecH18zgcvO6O91SSHJuvNNNbMdwAXOuQNmNgB4xzk3+iTLXgD8p3PuytOZ31ZeXp7Lz88/7bpFRD7NlpIqvvvnTWwoquL8URncf/V4BqUleFZPY3MLP1+6k0eW72FYvz78cu5kxucke1ZPMDKzNc65vM4s291TnYuBeb7X84CXu7KyLyyx1n51rgY2d7MeEZFuOyM7mT//y7nce9U41hRWcMmDy/nx69s4WFUX8FoKy2q47jcf8Jt39jD37MG8+vXPKPS6qbtHfOnAc8BgYB9wvXOuwszygDucc7f7lnuP1lOafYFy4Dbn3Jtm9jcgAzBgvW+dY6faro74RCRQSipr+fHr23ltYwmREcbsM3P46vnDGNM/qUe365zjxbXF3PvyZqIiI/jJtRO4fMKAHt1mKOvKEV+3gs8rCj4RCbR95cd54u97eXb1fmobmzl/VAbzz8vl3BHpfu0MuqXFsbqwgkUrClmy6SDThqXx4BcnkZ0S77dt9EYKPhGRHlJ5vIGnV37Ekx98RNmxesYOSOILeQM5e2gaYwckEXkaPac459hSUs3iDSW8sqGEA1V1xEVHcNdnR3DnBSNO6zPDjYJPRKSH1TU28/L6Yh5/by+7DrdeoekbG8XkwSmcPTSNvKGpTBqUQkLMJ0d/a25x1DU2c6Cqjtc2HmDxhmL2lNYQFWHMHJXB7EnZXDw2iz6xGjmusxR8IiIBVFxZS35hBasLK8gvPMKOQ0dxDqIijIGp8TQ2twZdfVMLdY3NNLUbEX7asDRmT8rmc+MHqAeW09SV4NO/EyIi3ZSTEk/OpBzmTGod8bzqeCNr9x1hdWEF+4/UEhsVQWxUBHHRkf/0nBgXzQWjM3T9LsAUfCIifpacEM1nx2Ty2TGn05mV9DQNSyQiImFFwSciImFFwSciImFFwSciImFFwSciImFFwSciImFFwSciImFFwSciImElJLssM7NS4KMe3kw/oKyHtxFqtE86pv3ySdonHdN+6Zg/9ssQ51xGZxYMyeALBDPL72y/b+FC+6Rj2i+fpH3SMe2XjgV6v+hUp4iIhBUFn4iIhBUF38kt9LqAIKR90jHtl0/SPumY9kvHArpfdI1PRETCio74REQkrCj4REQkrCj4fMzsejPbYmYtZnbSZrVmNsvMdpjZbjO7J5A1BpqZpZnZUjPb5XtOPclyzWa23vdYHOg6A+VU33szizWzZ33zV5nZ0MBXGVid2Ce3mFlpm5+P272oM5DM7AkzO2xmm08y38zsl759ttHMpgS6Ri90Yr9cYGZVbX5Wvt9TtSj4/mEzcC3w7skWMLNI4GHgcmAccIOZjQtMeZ64B1jmnBsJLPO970itc26S7zE7cOUFTie/97cBR5xzI4AHgZ8GtsrA6sLvw7Ntfj4eD2iR3ngSmPUp8y8HRvoe84HfBKCmYPAkn75fAN5r87NyX08VouDzcc5tc87tOMViU4HdzrkC51wD8Awwp+er88wcYJHv9SLgag9r8Vpnvvdt99cLwEVmZgGsMdDC7fehU5xz7wIVn7LIHOAp12olkGJmAwJTnXc6sV8CRsHXNTnA/jbvi3zTeqss59wBAN9z5kmWizOzfDNbaWa9NRw7873/eBnnXBNQBaQHpDpvdPb34fO+U3ovmNmgwJQW1MLt70hXzDCzDWb2upmd0VMbieqpDw5GZvYW0L+DWd9zzr3cmY/oYFpI3w/yafukCx8z2DlXYma5wN/MbJNzbo9/Kgwanfne97qfj1PozNf7CvAn51y9md1B6xHxhT1eWXALt5+TzlpLa3+bx8zsc8BLtJ4O9ruwCj7n3MXd/IgioO1/rAOBkm5+pqc+bZ+Y2SEzG+CcO+A7FXP4JJ9R4nsuMLN3gMlAbwu+znzvTyxTZGZRQDJBcmqnh5xynzjnytu8fYxeft2zk3rd3xF/cM5Vt3m9xMx+bWb9nHN+79Rbpzq7ZjUw0syGmVkMMBfota0Yaf3a5vlezwM+cVRsZqlmFut73Q84F9gasAoDpzPf+7b76zrgb6539xBxyn3S7trVbGBbAOsLVouBL/tad04Hqk5cUghnZtb/xDVxM5tKaz6Vf/papyesjvg+jZldAzwEZACvmdl659xlZpYNPO6c+5xzrsnM7gLeBCKBJ5xzWzwsu6f9BHjOzG4D9gHXA/hu97jDOXc7MBZ41MxaaP1B/YlzrtcF38m+92Z2H5DvnFsM/Bb4vZntpvVIb653Ffe8Tu6Tu81sNtBE6z65xbOCA8TM/gRcAPQzsyLgXiAawDn3CLAE+BywGzgO3OpNpYHVif1yHXCnmTUBtcDcnvrHUV2WiYhIWNGpThERCSsKPhERCSsKPhERCSsKPhERCSsKPhERCSsKPhERCSsKPhERCSsKPpEQZGY/NbN/NbPNZrbazMZ6XZNIqFDwiYQYMzsHuBRYDxQD/wP8wtOiREKIgk8k9EylddQDAxqBN4CzPK1IJIQo+ERCT0f9DDYHvAqREKXgEwk97wFX4OvgF/i8b5qIdIJGZxAJMc65tWb2Aq3j26XROu7fTd5WJRI6NDqDSIgyswuA/3TOXel1LSKhRKc6RUQkrOiIT0REwoqO+EREJKwo+EREJKwo+EREJKwo+EREJKwo+EREJKz8/23ES97CcMA7AAAAAElFTkSuQmCC\n", "text/plain": [ "