{
"cells": [
{
"cell_type": "markdown",
"source": [
"# Post processing and visualization\n",
"\n",
"\n",
"\n",
"*Figure 1*: Heat flux computed from the solution to the heat equation on\n",
"the unit square, see previous example: Heat equation."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"## Introduction\n",
"\n",
"After running a simulation, we usually want to postprocess and visualize the results in\n",
"different ways. Ferrite provides several tools to facilitate these tasks:\n",
"\n",
" - L2 projection of (discrete) data onto FE interpolations using the `L2Projector`\n",
" - Evalutation of fields (solutions, projections, etc) at arbitrary, user-defined, points\n",
" using the `PointEvalHandler`\n",
" - Builtin functionality for exporting data (solutions, cell data, projections, etc) to\n",
" the VTK format\n",
" - [Makie.jl](https://docs.makie.org/) based plotting using\n",
" [FerriteViz.jl](https://ferrite-fem.github.io/FerriteViz.jl/)\n",
"\n",
"This how-to demonstrates the VTK export, the `L2Projector` for projecting discrete\n",
"quadrature point data onto a continuous FE interpolation, and the `PointEvalHandler` for\n",
"evaluating the FE solution, and the projection, along a user-defined cut line through the\n",
"domain."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"> **Custom visualization**\n",
">\n",
"> A common assumption is that the numbering of degrees of freedom matche the numbering\n",
"> of the nodes in the grid. This is *NOT* the case in Ferrite. If the available tools\n",
"> don't suit your needs and you decide to \"roll your own\" visualization you need to be\n",
"> aware of this and take it into account. For the specific case of evaluating the\n",
"> solution at the grid nodes you can use `evaluate_at_grid_nodes`.\n",
"\n",
"This example continues from the Heat equation example, where the temperature field was\n",
"determined on a square domain. In this example, we first compute the heat flux in each\n",
"integration point (based on the solved temperature field) and then we do an L2-projection\n",
"of the fluxes to the nodes of the mesh. By doing this, we can more easily visualize\n",
"integration points quantities. Finally, we visualize the temperature field and the heat fluxes along a cut-line.\n",
"\n",
"The L2-projection is defined as follows: Find projection $q(\\boldsymbol{x}) \\in U_h(\\Omega)$ such that\n",
"$$\n",
"\\int v q \\ \\mathrm{d}\\Omega = \\int v d \\ \\mathrm{d}\\Omega \\quad \\forall v \\in U_h(\\Omega),\n",
"$$\n",
"where $d$ is the quadrature data to project. Since the flux is a vector the projection function\n",
"will be solved with multiple right hand sides, e.g. with $d = q_x$ and $d = q_y$ for this 2D problem.\n",
"In this example, we use standard Lagrange interpolations, and the finite element space $U_h$ is then\n",
"a subset of the $H^1$ space (continuous functions).\n",
"\n",
"Ferrite has functionality for doing much of this automatically, as displayed in the code below.\n",
"In particular `L2Projector` for assembling the left hand side, and\n",
"`project` for assembling the right hand sides and solving for the projection."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"## Implementation\n",
"\n",
"Start by simply running the Heat equation example to solve the problem"
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"include(\"../tutorials/heat_equation.jl\");"
],
"metadata": {},
"execution_count": 1
},
{
"cell_type": "markdown",
"source": [
"Next we define a function that computes the heat flux for each integration point in the domain.\n",
"Fourier's law is adopted, where the conductivity tensor is assumed to be isotropic with unit\n",
"conductivity $\\lambda = 1 ⇒ q = - \\nabla u$, where $u$ is the temperature."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "compute_heat_fluxes (generic function with 1 method)"
},
"metadata": {},
"execution_count": 2
}
],
"cell_type": "code",
"source": [
"function compute_heat_fluxes(cellvalues::CellValues, dh::DofHandler, a::AbstractVector{T}) where {T}\n",
"\n",
" n = getnbasefunctions(cellvalues)\n",
" cell_dofs = zeros(Int, n)\n",
" nqp = getnquadpoints(cellvalues)\n",
"\n",
" # Allocate storage for the fluxes to store\n",
" q = [Vec{2, T}[] for _ in 1:getncells(dh.grid)]\n",
"\n",
" for (cell_num, cell) in enumerate(CellIterator(dh))\n",
" q_cell = q[cell_num]\n",
" celldofs!(cell_dofs, dh, cell_num)\n",
" aᵉ = a[cell_dofs]\n",
" reinit!(cellvalues, cell)\n",
"\n",
" for q_point in 1:nqp\n",
" q_qp = - function_gradient(cellvalues, q_point, aᵉ)\n",
" push!(q_cell, q_qp)\n",
" end\n",
" end\n",
" return q\n",
"end"
],
"metadata": {},
"execution_count": 2
},
{
"cell_type": "markdown",
"source": [
"Now call the function to get all the fluxes."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"q_gp = compute_heat_fluxes(cellvalues, dh, u);"
],
"metadata": {},
"execution_count": 3
},
{
"cell_type": "markdown",
"source": [
"Next, create an `L2Projector` using the same interpolation as was used to approximate the\n",
"temperature field. On instantiation, the projector assembles the coefficient matrix `M` and\n",
"computes the Cholesky factorization of it. By doing so, the projector can be reused without\n",
"having to invert `M` every time."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"projector = L2Projector(ip, grid);"
],
"metadata": {},
"execution_count": 4
},
{
"cell_type": "markdown",
"source": [
"Project the integration point values to the nodal values"
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"q_projected = project(projector, q_gp, qr);"
],
"metadata": {},
"execution_count": 5
},
{
"cell_type": "markdown",
"source": [
"## Exporting to VTK\n",
"To visualize the heat flux, we export the projected field `q_projected`\n",
"to a VTK-file, which can be viewed in e.g. [ParaView](https://www.paraview.org/).\n",
"The result is also visualized in *Figure 1*."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"VTKGridFile(\"heat_equation_flux\", grid) do vtk\n",
" write_projection(vtk, projector, q_projected, \"q\")\n",
"end;"
],
"metadata": {},
"execution_count": 6
},
{
"cell_type": "markdown",
"source": [
"## Point evaluation\n",
"\n",
"\n",
"*Figure 2*: Visualization of the cut line where we want to compute\n",
"the temperature and heat flux."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"Consider a cut-line through the domain like the black line in *Figure 2* above.\n",
"We will evaluate the temperature and the heat flux distribution along a horizontal line."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"points = [Vec((x, 0.75)) for x in range(-1.0, 1.0, length = 101)];"
],
"metadata": {},
"execution_count": 7
},
{
"cell_type": "markdown",
"source": [
"First, we need to generate a `PointEvalHandler`. This will find and store the cells\n",
"containing the input points."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"ph = PointEvalHandler(grid, points);"
],
"metadata": {},
"execution_count": 8
},
{
"cell_type": "markdown",
"source": [
"After the L2-Projection, the heat fluxes `q_projected` are stored in the DoF-ordering\n",
"determined by the projector's internal DoFHandler, so to evaluate the flux `q` at our points\n",
"we give the `PointEvalHandler`, the `L2Projector` and the values `q_projected`."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"q_points = evaluate_at_points(ph, projector, q_projected);"
],
"metadata": {},
"execution_count": 9
},
{
"cell_type": "markdown",
"source": [
"We can also extract the field values, here the temperature, right away from the result\n",
"vector of the simulation, that is stored in `u`. These values are stored in the order of\n",
"our initial DofHandler so the input is not the `PointEvalHandler`, the original `DofHandler`,\n",
"the dof-vector `u`, and (optionally for single-field problems) the name of the field.\n",
"From the `L2Projection`, the values are stored in the order of the degrees of freedom."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"u_points = evaluate_at_points(ph, dh, u, :u);"
],
"metadata": {},
"execution_count": 10
},
{
"cell_type": "markdown",
"source": [
"Now, we can plot the temperature and flux values with the help of any plotting library, e.g. Plots.jl.\n",
"To do this, we need to import the package:"
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"import Plots"
],
"metadata": {},
"execution_count": 11
},