{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Aula 4\n",
    "\n",
    "# Soluções Numéricas\n",
    "\n",
    "Frequentemente, em situações reais, os problemas de engenharia são muito complexos para que uma solução analítica (uma equação explícita). Nesses casos geralmente utilizamos soluções numéricas. De forma geral essas soluções consistem em dividir um problema real em pequenas partes (discretização) e aplicar os princípios físicos que regem o problema às essas pequenas partes. Nessa aula construiremos soluções numéricas para problemas de magnetismo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "mi0 = 4*np.pi*1e-7 #[H/m]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Trataremos do problema de um dipolo magnético, conforme a seção 8.4 do livro Física Básica - Eletromagnetismo de Alaor Chaves. Esse problema consiste no cálculo do campo magnético ao longo do eixo Z gerado por uma espira de raio $R$ onde circula uma corrente $I$, conforme figura abaixo:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"problema_aula_4.png\" alt=\"drawing\" style=\"width: 25%;\"/>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Para esse caso específico, uma solução analítica é possível e é dada por:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "B = \\frac{\\mu_0}{2} \\frac{IR^2}{(R^2+z^2)^{3/2}} \n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Que podemos implementar em Python facilmente:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def CalcularB(I, R, z):\n",
    "    return (mi0*I*np.power(R,2))/(2.0*np.power(np.power(R, 2.0) + np.power(z, 2.0),3.0/2.0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "E podemos testá-la com valores arbitrários"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.1848404026933503e-07"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "R = 5 # [m] Raio da espira condutora\n",
    "I = 1 # [A] corrente que circula na espira\n",
    "z_pos = 1 # [m] posicao em Z para o qual o valor de campo magnetico sera calculado\n",
    "\n",
    "CalcularB(I,R,z_pos)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lembrando que o valor que essa função retorna é referente a componente ao longo do eixo Z. AS componentes $î$ e $ĵ$ são nulas."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Para desenvolver um método numérico que resolve o mesmo tipo de problema, precisamos primeiramente discretizar o sistema físico. Começaremos representando o círculo por um conjunto de pontos igualmente espaçados ao longo de seu perímetro."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_pts = 30 # numero de pontos para discretizacao da espira\n",
    "circle = np.zeros((num_pts,3)) # matriz de coordenadas x,y,z dos pontos que representam a espira\n",
    "angles = np.linspace(start=0, stop=2*np.pi, num=num_pts) # os angulos para os quais calcularemos os pontos\n",
    "\n",
    "for idx, angle in enumerate(angles):\n",
    "    circle[idx,0] = R*np.cos(angle)\n",
    "    circle[idx,1] = R*np.sin(angle)\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "No código acima nos criamos uma matriz nula de tamanho $num\\_pts\\times 3$. Em seguida é gerada uma lista de ângulos para os quais queremos calcular a posição $(x,y)$ do cículo, que correspondem às colunas 0 e 1, respectivamente. Para isso utilizamos as seguintes equações: "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "    x = R cos(\\theta)\\\\\n",
    "    y = R sin(\\theta)\\\\\n",
    "    z = 0\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmMAAAJQCAYAAADCJmN+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3X9w7fld1/HX22yUwy8Ds0VMtsviiFHkh9FLqwMIAiUVaomjOIAgA8iCI9IONIVQhj8cBu1kBFQYsRSRmVaQGdLwYwZDK1RB6cpu0xKgBIGB0pMybAcz5ceZNk0//nGTy9692XuTu0k+58fjMXOnOZ+cm/PeZjp97vfz/VGttQAA0Mef6j0AAMAsE2MAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6eqD3ABfx4IMPtkceeaT3GAAA9/TEE0+8q7X2nHu9b6Ji7JFHHsnjjz/eewwAgHuqqt8+z/tsUwIAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjh7oPQDAJNneHWZzZz8Hh6MsLgyyvrqctZWl3mMBE0yMAZzT9u4wG1t7GR0dJ0mGh6NsbO0liSAD7pttSoBz2tzZvxVip0ZHx9nc2e80ETANHBkDLt20buUdHI4utA5wHo6MAZfqdCtveDhKy59s5W3vDnuP9qwtLgwutA5wHmIMuFTTvJW3vrqcwfzcbWuD+bmsry53mgiYBrYpgUs1zVt5p1ut07gFC/QjxoBLtbgwyPCM8JqWrby1lSXxBVwq25TApbKVB3AxjowBl8pW3mSa1itgYRKIMeDS2cqbLG5mC33ZpgSYcdN8BSxMAjEGMOOm+QpYmARiDGDGuZkt9CXGAGacK2ChLyfwA8w4V8BCX2IMAFfAQkdiDDpzfyeA2SbGoCP3dwLACfzQkfs7ASDGoCP3dwJAjEFH7u8EgBiDjtzfCQAn8ENH7u8EgBiDztzfCWC22aYEAOhIjAEAdCTGAAA6EmMAAB05gR+AqeJ5r0waMQbA1PC8VyaRbUoApobnvTKJxBgAU8PzXplEYgyAqeF5r0wiMQbA1PC8VyaRE/gBmBqe98okEmMATBXPe2XS2KYEAOhIjAEAdCTGAAA6EmMAAB2JMQCAjlxNycTzUGAAJpkYY6J5KDAAk842JRPNQ4EBmHRijInmocAATDoxxkTzUGAAJp0YY6J5KDAAk84J/Ew0DwUGYNJ1j7GqmkvyeJJha+1Fvedh8ngoMACTbBy2KV+S5G29hwAA6KFrjFXVQ0k+N8mre84BANBL7yNj35nk5Une/0xvqKpHq+rxqnr8ySefvL7JAACuQbcYq6oXJfm91toTd3tfa+1VrbUbrbUbz3nOc65pOgCA69HzyNgnJ3lxVf1Wkh9K8hlV9ZqO8wAAXLtuMdZa22itPdRaeyTJFyT56dbaF/eaBwCgh97njAEAzLTu9xlLktbaG5O8sfMYAADXzpExAICOxuLIGADMuu3doUe7zSgxBgCdbe8Os7G1l9HRcZJkeDjKxtZekgiyGWCbEgA629zZvxVip0ZHx9nc2e80EddJjAFAZweHowutM13EGAB0trgwuNA600WMAUBn66vLGczP3bY2mJ/L+upyp4m4Tk7gB4DOTk/SdzXlbBJjADAG1laWxNeMsk0JANCRGAMA6EiMAQB0JMYAADoSYwAAHYkxAICOxBgAQEdiDACgIzEGANCRGAMA6MjjkHhG27tDz0kDgCsmxjjT9u4wG1t7GR0dJ0mGh6NsbO0liSADgEtkm5Izbe7s3wqxU6Oj42zu7HeaCACmkxjjTAeHowutAwD3R4xxpsWFwYXWAYD7I8Y40/rqcgbzc7etDebnsr663GkiAJhOTuDnTKcn6buaEgCulhjjGa2tLIkvALhitikBADoSYwAAHYkxAICOxBgAQEdiDACgIzEGANCRGAMA6EiMAQB0JMYAADoSYwAAHYkxAICOxBgAQEceFA4AnGl7d5jNnf0cHI6yuDDI+upy1laWeo81dcQYAHCH7d1hNrb2Mjo6TpIMD0fZ2NpLEkF2yWxTAgB32NzZvxVip0ZHx9nc2e800fQSYwDAHQ4ORxda5/6JMQDgDosLgwutc//EGABwh/XV5Qzm525bG8zPZX11udNE08sJ/ADAHU5P0nc15dUTYwDAmdZWlsTXNbBNCQDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHD/QegDtt7w6zubOfg8NRFhcGWV9dztrKUu+xAIArIMbGzPbuMBtbexkdHSdJhoejbGztJYkgA4ApZJtyzGzu7N8KsVOjo+Ns7ux3mggAuEpibMwcHI4utA4ATDYxNmYWFwYXWgcAJpsYGzPrq8sZzM/dtjaYn8v66nKniQCAq+QE/jFzepK+qykBYDaIsTG0trIkvgBgRtimBADoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKCjB3oPAABwant3mM2d/RwcjrK4MMj66nLWVpZ6j3WlxBgAMBa2d4fZ2NrL6Og4STI8HGVjay9JpjrIbFMCAGNhc2f/VoidGh0dZ3Nnv9NE10OMAQBj4eBwdKH1aSHGAICxsLgwuND6tOgWY1X13Kr6map6W1X9clW9pNcsAEB/66vLGczP3bY2mJ/L+upyp4muR88T+N+X5Otba2+uqg9J8kRVvb619isdZwIAOjk9Sd/VlNektfbOJO88+foPquptSZaSiDEAmFFrK0tTH19PNxbnjFXVI0lWkjx2xvcerarHq+rxJ5988rpHAwC4Ut1jrKo+OMmPJHlpa+3dT/9+a+1VrbUbrbUbz3nOc65/QACAK9Q1xqpqPjdD7LWtta2eswAA9NDzaspK8n1J3tZa+/ZecwAA9NTzyNgnJ/mSJJ9RVW85+fM5HecBALh2Pa+m/Lkk1evzAQDGQfcT+AEAZpkYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDoSIwBAHTU7dmUk257d5jNnf0cHI6yuDDI+upy1laWeo8FAEwYMXYftneH2djay+joOEkyPBxlY2svSQQZAHAhtinvw+bO/q0QOzU6Os7mzn6niQCASSXG7sPB4ehC6wAAz0SM3YfFhcGF1gEAnokYuw/rq8sZzM/dtjaYn8v66nKniQCASeUE/vtwepK+qykBgGdLjN2ntZUl8QUAPGu2KQEAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDoSIwBAHQkxgAAOhJjAAAdPdB7AACAy7S9O8zmzn4ODkdZXBhkfXU5aytLvcd6RmIMAJga27vDbGztZXR0nCQZHo6ysbWXJGMbZLYpAYCpsbmzfyvETo2OjrO5s99ponsTYwDA1Dg4HF1ofRyIMQBgaiwuDC60Pg7EGAAwNdZXlzOYn7ttbTA/l/XV5U4T3ZsT+AGAqXF6kr6rKQEAOllbWRrr+Ho625QAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjh7oPcA42d4dZnNnPweHoywuDLK+upy1laXeYwEAU0yMndjeHWZjay+jo+MkyfBwlI2tvSQRZADAlbFNeWJzZ/9WiJ0aHR1nc2e/00QAwCwQYycODkcXWgcAuAxi7MTiwuBC6wAAl0GMnVhfXc5gfu62tcH8XNZXlztNBADMAifwnzg9Sd/VlADAdbprjFXVu+/x9yvJO1trf+nyRupnbWVJfAEA1+peR8Z+o7W2crc3VNXuJc4DADBT7nXO2D84x884z3sAADjDXY+MtdZ+86mvq+pDn/p3Wmu///T3AABwfuc6gb+qvirJv0wyStJOlluSv3BFcwEAzITzXk35siR/tbX2rqscBgBg1pz3PmO/keSPr3IQAIBZdN4jYxtJ/ndVPZbkPaeLrbWvvZKpAABmxHlj7D8m+ekke0nef3XjAADMlvPG2Ptaa193pZMAAMyg854z9jNV9WhV/fmq+vDTP1c6GQDADDjvkbEvOvnPjaesubUFAMCzdK4Ya6199FUPAgAwi+66TVlVf/1eP+A87wEA4Gz3OjL2/VX16UnqLu/5viR3fZg4AABnu1eM/dkkT+TuMfbk5Y0DADBb7vWg8EeuaQ4AgJl03ltbAABwBcQYAEBH97qa8kpjrapeWFX7VfXrVfWNV/lZAADj6F4n8L+5qv5Za+3nL/uDq2ouyXcneUGSdyT5har6sdbar1z2ZwEAfPP2Xn7wsd/JcWuZq8oXPv+5+da1j+891j23Kb8qyb+tqu+tqg+75M9+XpJfb639ZmvtvUl+KMnnXfJnAADkm7f38po3vT3HrSVJjlvLa9709nzz9l7nye4RY621x5I8P8mbkzxeVd9VVf/u9M+z/OylJL/zlNfvOFkDALhUP/jY71xo/Tqd53FIH57kk3LzfmJPJHn/JX32Wfcua3e8qerRJI8mycMPP3xJHw0AzJLTI2LnXb9Od42xqvrqJOtJNpN8RWuXOvE7kjz3Ka8fSnLw9De11l6V5FVJcuPGjf7/jQEAE2eu6szwmqu73df+etzrnLFPTfK3Wmvfc8khliS/kORjquqjq+pPJ/mCJD92yZ8BAJAvfP5zL7R+ne51B/5/fFUf3Fp7X1V9TZKdJHNJ/lNr7Zev6vMAgNl1etXkOF5NWZd/wOvq3Lhxoz3++OO9xwAAuKeqeqK1duNe73MHfgCAjs5zNWWSpKo+MTfPIUuSn22tvfVqRgIAmB3nOjJWVS9J8tokH3Hy5zVV9S+ucjAAgFlw3iNjX5Hk+a21P0qSqnplkp9P8u+vajAAgFlw3nPGKsnxU14f5+ybtgIAcAHnPTL2/Ukeq6rXnbxeS/J9VzMSAMDsOFeMtda+varemORTcvOI2Je11navcjAAgFlwrhirqs9qrb0hNx8Yfrr2pa21H7iyyQAAZsB5zxn7lqr6D1X1gVX156rqx5P8vascDABgFpw3xj4tyW8keWuSn0vyX1pr//DKpgIAmBHnjbEPS/L83Ayy9yT5qKoxeMw5AMCEO2+MvSnJT7bWXpjkk5IsJvlfVzYVAMCMOO+tLT6rtfb2JGmtjZJ8bVX97asbCwBgNtw1xqrqL7fWfjXJg1X14NO+/YdXNxYAwGy415Gxr0/ylUn+zRnfa0k+49InAgCYIXeNsdbaV57859+5nnEAAGbLXU/gr6qXP+Xrz3/a977tqoYCAJgV97qa8gue8vXG0773wkueBQBg5twrxuoZvj7rNQAAF3SvE/jbM3x91uuJt707zObOfg4OR1lcGGR9dTlrK0u9xwIApti9YuwTq+rduXkUbHDydU5ef8CVTnbNtneH2djay+joOEkyPBxlY2svSQQZAHBl7rpN2Vqba619aGvtQ1prD5x8ffp6/rqGvA6bO/u3QuzU6Og4mzv7nSYCAGbBeR+HNPUODkcXWgcAuAxi7MTiwuBC6wAAl0GMnVhfXc5gfu62tcH8XNZXlztNBADMgvM+KHzqnZ6k72pKAOA6ibGnWFtZEl8AwLWyTQkA0JEYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDo6IHeAwAAXKbt3WE2d/ZzcDjK4sIg66vLWVtZ6j3WMxJjAMDU2N4dZmNrL6Oj4yTJ8HCUja29JBnbILNNCQBMjc2d/Vshdmp0dJzNnf1OE92bGAMApsbB4ehC6+NAjAEAU2NxYXCh9XEgxgCAqbG+upzB/Nxta4P5uayvLnea6N6cwA8ATI3Tk/RdTQkA0MnaytJYx9fT2aYEAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdPRA7wEm1fbuMJs7+zk4HGVxYZD11eWsrSz1HgsAmDBi7D5s7w6zsbWX0dFxkmR4OMrG1l6SCDIA4EJsU96HzZ39WyF2anR0nM2d/U4TAQCTSozdh4PD0YXWAQCeiRi7D4sLgwutAwA8EzF2H9ZXlzOYn7ttbTA/l/XV5U4TAQCTygn89+H0JH1XUwIAz5YYu09rK0viCwB41mxTAgB0JMYAADoSYwAAHYkxAICOxBgAQEdiDACgIzEGANCRGAMA6EiMAQB0JMYAADoSYwAAHXWJsararKpfrapfrKrXVdVCjzkAAHrrdWTs9Uk+rrX2CUl+LclGpzkAALrqEmOttZ9qrb3v5OWbkjzUYw4AgN7G4ZyxL0/yk72HAADo4YGr+sFV9YYkH3nGt17RWvvRk/e8Isn7krz2Lj/n0SSPJsnDDz98BZMCAPRzZTHWWvusu32/qr40yYuSfGZrrd3l57wqyauS5MaNG8/4PgCASXRlMXY3VfXCJN+Q5NNaa3/cYwYAYPxs7w6zubOfg8NRFhcGWV9dztrKUu+xrlSXGEvyXUn+TJLXV1WSvKm19tWdZgEAxsD27jAbW3sZHR0nSYaHo2xs7SXJVAdZlxhrrf3FHp8LAIyvzZ39WyF2anR0nM2d/amOsXG4mhIAIAeHowutTwsxBgCMhcWFwYXWp4UYAwDGwvrqcgbzc7etDebnsr663Gmi69HrBH4AgNucnhfmakoAgE7WVpamPr6ezjYlAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB090HsA7rS9O8zmzn4ODkdZXBhkfXU5aytLvccCAK6AGBsz27vDbGztZXR0nCQZHo6ysbWXJIIMAKaQbcoxs7mzfyvETo2OjrO5s99pIgDgKomxMXNwOLrQOgAw2cTYmFlcGFxoHQCYbGJszKyvLmcwP3fb2mB+Luury50mAgCukhP4x8zpSfqupgSA2SDGxtDaypL4AoAZYZsSAKAjMQYA0JEYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMAKAjMQYA0JEYAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI4e6D0AADCetneH2dzZz8HhKIsLg6yvLmdtZan3WFNHjAEAd9jeHWZjay+jo+MkyfBwlI2tvSQRZJfMNiUAcIfNnf1bIXZqdHSczZ39ThNNLzEGANzh4HB0oXXunxgDAO6wuDC40Dr3T4wBAHdYX13OYH7utrXB/FzWV5c7TTS9nMAPANzh9CR9V1NePTEGAJxpbWVJfF0D25QAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMPCucZbe8Os7mzn4PDURYXBllfXfbAWAC4ZGKMM23vDrOxtZfR0XGSZHg4ysbWXpIIMgC4RLYpOdPmzv6tEDs1OjrO5s5+p4kAYDqJMc50cDi60DoAcH/EGGdaXBhcaB0AuD9ijDOtry5nMD9329pgfi7rq8udJgKA6eQEfs50epK+qykB4GqJMZ7R2sqS+AKAK2abEgCgIzEGANCRGAMA6EiMAQB0JMYAADoSYwAAHYkxAICOxBgAQEdiDACgIzEGANCRxyEBwBjY3h16HvCMEmMA0Nn27jAbW3sZHR0nSYaHo2xs7SWJIJsBtikBoLPNnf1bIXZqdHSczZ39ThNxncQYAHR2cDi60DrTRYwBQGeLC4MLrTNdxBgAdLa+upzB/Nxta4P5uayvLneaiOvkBH4A6Oz0JH1XU84mMQYAY2BtZUl8zSjblAAAHXWNsap6WVW1qnqw5xwAAL10i7Gqem6SFyR5e68ZAAB663lk7DuSvDxJ6zgDAEBXXWKsql6cZNhae2uPzwcAGBdXdjVlVb0hyUee8a1XJPmmJJ99zp/zaJJHk+Thhx++tPkAAMZBtXa9u4RV9fFJ/nuSPz5ZeijJQZLntdZ+925/98aNG+3xxx+/4gkBAJ69qnqitXbjXu+79vuMtdb2knzE6euq+q0kN1pr77ruWQAAenOfMQCAjrrfgb+19kjvGZhs27tDjxABYGJ1jzF4NrZ3h9nY2svo6DhJMjwcZWNrL0kEGQATwTYlE21zZ/9WiJ0aHR1nc2e/00QAcDFijIl2cDi60DoAjBsxxkRbXBhcaB0Axo0YY6Ktry5nMD9329pgfi7rq8udJgKAi3ECPxPt9CR9V1MCMKnEGBNvbWVJfAEwsWxTAgB0JMYAADoSYwAAHYkxAICOxBgAQEeupgRgqmzvDt3uhokixgCYGtu7w2xs7d16Zu3wcJSNrb0kEWSMLduUAEyNzZ39WyF2anR0nM2d/U4Twb2JMQCmxsHh6ELrMA7EGABTY3FhcKF1GAdiDICpsb66nMH83G1rg/m5rK8ud5oI7s0J/ABMjdOT9F1NySQRYwBMlbWVJfHFRLFNCQDQkRgDAOhIjAEAdCTGAAA6cgI/dOY5egCzTYxBR56jB4BtSujIc/QAEGPQkefoASDGoCPP0QNAjEFHnqMHgBP4oSPP0QNAjEFnnqMHMNvEGADudwcdiTGAGed+d9CXE/gBZpz73UFfYgxgxrnfHfQlxgBmnPvdQV9iDGDGud8d9OUEfuDSuTJvsrjfHfQlxoBL5cq8yeR+d9CPbUrgUrkyD+BiHBkDLtW0X5lnCxa4bI6MAZdqmq/MO92CHR6O0vInW7Dbu8PeowETTIwBl2qar8yzBQtcBduUwKWa5ivzpn0LFuhDjAGXblqvzFtcGGR4RnhNwxYs0I9tSoBzmuYtWKAfR8YAzmmat2CBfsQYwAVM6xYs0I9tSgCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQkRgDAOhIjAEAdCTGAAA6EmMAAB2JMQCAjsQYAEBHYgwAoCMxBgDQUbXWes9wblX1ZJLf7j3HhHgwybt6D8Ed/F7Gj9/JePJ7GT9+Jxf3Ua2159zrTRMVY5xfVT3eWrvRew5u5/cyfvxOxpPfy/jxO7k6tikBADoSYwAAHYmx6fWq3gNwJr+X8eN3Mp78XsaP38kVcc4YAEBHjowBAHQkxmZAVb2sqlpVPdh7lllXVZtV9atV9YtV9bqqWug90yyrqhdW1X5V/XpVfWPveWZdVT23qn6mqt5WVb9cVS/pPRN/oqrmqmq3qn6i9yzTRoxNuap6bpIXJHl771lIkrw+yce11j4hya8l2eg8z8yqqrkk353k7yb52CRfWFUf23eqmfe+JF/fWvsrSf5mkn/udzJWXpLkbb2HmEZibPp9R5KXJ3Fy4Bhorf1Ua+19Jy/flOShnvPMuOcl+fXW2m+21t6b5IeSfF7nmWZaa+2drbU3n3z9B7n5f/xLfaciSarqoSSfm+TVvWeZRmJsilXVi5MMW2tv7T0LZ/ryJD/Ze4gZtpTkd57y+h3xf/xjo6oeSbKS5LG+k3DiO3PzX+zf33uQafRA7wF4dqrqDUk+8oxvvSLJNyX57OudiLv9TlprP3rynlfk5pbMa69zNm5TZ6w5gjwGquqDk/xIkpe21t7de55ZV1UvSvJ7rbUnqurTe88zjcTYhGutfdZZ61X18Uk+Oslbqyq5uR325qp6Xmvtd69xxJnzTL+TU1X1pUlelOQzm3vL9PSOJM99yuuHkhx0moUTVTWfmyH22tbaVu95SJJ8cpIXV9XnJPmAJB9aVa9prX1x57mmhvuMzYiq+q0kN1prHvLaUVW9MMm3J/m01tqTveeZZVX1QG5eRPGZSYZJfiHJF7XWfrnrYDOsbv6b4w8k+f3W2kt7z8OdTo6Mvay19qLes0wT54zB9fquJB+S5PVV9Zaq+p7eA82qkwspvibJTm6eKP7DQqy7T07yJUk+4+R/H285ORoDU82RMQCAjhwZAwDoSIwBAHQkxgAAOhJjAAAdiTEAgI7EGABAR2IMGGtVdfyUe069paq+8WT91VX1sc/yZ39IVf1GVX3Myev5qtqrquffZY7FC/z8T62qX6mqX3o2cwLTzX3GgLFWVX/YWvvgK/z5/yjJP22tfXZVbSR5pLX2VZc1x8kDr3+itfZxz3pYYCo5MgZMpKp6Y1XdqKqPqqr/W1UPVtWfqqqfrarPPnnP11XVL538OfPxOq21H07y/qp6eZKvTrJxzs//w6p6ZVU9UVVvqKrnncz0m1X14kv7BwWmngeFA+NuUFVvecrrf9Va+6+nL1prv11Vr0zyPUkeS/IrrbWfqqq/keTLkjw/SSV5rKr+R2tt94zPeGluPhLp0dba759zrg9K8sbW2jdU1euSfGuSFyT52Nx8vuKPXewfE5hVYgwYd6PW2l+72xtaa6+uqs/PzSNbp+/9lCSva639UZJU1VaST01yVoy9MMk7k1xkK/G9Sf7bydd7Sd7TWjuqqr0kj1zg5wAzzjYlMPGq6gOTPHTy8vS8rjrn311M8rVJnpfkc6rqE875sUftT066fX+S9yRJa+398S+6wAWIMWAavDLJa5N8S5LvPVn7n0nWquoDq+qDkvz9JD97xt/9jiTf1lp7R5KvS/LdVXWukAO4DGIMGHeDp93a4l8/9ZtV9WlJPinJK1trr03y3qr6stbam5P85yT/JzfPJXv1088Xq6oXJHk4yfclSWvtx5P8vyT/5Kr/oQBOubUFwDm4tQVwVRwZAzifd9/PTV+T/HiSd13dWMCkc2QMAKAjR8YAADoSYwAAHYkxAICOxBgAQEeV8+wpAAAAC0lEQVRiDACgo/8PrtDp1Z0P1q8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7ff1b191a6d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,10))\n",
    "plt.scatter(circle[:,0], circle[:,1]) # Plot do tipo scatter para visualisarmos os pontos \n",
    "plt.xlabel('Eixo X [m]')\n",
    "plt.ylabel('Eixo Y [m]')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Com essa estrutura já podemos, por exemplo, calcular o perímetro da espira numericamente. Nesse caso devemos implementar a seguinte integral:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "        P = \\int_{C}|d\\overrightarrow{L}|\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "No domínio discreto, a integral acima vira um simples somatório:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "        P = \\sum_{n=1}^{k-1} |c[n+1]-c[n]|\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Onde $c[n]$ é o ponto na n-ésima linha da matriz de pontos do círculo, $k$ é o número de pontos utilizados para a discretização e $c[n+1]-c[n]$ é a aproximação discreta para  $d\\overrightarrow{L}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "O perimetro calculado eh 31.3545153 e o valor teorico eh de 31.4159265 \n"
     ]
    }
   ],
   "source": [
    "perimeter = 0 # Acumulador para o perimetro\n",
    "dL = np.zeros(circle.shape) # guardaremos os dLs para uso futuro\n",
    "for idx in range(num_pts-1):\n",
    "    dl = circle[idx+1,:] - circle[idx,:]\n",
    "    dL[idx,:] = dl # armazena na matriz de dLs\n",
    "    perimeter = perimeter + np.linalg.norm(dl) # acumula a norma de dl\n",
    "print(\"O perimetro calculado eh %.7f e o valor teorico eh de %.7f \" % (perimeter, 2*np.pi*R))    \n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lembrando que quanto maior o número de pontos utilizados na discretização, melhor será a aproximação. Altere a variável $num\\_pts$ e veja o impacto no valor da aproximação."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora que temos os pontos da espira e uma matriz com todos os $dL$s que precisamos, conforme figura abaixo, podemos voltar ao problema inicial."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlUAAAJCCAYAAADp1TKRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzs3Xd4VGXi/v/7pECGUBKqJgaR3l0ggIBYQIyAaAAFEVZXFESxSxELNr4Ug+zvI2BBdhUBC12REhEWAQUkmNCUrpQEEIGElp7z+4NlliEVcjJnMvN+XRfXcp7zzMy9qHDznDPPMUzTFAAAAIrHz+4AAAAA3oBSBQAAYAFKFQAAgAUoVQAAABagVAEAAFiAUgUAAGABShUAAIAFKFUAAAAWoFQBAABYIMCOD61atapZq1YtOz4aAADgimzevPkv0zSrFTbPllJVq1YtxcXF2fHRAAAAV8QwjANFmcflPwAAAAtQqgAAACxAqQIAALAApQoAAMAClCoAAAALUKoAAAAsQKkCAACwAKUKAADAApQqAAAAC1CqAAAALECpAgAAsAClCgAAwAKUKgAAAAtQqgAAACxAqQIAALAApQoAAMAClCoAAAALUKoAeI39+/dr9+7ddscA4KMoVQBKtRMnTuiDDz5Qhw4d1KhRI1WsWNHuSAB8VIDdAQDgSqWlpenbb7/VrFmztHTpUmVmZkqS/v73v+uaa66xOR0AX0WpAlAq5OTkaO3atZo5c6bmzZunlJSUXHOeffbZIr/fovhExcTuUlJyqsJCHBoe1UDRLcKtjAzAx1CqAHi0HTt2aNasWZo9e7YOHTqU77x69erp3LlzWrJkiU6fPp3njzNnzuj06dM6cOQv/XH0LxllgnXN3ycqMTlVoxZskySKFYCrRqkC4HGOHDmiL774QjNnzlRCQkKRXrNnzx7deuutV/Q5fo7/3X+VmpmtmNhdlCoAV41SBfiQVatWadu2bUW+TObOS2RnzpzRwoULNWvWLK1cuVI5OTkl8jmXykk/L9M0ZRiGJCkpObXEP/Nqmaapv/76S5UqVVKZMmXsjgMgD5QqwAf89ddfGjZsmGbMmKHp06cXOj8xMVHjP/5KM+Z+LbNseVW588kSuURmmqa+++47ffbZZ1q4cKFSU91canKypOxMKeBCSQkLcbj38y9z7tw5/f7779q/f79+//33XD9uvvlmLVmyxNaMAPJHqQK8mGmamjVrll544QX99ddfkqSWLVvmmpeWlqa1a9cqNjZWsbGx2r59u/OcX1AFVb7jcRl+/pZfIjMMQ23btlVWVpauu+46rVmzRnFxccrKyrqq9wsODlbFihUL/bE/JUff/HpKWf5lJePCzjKOQH8Nj2pgyf+v/GRmZurgwYMuRenSAnX8+PF8X1ujRg19+umn8vNjJxzAU1GqAC+1b98+DRkyRN9//71zLDAwUE2aNJFpmtq1a5eWL1+u2NhY/fDDD/muEuWknVHG0b0qG3ahcFh9iSwkJETdu3dX9+7dJV1YrdmwYYPWrFmjNWvWaMOGDUpLSyv0ffz9/bVx40Y1adKkSJ97Vwlf2jx16pQ++ugj7d6921maDh06dNWXNWfMmKEaNWpYlg+A9ShVgJfJzMzUxIkT9dZbb+UqIyEhIXrqqacUGxurgwcPFvk9U3//xVmqSvoSWXBwsDp37qzOnTtLktLT0xUXF6c1a9Zo7dq1Wrdunc6cOZPrddnZ2XrhhRe0fPly5z1SBYluEV6iN6WHhoaqb9++Gj9+vH788UfnXlpXo02bNgoNDdWZM2dUoUIFC1MCsJJhmqbbPzQyMtKMi4tz++cC3m7Dhg0aPHiwtm3bVuz3CioXrIDrmqnM9S3kqNNaAZWqyxHor3G9mtn6DbmsrCxt3brVuZK1Zs0anThxwnl+yZIl6tatm2358nLo0CG98847+vjjj5Wenl6s94qIiFDjxo1z/QgJCbEoLYDLGYax2TTNyELnUaqA0i8lJUUvv/yyPvjgAxXnv+mWLVsqKipKUVFRateunZbuOO7xG2Tm5ORo586dzoL1559/atmyZQoMDLQ7Wi5HjhzRxIkT9eGHH+r8+fOWvve1116rxo0bK6haTf2WVknp1Rvr+hvqeOQ/M6C0oVQBPsA0TS1cuFBPP/20kpKSrvj1NWrU0J133qmoqCh16dJF1atXL4GU7mWapkzT9Ogbuo8fP65JkyZpypQpOnv2bIl8RuWooarwt64esboIlHZFLVWe+7sOgAIdOnRI0dHR6t279xUXqjfeeEPx8fFKSkrSZ599pv79+3tFoZIufKPQkwuVJFWrVk3jxo3TgQMHNHr0aFWqVCnfuQ6HQ7Vq1brizwisEiHpf5uaAih5nv07D4BcsrOz9d5776lx48b65ptvruo9tm7dqmbNmnl8+fB2lStX1ptvvqkDBw5ozJgxqlKlSq45OTk52rZtm86ePavNmzdr5syZGjVqlO69917Vq1cv33+GgVVrOn/uyZuaAt6Ey39AKZKQkKDBgwdr06ZNxX6vJ554QlOnTi3SN+XgHmfPntWHH36oiRMn6tixY87xefPmqXfv3nm+Ji0tTbt379avv/6q4R8v0cnD+5V9LlnXDHjHOSc8xKEfX+pU4vkBb8U9VYCX2bhxoyZOnKjU1FRlZmYqIyNDmZmZzh9FOb78v/e33npLr732mk3/j5Cf1NRUffzxx5owYYKSkpI0YMAAzZw5s9DXLYpP1KgF25Same0c454qoPgoVQByyc7OdilcWVlZql69OqtVHio9PV2ffPKJpk2bpo0bNxbpG43ufF4j4CsoVQDgJTIzM2WaJg9SBmxS1FLFjuqARVghQEnxxD23AORGqQIssCg+US9Oj1VatqnAkGuUmJyqUQsu7GpOsQIA38D3qYFiSkpK0pAnntD+Dwfr1KrpznH2B4K3OXbsmHbs2KG//vrrqh8MDXgzVqqAfGzZskWHDx9W9+7d8zx/6tQpTZgwQe+9955SUy/sA5S6Z4PSE3eqbHhDSewPBO8SGhqqt956S++//778/f1VrVo1Va9eXTVq1Cjwf6tXr66yZcvaHR8ocZQqIA/x8fG64447NG7cuFznzp8/r/fee08TJkxQcnJyrvOn1sxQjQfGyjAMhYU43BEXcIsyZcpo6tSpioyM1BNPPKGjR4/q6NGjRXptpUqV8ixdYWFhevjhh7kJH16BUgVcZvPmzerSpYtOnTqljIwM53hmZqamT5+ut956K/8/SPwDVKZ6bSknW46gshoe1cBNqQH3eeSRR9SsWTP16tVLhw4dKtJrUlJSlJKSot27dzvHHA6H5syZQ6GC16BUAZfYtGmT7rzzTucKVEZGhnJycvTll19q9OjR2rdvX56v8/Pz0+097ldyw3t1QhX59h+8XmRkpDZv3qw+ffpo9erVV/z6qlWr6ttvv1Xbtm2tDwfYhFIF/NfGjRt155136vTp086xNWvW6LPPPtOWLVvyfV2vXr00ZswYNWrUyB0xAY9RrVo1rVixQiNHjtSkSZOK/Lry5cvrhx9+UOPGjUswHeB+fPsPkLR+/Xp16dLFpVBJ0tdff51voerUqZM2btyo+fPnU6jgswICAvTuu+/q888/l8NRtHsIz549q86dO+uVV17R77//XsIJAfehVMHnrVu3TnfeeafOnDlTpPmtWrXSihUrtHLlSrVp06aE0wGlQ79+/bRhwwbVrl27SPOPHj2qsWPHqk6dOoqKitL8+fOVmZlZwimBkkWpgk9bs2aN7rrrLp09e7bQuQ0aNNC8efO0adMm3XHHHW5IB5QuzZs316ZNmxQVFVXk15imqe+++0733XefIiIiNGrUqHzvXQQ8HaUKPmv16tXq2rWrzp07V+jcFi1a6JdfflHv3r15+DBQgMqVK2vJkiV65ZVX8jxfp06dfF977NgxjR8/XnXr1lWXLl00d+5cl2/gAp6OUgWftHLlSnXr1k3nz58v0vz4+Hjdcsst2rWLHdKBwvj7+2vMmDFasGCBypcv73JuxowZ2rx5s4YMGaIKFSrk+x7ff/+9+vTpo+DKNRTS/gG1H7dSi+ITSzo6UCyUKvicFStW6O6773bugl5UmzdvVsuWLfXxxx/LNM0SSgd4j549e+rnn39Wgwb/26/tzz//VMuWLfXBBx8oKSlJ06dPL/DexKxzyco8cUhJKWkatWAbxQoejVIFnxIbG6sePXooLS3til4XERGh6OhovfTSSwoPD+eSBFBEjRo10s8//6x7771XknT8+HHnufLly+vRRx/Vxo0bFR8fryeffFIVK1bM9R7lb7xwjxbP04SnM+z4G3dkZKQZFxfn9s+Fb1u6dKl69eql9PT0Auddf/31atmypVq1aqVWrVqpZcuWql69uptSAt4pJydHY8eOlZ+fn15++eV85507d05z5szRE69OUHrSLvlXqqHwxz+WYVxYAzAk/T4+7+dxAiXFMIzNpmlGFjaPzT/hE7799lv17t071wrTDTfckKtAVa1a1aaUgPfy8/PTq6++mufzMi8VHBysRx55RNOPXa/fd/+q7LMnnYVKEs/ThEejVMHrff311+rTp48iIiJyFajKlSvbHQ/wKSEhIUWaNzyqgUady1BqZrZzzBHoz/M04dEoVfBqWVlZqlKlio4eParQ0FC74wAooovPzYyJ3aWk5FSep4lSgXuqAAAAClDUe6r49h8AAIAFKFUAAAAWoFQBAABYgFIFAABgAUoVAACABdhSAR5jUXwiX58GAJRarFTBIyyKT9TIr+KUmJwqU1JicioPTwVguZ07dxb6qCrgalGqYLvs7Gw9+cxz+mPGcOVk/u9Bxzw8FYDVsrOzVadOHT3wwAP68ssvdfr0absjwYtQqmCrs2fPqlevXjqybp4yju7RiW8nyTRznOeTklNtTAfA2zRp0kTvv/++vvrqK/Xr109Vq1ZV165dNW3aNB09etTueCjlKFWwzeHDh9WxY0d98803zrHzu39S8prPnMc8PBWA1e655x69/fbbkqTMzEwtX75cjz/+uMLCwtS+fXvFxMRoz549NqdEaUSpQokxTVPZ2dl5ntu8ebPatm2rhIQE1xN+/gqodI0kHp4KoOS88soruv/++13GTNPU+vXrNWLECNWvX19NmzbVq6++qri4ONnxSDeUPjz7DyXm3//+t5o1a6bWrVu7jC9atEj9+/fX+fPnXcbLla+oWn1f0/mqjfj2H4ASd+7cOXXo0EFbtmwpdO51112n6Oho9ezZUx07dlRgYKAbEsJTFPXZf5aVKsMw/CXFSUo0TfPuguZSqrxfUlKSGjdurFdffVXDhg2TdOFvge+++65GjBiR6299tWvX1rfffqtGjRrZEReAj/rjjz/UunVr/fXXX0V+TWhoqO6++2717NlTd955p4KDg0swITyBHQ9UflbSbxa+H0op0zQ1dOhQpaSk6IcffpB04b6FwYMHa/jw4bkKVYcOHbRhwwYKFQC3q1WrlubNm6eAgKJv23jq1CnNnDlTAwcO1MyZM7k0CCdLSpVhGNdJ6i5puhXvh9Jt/vz5WrRokSRp7dq1OnHihO666y5Nn577X4/+/ftr5cqVqlatmrtjAoAk6dZbb9V7771X5PkBAQF65plntHfvXg0ZMkSGYZRgOpQmVu2o/v9JGiGpQn4TDMMYLGmwJNWsWdOij4WnOXnypIYOHeo8TklJUfPmzZWUlJRr7ltvvaVXX32V35AA2G7IkCFKSEjQtGnTCp07e/Zs9enTxw2pUNoUe6XKMIy7Jf1pmubmguaZpjnNNM1I0zQjWZXwXi+88IL+/PNPl7HLC1XZsmX1xRdf6LXXXqNQAfAIhmFo8uTJuvnmmwud279/f7399tvKzMx0QzKUJlZc/usg6R7DMP6Q9KWkToZhzLLgfVHKxMbGasaMGQXOqVatmv7zn//ogQcecFMqACiaMmXKaP78+YqIiChwXlZWlkaPHq327dvr119/dVM6lAbFLlWmaY4yTfM60zRrSXpA0irTNAcUOxlKlTNnzmjw4MEFzjEMQ59++qnatWvnplQAcGWqV6+ur7/+Wg5H4RsPx8XFqWXLloqJicl3Tz74Fjb/hCVeeeUVHTx4sMA5pmmqe/fuuuOOO/T111/zmxAAj9SiRQt98sknLmMjR47UiBEjct2ykJ6erhEjRuiWW25hF3ZYW6pM01xd2B5V8D4//vijpkyZUuT5K1euVHR0tJ577jllZWWVYDIAuDp9+/bVqFGjnMeNGzfWhAkTtG7dOtWtWzfX/J9++kk33nijpkyZopycnFzn4RtYqUKxpKWl6bHHHivyPi3lypXT0KFDtWvXLk2ePPmK9oYBAHcaM2aM7r77wjpB7dq1JUnt27dXQkKCnn766VzzU1NT9fTTT6tLly46cOCAW7PCM1CqUCxjxozRzp07C51Xs2ZNxcTE6PDhw5oyZYrq16/vhnQAcPX8/Pw0a9YsNWzY0FmqJCk4OFjvvfeeVq5cqeuvvz7X61atWqU6DRqravfn1WH8Ki2KT3RnbNiIUoWrlpCQoAkTJhQ4p0OHDpo7d6727dunYcOGKTQ01E3pAKD4KlWqpMWLF+uaa67Jda5Tp07aunWrHnvssVznstPPK/t8shKTUzVqwTaKlY+gVOGqZGVl6dFHH83znqiAgAD1799fP//8s9atW6f77ruPy3wASq26devKzy/vPy4rVqyojz/+WEuXLlVYWJhzvGxYQ1Vs3VOSlJqZrZjYXW7JCnvxJx2uyqRJk/TLL7+4jFWpUkWPP/64nnzySYWHh9uUDADcr2vXrtq+fbsiOvbW+V0/qUr352X4+TvPJyWn2pgO7kKpwhXbvXu3Xn/9dedx48aN9dxzz6l///4qV66cjckAwD6hoaG6ccCrOnDgDwVUquFyLiyk8H2vUPpx+Q9XJCcnR4MGDVJaWpq6deum7777Ttu3b9egQYMoVAB83vCoBqpQNcxlzBHor+FRDWxKBHdipQpXZPbs2WratKmmTZumBg34TQIALhXd4sKtDzGxu5SUnKqwEIeGRzVwjsO7GUXdX8hKkZGRZlxcnNs/F8WXnp6usmXL2h0DAAC3MQxjs2makYXN4/IfrgiFCgCAvFGqAAAALECpAgAAsAClCgAAwAKUKgAAAAtQqgAAACxAqQIAALAApQoAAMAClCoAAAALUKoAAAAswLP/fMCi+ESeQwUApcimTZuUnJysLl262B0FV4BS5eXmbzqgUXPilOEfJElKTE7VqAXbJIliBQAeqkWLFmrcuLEaNmyod999V/Xq1bM7EoqAy39e7tkRL+v36U8p49g+51hqZrZiYnfZmAoAUJCAgAC99tprWrx4sZo0aaJhw4YpJSXF7lgoBKXKiy1btkyJqz9XVvJRHZk5TGcSlss0TUlSUnKqzekAAAXp16+f6tevr8zMTOdq1bRp05SdnW13NOSDUuWlDh48qAEDBvxvIDtTJ797X1knDkuSwkIcNiUDABRFQECARo8e7Tw+fvy4Hn/8cbVq1UqrV6+2LxjyRanyQhkZGerbt69OnjzpMh5yy0MKrBohR6C/hkc1sCkdAKCoHnjgATVs2NBlbMuWLbr99tvVu3dv7d+/36ZkyAulyguNHDlSGzZscBkLbXiTKrXtpfAQh8b1asZN6gBQCvj7+7usVl1qwYIFatSokUaNGqUzZ864ORnyYly8x8adIiMjzbi4OLd/ri+YP3++7rvvPpex66+/Xr/88osqV65sUyoAwNXKzs5Ws2bN9Ntvv+U755prrtHYsWP18MMPy8+P9RKrGYax2TTNyMLm8StfSm3fvl3x8fEuY3v37tXAgQNdxgIDAzVnzhwKFQCUUv7+/nr99dcLnHP06FENHDhQrVu31rp169yUDJejVJVSn376qcaMGeM8Tk1N1f3336/Tp0+7zJs0aZLatGnj7ngAAAvdf//9atKkSaHzfvnlF3Xs2FF9+/bVgQMH3JAMl6JUlUKZmZmaNWuWFixYoG3bLmzk+dxzzykhIcFlXp8+fTR06FA7IgIALOTn51foatWl5syZo4YNG+q1117TuXPnSjAZLkWpKoViY2N17NgxSdKYMWM0a9YsTZs2zWVO/fr1NX36dBmGYUdEAIDFevfurWbNmhVpbs2aNfXoo4+qY8eOCgjg4Snuwq90KfTpp586fz537lwtXrzY5XxQUJDmzp2rChUquDkZAKCkXFytuvzLSJebO3euevfuzV+qbcBKVSlz4sQJffPNN85j0zSVmuq6O/r777+v5s2buzsaAKCE9ezZs9Df31988UUdOnTITYlwKUpVKfPFF18oMzMz3/Ndu3bVQw895MZEAAB38fPz0xtvvFHgnIMHD+qOO+5w3iYC96FUlTKXXvrLy7Jly1ShQgW1bdtWjz/+uObOneueYAAAt4iOjtbf/vY35/Hf/va3XKtXe/bsUZcuXXI9WQMli1JVimzbtk2bN28udF5qaqp+/vln7dq1S+3bt3dDMgCAuxiG4bJa9dhjj+m7775T/fr1XeZt27ZNXbt2Zbd1N6JUlSIzZswo0jx/f3+9/fbbWrlypcLDeRwNAHibe+65Ry1btlRAQID69u2rGjVq6Pvvv1fNmjVd5v3888+65557ct17i5JBqSolLu5NVZiaNWvqhx9+0Kuvvip/f383JAMAuNvF1apu3bqpatWqkqSIiAitXLlS11xzjcvc1atX67777lNGRoYdUX0KpaqUuHRvqvz07t1bCQkJ6tChg5tSAQDscvfdd7s8WUOS6tatqxUrVuR6NNnSpUs1YMAAZWVluTOiz6FUlRIF3aDucDg0bdo0zZ07V6Ghoe4LBQCwjWEYeW4G2rRpU8XGxubaq3Du3LkKbdFV7cd+r0Xxie6K6VMoVaXA5XtTXapp06aKi4vToEGD2OgNACBJioyM1JIlS+RwOFzGz27/Xtvn/59emr+VYlUCKFWlQH57Uz355JP6+eef1bhxYxtSAQA8WceOHbVw4UIZ/q4PT/GvUFVpWTmKid1lUzLvRakqBS6/9BcaGqqFCxdq6tSpuf4WAgDARVFRUaraY4RkXPjjvnKXJ1TppguPuUlK5huBVuPZfx7u8r2pOnbsqNmzZysiIsLGVACA0qJe287KyUyTZKp8087O8bAQ/lJuNUqVh7u4N5Wfn59Gjx7NVgkAgCsyPKqBRp3LUGpmtnPMEeiv4VENbEzlnShVHiwrK0uzZs1SRESEZs+erY4dO9odCQBQykS3uLAJdEzsLiUlpyosxKHhUQ2c47AOpcqDxcbGqn379po+fXquPUcAACiq6BbhlCg3oFR5sPr162v+/PlslQAAQClAqfJg9erVszsCAAAoIrZUAAAAsAClCgAAwAKUKgAAAAtQqgAAACxAqQIAALAApQoAAMAClCoAAAALUKoAAAAsQKkCAACwAKUKAADAApQqAAAAC1CqAABAntLS0uyOUKpQqgAAQJ4+/PBDbdy40e4YpQalqgQtik9Uh/GrdMNLS9Rh/Cotik+0OxIAAEXWvHlzderUScuXL7c7SqlAqSohi+ITNeyTlTp88qxMSYnJqRq1YBvFCgBQarRr105ZWVnq0aOHZs2aZXccj0epKiHvLN+pg1+9qcRpg3V682LlZKYpNTNbMbG77I4GAECROBwOtW/fXllZWfr73/+uSZMm2R3Jo1GqSsj+hPXK/HO/slOO6dT3Hynxg4HKOH5AScmpdkcDAKDIbr/9dufPX3zxRY0cOVKmadqYyHNRqkpI2uaFLsd+ZRwKrHKdwkIcNiUCAODKderUyeX4nXfe0SOPPKLMzEybEnkuSlUJ2Lx5s1L2/eIyVrFNT5UrW0bDoxrYlAoAgCvXpk0blStXzmVsxowZ6tmzp86fP29TKs9EqSoBMTExLsd+joqqd3MPjevVTNEtwm1KBQDAlStTpow6dOiQa3zJkiXq0qWLTp48aUMqz0Spsti+ffs0d+5cl7HXX3pRG0Z3o1ABAEqlyy8BXvTTTz+pY8eOOnz4sJsTeSZKlcUmTZqknJwc53G5cuU0dOhQGxMBAFA8l96sfrlff/1V7du312+//ebGRJ6JUmWhP//8U//+979dxh577DFVqVLFpkQAABRfq1atVKFChXzPHzp0SDfffLM2bNjgxlSeh1J1lXJycjR58mSXsSlTprg8J8nf31/PP/+8u6MBAGCpgIAA3XLLLQXOOXnypDp37qxly5a5KZXnoVRdpcOHD+u5557TihUrJElnz57VlClTXOb07dtXtWrVsiEdAADWyu++qkudP39e99xzj2bOnOmGRJ6HUnWV9uzZo5ycHPXt21f79+/Xv/71L506dcplzogRI2xKBwCAtQq6r+pSWVlZeuihh/Tuu++WcCLPE2B3gNJq7969kqRTp06pZ8+eSk5OdjkfFRWlG2+80Y5oAABY7sYbb1RoaGiuBYSLypYtq5tvvllVq1ZV1apVdebMGSUmJio83He++U6pukp79uxx/nzr1q25zo8cOdKdcQAAKFF+fn667bbbtHDhwjzPp6en680338xzTytfweW/q3RxpSovQUFBWrhwoR5//HH94x//UL9+/bR//343pgMAwHqFXQIcOnSosrKy3JTG81CqrtKlK1WXS0tL0+TJkzVt2jTNmDFDbdq0Ue3atd2YDgAA6116s3qjRo301ltvuZzfsmWLPvjgA3fH8hiGHU+ajoyMNOPi4tz+uVbJyclRuXLllJ6eXujc559/XpMmTXJDKgAASpZpmrrmmmv0559/au3atWrTpo1uvPFG7dy50zmnYsWK2r17t2rUqGFjUmsZhrHZNM3IwuaxUnUVDh8+XKRC1adPH02cONENiQAAKHmGYej222/XkCFDdPPNN6tMmTK5thM6ffq0z377nVJ1FQq69HfRLbfcohkzZsjPj19iAID36N+/v8aNG+c87ty5s/r27esy57PPPtPatWvdHc12/Il/FQq6SV2SGjdurEWLFikoKMhNiQAAcI8ePXooJCTEZezdd99V+fLlXcZ88aZ1StVVKGilKiwsTMuWLVNoaKgbEwEAYJ/w8HC98cYbLmPbtm3T1KlT7QlkE0rVVchvpapChQpaunSK7/E1AAAgAElEQVSpatas6eZEAADY65lnnlHjxo1dxkaPHq0jR47YlMj9KFVXIa+VqoCAAC1YsIBd1AEAPikwMDDXytTp06dVr0t/dRi/SoviE21K5j6UqiuUk5Ojffv25Rr/5JNPdMcdd9iQCAAAz3DbbbfpwQcfdBk7t+M/2rdtk0Yt2Ob1xYpSdYXy2k5h3LhxGjBggE2JAADwHDExMfIvW855XOaaujICg5Sama2Y2F02Jit5PPvvCl1+6e/JJ5/kOX8AAPxXWFiYKnZ4UCk/famQWx5S+RujZPj5S5KSklNtTleyKFVX6NKb1O+991699957MgzDxkQAAHiWBp3u1+Emt8u/XCWX8bAQh02J3KPYl/8Mw4gwDOM/hmH8ZhjGDsMwnrUimKe6uFJ100036fPPP5e/v7/NiQAA8CwjuzVR+UqVXcYcgf4aHtXApkTuYcVKVZakF03T/MUwjAqSNhuGscI0zV8teG+Ps3fvXtWrV0+LFy9WuXLlCn8BAAA+JrpFuCQpJnaXkpJTFRbi0PCoBs5xb1XsUmWa5hFJR/778zOGYfwmKVySV5aqM2fOaPny5apatardUQAA8FjRLcK9vkRdztJv/xmGUUtSC0kb8zg32DCMOMMw4o4fP27lx7rVe++9p9q1a9sdAwAAeBjLSpVhGOUlzZf0nGmapy8/b5rmNNM0I03TjKxWrZpVH+t2TZo0sTsCAADwQJaUKsMwAnWhUM02TXOBFe8JAABQmljx7T9D0r8k/Waa5qTiRwIAACh9rFip6iDp75I6GYaR8N8f3Sx4XwAAgFLDim//rZPE7pcAAMCn8ew/AAAAC1CqAAAALECpAgAAsAClCgAAwAKUKgAAAAtQqgAAACxAqQIAALAApQoAAMAClCoAAAALUKoAAIAtjh49ancES1GqAACALZ555hmlpqbaHcMyPl+qFsUnqsP4VbrhpSXqMH6VFsUn2h0JAACfsG/fPo0cOdLuGJbx6VK1KD5RoxZsU2JyqkxJicmpGrVgG8UKAAA3CAoK0uTJk7V8+XK7o1jCp0tVTOwuJR/4VVkpx2SapiQpNTNbMbG7bE4GAID3czgckqR//OMfOn78uM1pii/A7gB2SkpO1V9L/qmsk4flFxyismENValdHyWpvt3RAADwekFBQZKkY8eO6dFHH9XXX38twzBsTnX1fHqlqnrZLGWdPCxJyjmXrNQ9G2RmZyssxGFzMgAAvN/FlSpJWrx4saZNm2ZjmuLz6VIVVe2M64BfgCpF1NfwqAb2BAIAwIdcXKm66Pnnn9fOnTttSlN8Pl2q/P7a63JcPqyuJvRppegW4TYlAgDAd1y6UiVJqamp6t+/vzIyMmxKVDw+Xao2bNjgcjyw150UKgAA3OTyUiVJv/zyi15//XUb0hSfz5aqnJycXKWqXbt2NqUBAMD3XH7576IJEyZo9erV7g1jAZ8tVbt27VJKSorL2E033WRTGgAAfE9eK1WSZJqmHnroIZ06dcrNiYrHZ0vV5atUNWrU0PXXX29TGgAAfE9+K1WSdOjQIT3xxBPOfSRLA58tVevXr3c5vummm0r13hgAAJQ2+a1UXfTVV19p1qxZbkpTfD5Rqs6fP6/Tp0+7jHE/FQAA9ipopeqioUOH6vfff3dDmuLziVKVnp6udu3aae/eC1sonDlzRtu3b3eZw/1UAAC4V2ErVdKFP7MHDBigrKwsNyQqHp8oVZUqVdLOnTvVunVrLV++XJs2bXK5Ruvv76/IyEgbEwIA4HuKslIlST/99JPGjRtXwmmKzydKlZ+fn6pUqaLk5GR1795dI0eOdDnfvHlzBQcH68yZM1q+fLnGjh1bKhoxAAClWWErVe3atVODBg1UpkwZvfnmm9q4caObkl0dnyhVklSlShVJF/aniouLczmXkpKidu3aKTQ0VF27dpWfn58CAnz6WdMAAJS4wkrVgAEDtHPnTp0/f1779u1TpUqV3JTs6vhcqcrL/v37tWHDBmVnZ6tmzZp69tln3ZgMAADfdOnlvwcffFC9evVyOT9nzhxJF27Tuf7669WwYUO35rtSlKrLjB07tkg3zgEAgOK5+Odtv379NGPGDPXr18/l/Jo1a3TkyBE7ol0VStUlIiMjc/0DBQAAJSMoKEgPPPCAPvvsMwUEBKhbt24qV66c87xpmpo3b56NCa8MpeoS7777rvz8fOaXBAAAW9WuXVszZ8503sdcrlw59ejRw2XOxUuApYHPNIiqVasWeD46Olq33HKLm9IAAACHw5Hri2F9+/Z1OV63bp0SExPdGeuq+UypKmilKiAgQBMmTHBjGgAAkJe77rpL5cuXdxkrLZcAKVWSnnjiCdWvX9+NaQAAQF4cDofuuecel7GvvvrKpjRXxudLVaVKlTR69Gg3pwEAAPnp06ePy/H69et18OBBm9IUnc+XqldeeaXQ+60AAID7REVFqWLFii5jpeESoM+UqryKU61atfT000/bkAYAAOQnKChI9957r8vY6//fdC2K9+wb1n2mVFWuXDnX2Lhx44r8MEcAAOA+Ea06uxyfPfSbXvz39x5drHymVAUGBrosJbZt2zbX1zYBAIBnWHWmhoyywZKkwGq1FNLx78pQgGJid9mcLH8+9dTgKlWq6PTp05IubPRpGIbNiQAAQF6Ons1W1e7PKzA0XIFVI5zjScmpNqYqmM+sVEn/u1m9d+/e6tChg81pAABAfsJCHCpX7yaXQnVx3FP5VKmqWrWqAgMDNX78eLujAACAAgyPaiBHoL/LmCPQX8OjGtiUqHA+d/lv6NChqlu3rt1RAABAAaJbhEuSYmJ3KSk5VWEhDg2PauAc90Q+Varq1avHFgoAAJQS0S3CPbpEXc6nStULL7ygChUq2B0DAAB4IZ+6p4pCBQAASopPlSoAAICSQqkCAACwAKUKAADAApQqAAAAC1CqAAAALECpAgAAsAClCgAAwAKUKgAAAAtQqgAAACxAqQIAALAApQoAAMAClCoAAAALUKoAAAAsQKkCAACwAKUKAACUClu3blVKSordMfJFqQIAAKXCjz/+qBkzZtgdI19eWaoWxSeqw/hVuuGlJeowfpUWxSfaHQkAABRTQkKCpk6dqpycHLuj5MnrStWi+ESNWrBNicmpMiUlJqdq1IJtFCsAAEq5LVu2aPfu3fr+++/tjpInrytVMbG7lJqZ7TKWmpmtmNhdNiUCAADFlZ2dra1bt0qSpkyZYnOavAXYHcBqu9d9q7M7VjuPy4Y3UsjNDyopOdW+UAAAoFj27t2r1NQLf5Z/++232r9/v2rXrm1zKldet1IVlPqX0v6Id/7IOP67JCksxGFzMgAAcLUSEhKcPzdNUx988IGNafLmdaWqfZ0qucYcgf4aHtXAhjQAAMAKW7ZscTn+17/+pfPnz9uUJm9eV6oaXlvR5dgR6K9xvZopukW4TYkAAEBxXV6qTp06pS+++MKmNHnzulJ1uVvrV6NQAQBQyl16+e+iKVOmyDRNG9LkzetLFQAAKN2OHz+upKSkXOMJCQn66aefbEiUN0oVAADwaJdf+ruUJ22vQKkCAAAeraBSNW/ePB05csSNafJHqQIAAB4tr/upLsrKytJHH33kxjT5o1QBAACPVtBKlSR99NFHysjIcFOa/FGqAACAx0pPT9dvv/1W4JyjR49qwYIFbkqUv1JfqtatW6eUlJQizc3IyNCKFStKOBEAALDKr7/+qqysrELnecIN66W+VO3bt0+1atXS22+/nW+5ysjI0EcffaR69epRqgAAKEUKu/R30Y8//qj4+PgSTlOwUl+qunXrppSUFI0ePVq1atXSN99843J+y5YtqlevnoYMGaKDBw+qR48eNiUFAABXqqCb1CXpxhtvdP586tSpJR2nQKW+VFWrVk3t2rWTJCUnJ+f6xd+/f78OHjwoSapcubJzLgAA8HwXV6puuOGGPEvTmjVrlJCQoOeff14rV67UyZMn3R3RqdSXKklFXn3q1q2bAgICSjgNAACwgmma2r59u1588UVt27ZNDz30UK45hw8f1o033qhJkyZp9+7dCgwMtCHpBV5Rqu6++25L5wEAAPudP39eS5cu1cSJExUcHKzy5curUqVKLnMOHz7s/HlgYKAqVKjg7phOXlGqmjRpolq1ahU4JyAgQHfddZd7AgEAgGILDg5W69atXcauu+46l+NLS5XdvKJUGYZR6CXAW265JVe7BQAApQulyg0KK1V86w8AgNKPUuUGt956q8qXL5/veUoVAAClH6XKDcqUKaOoqKg8zzVs2FB16tRxcyIAAGC1y0tVYmKiTUly85pSJeW/GsUqFQAA3sHrV6oMw7jLMIxdhmHsNQzjJSve82p069ZNhmHkGqdUAQDgHS4vVSdPntT58+dtSuOq2KXKMAx/SVMldZXUWFI/wzAaF/d9r8alu6tfxC7qAAB4j/Dw8Fxjwz5dZUOS3KxYqWojaa9pmvtN08yQ9KWkey1436ty+aoUu6gDAOA9Jq4+JAWUcRlbsHabXl20zaZE/2NFqQqXdOiS48P/HbPF5bums4s6AADe48ufD8svsKzLWGbKUX2x8VA+r3AfK0pV7puYJDPXJMMYbBhGnGEYccePH7fgY/N26e7q7KIOAIB3yTZNBTe+7X8D/gEKimiqbDNX9XA7K0rVYUkRlxxfJynp8kmmaU4zTTPSNM3IatWqWfCxebt0d3V2UQcAwLv4G4ZCbv2HykY0VWD1GxTSvp8CQ8Pkn8cX1dzNilK1SVI9wzBuMAyjjKQHJH1jwftetYulim/9AQDgXfq1jZBfYFld8+B4hT0yWZXa93WO263Yd3CbppllGMZTkmIl+Uv6t2maO4qdrBgu7q5OqQIAwLuMiW4mSfpi4yFlm6b8DUP92kY4x+1kydfiTNNcKmmpFe9lhTJlyuj5559nF3UAALzQmOhmHlGiLudVO6pf6qWXbNuDFAAA+CCvLVXlypWzOwIAAPAhXluqAAAA3IlSBQAAYAFKFQAAgAUoVQAAABagVAEAAFiAUgUAAGABShUAAIAFKFUAAAAWoFQBAABYgFIFAABgAUoVAACABShVAAAAFqBUAQAAWIBSBQAAYAFKFQAAKDW2b99ud4R8UaoAAECp8eqrr+rkyZN2x8gTpQoAAJQKf/zxhxYvXqwtW7bYHSVPXlmqFsUnqsP4VbrhpSXqMH6VFsUn2h0JAAAU0wcffKCcnBwlJCTYHSVPAXYHsNqi+ESNWrBNqZnZkqTE5FSNWrBNkhTdItzOaAAA4CqlpqZq+vTpksRKlbvExO7SubNnlHpgi1J++kop6+coNTNbMbG77I4GAACu0pdffum8l8pTS5XXrVTt2bBCx795RzJzJEn+5Sur4k33Kyk51eZkAADgapimqcmTJzuPd+zYoYyMDJUpU8bGVLl53UpV+A31nIVKkrLPnlT2meMKC3HYmAoAAFytDRs2KD4+3nmcmZmpnTt32pgob15Xql4b0EV+Zcu5jJlHd2t4VAObEgEAgOK4dJXqIk+8Wd3rSlWvVhFq1iLSZayV4y9uUgcAoBQ6cuSI5s6dm2vcE++r8rpSJUn3dLnV5fj4fs/dfRUAAOTv448/VlZWVq5xSpWbtGvXzuX4l19+UXp6uk1pAADA1cjMzNSHH36Y57mEhASZpunmRAXzylLVpk0bl+OMjAyXG9wAAIDnW7hwoY4cOZLnuRMnTigpKcnNiQrmlaWqSpUqql+/vsvY+vXrbUoDAACuxpQpUwo872mXAL2yVEm5LwFu2LDBpiQAAOBKbdmyRWvXri1wjqd9A9BrS9VNN93kcnx5qfr++++VkpLizkgAAKCIClulklipcpvLV6oOHjyopKQkmaapf/7zn7rrrrvkcLAhKAAAnubkyZOaPXt2ofM8baXKax5Tc/EbAIZhSJKaNGmi4OBgnTt3zjln9erVio2N1WeffaYKFSp43Pb2AABA+uSTT5SaWvjj5fbs2aNz584pODjYDakK5zWlKjs7Wz169FD16tV122236fbbb1fr1q21evVq55xHHnlEGRkZki7czA4AADxLdna23n//fedxlSpVdOLEiTznmqap7du3q23btu6KVyCvufwXEBCgp556Sp999pkGDhyoG264IdcNbhcLlUSpAgDAEy1btkzJycl66qmntGnTJv3zn/90OV+zZk0NGjTIeexJlwC9plRJUrdu3dS5c2fncXZ2dr5zq1at6o5IAADgCjRs2FBHjhzR5MmTFRkZqcTERJfzderU0bRp07Rq1SrVrVvXo25W96pSZRiGJk6c6LyvqiCsVAEA4Hnq1q3rcs/z4cOHXc5fd911kqTbb79dW7du1a23uj6azk5eVaok6W9/+5sefvjhQudRqgAA8Hz5lSpJcjgc6tu3r7sj5cvrSpUkjRkzptDtEihVAAB4voJKlafxylIVHh6uYcOGFTiHUgUAgOejVHmAESNGqEaNGvme50Z1AAA8W0ZGhv7880+XMUqVDcqXL6+333473/OsVAEA4NmOHDni3Nz7ovDwcJvSFM5rS5UkDRw4UE2bNs3zHKUKAADPdvmlv8DAQFWrVs2mNIXz6lLl7++vmJiYPM9RqgAA8GyXl6rw8HD5+XludfHcZBa56667dOedd+Yap1QBAODZlm/c4XLsCKluU5Ki8fpSJUkxMTEuG4KWLVvWYx6+CAAAclsUn6ivf9ruMnY4w6FF8Yn5vMJ+PlGqmjdvroEDBzqPq1SpUqRd1wEAgD1iYncp41yKy5gRXFkxsbtsSlS4ALsDuMvbb7+tL774QufPn+fSHwAAHi4pOVVVu7+gyl2eVPbZE8o685f8g0OVlJxqd7R8+cRKlSRde+21GjFihCTupwIAwNOFhVx4MopfmSAFVg6X4/obVaZqTee4J/KZUiVJw4YN07XXXkupAgDAww2PaiBHoL/LmCPQX8OjGtiUqHA+VaqCg4M1ZswYShUAAB4uukW4xvVqpvAQhwxJ4SEOjevVTNEtPHfzT5+5p+qihx9+WOXLl7c7BgAAKER0i3CPLlGX86mVKunChqD333+/3TEAAICX8blSJYntFAAAgOV8slQBAABYjVIFAABgAUoVAACABShVAAAAFqBUAQAAWIBSBQAAYAFKFQAAgAUoVQAAABagVAEAAFiAUgUAAGABShUAAIAFKFUAAAAWoFQBAABYgFIFAABsdfToUaWlpdkdo9goVQAAwDamaeqRRx6RYRh2Ryk2ShUAALDNlClTFBsbqzJlytgdpdgC7A5gt0XxiYqJ3aWk5FSFhTg0PKqBoluE2x0LAACvt2PHDg0fPlxBQUFesVLl06VqUXyiRi3YptTMbElSYnKqRi3YJkkUKwAASlB6eroefPBBpaenKzQ01O44lvDpUhUTu8tZqDJPHdHpTQsVUKGqYoLLUKoAAChBL7/8srZu3SpJcjgcNqexhk+XqqTkVGWeTFTympk6v/snycyRX1AFHW51j93RAADwWitWrNCkSZOcx0FBQTamsY5P36geFuKQmZWu87vWSWaOJCkn7Yz89/zH5mQAAHinEydO6B//+IfLmLesVPl0qRoe1UCVwusp6IaWLuNn4xYpMzPTplQAAHgn0zQ1aNAgJSUluYxTqrxAdItwjevVTHU6P+gyfvxooubMmWNTKgAAvNO///1vLVy4MNc4l/+8RHSLcG2b9oIiIyNdxt955x2ZpmlTKgAAvMuePXv07LPP5nmOlSovYhiGRo4c6TK2detWxcbG2pQIAADvkZmZqQEDBujcuXN5nmelysv07NlTderUcRmbMGGCTWkAAPAeb731ln7++ed8z7NS5WX8/f01bNgwl7HVq1fn+pdg+/bt7owFAECptm7dOo0dO7bAOaxUeaGHH35Y1atXdxl75513nD//7bff1L17d+61AgCgCFJSUjRgwADl5OQUOI+VKi/kcDj0zDPPuIwtWLBAe/bsUUpKiqKjo3Xw4EGdPHnSpoQAAJQeTz31lA4cOFDoPEqVl3ryyScVHBzsPDZNU++8844GDBig3bt3S7rwDQYAAJC/L774QrNmzSrSXC7/eanQ0FANHjzYZWz69On69ttvncd79+51dywAAEqNAwcO6IknnlClSpXUqlUr3XfffQXO95aVKp9+9t9Fu3bt0t69exUUFKSyZcuqY8eO+uc//5nvfFaqAADIn8Ph0L59+1S5cmUZhqE1a9Zo3rx5+c73lpUqSpWkGjVqqHfv3tqxY0eR5rNSBQBA/i7/0tflTylp1KiRXnjhBQ0aNEiS96xUcflPUkhIiJYtW6awsLAizWelCgCAosnOzs61StWnTx899thjmjZtmiRWqrxORESEli1bpo4dO+r06dMFzt2zZ49M05RhGG5KBwBA6bRmzRodO3bMZaxPnz6SpEGDBsk0TZUrV86OaJYr1kqVYRgxhmHsNAxjq2EYCw3DCLEqmB2aN2+uhQsXKjAwsMB5ycnJbKsAAEARXH7pr2nTpmrcuLHzePDgwYXeyF5aFPfy3wpJTU3TbC5pt6RRxY9kr06dOunTTz8tdB6XAAEAKFhWVpbmz5/vMnZxlepSZcuWdVekElWsUmWa5nemaWb993CDpOuKH8l+Dz74oMaPH1/gHG5WBwCgYKtXr9bx48ddxu6//36b0pQ8K29UHyhpmYXvZ6sRI0Zo6NCh+Z5npQoAgIJdfumvefPmatiwoU1pSl6hpcowjO8Nw9iex497L5nziqQsSbMLeJ/BhmHEGYYRd3lr9USGYej//u//FB0dned5VqoAAMhfZmamFixY4DLWt29fm9K4R6Hf/jNN846CzhuG8bCkuyV1Ngt40rBpmtMkTZOkyMjIUvFEYn9/f33++efq3Lmz1q9f73KOlSoAAPK3atUqnThxwmXMmy/9ScX/9t9dkkZKusc0zfPWRPIsDodDixcvVv369V3GWakCACB/MR984nJcu2FT1atXz6Y07lHce6qmSKogaYVhGAmGYXxoQSaPU6VKFS1fvtxlh9hTp07lauAAAECa+/Pv+k/sEpexM2GttSg+0aZE7lHcb//VNU0zwjTNv/33xxCrgnmaG264QUuXLlVwcLBzjNUqAABym/TdblW+80k56reT/C/s/RhYr4NiYnfZnKxksaP6FWjVqpXmzp2rHj16KDs7W3v27FHbtm3tjgUAgEc5ejZLwY1uUXCjW5STfl5ph7YrMOQaJSWn2h2tRPHsvyvUtWtX57OKWKkCACC3sJD/PSDZr2w5lavbJte4N2Kl6ioMHDhQhw4d0u7du+2OAgCAxxke1UCjFmxTama2c8wR6K/hUQ1sTFXyWKm6SqNHj9Ztt91mdwwAADxOdItwjevVTOEhDhmSwkMcGtermaJbhNsdrUQZBWwtVWIiIyPNuLg4t3+u1UzTlGEYdscAAAAlyDCMzaZpRhY2j5WqYqBQAQCAiyhVAAAAFqBUAQAAWIBSBQAAYAFKFQAAgAUoVQAAABagVAEAAFiAUgUAAGABShUAAIAFKFUAAAAWoFQBAABYgFIFAABgAUoVAACABShVAACgSObMmSPTNO2O4bEoVQAAoECmaer111/XV199JcMw7I7jsQLsDuDNFsUnKiZ2l5KSUxUW4tDwqAaKbhFudywAAIosOztbQ4cO1UcffaQpU6bYHcejUapKyKL4RI1asE2pmdmSpMTkVI1asE2SKFYAgFIhLS1N/fv314IFCyRJt99+u82JPBulqoTExO5yFipJOvfrap3NzlZMcBlKFQDA46WkpOjee+/VDz/8IEmqUaOGGjVqZHMqz0apKiFJyamSJDMrQydXfqyzCctkBJTRHzVukNTJ3nAAABTgyJEj6tq1q7Zs2eIcu/3227mfqhDcqF5CwkIcyslI09FZw3U2YZmk/xasbyYoJSXF5nQAAORt79696tChg0uhkrj0VxSUqhIyPKqBgoODVebaei7jaScS9dhjj/GVVACAx9m8ebPat2+v33//Pde5Tp24ylIYSlUJiW4RrnG9mqlpr2dVpkYdl3Pz5s3T5MmTbUoGAEBuK1eu1G233abjx4/nOnfdddepTp06ebwKl6JUlaDoFuFa/9pd+u2n71SpUiWXc8OGDdPGjRttSgYAwP/MmTNH3bp109mzZ/M836lTJ+6nKgJKlRvUrl1bn3zyictYZmam+vTpoxMnTriMm6apxMREd8YDAPiwqVOn6oEHHlBGRka+c7ifqmgoVW7Ss2dPvfDCCy5jBw8e1EMPPaScnBzn2HfffaexY8e6Ox4AwMeYpqnRo0frqaeeKvQ+X0pV0VCq3Gj8+PFq166dy9jSpUv1zjvvSLrwL/ibb76pL774QmlpaXZEBAD4gOzsbA0ZMkRvv/12oXNr166t66+/3g2pSj9KlRsFBgbqq6++UpUqVVzGX3nlFf3www9atWqV1q9fr1OnTmnx4sU2pQQAeLO0tDTdf//9mjZtWpHm862/oqNUuVlERIRmz57tcsNfTk6OHnjgAY0aNco59umnn9qQDgDgzbKzs/Xss89q9erVRX4Nl/6KjlJlg6ioKL366qsuY0ePHtWmTZucx8uXL9eRI0fcHQ0A4MX8/f310Ucf6eTJk8rMzNSxY8fUsWPHAl9DqSo6SpVNXn/99QKXVHNycjRr1iw3JgIA+JKAgADFx8dr7dq1+c5p2LChrr32WjemKt0oVW6Uk5OjvXv3at68eXr99deVnp5e4PwZM2aw8zoAoESkp6fr6aefdhkLDQ1V69atncesUl0ZHqjsRitWrFC/fv106tSpIs3fsWOHNm/erMjIyBJOBgDwNe+++6727NnjMjZ27Fh17NhRLVq0UGZmJjepXyFWqtwoKipKW7ZsKfT69aW4YR0AYLUDBw5ozJgxLmOtWrXSoEGD1KRJE7300kuSpNtuu82GdKUXpcrNIiIi9J///EdvvPGG/PwK/+X//Nj86LMAABykSURBVPPPC71MCADAlXj++eeVmprqPDYMQ++//778/f0lSS+//LL69OmjqlWr2hWxVKJU2cDf31+vv/66Vq9erYiIiALnsmcVAMBKy5Yt08KFC13GHnvsMbVp08Z5HBQUxJWSq0CpslHHjh2VkJCgnj17FjiPf7EBAFZIS0vLdXN65cqV83w8msPhcFcsr0GpslnlypU1f/58ffDBBwoKCspzDntWAQCsMHHiRO3bt89lbNy4cVzmswilygMYhqEhQ4Zo06ZNaty4ca7z2dnZmj17tg3JAADe4o8//tD/+3//z2WsdevWevTRR21K5H0oVR6kadOm2rRpkx5//PFc5z799FP2rAIAXJVF8Ylq1b2/0tLSnGOGYWjq1KnOm9NRfJQqD1OuXDl9+OGHmjdvnkJCQpzjF/esAgDgSiyKT9TwGT8o+fftLuNdevV32egTxUep8lC9e/fWli1b1KFDB+cYN6wDAK5UTOwuZQZVUtigj1S+RTdJhvwcFXWqUS+7o3kdSpUHq1mzplavXq3Ro0fLz8+PPasAAFcsKfnCflT+jgqqcueTuuahSarS9VkdzyxjczLvQ6nycAEBAXrzzTe1atUqlStXjj2rAABXJCzEdWuEstfWU7l6bXONo/goVaXErbfeqi1bthT5uYEAAEjS8KgGcgS63ozuCPTX8KgGNiXyXjxQuRSpUqWKBg0aZHcMAEApEt0iXNKFe6uSklMVFuLQ8KgGznFYh1IFAICXi24RTolyAy7/AQAAWIBSBQAAYAFKFQAAgAUoVQAAABagVAEAAFiAUgUAAGABShUAAP9/e3ceFdWVaA18H5lnUBwAUSFxgNYohphWgi8rxtbGIYAfxgk7UaNETRANJuqzh8Sv0We3pmOwlbz0Wh0cYgYBxYFg54VoTPrJ4NhKGwMokCiiCIQSCuq8PxA61ZQyeKtuFbV/a7FC3aFqs26EXbfOPZdIASxVRERERApgqSIiIiJSAEsVERERkQJYqoiIiCzA2bNnkZaWpnYMegCWKiIiIjNWUVGBJUuWYNy4cXjyySfVjkMPwBsqW4H0gjLenZyIyMI0NDRg27ZtePPNN1FdXY34+Hj4+vqqHYsegKWqm0svKMOa/eeg0TahsfoGvi0ux5ofGwCAxYqIyAxJKZGZmYlVq1bh8uXLAABHR0e8/vrrKiej9rBUdXObswqh0TZBV1+HG5+8Ce3Nq2icGIfNLvYsVUREZubChQtISEhAdna23vKlS5fCx8dHpVTUURxT1c2VV2kgdU24eXAztBXFgNTh1mfbcf7Td9DU1KR2PCIiAlBZWYnly5dj5MiRbQqVk5MTVq9erVIy6gyWqm7O19MJtef+Bs2VU3rLq3MzEBUVhdraWpWSERGRVqvFO++8g8GDByM5Odngm91ly5ahb9++KqSjzmKp6uYSJw1F79G/gPuY6DbrDh48iPDwcJSWlhrc99KlS8aOR0RktY4ePYqRI0ciPj4et2/fNriNs7MzEhMTTZyMuoqlqpuLDPHDxv83CsOjlqHXpOUQPWz01p8+fRpjxoxBXl6e3vL6+npMnjwZlZWVpoxLRNTtFRYWYsqUKfjlL3+JixcvPnDb5cuXo0+fPiZKRg+LpcoKRIb44as3nsHNo9uQ/VkWPDw89NZ///33GD9+PNLT01uXnTp1CiUlJUhISDB1XCKibun27dtISEjA8OHDcfjw4Xa3d3Fx4VkqC8NSZWUmTJiAb775BoGBgXrL6+rqEB0djc2bN0NKiS+//BIAkJqaiiNHjqgRlYioW2hsbMSf//xnDB48GG+//TYaGxs7tN8rr7wCb29vI6cjJQkppclfNDQ0VObm5pr8delfbt68iaioKJw4caLNuoULF6KoqAiff/45AMDf3x8XLlyAm5ubqWMSEVk8KSVOnz6NAwcO4ODBg22GWxji6uqK4uJi9OrVywQJqT1CiDwpZWh72/FMlZXy9vbGsWPHMG/evDbr3n///dZCBQDXrl3D2rVrTRmPiKjbEEIgJCQEv/nNb5Cbm4vz58+3+yb11VdfZaGyQCxVVszBwQEffPAB3nzzzXa3TU5ONnhWi4iIOu7HH3/EkiVLUFNTc99t3NzcsGrVKhOmIqWwVFk5IQTWr1+PvXv3wsHB4b7bSSmxaNEi3L1714TpiIi6j7t37yIqKgpfffXVA7eLj49Hz549TZSKlMRSRbh9+zbKysrg5eX1wO0KCwvx1ltvmSgVEVH3odVqMWvWrDazpffu3RvTp09vfezu7o6VK1eaOh4phPf+s2IXLlzAtm3bkJqairq6ug7ts2nTJsTExGDUqFFGTkdE1D3odDq8+OKLyMjI0Fvu4eGBzz77DFqtFgcOHAAArFixot03uGS+eKbKitnb28POzg5CiA7v09TUhAULFnT4kmAiImsmpcTSpUuxe/duveUuLi44cuQIRo0ahdDQUAwbNgweHh6cG9DCsVRZscGDB2Pbtm0oLS3FH/7wBwwcOLBD+xUUFOCPf/yjkdMREVk2KSVWr16NnTt36i13cHBARkYGxo4dC6B5bGtsbCwSEhLg6empRlRSCOepolaNjY3IyMjA22+/3e6Vfg4ODjh79iyGDBlionRERJblrbfewq9//Wu9Zba2tti/fz+mTZumt7y0tBSurq4sVWaK81RRp9na2mLGjBk4fvw4Tp06hblz58LW1vCwu/r6eixatAg6nc7EKYmIzFt6QRkCpi5rU6iEEEhNTW1TqACgf//+LFTdAEsVGRQaGopdu3ahpKQE69atMzgJ3fHjx9uc1iYismbpBWV4+dd/QPGh7W3WpaSkYNasWSqkIlNhqaIH8vX1xYYNG3Dt2jW89957+NnPfqa3fvXq1bh69apK6YiIzMvmrELonLwgbPXn/Rs0ZSkWLVqkUioyFZYq6hAnJycsWrQI586dQ3Z2NqZMmQIAqK2tRVxcHNQYm0dEZG7KqzRwCnwcfWb+DsLeGQDg8dRcYHiEysnIFFiqqFOEEHj22WeRmZmJwsJCLFu2DDk5OdizZ4/a0YiIVOfr6QQAcPQfjr6zfw+PsDnwGDerdTl1byxV1GVDhgzBu+++i9LSUmi1Ws5dRURWL3HSUDjZ2QAAHPo9Cs+n5sDZ3haJk4aqnIxMgTOq00Pz8vLCCy+8oHYMIiLVRYb4AWgeW1VepYGvpxMSJw1tXU7dG0sVERGRgiJD/FiirBQ//iMiIiJSgCKlSgjxmhBCCiG8lXg+IiIiIkvz0KVKCOEPYCIATlZEREREVkuJM1VbAawGwImKiIiIyGo9VKkSQkwHUCalPKNQHiIiIiKL1O7Vf0KIYwD6GVi1DsBaAL/oyAsJIRYDWAwAAwYM6EREIiIiIvMnunp7ESHECAB/A1B3b1F/AOUAxkgpf3jQvqGhoTI3N7dLr0tERERkSkKIPCllaHvbdXmeKinlOQB9fvKCxQBCpZQ3u/qcREREpqbT6ZCbm4vGxkaMGzdO7ThkwTj5JxERWZ2Ghgbk5OQgLS0NGRkZ6NmzJ77++mu1Y5GFU6xUSSkHKfVcZJ3SC8qwOasQ18q+R3/fflg9eRhnJSYixdTW1uLo0aNIT09HZmYm7ty5AwDo2bMnMjIy4OrqqnJCsnQ8U0VmIb2gDGv2n0NdgxY/7PtPXLdzQPzFBUDCXBYrIuqyiooKHDhwAOnp6cjOzkZ9fb3eehsbG3z88ccIDAxUKSF1JyxVZBY2ZxVCo23Cj//IgbaiGABwNfV1vPj3TzDww50ICQm5775ZWVl49tlnYWNjY6K0RGTOioqKkJ6ejrS0NHz11VfQ6XT33Xbr1q145plnTJiOujPe+4/MQnmVBrJJizsndustr/rnKYwePRqzZs3C5cuXDe6bmpqKF198EU1NTaaISkRmRkqJ06dP47e//S1GjhyJwMBArFy5EsePH39goVq4cCGWL19uwqTU3bFUkVnw9XRCY9UPkDrDxWjfvn0ICgpCXFwcysvL9dY1NDQgNTUV8+fPR2NjoyniEpHKpJT48ssvkZCQgMDAQISEhOB3v/sdzp4926H9x40bh+TkZAghjJyUrAlLFZmFxElD4d5vEPxeSoHXhJfQw8m9zTZNTU3YuXMnHn30Ubzxxhu4ffs2ALSOkdizZw/mzZvHYkVkBYQQaGxsxBdffIHi4uJO7evn54dPP/0UDg4OxglHVoulisxCZIgfkqJHoL+3OzxCn8MTr+/GrLhVcHNza7OtRqPBpk2bEBAQgKSkJFRVVbWu27dvH2bPng2tVmvK+ESkgmeeeQa5ubn4y1/+Ah8fnw7t4+joiPT0dPTrZ+hGIUQPp8szqj8MzqhOHVVRUYGkpCQkJyejoaGhw/tFRUXhww8/hL29vRHTEZE50Ol02LRpE9auXdvutrt27cLcuXNNkIq6k47OqM4zVWTWevfujS1btuDy5ctYsGABevTo2P+yaWlpmDlzZqeKGBFZnuLiYkyYMKFDhSoxMZGFioyKpYoswoABA/D+++/j/PnzmDFjRof2ycjIwIwZM9rMS0NElk9Kiffeew8jRozAF1980e72kydPRlJSkvGDkVVjqSKLEhQUhE8++QQnT56Eu3vbwez/LjMzE9HR0bh7964J0hGRKZSVlSEiIgKLFy9GbW1tm/VeXl56j4cMGYK9e/dyLjsyOpYqsjhFRUV47bXXUF1d3aHtDx8+jMjISGg0GiMnIyJjklJi165dGD58OI4ePdpmvZ+fH44cOYLFixe3LnN3d0dGRgY8PT1NGZWsFEsVWZTdu3dj1KhROHnyZKf2y8rKwnPPPYe6ujojJSMiY7px4wZmzJiB2NhYvSt+W8yfPx/nz5/H5MmTW285I4TA3r17MWzYMFPHJSvFUkUW49q1a8jPz0dISIjBqRbak52djWnTprFYEVmI9IIyhG38HH0i18AvcAjS0tLabNOnTx+kp6fjr3/9a+vZqJZSlZSUhIiICJNmJuvGKRXIIul0Onz77bfIy8tDXl4e8vPzkZ+f33rX+Qd5+umnkZmZCRcXFxMkJaKuaLnJeumhbagtOGxwm5iYGGzfvh3e3t56y69cuYL169dj9+7dnDGdFNHRKRVYqqjb0Ol0+O677/SKVl5ensGPCsLDw3H48GG4urqqkJSI2hO28XOUVWlQe+5vqDy8VW9dz549sX37djz//PMG921sbERDQwOcnZ1NEZWsAEsVEZoHthYVFbUpWrdu3UJYWBiOHDnSpY8Sici4At44BInmf8MV+9+C5tv/BQA4PfIEvjtxgDOik0l1tFTZmiIMkVqEEAgMDERgYCBiYmIANP+SLikpQV5eHnJycjB16lSVUxLRv/P1dEJZlQZCCPSa9Ap+uJkIj3GzMOSpqSxUZLZYqsjqCCEwaNAgDBo0SO0oRHQfiZOGYs3+c9Bom2Dj6gXfl3bC2cEeqyfzSj4yXyxVRERkdiJD/AAAm7MKUV6lgV9PVyROGtq6nMgcsVQREZFZigzxY4kii8J5qoiIiIgUwFJFREREpACWKiIiIiIFcEwVkcquXbuGvn37wt7eXu0oRO368ccfcf36ddy4cQPXr1/X+z48PBwzZ85UOyKRaliqiBSSXlCGzVmFKCkphreDDr+ZPxlRo/u3u9+VK1cwefJkpKSkICwszARJif5Fp9Ph1q1brcWovf8aunemEAK///3vW+eCI7JWnFGdSAEt9ynTaJtw+3/+gur/3Q9bt14Y//QEvDQnChMnTkSvXr0M7ltVVQUvLy8AwOLFi7Fp06bWG8MSGcuBAwcQFxeHGzduoKmpqcvP4+Xlhb1792LSpEkKpiMyLx2dUZ1jqogUsDmrEBpt8x8mTVE+AKCxphKfH/wIs2fPRu/evTFmzBisX78eJ06cQGNjY+u+np6eeOSRRwAAKSkpCAoKwkcffQQ13vCQ9Zg+fTo2btz4UDccHjlyJHJzc1moiO5hqSJSQHmVBkBzkdJWFLdZL6XEqVOnsGHDBoSHh6NXr16Ijo7Gzp07UVxcjNGjR7du+8MPP+D555/HtGnTUFJSYqofgazQ/PnzkZmZCRcXl07vO2fOHJw8eRKBgYFGSEZkmViqiBTg6+kEAGi8cx02rj3b3b66uhppaWmIi4tDQEAADh482GabQ4cOITg4GFu2bNE7s0WkpN69e2Ps2LEd3t7GxgZbt27Frl274OzsbMRkRJaHY6qIFPDTMVVSSmhvlqDx6mn41v4Tlwr+jvr6+od6/tGjRyMlJQWPP/64QonJmtXU1GDv3r1ISUlBXl5eh/fr3bs3PvroIzz99NPGC0dkhjimisiEIkP8kBQ9An6eTughBAIGByHlv36L01/n4NatWzhy5Aji4+MxbFjXbgabn5+PMWPGYOXKlaitrVU4PVkDKSVyc3OxePFi+Pj4YMmSJZ0qVE888QTy8vJYqIgegGeqiEzs6tWryMrKQlZWFo4dO4Y7d+50av8BAwYgOTkZU6dONVJC6k6qq6uxZ88epKSkoKCgoEvPsWDBAiQnJ8PR0VHhdESWoaNnqliqiFRUXl6Oxx57DJWVlZ3eNyYmBn/605/g4+NjhGRkqdILyvBfRy+h+OIZNF08hjvnv0D9Xc0D9/Hx8cGCBQvw6aef4tKlS63L7ezs8M4772DJkiUPdZUgkaXjx39EZq6mpgbTpk3rUqECgI8//hhBQUHYsWMHdDqdwunIErWM7Sv8Yj++T12FG7lH7luohBCIiIhAWloaSkpKsGHDBr0LInx8fJCTk4O4uDgWKqIOYqkiUkFDQwOio6ORn5//UM9z584dvPzyywgPD8eFCxcUSkeWqmW+NKdHnwSE4V/vvr6+WL9+PYqKinDo0CFERkbCzs4OAFBRUQEACAsLQ15eXqeuCiQi3qaGSBWZmZnw9fXFwoULYW9vDzs7O9jZ2d33+weta/m+oaEBUkrVzyqcP38eJ06cQExMzH1nkbc0UkqUl5fDx8cHPXqY73vRlvnSbN16wenRMdBc/qZ5heiBqVMi8NJLLyEiIgK2tm1/9dfX1+POnTtYtmwZtmzZwntREnUBx1QRkaK0Wi0CAwNx/fp1REREIDY2FlOmTLH4Qc7x8fHYsWMHBg4ciICAAAQEBCAwMLD1+4CAAPTs2VPVUhu28XOU3StWmiu5qMx6F66P/QKPPjUNuRtnPXDfiooKHDp0CC+88IIJkhJZFg5UJyLVbNy4EWvWrGl97OHhgZiYGMybNw/h4eFmfbbnfu7evYuf//znOHPmzH23cXNzM1i2Wr6MPVmm/nxpOkBKODvYIyl6BCJD/Iz62kTdGUsVET209IIybM4qxNWrJejn6YpVU0fh+bGPtns2prKyEv7+/tBo2g6S9vf3x9y5cxEbG4vg4GBjRTeKixcv4vHHHzf4c3VE37599c5yLV68GAMGDFA0Y8sxK6/SwNfTCYmThrJQET0klioieig/PetR/t9Loa28CgDo0aMH3N3dW7/c3Nz0Hrd87du3T+/yfENCQkIwb948zJ492+ynhqirq8OlS5ewYcMGpKWldfl5hBCYOXMm1q1bhxEjRiiYkIiMhaWKiB7KT8fnlG5/AU01N432Wj169MCECRMQGxuLqKgouLq6Gu212lNTU4OLFy/iH//4h95XcXExHub3pY2NDebMmYO1a9d2eWZ9IlJHR0sVr/4jIoNariQDAF19nVFfS6fTITs7G9nZ2XB2dkZkZCTmzZuHiRMnGrxSTQnpBWX4/x9/haKCE7CvKUc/3U3cLP0O165dU/R1bG1t8atf/Qpr1qzBI488ouhzE5F5YakiIoN8PZ1QVqWBlBKyoWtjiLqirq4Oe/bswZ49e+Dj44OkpCTExsYqOri95aPN20VXUHl0GwDge8WevZm9vT0WLlyI119/HQMHDlT42YnIHLFUEZFBiZOGYs3+c6hraIR/wsfQNdTBQVePl8N8MdrHEdXV1XpfNTU1eo+zsrJQX1/fqde0tbXFE088gfHjxyM8PBxhYWHw9PRU/GdrmSTTztu/0/v6+/ujtrYWt2/fNrje0dERS5YsQWJiIvz8OECcyJqwVBGRQS1XjDVfSSbg38erw1eS5eTk4MCBA+1u5+joiLFjx2L8+PEYP348nnzySbi4uDx09va0fLTZw9kTPRzdoLtbo7deCIGAgAAEBwcjKCgIwcHBCA4OxrBhw1qnTfj3UuXi4oKlS5di1apV6Nu3r9F/BiIyPyxVRHRfkSF+nb4cv6mpCQkJCQbXubm54amnnmotUaGhoarM3N3y0aYQAk6Bj0OnvQs77wHo6/8I9r0Rg6FDh953TqkzZ86gpKSk9bG7uzteeeUVrFixAt7e3qb6EYjIDLFUEZGiPvjgAxQUFAAAvL29ER4e3lqiRo4cCRsbG5UT/uujTY22Cd7TXgMAONnZICl6BELaKZHp6ekAAC8vL6xYsQKvvvqqUT6iJCLLw1JFRIppamrChQsXsGPHDowfPx7Dhg1T/V6Ehuh/tNm5STK//vprbNy4ES+//DLc3d2NHZWILAjnqSIi6iApJTQajdFvN0NE5qWj81RZ3g24iIhUIoRgoSKi+2KpIiIiIlIASxURERGRAliqiIiIiBTAUkVERESkAJYqIiIiIgWwVBEREREpgKWKiIiISAEsVUREREQKYKkiIiIiUgBLFREREZECWKqIiIiIFMBSRURERKQAlioiIiIiBbBUERERESmApYqIiIhIASxVRERERApgqSIiIiJSAEsVERERkQJYqoiIiIgUIKSUpn9RISoAlJj8hS2TN4CbaoegNnhczA+PiXnicTE/PCadN1BK2bu9jVQpVdRxQohcKWWo2jlIH4+L+eExMU88LuaHx8R4+PEfERERkQJYqoiIiIgUwFJl/lLUDkAG8biYHx4T88TjYn54TIyEY6qIiIiIFMAzVUREREQKYKmyIEKI14QQUgjhrXYWayeE2CyEuCSEOCuESBNCeKqdyZoJISYLIQqFEN8KId5QO4+1E0L4CyH+RwhxUQhxQQgRr3Ym+hchhI0QokAIkal2lu6GpcpCCCH8AUwEcFXtLAQAyAYwXEr5GIB/Alijch6rJYSwAZAM4JcAggHMFkIEq5vK6jUCWCWlDALwcwDLeEzMSjyAi2qH6I5YqizHVgCrAXAQnBmQUn4mpWy89/AbAP3VzGPlxgD4Vkr5nZSyAcCHAJ5TOZNVk1J+L6XMv/d9DZr/gPupm4oAQAjRH8AUAP+tdpbuiKXKAgghpgMok1KeUTsLGbQAwBG1Q1gxPwDXfvK4FPwDbjaEEIMAhAD4u7pJ6J630fwGXad2kO7IVu0A1EwIcQxAPwOr1gFYC+AXpk1EDzomUsqMe9usQ/NHHbtNmY30CAPLeEbXDAghXAF8CmCFlLJa7TzWTggxFcANKWWeEOJptfN0RyxVZkJK+ayh5UKIEQACAJwRQgDNHzPlCyHGSCl/MGFEq3O/Y9JCCPErAFMBTJCcm0RNpQD8f/K4P4BylbLQPUIIOzQXqt1Syv1q5yEAQBiA6UKICACOANyFELuklPNUztVtcJ4qCyOEKAYQKqXkzTBVJISYDGALgP+QUlaonceaCSFs0XyxwAQAZQBOAZgjpbygajArJprfAf4VwC0p5Qq181Bb985UvSalnKp2lu6EY6qIuuZdAG4AsoUQp4UQO9QOZK3uXTCwHEAWmgdEf8RCpbowALEAnrn37+P0vbMjRN0az1QRERERKYBnqoiIiIgUwFJFREREpACWKiIiIiIFsFQRERERKYClioiIiEgBLFVERERECmCpIiIiIlIASxURERGRAv4P6w76yyUyC7AAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7ff1b191a4d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,10))\n",
    "plt.scatter(circle[:,0], circle[:,1])\n",
    "plt.quiver(circle[:,0], circle[:,1], dL[:,0], dL[:,1],units='inches')\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A integral que deve ser resolvida é:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "    \\overrightarrow{B} = \\frac{\\mu_0 I}{4\\pi} \\int \\frac{d\\overrightarrow{L} \\times \\overrightarrow{r}}{|\\overrightarrow{r}|^3}\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "O valor teorico do campo na direcao Z eh 1.18484E-07 e o valor calculado numericamente eh 1.17559E-07\n",
      "O erro da aproximacao numerica eh 0.780536%\n"
     ]
    }
   ],
   "source": [
    "B_resultante = np.zeros((1,3))\n",
    "z_vector = np.array([0,0,z_pos])\n",
    "K = (mi0/(4*np.pi))*I # Representamos todos os termos constantes apenas por K para facilidade\n",
    "\n",
    "for c, dL_vec in zip(circle, dL): # Loop para percorrer pontos e seus dLs associados\n",
    "    r_vec = z_vector - c\n",
    "    dem = np.linalg.norm(r_vec)**3.0\n",
    "    B_resultante = B_resultante + np.cross(dL_vec, r_vec)/dem\n",
    "\n",
    "B_resultante = B_resultante*K\n",
    "valor_analitico = CalcularB(I,R,z_pos)\n",
    "print (\"O valor teorico do campo na direcao Z eh %.5E e o valor calculado numericamente eh %.5E\" % (valor_analitico, B_resultante[0,2]))\n",
    "print (\"O erro da aproximacao numerica eh %f%%\" % ((valor_analitico - B_resultante[0,2])*100/valor_analitico))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Utilizando a solução numérica temos acesso ainda a todas as componentes de $\\overrightarrow{B}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-3.46944695e-25 -1.01915004e-24  1.17559229e-07]]\n"
     ]
    }
   ],
   "source": [
    "print(B_resultante)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pode-se notar que as componentes X e Y não são exatamente 0, porém muito próximos, por conta da precisão da simulação. Novamente, quanto mais pontos forem utilizados, mais precisa é a simulação numérica."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercício 1\n",
    "\n",
    "Refazer a simulação substituindo o círculo de raio $R$ por um quadrado de lado $L=1m$ centrado na origem, mantendo o cálculo do campo no mesmo ponto. \n",
    "\n",
    "Dica: basta gerar a lista de pontos que representa o quadrado, o restante do problema é identico ao problema do círculo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercício 2\n",
    "\n",
    "\n",
    "Dada a simulação do Ex. 1, agora calcular o módulo do campo magnético $\\overrightarrow{B}$ ao longo de uma reta que passa pelo ponto $(0,0,0.1)$ e tem direção $î$ e plotar o resultado no intervalo de $-2 m$ à $2 m$.\n",
    "\n",
    "Dica: Você terá que dividir a reta em diversos pontos e calcular o campo em cada um deles. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}