
{
"cell_type": "markdown",
"source": [
"Firstly, we are going to plot the temperature values along the given line."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=1}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdeUAU5f8H8GdmlmO5QeQ+5PZERUXNE/IEUREzz9IsNa/6amn1VSurb2Z2aZm3+bMyTU1FUFG8EEnxQDw4FJBDbrmXBXZn5vfHFhEuCrj3vl9/scOzy4dhed77zDzzDMXzPAEAANBXtLoLAAAAUCcEIQAA6DUEIQAA6DUEIQAA6DUEIQAA6DUEIQAA6DUEIQAA6DUEIQAA6DUEIQAA6DUEIQAA6DVdC8KDBw9evny5lY1ZllVqMdAUx3HqLkGP8DyPHa5K6ExUSeF7W9eC8OLFi9evX29l49raWqUWA02JxWJ0zSojlUobGhrUXYUeQWeiSgrf27oWhAAAAG2CIAQAAL2GIAQAAL2GIAQAAL2GIAQAAL2GIAQAAL2GIAQAAL2GIAQAAL0mUHcBANBOYikpEPMFtaRIzOfXkmIxnycixWI+v5YUignH0fYmUmdTYmdMOZkQBxPKQUgcTSh7IXE2pUzxrw/wN/w3AGiukjpSJOYLakmhmC+sJfm1fLGY5In44jqSL+IlPHEUUo4mxF5IOZkQOyE12IHYC2lHIbE1lPIcX8kxj0SkuI7PryUZVXx8Ecmv5YrEJF/EE0KcTSk7IXEyoRyExMGEcjIhdsaUsymxE1L2QkKp+3cHUBkEIYBGEEnJmUdcbD7/sPqfEZ6lIbEXUo4mxEFIOZgQdzOqX0fibELbCYmzKWVh0OKrSSSEZXl3Y8rfhhB5oSaSkkcivkhMCmr5QjEprOXPVZDiOu6RiBSJ+bJ60tGYcjAhjkLiYkoFOVGjnGlrI+X99gDqhCAEUKcHVXx0Lh+VwyUU84EdqbGu9AgnYi+kZSMzQ6WdxDcVEF9LyteSyI1JKUeK6/iCWlJQSx7W8D8/4ObFsb06UCGudIgr1cMGw0XQKQhCAFWrZ0lcIR+dy0Xl8jUSMtaVmt+FPjiCNm95hKdiApo4mVBOJrJH1OKutFhKLhTyx3O4iac5KU9CXKkQV+pFJ9oEXQhoP7yLAVTkkYiPzuWjc/lzBVxXKyrUjf4tmO7VgdKK4ZVQQMa4UGNcGEJIagUflct/e4ebcY59wZ4KdaVDXCkvC634PQDkQBACKBHLkz+L+agc7kQen1vDj3ahX/Kkdgw16KDN59s6W1GdrajlPegqCTmdx53I4z+/xVoYUCGuVIgrPdRRiUd0AZQBQQigeKV15GQeF53Lx+RxbmbUWFfqhxeY/nYUo1ujJgsDEuFBR3gQnjA3S/noXH7NdfZeBR/sRIe4UmNdKGdT3fqFQUchCAEUgydEFgZRuVzK32Gwob/AyUT3w4AiJMCWCrClVvWmS+vIqTwuKpd/7yrrakaFuFKhrrTufQgAXYIgBHheCcX8zjQuOpczN6BCXalP+jD6fHjQ1pjM8KZneBOWZ2SHhRddZmWHhV/1pUc5Iw9B4yAIAdovvohfe4NNryJLu9Hv9xRgwkhTDEUG2VOD7Jn/9SN5Ij4ql1/2J2tmQNb0ZkJcsaNAgyAIAdojvohfd4u9XUb+050+1oU2YtRdkGZzMaXmd6be8KOjcrnV19kPEsl/e9OTPWjkIWgCBCFA21wo4D++weaIyH970YdH0Ab6egi0HWiKhLnRoa70Hw+5T25yX9zi1vSmw9wRh6BmCEKA1rpUyH90g82sJu/1pF/zpQWIwHahKRLhQU/yoI/ncGtvcKuuc+/0oGd608hDUBcEIcCznX7Er73JlojJf3vT071oTIB8fhQhYW70ODc6Mptbe5P76ja3pjcd3glxCGqAIAR4mlN5/NqbbHk9+W8veioiUNEoQsa702HudFQOv/Ym+/FNblUverIH4hBUCkEIIF90Lv/JTbZaQlb1oqd4omtWIoqQcW7UODdBdC6/9ia79iaHfQ6qhCAEaC4ql//4BlvHEoxOVCzElQpxFchG4Wtvcqt60y97YhQOSocgBPjHmUf8f6+xtVLyrj89wxtdsHqMdqFGuwguFfIf3mBXX+NWYmoSKBmCEIDwhBzL5tbe5FiOfBhAT+yEQaD6DXagYkMEsfn8xzfYr25z/+1FT/dCHIJSIAhBr/GEyK5poymyujc9Ade0aZgXnagXnQTnCvi1N9hPbnL/7UXP9EYcgoIhCEF/HcvmVl/nDGmytg89zg0RqLmCHKmgUMHFQv7jG+wnN7k1AfQrPvh7gcIgCEEf1UjIkgT2chH/9QAmFOteaomhDlRsiOBSIb/sCnswi9s1VNDRWN01gU7AIQbQO9dK+YAjUpYj1ycKkIJaZ7ADdTlM0MeW6nlYEpXLq7sc0AUIQtAjPCHf3eFCTko/DqD/bzhjZqDugqBdBDT5KID5NUjw5iX2rQS2nlV3QaDlEISgL/JE/ItR0iPZ3PVwwTQvvPO13nBHKmmSIE9EBhyTplRgaAjth+4A9MIfD7m+R6RDHakzIQJXUxwO1RE2RuTQCOY/3emhx6Xf3eEQhtA+mCwDOk4sJe8lssdy+EMjBIPsEYE66BUfur8dNf0ce76A3zGU6WCk7oJA2yhxRFhXV3f58uW0tLSntCkvL8/JyWm6heO4lJSUhISE8vLyxo21tbWZTdTW1iqraNAtNx/zAUekFQ0keRJSUJf5WVIJ4wXeFqTXYenZfIwMoW2UFYSpqak+Pj7vv//+6NGjZ82axfPN35pnz5719vbu0KHD0KFDGzdmZGQ4OzuHh4evXLmyU6dOW7ZskW2PjY3t1q3byL9dvXpVSWWDzuAJ+eo2N+akdE1ves8wxhzzYnSdIU2+7M/sGsq8coF9L5GVcOouCLSHsoJw9erV06ZNu3DhQnJy8sWLF2NjY5s18PLy+u23344cOdJ0o7m5+cmTJ1NTUy9evHj48OElS5ZUVVXJvtWnT5+Mvw0fPlxJZYNuKBaTcaekv2dxl8MwL0a/jHSmksIFKRVk4DFpeiWGhtAqSukjGhoajh49OmfOHEKIhYVFeHj4wYMHm7Vxd3fv27evoaFh0412dnY9e/aUfd2nTx+pVFpRUSF7yHFcRkZGWVmZMgoGXXIqj+9zRNrDhoobJ/CywOFQvWNrTI6OZJZ2owdFSr+7g4EhPJtSJssUFhZKJJJOnTrJHnbq1On06dNtfZHNmzcHBga6ubnJHt65c2f8+PG5ubmBgYH79u3r2LGj3GfV1NSkpqY2/jiGYYYOHUrT8vOe4ziOw/+Jiqhgb4ul5N1E7mQe2TecesGeIkR//7zc39RdiNrM9CK9beiZF7hLhdyWQbS1kmfQ6PneVrE27e2W+v+mlBKEtbW1FEU1jvaMjY1FIlGbXiEyMvK77747f/687OGwYcNKSkqMjIxqamomT578zjvv7NmzR+4T8/PzMzMzU1JSZA+FQmGXLl3Mzc3lNhaLxQzDtKkwaDexWMzzfGvelO1zp4J6LcGghxV/aZTUwoDX8wlVEomEZVk975o9jEjsi2TNLUHAEWbbAMmgjkrcG+hMVKlNe9vY2FggeEbSKSUIHR0deZ4vKyuTjdtKS0sdHR1b//SYmJjXX389MjKyS5cusi0WFhayL8zMzBYtWvTWW2+19FxfX9/Q0NAlS5a05gfxPG9mZtb6wuB5UBQlFAqVEYQ8IdtTuTXX2bV9mHmdaUIwff6vIDQ21ve1OM0I2TyUnHnEz75IRXSiNvRnDJTzSQydiSopfG8r5U1haWnp5+d38eJF2cO4uLjAwMBWPvfSpUuzZs36/fffW3pKdnZ2hw4dFFMoaL8iMQk9Jd2VzsWHCeZ1xrwYkGOEM3V9ouB+FT/0uDSzGjNooDllXVC/bNmyd955hxBy9+7dmzdv7tu3jxCSnp4eGBiYnZ1taWlZVFT0008/paenV1ZWfvHFF46Ojq+88sr9+/fHjBkTHByckJCQkJBACJkxY4aLi8vnn39uYmLi7u6ekpKybt26zZs3K6ls0C4ncvnX49g5vtRHAQzuUQdPYS8kUaMFG+9wA45KvxrAzPLG2wX+oawgnDdvnqmp6aFDh6ysrOLi4mxsbAgh1tbWCxYsMDIyIoSwLFteXt6xY8f58+eXl5fLxrkMwyxevJgQ0ng1vVQqJYT07dv3jz/+iI+Pt7e3j4qKGjx4sJLKBm1Rz5IPb7C/ZfC/BjHDHDE1FJ6NIuSt7vRIF2raWfZYNr9tMKPsGTSgLagnL3XXakuXLvXx8WnlOcLq6uqW5tGAwolEIkWdI0yp4KefYz3Nqe1DGBv0ZfLgHOFT1LFk5VX2aDa/dzgzxEExn6LQmaiSwvc2jg+AlvnjITfsuHRJN/rQCKQgtIcxQ74byGx6gZ4SK/0xRa8n1oIMFt0GbbIvg1t+hT01VtC7Aw6HwnMJc6MTxlMvRrNiKVnWA0MCvYY/P2iNXx5w717lTo1BCoJidDKn4scLdqVz7yXi3r56DUEI2mFbKvffa9y5EKaHDVIQFMZBSM6GCE7k8shCfYYgBC3wwz3u81tcbAjjY4kUBAWzE5LzoYJz+fyiy6xOTR2EVkMQgqb7XxK38S4XN47BCtqgJNZG5NRYwfVSflE8slAfIQhBo310g/35AXculHExRQqCElkZktgQQVolP/McK8VMUj2DIAQNxROy/Ap7LJu/ME7gZIIUBKUzFZDIUYLSOn7medzXV78gCEET8YS8ncBeKOBPhwg64qJwUBUTAYkcLahnyaQz0jrMntEbCELQOCxP5l5kbz7mz4UKOuCSeVAtQ5rsf5ExYqgIZKHeQBCCZpFyZNZ5NlfEnxgjMDdQdzWglwxp8lsQY2VIjTslFUnVXQ0oH4IQNEgDR6aeYysa+GMjBaZY9QjUR0CTvcOZTubUmJPSKom6qwElQxCCpqhnyZRYtoElf4wQCJGCoG40RbYPYfrYUsFR0sf16q4GlAlBCBqhVkrCYqTGDDk0gjFi1F0NACGEEIqQbwcwY12pEdHSkjp1VwNKgyAE9auWkDEnpS6m1C9BjAHekqBhPunDjHOjgqKkhWJ1lwLKgV4H1KyigYw8Ie1uTe0cyjC4XBA00id9mBne9LDj0lwRVp7RQQhCUKfyejLqhLSvLfXDIIQgaLT3e9JvdaOHHmczqpCFugZBCGpTJCbDoqTBTtT3LyAFQQss7Eq/35Meepy9V4Es1CmYnAfqUVBLRkRLp3jSHwbg0xhojXmdaTMDMuoEe2I07gimOxCEoAZZ1fyIaHZRVxp3BgetM92LpgkZfVIaNRr3iNYR6IZA1dIr+eFR7FvdkYKgraZ60dsGC8aelCYU4xipLkBPBCqVWsGPiGY/DKCXdsN7D7TYODfqp2GC8THSs/nIQq2HzghUJ+kxHxQl/bQv/Zov3nig9ca4UAdfFEw/Jz39CFmo3dAfgYpce0yFxHCbBzGv+OBdBzpimCN18EXBzPPSmAKsh6TF0CWBKmTX8FPjBDuGMOGd8JYDnTLYgYocJVhwRXDzMcaF2gq9EihdHUsmn2H/04ULcVF3KQBKENiR+q6vZNIZFuuRaikEISjdwnjWz4pa5IebnILOCnPhpnpSU89KpZy6S4G2QxCCcn1zh7v5mN82GGdQQMd91pcRMmTFVXzg0z4IQlCis/n8F7fYwyMYE6zcALqOpsje4YLIHH53OkaFWgZBCMqSU8PPPC/dFyzwMMfqG6AXrI1I5CjmvUT2WikmzmgTBCEohVhKws+wK/yZIEekIOiRzlbU5heYl2LZUkyc0R4IQlCKhZfZLlbU293xBgO9E+FBT/eipsRi4ozWQD8Fivf1bS4JE2RAj33ShzERYOKM1kAQgoKdzee/SMYEGdBrmDijXRCEoEiyCTK/YYIM6D1MnNEiCEJQmDqWRJxhV2KCDAAhhJDOVtTWwcxkrDij8RCEoDCvx7FdrKi3MEEG4G8T3enpXlhxRtOhzwLF+Po2d7ec34IJMgD/9ilWnNF4CEJQgLP5/HpMkAGQBxNnNB+CEJ5XNlaQAXgqTJzRcAhCeC5iKYk4w77XExNkAJ4GE2c0GYIQnotsBZml3fBGAniGie70DG9q0hmpBIdINQz6L2i/r29zt8uwggxAa33ShzE3wMQZjYMghHaKzefXJ7O/v8gIMUEGoHVoivwWLDiVh4kzmgVBCO2RXcPPwgoyAG1nYUAOj2DeS2QTSzBxRlMgCKHNGifIDMcEGYC262xFbRvMvBTLFovVXQoQQhCE0A6YIAPwnCa40zO9qYhYaQMOkWoA9GXQNl/d5pIe81sxQQbg+aztw1hi4oxmUFYQ5uXlRUREeHl5jR07Ni0t7ckG8fHx7777bkhIyPr165tuT05OHjFihJeX1/Tp04uLi2UbeZ7//PPPu3bt2rNnz507dyqpZnim2Hz+q9vs0ZFYQQbgedEU+TlIEJ3L732AUaGaKSsIZ8yY4eTkdPHixYEDB06YMIHjmv+lb926ZWhoyDDMvXv3GjdKJJLQ0NCxY8devHjR2Nh47ty5su179+7duXPn77//vnXr1vfee+/ixYtKKhueQjZB5pcggZsZTg0CKICVITk2knnnCibOqJlSgvDevXuJiYnr1q1zdnZetWpVZWXl+fPnm7VZuHDhZ5991qNHj6Ybo6KijIyMli9f7uzsvH79+lOnTuXm5hJCtmzZsmLFim7dug0YMGDevHlbt25VRtnwFFhBBkAZMHFGEygrCP38/ExNTQkhNE337t37zp07rXni3bt3+/TpI/va1tbW1dVVNl68e/duQECAbHtAQMDdu3eVUTY8BSbIACgJJs6onVJO9ZSWllpYWDQ+tLKyKikpafcTGxoaqqqqGrc//dWSkpJ27NixZs0a2UNTU9MrV65YWlrKbSwSiSgK45tn25jK3ChhTr/YUFPT/hepra1lWZamEaWqIJFIWJaVSqXqLkRfPGdnssKPXCsSLIuXruuNP9mztWlvGxsbCwTPSDqlBKGVlZVIJGp8WFVVZWNj08onFhUVNXuioaGhqalp4wtWV1c/5dX8/f1DQ0PnzZsne8gwTNNkbYbneTMzs9YUps9i8/mNadKE8QI7c8PneR2KooRCIYJQNWRBaGxsrO5C9MXzdyYHRpEBR6W/5xvO8cX/yDMovOtWyh739PR88OCBRCKRPUxNTfXy8mrlE1NSUmRfi0SivLw8T09P2fbU1NTGV5NtlIumaRMTE+u/PSUFoTXya/mZ57CCDIDSWRiQ319k3ktk75Rj4oyqKSUI+/Xr5+TkJJvScvjw4erq6jFjxhBC4uPjv/zyS1mb+vr68vJysVjc+AUhZOLEiZmZmTExMYSQjRs3+vv7d+7cmRDyyiuvbNq0SSwWP378eOfOnbNmzVJG2fCktxK4eZ1prCADoALdrKnP+jLz4lgOUahaSglCiqL27t27adMmW1vbt99+e9++fYaGhoSQjIwMWcgRQnbv3u3l5bVnz55Tp055eXl98cUXhBBzc/Nffvll9uzZtra2+/bt2717t6zxkiVLPD09HR0dvby8wsLCIiIilFE2NHMil08u49/vhWvnAVRkrh8tFJCtqZg2o1IUzyvxs4dIJJLNHVXIE+vr6xmGefppz6VLl/r4+CxZsqQ1P6W6utrc3Lwd5emDWinpfki6YwgT7KSY4aBIJMI5QpXBOUIVU2Bnkl7JDzkuvTFR4GyKIzHyKbzrVm6v1L4UbOmJRkZGz5z8A4ry32tskCOlqBQEgFbytaQWdKbf/hODQtV5WhA2NDTk5+cXFxc/uS4M6LZbZfxvGdwXgTgoCqAGH/Ri7pbzx7LR8aqInCC8du3a0qVLu3XrZmRk5OzsbG9vb2RkNHDgwA8//PDhw4cqrxBUjeXJ3Ivshv6MLY6rAaiDEUO2DGaWJnA1EnWXoh/+daQxLi5uxYoVf/75p7Oz88CBAydNmmRjYyOVSouLi+/cubNx48ZPP/00PDx83bp13t7e6qoYlO27O5yFAZnujZN5AGoz1IEKdqLWXGe/HoADM0r3TxAePXp06tSps2fP/uabbwYMGPBkU4lEEhMTs3379u7du6enp7u5uamwTlCR7Bp+3S02YbwA5wYB1OvL/kyPQ5IZ3nQfW/w7Ktc/Qejn53f//n0XF5eWmhoYGISGhoaGht68eRMLsuiqpQncsh6MlwX+8QDUrIMRWR/IvHaRvTZRYIADNMr0z97t3LnzU1Kwqd69e7dyyTTQLr9ncZlV/PIe+J8D0AgzvWl7Idl0F7NmlOtpXV5lZWViYmJ0dLTKqgE1qpKQ5X9yWwcz+OwJoDm2DWbW3WKzqrHYjBLJ7/MaGhoWLlzYsWPHwMDAxgWsFyxYMHXqVBXWBiq18iob5k69YI+DogAapJM59Z8ezKLLrLoL0WXyg3DFihV79+793//+t3nz5saN4eHhx44dq6urU1VtoDpXivnIHP6zvpifBqBx3u1BF9SS37NwgFRZ5ARhQ0PD9u3bv/rqq3feeadr166N23v06CEWi2W3jAddIuXI/EvsdwNoq+e6zxIAKIWAJlsHM28ncOX16i5FR8kJwpKSktra2iFDhjTbLpspWlFRoYq6QIXWJ3NOpiTCA+cGATRUYEdqgjv1wTUcIFUKOX2flZUVwzBPLiJz48YNQkgrZ5aCtnhQxX97h90yCAdFATTa5/2Y4zn85SLMmlE8OUFoamo6YsSI1atXl5SUUNRfUycKCgqWL18+YMAAR0dH1VYIyrUwnn2/F+NmhjkyABrN0pB8M4Ced4ltwLlCRZN/M4eNGzcOGTLE19fXz8+voqJi3LhxcXFxhJDz58+rtDpQsr0PuGIxWdIVB0UBtMBkD3rvA35DMvdBL/zPKpL8venr65uUlPTqq69WV1cbGhreu3cvIiLi+vXrvXv3VnF9oDxl9WTlVXbnUEaA/ykALbFpIP31bTatEgdIFUnOiFAkEn3yySdTpkz59ttvVV8QqMzyK+w0LyxjCKBN3Myo//ZmFlxiz4ZiQWCFkTMWqK6u/uKLLyQS3P9Dl10o4M/l8x8HYI4MgJZZ2o2ulpCfH+BUocLICUI7Ozt7e/vMzEzVVwOqUc+SN+PZjQNpMwN1lwIAbcRQZOdQ5t0rbAlWN1EQOUFI0/T69etXr159584d1RcEKvBZEtvdmhrvjnODAFqppw013YtecRWXFSqG/Fmjx48fr66u7tmzp4eHh6urq0DwT7PTp0+rqjZQinsV/NZULikcg0EALba2D9P9kPRcAR/kiHOFz0t+EBJC/P39VVkHqAZPyOJ49uMAxtFE3aUAwHMwMyA/DmLmX2KTJwmMca7/+cgPwgMHDqi4DlCN7alcLUvmdcZBUQCtN9aV2pVO/S+JXdsHSfhc0CHqkSIxWXWN3TqYoXEoBUAnbHqB2ZrK3avAZYXPRf6IMDIysqXbLb300kvKrAeU6K0E9o3OdE8bxCCAjnAQko8CmPmX2IvjcFlh+8kPwjfeeKOoqEjut3geHz200sk8/nopv3toi2eFAUAbze9M//KA25HKvYFTHu0lv1u8cuUKy/4zMbeqqurChQtfffXVpk2bVFUYKFKtlCyKZ7cMZoTIQQDdQlNky2AmKEoa6kY5mWBY2B7y+0V3d/dmW3r16mVjY/Puu++GhYXRND53aJmPbrCD7KmRzvgnAdBB3a2peZ3p5Ve4fUGYNdMebYi0kSNH3r9/PyUlRXnVgDIkPeb/7z731QD8hwDorFW9mGsl/IlcnLpqjzYEYVpaGiHEyMhIacWA4nE8WRjPruvHdDRWdykAoDRCAdk+hFl4ma3BKtFt16pZo1KpNCMjY/PmzZ06dfL09FRVbaAA39/jDBnyqi+OZgPouOGO1GB7au1Ndn0gDv+0TWtnjdI0HRwc/O233+IEoRYpqCWfJbHncbsWAP3w7UCm+0HJNC+6dwf807dBq2aNCgQCBwcHQ0NDVVUFirEwnl3clelihX8JAL3QwYj8rx8z/xKbMF7A4P++1eQHoZWVlampadO1tgkhUqm0urra2tpaJYXB84rK5VMq+N+CcZAEQI/M9qV/ecBtvsct6Yajd60lf0/5+fldvXq12cbExEQbGxvllwQKUC0hb15ifxzMGCEHAfQJRciPg5hPbrJ5Iswgba02fGSQSqUGBrh3j3ZYe4Md7ULh/iwAesjHklrcjVl+Bbewb61/Hfysq6sTi8WEEI7jqqury8vLG78lEokOHTrk7Oys6gKh7YrEZFc6lzwJq8gA6KkV/rT3AWlyGe+PtYVb4V995Y8//rhs2TLZ12PGjHmy9RdffKGKouD5rE9mZ/nQzqb4BwDQU8YM+U93+tOb3IEXcXbk2f4VhCNGjNi6dSsh5N13312wYIGXl1fjt8zNzbt3796jRw9VFwhtVFpH9qRztzAcBNBvb3ahNyRLbpfRPTAofJZ/dZc9evSQRV1DQ0N4eDgOhGqjdbfY6d4YDgLoOxMBWdaD+SyJw9TxZ5I/bli8eLGK6wCFKK0je+5zN8MxHAQAsrAr7bUfg8Jna7HHzM3NPXbsWGZmZk1NTdPtsmOnoJnWJ7NTPWkXDAcBgBBTAflPd+Z/t3BXimeQH4SxsbETJkyQSqWGhoYMw3AcV1VVJRQKHR0dVVwftF5pHdmdjuEgAPxjcTfaa7/kTjnd3Rqfj1sk/zrCt99+u3///qWlpREREYsWLaqsrIyLi3N2dl61apWK64PW+zKZnYLhIAA0YSogb3VnPk/CNYVPIycI6+vrU1JS1qxZY2ZmRgiRSqWEkMGDB2/btu2tt96SXWgImuZxPdmRxq3wx6JKAPAvS7rSZ/O51AosNNMiOf1mRUUFy7Kyo6BWVlaNl9UPGDCgurpadldC0DQbktkpnrS7GYaDAPAvZgZkSTfmUwwKWyYnCDt27GhkZPTo0SNCiIeHx8WLF2V3okhKSiKEyIaJoFEe15PtqdxKDAcBQJ6l3egzj7i0SgwK5ZPTddI0PXTo0OjoaELI1KlTs7OzBw8ePH/+/PHjx/v7++PGvBro69vsZA+6kzmGgwAgh5kBWdyV+ewmBoXyyZ9h+MMPP0rym1IAACAASURBVMiOiNrZ2R07duzzzz8/ffr00KFDN2zY0Pob816/fj06OtrW1nb69OmWlpZPNiguLv71119ra2vDw8O7dOlCCKmsrNy/f3/TNoMHD+7atWt2dvapU6caN44dO9bV1bWVZei8snqyNYVLnIjJogDQoqXdae/9krRK2s8Sn5ibk5NqPM/b2tp27dpV9nDEiBGxsbGZmZmHDh3y8PBo5etGRUWNHDmS47hz584NHDjwySk2paWlAQEBycnJtbW1AwcOlN31qa6u7vrfLl++PH/+/IKCAkJIcnLyhx9+2PitysrK9v/GOufr2+wkD9oDw0EAaJmFAVncDdNH5ZMzjCgoKHB2do6JiRk5cmS7X/fTTz9dt27dvHnzeJ4PDAw8cODAq6++2rTB9u3be/XqtWvXLkKIkZHRunXrDh8+bG9v33jB/oEDBy5cuBAUFCR76OXlhWv5n1TRQLalclcmYDgIAM/wdnfaa78kvZL2xaDw3+SMCK2trQUCAcO0fyUCkUj0559/hoSEEEIoiho7duyZM2eatYmNjR07dqzs67Fjx8bGxjZrsHPnzjlz5jQeiS0vL9++ffvvv//++PHjdheme76+zY53x3AQAJ7NwoAs7sqsu4VBYXNyRhJCoXDKlCk//fRTcHBw+15UdjzTzs5O9tDBweHy5ctPtrG3t29sUFVVVVNT0zglNS8v79y5c9u2bZM9NDIycnd3T05OTk1NffPNN2NiYgICAuT+6Ly8vLt37+bk5Pz16wkEK1asMDExkdu4vr7e0NCwfb+jJqhsIJvvURfG8PX1UnXX8mz19fU0Tbf+HDM8D4lEwrIsReETkopoS2fypg/peoS6WyL1tlB3Kc+hTXu7NeM6+YfUQkNDly9fPnTo0PHjxzs7OwsE/zR76aWXnvmDZZ0dz/81VZdl2SfroGma4/76YCL7ommbnTt3BgUFubu7yx6OGjVq1KhRsq+XL1++cuXK06dPy/3RhoaGJiYm1tbWsocmJiaGhoYtdb7a3i9/n0bGuRJfK+3o7Oi/qbsQvUDTNM/z2Nsqoy3vbWtj8mZnsuEutW2Qukt5Dm3a2635OCg/CJctW1ZUVFRYWBgXF9fsW43x9hSOjo4URRUUFHTq1IkQUlBQ8OQipU5OToWFhbKvCwoKrK2thUJh44/Yu3fvZ599JvfFR4wYcfDgwZZ+tJ2dnY+Pz5IlS55ZJCHEwMDAwMCgNS01UGUD+TFVEh8mMDDQjiCU7W2t6Cx0A03T2vv21jpa1Jm805P4HJDk1gk8tfaUisL3tvxe6cqVKxktaM2LCoXC4cOHHzlyhBAilUqPHz8uO18omxQquzx/7NixR48elcXqkSNHGs8XEkJiY2PLy8snTJjQuEX2FJkzZ874+fm153fVLd/d5UJdcdIbANrG0pAs6EL/D9NHm5A/Imw8JtluH374YXh4+P3791NTU01MTCZOnEgIycrK6tu3b3l5uZWV1Zw5c7Zs2TJhwgR7e/vDhw9fuHCh8bm7du2aOXOmsbFx45Y5c+ZUVla6ubmlpKQkJyc3vaZQP1VJyPd32UthmCwKAG22rAfje0CS1Qvz7P5CtXSok2XZc+fO3blzp7a29oMPPiCEpKSkmJmZtf5K9qysrNjYWGtr63HjxhkZGRFCRCJRfHx8cHCw7KRjTU3N8ePHa2trx44d2/TYaVxcnJ+fX+NcG0JIQUFBfHz848eP7e3tg4ODLSxaPM+7dOnS1h8ara6uNjc3b+Wvo1E+uck9qOL3DNOme4yJRCKhUIhDo6ohmyzT9NMkKJXWdSarrrHFdWTbYG3qQxopfG/LD8KioqKQkJAbN24IBAJ7e/u8vDxCyH/+85/ExMRLly4p8McrnD4EYZWEeO+XxIUJtGuFCAShKiEIVUzrOpOyeuJ7QJI4UaCNg0KF7235vdIbb7xRXl5++fLlmJiYxo0vv/xyQkJCVVWVAn88tMOmu9wYF6yTBADtZ2NE5nehv8A1hYQQuUEoEomio6O/++67gQMHNv387uPjw3Fcbm6uCsuD5mokZNNd9r+9Ma4CgOeyrAdzMIt7WI1bUsgLwsrKSpZlvb29m22XXe1XX1+virqgBRvvciOcMRwEgOfVwYi80Zn+IhmDQnlBaGtra2Jicu3atWbbY2NjGYbx8vJSSWEgh0hKNt1lV/XCcBAAFOAdf+ZAJpddo++DQjldqqGh4dSpU1euXHnhwoXGa/JPnz69bNmy8PBwuTdUAtXYdJcLdqI7a8lSMgCg4WSDwvV6PyiUfyHaN998k5aWNnz4cFNT04aGhg4dOpSVlXXv3v2HH35QcX3QqEZCvr3DngvFtYMAoDDLezBdfpe835N2MdXfT9jye1ULC4vz588fPnw4JiamsLDQ2tp62LBhzS5yBxXbnMINd6S7YDgIAIrT0Zi85kevu8V9/4JWXlOoEC0OLwQCwZQpU6ZMmaLKaqAlIin5+jZ7eiyGgwCgYCv8mc6/S97T40Hh0zrWlJSUW7duPXr0yN7evnv37r169VJZWdDM5nvcMEe6h42evk0BQHlsjclsX3p9MrdxoJ4OCuUHYXV19ezZsw8fPtx0Y1BQ0L59+xpvIggqUyslX99mT2E4CADK8a4/0/WgZKU/7ayXg0L5E/Hnzp0bHR396aef3rt3r6ysLD09/bvvvktKSoqIiFBxfUAI+TGFG+JA+2M4CADKYS8kr/rQX97W0+mjcgYZ1dXVhw8f3rx587x582RbrK2tfXx8PD09w8LCHjx48OS19qA8dSz55g4XPVpPD1kAgGqs6Ml0OyhZ4U87mejdZ245I0KxWMyy7LBhw5ptHz58OCGkurpaBWVBo833uAF2FIaDAKBUDkIyy4feoJfXFMoJwo4dO3p7eyckJDTbnpCQYGlp2blzZ5UUBoQQUseSr+9wq7GyKAAo33s9mT33ufxavVtoRs6hUYqidu3aNWPGjOrq6vDwcEdHx+Li4pMnT65du/ann34SCoWqr1JvbUnhAjtSPTEcBADlcxCSmd7017e5Df3161yM/KHGSy+9lJubu3TpUldXV4FA4OTk9Nprrz18+DA8PJz6244dO1Rcq74RS8mXydwaDAcBQFVW9qR3p3NFYnXXoVryZ+SvWrVKJBI9/ZmBgYFKqAf+sT2N69eR6tUBw0EAUBEnE2q6F70hmf1SnwaF8oNw8eLFKq4DmqlnyZfJ3B8j9ei9CACa4INeTI9Dknf9GTu9OQ+Gw24aalsqF2BL9bXFcBAAVMrRhLzsRX91m1V3IarT4mIliYmJx44dy87Orqura7r9wIEDyq9K39WzZH0yd3gEhoMAoAbv9aR7HZYu76Evg0L5Qfjll1+uWLHC1NTUzc3NxMRExTXBjjSuVwfSryOGgwCgBq6m1Mue9Dd32M/76cXHcTlBKJVKP/zww9dee+3777/HxRKq18CRL5O5Ay/qxfsPADTT+73oXoely3owHfXg5ntyzhGWlZWJxeKFCxciBdViZxrXzZoEYjgIAOrjakq95EF/e0cvzhTKX1nG3d09MzNT9dWAhCPrk7k1vTEcBAA1+6AXvSWFK6l7dkttJycIKYravHnz6tWrr1y5ovqC9NyudK6zJelvh+EgAKiZmxkV4UF/pweDQvmTZYKDg/v16zdgwAALCwtbW9um38rIyFBJYfpIypF1t7h9QRgOAoBG+KAn3feI9B1/xspQ3aUok/wgnDNnzv79+4cOHert7S0Q4H6wKnIij3M2IQMwHAQAzdDJnBrpQv/ygFvUVZcvOpcTclVVVfv379+wYcOyZctUX5A+25HGz/XT5XcbAGidub708iusbgehnN9NIpHwPC+7+yCoTKGYxBVyL3no8rsNALTOi86USEqul+ryvZnkdLsdOnQIDAyMj49XfTX6bFcaN8WDNjNQdx0AAE1QhMz2oXem6fINe+Wf//vss8/mzp3b0NAwatQoU1PTpt/y9PRUSWH6hSdkdzr3K6bJAIDmmeNL+R9mN/RnTHR0xoj8X2vmzJlFRUXvvPPOk9/ieV0eIKvLuXxeKMCaagCgiZxNqRfsqd+zuFd9dPPcjfwg3L59e7O1tkGpdqZxb2CaDABoqrl+9Ne39SwIw8LCVFyHPqtoICfzuE0v4PQgAGioca70ongupYLvYqWDB66eFu81NTU3b96Mi4tTWTX6ae99bqwrbWOk7joAAFogoMkrPtRP6bo5ZUZ+EEokkrfeeqtDhw4BAQHTpk2TbZw3b9706dNVWJu+2JXO4fJBANBwb3Sm/+8+J9HFKJTf/65cuXLnzp0fffTRpk2bGjdOmjTpyJEjOHeoWIklfGUDGeagg0cbAECXeJpTflZUZI4OJqGcIGxoaNi2bduXX375/vvv9+jRo3G7v7+/WCzOzc1VYXm6b2ca97ofTSMHAUDjzfXTzQsK5QRhaWmpSCR6cmUZMzMzQkhFRYUKytITIin5PYub7YsYBAAtMLkTfaWYzxXp2kV0coLQ0tKSYZicnJxm25OSkgghzs7OqqhLPxzI5AY70E4mCEIA0AJCAZnqRe9O14MgNDU1DQ4OXr169ePHjynqrz5adn19YGCgk5OTaivUZTvTuLkYDgKA9pjXmd6RyrG6FYXyryPcuHHjkCFDfH19u3TpUlFRER4efv78eZZlz507p+L6dFhaJZ9ZzY91xXxRANAa/jaUnZCczedHOuvOh3j5vXDnzp2TkpKmTZtWVFRkYGBw/fr1sLCwxMTEPn36qLg+HbY9lZvjSxsgBwFAq+jelJkWl1B1dnb+/vvvVVmKXmngyM8PuEthOrqELQDorule9AeJkpI6pqOxuktREPnjkSlTpty7d6/Zxnv37o0cOVL5JemFY9lcN2vK20J3ji0AgJ6wNCTj3emfH+jOoFB+EF68ePHJyyQqKipiY2OVX5Je2JmG1WQAQFvN9aO3pep6EMqVnZ1ta2urvFL0R56ITyzhw90RhACglYY6UBQhfxbryOTRf52jOnTo0JYtWwgh5eXlb7/9tqWlZeO36urqrl+/HhIS0vqXlkgkubm5Tk5OxsYtHkjOy8szMTGxsbFpzQvm5+cbGBh07Nix9TVoph1p3HRvWojzgwCgtWb70jvTuAF2unA78dYOSqytrZctWyaLydY4d+6cm5tbaGioi4vL4cOHn2xQWlrav3//wYMH+/r6Ll68WHa/36KiIqqJjz76SNa4pqZmxIgRgYGB3bt3nzFjhlQqbWUZGojjyU/p/Gu+GA4CgBZ71Yc+9JCrlqi7DkX416gkIiIiIiKCEDJq1KgNGzb4+/u370U5jps7d+6GDRtmzJhx9uzZyZMnjxkzxsTEpGmbTz75xN3d/c8//ywvLw8ICBg3btyYMWMIIQYGBg0NDc1e8JtvvqEoKjs7u76+vn///vv27Zs1a1b7alO7mEe8nZD06oBpMgCgxeyFZLgjvT+Te137pzvI/wViYmLanYKEkISEhKqqqqlTpxJCgoODHRwcoqOjm7X59ddfFy1aRFGUjY3NtGnT9u3b1/it+vr6ZmO+X3/9dcGCBQzDmJiYzJkzp2ljrYNpMgCgG173o3foxAWF//TIDx8+bOWC2llZWZWVlU9p8PDhQ09PT4b569ixj4/Pw4cPmzaoqakpLS318fGRPfT29s7KypJ9LZFIHBwczM3NR44cmZmZKduYnZ3t7e39ZOMnNTQ05ObmXv9bSkpKa34jlSmpI7H53DQvBCEAaL3RLlS+iNwp1/opM/8cGr1y5cqbb765ePHi2bNne3p6ym199erV7du37927NyUlpelUmmaqq6uFQmHjQxMTk6qqqqYNZA8b25iamsq2WFpa3r9/39vbu6amZuHChVOnTr169apEIhGLxU82luv+/fvp6emnT5/+69cTCP7444+WSq2pqWnpdZRke5og1Imi6uqq9e+ujrW1tVKplKbxIUAVJBIJy7ISiU6cwNEGqu9MNMT0ToIfb0vW9VbpvI027W1jY2MDA4Ont/knCF9++WUDA4P33nvv008/7dmz58CBA/38/GxsbKRSaUlJyZ07d+Lj4zMzM4cNG3bp0iUPD4+nvKi9vX15eXnjw/Lycnt7+6YNOnbsSNN0RUWFtbU1IaSsrEzWwNjYWDbyMzMzW7dunbOzc1lZmY2NjbW1deNotbGxXN26dZs4ceKSJUue/ms3Mjc3b2VLhfjloXTbEMbcXFfWY2gLmqaFQiGCUDVkQfiUCdugcCruTDTEgu58v6PSrwYJjVU7e1Sxe/tfk2UmTZo0ceLEkydP7tmz59ChQ8XFxY3f8vLyGjly5BtvvNGa5Ub9/f3T09MrKystLS1Zlr127drq1aubNjAwMOjateuVK1dkgXrlypVevXo1e5GysjKGYWQDwZ49e165ciUwMLClxlohvohnefKCPabJAICO6GRO9e5AHc3mXvbU4s+4za9lo2k6JCQkJCSE5/nCwsKSkhKBQODg4NDKS/1kvLy8goODFy5c+O677+7atcvd3X3QoEGEkAMHDhw4cODgwYOEkMWLF69Zs8bd3T0rK+vIkSPXr18nhERFRRUXF3fp0qWkpGT16tXTpk2TBeGSJUuWLFnSo0ePqqqqnTt3xsTEKGwHqNDONO71zrgXPQDoFNka3DoVhI0oinJ0dHR0dGzf6/7888+rVq168803u3TpcuzYMdnGDh06+Pn5yb6eN29eXV3dypUrzc3Njxw5IjsiamNjs3fv3u3bt1tZWb3yyisLFy6UNZ40aVJFRcXHH39saGi4Z8+efv36ta8qNaqRkKPZ3Lp+zzhUDQCgXcI70UsT2Iwq3ktrF0+mZFey64ylS5f6+Pi08hxhdXW1yg7rb03lzjzif39RF1ZhaB+RSIRzhCqDc4QqpsrORAMt+5M1MyBr+6iof1P43kavpCK4fBAAdNW8zvTudF57b1uPrlkVbpfxhbVEl27oDADQqLMV5WJKTuZpaxIiCFVhRxr3mh/FIAcBQEdp9W3rEYRKV8+S3zK5OVhlGwB011RP+nwBV1Cr7jraBb2z0h16yPXuQLmbYTwIADrLzIBEdKL/775WDgrlXz4xZcqUpkvDNNW4ehm00s40bkEXfOAAAB33uh/9ygV2RU/tu1q6VR10ZWVlQkLC1atXlV2N7smq5m+X8ePdEIQAoOP621FChsQVat+UGfkjwgMHDjTbUlJSMnny5FGjRim/JJ2yI417xYc20t+rBwFAj8huWz/UQcu6vNaOVDp27Pj111+vWbPm6TdggqakHNlzn38Nlw8CgH6Y5UMfy+bK69VdRxu1oY92cnKqq6vLyMhQXjU65kQe18mMdLXSugPmAADt0cGIjHah92Vo2ZSZ1gYhz/Pbt28nhHTq1EmJ5eiWnWk8VpMBAL0y14/emqplQdiqWaNSqTQrKys7O3vmzJltug2FPisUk4uF3M/Dsco2AOiREc6USEpulPIBtlpzMKzFu080ZWRk9OKLL44ZMyYiIkLZBemM3encSx60GXIQAPQJRcirPvTOdC7AVmumzLR21ii0CU/I7nTu5+Fa8z4AAFCU13wp/8Psl4GMSauGWuqHM1hKcb6AN2ZIYEetOTIAAKAozqbUADvq0EOtOVOIIFSKnWnc65gmAwD6SrvW4EZnrXgVDSQqh5vuhX0LAHpqvBudXsmnVmjHKjPorBXv5wfcWFfaFvcGBwB9JaDJLG96d7p2DAoRhIqH46IAAHP96P+7z0m0IQrRXyvYtVK+soEMd8Q0GQDQa76WlK8lFZWrBUmIIFQw2XBQ+25DAgCgaNoyZQZBqEhiKTmYxc32RQwCAJCXPOiEIj5XpOlTZhCEirQ/kxtoRzuZIAgBAIhQQKZ40nvSEYT6ZGc6N9cPKQgA8Je5fvSudI7T7ChEECpMWiX/oJIPccUuBQD4Sx9bysqQnM3X6CREr60wO1K5Ob60AfYoAEATc/3onZp9QSG6bcWQcuSXDG6OL/YnAMC/zPSmT+ZypXXqrqNl6LgV42g219mS8rHECUIAgH+xNCRhbvQvDzR3UIggVIyd6RxuRg8AINdcP3qHBl9QiL5bAfJE/JViPrwTdiYAgBzDHCkpT64Ua+iUGfTdCrD3AT/Vi9aWW1ACAKjeqz70nvsaOihEECrA0Wxusgf2JABAi17yoI5m85o5JET3/byKxCS9kh9sj2kyAAAt8rKgLAzJjVJNjEIE4fOKzOHGuuDyQQCAZxjvRkXmaOLRUfTfzysyhw9zw3AQAOAZwtzoY9kYEeocsZRcKODGYFk1AIBnecGeyq/ls2s0LgvRgz+XM/lcX1vKylDddQAAaDyaImNc6KgcBKFuiczhw9ywDwEAWiVMI08TohNvP56Q6Fw+FCcIAQBaZ7QLfbmIr5aou45/QxC2X2IJb21IvC0QhAAArWJmQF6wp2LyNGtQiCBsv8gcbrw7UhAAoA3C3OhIDTtNiCBsv2PZOEEIANA2E9ypqFyO1aQoRD/eTjk1fKGYD+yIESEAQBs4m1JuZlRCkQYlIYKwnY5m82FuNI0cBABoI02bO4ogbKfIHA4LygAAtEOYG31Mk04TIgjbo0pCrpbwI5yx9wAA2izAlqqRkPRKTclCdOXtcSqPG2RPmeIGhAAAbUcRMs6N0py5owjC9jiWzU9wx64DAGin8W605pwmVNaghuO4b775JjIy0tbWdsWKFYGBgU+2OXPmzMaNG0Ui0bRp015//XVCSE1NzY4dOy5cuFBZWdm7d+8VK1bY29sTQpKTk3/44YfGJy5ZsqR79+5KqvyZWJ6czOM+74fxIABAOwU7UdPP8Y/rSQcjdZeivCD89ttvd+/evWXLljt37owePTolJcXBwaFpg5SUlEmTJv34449OTk6vvvqqqanptGnTsrKyEhISZsyYYWdnt3HjxpEjRyYlJdE0nZ2dfeHChU8++UT23A4dOiip7Na4VMh7mFMuppgpAwDQTkYMCXKiT+RyM73Vf3RNWUG4cePG77//fvDgwYMHD46Ojv7pp5/ee++9pg1+/PHHqVOnzpgxgxDy4Ycfbty4cdq0aT169Ni/f7+sQa9evSwtLR8+fOjp6UkIsbW1femll5RUbZtE5nC4jh4A4DmFuVGROfxMb3XXoaRzhOXl5dnZ2QMHDpQ9HDhw4M2bN5u1SUpKGjBgQGODpKQknv/XidOUlBShUNg4jszKypo6deqCBQvOnj2rjJpbD3fiBQB4fuPc6NOPuAYNOFGolBFhYWEhRVFWVlayhzY2NoWFhc3aFBUVWVtby762trauq6urrKxsfEpNTc0bb7zx8ccfm5iYEEKcnZ2XLVvm5eWVmpoaHh6+efNm2VDySXfv3v3jjz92794te8gwzJEjRywtLeU2FolEFNW2SEuvosRSQy8jUU1Nm54HpLa2lmVZmsZgWhUkEgnLslKpVN2F6It2dCYgJMTX3OBkZm2wQ9vCsE1729jYWCB4RtIpJQgtLCx4nq+rqzM1NSWEiESiJ6PI3NxcLBbLvq6traVpWtaYECIWi8PCwvr37//OO+/ItgQEBAQEBMi+NjMz+/bbb1sKQm9v7z59+rz88suyhwYGBs7Ozi3VyfO8mZlZm361M5ncBHfevI3PAkIIRVFCoRBBqBqyIDQ2NlZ3IfqiHZ0JEEImenBnSgTjvZk2PUvhe1spQWhvb29sbJyRkeHv708IycjIcHNza9bGzc0tIyND9nVGRoaTk5OBgQEhpK6ubuLEia6urlu3bpWb+e7u7o8fP27pRxsZGbm6uvbp00dhv8y/ReZwq3u37W8GAAByjXenxp7kNg5UcxlK+XguEAgmT578448/EkKKi4sPHz48bdo0QkhZWdn69evr6+sJIdOmTfvll19qamp4nt+6dausQUNDw5QpUywsLHbt2tV06JCamio7yFNTU/PDDz8MGTJEGWU/0+N6cqeMH+aAAyAAAArQ1YoyYkhymZqvrFfWcarPPvvs4sWLXbt27dat2/Tp0wcNGkQIKSoqWrlypeyI6KRJkwIDA729vX19ffPz82VzSi9fvhwZGXnw4EEDAwOKoiiKunTpEiFk48aNHTp06Nq1q6OjI0VRGzZsUFLZT3c8hxvpTBthQAgAoCDjXKlj2WoOQqrZXE0F4nk+IyPDysrK1ta2cSPLsgzzT5IUFRWJRCLZBRJPV1ZWVlpa6uDgYGFh8ZRmS5cu9fHxWbJkSWsqrK6uNjc3b01Lmcmx7AR3apYGXPWijUQiEc4RqgzOEapYWzsTaHSugH/vKntlQhvO0yl8bytxeRSKory9m18h0jQFCSGyhWNaw8bGxsbGRjGVtUs9S2IfcVsGGaixBgAAHTPUgcqs5h+JeGf1rVKCj+etda6A97ehbPEJGwBAcRiKjHKmT+Sp8+gogrC1InO4MCy0DQCgaGHqvhMFevbWisrhx2NBGQAARRvrSl8o4GrVt/wDgrBVkh7zRgzxtUQQAgAomKUh6WtLxearbbE1BGGrHM3mJ7gjBQEAlCLMjVbj0VEEYavgjhMAAMozsRN1LJvj1BSF6NyfLb+Wf1jND7TDiBAAQCnczaiOxlRiiXqSEEH4bMey+VA3WoBdBQCgNOPdqcgc9ZwmRO/+bJE5HG5ACACgVGo8TYggfAaRlMQX8aNdsKMAAJQosCNVLOazqtWQhejfnyEmjxtgR5ljYTUAAGWiKRLqRh9Xx6AQQfgMkTk85osCAKhAmJt6ThOii38ajifRuVyoK04QAgAo3Shn+moJX9Gg6p+LIHyaP4t5RxOqkzmCEABA6YQCMsSBislT9aAQQfg0kTkc1hcFAFAZtcwdRRA+zbFsHnecAABQmTA3+kQuJ1XtmBC9fIsyqviKBtLHFiNCAAAVcTQhnhbUpSKVDgoRhC06lsOHuVGIQQAAVQpzo1U8dxRB2KLIbNyJFwBA1ca7U0ezMSLUABUN5MZjPtgRA0IAAJXqaUNJOZJSobosTJToewAAGCpJREFURBDKF53LBTnSQoG66wAA0D9hbpQq544iCOWLzOHDcCdeAAB1CHNX6WlCBKEcEo7E5HEhrtg5AABqEORI3S3ni8Qq+nHo6+W4WMj7WVIOQnXXAQCglwxoMsKJPpGrokEhglCOyBzMFwUAUKcwd9WdJkR3L0dUDo878QIAqFGoK302n6tjVfGzEITN3SnnpTzpbo0gBABQGxsj4m9DnctXxaAQQdjcsWx+AuaLAgCom8rmjiIIm4vM4XAnXgAAtZvoTh3L4VUwJESP/y/FYpJWyQ9xwIgQAEDNvC0oMwG5War0KEQQ/ktkDjfahTbEXgEA0ADj3aljyj86ii7/XyIxXxQAQGOo5j69CMJ/iKXkfAE31gX7BABAI7xgT+XW8Hki5WYhOv1/xObzAR0oayN11wEAAIQQQhiKjHVV+qAQQfgPLCgDAKBpwtwoZV9EgX7/LzwhUbk4QQgAoFlGu9DxhXyNRIk/AkH4l+ulvIUB8bZAEAIAaBBzAzLQnop5pMRBIYLwL8eyOSwoAwCggZQ9dxRB+Jdj2TwWlAEA0EDj3anoXI5VWhSi6yeEkJwavkDM97fDiBAAQOO4mlJOJtSfxcpKQgQhIYQcy+FDXWkGOQgAoJHGK3PuKIKQEEIisznMFwUA0Fhh7nRkNkaESlMjIVdK+BHO2BUAABqqjy1VJSH3K5WShej9yYk8bpA9ZW6g7joAAKAFFCGhrtTxXAShckRivigAgMYLc6cjs5VymlDfA4Dlyck8bhxOEAIAaLYRTtTNx/zjesW/sr4HYXwR725GuZgiCAEANJoRQ4Y70qfyFD8oFCj8FRtdv379xo0bXbp0GTx4sNwG+fn5p0+fNjMzCwkJEQqFjdvPnj2bmZkZGBjo7+/fuLG0tPTkyZMGBgYhISHm5uaKKjIyGwttAwBohzB3KjKHD7NT8MsqKwO+/vrrCRMm3Lp1a/bs2StXrnyyQVJSUvfu3S9cuLBly5ZBgwaJxWLZ9gULFixatCgpKWn06NE7duyQbczIyOjWrdvJkyd//fXXgICAsrIyRdV5DHfiBQDQEuNc6VN5XIPCx4S8EtTU1FhaWiYmJvI8n52dbWxs/OjRo2ZtIiIiVq1axfM8y7L9+/ffsWMHz/P37983MTEpKirief7s2bMODg719fU8z8+fP3/BggWyJ44dO/bzzz9v6UcvWbJk48aNrazz+qNq518lXJt/P2iPmpoalmXVXYW+aGhoEIvF6q5Cj1RVVam7BH0x8Kjk2P0axb6mUkaEly5dsrKy6tu3LyHEzc2tV69eMTExzdI3KioqIiKCEELT9MSJE6Oioggh0dHRgwYNsrOzI4QMHz5cIpFcv36dEBIVFTVp0iTZcyMiImSNn1/0I2aCO4XxIACAtghzp0/kKzi5lHKOMD8/39nZufGhs7Pzo0ePmjYoKyurq6trbOPk5JSfn9/siRRFOTo6Pnr0iOO4wsLCJxvLVVJSkpWVVV1dLXtoYmLyxhtvGBoaym184hH1fi9OIlHuLR9BRiKRCAQCmsYZWVWQSCQsyzIMo+5C9IVEIpFIlHnHPPhbiBOZcI9p/d5mGOaZ3Y5SgpBlWarJQIumaZZlmzUghDS2YRhGKpW29ESO4ziOe7KxXPX19bW1tY0nEauqqurr61vqDsKdG4baGf67NFAWlmVlh0bVXYheYP+m7kL0Bfa2yviZk/meDSwrf3jzpNZ8+FZKEDo6OhYXFzc+LCoqGjVqVNMGtra2BgYGxcXFtra2sgZOTk6yJ6akpDR9opOTk0AgsLOzKy4u7tKlCyGksLBQ1lguFxeXoKCgJUuWtKbO+Z0l5ibGbfzloJ1YljU2NsaIUDUYhpHtcHUXoi8kEgn2tsos7qrgva2UXmngwIF5eXkPHjwghJSXlycmJg4fPpwQIhaLZWM1mqaHDx9+8uRJWfuTJ08GBwcTQoKCguLi4kQiESEkKSmprq4uICBAtv3UqVOyxqdOnZI1BgAAeH5KGRHa2NgsXLgwPDx89uzZBw8enDRpkre3NyHkp59+2rp1a1JSEiHk/fffDw8PF4lE2dnZaWlp+/fvJ4T07t172LBhoaGhYWFh27ZtW758uampKSHk3XffHT58OMMwVVVVp0+fvnHjhjLKBgAAPaSs41RffvnlRx99VFZWtnjx4j179sg2BgcHr127VvZ1UFDQ+fPnOY7r2rVrYmKitbW1bPvBgwdnz55dVla2YcOG1atXyzb27t37ypUrRkZGTk5ON27ccHFxUUiRMTExHIeZMipy9erV0tJSdVehL7Kzs+/evavuKvRFbW3thQsX1F2FHmk8mqgolI5NXli6dKmPj08rzxE6ODgkJSU5ODgouyoghEycOHHWrFmya2ZA2b766qvs7OyNGzequxC9cPXq1YULF167dk3dhegLQ0NDkUhkYKCwewZh5gKojo596tJk2NUArYcgBAAAvYYgBAAAvaZr5whDQkIyMzNdXV1b0zguLq5///4trTsDipWcnOzg4CBbPw+ULScnp66uztfXV92F6IWqqqrU1NTAwEB1F6Ivzp49GxQURLVufczw8PCFCxc+vY2uBWFiYmJeXl4r79OUlZXl4eGh7JJAJj8/38bGBhcdq0ZVVVVDQ4NswQpQNpZlHz165Obmpu5C9EWbum4PDw8vL6+nt9G1IAQAAGgTnCMEAAC9hiAEAAC9hiAEAAC9hiAEAAC9ppRFtzVfUVFRSkqKp6dnSxO9KioqTp48yTDMmDFjWjkHFVpSXV198uRJjuPGjBljaWnZ7LsikSghIaHxoZ+fXyuvfoFGcXFxDx486NevX/fu3eU2yMjIiIuLc3Z2fvHFF3EnrOeUk5Nz7tw5e3v7ESNGCATNu9Dc3Ny0tLTGhwMHDpTdOQDaJz8/PzU1tXPnzi3dfa+0tDQmJsbY2HjMmDEmJibt/DG8/hk1apRQKDQ1Nf3mm2/kNsjNzXV2do6IiBg/frynp2dxcbGKK9QlRUVFHh4e48ePj4iIcHFxyc3NbdYgJSXF0NBwxN8OHz6sljq116JFi3x8fObPn29vb79z584nGxw/frxDhw6vv/567969J06cqPoKdcnZs2dtbGzmzp0bGBg4atQo2b2mm9q0aZODg0Pj+/nhw4dqqVM3BAYGmpmZCYXC3bt3y22Qnp5uZ2c3ffr00aNHd+vWrbKysn0/SB+DMCsrSyKRjBkzpqUgXLZs2axZs2Rfh4eHf/TRRyqsTtesWbMmPDxc9vWsWbOWL1/erEFKSoqdnZ3K69IRGRkZQqGwoKCA5/nY2FgHB4eGhoZmbXr27PnTTz/xPF9dXe3o6BgXF6eGQnXFoEGDvv/+e57nxWKxp6fniRMnmjXYtGnT9OnT1VGaDsrMzJRKpQMGDGgpCF977bUlS5bwPM9xXFBQUEtd+jPp40GSTp06PXlAo6njx49PnjxZ9nVERMTx48dVUpduas3OZFk2NjY2Pj5edk9maL3o6OgXXnhBdgeVoKAgiUSSmJjYtEFOTk5ycrLsph9mZmajR4/G+7ndysrK4uPjZTvT2Ng4NDRU7s4sLy8/ceLEzZs3cZe35+Th4cEwzFMaHD9+XPbn+P/27j4oqqoNAPjZD5bvBZZguYCAS8sI8aVAQFhQYg3GgIGpxFROYsUIjGnaFOCAkWEwKoOZlTISSuQSKkQpaIpCTATxDYFIAUIssLvsALuwy+6+f5zpzg4siCLwAs/vr3ufc+69Z8/u3Gf33rP3UCiU+ZyrV2MifKi+vj4rKyu8bGVl1dvbu7TtWdZ6e3sf2plMJjMjIyMmJobL5arfLwQP1dvbS07PSaFQCIKY0sN9fX3GxsYGBgZ4FT7P89HX10en09lsNl7V2JkUCoXP53/11Vdbt2718fERCASL3szVQi6XDw0NkZ//+Xy2V+ZgmV27dtXV1U0Jbtmy5ejRo3PZfHJyknyKHY1Gm5ycfMLtW1kuXbqksWPxW6BQKMjRGRo7k8vldnZ24uXk5OQ9e/Y0NTUtZHtXFIVCof7ERTqdPqWHp1SAz/N84A/z7CeH999/f+/evQghmUwWFBR05MiRjIyMxW7o6qBQKJRK5RM5V6/MRHj48GGJRDIlaGxsPMfNCYIYHBzEy3w+f6bRSgALDAx0cnKaqZQgiIGBAbyssTPVL31EREQkJydPTExoa2svRFNXHoIg1Gein97DFhYWw8PDMpkMP1yez+cTBLHYrVwp8C3Y4eFhfDLR2Jnk55nBYGzbti0vL2+xW7lq6OjomJiYDA4OcjgcNL9z9cq8NMrhcJynIX9BaySXy0UiEV4OCAi4fv06Xi4pKQkICFjoBi9rLBZrem+T4/hn6kyRSCSXy6fs6s8//zQ3N4csOHcBAQHl5eX4a19DQ8PY2JiHhwdCSCKRjI6OIoTWrl1rY2Nz48YN9N+92BdffHFp27x8sdlsJyenkpIShJBKpSotLcWdqVAoNF4Crampgf8CPXH4uwheflLn6pX5i3B258+fr6ysbG5uFovFra2tUVFRXl5eJSUlb775plAoRAjt37/fz89PX19fJpPl5+dPGX0AHklcXJyXl5eJiQmDwcjKyqqoqMBxDoeTm5sbFBSUnp7e1tbm4ODQ29ublZV14sSJpW3w8rJ+/frnn38+ODg4JCTkm2++2bdvH74dmJCQ8ODBg0uXLlGp1I8//vjdd9/dv39/eXk5k8ncsmXLUrd6GYuPj4+Li+vu7q6pqZHJZGFhYQihuro6T09PmUympaUVGRlpaWnJZrOrqqpKSkrKy8uXusnL2OnTp+vr6zs7O7OzsysrK2NjY52dnfPz8xMSEvD9lEOHDr3yyisUCkUgEJSVlZ08efLxDkRLSkp6kg1fDoaHh7W1tV944QV3d3dLS8tnnnmGxWLp6uo6Ojq6ubkhhMzNzcPCwlpbW3V0dDIzMx86hQeYBYvF2r59e1tbG5VKzcjIcHR0xHFLS0sfHx8mk0kQxMjISH9/P5vNPnr0aFBQ0NI2eNkJDw9XKpUPHjx46623oqOjcdDIyMjV1RVfMvLw8HB1dW1paXFxcTl58iTMhDUfLi4uzz77bHNzM5fLPXXqFP7awWAwuFyup6cnhUKxtLQcHBwUiUSurq5nzpyBid7mQygU6unpvfTSSy4uLpaWls7OzkZGRnp6ek5OTviak7W19auvvtrS0mJkZPTll1+S4/IeFUzDBAAAYFVbmfcIAQAAgDmCRAgAAGBVg0QIAABgVYNECAAAYFWDRAgAAGBVg0QIAABgVYNECMCCEwqFPB5PJpMtdUNQQUHBzZs38XJXV1d2djZ+AM0iE4lEPB5vYmJi8Q8NwHSQCAFYcPHx8d9++y1+2ufSSklJOXv2LF7+448/du3axefz57/b6urqrKysuddnMpmJiYmnTp2a/6EBmD9IhAAsrPb29rNnzyYmJi51Q6Zydnb+9NNPTU1N57+rwsLC2NjYuden0WiHDh1KSUkRi8XzPzoA8wSJEICFdebMGVtb240bN06JKxQKPp+vcS5iqVTa398/06XU2UvJCuRDo8bHx/l8/vRJYtetW5eQkDB9VpbJyUly9pUplErlwMDAI012I5PJ+vv7pVLplPiOHTsmJycvXLgw910BsEAgEQLwaGJjY62trTs6OvDq5OTkpk2bNmzYoPFmm1KpzMnJee2119QnBRwfHz9w4ICZmZmFhYWBgYG9vT2eHQIhxOfzw8PDjYyMCIIwMTGJjo5Wn1Csr68vNDQUP6CVxWLFxsaOj4/jor///pvFYuXk5EREROAKQqFQoVB8+OGHLBbLwsLCysrq+++/V29bcXExQRBdXV141dPTMy4u7vjx42ZmZubm5iwWKzMzk6xcVVUVEBCgo6PDZrP19fX9/Pzq6+tx0cGDB9PS0qRSKYvFwsfCcbFYHBUVZWJiQhAEk8kMDw9Xz6/6+vqbN2/Ozs5+nPcAgCdLBQB4FKOjo46Ojq6urlKpVKVSffTRR3Q6vaKiQmPl2tpahFBBQYF6MDg4mMFgJCcnV1dXV1VVZWRkXLlyRaVSTUxMuLm5GRsbZ2dn19fXp6ena2lphYWF4a2kUqmTkxOLxbp48WJ9fX1qaiqdTn/jjTdw6b179xBC5ubmb7/99o0bN37++WeJRHL48GEKhZKYmNjY2HjlyhUbGxtjY+OdO3fiTXg8HkKoo6MDr3K5XAsLi40bN167du23334LDQ2l0WgtLS24tKioKCEh4fbt262trcXFxe7u7hYWFqOjoyqVqrW1NTIyUltbu7S0tLS09ObNmyqVSi6X+/r6EgSRnZ3d3Nx8+fJlW1tbb29vhUJB9sMXX3xBo9FEItGTeWMAeFyQCAF4ZI2Njbq6unFxcb/88guVSk1LS5up5nfffYcQampqIiN40GZ6evr0yjgz5ebmkpGEhASEUH19vUqlysnJQQj9+OOPZOnBgwcpFEpbW5vqv0QYGBhIlkqlUkNDw4iICDKCf3fOkgjNzc3FYjFeFQqFWlpaqampGl8XPlxRURFeTUxM1NPTU6+Qm5uLELp79y4ZKSsrQwjdunWLjFy9ehUhVF5ervEQACya1TgfIQDz5OzsfOLEiejo6PPnzwcFBR04cGCmmvhiIIvFIiOlpaUIod27d0+vXFdXR6PRtm3bRkZ27tyZkpJy584dV1fXuro6bW3trVu3kqU7duxIS0u7c+eOg4MDjoSEhJCl9+/fHxkZCQ8PJyObNm2afkdQnZ+fH5PJxMsmJiYWFhY9PT1kaU9PD4/H6+7uxjf8qFQqeX14uuvXrxsaGo6Pj5NXfRUKBY1Ga2pqImdPxeN0BgYGZmkSAIsAEiEAjyMyMvKTTz4RCoVJSUnq9/+moNPpCCH10SVDQ0P6+voaE1J3d/dTTz2lpaVFRiwtLRFCePbz7u5uNptNpVI1lmLk/TmEUG9v75QIQoggiFlelHrCRghpa2uTQ3IuXLjwzjvvODo6ent7m5iYUCgUKpU6y5hPPp8vkUi2b9+uHmQymSKRiFyVy+UIIfXXC8CSgEQIwOOIiYmRy+V2dnYxMTF3796d6WzOZrMRQgKBYM2aNThibGw8NjY2OjqK53RVp6+vLxQKlUolme3wD0ojIyNcOjQ0pF5fvRRTT8nm5uYIIY2bPIYjR44EBwcXFBTgVbFYfOzYsVnqGxkZmZqazv4nRZzCcRcBsIRg1CgAjywvLy87O/v06dM8Hq+2tjY+Pn6mml5eXgihxsZGMuLv748QwvfnpvD29pbL5SUlJWSkqKgIIeTr64sQ8vHxkUgkt27d0lg6HZfL1dHRuXbtGhmprq6ekhfnSKVS/fPPPxs2bCAjxcXF6hUMDAwmJibU/6Hh7+8/MDBw+/btWXbb0NCgo6Pj4uLyGE0C4Ela6puUACwz7e3thoaGUVFRePX48eMUCuXq1asz1edwOO+99x65qlAofH198dDQoaGhwcHBwsLCX3/9VaVSjY6O2tnZ2dnZlZWVicXi/Px8JpPp7++PNxSLxdbW1vb29hUVFcPDw3l5eQYGBi+//DIuxaNXeDye+qH37t3LYDDOnTs3MjJSW1vr4uKip6c3y2CZ3bt3q2/+9NNP79mzBy97eHhwudympiaJRFJYWGhra0uj0ZKSknDpTz/9hBBKTU2trKysqalRqVRjY2OOjo5WVlY8Hk8gEAiFwt9///2DDz5ob28n9x8YGLh58+ZH6HoAFgYkQgAegVQqdXd3d3JyGhsbwxGlUhkaGspisbq6ujRu8tlnn5mZmclkMjIiEAjCwsLI65+GhoaXL1/GRX/99ZeHhwf5PTUoKGhgYIDcsKmpyc3NDRdRKJSQkBCBQICLNCZCiUTy+uuv4/oMBuPYsWPr169/vERYXV1tZ2eHd2VqalpYWKirq0smQqVSuW/fPoIgKBQKg8HAQfyvR/JlUqnU5557rqenB5fy+Xw6nf7DDz/Mue8BWCgU1X+PnwAALITBwUF7e/tz586ROQkbGhrq7Ow0MDDgcDg6OjrqRffv3xcIBNbW1ng4zBQdHR1CoXDNmjWzj3whdXd39/f3Ozg4TB+hg0dyzvGFyOXytrY2pVK5bt26uT83VSQS3bt3T09Pz9raWr0BqampX3/9dVtb2//DI1jBKgeJEIAF9/nnn1+8eLGhoUF9zOdqNjY2xuFwMjMzpwwrBWBJQCIEYMHJ5fKenh4bGxv8bwoglUr//ffftWvXzvLPEwAWDSRCAAAAqxpcqAEAALCqQSIEAACwqkEiBAAAsKpBIgQAALCq/Q+q9CxlurG6fAAAAABJRU5ErkJggg==",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 12
}
],
"cell_type": "code",
"source": [
"Plots.plot(getindex.(points, 1), u_points, xlabel = \"x (coordinate)\", ylabel = \"u (temperature)\", label = nothing)"
],
"metadata": {},
"execution_count": 12
},
{
"cell_type": "markdown",
"source": [
"*Figure 3*: Temperature along the cut line from *Figure 2*."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"Secondly, the horizontal heat flux (i.e. the first component of the heat flux vector) is plotted."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=1}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdeXxMV/8H8DP3zmTf98WaIEEoUhoiCGLfIqhdG081bZ8+WqVapeqptbRFq4q2llDUvstCRO0iREJIJEEE2SP7ZOae8/tjnt/QZEJkmcnMfN5/eM3cOZn7zWTMZ865554rYowRAAAAfcVpugAAAABNQhACAIBeQxACAIBeQxACAIBeQxACAIBeQxACAIBeQxACAIBeQxACAIBeQxACAIBeQxACAIBe0+4gjIuL27p1aw0bM8awnpyaCYKg6RL0DqVU0yXoF7zg6lfvHyzaHYTx8fGRkZE1bFxRUSGXyxu0HqiktLRU0yXoF0ppeXm5pqvQL3iTq1+9v+baHYQAAAB1hCAEAAC9hiAEAAC9hiAEAAC9hiAEAAC9hiAEAAC9hiAEAAC9hiAEAAC9hiAEAACtsTOFTjxnUL/PKa7fpwMAAGggvybSb6/TvX4yQgzr8WkRhAAAoAVWxNENd2j0MN5RVM+rRiMIAQCgUWOEzL0iHHvI/h7Gu5qKiorq+fkRhAAA0HgJjMz4W0gsYOeGi63rc0D0OQQhAAA0UlKBTDojFFWwiCFi0wbLK8waBQCAxqhETkZGyOWUHBrQgClIEIQAANAI5UvJgBNyR2PR3n68Ed+w+0IQAgBA4/K0jPQ5Jn/TTrSlNy9u+JhCEAIAQCOSUsh8D8snunNruvMitewRk2UAAKCxiM9jg8OErztzMzzV109DEAIAQKNwNZuNjJCv9uHHual1tBJBCAAAmhf1hE04Ld/SWzyoiXoGRJ9DEAIAgIYdekBDzgn7+ot9HdWdggRBCAAAmvXbXfpNLD05WPyGjQZSkCAIAQBAg9beoqsTaNQQvrWlZlKQIAgBAEBTVsTRrcn07DC+ianGUpAgCAEAQP0oIx+eF67lsLPDxHZGGi4GQQgAAGpVQcnUM0JmGTs1VGwh0XQ1CEIAAFCnUjkZc0ouFolODBI39CKiNYQl1gAAQE0KKsjAk3JbQ9H+/g2+lHbNIQgBAEAdMsuI/zF5JxvRVrUspV1zjakWAADQUQ+Kmd9R+cAmop968Jwmp4iqoLFjhCUlJYQQU1PT6hoUFxfzPG9sbKzGogAAoP7lSUmfY8KcDtyH7Rpj70sDNcnl8mnTprm6ujZp0mTq1KkymaxSg0uXLrm4uDRr1szR0bFLly5xcXHqLxIAAOrLh+eF0S1EjTMFiUaCcPPmzXFxcRkZGRkZGbdu3frjjz8qNXB3d7969WpeXl5+fr6/v/+MGTPUXyQAANSL9Yk0uZAtfbPRzI2pQgNBGBoa+v7775uampqYmISEhISGhlZqYG9v7+rqSgjhed7f3z87O1v9RQIAQN3dLmALrwk7+vCGjTcHNXGMMCUlxdPTU3Hbw8MjNTW1apuKiootW7bk5+fv2rVr2bJl1T0VY6y4uFj5DBKJxNXVleMaae8bAECvlAtk4mlhRTfe06qRTY/5Jw0EYWFhoXKOjJmZWUFBQdU2jLHU1NRHjx6VlZWJRNW+gvfu3YuMjOzXr59yy6ZNm3x8fFQ2lkqlHMdJJI1gGQO9UVJS8pI/H9Q7Sml5eTmlVNOF6BHFvD9QadY1sZupaKxLeXFxfT7ta32wmJiYvLJ3pIEgtLe3V4ZfQUGBg4ND1TaGhobLly8nhJw7d27QoEHDhw9XOX20devWo0aNqjq4qpJEIkEQqhljzMzMTNNV6BFKqVgsNjEx0XQh+gVvcpWOp7PwJ8L1QLGZoWH9PnO9f7BoYBTRy8srJiZGcTsmJsbLy+sljd3d3UtKSkpLS9VSGgAA1INHJSz4rHxXX966nkOwQWigR/jBBx8EBwf37dtXJBKtXr1606ZNiu0BAQFLly7t2rXrgQMHxGJxmzZt8vLylixZ4u/vb2trq/46AQCgFigj06KF/7TnfRy048iIBoJw8ODBixYt+vjjjxlj33zzzdChQxXbra2tFeOWxsbGP/30U1pamrW1tZ+f3+eff67+IgEAoHaW3KACI3Pf0Jp5iyLGmKZrqL3t27eHhYXV8BghJsuoX1FRkbm5uaar0COKyTI4RqhOxcXFOEb4oivZbHi4/OpIcTOzhuoO1vsHi9YkNgAANHIFFWT8aWFjT77hUrAhIAgBAKB+fHheGN5MNLK5liULLswLAAD14Le79HY+uzRS+2JF+yoGAIDGJrGAzbsqRA9rLBedfy1a1oEFAIDGRiqQiVHC0q5828a9lFp1EIQAAFAnc64I7haif3loa6BgaBQAAGrvRDo79IBdD9TiNNHi0gEAQLMyy8h754TdfXkbbVhKrTra2pMFAADNooxMPiN/35PzddTKQ4NKCEIAAKiNFTepjJJ5nbQ+RzA0CgAAr+1SFluTIMSMEvPa3RskBD1CAAB4XcUyMjVa+LUn38RU+2MQQQgAAK8r5LwQ4CoapW1LqVUHQ6MAAPAaNifRG7nsqhYupVYd3flNAACgod0rZHOvCJFDxMY6lB460rEFAICGJqNkUpTwX2++o40uHBpUQhACAECNfHFVcDUVhbTVteDQoc4tAAA0mBPpbG+adi+lVh0d/JUAAKB+ZZWRf/0thPbR7qXUqqNrPVwAAKhfjJDgs/L3PEV9XXTq0KASghAAAF5m5U36TEbmd9LCS+7WDIZGAQCgWtdy2PfxwuURYrHu9pt09zcDAIC6KZaRSVHCGh++hbluDooqIAgBAECFEjl5+7S8t7NovLuOJ4WO/3oAAFALT0pJ76NyJ2PRzz109tCgEoIQAAD+IT6PdT8sH9ZM9HsvXqIHKYHJMgAA8Fx4BptyRr7ah5+g6yOiSghCAAD4n/WJdPF1eihA7OOgy7NjKkEQAgAAoYzMuSKcSGfnhvMtdXqOaFUIQgAAfVcukHeihexydmGE2MpA09Wonb4MAQMAgEpPy0ivo3JjMTkxSB9TkCAIAQD02a181v2wfEhT0eZevIG+BgKGRgEA9FRkBpt0Rv7DW/ykVvqagYQQBCEAgH76I4l+dVXY00/cy0m/psZUpbEgTEhIyMnJ8fb2Njc3r/qoVCqNj48vKyvz8vKytrZWf3kAALqKEbIoVthxj50ZJvaw1PcUJC8JQrlcnpaWlp2dLRaLnZ2dmzZtWo97nTZt2tmzZ9u0aRMXF3fixInOnTu/+OilS5eGDBnSsmVLMzOzuLi49evXT5gwoR73DgCgt8rkZGq0kF3OLo8U6+RVdmuhchDKZLKDBw9u27YtOjq6qKhIud3JySkgIGD69Om9e/eu4y7//vvviIiI27dvW1lZLV269Msvvzx58uSLDVxcXGJjY1u0aEEI2b179/vvv//2229znF4PYQMA1F1WGRkRIW9jIQofLNbbqTFVPQ9CSumOHTvmzZv35MkTHx+f999/39PT09raWi6X5+TkxMfHnzt3rk+fPt7e3t9//31d4nDfvn3Dhw+3srIihEydOnXBggWFhYUWFhbKBs2aNVPe7tixY3FxsVQqNTY2rvUeAQAg+RkbFi5McBct7MJjPPRFz4Nw9+7dc+bMmTVr1uTJk11cXFS2vnXr1m+//TZ48OD4+Hh3d/fa7TI9Pd3b21tx29XVleO4jIyMF4PwRWvWrBkxYkR1KVhWVvbw4cO//vrrf7+MWNy3b9/qnopSqvwX1INSihdcnej/03QhekRbXvDTT8jkM8J33bjJ7iJGKdN0PXXxWq95TUYTnwehn59fSkqKqanpS1q3b9/+xx9//Pzzz1/e7OXKy8slEonitkgkkkgkZWVlKltu2LAhLCzswoUL1T1VTk5OWlrarl27lFtcXFy8vLxUNpZKpRzHKXcNalBWVsbzun8Nl8aDUlpeXq7pKvRLWVlZ4z9wsyON/zqO39pD3tOBlpZqupo6e60PFhMTk1f+gZ4HYZMmTWr4vM7OzjVsqZKTk1Nubq7idklJSVlZmcon3Lp16+LFi6Oiol6yu6ZNm/bu3Ts0NLQm+5VIJAhCNWOMmZmZaboKPUIpFYvFJiYmmi5EvzTmNzkjZH6MsCeN/T2cb2OpI8vG1PsHiwa+yPj4+Jw9e1ZxOzo6umXLlo6OjpXa7NmzZ968eWFhYa1atVJ7gQAAukAqkClnhKjH7PxwcRucJlE91UFYWlq6fPlyHx8fFxcXm3+q+y4nTJjw8OHD2bNn//XXXzNnzpw1a5ai3xoYGLh48WJCyIULFyZMmNC/f/8jR46sWLFixYoV+fn5dd8vAID+yJWSASfkFQI5NURsb6Tpaho31ecRTpky5cCBA/369Rs3bpyRUT2/hGZmZufPn1+7du2RI0cWLVo0ceJExfagoCDFJB1TU9PZs2cTQpT5pxXHogEAGol7hWxomDCoiehHH55DV/BVRIxVnj1UUVFhZma2dOlSRRo1Ztu3bw8LC6vhMUJMllG/oqIilSsHQQNRTJbBMUJ1Ki4ubmzHCE+ks+Cz8uXd+GmtG/ssntqp9w8WFT3CgoICmUzWr1+/etwNAAA0tApKvrwq7E1jf/UT++n9CqI1p+L7gr29fZs2beLj49VfDQAA1E5aEet9VJ5aSGIDkYKvR0UQikSiP/74Y+nSpREREXK5XP01AQDAa9mWTLsdko934/YH8LZYQfQ1qZ4sExwc/ODBgwEDBkgkkkrD33l5eWopDAAAXq1QRkLOCfF5LGqo2MsaHcHaqHbWaHFxsZpLAQCA13I1m02MEro7iC6PFJvg8rK1pfqVmz9/vprrAACAmmOErE2gS24IP/Xg33bTzdmhavPqrxCNcHIwAIA+e1JKppyRyyi5Fihuaorh0Lqq9nvEqVOn/P39ra2tzc3NLSwsevTosX//fnVWBgAAVR1LZ94HZb2cudNDkYL1Q3WP8ODBg0FBQQ4ODuPHj3dycsrJyTl27FhQUND69etDQkLUXCIAABBCpAJZGCvsSmG7++IEifqkOghnz57dv3//Q4cOKddX+/HHH4ODg7/88svg4GADAx1ZwhwAQFvcfcYmnBZamItiA8U2OEGiXqkYGs3Ozk5JSVm4cOGLq4yKxeJvv/22oKAgMTFRjeUBAAD5/S71OyIPacvt788jBetdtZNlRKLK/e6qWwAAoEEVysgH54QbuezUEHEHG3wINwjVS6y5ubktXry4oqJCuZFS+t///tfS0tLT01ON5QEA6K8r2azLATkvIldHIQUbkOoe4XfffTd27Fh3d/dRo0a5uLhkZWUdP348KSnpp59+MjREtxwAoGFRRpbF0Z9vCRt68iOa4zTBhqU6CIOCgk6cOPHtt9/++uuvcrmc47jOnTvv3r173Lhxaq4PAEDfZJSwKWcEQkjMKLErTpBoeNUeIxw4cODAgQOlUmlpaamRkZGxsbE6ywIA0E+HH9AZ54SQttyCzjyPEFSLV6wsY2hoiLFQAAA1KBfI3CvCoQdsbz9xT5wmqEbPg/DmzZtHjx719/fv3r37mjVrSkpKVP7AvHnz1FUbAIC+OJ/JPjgvtLMSxY0WW+JUbfV6HoTXr19fsGDB0qVLu3fvvmLFiszMTJU/gCAEAKhHGSVs7lUa/YSt7MaNd8e8GA14HoTTpk2bNm2a4vbjx481VA8AgL6QUfLLbfrtdWFSK+72GLG5RNMF6SvVxwifPXtmamoqFv/jUblcXlRUZG1trZbCAAB02aEH9LPLtIO16MoosZs5jghqkupuuIeHx5UrVyptvHr1qo2NTcOXBACgyxIL2MAT8nlX6Xpf/kAAjxTUuNe4pLFcLpdI0HUHAKilYhlZFS+sT6SzvPhPO3AGOCDYOPwjCMvLy8vKygghlNKioqL8/HzlQ6Wlpfv27XN1dVV3gQAA2o8ysjmJzo8RhjXj4kdLHHBidmPyjyBcv379rFmzFLcHDRpUtfXy5cvVURQAgA65ms1mXhTkjBwIEPs4YCC00flHEPbv33/Dhg2EkDlz5oSEhLi7uysfMjMz69ChQ4cOHdRdIACA1npSSr6JFY6nsyVvclNac8jAxukfQaiMuoqKisDAQAyEAgDUjlQgqxPoqnjhPQ8ucYzYDPMrGjHVk2WCg4MFQai0sbi4WCQSmZqaNnxVAABaLDKD/eei4G5BLo0Qu1ugH9jYqZ60NHDgwK+//rrSxl9//fWNN95o+JIAALTVnQI26KT8PxeFNd35IwOQgtpBRRBWVFRcuHBh1KhRlbYHBgampKQ8fPhQLYUBAGiTZxXks8tCr6PyQU24uNHiAFdEoNZQEYQ5OTmUUkdHx0rbHRwcCCHVrUEKAKCfGCHbkmnbvbLsMhIfJPnEi5PgBEGtouIYoZWVlVgsTkhI8PT0fHF7fHw8IcTW1lZNpQEANHqXc7gvIuUGPDkyQOxth16gVlLxvcXExKRv375z5869f/++cmNWVtbMmTPbtWvn5uamvuoAABqrK9lsSJg8+KJkphd3bjhSUIupnjX6448/9uzZ09PTs2/fvq6urk+fPj1z5gylNDw8XM31AQA0NnF5bMl1ejmbzfLitvmU2lmaaboiqBPVQdiuXbvY2NilS5dGREScP3/e0tIyMDDwyy+/bNu2bb3stbCwcM+ePUVFRUOHDm3dunXVBlKpNC4u7sGDB/7+/nZ2dvWyUwCAOrqZxxYrI7APb8ST4mJN1wR1Vu2i2y1atNi4cWND7LK4uLhbt27t2rVr2bJl165dT5w40b179xcbMMYsLCyaN29+//79qKgoBCEAaFzVCASd8bKrT1RUVKSmpmZnZ/v5+dXjLkNDQ+3s7Pbt2ycSiezt7ZcsWXL06NEXG4hEoqysLEtLS0zMAQCNu5bDvokVbuaRL97gtvvzuGSE7lH9JxUEYe7cuVZWVm3btp0wYYJi4/Tp0ydPnlz3XYaHhw8bNkwkEhFChg8fHhERQSmt1MbS0rLuOwIAqIvYHDYiXAiMEAY34ZLGij9oiwsn6SbVPcKvvvpq7dq1c+fOtbCw+OGHHxQbR48e/fbbb0ulUkNDw7rsMiMjw9nZWXHbxcWloqIiJydHcZLi68rNzb1+/fqcOXMUd3menz59erNmzVQ2lkqlHMdVDV1oOFKp1MDAQNNV6BFKqVQq5XkM29VVQj5ZFs+dy2QftyWhPZkRLxA5kcpVtJRKpbhQq5q91geLgYGBot/1EiqCUCaTrVu37rvvvvv444+jo6OV2zt16lRSUpKent6qVauaV1yVSCRijCluK268ssqXkEgkVlZWLz55XWoDAD13I48sucnF5pLZXuwPX2aILxV6QEUQZmdnFxcX9+vXr9J2CwsLQsiLV+utHRcXF+XyNE+fPjUwMKj1sUBbW1svL6+vvvqqhu05jsN3N3WqqKio4/gBvBZKKWMMr3nt3Mhli2JpTA77vCP3V3+uhtNhZDIZXnA1q/cPFhUD3hYWFhzHZWRkVNp+8+ZNQoiLi0sddzlgwICjR48q+oJHjhwJCAjgOI4Qcv/+/YKCgjo+OQDA64rLY6MjhWHhQl8XUfI48cfta5qCoBtUBKGZmVnv3r2/+eabZ8+eKUca8/Ly5syZ4+3tXfeLFE6ePDk7O3vMmDGzZ89evny5sj8XGBj4559/Km7PmjVr3LhxxcXFCxYsGDdu3NOnT+u4UwCAqm7msaBIYWiY0MdZdA8RqK9UT5ZZu3atn5+fh4dHu3btCgsLJ06cGBERUVpaeurUqbrv0tzc/MqVK3v37i0qKrp69aryiOMPP/zQokULxe2BAwcWFhaOHTtWcdfMDAs3AEB9upXPVsTRiAz6iRe/vQ9v/LJTyUDHqf7je3l5xcbGLl68WJl//v7+X3/9tZeXV73s1cLCIjg4uNJGf39/5e2BAwfWy44AACqJesJWJ9CYbDanI7ehpwQRCCreAmVlZWvXrh06dOjvv/+u/oIAABqCVCC7UunqBFohkJle3C5/9ALhf1S8EQoKCr744os+ffqovRgAgPqXXU5+TaTrE4WONqLlXfkBTXCWFfyDiskyDg4OdnZ2uBI9AGi7W/nsvb8Fjz2y9BIWOUR8cpB4IFIQqlDRI+R5fvHixV9//bW3tzeuPggAWocREvaI/RgvxOezD9vySWMldkaargkaMdVj5KdOncrKyvL09Gzfvn2liz9ERESopTAAgNdWJieh9+iaBGrAk0+9uPHuWB0UXq3ag8VdunRRZx0AAHXxpJT8kihsvEN9HLiffXl/Z4yAQk2pDsK//vpLzXUAANTO9Vy2OoEefUgnuHPnholbWyIC4fVg+jAAaCXKyJGHdHUCTSkk/27PrfaRWGPJT6iV50GYlZWVnJzcokULV1fXq1evVlRUqPwBX19fddUGAKBCsYxsTqJrb1FbI/KpFxfUghPjQCDUwfMgPHbsWHBw8NKlS7/88sthw4ZlZWWp/AHlFZQAANTsSSnZcEdYd5v6OIjW+/L9XTEKCvXgeRAOHTr077//Vqz2eeTIkep6hAAAakYZCXvENt6hfz+l77ThYkaJm5shAqHePA9CBwcH5WXiu3XrpqF6AACee1pGtibRDXeojSGZ4cmF9pGY4YqiUN8wWQYAGh3KSEQG23iHRj2h41pye/vxXezQBYSG8jwIw8LCVq1a9cofwAn1ANBwMsvIliS68Q414snU1txvfpgLCg2u2h5hTEzMs2fPPD09XVxcsrKyEhMTJRIJpowCQENghJzKYBvu0FOPaVALbnc//k10AUFdngfhwIEDlVcBXLNmzaNHj86dO9e+fXvFlvv370+cOLFz584aqBEAdFe+lOxJo2sSKCNkWmtuQ0+JDbqAoF4qeoQymeybb77Zu3evMgUJIS1atPj99987duz4+eefV1p9FADgdTFCoh6zDXdoRAYNbM790Yt/ywFdQNAMFUGYk5NTUFDg5ORUabuTk5NcLr9//z6CEABqraCC/JVKf7pF5Yy805pb74suIGiYivUYbGxszMzMNm7cWGn7pk2bOI5r1qyZWgoDAJ3CCDnzhE2MElrukl3IZBv9+MQx4rlvcEhB0DgVPUJDQ8Mvvvhi/vz5CQkJgYGBTk5OOTk5J06cOHr06AcffKA81xAAoCYUXcCfb9MKgbzbhvuph8QW4QeNiepZo1999ZWtre2yZcs+/vhjxRZHR8fFixd//vnnaqwNALSYwEjYI7YliUZk0BHNufW+vK8jjgJCY1Tt6RMhISEhISF5eXmPHj1ycnJCRxAAaiixgG1NpqHJrJkZeacNt9FPYmWg6ZoAqveKlWWsrKxkMpmNjY16qgEA7fWsguxKpVuS6MNiMqW16NQQ3tMKXUDQAq8IwkePHjVv3jw2NhZnEAKASpSRU4/ZliR6PJ0GuHILOvMDm4h4JCBoD6w1CgC1lPyMbU2m25KZgzF5pw33Uw+cCAFaCUEIAK+nSEb2pNHNSTT5GZvUijs2kO9ggw4gaDEEIQDUCGUk+inbkkQPP6B9nLk5HbjBTTkJLg0P2u8VQejs7BwTE+Pp6ameagCgEUorYluT6dZkZikh77bhVr0lsTfSdE0A9Ud1EBYUFFhZWRFCJBKJt7e3cvvdu3c9PDzUVBoAaFSJnOxLo5uT6K18NsGd29+f72yLIVDQQarHNfz8/E6ePFlp4969e996662GLwkANIkycuYJCz4rNN0p25vGPm7PPZooWdMdKQg6S3UQtmrVaujQoQsWLBAEgRBSWloaHBw8duzYwMBA9ZYHAOqTWMC+ihFa7pZ/eknoYCO6PUZyeAA/ugVngAOBoNNUD43u27fvu+++W7BgwdmzZ+fPnz9z5sz09PRt27ZNmTJFzfUBQEPLLie7UmjoPfq4lEx0Fx0byHtZo/MHekR1EHIc98UXX/j5+Q0ePHjAgAEeHh4xMTE4OgigS6QCOfKQht5jZ5/Q4c24JW/y/VxEHBIQ9E+1s0aLiorWr19fVFTUtGnTBw8enD59GkEIoAMYIeefstB7dF8a7WwnmtKK29FHYibRdFkAmqN67P/27dvdu3c/cuTIzp07U1JSZs6c+dFHHwUGBubn59fLXg8ePDht2rSPP/74zp07KhskJiZ+/PHH06ZNO3z4cL3sEQDuFbJvYoVWu+Uh5wV3C9GN0eKIweKprTmkIOg51UE4YMAAIyOj2NjY8ePHSySS5cuX79+/Pzo6+sVTKWrtr7/+CgkJGTBggJ2dXc+ePTMzMys1ePr0ac+ePe3t7QcMGDBjxoy9e/fWfacAeitPSn5NpL5H5H5H5M8qyN7+fEKQ+POOXBNTDIMCEFLd0GhISMicOXMMDZ+vGzhq1KjOnTtPmDCh7rtcuXLlsmXLJk2aRAiJjY39/fff582b92KD3377zc/P7+uvvyaESKXSlStXjhkzpu77BdArFZScSKfbktnpx3RQU+6rTvwAV5EY8z8BqlAdhPPnz6+6sXnz5tHR0XXcn0wmi42N7d27t+Junz59oqKiKrW5dOlS//79Fbd79+793nvvyeVysRirwQHUyOUsFnqP/pVK21mLprbmNveWWGDwE6B6r5cuEkld/z9lZWVRSm1tbRV37e3tnzx5UqnN06dPlQ3s7OwopZmZma6urlWfLT09PTo6evTo0cotc+bM6dChg8pdS6VSjuPq/itAzZWUlIhEGH9Tk4clZGcav+s+LxLJJrSkZwLkzUwJIYRISbFUw7XpsJKSEk2XoHde64PFxMSE414xEqLubpaxsTEhpKKiQnG3vLzcxMSkahtlA6lUSgip2kbBzs6uZcuW48ePV9zleb5t27bVNeZ5HkGoZoIgVPfngPryrILsu89C77HEAja2BdnUXdarqWIlUFwVXh0opXiTq9lrfbC8MgWJ+oPQ2traxMTk4cOH9vb2hJCHDx9W7eq5uro+fPhQcfvhw4empqaKhU+rMjY2btas2bhx42qya+7/1aF8eD14wRuOnJKwDBaaTMMe0X6u3GcduEFNOTGh5eUyvObqhDe5+tX7a67uv59IJBo9enRoaCghpKysbM+ePUFBQYSQ0tLSnTt3lpWVEUKCgoL27NmjuB0aGhoUFDePT40AACAASURBVIThNQClaznsk0tCk52yZTcEfxdR6tuSvf34Ec2xEBpALWlgBsrChQv79u178+bNx48ft2nTZuTIkYSQ3NzciRMnPnr0yNXVNTAwcOvWrV26dHF2dk5JSTl9+rT6iwRobB6VsO33WGgylVIypRV3frjY3QJfEAHqgQaCsFWrVnfv3o2JiTE3N+/UqZNio4uLy/37952cnAghYrH4yJEjcXFxRUVFXbt2NTLCpc9AfxXJyP77dFsyjctlY924jX58D0eMkADUJ9VBSCndtWvXoUOH0tPTlfNWFGJiYuq+V2NjYz8/vxe38DzfvHlz5V2RSKTMSAA9JDASmcFC79FjD2lvZ+6jdtzQppwhr+myAHSR6iD84IMPNm7c2KZNm7Zt2xoYYO4ZgPrczGOhyfTPFNbUjExuxa3pLrE1fPVPAUCtqQhCuVy+bdu2OXPmrFixArNUANTjSSnZmUK3JdOCCjK5lej0UN7DEv/7ANRBRRDm5uaWl5ePHz8eKQjQ0MoFcuA+3ZZMr2SzUc25Nd35Xs74jwegViqC0N7e3tXVNS0trUuXLuovCEBPXM5iW5LpX6m0q71oWmtuf3/OGMsIAmiCiv95HMetW7du3rx5bdu2bdeunfprAtBhT0pJ6D26JYkKjExrzcWNFuMqEACapfor6Lp1654+fdqhQ4dmzZopl/1UqJdZowD6RiqQww/p1iR6IYsFteA2+fG+jsg/gEZBdRA2b968ulXNAOC1XMthW5LorlT6ho3onTbcX/04EwyBAjQmqv9Hbtq0Sc11AOiYrDKyI4VuTqIlMjKtDRczStzcDF1AgMYIX00B6pOMkuPpdHMSO/uUjmjG/YRZoACN3vMgfPDgwdWrVzt06ODh4XH06NHy8nKVP4CLxQOodDOPbU6if6ZQT0vRO2247X0kZrjkF4A2eB6Ep0+fDg4OXrp06Zdffjl9+vSsrCyVP8AYU1dtAFogV0r+vEc3J9FcKZnWWnQBa2EDaJvnQThmzJhevXop5oheuXJFLpdrriqAxk5g5OQjtjmJnsqgw5pxK9/i/Z1FHBIQQAs9D0Jzc3Nzc3PF7RfXvwaAF2WUsN+T2G93qIspmd6G+91PYonleAG0GSbLANQIZeTkI7bxDv37KR3vzh0ZyL9hgw4ggC5AEAK8wuNS9sdd9ttd6mhMZnhyO/wlpvh/A6BD8B8aQDXKSEQG23CHRj+h49y4AwF8Z1t0AQF0EIIQoLKnZeSPu/S3u9TGkMzw5Lb1xokQALoMQQjwP5SRU4/Zhjv09GM6piW3px/vbYcuIIDuUx2Eqampbm5uVbdHRUX5+/s3cEkA6pZZRrYk0U13qYWEzPDkNveSmKMLCKA3OJVb+/btu2bNmhe3UEpXrFgREBCglqoA1IERcuoxe/u00HavLLmQ/enPxwaKQ9pySEEAvaK6Rzhq1KhPPvnk0qVLGzZssLCwePz48aRJk86dO7do0SI11wfQELLL/9cFNObJ+57cxp44FxBAf6kOwtWrV/fs2fNf//qXt7f3rFmzFi5caGhoePr0aT8/PzXXB1C/zj1lvyTSE+l0VAtuW2/exwFHAQH0XbWTZcaMGdOxY8c333zzww8/7Nix46lTp+zs7NRZGUA9klGyJ43+GE8LZeTf7bhffCVW6AICACGkumOEhJD09PTg4GCpVOrr6xsfH79s2bKKigp1VgZQL3KlZFkcbblb/sddurALnzhG/HF7DikIAEqqg/DIkSOdO3dOT0+Pioo6d+7cli1b1q9f36NHj5SUFDXXB1Brd5+xD84Lrf+SJT9jxwfykUPEw5phXWwAqEx1EL733nu9evWKi4vr0aMHIWTq1Knnzp179uyZt7e3essDqI1zT9m4U0Kfo3JLAxI/WvxHL74j1gUFgGqoPka4atWqSZMmiV64sHaXLl2uXbv2/vvvq6swgNdWLpA/U+jqBEoImdme29ZHYsRruiYAaPRUB+HkyZOrbrSwsNi5c2cD1wNQG1llZHMS/fk29bQiS97khjXDCCgA1BSWWAPtdjOPrU6gBx/Qt924iMG8pxUSEABej+ogdHJyyszMVPkQY6wh6wGoEcrIiUfsx3jhzjPyUTsueZzE1lDTNQGAdlIdhPPnzy8pKVHeLS8vj46Ovnr16qeffqquwgBUK5GTbcl0TQI1k5BPvbhxbpyk2pOAAABeTXUQ/vvf/6668bPPPouPj2/gegCq9aiErbtNf7tLezlxm/x4PyeMggJAPXiNY4QfffSRu7t7enp606ZNG64ggKri89iyOBr2iE5tzV0eKXYzRwQCQL15jSCklBJC8vLyEISgNgn5bFEsPZ9JP+vA/9pTYoHrQgBAfavR0ZXy8vK4uLgPP/zQ2NjYw8Oj7ns9fPiwr69vhw4dFi1apMjXSjZt2hQSEhIQEHDr1q267w600e0CNvWMEHBc/qad6N44yWcdOKQgADSE15g1ampq+vPPPxsZGdVxl3fu3Jk8efKOHTtatWo1fvx4KyurmTNnvtiAMRYeHu7t7b1jx46CgoI67g60zq189u11Gv2EzunI/9pTYoJzfACgIdVo1qiRkVHz5s19fX3t7e3rvstNmzYFBQUNHz6cELJw4cJ58+ZVCkKRSLRnzx5CyMqVK+u+O9AiiQVs2Q0ankFD2nIb/TAQCgDq8BqzRutLQkLCiBEjFLfffPPNpKQkqVRqaIizwPTaixGY5IsIBAD1aZBRp5ycnGvXrlXd3qtXL2Nj46ysLCsrK8UWa2trxlh2dnaTJk1qsaOkpKSDBw+2bNlSuWXjxo0+Pj4qG0ulUo7jJBJ8xKpPcXHxK9vcLRR9f1t8OpOb3kqIHSKYixkpJ0XlaqhOB1FKpVKpIAiaLkSPlJSUYJkRNavJB4uSiYkJz79i0eHnQXj27NkNGza88kl37Njxyjb3799fu3Zt1e2dOnUyNja2srJS/hpFRUWEEGUuvq7WrVv379//+++/V9zleb5p06Ycp3oGkIGBAYJQ/czNzat7KLWIrYijhx7QkLbchj48eoF1RymVSCQmJiaaLkSPiEQiMzMzTVehd17ywVILz4MwPz8/ISGhXp70zTffPHbsWHWPurm5JSUlKW4nJSU5ODjU+m2keAu6ubnV7sdBU+4+Y99ep5EZdFYH/t44iRkiEAA053kQjhw5cuTIkWrY5eTJkydNmjR79mxHR8c1a9ZMmTJFsX3NmjXdunXr3r07IaSwsFAQBMZYUVFRfn6+paVldf080C5pRWz5//cC1/WQWOJK8QCgac/T5c8//1SEECFk8eLFDx48aKBd+vv7BwcHt23b1sHBobS0dP78+YrtBw4cSExMVNweOHCgu7s7IWTixInu7u4NVwyoTdIzNvWM4HNY3sJcdG+c5JsuPFIQABqD5z1CxlheXp7i9s8//9y3b9/mzZs30F4XL168YMECmUz24qDomTNnlLcvXrzYQLsG9Ut+xr69TsMy6Mz2/DpfiTkGQgGgMXkehG5ubmlpaXv37u3YsaMgCI8fP05NTa36A/V1QM7Q0BCnTOi8lCLRj7HCyUf0P4hAAGisngehj4/P0KFDx44dq7irvFEJJgpDTRRUkG9ihR3JBp90EP3UA+cFAkDj9TwIRSLRgQMH4uPjU1NTp02bNnv27HpZVhT0DWVk+z0694oQ4MpdHlThZo8jgQDQqFU+ob5Dhw4dOnQYMGDA6NGj27Vrp5GaQHudz2T/uSiYiMnxQeLOtqKiIowfAEBjp3plmb/++kvNdYC2yyhhc6/S6CdsRTdugjuHCwYCgLZ4fvrEjRs3anhCfVhYWFZWVoOVBFpGRsmaBPrGfrmtIbk9RjwRKQgAWuV5EGZmZnp7ewcFBR09elQul1dtWlRUtGXLlu7du48bN04mk6mxSGi8Dj2gbffKzzxhV0aJ13TnMS8UALTO86HRgQMHxsbGfvnllyNHjjQxMfH29m7Tpo2NjY1cLs/Ly7t58+bNmzd5nn/33XcPHTrk4OCgwaKhMUgsYJ9cFB6VkPW+fIArOoEAoK3+cYywffv2hw8fTk1NDQ0NjYqK2rdvX15eHs/z9vb2Xbt2nTp16uTJk21sbDRVKzQSBRVkeZywJYl+6sV/2oEzwOJ3AKDNVEyWcXNzW7hw4cKFCwkhMplMLBaLRPi+D4QQQhn5I4kuiBFGNOfigyT2RpouCACgzl5xPUJctAiUzmeymRcFI/5/p0ZouhwAgPrRIBfmBR3zuJR9cYWefsKWvslNaY1JoQCgU3B4B15GKpBlcbTjPnlTM3J3rHgqUhAAdA56hFCtww/orMvUy1p0eaTY3QIJCAC6CUEIKiQ9Y59eEu4VkrXd+SFNEYEAoMswNAr/UEHJl1eFnkfk/V25hCAxUhAAdJ7qIDxy5EjVjZTSFStWNHA9oEn3CpnvYfntApIQJPnUi5PgaxIA6AHVH3Xvvffeu+++W1JSotySlZU1ZMiQefPmqaswULc9abTnEfnkVtzBAN7BWNPVAACoi+og/P777/fu3du1a1fFMtxRUVGdOnWKi4sLCwtTb3mgDkUyMvWM8PU1enKQeKYX5oUCgH5RHYSTJk2KiYkxMDDo1q3blClTAgIC2rdvf+PGjf79+6u5PmhoV7NZlwNyYzG5NkrcCafJA4D+qfYokIeHx4EDB8Ri8fbt2728vI4dO+bo6KjOyqChMULWJNChYfJvvbkNPXkTzCAGAL1UbRCGh4f7+PhYWlrOnj379u3b/fr1e/TokTorgwaVVUaGhsl3ptJLI8Xj3TErBgD0l+pPwIULFw4ePLhr167Xr19fuXLl2bNn09PTO3XqdPToUTXXBw3hRDrrdEDW1U50bpjYzRzDoQCg11QH4R9//LFixYojR47Y2dkRQnx8fK5fv+7n5zdixAj1lgf1rIKSzy4LIeeFXX3Fi7x5MbqCAKD3VB8XOnHihJeX14tbrK2t9+/fv3btWrVUBQ0i6RmbECU0MxXFBoptDTVdDQBA46C6R1ApBRVEItHMmTMVt0+fPu3v79+AdUF925ZMfY/Ip7biDgTwSEEAAKVazhQsKyt78uRJ/ZYCDaRQRkLOCfF57MxQcXtrHBEEAPgHHCPScVezmfcBuVhELo9ECgIAqIBzx3QWI2RtAl1yQ/ipB/+2G77xAACohiDUTY9L2dQzgoySa4HipqboCAIAVAsdBR10LJ15H5D3cuZOD0UKAgC8AnqEOkUqkIWxwq4U9lc/sZ8TIhAA4NUQhLrjTgGbGCW0MBfFBoptcIIEAEDN1HJotEWLFpMnT67fUqAudqXQXkflH7Tl9vfnkYIAADWnukf47NkzSulLfszFxWXOnDl13Hd5ebmRkVF1jzLGZDKZgYFBHfeiD5bcoL/dpaeHir1wggQAwGtS3SP08PCweZXQ0NBa7/XXX3+1tbV1cnLq169fVlZWpUfv3r3br18/U1NTGxubt9566/r167Xekc6TUxJyTjj0gF4cgRQEAKgN1T3CVatWffbZZ66urqNHj3Z0dMzNzT1+/HhcXNzixYudnJwUbd58883a7TIpKenzzz+/cOFC27Ztp0+fPnfu3M2bN7/YoLy8/L333jt8+LCxsfH8+fPHjBmTkpJSu33ptnwpCYqUWxmKzgwV42qCAAC1I2KMVd06duxYc3PzP/7448WNCxYsiIyMvHjxYh13OX/+/Hv37u3atYsQcuvWra5du+bl5VU3Rnrv3r3WrVsXFxebmppWfXT79u1hYWE17JtKpVKO4yQSSV2KbzxSi9iwMGFIU9F33XiusXYFi4qKzM3NNV2FHqGUlpeXm5iYaLoQPVJcXGxmZqbpKvRLvX+wqBgaLS4u3r9//yeffFJp+8yZMy9dunTv3r067jIlJaVt27aK256enlKpNCMjo7rGhw4d8vb2VpmChBDGWEVFRf4LVOa67rmUxXoekf+nPbfqrcabggAAWkHFgFpxcTGl9NmzZ5W2K7YUFha+8kkzMzMrjXYqTJ061cXFpaCgQPkFiud5Y2Pj/Px8lc9z4cKFJUuWREREVLej5OTkgwcPhoeH/++XEYv//PPP7t27q2ysMz3CQ+ncp9fEv3QTBrlIi4s1Xc1LlZSUiEQIavVR9AhfPtMN6ldJSYmmS9A7r/XBYmJiwnGvOD9CRRA6Ojq2adPmk08+OXjwYNOmTRUbc3JyPvjgA1tb23bt2r1yx4IgFKv6hBYEgRBiZ2enTFOZTFZaWurg4FC18bVr1wIDA//8809vb+/qdtSmTZtx48bVcGhUIpHoQBCuSaCr4mnYYL6zrRbMp2WMYdRInSilYrEYQ6Nqhje5mtX7B4uKIBSJRL///vuQIUPc3d27du3q7OyclZUVExNDKd29e/dLTnhQcnFxWbx4cXWPtmvX7sKFC4rb169ft7KycnZ2rtQmLi5u2LBhmzZtGjRo0Ov8OrpMRslHF4RrOezySN7FBN0sAID6obrD2LNnz1u3bs2aNcvMzCwxMdHAwOD999+Pi4sbOXJk3Xc5bdq06OjoAwcOZGRkzJ8//91331X00r766qtt27YRQpKSkvr37x8YGGhiYhIZGRkZGVlaWlr3/Wq1ZxVkaJj8aSmJHipGCgIA1KNqJ903bdp0+fLlDbFLFxeXvXv3LliwICcnZ9CgQcq+I6VUMdUlIyOjU6dOycnJK1asUDy0detWfR7tyShhw8KFbvaidT14MZZJBwCoV6pPn9AW+nD6xJVsFhghzH2D+0977ctAnD6hZjh9Qv1w+oT61fsHC07DbtQO3Kch54Xf/PjhzbQvBQEAtAKCsPFaFU/XJtCTg8SdbXFQEACgoSAIGyOBkZkXhbNP2bnhfDMzpCAAQANCEDY6xTIy/rRcRsm54WILLTugCQCgfXDkqXHJKGG9jspdTEXHBiIFAQDUAUHYiFzLYT6HhamtuY09cZoEAICaYGi0sTj8gL53TtjYkx/ZHBkIAKA+CMJGYU0CXRlPjw0Uv2mHqTEAAGqFINQwgZFPLwlRj9n54XxzTBAFAFA7BKEmSQUSFCmXM3J+BKbGAABoBoJQk764Kog50cF+mBoDAKAxCEKNOfmI7Utj10eLkYIAABqEINSMzDIy/ayw3Z+3NdR0KQAA+g2dEQ2gjEw5I5/hyfk7Y3YMAICGIQg1YOVNWiaQ+Z3x4gMAaB6GRtUtJof9kCBcGSnm0RsEAGgE0ClRq0IZGX9aWO+LUwYBABoLBKFafXhe6O8iGt0CLzsAQGOBoVH12ZZMr+ewq6PwmgMANCL4UFaTlEI2+7IQMURsgpccAKAxwRidOsgomXRGWNiFf8MGhwYBABoXBKE6zLsq2BqSD9vh1QYAaHQwTtfgwjPY7lQWGyhGZxAAoBFCEDaszDLybrSww5+3M9J0KQAAoAoG6xoQI+Rff8uDPUR9sJQaAEBjhSBsQD/E06wy8nVnXtOFAABAtTA02lBic9h3N4VLI8QSfNkAAGjE8CHdIErkZGKUsNqHb2mOQVEAgEYNQdggPjov+DmJJrjj5QUAaOwwNFr/9qTRC1nsGpZSAwDQBviwrmepReyj80LYYLG5RNOlAABADWDsrj7JKZkcJXzVie9si0ODAADaAUFYnxZcE6wMyX+88KoCAGgNDI3Wm+gnbGsyvR4oQWcQAECLaCYIGWM3b95kjHXs2JHjVPSfnj17lpyczPO8h4eHiYmJ+it8XdnlZPIZYVtvsaOxpksBAIDXoYEgLCoqGjBgQHFxsUgkMjExCQ8Pt7CweLHBqVOnxo4d6+npWV5enp6evn379oEDB6q/zppjhPzrb2FKK1F/V/QGAQC0jAaOZq1fv97Y2DguLu7GjRvm5ua//PJLpQY9evTIzs6+cOFCbGzsvHnzZs2apf4iX8vaBPq4hH3jjaXUAAC0jwaCcM+ePdOmTeM4juO4adOm7dmzp1IDY2Njnv9fqLi4uIhEjbqbdSOXLY0TdvfjDTBFBgBAC2lgaPThw4ctW7ZU3G7ZsmV6enrVNlKpdOHChVlZWTdv3tywYUN1TyWVSp88eRIZGanc4uPjY2ZmVu81V0exlNoPb/FuWEoNAEA7NUgQnjhxYsuWLZU28jz/559/EkLKysoMDAwUG42MjEpKSqo+g0gkcnNzMzMzO3/+/Pnz5319fVXu6OnTpwkJCUuXLlVu+fbbb9944w2VjaVSKcdxEkl9nuj+0RVxF2sy0kleXFyPz6o7SkpKGnmHXsdQSsvLyymlmi5Ej6j8BIMG9VofLCYmJiqnZL6oQYKwTZs248aNq7RRWYqTk1NeXp7idm5urrOzc9VnMDAwmDFjBiEkMDCwU6dOISEhlSbUKDRv3jwgICA0NLQmVUkkkvoNwn1p9EIOjQ0Um2ERmWowxtTZQQdKqVgs1oqJ1roEb3I1q/cPlgYJQnd3d3d39+oe7dq167lz54YMGUII+fvvv7t27fqSpzI2Nm6cX2/TS9hHF4RDAVhKDQBAu2ngGOF//vOfAQMGeHh4cBz3008/nTx5UrG9RYsWW7du7d2794YNG4qKijw9PbOzs9euXfv222+r7A5qkJySCaeFOR35txww7gcAoN00EIRvvfXWgQMHfv/9d0LIvn37fHx8FNsnTZrk4uJCCOnatev27dvPnz9vaWk5c+bMSZMmqb/Il/v6mmAuIbM6YJ4oAIDW08zKMn379u3bt2+ljUuWLFHc6NKlS5cuXdReVE0tvCbsv8+ih4nRGQQA0AFYa/Q1MEI+uyREP2Vnh4kdsJQaAIBOQBDWlMDIe38LKYXs9BCxpYGmqwEAgHqCIKwRqUAmRglSgZ0cJDbGawYAoEMw3ePVimVkWLjciCcHApCCAAC6BkH4CnlS0v+EvLWFKLQPL8GrBQCgc/DR/jKPS1nvo/I+zqJffHkOk0QBAHQRgrBaqUWs11HhnTbc8q64vhIAgM7CIS/VYnPYsHD5sq78tNb4rgAAoMsQhCqcfcrGnpJv6MmPao4UBADQcQjCyo6ls3ej5Tv8xQGuOCoIAKD7EIT/sDOFfnJJODRA3B2raQMA6AcE4XPrE+mSGzRysLiDDVIQAEBfIAj/Z/F1ujWZnhvGtzBHCgIA6BEE4f+W0j79mP09XOyEpbQBAPSMvgehwMj754T4PHZqqNjWUNPVAACA2ul1EFZQMilKyJOyU0PEZhJNVwMAAJqgv+fJlcjJ8DC5nJLjA5GCAAD6S0+DMF9KAo7LnU1Ee/rxhlhADQBAj+ljED4tI32OyXs5izb35sX6+AIAAMBzepcDaUXM74h8UitueVce50kAAIB+TZa5XUCGRwoLOnPveerdNwAAAFBJj4LwWi4ZHUV+7sGNaYkUBACA/9GjIKygoh29SL+mSEEAAHhOj4Kwuz3jcJl5AAD4J3SPAABAryEIAQBAr+lREN65cyc1NVXTVeiXsLAwxpimq9Aj2dnZV69e1XQV+uX8+fOFhYWarkKPyGSyU6dO1e9z6lEQ7tixY9++fZquQr+89957BQUFmq5Cj0RHR//000+arkK/LF++HF8+1On+/fuzZ8+u3+fUoyAkhKB3AroN73CAWtCvIAQAAKgEQQgAAHpNu88jfPTo0enTpwMCAmrSODk5WSKRnDlzpoGLgudKSkpGjx4tFmv320yLZGZmZmVl1fB/BNSLuLi4uXPnWltba7oQfVFWVvbkyZOav8kDAwM//PDDl7cRafVBhSdPnpw+fdrR0bEmjfPy8niet7S0bOiqQCktLa1ly5aarkKPlJeX5+fnOzs7a7oQPfLo0SNHR0eJBBc1VRPG2IMHD1q0aFHD9i1btnR3d395G+0OQgAAgDrCMUIAANBrCEIAANBrCEIAANBrCEIAANBrejGvvby8/ObNm4yxt956S2UDxlh0dHRaWlqPHj08PDzUXJ5OevjwYVRUlIODQ0BAQNXTJx49enTnzh3lXR8fHzMzM/UWqAtiYmLi4+O9vLy6du2qssHjx49PnTplbW09YMAAAwMDNZene+RyeURERFZWVt++fZs2bVq1QVRUlCAIittOTk5eXl7qLVAHpaSkpKWldevWzcLCQmWDBw8enDlzxtHRMSAggOf5Wu6G6bqtW7caGBjY2tp26tSpujZTp05t167djBkz7O3td+3apc7ydFJUVJSNjU1wcHC3bt0CAgIEQajUYN26dY6Ojv3/X2pqqkbq1GrLli1r2rRpSEhIs2bNvv3226oNLl++bGNj8+677/r6+vr6+kqlUvUXqUvkcnnfvn3feuut4OBgGxub6Ojoqm3Mzc179uypeFcvW7ZM/UXqkoqKCmtra2tra57nL1++rLJNZGSkjY3N9OnTu3btOmjQoKofNTWk+0GYlZVVUFCwc+fO6oLw5s2blpaWubm5jLGDBw+6ubnV+tUEhZ49e65du5YxVlZW5ubmdvz48UoN1q1bN378eE2UpiPy8/NNTU0TEhIYY3fu3DExMcnLy6vUZvDgwUuWLGGMVVRUdOjQAd/w6ujIkSPu7u7l5eWMsTVr1vTq1atqG3Nz8/v376u9NN0kCEJKSgpjzMLCorog7N69+y+//MIYKy0tbdGiheJyN7Wg+8cI7e3tX34S/dGjR/39/W1sbAghQ4YMefLkya1bt9RVnQ4qKCg4d+7cmDFjCCFGRkbDhg07evSoymYnTpy4fv06pVTtNWq9qKiopk2btm/fnhDi4eHh7u4eGRn5YgOZTBYeHh4UFEQIkUgkI0eOVPlXgJo7evTo8OHDDQ0NCSFjxow5e/asyqsvXb169dSpUzk5OWovUNdwHOfm5vaSBjk5ORcvXlS8yY2NjYcMGVLrN7nuB+ErZWRkNGnSRHFbIpE4ODhkZGRotiSt9vjxY7FY7OTkpLjr6upa9fUUiURZWVnr168fNWqUj49Pbm6u2svUbi++aYmqF/np06eCICjbqPwrwGvJyMhwdXVV3HZycuJ5vupLamtru3nz5v/+978tW7bcQHXtZAAADVxJREFUunWr2mvUL48fP1Z8Yivu1uVNrguTZfbv3//f//636vaYmJiarHIpCMKLyyOJxWK5XF6f9emi2bNnV+qCEEK6deu2ceNGQRBEIpFyI8/zVV/PGTNmfPDBB4SQioqKoUOHLlq0aO3atQ1dsy6p9CJXfdMqpmwo26j8K8BrEQSB4/7XcxCJRBzHVX1J7927p5ivcfz48aCgoOHDhyuGmqAhvPgXIXV7k+tCEPr7+7dp06bq9hrOIHJ2dr59+7biNqU0KyvLxcWlPuvTRTNnznznnXcqbVTM/HRycpLJZPn5+YqPgMzMzKpLXyr/NAYGBkFBQTt27GjwinWLs7NzVlaW8m5mZmalN62iR56dnd28eXOVDeB1vfia5+XlyWSyqi+p8o09ZMgQQ0PDxMREX19ftVapT5ycnKRSaWFhoWJCqcqPmhrShaFRa2trL1Ve/MpcVX5+vkwmI4T06dPnzJkzFRUVhJCLFy8aGRkpDr3ASzRt2rTqC65YBtfe3t7Lyys8PJwQwhiLiIjw9/cnhAiCoHIINDY2VuVMdHiJnj173rlz58mTJ4SQrKys+Ph4Pz8/Qkh5ebniwJWRkZGPj4/ir0AICQ8P79Onj+bq1QV9+vQJDw9njBFCwsPDO3bsaGtrSwgpKioqKyur1PjevXuFhYV4YzeEsrKyoqIiQoiTk5Onp6fiTU4pjYyMVHzU1EYtJ/Rojzt37syYMaNfv362trYzZsxYvXq1Yru9vf2hQ4cUt3v16jVo0KAffvjB3d195cqVmitWR/z5558ODg4rV64cP35827ZtFRPtrl27RghR3J44ceLs2bNXrlw5btw4S0tLxVme8FpmzJjh7e29evXqbt26TZ8+XbFx1apV3bt3V9w+fPiwjY3NihUr3nnnnZYtWxYWFmquWF1QVlbm4eExYcKElStXOjg47N69W7E9ICBg0aJFjLHjx4+PHTt26dKl8+fPd3V1/eCDDzRary6YP3/+jBkzDAwMRo0aNWPGjJycHMbYV199NXToUEWD0NBQJyenVatWjRs3rn379rU+R0j3rz7x+PHjF6cSNWnSZMiQIYSQ7du39+rVq1mzZoSQsrKyLVu2pKen+/r6Dh06VGO16pCzZ8+Gh4fb2tq+++67VlZWhJDc3Nz9+/dPnz6d47hLly5FR0cXFBQ0adJk7NixysPdUHOU0p07dypOqJ8wYYJiUC4+Pj4tLW3EiBGKNpcuXTp27JilpeU777xjZ2en0Xp1QX5+/ubNm/Py8gYNGtSzZ0/FxhMnTjg5OXXu3DkvL+/gwYOpqalGRkbdu3fv16+fZqvVAbt373727Jny7sSJE83MzGJiYnJzcwcOHKjYeObMmYiICHt7+3fffbfWV9nT/SAEAAB4CV04RggAAFBrCEIAANBrCEIAANBrCEIAANBrCEIAANBrCEIAANBrCEKAhlJWVrZnz54XT4TSlJMnTx4+fFhxOzs7e+vWrZmZmeovo7y8fM+ePQUFBerfNcBLIAgBGsr333+/aNEic3NzTRdC1q1bt2rVKsXte/fuvfPOO3fv3q370yYmJm7cuFEqldawvZGR0bp161QukQ+gQQhCgAaRm5v73XffzZ8//8UF8hsDxRXtFQvD1tHZs2fff//90tLSmv/IvHnz1q1bd//+/brvHaC+NK7/ogA6Y/PmzSKRaOTIkZW2M8aysrJUXtNVKpU+ffq0vLxc5RMqHn1J90vRQHEBJkJIRUVFZmZm1QvTuLq6zp8/X7G44IsEQcjOzq5uqamsrCzFwvQ1JJfLnz59WlJSUml7//79HR0dN23aVPOnAmhoCEKAGlm6dKm9vf2VK1cUdxlj48ePb9WqVXUH27Zt2zZ06FBjY2PlFkEQvv32WycnJ0dHR0tLyyZNmuzatUvxUGFh4TvvvGNlZeXs7GxlZTVhwoS8vDzlDxYUFEyePNnS0lLx6JQpU5SH2crLy21sbNasWRMSEqJooLim2PLly+3s7JycnBwcHCpd6/HatWvOzs6XL19W3B0+fPiYMWNCQ0NdXV0dHBwsLCwWLFigbJycnDxw4EBjY2NHR0dTU9POnTtHR0crHvrxxx9nzZpFCHFzc7OxsbGxsVFcll0qlc6aNcvOzs7Z2dnCwiIgIODF/h/HcSNHjsRFa6FxqY8lwgF0n0wm8/X1bdasWW5uLmNs7dq1IpHowIEDKhtnZWWJRKK1a9e+uHHGjBkcx82aNevSpUuxsbEbN27csmWL4qGAgAAjI6Off/755s2bGzduNDU17dGjhyAIjDFKaa9evUxMTNavX3/z5s1ffvnF2Ni4d+/elFLGmGJY0sHBYfjw4SdPnoyMjMzOzt64cSMh5MMPP7xx40Z4eLinp6etra2fn59iXxcuXCCEREdHK+726dPHwcHhjTfeOHTo0OXLl4ODgwkhp06dUjx68eLFWbNmRUZGJiYmRkZG9u7d28zMLCMjgzGWlpb2ySefEEIOHDgQERERERGhWPt/1KhRVlZWv/zyS0JCwsmTJ728vFq1alVaWqp8HRSXn7xz5069/W0A6gZBCFBT6enptra2w4cPv3HjhpGR0aefflpdy8jISELIyZMnlVsSExNFIpHKH1Ek04vX//r1118JIcePH2eMRUVFEUJezNTVq1crs0oRhO3atVOkpkKLFi169eqlvHvr1i2RSPSSIDQ2Nk5PT1fclUql9vb2//73v1X+Xs+ePTMwMFi/fv2Ldebl5SkbKPqLO3fuVG5JSkoSiUTbt29XbomJiSGEKC9jBKBxunCFegD1aNKkyebNm0eOHBkdHd2hQ4fly5dX11IxSGhjY6Pcosit6dOnV21848YNQsjbb7+t3PL222+HhIScPXt28ODBVR8dP378J598cvbs2b59+yq2DB8+XDklp6Cg4P79+59++qmyfbt27dq1a/eS38vLy6tJkyaK2wYGBq1atUpPT1c+mp2dvXv37tTUVMUBP0NDw3v37lX3VGFhYSKRyNTUVPFVQMHGxiYhIUF5V3E9W+XV3gE0DkEI8BoGDx7cqlWr5OTkr776ysDAoLpmYrGYEPLiRBVFNCrz5kUPHz4UiUTOzs7KLVZWViYmJrm5uYpHxWKxvb298lEHh/9r7/5e2fvjOIC/t9mHHcJkP9BEKPNzkkTb/ChpUeIPsAsWV4rIhQtXI0pKKSUSd5QSk6aoFaWE5mplNltrOyxmfmbnOJ+Ld9/TifH19f36zDevx9Xp/XqfH9suntvO65wjjYqKwlVMJpOxyx6PByEkl8u5u0hJSXmny4Yb2Aih6Ohoti/GbDY3NzenpaVptdqkpCQ+ny8QCMJ2+mAkSfJ4PL1e/2KcezElfluEQuFbGwHgD4MgBOAfGBwcdDqdubm5fX19tbW1b10jiJOJm1X46cQkSb5+dihBEAzD+P1+Nr3u7u7u7+/xzNjYWIqiAoGAWCzG1cvLS4qiuNvh8XjsMn7KMc5d1sXFRXx8/Cde7/DwcHFxscViwdHOMMzk5OQ78xMSEvh8vsfj4XYJvYCP7UVUAxBB0DUKwEdtb2+PjIwYjUaTyXR+ft7e3v7WTJVKJRQKj4+P2RGtVosQWlxcfD25oqICIbS2tsaOrK6usuPl5eUIIZPJFLb6mkQikclkGxsb7IjL5cKtpJ/gcDhUKhVOQYSQxWK5vb1lq/h7wMPDAztSVVVFUdTy8vI727RarQihsrKyzx0SAP+9CJ+jBOB/wufzpaSk6HQ63K65tLSEEJqenn5rvkajaWxs5I60tLSIRKKJiQmv13t1dbW5ubmyssIwDE3TpaWlycnJJpMpGAyazWa5XK5UKp+enhiGCYVChYWFONuCweD6+rpEIikqKqIoivmrWWZ8fJy7o6GhIT6fPzo6GggEbDabWq0mCOKdZpn6+nru6tXV1TqdDi83NTVJpdK9vb3Hx8etra2srCyRSGQwGHD18PAQIdTb27uzs7O/vx8KhWiaVqvVYrF4ZmaGJMlgMHhwcDAwMLC7u8tuX6/X5+XlfeIjAOCLQBAC8Pdomq6rq5PJZD6fjx3s6OiIiYk5OjoKu8rc3NyvX7/8fj87cnd319bWxv66EolEU1NTuOTxeGpqativp5WVlQ6Hg13R5XJpNBq2qtVqXS4XLoUNQoqiOjs7cfuMQCDo7u5uaGj4XBCenJzk5+fj/cbFxc3OzmZnZ7NByDCM0WhUKBQCgQAhdH5+zjBMIBBobW1lXyZCqKSkxGq14vkPDw+JiYljY2Mfet8B+CN4zBs3kgAA/BuPj485OTk9PT3cBk6E0PX1tc1mIwgiIyMjLi6OW3K73V6vVyqVhr3/2dnZGUmScrn89U1hwvJ6vW63OzMzk9tog9E0jaPrI2iattvtt7e3SqXynTN/L9zc3NhstqioKIVCgdtEsYWFha6urtPTU/aUJwARB0EIwFeZn5/v7++32+0EQUT6WL4FmqYLCgoMBgO+JQ0A3wQEIQBf5fn52el0pqamxsTERPpYvoVQKOR2u9PT07l/nAIQcRCEAAAAfjS4fAIAAMCPBkEIAADgR4MgBAAA8KNBEAIAAPjRfgOlTdPWogMVFgAAAABJRU5ErkJggg==",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 13
}
],
"cell_type": "code",
"source": [
"Plots.plot(getindex.(points, 1), getindex.(q_points, 1), xlabel = \"x (coordinate)\", ylabel = \"q_x (flux in x-direction)\", label = nothing)"
],
"metadata": {},
"execution_count": 13
},
{
"cell_type": "markdown",
"source": [
"*Figure 4*: $x$-component of the flux along the cut line from *Figure 2*."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"---\n",
"\n",
"*This notebook was generated using [Literate.jl](https://github.com/fredrikekre/Literate.jl).*"
],
"metadata": {}
}
],
"nbformat_minor": 3,
"metadata": {
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.11.6"
},
"kernelspec": {
"name": "julia-1.11",
"display_name": "Julia 1.11.6",
"language": "julia"
}
},
"nbformat": 4
}