diff --git a/code/+phys/luttinger.m b/code/+phys/luttinger.m index f225f1a..92ce238 100644 --- a/code/+phys/luttinger.m +++ b/code/+phys/luttinger.m @@ -1,10 +1,11 @@ -function [I,Lfig] = luttinger(w,sloc,gloc) +function I = luttinger(w,sloc,gloc) %% Computes the Luttinger sum-rule, as defined in PRB 90 075150 % % Input: -% w : real valued array, \omega domain -% sloc : complex valued array, \Sigma(\omega) function +% w : real valued array, \omega domain +% sloc : complex valued array, \Sigma(\omega) function % gloc : complex valued array, G_{loc}(\omega) function +% % Output: % I = \frac{1}{\pi}\Im\int_{-\infty}^{\infty}dwG_loc(w)d\Sigma(w)/dw % @@ -12,49 +13,43 @@ % % Copyright (c) 2022, Gabriele Bellomia % All rights reserved. - global DEBUG -% The original definition, as given by Logan, Tucker and Galpin, prescribes -% to integrate over the negative semiaxis: -% -% I = \frac{2}{\pi} \Im\int_{-\infty}^{0} dw G_loc(w) d\Sigma(w) / dw -% -% Nevertheless we find most useful to exploit particle-hole symmetry and -% thus integrate over the whole real-axis. Care has to be taken whenever -% the frequency domain contains the origin, for therein resides a nasty -% pole, arising from the derivative of the self-energy. Thus we divide -% the integrand in (w < 0) and (w > 0) parts, and integrate the union. -% Failing to exclude the (w = 0) point, if exist, leads to a diverging -% uncontrolled result. - -wl = w(w<0); % Build w < 0 patch -gl = gloc(w<0); -sl = sloc(w<0); -wl = wl(1:end-1); -gl = gl(1:end-1); -dl = diff(sl); - -wr = w(w>0); % Build w > 0 patch -gr = gloc(w>0); -sr = sloc(w>0); -wr = wr(1:end-1); -gr = gr(1:end-1); -dr = diff(sr); - -w = [wl,wr]; % Join the patches -ds = [dl,dr]; -g = [gl,gr]; - -integrand = imag(g.*ds); -I = 1/pi*sum(integrand); -I = abs(I); - if DEBUG -Lfig = figure("Name",'Luttinger integrand','Visible','off'); -plot(w,integrand); -xlabel('\omega'); -ylabel('Im[G(\omega)\partial\Sigma/\partial\omega]'); -ylim([-0.1,0.2]); - end + % The original definition, as given by Logan, Tucker and Galpin, prescribes + % to integrate over the negative semiaxis: + % + % I = \frac{2}{\pi} \Im\int_{-\infty}^{0} dw G_loc(w) d\Sigma(w) / dw + % + % Nevertheless we find most useful to exploit particle-hole symmetry and + % thus integrate over the whole real-axis. Care has to be taken whenever + % the frequency domain contains the origin, for therein resides a nasty + % pole, arising from the derivative of the self-energy. Thus we divide + % the integrand in (w < 0) and (w > 0) parts, and integrate the union. + % Failing to exclude the (w = 0) point, if exist, leads to a diverging + % uncontrolled result. + + wl = w(w<0); % Build w < 0 patch + gl = gloc(w<0); + sl = sloc(w<0); + wl = wl(1:end-1); + gl = gl(1:end-1); + dl = diff(sl); + + wr = w(w>0); % Build w > 0 patch + gr = gloc(w>0); + sr = sloc(w>0); + wr = wr(1:end-1); + gr = gr(1:end-1); + dr = diff(sr); + + w = [wl,wr]; % Join the patches + ds = [dl,dr]; + g = [gl,gr]; + + integrand = imag(g.*ds); + I = 1/pi*sum(integrand); + I = abs(I); + end + diff --git a/code/+phys/sopt.m b/code/+phys/sopt.m index c4297da..dd0a575 100644 --- a/code/+phys/sopt.m +++ b/code/+phys/sopt.m @@ -1,4 +1,4 @@ -function im_sloc = sopt(A,F,U) +function iS = sopt(A,F,U) %% SOPT stands for Second Order Perturbation Theory: % it computes the imag of the 2nd order diagram for the self-energy % @@ -58,7 +58,7 @@ D = D + flip(D); % flip{v(1:end)}=v(end:1) %% Imaginary part of the Self-Energy according to SOPT - im_sloc = -pi * U^2 * D; + iS = -pi * U^2 * D; end diff --git a/code/dmft_loop.m b/code/dmft_loop.m index 691d1a7..12ce62b 100644 --- a/code/dmft_loop.m +++ b/code/dmft_loop.m @@ -35,7 +35,7 @@ % %% BSD 3-Clause License % -% Copyright (c) 2020, Gabriele Bellomia +% Copyright (c) 2020-2022, Gabriele Bellomia % All rights reserved. if(~exist('pmode','var')) @@ -55,7 +55,7 @@ while LOOP % Weiss field from local Green's function - g0 = 1 ./ (w + eta - 0.25 .* gloc); + g0 = 1 ./ (w + eta - 0.25*D^2 .* gloc); % Spectral-function of Weiss field A0 = -imag(g0) ./ pi; diff --git a/code/main.m b/code/main.m index 9082b72..9e33ec1 100644 --- a/code/main.m +++ b/code/main.m @@ -11,8 +11,8 @@ end %% INPUT: Physical Parameters -D = 1; % Bandwidth -U = 0.1; % On-site Repulsion +D = 4.0; % Bandwidth +U = 4.0; % On-site Repulsion beta = inf; % Inverse Temperature %% INPUT: Boolean Flags @@ -20,7 +20,7 @@ ULINE = 1; % Takes and fixes the given beta value and performs a U-driven line TLINE = 0; % Takes and fixes the given U value and performs a T-driven line UTSCAN = 0; % Ignores both given U and beta values and builds a full phase diagram -SPECTRAL = 0; % Controls plotting of spectral functions +SPECTRAL = 1; % Controls plotting of spectral functions PLOT = 1; % Controls plotting of *all static* figures GIF = 0; % Controls plotting of *animated* figures UARRAY = 0; % Activates SLURM scaling of interaction values @@ -32,11 +32,11 @@ mloop = 1000; % Max number of DMFT iterations err = 1e-5; % Convergence threshold for self-consistency mix = 0.10; % Mixing parameter for DMFT iterations (=1 means full update) -wres = 2^15+1; % Energy resolution in real-frequency axis -wcut = 6.00; % Energy cutoff in real-frequency axis -Umin = 0.00; % Hubbard U minimum value for phase diagrams -Ustep = 0.10; % Hubbard U incremental step for phase diagrams -Umax = 6.00; % Hubbard U maximum value for phase diagrams +wres = 2^13+1; % Energy resolution in real-frequency axis +wcut = 6.00 * D; % Energy cutoff in real-frequency axis +Umin = 0.00 * D; % Hubbard U minimum value for phase diagrams +Ustep = 0.10 * D; % Hubbard U incremental step for phase diagrams +Umax = 6.00 * D; % Hubbard U maximum value for phase diagrams Tmin = 1e-3; % Temperature U minimum value for phase diagrams Tstep = 1e-3; % Temperature incremental step for phase diagrams Tmax = 5e-2; % Temperature U maximum value for phase diagrams @@ -94,8 +94,7 @@ fprintf('< U = %f\n',U); [gloc{i},sloc{i}] = dmft_loop(gloc_0,w,D,U,beta,mloop,mix,err,'quiet'); Z(i) = phys.zetaweight(w,sloc{i}); - [I(i),infig] = phys.luttinger(w,sloc{i},gloc{i}); - plot.push_frame('luttok.gif',i,mloop,dt,infig); + I(i) = phys.luttinger(w,sloc{i},gloc{i}); S(i) = phys.strcorrel(w,sloc{i}); end if(PLOT) diff --git a/legacy/PYTHON/.ipynb_checkpoints/real_ipt-text_v3-checkpoint.ipynb b/legacy/PYTHON/.ipynb_checkpoints/real_ipt-text_v3-checkpoint.ipynb index 1ca35e4..647a778 100644 --- a/legacy/PYTHON/.ipynb_checkpoints/real_ipt-text_v3-checkpoint.ipynb +++ b/legacy/PYTHON/.ipynb_checkpoints/real_ipt-text_v3-checkpoint.ipynb @@ -3,12 +3,10 @@ { "cell_type": "code", "execution_count": 1, - "metadata": { - "collapsed": true - }, + "metadata": {}, "outputs": [], "source": [ - "%matplotlib inline" + "%matplotlib notebook" ] }, { @@ -43,7 +41,7 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ @@ -65,7 +63,7 @@ " return 1. / (np.exp(exponent) + 1)\n", "\n", "\n", - "def semi_circle_hiltrans(zeta, D=1):\n", + "def semi_circle_hiltrans(zeta, D):\n", " \"\"\"Calculate the Hilbert transform with a semicircular DOS \"\"\"\n", " sqr = np.sqrt(zeta**2 - D**2)\n", " sqr = np.sign(sqr.imag) * sqr\n", @@ -84,7 +82,7 @@ " return -np.pi * U**2 * (Appp + Appp[::-1])\n", "\n", "\n", - "def dmft_loop(gloc, w, U, beta, loops):\n", + "def dmft_loop(gloc, w, U, beta, D, loops):\n", " \"\"\"DMFT Loop for the single band Hubbard Model at Half-Filling\n", "\n", "\n", @@ -98,6 +96,8 @@ " On site interaction, Hubbard U\n", " beta : float\n", " Inverse temperature\n", + " D : float\n", + " Half-bandwidth of the noninteracting system\n", " loops : int\n", " Amount of DMFT loops to perform\n", "\n", @@ -130,7 +130,7 @@ " sigma = hsi + 1j * isi\n", "\n", " # Semi-circle Hilbert Transform\n", - " gloc = semi_circle_hiltrans(w - sigma)\n", + " gloc = semi_circle_hiltrans(w - sigma, D)\n", "\n", " return gloc, sigma" ] @@ -144,38 +144,1011 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 5, "metadata": { - "scrolled": true + "scrolled": false }, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAEUCAYAAAAvLpGtAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl4XOV1+PHv0W7tsixZsuR9w7swDgYCZt+cgAlLMJACBUrJ0j78kraBkIBDQ8FtCgklTXGBBggQwIFgEhPH2EAgYRPGtiyv8ip502bt+8z5/XGv7JEYbdZoRtKcz/PMMzP3vjP3zNVozn2Xe19RVYwxxph2EaEOwBhjzOBiicEYY0wHlhiMMcZ0YInBGGNMB5YYjDHGdGCJwRhjTAeWGIwxxnRgicEYY0wHlhiMMcZ0EBXqAE7GqFGjdMKECaEOwxhjhpTPPvusXFUzeio3JBPDhAkTyM/PD3UYxhgzpIjI/t6Us6YkY4wxHVhiMMYY04ElBmOMMR0MyT4Gf1pbWykpKaGpqSnUoQw6cXFx5ObmEh0dHepQjDFDwLBJDCUlJSQlJTFhwgREJNThDBqqSkVFBSUlJUycODHU4RhjhoBh05TU1NREenq6JYVORIT09HSrSRljem3YJAbAkkIXbL8YY/piWCUGY4wZrj7eU8Gjf9pBq8c74NuyxGCMMUNA/v5jPL6+CK/qgG/LEoMxxgwB6iaEiCA0DVtiCLAnn3ySrKws5s2bx+TJk3nuued6fE1UVBTTpk0jLy+PvLw8HnzwQQAaGxs599xz8Xg8Pb5HS0sLixYtoq2trd+fwRgz+HjdikIwegwtMQTY5s2bWbZsGZs2beKll17iu9/9bo+vWb58ObNmzeLzzz9n48aN3H///QA888wzXH311URGRvb4HjExMVx44YW8/PLL/f4MxpjBp70FKRg1hmFzHoOvH79ZyNZDNQF9z5ljknngilk9lisoKODaa68FYOLEicTExACwd+9e7r77bg4ePEhERATPP/8806dPB+B73/seGzZs4Cc/+Qk/+tGPjr/XCy+8wIsvvnj8+apVq3juuedYuXLl8WW//OUv2bZtG48//jhXXXUV9957LzfddFNAPrMxZvBo71sIxiDDgNQYROQyEdkhIkUico+f9bEi8rK7/mMRmeCzbq6IfCgihSJSICJxgYgpVAoKCpg+fTqqyhNPPMFDDz1Ea2srd9xxB48++ij5+fksW7aMRx55pMPrnnrqKZ5//nnef/99wGka2rNnD76XF7/vvvtYtmxZh9dNnjyZrVu3AjB79mw+/fTTAf18xpjQUFVEgjP8vN81BhGJBH4BXAyUAJ+KyCpV3epT7HbgmKpOEZGlwHLgehGJAn4N/I2qbhKRdKC1vzH15sh+IBQXF1NbW8vixYs5ePAgc+fOZdmyZaxcuZLCwkKuueYaANra2jjnnHM6vHbbtm14PB7mzZsHQHl5OampqcfXb9q0Ca/Xy+zZs9m/fz+rV6/mm9/8Jq2trce/KJGRkcTExFBbW0tSUlKQPrUxJhi8Gpz+BQhMU9LpQJGq7gEQkd8ASwDfxLAEWOY+Xgk8Ic6v2SXAZlXdBKCqFQGIJ2Q2b97MokWLWL9+PceOHWP27Nl8+OGHbNq0iYceeojbb7/d7+uam5u59dZb+eUvf0lycjIAI0aM6HC28saNGznttNMAWLt2Lbt27QJg69atx5NJ+3vFxQ3pSpcxxg9Fg9K/AIFpSsoBin2el7jL/JZR1TagGkgHpgEqImtEZIOI/EtXGxGRO0UkX0Tyy8rKAhB24BUUFHDqqacCkJaWxo033sgf/vAHsrOzWbNmDV6v93g59RmLfP/993PGGWdwySWXHF+WlpaGx+M5nhy8Xi91dXV4PB5ee+01amtraWxs5Fe/+hU33ngjABUVFWRkZNjF8owZhrwanI5nCExi8Bdp5zMwuioTBZwN3OTef01ELvS3EVVdoaoLVHVBRkaPM9OFhG9iALjiiitYvXo1t912G16vlxkzZpCXl8fy5cs7tBP+9Kc/5b333js+XPUf//EfAbjkkkv44IMPAFi8eDF79uwhLy+Pu+66i8LCQhYsWMCdd97J/PnzAXjnnXdYvHhxED+xMSZYvBq8tqRANCWVAGN9nucCh7ooU+L2K6QAle7y91S1HEBEVgPzgXUBiCvoXnjhhQ7PFy1axOeffw7QYSRRZ12dp/Cd73yHRx99lIsuuojRo0ezcePG4+uuvPLKL5R/8cUXefjhh08mdGPMYKcQEaTEEIgaw6fAVBGZKCIxwFJgVacyq4Bb3MfXAuvVaUtZA8wVkXg3YZxLx76JsHbqqady/vnn9/oEt6uuuur4EFhjzPDi1eD1MfS7xqCqbSLyHZwf+UjgGVUtFJEHgXxVXQU8DTwvIkU4NYWl7muPicijOMlFgdWq+of+xjSc3Hbbbb0qFxMTw8033zzA0RhjQmWojUpCVVcDqzstu9/ncRNwXRev/TXOkFVjjDFd0CHW+WyMMWaAed0T3ILBEoMxxgwBzpnPVmMwxhjjUobWqCRjjDEDLJijkiwxGGPMEODV4M3fbonBGGOGALXO56HrySef5K677up1+a5mbwObwc0Yc4IG8cznYTlRTyht3ryZOXPm9Lr88uXL+eCDD3jttde+UE082RncbKIeY4afIXXm86D01j1wpCCw75k1By5/pMdiBQUFLF26lOuuu+749Y2Ki4t54YUXWLFiBR999BHnnHMOTz/9NND17G3QcQa3nmZvA2wGN2OGsWCe+WxNSQG2ZcsW5syZQ0FBAZMmTeKDDz7glltu4fbbb2f58uVs2bKF1157jebm5uOv6Tx7G3xxBreeZm8Dm8HNmOFMg9j5PDxrDL04sh8IxcXFJCYmEhcXR1VVFXfffTfgTLpz++23k52dDUB8fPzxuaDhi7O3QccZ3HozexvYDG7GDGeqSkSQDuWtxhBA7f0LhYWFzJ8/nwj3r7hp0yYWLlwIQElJCWPGjDn+g+5v9jboOINbb2dva38/m8HNmOHHq4oEqTHJEkMAFRQUMHfuXAoKCjr8YG/evJm5c+cCTpJofwz+Z2+DjjO49Wb2NrAZ3IwZzuzM5yGqoKDgeP9CXl4eAE1NTTQ2NpKWlgZ0TBLQ9extcGIGt97M3gY2g5sxw1kwp/YU37mHh4oFCxZofn5+h2Xbtm1jxowZIYpoYHz++ec8+uijPP/8870qf/XVV/Pwww/7naxnOO4fY8LJt1/cwLbDNaz/3nkn/R4i8pmqLuipnNUYBjGbwc0Yc1wQawzDc1TSMGIzuBljoP0Et+Bsy2oMxhgzBNjVVU/SUOwvCQbbL8YMfd4g/hsPm8QQFxdHRUWF/Qh2oqpUVFTYuQ3GDHFHa5rszOe+ys3NpaSkhLKyslCHMujExcWRm5sb6jCMMf0gQFx0cI7lA5IYROQy4OdAJPCUqj7SaX0s8BxwGlABXK+q+3zWjwO2AstU9acnE0N0dDQTJ048uQ9gjDGDXItHyUkKTs2/3+lHRCKBXwCXAzOBG0RkZqditwPHVHUK8BiwvNP6x4C3+huLMcYMVy1tHmKjglNjCMRWTgeKVHWPqrYAvwGWdCqzBHjWfbwSuFDcxjIRuQrYAxQGIBZjjBmWWj1KdOTQGZWUAxT7PC9xl/kto6ptQDWQLiIJwPeBHwcgDmOMGbZa2rzEDKEag78U1nloUFdlfgw8pqp1PW5E5E4RyReRfOtgNsaEm1aPl+jIodP5XAKM9XmeCxzqokyJiEQBKUAlsBC4VkT+HUgFvCLSpKpPdN6Iqq4AVoBzraQAxG2MMUNGMGsMgUgMnwJTRWQicBBYCtzYqcwq4BbgQ+BaYL06Jxyc015ARJYBdf6SgjHGhLvGVs/QSQyq2iYi3wHW4AxXfUZVC0XkQSBfVVcBTwPPi0gRTk1haX+3a4wx4aKhpY22IJ76HJDzGFR1NbC607L7fR43Adf18B7LAhGLMcYMN9WNrQDkpI4IyvaGzSUxjDFmuGpocS69nxwXnNkZLTEYY8wg1+gmhhExkUHZniUGY4wZ5NprDPGWGIwxxgCU1zUDMCLaEoMxxhigor4FgJw063w2xhgD1LijktLiY4KyPUsMxhgzyFU3thITGTGkrq5qjDFmAB2saiQ2OiJoM7hZYjDGmEGusq6FyRmJQdueJQZjjBnkyuuayUiKDdr2LDEYY8wgpqocqGxgVKIlBmOMMUBNUxvNbV4mjooP2jYtMRhjzCBWcqwBgMykuKBt0xKDMcYMYvsrnMQwJdM6n40xxgC7S52Zj3ODdNYzWGIwxphBbW9FPekJMaQG6axnsMRgjDGD2tZDNUwOYjMSWGIwxphBy+tVikrrmJmdHNTtWmIwxphBamdpLW1e5ZSspKBu1xKDMcYMUhsPVAEwNzc1qNu1xGCMMYPU5weqiI2KYEb2EKwxiMhlIrJDRIpE5B4/62NF5GV3/cciMsFdfrGIfCYiBe79BYGIxxhjhoNP91Vy6rjUoF1VtV2/E4OIRAK/AC4HZgI3iMjMTsVuB46p6hTgMWC5u7wcuEJV5wC3AM/3Nx5jjBkOyuua2VNez+kT04O+7UDUGE4HilR1j6q2AL8BlnQqswR41n28ErhQRERVP1fVQ+7yQiBORIJ3pShjjBmk3ttRBsC500YFfduBSAw5QLHP8xJ3md8yqtoGVAOd0+A1wOeq2hyAmIwxZkh7Z0cpCTGR5I1NC/q2owLwHv4av7QvZURkFk7z0iVdbkTkTuBOgHHjxvU9SmOMGSI8XmX99lLOm55BZERw+xcgMDWGEmCsz/Nc4FBXZUQkCkgBKt3nucDrwM2qururjajqClVdoKoLMjIyAhC2McYMTn/eVUZDi4dLZ2WFZPuBSAyfAlNFZKKIxABLgVWdyqzC6VwGuBZYr6oqIqnAH4B7VfUvAYjFGGOGvDc+P0hMZAQXzxwdku33OzG4fQbfAdYA24BXVLVQRB4UkSvdYk8D6SJSBHwXaB/S+h1gCvAjEdno3jL7G5MxxgxV9c1trN5yhEtnZxEfE4jW/r4LyFZVdTWwutOy+30eNwHX+XndT4CfBCIGY4wZDl7//CAtbV6uXzC258IDxM58NsaYQUJVef7D/eSmjeDLU4J//kI7SwzGGDNI/KWogh1Ha7n1rAlBP9vZlyUGY4wZJB5fv4uk2ChuWjg+pHFYYjDGmEHgoz0VfLK3ktvPmciImMiQxmKJwRhjQkxV+bfV20iOi+KOcyaFOhxLDMYYE2q/23iQzSXV3H3RNBJjQzNE1ZclBmOMCaHaplYefHMrk0YlcPOZoe1baBf61GSMMWHswTe3cqyhlf/5xmlERQ6OY/XBEYUxxoShd7aX8upnJdxw+lgWTgrdeQudWWIwxpgQKK1t4ruvbCQ3bQQ//Ernuc1CyxKDMcYEmcerfOvXG6hubOWJG+eTMAg6nH0NrmiMMSYMPLBqC/n7j/HgklnkjU0NdThfYDUGY4wJomc+2MuvPzrA1xfkcvOZE0Idjl9WYwhjqkpdcxs1TW3UNbVR19xGU6uHplYPLW1eWjxevKp4vc50ewJERgiREUJ0ZASxUc5tREwk8TFRJMRGkhQbTVJcFBEhmHXKmMHujY0HefD3WzlzUjoPfW1OqMPpkiWGYcrjVQ5VNbK/ooGSYw0crGrkcHUTR2uaKKttpryuhaqGFtq8nWdhDYykuChS46NJi49hZIJzy0iMJSMplszkOEYnxZKVEkdWShyxUaE9/d+YYPjjlsPc/fJGZo1J5n9vWUD0IBma6o8lhiFOVTlY1UjhoRq2H65l59FadpXWsq+8gRaP93i5CIHMpDhGJ8eSmxZP3thURibEkBofTXJcNElx0cTHRhIfHUlcdCQxURFERwqRERFECAiConi8zq3F46WlzUtTq5emVg8NLR7qm9uobW6jprGV6sZWqhpaONbQSkVdCzuP1FJe30JLm/cLn2FUYixjUuPITRvB2LR4ckfGM25kPONHxpObNmLQjO025mS9VXCYb7+4gWmZSfz69oWD4uzm7gzu6MwXNLV62FhcRf6+SjYcqGJTcRUV9S0AiMC4kfFMzUzk/OmZTByVwLj0eMamxZOdEhfyH1hVpaaxjaO1Ts3lSHUTh6ubOFTVyMGqRrYfruXtraUdElpUhDB2ZDwT0uOZlJHI5IxEJmckMDkzkVGJsSH8NMb0ziv5xdzz281Mz0rmxTsWkpYQE+qQemSJYZDzepXCQzX8eVcZ7+8qY8OBquNH3VMyEzn/lEzm5qYwOyeF6aOTBt2wN18iQkp8NCnx0UwbneS3jNerHK1t4kBFA/srG9hfUc/e8nr2ljfw4Z4KmlpPJI20+Gimjk5i2uhEpo1OYtroJGZkJZMSHx2sj2RMl1SVx97exePrdnHW5HRW3Lxg0NcU2g2NKMNMS5uXv+4uZ03hUd7edpSy2mYATslK4uYzxnPGpHQWTEgjNX7wH3n0VUSEkJ0yguyUEV84E9TrVQ5VN7K7rJ6i0jqKSmvZebSONzYeorap7Xi5rOQ4TslOYkZ2MjOyk5mZnczEUQlEWoe4CZL65jb+ZeVm/lBwmKvn5/DI1XOJiRo6TaKWGAYJr1f5eG8lb2w8yOqCw9Q0tZEQE8m50zO44JTRnDstg4yk8G46iYgQctPiyU2L59xpGceXqypHaprYcaSW7Udq2X64hu1HavlLUTmtHqdzPS46ghnZycwek8KsMcnMzklh2uikIfXPaoaGXUdr+dYLG9hdVsf3LzuFu86dFNLZ2E6GqA7MqJSBtGDBAs3Pzw91GAFxpLqJV/KLeSW/mJJjjcTHRHLprCwWz8lm0bRRNmKnH1ravOwqrWXroRq2Hq5h66EaCg/VUNfs1C6iI4VTspwkMTc3hTmWLEw/qCovfVLMg78vJCEmip8tzeOcqRk9vzCIROQzVV3QY7lAJAYRuQz4ORAJPKWqj3RaHws8B5wGVADXq+o+d929wO2AB/hHVV3T0/aGQ2L4dF8l//eXvfxxyxG8Cl+eks61p+Vy2azskM/eNJx5vcr+ygYKD1VTcLCagpJqthyspsZtioqJimBmdjLzclOYm5vKvLEpTBqVaOdlmG4dqW7iB68XsH57KWdNTuex6/MYnRwX6rC+IGiJQUQigZ3AxUAJ8Clwg6pu9SnzLWCuqt4lIkuBr6nq9SIyE3gJOB0YA7wNTFNVT3fbHKqJwetV3t52lF+8u5tNxVWkjIjm+i+N5aaF4xifnhDq8MKWqnKgsoFNJdVsLq5i80EnWTS0OF/DpNgoZuekkDculbyxqczLTSUrZfD905vg83iVFz/ez/I/7qDN6+VfLj2FW8+aMGgPJHqbGALRx3A6UKSqe9wN/wZYAmz1KbMEWOY+Xgk8IU6j2xLgN6raDOwVkSL3/T4MQFyDhqqypvAoP3t7J9uP1DJ25AgeXDKLa0/LJT7GunlCTUQYn57A+PQErpw3BnD+4YtK69hcUsXG4io2lVTx1Pt7jvdZZCXHkTc2lbxxTqKYm5syqEeEmcD7bH8lP35zK5tLqjlzUjoPXz2HCaOGxwFeIL7JOUCxz/MSYGFXZVS1TUSqgXR3+UedXpsTgJgGjY/3VPBvb21nU3EVk0Yl8J/XzWNJ3piQn1NgutDWAo2VRDZUMr21jumptVwXXw8T62htqOFoRQVllVWUV9dRvr+Bhh3N7MDDHvEwckQkoxKiSUuKJz0pnpTEeCIioyDCvUXFQUwixCZCTALEJDn3sYkQlwoJoyB6RKj3gOnBvvJ6/uNPO/jD5sNkJsXys+vzWJI3Zsh1MHcnEInB397o3D7VVZnevNZ5A5E7gTsBxo0b15f4QqLkWAP/tnobqwuOkJUcx/Jr5nD1/NxBfRr8sNZcC9UHoabEvT8INYegvhwaKqChHOoroLm6y7eIBnLdGxHREBGFNz4SD1G0aATNHqHtmBep9AIe6vEQHeElCi+R2ob4/2p32kgCJKRDfDrEj3KSRVI2pORAcq57nwMj0pwzGk3QFFc28N/vFvFqfgnRkRH8wwVTuOvcycOyphiIT1QCjPV5ngsc6qJMiYhEASlAZS9fC4CqrgBWgNPHEIC4B0Sbx8vTH+zlZ2/vQlHuvmgqf79osnUoDzRVqDsKFUVQsRsqdzv3x/ZBdTE0df7BF0gcDQkZzg9x6qknfojj0yF+JMQmOUf4Me1H+O7RflTc8R/lCPcWDSTgNBvuLa93mp+Kq/i8uIpth2to9SgxtDI+SVmQFc280VHMSI9kSiokaCM0HjuRnI4nqjIo3QZ1R8Db1jH86HgnQYycCCMnQ/pkGDnJuU8ZCxH2fQuUHUdqefLPu3lj4yEiBG5cOI7vnD+FzEHYuRwogeh8jsLpfL4QOIjT+Xyjqhb6lPk2MMen8/lqVf26iMwCXuRE5/M6YOpQ7XzecaSW7726kS0Ha7hoRiYPXDGLsSPjQx3W8KIKtYfhaKFzK93q/HhW7oGWuhPlIqLdH81Jzg9l5yPu5DEQGZwzpJtaPWw9XMOmYqe/YnNJNXvL64+vn5SRwLzcVObkpDBvbAozs1M6Hkh4PVBX6tRyqkvc+4NQfQAq9zlJsLXhRPnIGEibCBnTYfQsyJzp3KdNhAirsfaGx6u8s72UZz/cx/u7yomLjmDpl8Zx56JJjEkdus19wR6uuhj4Gc5w1WdU9SEReRDIV9VVIhIHPA+cilNTWOrTWX0fcBvQBtytqm/1tL3BlhhUlac/2Mu/r9lBYmwUDy6ZxVfnjgl1WEOf1wNlO+DQBji8CY5uhaNboKnqRJmkMZA5A0ZNdY+cJzn3KWMhcvBW8asaWo6PgtpUUs3mkipK3TPcIyOEqZmJzMlJYY57uZMZWcld1zpVofaIU1tqrylV7IaybVC5l+Ots9HxkHEKjJ4JWfMgZz6Mng3Rw/fIt6+KKxtY+VkJKz8r4WBVI6OTY/nGwvHcdMZ4Rg6Baxz1JKiJIdgGU2Koamjhu69sYv32Ui6akckj18y1i7udDFWn2efQBjjo3g5vglb3yDom0T3ynQmZs9z7mU6TzzBxtKaJTW6NosAdMtt+gcQIca6NNXtMCjPHJDNzTDKzslN6vi5USz2UbofSQiexlro1rYYKZ31EtFObyJkPY+ZDzmlOTSOMmqLKapv545bDvLnpMJ/sqwTg7CmjuHHhOC6eOXpY9QtaYgiCrYdquPP5fI7WNPGDxTO49awJw2pkwoDytMGRzXDgQ9j/VzjwkdOuDhAZC9lz3R8q9wcrfUrYNYOoKoeqmygoqWbroWq2HKqh8FA1R2uaj5fJSR3BjOwkTslK5hT3fnx6fPc/ZqpOc9TBDXDwMycZH9oIzTXO+phEyP0SjD8Lxp0JuQuG1WgpVWVPeT3rt5WydttR8vdV4lWYnJHA107N4apTc8hNG55NwJYYBtiawiP8v5c3khQXxS+/cRrzx6WFNJ5Br7UJDubD/g/hwF+h+JMTfQKp450fobGnO0kgcyZEDf1q+0Apq20+fomPrYdr2H64hj3l9XjcSZdioiKYNCqB6VnOFWenZCYyJTOR8SPjux4m7fU6TVGHNkDJp87fqXQroE6tYsypMP5MGHeWcx+XErwPHAClNU18vLeSv+4u5/1d5ZQcawRg+ugkLp2dxeWzszglK2nYH9hZYhhAz324jwdWFTInJ4Un/+Y0slOGz9FUwHi9Tn/Anndg9ztOzaCtCRDnh3/8mc7R6PiznI5g0y9NrR6KSuvYfsSZrGnHkVp2Ha3lUHXT8TJREcKEUQlMHJXAJPfeObEvnqzkuC+erdt4DA587CTy/R/Coc/B2woS6TQ5TToPJp8POQsGVSJvavWw7XANm0uq2VRcxWcHjrG/wumcT4qNYuGkdM6dnsF50zLCbnCIJYYBoKo8tnYnj68v4sJTMvmvG0+1M5d9VR88kQj2vHuiaSjjFJh0Pkw6F8ad4YzBN0FR29RKUWkdu8vq2V1WR1FpHXvL6zlQ0XGGv5jICMaOHMHYkfHkpI5gTOoIctOcy59nJceRmRxLnDY7tb69f3b+xoc2gHqdcy8mfNn9G5/nDAYY4CPv9kmf9lc683XsLqtnd2kd24/UsLe8nvYZa0clxnLquFS+NCGNhRPTmTUmOaxPLrXEMAAefmsbT763hyvmjeHRr88bVp1SJ8XrcZoddv4RdvzRGQUDzvkBk8478UORnB26GI1f7XOC76uoZ19FA8WVDRyoaKCkqoGSY41UNbR+4TWp8dFkJMYyKjGWkYkxjI1rYVbLJqbWf0ZO5cck1e8DoCVhDM2TLsYz+RJk0iJiRyQQExmBCH6bato8Xlo8XppbvTS0emhwp4itbmylprGVyvoWKupaKKtt5mhtE4ernFn/aptPnNsRITB2ZLw7WVMSM8c4I7rGpMQN++ahvrDEEGCPrd3Jz9ft4op5Y/jZ9XnhO+lLYxXsXgc718CutdBY6VzuYfxZMPUSmHyB01Rk/4xDWkNLG4eqGjlU5UzBWlrbxNGaZkprmyiva6GyvoXyuuYOEyTlUMbZkVu4IOJzzo4oIEGaadQY/uKdxXrvfN7x5lEmo45/NVTBo0pvfoIiI4SRCTGMTo4lK3kEY1LjGJsWz9iR8W6TWDxx0eEzkupkBfMiesPe8x/u4+frdnHJzNHhmRSO7Ydtbzo1gwMfOmfhxqfDtEud2+QLhlxnpOlefEwUUzKTmJLpfwrWdm0eL1XukX1NUxu1TV+hvrmNPzU0klz6MVlH3uVL5e9zUdPTAByNn8autHPYMfJ8KuKnEBkZQVREBDFREcRFRxAfE0l8TBSJsVEkj4giZUQMIxNiSBkRHX7/dyFkNYYevLO9lL/91ad8aUIaL9xxRvhM4lJeBNvegK2r4PBGZ1nmLCcRTL/c6XwMo7Huph9UoXzniSbH4o+cvomRk2DGFTBjiTMs2WqZA86akgJgT1kdX/2vD8hMiuWN75xNyohhPMm8qnNpia1vwLZV7lBFnBEnM5c4/8AjJ4Y2RjM81JXB9t87tdC97zk10OQcOOWrMPNKZ7SaHXQMCEsM/dTU6mHx4+9zpLqJ3//D2UzKSBzQ7YWEqnMW7JbfOsmgoggQp79gxpUw46uQkhvqKM31KRUsAAAWhUlEQVRw1njMqUVsexOK3gZPs3Mxw5lLYM61MPaMsDuxcSBZH0M//fjNQvaU1fPfN80ffkmhci9sWQkFv3VGEkkkTDwHzviWc9SWNDrUEZpwMSIN8m5wbs11sOtPTq1144uQ/7RTk5h9Ncy+FrLnWXNTkFhi8OOdHaW89EkxS780lsVzhslQy7pSKHwdCl51hpiCU2X/yn/CzKucy00bE0qxiW4SuNpJEjvecg5gPvol/PW/nMuizL4WZl8DGdNCHe2wZk1JndQ3t3Huf7xLTKSw7nvnDe15FJqqYdvvnWSw9z2nw2/0HKeKPvtqSB38Ex4ZQ0Ol09RZsBL2fQAoZM2BOdc5NztzvtesKekk/ceaHZTXNfPCHQuHZlLwtELROtj4gnOugacZ0ibA2d91EkLmjFBHaEzfxI+E0251brVH3JrvSlh7P6x9wDmJct4NTp9YzPCYcznULDH4KCqt41d/3cdX5mbz5SlDrGnlaKHTLrv5FagvdTrwTrsV5n7dGVpqbbNmOEjKgjO+6dwqdjvf900vwet3wu8TnE7reUthwjnWad0Plhh8PLx6G1ERwg+/MkSOqusrnDbYjS84cxdERMG0yyDvJph6cdBmKDMmJNInw/n3wnn3OCdebnoJCn8Hm150Zuub+3WnJmH9EX1micG1sbiKddtLuXPRpMF9tVRPqzOsb+MLzjA/b6szWuOy5U5TkXUim3Aj7hDr8WfB5f8OO1bDpt/AX34OHzzq1Jjn3eB0Wg+jiZ0GkiUG1+PrdhEbFcG3z5sS6lD8O7LFaSoqeMWZJD4hAxb+vfOFz5od6uiMGRyiRzgJYPY1UHvUqVFveglW/xP88V7nzP28G2HKxYPqUuGDjSUGoKi0lvXbS7nj7Ik9T5UYTPUVzoiijS84s51FRMN0t6loykXWVGRMd5JGw5nfdm5HtjgJYvMrzlnX8enOiKZ5N9j5EX5YYgCe/mAfAH+3aFJoAwGnqWjX2hOjirytkJ3nVJFnXwsJ6aGO0JihJ2s2ZD0EF/3YuTrwxhch/xn4+H+ca4Dl3QBzvm4nd7rCPjHUN7fx2oYSLp+dxejkuNAFcqTgxKiihnJIyHSaivJudCZrN8b0X2TUiasCN1RC4Wuw8SX40w+doa9TLnT+56ZdDtEh/D0IsbBPDKs2HaK5zcs3zhgf/I3Xl/s0FRW4TUWXu01FF1pTkTEDKX4kfOkO51a20xnNtOllePVW5zLys6+BeTdC7oKwa2rqV2IQkZHAy8AEYB/wdVU95qfcLcAP3ac/UdVnRSQeeBWYDHiAN1X1nv7EczJe//wgmUmxnDU5SE00nlbnejAbX3QuQ+xtc5qKFv/URk0YEyoZ0+CiZXDBj5yrBGx80alJ5D8D6VOdpqa5SyElJ9SRBkV/awz3AOtU9RERucd9/n3fAm7yeABYACjwmYisApqBn6rqOyISA6wTkctV9a1+xtRrpbVNfLK3kr87Z+LAT//nr6nojG86RySjZw7sto0xvRMR6Uw8NfkCaKqBrb9zEsS6B2Hdvzrzlufd5FxsMiY+1NEOmP4mhiXAee7jZ4F36ZQYgEuBtapaCSAia4HLVPUl4B0AVW0RkQ1AUK/x/MctRwBYkjdARwGdm4oiY040FU2+0GnvNMYMTnHJMP9m51a5xzk3YtNL8NrfQUwSzFriHNiNP2vYNTX195dptKoeBlDVwyKS6adMDlDs87zEXXaciKQCVwA/72pDInIncCfAuHGBufjb29tKyUyKZdaY5IC8HwBtLSeainatcZqKxpxqTUXGDGUjJ8H5P4Bz74H9f3ESxJbX4fNfO9cim3eDcymOtAmhjjQgekwMIvI2kOVn1X293Ia/VHr8kq4iEgW8BDyuqnu6ehNVXQGsAOfqqr3cdpdaPV7+WlTO1fNz+t+MpAol+c7JZ1t+Cw0V1lRkzHAUEeHMXTLxHFj8H84EQxtfgHcfgXcfhvFnO/0RM5dAbPfzZQ9mPSYGVb2oq3UiclREst3aQjZQ6qdYCSeam8BpLnrX5/kKYJeq/qxXEQfIJ3srafMq507zV8nppbKdTjIoeBWO7YOoOPdaRTdaU5Exw11MglNLmLcUqoph82+cloI3vg2r/9mZDjfvRpiwaMhd0K+/v1yrgFuAR9z7N/yUWQP8m4ikuc8vAe4FEJGfACnAHf2Mo88+3lsJwJl9HY1Uc8ipFRS86ly4TiJg4iJY9C/OFyEugM1SxpihIXUsLPpnOOefoPgTZ+jrltdh88vOBf3mXe80N42aGupIe6VfE/WISDrwCjAOOABcp6qVIrIAuEtV73DL3Qb8wH3ZQ6r6fyKSi9P3sB1nhBLAE6r6VE/bDcREPTes+IjS2ibWfe+8ngu3T15e+BrsfR9Qp99gznVOv0GSv5Y2Y0xYa22E7X9w+iN2r3cnyprtzJg466qQJIneTtQTljO4eb3KnGVruHxONj+9bp7/QrVHnPbDrW84nU3qdTqg2meNGiKZ3xgzCNQcdg4sC38HJZ84yzJnOQli5lVBuzS4zeDWjb0V9dS3eJiXm9JxxbH9ziV7t74BBz4CFEZNc6qHM5c4l6YYZsPSjDFBkJx94oJ+1QedqUoLfwfvPOTcMmc6vzHTLoWseSHvkwjLxLDlYDUAc8YkOHPI7lzjDDEt2+4UyJwF593r/KEyTwlhpMaYYScl58QsdDWHYOsq50S69pFNiaOdibamXgqTzw/J6Kbwakra+2eo2M3HGzbgLc7njLh9SGuDc42i8We5F9e6zJkZyhhjgqmuDIrWOgequ9dDc43z25Q125lsKHMGJOc4Ix5Pci4Ja0ryZ/U/Q9l2TiOSoqgJyPybYfyXncnEbTSRMSaUEjOc4a15NzrXVDvwkXOJ8JJ856zrljqn3H1HBjyU8EoM1/4fxCZx2VM7mZCRwlOX95g4jTEm+CKjT5xIB+D1QF0p1B52ZqkbYEPrrIv+Gj2TtqQc9lY2k5s2iOd1NsYYXxGRTgd2zvzgbC4oWxlEKupb8HiVKZmJoQ7FGGMGpbBLDEWlTjtddkr4zs5kjDHdCbvEcLSmCYBJGVZjMMYYf8IuMRyudhLDyISTG+5ljDHDXdglhoNVjYyIjiQ5LrwGZBljTG+FXWKobmhlTGrcwE/laYwxQ1TYJYay2mZSRkSHOgxjjBm0wi4xFB9rIC3e+heMMaYrYZcYmlo9ZKfaUFVjjOlK2CWG+mYPSXHWlGSMMV0Jq8TQ3OahxeMlMdZGJBljTFfCKjGU1jgziMbHRIY4EmOMGbzCKjFUN7YCkGxNScYY06WwSgzNbR4ARiXFhjgSY4wZvMIqMTS1egGIiwqrj22MMX3Sr19IERkpImtFZJd7n9ZFuVvcMrtE5BY/61eJyJb+xNIbTa1OjSE22voYjDGmK/09dL4HWKeqU4F17vMORGQk8ACwEDgdeMA3gYjI1UBdP+PoleY2t8YQbTUGY4zpSn9/IZcAz7qPnwWu8lPmUmCtqlaq6jFgLXAZgIgkAt8FftLPOHqlor4FgLgoqzEYY0xX+psYRqvqYQD3PtNPmRyg2Od5ibsM4F+B/wQa+hlHr9S0j0qyayUZY0yXejzTS0TeBrL8rLqvl9vwdxlTFZE8YIqq/j8RmdCLOO4E7gQYN25cLzfdkderACTZJbeNMaZLPf5CqupFXa0TkaMikq2qh0UkGyj1U6wEOM/neS7wLnAmcJqI7HPjyBSRd1X1PPxQ1RXACoAFCxZoT3H70+YmhqgIu+S2McZ0pb9NSauA9lFGtwBv+CmzBrhERNLcTudLgDWq+ktVHaOqE4CzgZ1dJYVA8XiVyAixuRiMMaYb/U0MjwAXi8gu4GL3OSKyQESeAlDVSpy+hE/d24PusqBrcxODMcaYrvWrsV1VK4AL/SzPB+7wef4M8Ew377MPmN2fWHqjzeO1ZiRjjOlBWA3ob/OqJQZjjOlBWCUGj1eJigyrj2yMMX0WVr+S1sdgjDE9C6/EYH0MxhjTo7BKDB6rMRhjTI/CKjG0eZVo62MwxphuhdWvpNUYjDGmZ2GVGNq8XiLtrGdjjOlWWCUGjxcirMZgjDHdCqvEAIrlBWOM6V5YJQavgrUkGWNM98IqMagqEZYZjDGmW2GVGLzqf9YgY4wxJ4RVYlCwuRiMMaYH4ZUYVK2PwRhjehBmiQHrYzDGmB6EVWLwqlofgzHG9CCsEoPVGIwxpmdhlRi8asOSjDGmJ2GVGBTszGdjjOlBeCUGVcSqDMYY061+JQYRGSkia0Vkl3uf1kW5W9wyu0TkFp/lMSKyQkR2ish2EbmmP/H0RBUiwioVGmNM3/X3Z/IeYJ2qTgXWuc87EJGRwAPAQuB04AGfBHIfUKqq04CZwHv9jKdbXqsxGGNMj/qbGJYAz7qPnwWu8lPmUmCtqlaq6jFgLXCZu+424GEAVfWqank/4+mWc+bzQG7BGGOGvv4mhtGqehjAvc/0UyYHKPZ5XgLkiEiq+/xfRWSDiLwqIqP7GU+3nKurWmYwxpju9JgYRORtEdni57akl9vw90usQBSQC/xFVecDHwI/7SaOO0UkX0Tyy8rKernpzlu1+RiMMaYnUT0VUNWLulonIkdFJFtVD4tINlDqp1gJcJ7P81zgXaACaABed5e/CtzeTRwrgBUACxYs0J7i9seurmqMMT3rb1PSKqB9lNEtwBt+yqwBLhGRNLfT+RJgjaoq8CYnksaFwNZ+xtMtxeZjMMaYnvQ3MTwCXCwiu4CL3eeIyAIReQpAVSuBfwU+dW8PussAvg8sE5HNwN8A3+tnPN3yeq3z2RhjetJjU1J3VLUC50i/8/J84A6f588Az/gptx9Y1J8Y+sLmYzDGmJ6F1elealdXNcaYHoVZYrCrqxpjTE/CKjF4bQY3Y4zpUb/6GIaaRdMyyE6JC3UYxhgzqIVVYvjRV2eGOgRjjBn0wqopyRhjTM8sMRhjjOnAEoMxxpgOLDEYY4zpwBKDMcaYDiwxGGOM6cASgzHGmA4sMRhjjOlAnGkRhhYRKQP2n+TLRwEDOrf0SbK4+sbi6huLq2+Ga1zjVTWjp0JDMjH0h4jkq+qCUMfRmcXVNxZX31hcfRPucVlTkjHGmA4sMRhjjOkgHBPDilAH0AWLq28srr6xuPomrOMKuz4GY4wx3QvHGoMxxphuDMvEICLXiUihiHhFpMsefBG5TER2iEiRiNzjs3yiiHwsIrtE5GURiQlQXCNFZK37vmtFJM1PmfNFZKPPrUlErnLX/UpE9vqsywtWXG45j8+2V/ksD+X+yhORD92/92YRud5nXUD3V1ffF5/1se7nL3L3xwSfdfe6y3eIyKX9ieMk4vquiGx19886ERnvs87v3zRIcd0qImU+27/DZ90t7t99l4jcEuS4HvOJaaeIVPmsG5D9JSLPiEipiGzpYr2IyONuzJtFZL7PusDvK1UddjdgBjAdeBdY0EWZSGA3MAmIATYBM911rwBL3cf/A3wzQHH9O3CP+/geYHkP5UcClUC8+/xXwLUDsL96FRdQ18XykO0vYBow1X08BjgMpAZ6f3X3ffEp8y3gf9zHS4GX3ccz3fKxwET3fSKDGNf5Pt+hb7bH1d3fNEhx3Qo84ee1I4E97n2a+zgtWHF1Kv8PwDNB2F+LgPnAli7WLwbeAgQ4A/h4IPfVsKwxqOo2Vd3RQ7HTgSJV3aOqLcBvgCUiIsAFwEq33LPAVQEKbYn7fr1932uBt1S1IUDb70pf4zou1PtLVXeq6i738SGgFOjxBJ6T4Pf70k28K4EL3f2zBPiNqjar6l6gyH2/oMSlqu/4fIc+AnIDtO1+xdWNS4G1qlqpqseAtcBlIYrrBuClAG27S6r6Z5yDwK4sAZ5Tx0dAqohkM0D7algmhl7KAYp9npe4y9KBKlVt67Q8EEar6mEA9z6zh/JL+eKX8iG3KvmYiMQGOa44EckXkY/am7cYRPtLRE7HOQrc7bM4UPurq++L3zLu/qjG2T+9ee1AxuXrdpwjz3b+/qbBjOsa9++zUkTG9vG1AxkXbpPbRGC9z+KB2l896SruAdlXQ3bOZxF5G8jys+o+VX2jN2/hZ5l2s7zfcfX2Pdz3yQbmAGt8Ft8LHMH58VsBfB94MIhxjVPVQyIyCVgvIgVAjZ9yodpfzwO3qKrXXXzS+8vfJvws6/w5B+Q71YNev7eIfANYAJzrs/gLf1NV3e3v9QMQ15vAS6raLCJ34dS2LujlawcyrnZLgZWq6vFZNlD7qydB/W4N2cSgqhf18y1KgLE+z3OBQzjXIUkVkSj3qK99eb/jEpGjIpKtqofdH7LSbt7q68Drqtrq896H3YfNIvJ/wD8FMy63qQZV3SMi7wKnAr8lxPtLRJKBPwA/dKvZ7e990vvLj66+L/7KlIhIFJCC0zzQm9cOZFyIyEU4yfZcVW1uX97F3zQQP3Q9xqWqFT5P/xdY7vPa8zq99t0AxNSruHwsBb7tu2AA91dPuop7QPZVODclfQpMFWdETQzOl2CVOj067+C07wPcAvSmBtIbq9z36837fqFt0/1xbG/XvwrwO4JhIOISkbT2phgRGQV8Gdga6v3l/u1ex2l/fbXTukDuL7/fl27ivRZY7+6fVcBScUYtTQSmAp/0I5Y+xSUipwJPAleqaqnPcr9/0yDGle3z9Epgm/t4DXCJG18acAkda84DGpcb23ScztwPfZYN5P7qySrgZnd00hlAtXvgMzD7aiB62EN9A76Gk0mbgaPAGnf5GGC1T7nFwE6cjH+fz/JJOP+4RcCrQGyA4koH1gG73PuR7vIFwFM+5SYAB4GITq9fDxTg/MD9GkgMVlzAWe62N7n3tw+G/QV8A2gFNvrc8gZif/n7vuA0TV3pPo5zP3+Ruz8m+bz2Pvd1O4DLA/x97ymut93/g/b9s6qnv2mQ4noYKHS3/w5wis9rb3P3YxHwt8GMy32+DHik0+sGbH/hHAQedr/LJTh9QXcBd7nrBfiFG3MBPqMtB2Jf2ZnPxhhjOgjnpiRjjDF+WGIwxhjTgSUGY4wxHVhiMMYY04ElBmOMMR1YYjDGGNOBJQZjjDEdWGIwxhjTgSUGY4wxHVhiMOYkiEimiJSLSIuIzOyh7M9FREXkv4IVnzH9YYnBmJPzM5xrOT2qqj1dSO1T9/7cbksZM0jYtZKM6SMRWYgzE9oxYKKqVvdQ/izgL0ALzhSbnu7KGxNqVmMwpu/ud+//u6ek4Cpz72NwahnGDGqWGIzpA3e6x8vdp//bad04ETlbRDr/+Pv+n0UOZHzGBIIlBmP6ZjHOtfELVXV/p3UPA+/jzBPsa4x730b3s/YZMyhYYjCmb05z7z/1s+4c935Hp+WnuvcbrX/BDAWWGIzpm/HufbHvQhGZjzMnb6Oq1nZ6zWL3/q0Bjs2YgLDEYEzfRLv3nf93vuned+hDEJFTgPMBL/DcwIZmTGBYYjCmb9r7FS4UkQgAETkf+FvgdSBGRE53lycCz+D8nz2jqkUhiNeYPrPzGIzpAxG5AFjnPv0Mp0lpMfAu8AMgHzgKvAecjdPx/AlwgarWBzteY06G1RiM6QNVXQ/cDGwFZuN0Rj8OXKOqnwH/jHMi2xKgFuech/MsKZihxGoMxhhjOrAagzHGmA4sMRhjjOnAEoMxxpgOLDEYY4zpwBKDMcaYDiwxGGOM6cASgzHGmA4sMRhjjOnAEoMxxpgOLDEYY4zp4P8Dq4DJClDQdO4AAAAASUVORK5CYII=\n", + "application/javascript": [ + "/* Put everything inside the global mpl namespace */\n", + "/* global mpl */\n", + "window.mpl = {};\n", + "\n", + "mpl.get_websocket_type = function () {\n", + " if (typeof WebSocket !== 'undefined') {\n", + " return WebSocket;\n", + " } else if (typeof MozWebSocket !== 'undefined') {\n", + " return MozWebSocket;\n", + " } else {\n", + " alert(\n", + " 'Your browser does not have WebSocket support. ' +\n", + " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", + " 'Firefox 4 and 5 are also supported but you ' +\n", + " 'have to enable WebSockets in about:config.'\n", + " );\n", + " }\n", + "};\n", + "\n", + "mpl.figure = function (figure_id, websocket, ondownload, parent_element) {\n", + " this.id = figure_id;\n", + "\n", + " this.ws = websocket;\n", + "\n", + " this.supports_binary = this.ws.binaryType !== undefined;\n", + "\n", + " if (!this.supports_binary) {\n", + " var warnings = document.getElementById('mpl-warnings');\n", + " if (warnings) {\n", + " warnings.style.display = 'block';\n", + " warnings.textContent =\n", + " 'This browser does not support binary websocket messages. ' +\n", + " 'Performance may be slow.';\n", + " }\n", + " }\n", + "\n", + " this.imageObj = new Image();\n", + "\n", + " this.context = undefined;\n", + " this.message = undefined;\n", + " this.canvas = undefined;\n", + " this.rubberband_canvas = undefined;\n", + " this.rubberband_context = undefined;\n", + " this.format_dropdown = undefined;\n", + "\n", + " this.image_mode = 'full';\n", + "\n", + " this.root = document.createElement('div');\n", + " this.root.setAttribute('style', 'display: inline-block');\n", + " this._root_extra_style(this.root);\n", + "\n", + " parent_element.appendChild(this.root);\n", + "\n", + " this._init_header(this);\n", + " this._init_canvas(this);\n", + " this._init_toolbar(this);\n", + "\n", + " var fig = this;\n", + "\n", + " this.waiting = false;\n", + "\n", + " this.ws.onopen = function () {\n", + " fig.send_message('supports_binary', { value: fig.supports_binary });\n", + " fig.send_message('send_image_mode', {});\n", + " if (fig.ratio !== 1) {\n", + " fig.send_message('set_dpi_ratio', { dpi_ratio: fig.ratio });\n", + " }\n", + " fig.send_message('refresh', {});\n", + " };\n", + "\n", + " this.imageObj.onload = function () {\n", + " if (fig.image_mode === 'full') {\n", + " // Full images could contain transparency (where diff images\n", + " // almost always do), so we need to clear the canvas so that\n", + " // there is no ghosting.\n", + " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", + " }\n", + " fig.context.drawImage(fig.imageObj, 0, 0);\n", + " };\n", + "\n", + " this.imageObj.onunload = function () {\n", + " fig.ws.close();\n", + " };\n", + "\n", + " this.ws.onmessage = this._make_on_message_function(this);\n", + "\n", + " this.ondownload = ondownload;\n", + "};\n", + "\n", + "mpl.figure.prototype._init_header = function () {\n", + " var titlebar = document.createElement('div');\n", + " titlebar.classList =\n", + " 'ui-dialog-titlebar ui-widget-header ui-corner-all ui-helper-clearfix';\n", + " var titletext = document.createElement('div');\n", + " titletext.classList = 'ui-dialog-title';\n", + " titletext.setAttribute(\n", + " 'style',\n", + " 'width: 100%; text-align: center; padding: 3px;'\n", + " );\n", + " titlebar.appendChild(titletext);\n", + " this.root.appendChild(titlebar);\n", + " this.header = titletext;\n", + "};\n", + "\n", + "mpl.figure.prototype._canvas_extra_style = function (_canvas_div) {};\n", + "\n", + "mpl.figure.prototype._root_extra_style = function (_canvas_div) {};\n", + "\n", + "mpl.figure.prototype._init_canvas = function () {\n", + " var fig = this;\n", + "\n", + " var canvas_div = (this.canvas_div = document.createElement('div'));\n", + " canvas_div.setAttribute(\n", + " 'style',\n", + " 'border: 1px solid #ddd;' +\n", + " 'box-sizing: content-box;' +\n", + " 'clear: both;' +\n", + " 'min-height: 1px;' +\n", + " 'min-width: 1px;' +\n", + " 'outline: 0;' +\n", + " 'overflow: hidden;' +\n", + " 'position: relative;' +\n", + " 'resize: both;'\n", + " );\n", + "\n", + " function on_keyboard_event_closure(name) {\n", + " return function (event) {\n", + " return fig.key_event(event, name);\n", + " };\n", + " }\n", + "\n", + " canvas_div.addEventListener(\n", + " 'keydown',\n", + " on_keyboard_event_closure('key_press')\n", + " );\n", + " canvas_div.addEventListener(\n", + " 'keyup',\n", + " on_keyboard_event_closure('key_release')\n", + " );\n", + "\n", + " this._canvas_extra_style(canvas_div);\n", + " this.root.appendChild(canvas_div);\n", + "\n", + " var canvas = (this.canvas = document.createElement('canvas'));\n", + " canvas.classList.add('mpl-canvas');\n", + " canvas.setAttribute('style', 'box-sizing: content-box;');\n", + "\n", + " this.context = canvas.getContext('2d');\n", + "\n", + " var backingStore =\n", + " this.context.backingStorePixelRatio ||\n", + " this.context.webkitBackingStorePixelRatio ||\n", + " this.context.mozBackingStorePixelRatio ||\n", + " this.context.msBackingStorePixelRatio ||\n", + " this.context.oBackingStorePixelRatio ||\n", + " this.context.backingStorePixelRatio ||\n", + " 1;\n", + "\n", + " this.ratio = (window.devicePixelRatio || 1) / backingStore;\n", + " if (this.ratio !== 1) {\n", + " fig.send_message('set_dpi_ratio', { dpi_ratio: this.ratio });\n", + " }\n", + "\n", + " var rubberband_canvas = (this.rubberband_canvas = document.createElement(\n", + " 'canvas'\n", + " ));\n", + " rubberband_canvas.setAttribute(\n", + " 'style',\n", + " 'box-sizing: content-box; position: absolute; left: 0; top: 0; z-index: 1;'\n", + " );\n", + "\n", + " // Apply a ponyfill if ResizeObserver is not implemented by browser.\n", + " if (this.ResizeObserver === undefined) {\n", + " if (window.ResizeObserver !== undefined) {\n", + " this.ResizeObserver = window.ResizeObserver;\n", + " } else {\n", + " var obs = _JSXTOOLS_RESIZE_OBSERVER({});\n", + " this.ResizeObserver = obs.ResizeObserver;\n", + " }\n", + " }\n", + "\n", + " this.resizeObserverInstance = new this.ResizeObserver(function (entries) {\n", + " var nentries = entries.length;\n", + " for (var i = 0; i < nentries; i++) {\n", + " var entry = entries[i];\n", + " var width, height;\n", + " if (entry.contentBoxSize) {\n", + " if (entry.contentBoxSize instanceof Array) {\n", + " // Chrome 84 implements new version of spec.\n", + " width = entry.contentBoxSize[0].inlineSize;\n", + " height = entry.contentBoxSize[0].blockSize;\n", + " } else {\n", + " // Firefox implements old version of spec.\n", + " width = entry.contentBoxSize.inlineSize;\n", + " height = entry.contentBoxSize.blockSize;\n", + " }\n", + " } else {\n", + " // Chrome <84 implements even older version of spec.\n", + " width = entry.contentRect.width;\n", + " height = entry.contentRect.height;\n", + " }\n", + "\n", + " // Keep the size of the canvas and rubber band canvas in sync with\n", + " // the canvas container.\n", + " if (entry.devicePixelContentBoxSize) {\n", + " // Chrome 84 implements new version of spec.\n", + " canvas.setAttribute(\n", + " 'width',\n", + " entry.devicePixelContentBoxSize[0].inlineSize\n", + " );\n", + " canvas.setAttribute(\n", + " 'height',\n", + " entry.devicePixelContentBoxSize[0].blockSize\n", + " );\n", + " } else {\n", + " canvas.setAttribute('width', width * fig.ratio);\n", + " canvas.setAttribute('height', height * fig.ratio);\n", + " }\n", + " canvas.setAttribute(\n", + " 'style',\n", + " 'width: ' + width + 'px; height: ' + height + 'px;'\n", + " );\n", + "\n", + " rubberband_canvas.setAttribute('width', width);\n", + " rubberband_canvas.setAttribute('height', height);\n", + "\n", + " // And update the size in Python. We ignore the initial 0/0 size\n", + " // that occurs as the element is placed into the DOM, which should\n", + " // otherwise not happen due to the minimum size styling.\n", + " if (width != 0 && height != 0) {\n", + " fig.request_resize(width, height);\n", + " }\n", + " }\n", + " });\n", + " this.resizeObserverInstance.observe(canvas_div);\n", + "\n", + " function on_mouse_event_closure(name) {\n", + " return function (event) {\n", + " return fig.mouse_event(event, name);\n", + " };\n", + " }\n", + "\n", + " rubberband_canvas.addEventListener(\n", + " 'mousedown',\n", + " on_mouse_event_closure('button_press')\n", + " );\n", + " rubberband_canvas.addEventListener(\n", + " 'mouseup',\n", + " on_mouse_event_closure('button_release')\n", + " );\n", + " // Throttle sequential mouse events to 1 every 20ms.\n", + " rubberband_canvas.addEventListener(\n", + " 'mousemove',\n", + " on_mouse_event_closure('motion_notify')\n", + " );\n", + "\n", + " rubberband_canvas.addEventListener(\n", + " 'mouseenter',\n", + " on_mouse_event_closure('figure_enter')\n", + " );\n", + " rubberband_canvas.addEventListener(\n", + " 'mouseleave',\n", + " on_mouse_event_closure('figure_leave')\n", + " );\n", + "\n", + " canvas_div.addEventListener('wheel', function (event) {\n", + " if (event.deltaY < 0) {\n", + " event.step = 1;\n", + " } else {\n", + " event.step = -1;\n", + " }\n", + " on_mouse_event_closure('scroll')(event);\n", + " });\n", + "\n", + " canvas_div.appendChild(canvas);\n", + " canvas_div.appendChild(rubberband_canvas);\n", + "\n", + " this.rubberband_context = rubberband_canvas.getContext('2d');\n", + " this.rubberband_context.strokeStyle = '#000000';\n", + "\n", + " this._resize_canvas = function (width, height, forward) {\n", + " if (forward) {\n", + " canvas_div.style.width = width + 'px';\n", + " canvas_div.style.height = height + 'px';\n", + " }\n", + " };\n", + "\n", + " // Disable right mouse context menu.\n", + " this.rubberband_canvas.addEventListener('contextmenu', function (_e) {\n", + " event.preventDefault();\n", + " return false;\n", + " });\n", + "\n", + " function set_focus() {\n", + " canvas.focus();\n", + " canvas_div.focus();\n", + " }\n", + "\n", + " window.setTimeout(set_focus, 100);\n", + "};\n", + "\n", + "mpl.figure.prototype._init_toolbar = function () {\n", + " var fig = this;\n", + "\n", + " var toolbar = document.createElement('div');\n", + " toolbar.classList = 'mpl-toolbar';\n", + " this.root.appendChild(toolbar);\n", + "\n", + " function on_click_closure(name) {\n", + " return function (_event) {\n", + " return fig.toolbar_button_onclick(name);\n", + " };\n", + " }\n", + "\n", + " function on_mouseover_closure(tooltip) {\n", + " return function (event) {\n", + " if (!event.currentTarget.disabled) {\n", + " return fig.toolbar_button_onmouseover(tooltip);\n", + " }\n", + " };\n", + " }\n", + "\n", + " fig.buttons = {};\n", + " var buttonGroup = document.createElement('div');\n", + " buttonGroup.classList = 'mpl-button-group';\n", + " for (var toolbar_ind in mpl.toolbar_items) {\n", + " var name = mpl.toolbar_items[toolbar_ind][0];\n", + " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", + " var image = mpl.toolbar_items[toolbar_ind][2];\n", + " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", + "\n", + " if (!name) {\n", + " /* Instead of a spacer, we start a new button group. */\n", + " if (buttonGroup.hasChildNodes()) {\n", + " toolbar.appendChild(buttonGroup);\n", + " }\n", + " buttonGroup = document.createElement('div');\n", + " buttonGroup.classList = 'mpl-button-group';\n", + " continue;\n", + " }\n", + "\n", + " var button = (fig.buttons[name] = document.createElement('button'));\n", + " button.classList = 'mpl-widget';\n", + " button.setAttribute('role', 'button');\n", + " button.setAttribute('aria-disabled', 'false');\n", + " button.addEventListener('click', on_click_closure(method_name));\n", + " button.addEventListener('mouseover', on_mouseover_closure(tooltip));\n", + "\n", + " var icon_img = document.createElement('img');\n", + " icon_img.src = '_images/' + image + '.png';\n", + " icon_img.srcset = '_images/' + image + '_large.png 2x';\n", + " icon_img.alt = tooltip;\n", + " button.appendChild(icon_img);\n", + "\n", + " buttonGroup.appendChild(button);\n", + " }\n", + "\n", + " if (buttonGroup.hasChildNodes()) {\n", + " toolbar.appendChild(buttonGroup);\n", + " }\n", + "\n", + " var fmt_picker = document.createElement('select');\n", + " fmt_picker.classList = 'mpl-widget';\n", + " toolbar.appendChild(fmt_picker);\n", + " this.format_dropdown = fmt_picker;\n", + "\n", + " for (var ind in mpl.extensions) {\n", + " var fmt = mpl.extensions[ind];\n", + " var option = document.createElement('option');\n", + " option.selected = fmt === mpl.default_extension;\n", + " option.innerHTML = fmt;\n", + " fmt_picker.appendChild(option);\n", + " }\n", + "\n", + " var status_bar = document.createElement('span');\n", + " status_bar.classList = 'mpl-message';\n", + " toolbar.appendChild(status_bar);\n", + " this.message = status_bar;\n", + "};\n", + "\n", + "mpl.figure.prototype.request_resize = function (x_pixels, y_pixels) {\n", + " // Request matplotlib to resize the figure. Matplotlib will then trigger a resize in the client,\n", + " // which will in turn request a refresh of the image.\n", + " this.send_message('resize', { width: x_pixels, height: y_pixels });\n", + "};\n", + "\n", + "mpl.figure.prototype.send_message = function (type, properties) {\n", + " properties['type'] = type;\n", + " properties['figure_id'] = this.id;\n", + " this.ws.send(JSON.stringify(properties));\n", + "};\n", + "\n", + "mpl.figure.prototype.send_draw_message = function () {\n", + " if (!this.waiting) {\n", + " this.waiting = true;\n", + " this.ws.send(JSON.stringify({ type: 'draw', figure_id: this.id }));\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_save = function (fig, _msg) {\n", + " var format_dropdown = fig.format_dropdown;\n", + " var format = format_dropdown.options[format_dropdown.selectedIndex].value;\n", + " fig.ondownload(fig, format);\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_resize = function (fig, msg) {\n", + " var size = msg['size'];\n", + " if (size[0] !== fig.canvas.width || size[1] !== fig.canvas.height) {\n", + " fig._resize_canvas(size[0], size[1], msg['forward']);\n", + " fig.send_message('refresh', {});\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_rubberband = function (fig, msg) {\n", + " var x0 = msg['x0'] / fig.ratio;\n", + " var y0 = (fig.canvas.height - msg['y0']) / fig.ratio;\n", + " var x1 = msg['x1'] / fig.ratio;\n", + " var y1 = (fig.canvas.height - msg['y1']) / fig.ratio;\n", + " x0 = Math.floor(x0) + 0.5;\n", + " y0 = Math.floor(y0) + 0.5;\n", + " x1 = Math.floor(x1) + 0.5;\n", + " y1 = Math.floor(y1) + 0.5;\n", + " var min_x = Math.min(x0, x1);\n", + " var min_y = Math.min(y0, y1);\n", + " var width = Math.abs(x1 - x0);\n", + " var height = Math.abs(y1 - y0);\n", + "\n", + " fig.rubberband_context.clearRect(\n", + " 0,\n", + " 0,\n", + " fig.canvas.width / fig.ratio,\n", + " fig.canvas.height / fig.ratio\n", + " );\n", + "\n", + " fig.rubberband_context.strokeRect(min_x, min_y, width, height);\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_figure_label = function (fig, msg) {\n", + " // Updates the figure title.\n", + " fig.header.textContent = msg['label'];\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_cursor = function (fig, msg) {\n", + " var cursor = msg['cursor'];\n", + " switch (cursor) {\n", + " case 0:\n", + " cursor = 'pointer';\n", + " break;\n", + " case 1:\n", + " cursor = 'default';\n", + " break;\n", + " case 2:\n", + " cursor = 'crosshair';\n", + " break;\n", + " case 3:\n", + " cursor = 'move';\n", + " break;\n", + " }\n", + " fig.rubberband_canvas.style.cursor = cursor;\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_message = function (fig, msg) {\n", + " fig.message.textContent = msg['message'];\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_draw = function (fig, _msg) {\n", + " // Request the server to send over a new figure.\n", + " fig.send_draw_message();\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_image_mode = function (fig, msg) {\n", + " fig.image_mode = msg['mode'];\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_history_buttons = function (fig, msg) {\n", + " for (var key in msg) {\n", + " if (!(key in fig.buttons)) {\n", + " continue;\n", + " }\n", + " fig.buttons[key].disabled = !msg[key];\n", + " fig.buttons[key].setAttribute('aria-disabled', !msg[key]);\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_navigate_mode = function (fig, msg) {\n", + " if (msg['mode'] === 'PAN') {\n", + " fig.buttons['Pan'].classList.add('active');\n", + " fig.buttons['Zoom'].classList.remove('active');\n", + " } else if (msg['mode'] === 'ZOOM') {\n", + " fig.buttons['Pan'].classList.remove('active');\n", + " fig.buttons['Zoom'].classList.add('active');\n", + " } else {\n", + " fig.buttons['Pan'].classList.remove('active');\n", + " fig.buttons['Zoom'].classList.remove('active');\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype.updated_canvas_event = function () {\n", + " // Called whenever the canvas gets updated.\n", + " this.send_message('ack', {});\n", + "};\n", + "\n", + "// A function to construct a web socket function for onmessage handling.\n", + "// Called in the figure constructor.\n", + "mpl.figure.prototype._make_on_message_function = function (fig) {\n", + " return function socket_on_message(evt) {\n", + " if (evt.data instanceof Blob) {\n", + " /* FIXME: We get \"Resource interpreted as Image but\n", + " * transferred with MIME type text/plain:\" errors on\n", + " * Chrome. But how to set the MIME type? It doesn't seem\n", + " * to be part of the websocket stream */\n", + " evt.data.type = 'image/png';\n", + "\n", + " /* Free the memory for the previous frames */\n", + " if (fig.imageObj.src) {\n", + " (window.URL || window.webkitURL).revokeObjectURL(\n", + " fig.imageObj.src\n", + " );\n", + " }\n", + "\n", + " fig.imageObj.src = (window.URL || window.webkitURL).createObjectURL(\n", + " evt.data\n", + " );\n", + " fig.updated_canvas_event();\n", + " fig.waiting = false;\n", + " return;\n", + " } else if (\n", + " typeof evt.data === 'string' &&\n", + " evt.data.slice(0, 21) === 'data:image/png;base64'\n", + " ) {\n", + " fig.imageObj.src = evt.data;\n", + " fig.updated_canvas_event();\n", + " fig.waiting = false;\n", + " return;\n", + " }\n", + "\n", + " var msg = JSON.parse(evt.data);\n", + " var msg_type = msg['type'];\n", + "\n", + " // Call the \"handle_{type}\" callback, which takes\n", + " // the figure and JSON message as its only arguments.\n", + " try {\n", + " var callback = fig['handle_' + msg_type];\n", + " } catch (e) {\n", + " console.log(\n", + " \"No handler for the '\" + msg_type + \"' message type: \",\n", + " msg\n", + " );\n", + " return;\n", + " }\n", + "\n", + " if (callback) {\n", + " try {\n", + " // console.log(\"Handling '\" + msg_type + \"' message: \", msg);\n", + " callback(fig, msg);\n", + " } catch (e) {\n", + " console.log(\n", + " \"Exception inside the 'handler_\" + msg_type + \"' callback:\",\n", + " e,\n", + " e.stack,\n", + " msg\n", + " );\n", + " }\n", + " }\n", + " };\n", + "};\n", + "\n", + "// from http://stackoverflow.com/questions/1114465/getting-mouse-location-in-canvas\n", + "mpl.findpos = function (e) {\n", + " //this section is from http://www.quirksmode.org/js/events_properties.html\n", + " var targ;\n", + " if (!e) {\n", + " e = window.event;\n", + " }\n", + " if (e.target) {\n", + " targ = e.target;\n", + " } else if (e.srcElement) {\n", + " targ = e.srcElement;\n", + " }\n", + " if (targ.nodeType === 3) {\n", + " // defeat Safari bug\n", + " targ = targ.parentNode;\n", + " }\n", + "\n", + " // pageX,Y are the mouse positions relative to the document\n", + " var boundingRect = targ.getBoundingClientRect();\n", + " var x = e.pageX - (boundingRect.left + document.body.scrollLeft);\n", + " var y = e.pageY - (boundingRect.top + document.body.scrollTop);\n", + "\n", + " return { x: x, y: y };\n", + "};\n", + "\n", + "/*\n", + " * return a copy of an object with only non-object keys\n", + " * we need this to avoid circular references\n", + " * http://stackoverflow.com/a/24161582/3208463\n", + " */\n", + "function simpleKeys(original) {\n", + " return Object.keys(original).reduce(function (obj, key) {\n", + " if (typeof original[key] !== 'object') {\n", + " obj[key] = original[key];\n", + " }\n", + " return obj;\n", + " }, {});\n", + "}\n", + "\n", + "mpl.figure.prototype.mouse_event = function (event, name) {\n", + " var canvas_pos = mpl.findpos(event);\n", + "\n", + " if (name === 'button_press') {\n", + " this.canvas.focus();\n", + " this.canvas_div.focus();\n", + " }\n", + "\n", + " var x = canvas_pos.x * this.ratio;\n", + " var y = canvas_pos.y * this.ratio;\n", + "\n", + " this.send_message(name, {\n", + " x: x,\n", + " y: y,\n", + " button: event.button,\n", + " step: event.step,\n", + " guiEvent: simpleKeys(event),\n", + " });\n", + "\n", + " /* This prevents the web browser from automatically changing to\n", + " * the text insertion cursor when the button is pressed. We want\n", + " * to control all of the cursor setting manually through the\n", + " * 'cursor' event from matplotlib */\n", + " event.preventDefault();\n", + " return false;\n", + "};\n", + "\n", + "mpl.figure.prototype._key_event_extra = function (_event, _name) {\n", + " // Handle any extra behaviour associated with a key event\n", + "};\n", + "\n", + "mpl.figure.prototype.key_event = function (event, name) {\n", + " // Prevent repeat events\n", + " if (name === 'key_press') {\n", + " if (event.which === this._key) {\n", + " return;\n", + " } else {\n", + " this._key = event.which;\n", + " }\n", + " }\n", + " if (name === 'key_release') {\n", + " this._key = null;\n", + " }\n", + "\n", + " var value = '';\n", + " if (event.ctrlKey && event.which !== 17) {\n", + " value += 'ctrl+';\n", + " }\n", + " if (event.altKey && event.which !== 18) {\n", + " value += 'alt+';\n", + " }\n", + " if (event.shiftKey && event.which !== 16) {\n", + " value += 'shift+';\n", + " }\n", + "\n", + " value += 'k';\n", + " value += event.which.toString();\n", + "\n", + " this._key_event_extra(event, name);\n", + "\n", + " this.send_message(name, { key: value, guiEvent: simpleKeys(event) });\n", + " return false;\n", + "};\n", + "\n", + "mpl.figure.prototype.toolbar_button_onclick = function (name) {\n", + " if (name === 'download') {\n", + " this.handle_save(this, null);\n", + " } else {\n", + " this.send_message('toolbar_button', { name: name });\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype.toolbar_button_onmouseover = function (tooltip) {\n", + " this.message.textContent = tooltip;\n", + "};\n", + "\n", + "///////////////// REMAINING CONTENT GENERATED BY embed_js.py /////////////////\n", + "// prettier-ignore\n", + "var _JSXTOOLS_RESIZE_OBSERVER=function(A){var t,i=new WeakMap,n=new WeakMap,a=new WeakMap,r=new WeakMap,o=new Set;function s(e){if(!(this instanceof s))throw new TypeError(\"Constructor requires 'new' operator\");i.set(this,e)}function h(){throw new TypeError(\"Function is not a constructor\")}function c(e,t,i,n){e=0 in arguments?Number(arguments[0]):0,t=1 in arguments?Number(arguments[1]):0,i=2 in arguments?Number(arguments[2]):0,n=3 in arguments?Number(arguments[3]):0,this.right=(this.x=this.left=e)+(this.width=i),this.bottom=(this.y=this.top=t)+(this.height=n),Object.freeze(this)}function d(){t=requestAnimationFrame(d);var s=new WeakMap,p=new Set;o.forEach((function(t){r.get(t).forEach((function(i){var r=t instanceof window.SVGElement,o=a.get(t),d=r?0:parseFloat(o.paddingTop),f=r?0:parseFloat(o.paddingRight),l=r?0:parseFloat(o.paddingBottom),u=r?0:parseFloat(o.paddingLeft),g=r?0:parseFloat(o.borderTopWidth),m=r?0:parseFloat(o.borderRightWidth),w=r?0:parseFloat(o.borderBottomWidth),b=u+f,F=d+l,v=(r?0:parseFloat(o.borderLeftWidth))+m,W=g+w,y=r?0:t.offsetHeight-W-t.clientHeight,E=r?0:t.offsetWidth-v-t.clientWidth,R=b+v,z=F+W,M=r?t.width:parseFloat(o.width)-R-E,O=r?t.height:parseFloat(o.height)-z-y;if(n.has(t)){var k=n.get(t);if(k[0]===M&&k[1]===O)return}n.set(t,[M,O]);var S=Object.create(h.prototype);S.target=t,S.contentRect=new c(u,d,M,O),s.has(i)||(s.set(i,[]),p.add(i)),s.get(i).push(S)}))})),p.forEach((function(e){i.get(e).call(e,s.get(e),e)}))}return s.prototype.observe=function(i){if(i instanceof window.Element){r.has(i)||(r.set(i,new Set),o.add(i),a.set(i,window.getComputedStyle(i)));var n=r.get(i);n.has(this)||n.add(this),cancelAnimationFrame(t),t=requestAnimationFrame(d)}},s.prototype.unobserve=function(i){if(i instanceof window.Element&&r.has(i)){var n=r.get(i);n.has(this)&&(n.delete(this),n.size||(r.delete(i),o.delete(i))),n.size||r.delete(i),o.size||cancelAnimationFrame(t)}},A.DOMRectReadOnly=c,A.ResizeObserver=s,A.ResizeObserverEntry=h,A}; // eslint-disable-line\n", + "mpl.toolbar_items = [[\"Home\", \"Reset original view\", \"fa fa-home icon-home\", \"home\"], [\"Back\", \"Back to previous view\", \"fa fa-arrow-left icon-arrow-left\", \"back\"], [\"Forward\", \"Forward to next view\", \"fa fa-arrow-right icon-arrow-right\", \"forward\"], [\"\", \"\", \"\", \"\"], [\"Pan\", \"Left button pans, Right button zooms\\nx/y fixes axis, CTRL fixes aspect\", \"fa fa-arrows icon-move\", \"pan\"], [\"Zoom\", \"Zoom to rectangle\\nx/y fixes axis, CTRL fixes aspect\", \"fa fa-square-o icon-check-empty\", \"zoom\"], [\"\", \"\", \"\", \"\"], [\"Download\", \"Download plot\", \"fa fa-floppy-o icon-save\", \"download\"]];\n", + "\n", + "mpl.extensions = [\"eps\", \"jpeg\", \"pdf\", \"png\", \"ps\", \"raw\", \"svg\", \"tif\"];\n", + "\n", + "mpl.default_extension = \"png\";/* global mpl */\n", + "\n", + "var comm_websocket_adapter = function (comm) {\n", + " // Create a \"websocket\"-like object which calls the given IPython comm\n", + " // object with the appropriate methods. Currently this is a non binary\n", + " // socket, so there is still some room for performance tuning.\n", + " var ws = {};\n", + "\n", + " ws.close = function () {\n", + " comm.close();\n", + " };\n", + " ws.send = function (m) {\n", + " //console.log('sending', m);\n", + " comm.send(m);\n", + " };\n", + " // Register the callback with on_msg.\n", + " comm.on_msg(function (msg) {\n", + " //console.log('receiving', msg['content']['data'], msg);\n", + " // Pass the mpl event to the overridden (by mpl) onmessage function.\n", + " ws.onmessage(msg['content']['data']);\n", + " });\n", + " return ws;\n", + "};\n", + "\n", + "mpl.mpl_figure_comm = function (comm, msg) {\n", + " // This is the function which gets called when the mpl process\n", + " // starts-up an IPython Comm through the \"matplotlib\" channel.\n", + "\n", + " var id = msg.content.data.id;\n", + " // Get hold of the div created by the display call when the Comm\n", + " // socket was opened in Python.\n", + " var element = document.getElementById(id);\n", + " var ws_proxy = comm_websocket_adapter(comm);\n", + "\n", + " function ondownload(figure, _format) {\n", + " window.open(figure.canvas.toDataURL());\n", + " }\n", + "\n", + " var fig = new mpl.figure(id, ws_proxy, ondownload, element);\n", + "\n", + " // Call onopen now - mpl needs it, as it is assuming we've passed it a real\n", + " // web socket which is closed, not our websocket->open comm proxy.\n", + " ws_proxy.onopen();\n", + "\n", + " fig.parent_element = element;\n", + " fig.cell_info = mpl.find_output_cell(\"
\");\n", + " if (!fig.cell_info) {\n", + " console.error('Failed to find cell for figure', id, fig);\n", + " return;\n", + " }\n", + " fig.cell_info[0].output_area.element.on(\n", + " 'cleared',\n", + " { fig: fig },\n", + " fig._remove_fig_handler\n", + " );\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_close = function (fig, msg) {\n", + " var width = fig.canvas.width / fig.ratio;\n", + " fig.cell_info[0].output_area.element.off(\n", + " 'cleared',\n", + " fig._remove_fig_handler\n", + " );\n", + " fig.resizeObserverInstance.unobserve(fig.canvas_div);\n", + "\n", + " // Update the output cell to use the data from the current canvas.\n", + " fig.push_to_output();\n", + " var dataURL = fig.canvas.toDataURL();\n", + " // Re-enable the keyboard manager in IPython - without this line, in FF,\n", + " // the notebook keyboard shortcuts fail.\n", + " IPython.keyboard_manager.enable();\n", + " fig.parent_element.innerHTML =\n", + " '';\n", + " fig.close_ws(fig, msg);\n", + "};\n", + "\n", + "mpl.figure.prototype.close_ws = function (fig, msg) {\n", + " fig.send_message('closing', msg);\n", + " // fig.ws.close()\n", + "};\n", + "\n", + "mpl.figure.prototype.push_to_output = function (_remove_interactive) {\n", + " // Turn the data on the canvas into data in the output cell.\n", + " var width = this.canvas.width / this.ratio;\n", + " var dataURL = this.canvas.toDataURL();\n", + " this.cell_info[1]['text/html'] =\n", + " '';\n", + "};\n", + "\n", + "mpl.figure.prototype.updated_canvas_event = function () {\n", + " // Tell IPython that the notebook contents must change.\n", + " IPython.notebook.set_dirty(true);\n", + " this.send_message('ack', {});\n", + " var fig = this;\n", + " // Wait a second, then push the new image to the DOM so\n", + " // that it is saved nicely (might be nice to debounce this).\n", + " setTimeout(function () {\n", + " fig.push_to_output();\n", + " }, 1000);\n", + "};\n", + "\n", + "mpl.figure.prototype._init_toolbar = function () {\n", + " var fig = this;\n", + "\n", + " var toolbar = document.createElement('div');\n", + " toolbar.classList = 'btn-toolbar';\n", + " this.root.appendChild(toolbar);\n", + "\n", + " function on_click_closure(name) {\n", + " return function (_event) {\n", + " return fig.toolbar_button_onclick(name);\n", + " };\n", + " }\n", + "\n", + " function on_mouseover_closure(tooltip) {\n", + " return function (event) {\n", + " if (!event.currentTarget.disabled) {\n", + " return fig.toolbar_button_onmouseover(tooltip);\n", + " }\n", + " };\n", + " }\n", + "\n", + " fig.buttons = {};\n", + " var buttonGroup = document.createElement('div');\n", + " buttonGroup.classList = 'btn-group';\n", + " var button;\n", + " for (var toolbar_ind in mpl.toolbar_items) {\n", + " var name = mpl.toolbar_items[toolbar_ind][0];\n", + " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", + " var image = mpl.toolbar_items[toolbar_ind][2];\n", + " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", + "\n", + " if (!name) {\n", + " /* Instead of a spacer, we start a new button group. */\n", + " if (buttonGroup.hasChildNodes()) {\n", + " toolbar.appendChild(buttonGroup);\n", + " }\n", + " buttonGroup = document.createElement('div');\n", + " buttonGroup.classList = 'btn-group';\n", + " continue;\n", + " }\n", + "\n", + " button = fig.buttons[name] = document.createElement('button');\n", + " button.classList = 'btn btn-default';\n", + " button.href = '#';\n", + " button.title = name;\n", + " button.innerHTML = '';\n", + " button.addEventListener('click', on_click_closure(method_name));\n", + " button.addEventListener('mouseover', on_mouseover_closure(tooltip));\n", + " buttonGroup.appendChild(button);\n", + " }\n", + "\n", + " if (buttonGroup.hasChildNodes()) {\n", + " toolbar.appendChild(buttonGroup);\n", + " }\n", + "\n", + " // Add the status bar.\n", + " var status_bar = document.createElement('span');\n", + " status_bar.classList = 'mpl-message pull-right';\n", + " toolbar.appendChild(status_bar);\n", + " this.message = status_bar;\n", + "\n", + " // Add the close button to the window.\n", + " var buttongrp = document.createElement('div');\n", + " buttongrp.classList = 'btn-group inline pull-right';\n", + " button = document.createElement('button');\n", + " button.classList = 'btn btn-mini btn-primary';\n", + " button.href = '#';\n", + " button.title = 'Stop Interaction';\n", + " button.innerHTML = '';\n", + " button.addEventListener('click', function (_evt) {\n", + " fig.handle_close(fig, {});\n", + " });\n", + " button.addEventListener(\n", + " 'mouseover',\n", + " on_mouseover_closure('Stop Interaction')\n", + " );\n", + " buttongrp.appendChild(button);\n", + " var titlebar = this.root.querySelector('.ui-dialog-titlebar');\n", + " titlebar.insertBefore(buttongrp, titlebar.firstChild);\n", + "};\n", + "\n", + "mpl.figure.prototype._remove_fig_handler = function (event) {\n", + " var fig = event.data.fig;\n", + " if (event.target !== this) {\n", + " // Ignore bubbled events from children.\n", + " return;\n", + " }\n", + " fig.close_ws(fig, {});\n", + "};\n", + "\n", + "mpl.figure.prototype._root_extra_style = function (el) {\n", + " el.style.boxSizing = 'content-box'; // override notebook setting of border-box.\n", + "};\n", + "\n", + "mpl.figure.prototype._canvas_extra_style = function (el) {\n", + " // this is important to make the div 'focusable\n", + " el.setAttribute('tabindex', 0);\n", + " // reach out to IPython and tell the keyboard manager to turn it's self\n", + " // off when our div gets focus\n", + "\n", + " // location in version 3\n", + " if (IPython.notebook.keyboard_manager) {\n", + " IPython.notebook.keyboard_manager.register_events(el);\n", + " } else {\n", + " // location in version 2\n", + " IPython.keyboard_manager.register_events(el);\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype._key_event_extra = function (event, _name) {\n", + " var manager = IPython.notebook.keyboard_manager;\n", + " if (!manager) {\n", + " manager = IPython.keyboard_manager;\n", + " }\n", + "\n", + " // Check for shift+enter\n", + " if (event.shiftKey && event.which === 13) {\n", + " this.canvas_div.blur();\n", + " // select the cell after this one\n", + " var index = IPython.notebook.find_cell_index(this.cell_info[0]);\n", + " IPython.notebook.select(index + 1);\n", + " }\n", + "};\n", + "\n", + "mpl.figure.prototype.handle_save = function (fig, _msg) {\n", + " fig.ondownload(fig, null);\n", + "};\n", + "\n", + "mpl.find_output_cell = function (html_output) {\n", + " // Return the cell and output element which can be found *uniquely* in the notebook.\n", + " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", + " // IPython event is triggered only after the cells have been serialised, which for\n", + " // our purposes (turning an active figure into a static one), is too late.\n", + " var cells = IPython.notebook.get_cells();\n", + " var ncells = cells.length;\n", + " for (var i = 0; i < ncells; i++) {\n", + " var cell = cells[i];\n", + " if (cell.cell_type === 'code') {\n", + " for (var j = 0; j < cell.output_area.outputs.length; j++) {\n", + " var data = cell.output_area.outputs[j];\n", + " if (data.data) {\n", + " // IPython >= 3 moved mimebundle to data attribute of output\n", + " data = data.data;\n", + " }\n", + " if (data['text/html'] === html_output) {\n", + " return [cell, data, j];\n", + " }\n", + " }\n", + " }\n", + " }\n", + "};\n", + "\n", + "// Register the function which deals with the matplotlib target/channel.\n", + "// The kernel may be null if the page has been refreshed.\n", + "if (IPython.notebook.kernel !== null) {\n", + " IPython.notebook.kernel.comm_manager.register_target(\n", + " 'matplotlib',\n", + " mpl.mpl_figure_comm\n", + " );\n", + "}\n" + ], "text/plain": [ - "
" + "" ] }, - "metadata": { - "needs_background": "light" + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "text/html": [ + "" + ], + "text/plain": [ + "" + ] }, + "metadata": {}, "output_type": "display_data" + }, + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" } ], "source": [ - "w = np.linspace(-6, 6, 2**12)\n", - "gloc = semi_circle_hiltrans(w + 1e-3j)\n", + "U = 4.0\n", + "D = 4.0\n", + "beta = 10000000000\n", + "w = np.linspace(-6*D, 6*D, 2**12*int(D))\n", + "gloc = semi_circle_hiltrans(w + 1e-3j,D)\n", "#gloc= 0.0\n", - "gloc, sigma_loc = dmft_loop(gloc, w, 0.3, 100, 100)\n", - "\n", - "#plt.plot(w, -gloc.imag, lw=3, label=r'$\\pi A(\\omega)$')\n", - "plt.plot(w, sigma_loc.real, '-', label=r'$Re \\Sigma(\\omega)$')\n", - "plt.plot(w, sigma_loc.imag, '-', label=r'$Im \\Sigma(\\omega)$')\n", - "plt.xlabel(r'$\\omega$')\n", + "gloc, sigma_loc = dmft_loop(gloc, w, U, beta, D, 100)\n", + "fig,ax = plt.subplots(2,1,sharex=True)\n", + "ax[0].plot(w, -gloc.imag, lw=3, label=r'$\\pi A(\\omega)$')\n", + "ax[1].plot(w, sigma_loc.real, '-', label=r'$Re \\Sigma(\\omega)$')\n", + "ax[1].plot(w, sigma_loc.imag, '-', label=r'$Im \\Sigma(\\omega)$')\n", + "ax[1].set_xlabel(r'$\\omega$')\n", "plt.legend(loc=0)\n", "#plt.ylim([0, 0.1])\n", - "#plt.xlim([-1, 1])\n", - "plt.show()" + "#plt.xlim([-1, 1])\n" ] }, { @@ -296,7 +1269,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.3" + "version": "3.8.10" } }, "nbformat": 4, diff --git a/legacy/PYTHON/real_ipt-text_v3.ipynb b/legacy/PYTHON/real_ipt-text_v3.ipynb index eb25669..c55107c 100644 --- a/legacy/PYTHON/real_ipt-text_v3.ipynb +++ b/legacy/PYTHON/real_ipt-text_v3.ipynb @@ -3,9 +3,7 @@ { "cell_type": "code", "execution_count": 1, - "metadata": { - "collapsed": true - }, + "metadata": {}, "outputs": [], "source": [ "%matplotlib inline" @@ -43,7 +41,7 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ @@ -65,7 +63,7 @@ " return 1. / (np.exp(exponent) + 1)\n", "\n", "\n", - "def semi_circle_hiltrans(zeta, D=1):\n", + "def semi_circle_hiltrans(zeta, D):\n", " \"\"\"Calculate the Hilbert transform with a semicircular DOS \"\"\"\n", " sqr = np.sqrt(zeta**2 - D**2)\n", " sqr = np.sign(sqr.imag) * sqr\n", @@ -84,7 +82,7 @@ " return -np.pi * U**2 * (Appp + Appp[::-1])\n", "\n", "\n", - "def dmft_loop(gloc, w, U, beta, loops):\n", + "def dmft_loop(gloc, w, U, beta, D, loops):\n", " \"\"\"DMFT Loop for the single band Hubbard Model at Half-Filling\n", "\n", "\n", @@ -98,6 +96,8 @@ " On site interaction, Hubbard U\n", " beta : float\n", " Inverse temperature\n", + " D : float\n", + " Half-bandwidth of the noninteracting system\n", " loops : int\n", " Amount of DMFT loops to perform\n", "\n", @@ -116,7 +116,7 @@ "\n", " for i in range(loops):\n", " # Self-consistency\n", - " g0 = 1 / (w + eta - .25 * gloc)\n", + " g0 = 1 / (w + eta - .25 * D**2 * gloc)\n", " # Spectral-function of Weiss field\n", " A0 = -g0.imag / np.pi\n", " # Clean for PH and Half-fill\n", @@ -130,7 +130,7 @@ " sigma = hsi + 1j * isi\n", "\n", " # Semi-circle Hilbert Transform\n", - " gloc = semi_circle_hiltrans(w - sigma)\n", + " gloc = semi_circle_hiltrans(w - sigma, D)\n", "\n", " return gloc, sigma" ] @@ -144,16 +144,26 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": 7, "metadata": { - "scrolled": true + "scrolled": false }, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAEUCAYAAAAoQI39AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XmcXFWd9/HPr6q6uzpJdxKyJ52NQCIhhCSERRAQQcQoy6PiMKCCwUEYUXnUxxF5ZBhHHkEdVEZlyACjYFQUM8ooiiDLkJGdLE0SAzEs6bBkMSGd9FbLef64t6pvJb1UdVWl63Z9369XU3W3uqfCre89derce8w5h4iIDH2RwS6AiIgcHAp8EZEqocAXEakSCnwRkSqhwBcRqRIKfBGRKqHAFxGpEiULfDOLmtkqM/uNPz3TzJ40sxfN7G4zqy3VvkREpHClrOF/FtgQmL4R+LZz7nBgF3BpCfclIiIFslJcaWtmTcCPgOuBzwFnA9uBic65pJm9HbjOOfeevl5n7NixbsaMGUWXR0Skmjz77LM7nHPj+lsvVqL9fQf4ItDgT48Bdjvnkv50CzClvxeZMWMGzzzzTImKJCJSHczslXzWK7pJx8zeD2xzzj0bnN3Dqj1+lTCzy8zsGTN7Zvv27cUWR0REelGKNvyTgHPM7GXgZ8C78Gr8o8ws8w2iCXitp42dc8ucc4udc4vHjev3G4mIiAxQ0YHvnLvaOdfknJsBXAA85Jy7CHgY+JC/2sXAr4vdl4iIDFyp2vB78g/Az8zsa8Aq4PYy7ktEhoBEIkFLSwsdHR2DXZSKFI/HaWpqoqamZkDblzTwnXOPAI/4zzcDx5Xy9UVkaGtpaaGhoYEZM2Zg1tNPgdXLOcfOnTtpaWlh5syZA3oNXWkrIhWjo6ODMWPGKOx7YGaMGTOmqG8/CnwRqSgK+94V+2+jwBcJ2N7aye+aXx/sYoiUhQJfJOCKHz/LFcufY9e+rsEuikjJKfBFAl7b3Q7A3s5kP2uKhI8CXySg+DtLyVBw6623MnHiRI4++mhmzZrFnXfe2e82sViM2bNns2DBAhYsWMBXv/pVANrb2zn11FNJpVL9vkZXVxennHIKyWR5KhwKfJGAzL0EU2lFfzVbu3Yt1113HWvWrOGnP/0pn/vc5/rd5sYbb+TII49k1apVrF69mmuvvRaAO+64gw984ANEo9F+X6O2tpbTTz+du+++u+j30JNyXnglEjqZThBJBf6g+6f/Wsf61/aU9DXnTm7kH88+st/1mpub+dCHvBsFzJw5k9pabziPl156iauuuoqtW7cSiUS46667mDNnDgCf//znee655/ja177GV77ylexrLV++nJ/85CfZ6XvvvZc777yTe+65JzvvlltuYcOGDdx8882cd955XH311Vx00UUlec9BquGLBET8xFcNv7o1NzczZ84cnHN873vf4/rrryeRSPCJT3yCm266iWeeeYbrrruOG264IWe72267jbvuuovHHnsM8JpoNm/eTPC279dccw3XXXddznazZs1i/fr1AMybN4+nn366LO9LNXyRgGjEC/xkOp3X+utf28Nh40dQG1PdqdTyqYmXw5YtW2htbWXJkiVs3bqV+fPnc91113HPPfewbt06PvjBDwKQTCY5+eSTc7bdsGEDqVSKo48+GoAdO3YwatSo7PI1a9aQTqeZN28er7zyCvfddx9XXHEFiUQi28c+Go1SW1tLa2srDQ0NlJICXyTAz/u8avgtu9pYcvNjXHT8NK7/X0eVuWRysKxdu5ZTTjmFhx56iF27djFv3jwef/xx1qxZw/XXX8+ll/Y8eF9nZyeXXHIJt9xyC42NjQDU19fnXBm7evVqjjnmGAAeeOABXnzxRQDWr1+fPUlkXisej5f8valaIhJQSJPO7rYEAM++squsZZKDq7m5mYULFwIwevRoLrzwQn77298yadIk7r//ftL+t7/m5maCIwZee+21nHDCCZx55pnZeaNHjyaVSmVDP51Os3fvXlKpFCtWrKC1tZX29nZ++MMfcuGFFwKwc+dOxo0bN+AbpPVFgS8SEInkH/i6A8DQFAx8gLPPPpv77ruPpUuXkk6nOeKII1iwYAE33nhjzq0OvvWtb/Hoo49mu2V+5jOfAeDMM89k5cqVACxZsoTNmzezYMECLr/8ctatW8fixYu57LLLWLRoEQAPP/wwS5YsKct7U5OOSEDUMm34/Qd+CYaDlgq0fPnynOlTTjmFVatWAeT0rNlfb/3sr7zySm666SbOOOMMJkyYwOrVq7PLzjnnnAPW/8lPfsLXv/71gRS9X6rhiwRYAW34nUnvq72CX/qycOFCTjvttLwvvDrvvPOyXT1LTTV8kYBoAU06iZQX+GklvvRj6dKlea1XW1vLxz72sbKVQzV8kYBCAj+T87pIS8JCgS8SECmkDd+/805XMr8++yKDTYEvEtBdw88jxP1zQqZpR6TSKfBFArovvOp/3cx3ADXpSFgo8EUCupt0+k/8TBt+Qk06EhIKfJGAgn60zbThq0lHQkKBLxLQffO0/HvpqA1fwkKBLxJgBdxLJ7OGmvCHnltvvZXLL7887/V7G+0KKmvEK114JRIQLWAAFKcLroastWvXctRR+d8B9cYbb2TlypWsWLEi5/46MPARr8oxAIoCXyQg24afRzON4r7MfvcleKO5tK858Sh47w39rtbc3MwFF1zA+eefn73/zZYtW1i+fDnLli3jiSee4OSTT+b2228Heh/tCnJHvOpvtCtAI16JHCyFtOEr8Yeu559/nqOOOorm5mYOPfRQVq5cycUXX8yll17KjTfeyPPPP8+KFSvo7OzMbrP/aFdw4IhX/Y12BRrxSuSgiUW8OlAilX8vHSmTPGri5bBlyxZGjBhBPB5n9+7dXHXVVYA3mMmll17KpEmTABg2bFh2rFs4cLQryB3xKp/RrqC8I16phi8SkK3h59Oko7wfkjLt9+vWrWPRokVE/ErAmjVrOP744wFoaWlh8uTJ2aDuabQryB3xKt/RrjKvpxGvRMpsIN0yZWhpbm5m/vz5NDc35wTx2rVrmT9/PuCFf+Y59DzaFeSOeJXPaFegEa9EDppMz5u8rrQtd2FkUDQ3N2fb7xcsWABAR0cH7e3tjB49GsgNf+h9tCvoHvEqn9GuoMJHvDKzqcCdwEQgDSxzzn3XzA4B7gZmAC8DH3bOafBPqWjZ++Pk04YfqOI75w7ojifhlBnxKljrjsfjvPTSS9npq6++OmebvvrYZ0a8uuuuu/od7Qoqf8SrJPB559wRwAnAp8xsLvAl4I/OucOBP/rTIhWtkHvcB9fQDdSkN5U04lXRge+ce90595z/vBXYAEwBzgV+5K/2I+C8YvclUm7dNfzCfrTVPfGlL0uXLs37wqvQjHhlZjOAhcCTwATn3OvgnRSA8b1sc5mZPWNmz2zfvr2UxREpWKaZJpFXjb17HQW+hEHJAt/MRgC/BK5yzu3Jdzvn3DLn3GLn3OJx48aVqjgiAzLgGr5uoFYyumVF74r9tylJ4JtZDV7YL3fOrfBnv2lmk/zlk4BtpdiXSFkNsA1fNfzSiMfj7Ny5U6HfA+ccO3fuLKp/fil66RhwO7DBOXdTYNG9wMXADf7jr4vdl0i5Za6eza+XTvdz1fBLo6mpiZaWFtS827N4PE5TU9OAty/FrRVOAj4KNJtZps/Rl/GC/udmdinwKnB+CfYlUlbdvXTy6YevNvxSq6mpYebMmYNdjCGr6MB3zq0EeuuAfHqxry9yMKXdwGr4GgRFwkBX2ooEDLQffntX/32sRQabAl8kIBPi+dTYgz8stnaUZ4QikVJS4IsEZGv4eTTpBLV2JspQGpHSUuCL5PCCPq8xbQOrqIYvYaDAFwnIhHiiwF46CnwJAwW+SEBhd8vsfq7AlzBQ4IsEZO+lU+CtFd5qVxu+VD4FvkhAJsPzasP3HxvjMV5/q71sZRIpFQW+SEC2Db+AbplTDxnGlr+2lbNYIiWhwBcJyNTa2xP9X0iVWXf6mGFs2dWuq22l4inwRQIytfa2zjyunPUT/6gpo+hKptn4RmsZSyZSPAW+SA/aEql+b9Gb6ZZ5zHRvYOtVr2rIZqlsCnyRgEzGp9KOzn7ugJlZt2l0PeMa6lj16u4yl06kOAp8kYB0oFbfvPWtPtfNrGkG8yY3skFNOlLhFPgiAcFWnPP/7fG81jWMwyc08Jfte/PqzikyWBT4IgGO/AM7s64ZHDZ+BF3JtLpnSkVT4IsEOAfHzTyEi98+HYCOPrpndtfwYfaEBgBeeFPNOlK5FPgiAQ4vwOdNGQnAG2919LkueBvMmdBAbTTCs6+op45ULgW+SJDzmmgmjowD8Oae3gM/U8U3jPraKIumj+KxF3ccjFKKDIgCXyTA4TCMCY1+4Ld29rGux/wRnU+dPZ71r+/h1Z1qx5fKpMAXCXB+DT8T+Nv6qOEH2/ABzl0wGTNYsaqlzKUUGRgFvkiAwwv8xniMeE2k7zb8TJOOX8WfPKqeE2eN4RfPtNCZ1KDmUnkU+CIBznlNOmZes05eTTqBeZ88ZRZbd7dzx8qXy1lMkQFR4IsEZGr44DXr9PWjbbZJJ5D4Jx8+lvccOYF/+cNGntO9daTCKPBFAoJX2k4aGWfrrt4HNumu4XcnvpnxjQ8ezaRRca5c/pwGRpGKosAXCfBq+F6Av21iI1t3t7O7ravndff/1dY3clgNt1x0DHs6kvztsid4ece+MpZYJH8KfJH9ZPJ7fpN38VV/N1EzO3DevCkj+dHS49jdnuDsf13JL59t6fd2yyLlpsAXCQqE8lFNI4lGjCc3/7XPVXvIe8C7T/5/XfkODp8wgs//Yg0fvvVx1mzRLZRl8CjwRQKCP9o2xmtYNG0Uj7ywrZd1c7tl9mTqIcO45/IT+cYH5/OX7fs49/v/w0W3PcFjL25XjV8OuthgF0CkkjiXW2N/55zxfPP+jWxr7WB8Q/yAdaH3Gn5GJGJ8+NipvPeoifzkyVe5feVLfPT2p5g7qZFPnnoo7ztqErHogXWvP/1lB7c99hIRg7mTGjll9jgWTB3V47oi+Sh74JvZWcB3gShwm3PuhnLvU2SgHC6nxn7q7HF88/6N/PcLO/jQMU37revpo4KfoyFewydPncUlJ83gV6u2suy/N/PZn63mG7/fyKXvmMnfHDuV4XUxEqk0N//xRb738CbGN9Qxsr6Ghzdu5+aHNjGiLsbCaaOYNW4EdTURdu3rYldbgrfaE6TTjoZ4jBHxGuKxCPW1UepiESIRI2JGxCBq3jUGETOiEXKeR/znNbEItVGjJhqhNhbxHgPPa6JGXSxKvMbbR31NlHiNt6++vu3I4Ctr4JtZFPg+8G6gBXjazO51zq0v535FBmr/Gv6RkxsZ11DHIxu3HRj4gQFQClEXi/I3x07j/GOm8tCft7Hsvzfz1d+s5zsPvsCZR05kbctuXnhzL+cf08Q/nXskw2pjvNWeYOWLO3h88w6eeXkXq7fspjORZvTwGkYPq6WxvoaaaISd+7p4eWcbHYmU/5cm7Zz/543oVa6WJDOIx6LU10aJxyLEAyeD7GNmWU00e/KIRY1YJPM8QizinWxiUaMm4j3GohFqIv5yf340Yj1uE4tkTmLeX/a5GZEI2flRf341naTKXcM/DtjknNsMYGY/A84FFPhSkTL30skwM06dPY4H1r9JMpXOaU4JDoAyEJGIccbcCZwxdwLPvbqL21e+xAPr32TaIcP4t48cw1nzJmbXHVlfw/vmT+J98ycNbGcBLhD+aedIp7ufp9KORMrRlUqTSKZJpNJ0+o+JlKMrO887mbQnUrR3pehIpujoStGeCMxPpOj0H9u6kuzc15Wdbk+kSKYciVSaZNoN6khh5n/ziQROAhHjgBNG8CTSfQLxvh1lnme+SZn/mPnWZAfMy502gzPnTuS8hVPK+l7LHfhTgC2B6Rbg+JLvZfOj8PD/855nP305n9rceft/Qntd3tP6/a0zkNcoUTnMwCJgUYhkHqP7PeY5P1YLtSOgdjjUDPOfD/Onh3fPj1RIe3LHW7BjE7TvglQnpBIQi0OsDmrqvfI3TIJhh/SZ0JkbHgPwVgu88jj/Z+ePOSe5kz13/ZBDZh0Db3s/jJuTX03ZOejcAx17ch9TXd3/3vWHsOiQGSy6cFEJ/iH6Z2ZEDaK9fTNxDpKdkExAoh0Sbf5fO7g2SLcDbd5zl3nuP8bSEE1DXRpc5s8Fngemcdlp7ySU9h7TaW8V5307cWmHw5FOe+t5873lzrnuv3TgOc573cyjw9+ft5ye1sm8d+eA/bbxlzkHpMAlvRN+5hhw2f9kKgIuMN39z5r9Jw78JzN7t50HC79UxP/Z/pU78Hs6onI+JmZ2GXAZwLRp0wa0k+da9pDa2obhsMBOzchOG9539YhfKvNvg2u47HoAEX+97u0ITHvDY+y/TU/TmTdvmX2ZZV/De+72K6fbb3/d5cPAHIF99LY/h7k05tKQToFLBR7T+037j8WqGZZ7Uqip9/5icaiJ+6Eb9+fVQazemx+t86ajtd728ZEQHwX1o7yw/utm2PMatP/V23bsbJhxMtSNyN1/6xvw28/Dxt/l936itTB6Bkw9HuaeC4eeBtHuj4Fz3r8tqSTcegq07WRcYxOvWT28sRZe/i388avwyceAYQe+vnPw5jpo/gX85Y+wczMk8rzwqulYOOsGaFqcOz+VhC1PwuurvRNbXSMMHwcNE72TWE097N0G+7bBvh3ev1miHZIdkOjwHjN/iQ5ItnuBnvAfc6b99QYiFodIzK90+EezRXr5s5x1DK+WnK245FRoMvPoYznd05E8t8++Tn/7hL4qCb0rcJs5xX9760+5A78FmBqYbgJeC67gnFsGLANYvHjxgL7XRWaezPI3pvlnanLaKzNfWVPOqyGknSPlvA92Kp37lTaV2S4dWBaYl10nHZje7+uwc2S3TaUdyUH4qhqNGHWxiP8XpTYWoa4uQl2NN51ZVhuLEI8a8RgMr4H6mFEfddTVRBgRSTLcuhhmHd6f6yBOB3HXQZ3roDbVTizVTizZRiS5j2iynUhiH5FkG5bsxNp29h4yBYwbm6O2AY77Ozjps2ztrGP5j+/g0u03MNw62XL4xzl00WlER4z3TiSRmB9mfhk63oLWN6H1Ndi+EdbfC6vugsYp8IFlMOMd2d0YwN43oW0nLF5K5KwbuelHq9i6q52HPl4H/3EWvPI/eD9NBbLgry/BLy7xgtmi3msuegeMnOKd0OoaId4IdSMhWtN9Im7bAW8+D0/fDrefCadfCyd+xnvh9b+CP3wF3tpCwSI1gROvf8INnoQz/1bZZfudlDMn6prh/uMw/9vSsO7nmcdYfeV825NelTvwnwYON7OZwFbgAuDCUu9kwdRRLJg6qtQvWzJe8Kez7aOptCMZaLtMpHKXJfx1kylvu6T/PBV4ngy8RvYx7bWxdiZT/mOazoQ33elPZ5bv7Uxml3Uk0l4brN/+2rO4/5cfM7I/uEUzP6hFjJo6oz6aIh5JMsxSxKNJhtNBg9tLfaqV2sQe2pPGptQ43qqZwPSmJo6fWs/Jw1/l0FfuwVZ+m8STt/NyYgZfZA1v1h/K30c/x0PNhzDrjeFc+a7JvH/+ZGr667qY7IQX/+DV1n/8QVj6e5i8MFt29vj1ktlnQayWd8+dwLW/Xsem+DEcVjcSdv4FV3+Gtz7m1az/Y4lXm1/yLZh7HowYl/e/F7PfA8d+Au79DDz4j/Dn33g14S1PwoR5cOY/w4xToH601yS0bwe0vu59w0m0wYgJ3v6GjfXWqR3uNReJBJQ18J1zSTO7Ergfr1vmHc65deXcZyXyfuQJx4fPOUdnMp0N/45E5oe4wInB/5EukfRPPGnvB73uE1AP89Jp/0c678SVyJyoUo5252hz3r/TsNoow2qjvK02xq62LtZs2c3vN+wAjGG1F3DCsFP46L4fMq/mNd5a+GkmnHUN/x6N87vnX+d7D23if9+9hm/d/wJL/W6OI+p6OcRjdXDE2TDt7fCDt8MjN8CFd3e3s+7Z6j02TgbgjCO8wH9gwzYOGzEe9m3H+ec/M7zaeevrcNkjMHnBwP7x4yPh/B/Ccz+Cx7/vtW+/719g0SU5zU7U+01fYw8b2H6kapW9H75z7j7gvnLvR0rDzIj7XegqxbbWDp5+aRdPvbSTN/aM5dVZd3Ly8dOyPWaiwPvnT2bJvEleN8fHNvPPv1nPdx98gYtOmM4lJ87IjmB1gOFjYe45sOZuSKezQxzS+rq3vNHrNTF5VD1zJjTwp7/s4IrhY2HfDtwYbxUDr4ln8oKBh32GGRxzifcnUmK60lYq3viGeF5dEoPdHFe9uot/f2wztz76F257bDMXHT+dL7xnTs81/onz4enbYPfLOAdTujbDquXeD8v1o7OrHX/oIdzzbAvpI8YQ2fgbpo34EzDW68e9fSMcdnqJ37lIaelXFhmSFk4bzQ8uOoaHv/BOPnRMEz96/GXefdOjPPbi9gNXHj/Xe9z2ZxzwlS2fgDebveacQO+M42eOoa0rxV9pAOC8dZ8BwNp3w943YNycMr8rkeIo8GVImz5mOF//wHxWXHEiI+piXHzHU/z4iVdyV8oE9fYNuTc0G577o+uxM73a/qsd3d0x43Rib6z2JiYcWfLyi5SSAl+qwsJpo/nVp07inXPG839/9Ty/fLale2G8ERqbsjX8rGhNzmuMb4hz6NjhbNg7PDuvgXbs5ZVeN8yppb+mUKSUFPhSNYbXxbjlI4s46bAxXL2imXWvBQY2mXQ0bHkCgtdN2IEfj8PGj2BjW0N2us664OWVXpfOuoYD1hepJAp8qSp1sSj/+reLGDWshk//dBVtXUlvweFnwO5XmeYCFzgt/OgB288YO5zHW8dmpxusA7Y+C9PfXu6iixRNgS9V55DhtXznggVs3r6PWx/d7M2c7l1pOzu5iSQxWPgRmP/hA7adPmYYLyYnsPv0bwEwx7ZAOgHj1X4vlU+BL1XpxFljOXPuBH74p5fpSKRgpHfr4/FuOzGSXpt+D/dPmTKqHoBtEa+WPyvi99c/5NCDU3CRIijwpWp95ITpvNWe4JGN27z7w9Q1MsH53TZjtT1uM9kP/B1t3o3aRtPqLaiv3Ft7iGQo8KVqnThrDGNH1PGrVf59c2pH0Oj8AI/1fGXuxJHe/J1+4I80/06Y+sFWQkCBL1UrFo1w5pETWLlpB8lUGuoCgR/tuYbfGK9hRF2M7fu8H3tHkgn8xoNRZJGiKPClqr390DHs7Uzy/Gt7oHYEI5wf4H3c7G7SyDjb92Vq+HsB8+5OKVLhFPhS1Y6feQgAz72yC2qHM5w2b0EPffAzJo6M8+c3vfVGss9rzqmicVElvBT4UtXGNdQxalgNm7bvhVicWpfwl/Qe4MNqo2zzm3Tq6fRutSwSAgp8qWpmxuHjR7Dpzb0QraUWP/D7qOFPO2QYKbwmnzq6vJGlREJAgS9V77DxDbywrRUXraGWLm9mH000f//Ow0j5H504XbmDk4hUMAW+VL0jJjWwuy1BRzpCDf6tFvqo4Y8eXsvhE0YCELeEavgSGgp8qXpvm+h1qdzdCVEyY/r2/SNsJBao1UcV+BIOCnypenMmNmAGO9r7vlNmUDQY+KrhS0go8KXqjayv4cjJjWxtTXXP7KebZTRYq1cbvoSEAl8Er1nnrwXV8AOBrxq+hIQCXwTvLpi7u/K/eCoWDHy14UtIKPBF8K6eTRC4nUJBNXw16Ug4KPBF8Nrx0y7wceinDT8WC5wcVMOXkFDgi+DdBTOV83HoO/BrcwK/5ztrilQaBb4I0FgfIx0M+X6adGqCga8mHQkJBb4IXg0/TQFNOlE16Uj4KPBFgOF1sdwmnX5q+LFoYLlq+BISCnwRIF4Tya3h99OGH42qSUfCR4EvAsRrogXV8KM5I2Jp8BMJBwW+CFATjeSGfH+3Vojkv65IpSgq8M3sm2b2ZzNba2b/aWajAsuuNrNNZrbRzN5TfFFFyssiBVx4Fc3/24BIpSj2SH0AmOecmw+8AFwNYGZzgQuAI4GzgB+YWe+jQotUgGjOTdAKaMNX4EtIFHWkOuf+4JzzR4zgCaDJf34u8DPnXKdz7iVgE3BcMfsSKbdITogX0C1TgS8hUcojdSnwO//5FGBLYFmLP0+kYkULCfyIkXaW17oilaLf/mRm9iAwsYdF1zjnfu2vcw2QBJZnNuthfdfDPMzsMuAygGnTpuVRZJHyiASbdPrtpWPdB7Rq+BIS/Qa+c+6Mvpab2cXA+4HTnXOZz0ALMDWwWhPwWi+vvwxYBrB48eIeTwoiB0MhbfixqJEmQpSUAl9Co9heOmcB/wCc45xrCyy6F7jAzOrMbCZwOPBUMfsSKbdorIBeOqrhSwgVe4ng94A64AHz2jGfcM5d7pxbZ2Y/B9bjNfV8yjmX6uN1RAZdbpNO/234jgiQQhdeSVgUFfjOucP6WHY9cH0xry9yMEUL6YcfiQRq+Ap8CQd9FxXx5fTS6a8NP2Ld995Rk46EhI5UEV8kol46MrTpSBXxFXThVU4NX006Eg4KfBFfIbdLUA1fwkhHqoivoDb8qOEy6yjwJSR0pIr4Cru1QkSBL6GjI1XEFyngHvdek05mHbXhSzgo8EV8uYHf35i2gZBXDV9CQkeqiM+CgZ9HP/zAhuUpkEiJ6UgV8UULHOIw209H3TIlJBT4Ir5ChjhUDV/CSEeqiC8SLaBJJ6cNXzV8CQcFvogvYgX8aFvAD7wilUJHqojPCuiWWRO1QBu+PkYSDjpSRXyRAm+tkO+6IpVCR6qIL6eXTj9t+DUFtPeLVAoFvoivkCYd9dKRMNKRKuIr9NYK2TUU+BISOlJFfBMa67sn+glxM3XLlPBR4Iv4jpo6OjBVQIgr8CUkFPgiGQX0w+91O5EKpiNVJKvQZhr1w5dw0ZEqklFgDV8/2krY6EgVySigH763hmr4Ei46UkUybKB96/WjrYSDAl8kY6BdLVXDl5DQkSqSUWAbfn1NJO91RSqBjlSRjALb8LO3xFfgS0joSBXJGmiTjtrwJRwU+CIZBYxpm7udAl/CoSSBb2ZfMDNnZmNwB5iFAAAJ0UlEQVT9aTOzm81sk5mtNbNFpdiPSFkVeqWtU7dMCZeij1Qzmwq8G3g1MPu9wOH+32XALcXuR6TscmrqutJWhp5SHKnfBr5I9ugH4FzgTud5AhhlZpNKsC+R8hlok4764UtIFBX4ZnYOsNU5t2a/RVOALYHpFn+eSOXSzdNkiIv1t4KZPQhM7GHRNcCXgTN72qyHea6HeZjZZXjNPkybNq2/4oiUUYFNOmrDl5DpN/Cdc2f0NN/MjgJmAmv8wSCagOfM7Di8Gv3UwOpNwGu9vP4yYBnA4sWLezwpiBwUBd9aQYEv4TLgI9U51+ycG++cm+Gcm4EX8oucc28A9wIf83vrnAC85Zx7vTRFFikTdcuUIa7fGv4A3QcsATYBbcDHy7QfkdIZ6M3TVMOXkChZ4Pu1/MxzB3yqVK8tclAUeGuF7jZ81fAlHFQ1EckouJeO2vAlXHSkimQMuA1fHyMJBx2pIlkDvHmaLrySkFDgi2QUWlNXP3wJGR2pIhkDDW79aCshocAXyRhwcCvwJRwU+CIZA67h62Mk4aAjVSRjoDV8NelISCjwRQZMP9pKuOhIFSmaavgSDgp8kWKphi8hoSNVpFiq4EtIKPBFBkoXXknI6EgVKZqq+BIOCnyRYqmGLyGhI1WkWOqHLyGhwBcpmgJfwkGBLzJg+tFWwkVHqkix1KQjIaHAFymWavgSEjpSRYqmGr6EgwJfpFiq4UtI6EgVKZba8CUkFPgixVLgS0go8EWKpsCXcFDgixRLNXwJCQW+SLH0o62EhI5UkaKphi/hoMAXKZZq+BISOlJFiqU2fAkJBb5IsVTDl5Ao+kg1s0+b2UYzW2dm3wjMv9rMNvnL3lPsfkQql2r4Eg6xYjY2s9OAc4H5zrlOMxvvz58LXAAcCUwGHjSz2c65VLEFFqk4quFLSBR7pF4B3OCc6wRwzm3z558L/Mw51+mcewnYBBxX5L5EKpPa8CUkig382cDJZvakmT1qZsf686cAWwLrtfjzRIYgBb6EQ79NOmb2IDCxh0XX+NuPBk4AjgV+bmaH0vMnwPXy+pcBlwFMmzYtv1KLVBLV8CUk+g1859wZvS0zsyuAFc45BzxlZmlgLF6Nfmpg1SbgtV5efxmwDGDx4sU9nhREKpoCX0Ki2CadXwHvAjCz2UAtsAO4F7jAzOrMbCZwOPBUkfsSqUz60VZCoqheOsAdwB1m9jzQBVzs1/bXmdnPgfVAEviUeujI0KUavoRDUYHvnOsCPtLLsuuB64t5fZFQUA1fQkJHqkix1IYvIaHAFymaAl/CQYEvMlARv0U0Fh/ccojkqdgfbUWGlr/5cXeQ9+fvHoJX/gQR1ZskHBT4IkFHnJ3/upOO9v5EQkJVExGRKqHAFxGpEgp8EZEqocAXEakSCnwRkSqhwBcRqRIKfBGRKqHAFxGpEubdzbgymNl24JUBbj4W7178Q4HeS2UaKu9lqLwP0HvJmO6cG9ffShUV+MUws2ecc4sHuxyloPdSmYbKexkq7wP0XgqlJh0RkSqhwBcRqRJDKfCXDXYBSkjvpTINlfcyVN4H6L0UZMi04YuISN+GUg1fRET6MOQC38w+bWYbzWydmX1jsMtTLDP7gpk5Mxs72GUZKDP7ppn92czWmtl/mtmowS5TIczsLP+Y2mRmXxrs8gyUmU01s4fNbIP/+fjsYJepGGYWNbNVZvabwS5LMcxslJnd439GNpjZ28u1ryEV+GZ2GnAuMN85dyTwrUEuUlHMbCrwbuDVwS5LkR4A5jnn5gMvAFcPcnnyZmZR4PvAe4G5wN+a2dzBLdWAJYHPO+eOAE4APhXi9wLwWWDDYBeiBL4L/N459zbgaMr4noZU4ANXADc45zoBnHPbBrk8xfo28EUg1D+0OOf+4JxL+pNPAE2DWZ4CHQdscs5tds51AT/Dq1SEjnPudefcc/7zVrxgmTK4pRoYM2sC3gfcNthlKYaZNQKnALcDOOe6nHO7y7W/oRb4s4GTzexJM3vUzI4d7AINlJmdA2x1zq0Z7LKU2FLgd4NdiAJMAbYEplsIaUgGmdkMYCHw5OCWZMC+g1cZSg92QYp0KLAd+A+/eeo2Mxterp2FbkxbM3sQmNjDomvw3s9ovK+rxwI/N7NDXYV2RernvXwZOPPglmjg+novzrlf++tcg9essPxglq1I1sO8ijye8mVmI4BfAlc55/YMdnkKZWbvB7Y55541s3cOdnmKFAMWAZ92zj1pZt8FvgR8pVw7CxXn3Bm9LTOzK4AVfsA/ZWZpvPtTbD9Y5StEb+/FzI4CZgJrzAy8JpDnzOw459wbB7GIeevr/wuAmV0MvB84vVJPwL1oAaYGppuA1wapLEUzsxq8sF/unFsx2OUZoJOAc8xsCRAHGs3sx865jwxyuQaiBWhxzmW+ad2DF/hlMdSadH4FvAvAzGYDtYTwxkrOuWbn3Hjn3Azn3Ay8g2JRpYZ9f8zsLOAfgHOcc22DXZ4CPQ0cbmYzzawWuAC4d5DLNCDm1R5uBzY4524a7PIMlHPuaudck//ZuAB4KKRhj/+Z3mJmc/xZpwPry7W/0NXw+3EHcIeZPQ90AReHrDY5VH0PqAMe8L+xPOGcu3xwi5Qf51zSzK4E7geiwB3OuXWDXKyBOgn4KNBsZqv9eV92zt03iGUS+DSw3K9QbAY+Xq4d6UpbEZEqMdSadEREpBcKfBGRKqHAFxGpEgp8EZEqocAXEakSCnwRkSqhwBcRqRIKfBGRKqHAFxGpEgp8kf2Y2Xgz22FmXf0NEGJm3/VHJPvXg1U+kYFS4Isc6DvAGOAm51x/N7J62n88tbxFEime7qUjEmBmx+ONyrULmOmce6uf9U8E/gfvZn3DnHOp8pdSZGBUwxfJda3/+IP+wt6XGWuhFu9bgUjFUuCL+MxsOt5g5QD/vt+yaWb2DjPbP9SDn6FoOcsnUiwFvki3JXhDGq5zzr2y37KvA4/hjUQWNNl/TALbyls8keIo8EW6HeM/Pt3DspP9x437zV/oP65W+71UOgW+SLfp/uOW4EwzW4Q3rm27c651v22W+I+/K3PZRIqmwBfpVuM/7v+5uMJ/zGmjN7O3AacBaeDO8hZNpHgKfJFumXb7080sAmBmp+GNMfqfQK2ZHefPH4E3hnIEb5zbTYNQXpGCqB++iM/M3gX80Z98Fq9pZwnwCPBl4BngTeBR4B14P9g+BbzLObfvYJdXpFCq4Yv4nHMPAR8D1gPz8H7EvRn4oHPuWeD/4F1gdS7Qitdn/50KewkL1fBFRKqEavgiIlVCgS8iUiUU+CIiVUKBLyJSJRT4IiJVQoEvIlIlFPgiIlVCgS8iUiUU+CIiVUKBLyJSJf4/WP+nv5xDhnQAAAAASUVORK5CYII=\n", "text/plain": [ - "
" + "" + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAEQCAYAAACz0c/rAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAABGRElEQVR4nO3dd3xUVdrA8d9JnfRKQiCk0HsJSFEEREQEleKioq4NF/TVdX3XLbLW9V131XV1XXWxICtid0VFRSxUQVoCgZAQSAglgUAaCQnpmfP+cScYMCFtJjOTeb6fz3xm5pZzn5vyzJ1zT1Faa4QQQnR+bvYOQAghRMeQhC+EEC5CEr4QQrgISfhCCOEiJOELIYSL8LB3AE0JDw/XcXFx9g5DCCGcSlJSUoHWuktj6xw24cfFxZGYmGjvMIQQwqkopY40tc4qVTpKqaVKqTyl1N4m1iul1L+UUplKqT1KqQRrHFcIIUTLWasO/y1g2gXWXwX0sTwWAIutdFwhhBAtZJUqHa31RqVU3AU2mQm8rY1uvVuVUsFKqSitda41ji9ER8gtqeBfazLZmlWIr5c7Vw/txh2XxGHydLd3aEK0SEfV4XcHshu8z7EsOyfhK6UWYHwDICYmpoNCE6J5SUeKuP0/O6iuNTOxbxeKy2t4ZnU6q1NPsPS2UYT5e9s7RJdSU1NDTk4OlZWV9g7FbkwmE9HR0Xh6erZ4H4e6aau1fh14HWDUqFEyyI9wCFn5Zdy2dAddArx5646LiA3zA2D13hP85oNdLFyexHu/GouXh7Ry7ig5OTkEBAQQFxeHUsre4XQ4rTWFhYXk5OQQHx/f4v066i/0GNCjwftoyzIhHFp1rZl739uFp7vi3bvGnE32ANMGd+W5ucNIPHKKF74/YMcoXU9lZSVhYWEumewBlFKEhYW1+htORyX8lcCtltY6Y4ESqb8XzuDtLYfZl3uaZ64bSrdgn5+tv2ZYN+aOjOaNjVnsP1Fqhwhdl6sm+3ptOX9rNct8H9gC9FNK5Sil5iul7lZK3W3ZZBWQBWQCbwD/Y43jCmFLxeXVvLgmg4l9u3DFwMgmt1s0fQB+3h78/Zv9HRidEK1nrVY685pZr4F7rXEsITrK8i1HKK2sZdH0/he8mgr182L++Hie/+4Ae4+VMLh7UAdGKUTLyV0mIRpRWVPHsi2HmdSvC/27Bja7/e2XxBFg8uC1jVkdEJ0QbSMJX4hGfJ58jIKyahZM6Nmi7QNNnswd2YPVe3PJL62ycXTCWfz6178mNjb27PuKigomTpxIXV1di/avrq5mwoQJ1NbWWiUeSfhCNOLjxBx6dfFjXM+wFu9z05gYauo0HydlN7+x6PQOHz7MunXrqK6uprTUuKG/dOlS5syZg7t7yzrreXl5cfnll/Phhx9aJSaHaocvhCM4XHCGxCOn+MO0fq1qCdE7wp8x8aF8tCObeyb2cvlWJB3lz1+kknb8tFXLHNgtkMevGXTBbQ4ePMjYsWPx8/MjODiYo0ePEhISwq5duwgMDOTxxx/nkUce4Y033iA1NZWxY8fy7rvv8t57750tY+XKlSxbtoxPPvnk7LLFixeTlpbGSy+9BMCsWbNYtGgRN998c7vPS67whTjPip05KAWzR3Rv9b6zRnTncGE5qVZOQMLx9OrVi/Hjx7N8+XKSk5MZOnQon332GYGBgaSmprJ3715uuOEGBgwYwN69e6muriYrK4uGw74//PDD/PnPf/5Zufv27Tv7fvDgwezYscMqMcsVvhANaK35Yk8uF/cKIyro5+3umzNtUFce/WwvX+w5Lq11OkhzV+K2lJqayuDBgwHYt28f/fr1A+CRRx7hySefRCnFgAEDSE1NpaCggODg4LP77t69G7PZzODBgzly5AirVq3innvuoaam5pxvh+7u7nh5eVFaWkpAQEC74pWEL0QDB/PLOFRwhjsviWvT/iF+XozvE86Xu3N5aNqFm3MK51ZRUUFlZSUhISFkZ2cTHh6Ol5cX27ZtY/Xq1ezatYt7772XyspKhgwZgo+Pzzk9Y5OTkxk5ciQA3333HRkZGQCkpaUxbNiwc45VVVWFyWRqd8xSpSNEA9+mnQRgygU6WjVn+uAojhVXsC9Xet52ZmlpaQwYMAAwru7rX//pT3/iiy++4PDhwxw+fJjdu3eTmppKSEgIdXV1Z5O+2WymrKyMuro6VqxYQWlpKRUVFbz11lvcdNNNZ49TWFhIeHh4qwZJa4okfCEa+C7tJEOjg9pUnVNvUj9jdrl1+/OsFZZwQA2rc3x8fNi5cycvv/wy1dXVTJky5ex2kZGRlJWVUVRUxNSpU9m0aRMA06dPJysri+HDh3P33XeTmprKqFGjWLBgAQkJP80RtW7dOmbMmGGVmJXRCdbxjBo1SssUh6Ij5ZVWMvqpNTx4RV9+fXmfdpV19Us/4OPpzsd3X2yl6ERDDa+oncnOnTt54YUXWL58eYv3mTNnDk8//TR9+/b92brGfg5KqSSt9ajGypIrfCEsfjhQAMDkARHtLmtyvwiSjpyiuLy63WWJziMhIYHLLrusVR2vZs2a1WiybwtJ+EJYbD5YQKifFwNaMJRCcyb1j8CsYcOBfCtEJjqTO++8s1Udr2699VarHVsSvhAYzTE3ZxYwrlcYbm7tb1kzLDqYQJMHWw4WWiE6IaxDEr4QwMH8M5w8XcUlvcKtUp67m2J0fBhbsiTh24qj3n/sKG05f0n4QgA/HjTq7y/p3fKxc5ozrlcYRwrLOV5cYbUyhcFkMlFYWOiySb9+isPWts2XjldCAJszC+ge7ENMqK/VyqwfeG1rViFzEqKtVq6A6OhocnJyyM933Xsk9ZOYt4YkfOHyzGbNtkNFXDEg0qo9Y/t3DSDY11MSvg14enq2avJuYZAqHeHysgrOUFxew0VxoVYt181NMSY+VOrxhcOQhC9c3s4jpwBIiA2xetlj4sPILqrgREll8xsLYWOS8IXLSzpyimBfT3qG+1m97JGWD5GdR09ZvWwhWksSvnB5SUdPkRATYpX29+cb2C0Qk6cbSUck4Qv7k4QvXFpxeTWZeWVnr8StzdPdjaHdgyXhC4cgCV+4tPqqFlslfDDuDaQeL6GypmXjpwhhK5LwhUtLOnIKdzfFsOhgmx0jISaYmjrN3mMlNjuGEC0hCV+4tKQjpxjULRAfr5YNZtUWCXLjVjgISfjCZdXUmdmdXUJCjO2qcwDC/b2JDfOVenxhd5LwhctKzy2loqaOUXG2TfgAI2NCSDpS7LJjvwjHIAlfuKykI0WAbW/Y1kuIDaGgrIrsIhlITdiPJHzhspKOFhMVZGrX/LUtJR2whCOQhC9c1s4jp2wynEJj+kYG4O/tIfX4wq4k4QuXdPJ0JceKK2x+w7aeu5tiRIx0wBL2JQlfuKT6AdM6ov6+XkJMCOknTlNWVdthxxSiIUn4wiUlHTmFt4cbA6PaP2F5S42MDcGsYXd2cYcdU4iGJOELl7Tz6CmGdA/Cy6Pj/gWGxwSjFFKtI+xGEr5wOVW1dew9drrDbtjWCzR50i8yQBK+sBtJ+MLlpB4/TXWducNu2DaUEBvCzqOnMJulA5boeFZJ+EqpaUqp/UqpTKXUQ42sv10pla+USrY87rLGcYVoi59muAru8GOPjAmhtLKWzPyyDj+2EO1O+Eopd+AV4CpgIDBPKTWwkU0/1FoPtzyWtPe4QrRV0pFTRIf4EBFg6vBj17cKkmodYQ/WuMIfDWRqrbO01tXAB8BMK5QrhNWZzZpth4oYHW/dCctbKjbMlzA/L0n4wi6skfC7A9kN3udYlp3vOqXUHqXUf5VSPRorSCm1QCmVqJRKzM/Pt0JoQpwrI6+MojPVjOsZZpfjK6VIiA0h8XCRXY4vXFtH3bT9AojTWg8FvgOWNbaR1vp1rfUorfWoLl26dFBowpVsOVgAwFg7JXyAMfGhHC4s51ixDKQmOpY1Ev4xoOEVe7Rl2Vla60KtdZXl7RJgpBWOK0SrbckqJDrEhx6hvnaLYUJf42JmU4Z8ixUdyxoJfwfQRykVr5TyAm4EVjbcQCkV1eDttcA+KxxXiFapr7+3V3VOvT4R/kQGevNDRoFd4xCux6O9BWita5VS9wHfAO7AUq11qlLqSSBRa70SuF8pdS1QCxQBt7f3uEK0VnJOMcXlNYzvE27XOJRSjO/dhbXpJzGbNW5uyq7xCNfR7oQPoLVeBaw6b9ljDV4vAhZZ41hCtNWafSdxd1NM6hth71C4tE84n+zMYe/xEobacAJ1IRqSnrbCZazZl8eo2BCCfD3tHQqX9gnHTcF3aSftHYpwIZLwhUvILion/UQpUwZE2jsUAML8vRkTH8ZXKbkyz63oMJLwhUtYufs4AFcO6mrnSH4yfWgUWflnOHBShlkQHUMSvuj0tNas2JnDRXEhxITZrznm+aYN6oqbgpW7jzW/sRBWIAlfdHp7cko4mH+GOQnR9g7lHF0CvJnUL4KPEnOoqTPbOxzhAiThi07vrR8P4+flzvQhUc1v3MFuGRtDfmkV36bKzVthe5LwRaeWW1LBF7uPc/1FPQjysX/rnPNN7BtBdIgPb27Kkpu3wuYk4YtObfH6g2jgzkvi7R1Ko9zdFPdM6sXOo8WsPyBDLQjbkoQvOq3MvDLe3XaUeaN72HXsnObMHdmDHqE+/H31fmqlLl/YkCR80SmZzZrHV+7Fx9OdB6b0tXc4F+Tl4cYfp/UnLfc0/9l82N7hiE5MEr7olJZtOczmzEIWTe9PuL+3vcNp1owhUUwZEME/vttPxslSe4cjOilJ+KLT2X6oiL+tSmdy/whuGh1j73BaRCnFU7OH4O/tyYLlSZyurLF3SKITkoQvOpWs/DIWLk8kOsSH568fhlLOMxJlZKCJf9+cQHZROQvfTqKyps7eIYlORhK+6DQyTpZyw+tbcVOKpbdfRLCvl71DarXR8aH8fe5Qth4q5O53kqiqlaQvrEcSvugU9uWe5sbXtwLwwYKxxIX72Tmitps9Ipq/zh7C+v35/M87OyXpC6uRhC+c3u7sYua9sRUvDzc+WjiOPpEB9g6p3eaNjuEvswazJj2PhculekdYhyR84dS2ZRVy85JtBJg8+HDBOOKd+Mr+fLeMjeXpOUPYcCCfX72dKElftJskfOG01u/P49al2+kaZOLjhRc71EiY1nLj6BieuW4omzILmL9sBxXVkvRF20nCF05pc2YBC95Ook+kPx8tHEfXIJO9Q7KZ60f14B9zh7HlYCF3vLWdM1W19g5JOClJ+MLpJGcX86u3E+nZxY93548l1M/5WuO01pyEaF64YTjbDxVxx392UCZJX7SBJHzhVLKLyrnjP9sJ9/fm7TtHO8T8tB1l5vDu/GveCJKOnmLB24kyhr5oNUn4wmlU1tRxz7tJ1Jk1b985mojAzluN05Srh3bjmeuG8uPBQp76ap+9wxFOxsPeAQjRElprHv1sL3uPnWbp7aOcup19e/1iZDT7ck/z5qZDDIwK5PqLetg7JOEk5ApfOIUPdmTzcVIO91/eh8n9I+0djt0tuqo/l/YJ55HP9pKcXWzvcISTkIQvHN7u7GIe/zyViX278JvL+9g7HIfg4e7GS/NGEBnkzT3vJFFQVmXvkIQTkIQvHFrRmWrueSeJiEBvXrxxOO5uzjMYmq0F+3rx6i0jKTpTzX3v7ZTJU0SzJOELh1VbZ+b+93dRcKaaV28Z6ZSDodnaoG5BPH3dELZmFfH01+n2Dkc4OLlpKxzW37/Zz6bMAp79xVAGdw+ydzgOa/aIaHZnl7Bk0yGG9gjm2mHd7B2ScFByhS8c0pd7jvPaxixuGRvD9aOkFUpz/jR9ABfFhfDH/+4h/cRpe4cjHJQkfOFwkrOL+f3HexgZG8JjVw+ydzhOwcvDjVduTiDA5MHC5UmUlMuMWeLnJOELh5KZV2b0pA3wYvEtCXh5yJ9oS0UEmFh8y0iOF1dw73s7ZXRN8TPy3yQcRsbJUm5Zsg03pVh+5xgiAlyvJ217jYwN4W9zhrL5YAELZBx9cR5J+MIhJB4u4vrXtlBr1rxz1xiX7knbXr8YGc0zc4byQ0Y+Ny/ZRn6ptNEXBkn4wq7qzJpXNxzkhte3EuTjySf3jGNAVKC9w3J611/Ug5fnJZB2/DTXvryJzZkF9g5JOABJ+MJuth8q4tqXN/H01+lMHRjJyl+PJzZMruytZcbQKP57zzhMnu7cvGQbD360m5xT5fYOS9iR0lrbO4ZGjRo1SicmJto7DGFlZVW1fJ92kre3HGbn0WK6Bpp45OoBzBgShVLSi9YWKmvqeHFNBm/+cAiz1lw7rBu/GBXN2Pgw3KTncqejlErSWo9qdJ01Er5SahrwIuAOLNFaP33eem/gbWAkUAjcoLU+fKEyJeE7P601BWXVZJwsZefRU+w4fIotWYVU15rpEerDXeN7MndUNL5e0v+vI+SWVPDq+oOs2HmM0qpawv29ubhXGON6hTEgKpC+kf7yu+gEbJrwlVLuwAHgCiAH2AHM01qnNdjmf4ChWuu7lVI3ArO11jdcqNy2JnyzWVNuaZlw/rk1fHfOKt3wpW50m3P31U0sb7ycJl42uX1Ljnv+dk0tt9b5gKayxkxVrZmq2jqqaoznyhozlTV1lFTUUFRezakz1RSdqebE6UoOF5SfMzNTnwh/xvcJZ8aQKBJiQuTq0k4qquv4Nu0E69Lz2JRZeHbgNaUgKtBEZJCJroEmIgNNBPl4EmDywN/bA3+TB37eHni5u+Hp7oanu7I8//RaKVBKoSzluVleo0ChjPX8tI2bZUH9ckfgCN803ZXCx8u9TfteKOFb4+N8NJCptc6yHOwDYCaQ1mCbmcATltf/BV5WSiltg/qkU+XVjPzL99YuVrSAl7sboX5ehPh50SXAm5ExIcSF+9Gziz/Do4NdanYqR+bj5c7M4d2ZObw7WmuOFpWTfqKU/SdKOVJYzonTFRw4WcqmjAJKZSpFuxjeI5jP7r3E6uVaI+F3B7IbvM8BxjS1jda6VilVAoQB5zQdUEotABYAxMTEtCkYXy8PHp4+oEGZTW/b8JNcnbOcJparZrehRWU2X05T259/GdRkfE2V1YI4zim/wQpvDzdMnu54e7id89rk6U6Qjye+Xu4OcXUkWk4pRWyYH7Fhflw5qOvP1pvNmjPVtZRW1lJWZTxqas3UmjXVdWZqas3U1GlqzWaqa81obXyrNJ6Nb4lmbfmeaXnW2vhW2XC9o3CUULoEeNukXIeqsNNavw68DkaVTlvK8PFy51cTelo1LiFclZubIsDkSYBJvp11BtZolnkMaDi6VbRlWaPbKKU8gCCMm7dCCCE6iDUS/g6gj1IqXinlBdwIrDxvm5XAbZbXvwDW2qL+XgghRNOs1SxzOvBPjGaZS7XWTymlngQStdYrlVImYDkwAigCbqy/yXuBMvOBI+0OzrbCOe8+hAtx5XMH1z5/Vz53cPzzj9Vad2lshcN2vHIGSqnEppo/dXaufO7g2ufvyucOzn3+MrSCEEK4CEn4QgjhIiTht8/r9g7Ajlz53MG1z9+Vzx2c+PylDl8IIVyEXOELIYSLkIQvhBAuQhJ+Gyil5iqlUpVSZqXUqPPWLVJKZSql9iulrrRXjB1BKfWEUuqYUirZ8phu75hsTSk1zfK7zVRKPWTveDqaUuqwUirF8vvu1OOXK6WWKqXylFJ7GywLVUp9p5TKsDyH2DPG1pKE3zZ7gTnAxoYLlVIDMXoaDwKmAf+2DB/dmb2gtR5ueayydzC2ZPldvgJcBQwE5ll+567mMsvv2ynborfCWxj/xw09BKzRWvcB1ljeOw1J+G2gtd6ntd7fyKqZwAda6yqt9SEgE2P4aNE5nB0KXGtdDdQPBS46Ia31RoyRARqaCSyzvF4GzOrImNpLEr51NTZUdHc7xdJR7lNK7bF8/XWqr7dt4Iq/3/Np4FulVJJlOHNXE6m1zrW8PgFE2jOY1nKo4ZEdiVLqe+DnA4TDw1rrzzs6Hnu50M8BWAz8H0YS+D/gH8CdHRedsIPxWutjSqkI4DulVLrlStjlaK21Usqp2rVLwm+C1npKG3ZryVDRTqWlPwel1BvAlzYOx9463e+3tbTWxyzPeUqpTzGquVwp4Z9USkVprXOVUlFAnr0Dag2p0rGulcCNSilvpVQ80AfYbueYbMbyB19vNsbN7M6sJUOBd1pKKT+lVED9a2Aqnf93fr6GQ73fBjjVt325wm8DpdRs4CWgC/CVUipZa32l1jpVKfURxny+tcC9Wus6e8ZqY88qpYZjVOkcBhbaNRobs0zPeR/wDT8NBZ5q57A6UiTwqWUaSw/gPa31avuGZDtKqfeBSUC4UioHeBx4GvhIKTUfY/j26+0XYevJ0ApCCOEipEpHCCFchCR8IYRwEZLwhRDCRTjsTdvw8HAdFxdn7zCEEMKpJCUlFTQ1p63DJvy4uDgSEzv12ExCCGF1SqkjTa2TKh0hhHARDnuFL4Rov/zSKjLzytBa0y3Yh9gwXyzt6IULkoQvRCeUdKSIv3+zn61Z5w72GBnozawR3bllTCw9Qn3tFJ2wF0n4QnQiWmteWpvJ898dIDLQm99e0ZcRMcF4uruRlX+Gtel5LPnhEG/+cIibx8Tw68v7EO7vbe+wRQeRhC9EJ/LsN/tZvP4gs0d056nZg/H1+ulffGzPMG4aE0NuSQUvr83knW1H+Xz3cR6dMZA5Cd2lqscFyE1bITqJj3Zks3j9QeaNjuH564edk+wbigry4anZQ/jmgUvp3cWfBz/ezW3/2cGJksoOjlh0NEn4QnQCRwrP8MQXqVzcK4y/zBrcoqv13hEBfLRwHE/OHMSOQ0VMe3Ejq/ee6IBohb1IwhfCyWmteeiTFDzcFP+4fhjubi2vmnFzU9w6Lo6v7h9PjxBf7n4niUUrUiivrrVhxMJeJOEL4eTWpuexJauQB6f2IyrIp01l9Ozizyf3XMzdE3vxwY6jXP3SJvYeK7FypMLeJOEL4cTqzJq/fZ1Oz3A/bhoT066yvDzceOiq/rw7fwxnqmqZ/e/N/Ht9JnVmGUK9s5CEL4QT+zb1BJl5ZTw4tR+e7tb5d764dzirfzOBKwZG8uzq/dzw2haOFpZbpWxhX5LwhXBSWmte25hFbJgv0wY3Ns9824X4efHKTQm8cMMw9p8sZdqLG3l/+1FkwiTnJglfuJT80ioe/Wwv0/65kVuWbOObVOdtlZJ45BTJ2cXcNT6+VTdqW0opxewR0XzzwASG9whm0YoU7lqWKM03nZgkfOEyDhWc4ZqXNvHhjmy6BpnIPlXOwuVJPPb5XsxOWE/9ztYjBJo8+MXIHjY9TrdgH96ZP4bHrh7IpswCpjy/gf9sPiR1+05IetoKl1BaWcP8ZTuoqq3j03svZlC3IGrrzDz7zX5e35iFt4cbD88YaO8wW6ykvIav957gxot64OPlbvPjubkp7hwfz+UDInj081T+/EUaK3Ye46+zhzAkOsjmxxfWIVf4wiW88F0GhwrOsPiWkQzqZiQoD3c3Fl3Vn9vGxfLGD4dYlZJr5yhbbuWe41TXmrl+lG2v7s8XG+bHsjsu4qV5IzhxupJrX9nEbz9K5nhxRYfGIdpGEr7o9PafKGXZlsPMGx3D2J5h56xTSvHI1QMZFh3EohUpTlM//XFiNgOiAhnULbDDj62U4pph3fj+txNZcGlPvtyTy2XPreeZ1emUVNR0eDyi5SThi07vpbUZ+Hi684cr+zW63tPdjX/eOILKmjr+76u0Do6u9Q7ml7Enp4Tr7DzgWZCPJ4umD2DtgxO5anBXFq8/yPin1/Ls6nQKyqrsFpdomiR80akdKjjDqpRcbhkbS7CvV5PbxYf78T+TevPVnlw2ZxZ0YIStVz/ezYyhUXaOxBAd4ss/bxzBV/ePZ0LfLizecJDxz6zliZWpZOWX2Ts80YAkfNGpvbkpCw93N+4cH9fstgsn9iQm1JcnVqY6dAuUr/bkkhAT3OZhFGxlULcgXrk5ge9/O5Grh3bjna1HmPyPDdy8ZCtfp+RSU2e2d4guTxK+6LQqquv4fNdxrh4SRUSAqdntTZ7uLLqqPxl5ZXy661gHRNh6hwvOkJZ7mulDHOPqvjG9uvjz3Nxh/LhoMr+b2pfDBeXc8+5Oxv1tDU+sTCXpyCnpwGUn0ixTdFqrUnIprarlhota3pJl2uCuDOkexAvfHeCaYVF4e9i+yWNrfG2pzrnKgRN+vYgAE/dN7sM9k3qzLj2PFbtyeG/7Ud768TDRIT7MGBLF5P4RjIwNwcNKw0KIC5OELzqtD3dkEx/ux+j40Bbvo5Ti91f249al2/lgeza3XRxnuwDb4Ou9uQzrEUz3YMeqzrkQdzfFlIGRTBkYSWllDd+mnmTl7uMs3XyI1zZmEWjyYELfLlzWL4KxvcKc6tycjSR80SkdLSxn++Ei/jCtX6tbslzaJ5yxPUN5eV0mN1zUA5OnY1zl552uZE9OCb9vorWRMwgweXLdyGiuGxlNaWUNmzMLWLMvj3X78/lyj9EPIjrEhzHxYYyJD+Wi+FDiwnxl+kUrkYQvOqVVe43kcc3Qbq3eVynFby7vy7w3tvJxYja/HBdn5ejaZv2BfAAu6xdh50isI8DkybTBUUwbHIXZrNl34jTbDxWxLauIdfvz+GRnDgCBJg+GRAcxuHsQQ7sHMzQ6iOgQH/kQaANJ+KJT+joll6HRQfQI9W3T/mN7hpIQE8yrG7K4cXSM1YYebo/1+/OIDPRmQFSAvUOxOjc3xaBuQQzqFsQdl8SjtSYzr4zEI6dIOVZCSk4JSzcdoqbOuNnr7+1B7wh/+kb60ycigD6R/vSJDKBbkEk+CC5AEr7odHJOlbM7p4Q/Tuvf5jKUUtw3uTd3vpXI58nH+cXIaCtG2Ho1dWZ+yChgxpAol0hoSin6RAbQJzKAeZZlVbV1HDhRxp5jxew/UUrGyTLWpufxUWLO2f18vdyJCfUlJtSX2DBfYsL8jNehvnQP8XGID257koQvOp2vU4yWLNOHtG+M+Mv6RTAgKpB/r89k9ojuNhmCuKV2HjlFaWUtk/p1sVsM9ubt4c6Q6KCfDdZWdKaajJOlZOSVcTC/jKOF5RwqOMOGA/lU1f7U9t9NQWSgia5BJroF+dA1yERUkKnBsw+RAd6dusWQJHzR6azam8ugboHEhvm1qxylFPde1ov73tvFN6kn7Nr2fd3+fDzcFJf0DrdbDI4q1M+LMT3DGHPeOElmsyavtIqjReUcKTzD0aJyjhdXcuJ0BftOnGZteh4VNXXn7OOmjPLC/b0tDy/CGrw+uzzAizA/b7w8nOvDQRK+6FTyS6vYdbSYB6/oa5XyrhocRc/wAyxef5CrBne1W3XK+v15XBQXSoDJ0y7Hd0Zuboquliv4xprmaq05XVFL7ukKcksqyS2uJLekgoKyKvJLqykoq+JI0RkKSqt/9sFQz8/LnWBfL4J9PY2HjxdBvp4E+3gS4vvT67Pb+HgSYPLE5Olml78lSfiiU9lQ35Klv3Vasri7Ke66tCd/+jSFbYeKfjbaZkc4XlxB+olS/jS97fckxM8ppQjy9STI15P+XS886uiZqloKy6rJL6uioKyKwrJqCsuqOFVeQ3FFNSXlNRRX1JBecpqSihqKy2uovcDwHO5uCj8vdwJMnvh7e+Dn7Y6/yZMAy+veEf4smNDL2qcsCV90Luv25xER4G3VYYPnJHTnH9/uZ8kPWXZJ+Ov3d67mmM7Iz9sDP28PYsJa1upLa01ZVS3F5TVnPwCKK6o5VV5DWWUtZ6pqKauqpbSylrKqGs5U1VFSUcOxU+WUVdWSW1IpCV+IC6mtM7PxQL7Vq15Mnu7cMjaWF9dkkJlXRu8If6uV3RLr9+fRPdinw48r2k4pRYDJqL7p2ClqLsy57jgIcQE7jxZTWllrkyvhX46LxdvDjTc3HbJ62RdSVVvH5swCJvXr4hLNMYVtScIXnca6/XlGS5Y+1m/JEu7vzZyEaFbszOnQyT0SD5/iTHWdVOcIq7BKwldKTVNK7VdKZSqlHmpkvbdS6kPL+m1KqThrHFeIhtal5zEqLoRAG7VkmT8+nqpaM8u3HLFJ+Y1Zl56Hl7sbF/fu+HsHovNpd8JXSrkDrwBXAQOBeUqpgedtNh84pbXuDbwAPNPe4wrRUG6J0ZLFllfCvSP8mTIgguVbj1DZRDM9a1t/IJ8xPUPx9ZLbbaL9rPFXNBrI1FpnASilPgBmAg0nB50JPGF5/V/gZaWU0raYBaHyNHz1W8sbS52nUue9bmodP6274Haqme1aU0ZzMbVwnbun5eFleTR8bXl4eIF3EJjqH4Hg4U1ncLYli5WaYzaq8CCPhXzLjswt5Cx7n94jJsKgWcbP0gayi8rJzCtj3ugYm5Tf4cx1UFVqPKrLoKoMqk5D9Rmoq4a6GjDX/PS6rsFrbTYeaMtr3eC9Pu+9ufFtzmok7fwsFekLvv35+ub2b+X60Hi47E8/C7O9rJHwuwPZDd7nAGOa2kZrXauUKgHCgHMmD1VKLQAWAMTEtPGP3FwLx5Ia/AD1ua8bPDW+Tp/3upXb/WydtY7VgjLawsMEvmEQ2A0Cu0NQNATHQMQAiBgEfs5RlbAuPY9uQSb62KIlS+Vp+GYR7HqXGDQ+nuF45SRDzifw7SNw+WNw0V3nfohbQf3omE4znILZDEVZUJgBBRlQmAmnj0PZSSjLgzP5oNv6zchyAaXcLK/dmnjfyLL61+cU19jvSjWzTWvXt6P82spG4ms/h/qeqLV+HXgdYNSoUW3LYr6hcP8ua4blHMx1lqshyxVRbVWDKyXL8tpKI3lVnYbKEqgsNt6fKYDTOXAiBQ6sPvePLSAKeoyB+EshfiKE97HbKTalviXLzBHdrd+SpTgb3p4Jpw7BuHth3L1sOQT3v7+TD6/2ZkzWK7Dqd5C9DWb+2/gWZSXr0/OICfWlZ3j7hoiwGXMd5CTCwTWQvd240Ko6/dN63zDjAiKgK0QNA/9I4//Tyx+8A4yHlz94+RnfNOu/lbqd923VzQPcpH2JNVgj4R+Dc5qaRluWNbZNjlLKAwgCCq1wbFHPzR3cfMCznbMFaQ2lJyAvDU6mGh8CRzZD2mfG+i4DYPB1MHQuhMS1N2qrsFlLltO58NYMqCiG276EuEsAmD7YzDPBvjy318THCz+FTc/DmieN6oM5S6ySnCpr6vjxYCFzR0U7VnNMrSFnB+x8G/avgvJC4+o5chAM+QV0S4Au/SGsl5HchUOxRsLfAfRRSsVjJPYbgZvO22YlcBuwBfgFsNYm9fei/ZSCwCjj0ftyY5nWxlf1g2th7wpY9xdY9xQMuBrG/Rpizq/B61hnW7L0smL1U201fHSr8e3n9i+g+8izqzzc3Zg/Pp4nv0xjZ3YxCZc+CModvn8cwnpbpe51x+EiKmrqHKc6p64W9v4XNr9oXAx4+kH/6dB3GvSeAj7B9o5QtEC7E76lTv4+4BvAHViqtU5VSj0JJGqtVwJvAsuVUplAEcaHgnAWShlXbGG9YPSvoCQHEpcaj31fQN+r4IonoYt1BixrrfUH8hkdH4qftxVrKNc+CTnbYe5b5yT7ejdc1IMX12Tw+oYsXv3lSLjkN5C/Hzb+HeInQNz4dh1+XXo+Xh5ujOtp59ExtYb0L+G7x6HoIEQOhmteNL7leXe+iVg6O6tUjGmtV2mt+2qte2mtn7Ise8yS7NFaV2qt52qte2utR9e36BFOKijauFH5v2kw5QmjyuffY2HtU8aVcQfKOWW0ZLHqlfDxZNjyCiTcBoNmN7qJn7cHvxwbyzdpJ8jKLzM+FKf/3ajm+vxeqGnfTbf1B/IY2zMMHy87zqdbnA3v3wgf3mLUsd/wDiz8AUbeLsneSTnUTVtrKCmvYf6yHXi4K/y9PYyHyQN/b08CTB4E+3oS5PPzR4DJ064TXDglL18Y/78w4pdGa5WNzxo3fee+ZXwb6AD1zTGtlvDNZvjyAfDrYnxruYDbLo7j9R+yWLLpEH+dPQS8/WHG87B8Fvz4Ekz8fZtCOFpYTlb+GW4ZE9um/a0ifRV8drdRlTP1KRhzN7h3unThcjrlb9DLw43aOs3x4krKqoyR6Uora6muMze5j1IQ4O1hDJd63odBoGVs6zDLxAhhlkkRwvy8MHna8QrMUfiFw+xXYcA18Pl98MZkI+n3uszmh16/P5/uwT706mKl5pj7VsLxXTDr1WbrpbsEeHNdQjT/Tcrhf6f0pUuAt3HOA2fCD/+AkbeBf+tvJH+37yQAk23Zp6ApWsPav8APzxkta+a+BaE9Oz4OYROdLuEH+Xry3q/GNrqussYYgrT+UT90af3j9Nnl1ZRU1HCipJKSilpOV9Q0+WHh7+1hfAD41c+MY8yEExnoTdcgH6KCTEQGmgjz88Kts3+D6D/DaK3x/jx45zrjQ2Do9TY7XFVtHT8eLGC2tZpj1tUaN6O79G9x3L+6NJ4Pdhzl7S2HeXBqP2Ph5MeMexs//gum/qXVYXybeoJ+kQHEdXRzzLpa+OI3kPyO8a1t+nPgaerYGIRNdbqEfyEmT3dMnu5EBrbuj1hrzZnqOorKqik489PkB4VnjFlxis5UU1hWTXZROcnZxRSdqabuvMkPPN0VEQE/zaHZNdBEVLDP2cmWe4T42re+1lpC4mD+t0bSX7HA6A+Q8EubHCrx8CnKq+uYZK3mmGmfQcEBuP5to5lrC/Ts4s/UgZG8veUId0/sZdw4Du8NQ+bC9iVw8W/Av+XVTYVlVew4XMR9l/Vu40m0kbkOPplv/AwmPgSTHrJ6RzJhfy6V8NtKqZ/uB7RkAoQ6s6awrIoTpyvJLankREklJ05bnksqST1+mu/3naSy5txvDREB3sSE+hqPMOODoE9EAL26+DvXh4F3ANz0EXx4M6z8tTGEw8CZVj/M+v1WbI6ptXGjNqw39L+mVbsumNCLb1JP8lFiNndcEm8snPB72PMhJL5pJM8WWrMvD7OGqYPaNwF7q2htXNmnfQZX/B9ccn/HHVt0KEn4NuDupogINBERaGJodOPbaK05VV7D0aJy42GZZPloUTlbswr5NPnY2VEUlIIeIb70jfSnT2QAfSL86dc1gL6RAXi6O2gPRC9fuPE9WHaNcaUfEAU9Rlv1EOv3W7E5Zs4OOL7TqMZoZcepkbEhjIoNYckPh7hlbKzxOwnvY7RPT1wK43/b4h6436adoHuwj1Vn7GrW+qdh13LjQ0qSfacmCd9OlFKE+nkR6ufF8B7BP1tfVVtHdlE5GSfLOHCyjIy8UjJOlrHhQD41dcYngZeHGwOjAhkWHcSQ6GCGRQfRs4u/47Q28vSBeR/Am1cYVTx3bzI6dFlBzqlyMvLKuOEiK80ntHWxMQjasHlt2v1/LuvFnW8lsmJnDjdcZBkHavRCeG+ucSN4yC+aLaOsqpaNGQXcPCam43rXpn8FG56GYTfBZQ93zDGF3UjCd1DeHu70jgigd0QAVw35aXlNnZkjhWdIyy0lJaeY3TklfJyUwzLLGO0BJg/GxIcytmcYY3uGMSAq0L4fAH7hRtJ/fRJ8chfc+rlVmvd9n2bFlizlRUbnolHzjaaVbXBZvwiGRQfx0tpM5iREG1f5vacYLVwSl7Yo4a/ee4LqWjNXD+3WphharfAgrFgI3UbA1S9Inb0LkITvZDzd3c5+EFw7zEgMdWZNVn4Zu3NKSDxcxNasQr7flwdAoMmDi3uFc/mACCb3jyDM3w7DIXfpZ7RP/+xuo62+FYYe+Cb1JH0i/OlpjeaYqSuMweWGnz8iSMsppXhgSl/ueGsHnyTlcOPoGKNqaPhNRjPHU0cg5MLt6j/dlUNsmC8JMcFtjqPF6mrh04VGjDe8I61xXIRTJfyamhpycnKorLTN0KHObmSoiZlDB+Dp6UluSQXbsorYcrCQDQfyWZ16AqUgISaEKQMiuWZYFNEhzd+Atprh8yBrPWx8DvpfDVFD21zUqTPVbD9cxN0TrdQ+PPl9Y8iAdsQERuevYT2CeXmdcZXv5eEGQ28wEv6ejy7YEetESSU/Hizk/sl9OqY6Z/M/jfsW171p9JwWLsGpEn5OTg4BAQHExcU51giCDkBrTWFhITk5OcTHxxMV5MOsEd2ZNaI7WuuzLYO+33eSZ1an88zqdEbHhzJ7RHemD4kiyMc20wKeY9rfjAHYPr8XfrWuzVU7a9PzqDNrpg60QkuWggw4lmj0Jm0n4yq/D3f8Zwef7MwxJi4JjoHY8bDnA5jwuyarTT6z3KSfNaJ7u+NoVl66caN20GxjTBzhMhy0iUfjKisrCQsLk2TfCKUUYWFhjX77UUoxuHsQD0zpy5e/vpQf/nAZD17Rl4LSKhatSGH0U9/z+493s/dYiW2D9A2FGc/BiT2w9ZU2F/NN6gm6BpoYGm2Fmab2fGQM7ztkbvvLAib17cKImGBe+O4AZ6pqjYXDbjQmAzm2s9F96syad7cdYXRcKPG27mylNXz9e6MV1fTnpN7exThVwgck2V9AS382PUJ9+fXlfVjz4ERW3ncJ142M5ss9uVz90iauW/wjX+3JxWy20ejVA2caQ+pu+LsxC1IrlVfXsjEjn6mDIq3zt5D+FcRcDAGR7S8L43fwyIyB5JVW8eqGg8bCAdcYk3jUzylwnnXpeWQXVXDbxXFWieGCUj+FQxth8qPGDXXhUpwu4QvrUUoxNDqYv84ewtY/Xc6jVw+koKyKe9/byZX/3Mjnycd+1mPYKqY+BbUVRt12K32XZnRYmz7ECs07Tx2GvFRjXHcrGhkbwszh3Xh9YxY5p8qNMXniJxrDLTQyDcTSzYfoGmhi6iDrfOg0qabCGOSu6xAYdadtjyUckiR8AUCQjyfzx8ez9sFJvDRvBErBbz5I5ornN/B1Si5Wna8mvDdc9Cujs8+JlFbt+tmuY3QLMjE6zgqzKe3/2njud1X7yzrPH6f1Ryl49LO9xs9uwDXGNIkn956z3dasQn48WMj88fG270S3/Q04fQymPdPioSNE5yIJX5zD3U1xzbBurP7NBF69JQFPdzfueXcnN7y2lT05xdY70MQ/gHegMYZ+CxWUVbExo4Brh3e3zkB0+1cZUzbaYDTIbsE+/OHK/qzbn89HidlGyySUcZVvYTZrnl2dTkSAN78cZ+OhkCtLjKkYe085O1WjcD2S8Nvgtddeo2vXrgwbNoxevXrx9ttvN7uPh4cHffv2Zfjw4QwfPpwnnzTGWq+oqGDixInU1dU1W0Z1dTUTJkygtra23efQHDc3xbTBUXx1/3j+NmcIWQVlXPvyZh78aDeFZVXtP4BvKIy7Dw58bQxH3AJf7cmlzqyZNcIKHZMqTsHhzTa5uq93+8VxjOsZxp+/SGNviRfEXnxOwn97y2F2Hi3mD9P6236Y7S2vGOc8+RHbHkc4NEn4bZCSksITTzzB7t27ef/99/ntb3/b7D7PPPMMgwcPZteuXSQnJ/PYY48BsHTpUubMmYO7e/P/8F5eXlx++eV8+OGH7T6HlvJwd2Pe6BjW/W4S90zqxcrdx5jy/AY+ScppfzXPmIVgCob1zzS7qdaa97cfZUBUIP27WmGcmYzvQdcZQzrbiJub4sUbhxPi68Wdb+3gRPcrjPlgCzLZlFHAU6v2MalfF65LsHFTzPIiI+EPnGn0qhUuy6na4Tf05y9SSTt+2qplDuwWyOPXDGp2uz179nDddUb75fj4eLy8jIGxDh06xAMPPMCxY8dwc3Nj+fLl9OtnjJH+4IMPsnPnTv7yl7/w6KOPni3r3Xff5b333jv7fuXKlSxbtoxPPvnk7LLFixeTlpbGSy+9xKxZs1i0aBE333yzVc65pQJMnvxxWn9mj+jOohUpPPjxblbsyuGvs4cQG9bGpoSmQLj4PuPm7fFdF0xGO4+eIv1EqTGzlDXs/wr8I6FbgnXKa0JEoIm37riIX765nXkbw1jnAV9+vIQHsifQs4sfL94wwvYtz5LeguoymPAH2x5HODy5wm+DlJQU+vXrh9aal19+maeeeoqamhruuusunn/+eRITE3niiSd4+umnz9lvyZIlLF++nB9++AEwqmiysrKIi4s7u83DDz/Mn//853P269WrF/v27QNg8ODB7Nixw7YneAF9IwP4eOE4/jJrMHuyS7jqxR/4YPvRtl/tj7Zc5W987oKbvb3lCAHeHswcboXqnNoq4wq/77RWj4zZFn0iA1h53yWMGDKUfTqO7ifWMndUD/57z8UE+dq4w1tdDexYYrQS6jrYtscSDs9pr/BbciVuC9nZ2ZSWljJ9+nSOHTvG0KFDeeKJJ/jvf/9Lamrq2Sv/2tpaLr300nP23bdvH3V1dQwbNgyAgoICgoODz67fvXs3ZrOZwYMHc+TIEVatWsU999xDTU3N2atAd3d3vLy8KC0tJSDAPhNJu7kpbhkby+T+Efzu4908tCKF7/fl8fR1Qwhv7Vg9pkC46C5jSsDCg43OhZtbUsGqlFxuHhNrnaGQD2+C6lLoZ93mmBcSEWji+RuGo9fdCBueYcTUSDB1QO/mfSuNljkznrf9sYTDkyv8VkpJSWHChAkkJydz4MAB0tPT2bJlC7t37+app54iOTmZ5ORk9u7dy+LFi8/uV1VVxe23387ixYsJDDTqoH18fM7pGZucnMzIkSMB+O6778jIyAAgLS3t7IdEfVkmk/0Hu+oW7MM788fwyIwBbMzI58oXNp4dxbJVRi8Ad0/Y+u9GV7+2IQutYf74+HZGbLF/FXj6Qs+J1imvFdSAq1Hon5qE2trWxUYrpD5TO+Z4wqFJwm+lPXv2MGKEUdccEhLCTTfdxFdffUVUVBTffPMNZrMxi1VKSso51RyPPfYYY8eOZerUn/7xQkJCqKurO5v0zWYzZWVl1NXVsWLFCkpLS6moqOCtt97ippuMkRwLCwsJDw/H07MDrg5bwM1NcdelPfnivvFEBJq46+1EFq3Y89OwAi0REGnMIbvrXeMGYwMnT1fy3vajzEnoTo9QKwz2pi3JttdkY7z+jhY5GIJijB6+tpa9wxggbczdHVJ1JRyf/BW0UkpKytmED3DNNdewatUq7rzzTsxmMwMGDGD48OE888wz59yMe+6559iwYcPZZpn332/MLDR16lQ2bdoEwPTp08nKymL48OHcfffdpKamMmrUKBYsWEBCgnFzcd26dcyYYbuWJW3Vr2sAn917MXdP7MUHO7KZ8a8f2HX0VMsLGHef0ft2x5vnLH5mdTpaa+611hyvubuNKo4OrM45h1JGy6Cs9VBVZttjbVts9HVox7DPopPRWjvkY+TIkfp8aWlpP1vm7JKSkvQtt9zS4u1nz56t9+/f3+R6R/gZbTlYoC/+2xrdc9FX+vlv9+ua2rqW7bj8Oq2f7aV1dYXWWuttWYU69o9f6me+3me94Nb+VesngrUuy7dema2VtVHrxwO1Tv3MdscoztH6z6Far/6T7Y4hHBKQqJvIq3KFb2cJCQlcdtllLe54NWvWLPr27dsBkbXd2J5hfP3ApVw7rBsvrsngule3cKjgTPM7jrsXzuTD3k8oKa/hfz9MJjrEh/smW+nqHozmmD3G2HfgsJhx4BMC6atsd4wdS0CbYfSvbHcM4XQk4TuAO++8s8Udr2699dYOiKj9Ak2evHDDcF6aN4JD+WVMf/EHlm85fOHB2HpOgoiBmLe8wn3vJXHydCUvzRuBr5eVGpMVZxtj99iwd22LuHsYTUIPrDaaTVpbdTkk/ceoOgqJs375wmlJwhc2dc2wbnzzvxMYGRvCo5+nMuNfP7Apo6DxdvtKUZ6wALe8VGqzNvLX2UMYERNivWDODpbmAPdA+s+AymI48qP1y97zoTGMwph7rF+2cGqS8IXNRQX5sHz+aP59cwJlVbXc8uY2ZvxrE//ZfIh9uacpOlPNwfwy3tx0iCu/j6BIB/D36M1cf1EP6waS/iWE9TFG67S3XpPBw2Q0EbUmrWHbq9B1qDF2jxANOG3HK+FclFJMHxLF5P4RrNh5jLe3HObPX6T9bLuRsSFUR91BdPJLTXbEapOKU0aHq0vut0557eXlBz0vM5pnTnvaejNPZa2D/HSY9arMZiV+RhK+6FAmT3duGhPDTWNiyC4qZ+fRUxSWVRPo48nwHsH0jvCH0p6wZ7FxpTr979Y58IFvLYOlXWOd8qyh/wxjtNATKe2eQP2srYvBLwIGz7FOeaJTkYQv7KZHqG/jnakCusKQXxgdsS572Jgxqr3Sv4SAKMcaLbLvNEAZV/nWSPgFGZDxLUz6E3i0cogL4RLaVYevlApVSn2nlMqwPDd6h00pVaeUSrY8VrbnmMJFjL0Has7AzubnGmhWTQVkrjE6WzlSj1P/LhAz1nq9bre9Cu5eMOoO65QnOp32/vU/BKzRWvcB1ljeN6ZCaz3c8ri2nccUriBqGMSOh+2vQ107J3zJ2mB8eNhw7Ps26z8DTqbAqSPtK6fiFCS/B0Pmgn+EdWITnU57E/5MYJnl9TJgVjvLcwqvvfYaCxcubPH2Tc12BY4745VDGPc/UJIN6V80v+2FpH4KpiCIu7T5bTta/RAP7W2ts3M51JQb4+YI0YT2JvxIrXWu5fUJILKJ7UxKqUSl1Fal1KymClNKLbBsl5ifn9/O0GwnJSWFoUNbXufa1GxX4PgzXtlV32lGx6Gti5vdtEnV5Ub9/cCZ4OFltdCsJqyXMa9ue6p16mqMb0Kx461381d0Ss3etFVKfQ90bWTVww3faK21UqqpbpSxWutjSqmewFqlVIrW+uD5G2mtXwdeBxg1atSFZ9T4+iGjdYM1dR0CVz3d7GZ79uzh+uuvZ+7cuURGRpKcnEx2djbvvvsur732Gtu2bePSSy/lzTeNgcCamu0Kzp3xqrnZrgC7zXhlF27uRueh1X+EnCSIHtn6Mg58bcz2NGSu9eOzlv4zYNMLxkihvqGt3z/tc+ObkLVaNIlOq9krfK31FK314EYenwMnlVJRAJbnvCbKOGZ5zgLWAw7UVKL19u7dy9ChQ0lJSaFnz55s2rSJhQsXMn/+fJ599lnS0tL46quvqKr6abLv82e7gp/PeNXcbFdg/xmvOtyIm40RH5sYK79Zez6GgG4Qe4l147Km/jOMJqPpX7Z+X63hx5eMDmV9rrR+bKJTaW+zzJXAbcDTlufPz9/A0nKnXGtdpZQKBy4Bnm3ncVt0JW4L2dnZ+Pn5YTKZKC4u5oEHHgCMjkXz588nKioK+Glmqnrnz3YF58541ZLZrhqWa88ZrzqUdwAk3Gq0QCl5EoJaMeF3WT5kfm9Mlu7WfJWZ3XQbAeF9jXr4hFaOlXRkM+Qmw9UvOFYLJOGQ2vsX8jRwhVIqA5hieY9SapRSaollmwFAolJqN7AOeFpr/fMulk6ivv4+NTWVhIQE3Cz/ZLt372bMmDEA5OTk0K1bt7OJurHZruDcGa9aOttVfXmOMONVhxm9wBj5cfvrrdtv13Iw17Q+iXY0pWDk7ZCzHU6mtm7fH18C3zAYNs8moYnOpV0JX2tdqLW+XGvdx1L1U2RZnqi1vsvy+ket9RCt9TDL85sXLtWx7dmzhyFDhpCSknJOIt6zZ8/ZG7m7d+8+56ZuY7NdwbkzXrVktitwvBmvOkRILPS/GpLeguoWDLMMYK6DxP8YLXO69LNpeFYxbJ7Rhj5pWfPb1ju+yxhxc8zd9pm9Szgd+Q7YSvVX+CkpKQwfPhyAyspKKioqCAkx+p01TP7Q9GxX8NOMVy2Z7Qocd8Yrmxt3rzG65O73W7b9gW+g5KgxQboz8A2FgbOMtvQVLZwpbN3fwBQsTTFFyzU1M4q9HzLjVeOcYcYrmzCbtX59stbP9de6qqz5bV+bpPULg7Wure6Y+KwhN8WYCWvd081ve3Sbse3G52wfl3AqyIxXjqszznhlE0rBlU9B6XHY9M8Lb5vxLRzfCZf+DtydqOqr62CjI9bWf0NlSdPbmetg1e+MsYFGt7wDoBCS8B1AZ5zxyiZixsLgX8DmFyEvvfFtaqvgmz9BSLxz3sic9BBUnYa1f2l6mx1LjMnYr3wKvP07Ljbh9CThC+dy5V+NppqfzDcGRTvf+r9BYSbMeM4xe9Y2J2qYcd9hx5LGZ8M6kQLfPgq9p8AgGQJZtI7TJXzd2NR4AnCRn01AJMxabDRf/PAWY+iEerveMXqsJtxqJERnNfkRCO0JH90KeT91uqMgA96da9zglQlORBs41Xj4JpOJwsJCwsLCzumMJIxkX1hY6Brt8/tOhWv/BSt/Da9eYlTz5KfDvpXGROjTn7N3hO1jCoIb34dlV8OSKTDiFmP5zuXgaYLbvjSGVhailZSjXhWOGjVKJyYmnrOspqaGnJycs52VxLlMJhPR0dGu00Y/awN8/7jRHt0nBEbNh4l/6DyTfxRnG/cj6idf7z/dqNIKirZvXMKhKaWStNajGl3nTAlfiEbVVoObR+cdWqCuBpSbYw8PIRzGhRK+U1XpCNEoZ7w52xrO1LRUOLROekkkhBDifJLwhRDCRThsHb5SKh9o50SfNhcOFNg7CDtx5XMH1z5/Vz53cPzzj9VaN9qMy2ETvjNQSiU2dXOks3PlcwfXPn9XPndw7vOXKh0hhHARkvCFEMJFSMJvn1ZOwdSpuPK5g2ufvyufOzjx+UsdvhBCuAi5whdCCBchCV8IIVyEJPw2UErNVUqlKqXMSqlR561bpJTKVErtV0pdaa8YO4JS6gml1DGlVLLlMd3eMdmaUmqa5XebqZR6yN7xdDSl1GGlVIrl992pB7tSSi1VSuUppfY2WBaqlPpOKZVheQ6xZ4ytJQm/bfYCc4CNDRcqpQYCNwKDgGnAv5VSnX3Eqxe01sMtj1X2DsaWLL/LV4CrgIHAPMvv3NVcZvl9O2Vb9FZ4C+P/uKGHgDVa6z7AGst7pyEJvw201vu01vsbWTUT+EBrXaW1PgRkAqM7NjphQ6OBTK11lta6GvgA43cuOiGt9Uag6LzFM4FlltfLgFkdGVN7ScK3ru5AdoP3OZZlndl9Sqk9lq+/TvX1tg1c8fd7Pg18q5RKUkotsHcwdhCptc61vD4BRNozmNaS4ZGboJT6HujayKqHtdafd3Q89nKhnwOwGPg/jCTwf8A/gDs7LjphB+O11seUUhHAd0qpdMuVsMvRWmullFO1a5eE3wStdVsmRT0G9GjwPtqyzGm19OeglHoD+NLG4dhbp/v9tpbW+pjlOU8p9SlGNZcrJfyTSqkorXWuUioKyLN3QK0hVTrWtRK4USnlrZSKB/oA2+0ck81Y/uDrzca4md2Z7QD6KKXilVJeGDfoV9o5pg6jlPJTSgXUvwam0vl/5+dbCdxmeX0b4FTf9uUKvw2UUrOBl4AuwFdKqWSt9ZVa61Sl1EdAGlAL3Ku1rrNnrDb2rFJqOEaVzmFgoV2jsTGtda1S6j7gG8AdWKq1TrVzWB0pEvhUKQVG7nhPa73aviHZjlLqfWASEK6UygEeB54GPlJKzccYvv16+0XYejK0ghBCuAip0hFCCBchCV8IIVyEJHwhhHARkvCFEMJFSMIXQggXIQlfCCFchCR8IYRwEZLwhRDCRUjCF0IIFyEJX4gLUEpFKKUKlFLVzU12opR6USmllVIvd1R8QrSGJHwhLuyfQBjwvNY6rZltd1ieJ9g0IiHaSMbSEaIJSqkxwFbgFBCvtS5pZvuLgc1ANeDbyQfOE05IrvCFaNpjlud/N5fsLfItz14Y3wqEcCiS8IVohFIqFmOycoA3zlsXo5Qar5Q6P6k3/H/q7JPXCyckCV+Ixk0HFJCqtT5y3rq/AT8A8ect72Z5rsXJZkISrkESvhCNG2l53tHIukstz/vPWz7C8pws9ffCEUnCF6JxsZbn7IYLlVIJGPPaVmitS8/bZ7rl+WsbxyZEm0jCF6Jxnpbn8/9H7rE8n1NHr5TqD0wGzMDbtg1NiLaRhC9E4+rr7S9XSrkBKKUuA+4APgW8lFKjLcv9gaUYdf5LtdaZdohXiGZJO3whGqGUmgyssbxNwqjamQ6sB/4EJAIngQ3AeIwbttuByVrrMx0drxAtIVf4QjRCa70WuBVIAwZj3MT9F3Cd1joJ+D1GB6uZQClGm/1JkuyFI5MrfCGEcBFyhS+EEC5CEr4QQrgISfhCCOEiJOELIYSLkIQvhBAuQhK+EEK4CEn4QgjhIiThCyGEi5CEL4QQLuL/AbZh0yoHOiPmAAAAAElFTkSuQmCC", + "text/plain": [ + "
" ] }, "metadata": { @@ -163,19 +173,22 @@ } ], "source": [ - "w = np.linspace(-6, 6, 2**12)\n", - "gloc = semi_circle_hiltrans(w + 1e-3j)\n", + "U = 2.0\n", + "D = 2.0\n", + "beta = 10**6\n", + "w = np.linspace(-6*D, 6*D, 2**12)\n", + "gloc = semi_circle_hiltrans(w + 1e-3j,D)\n", "#gloc= 0.0\n", - "gloc, sigma_loc = dmft_loop(gloc, w, 3.3, 100, 100)\n", - "\n", - "#plt.plot(w, -gloc.imag, lw=3, label=r'$\\pi A(\\omega)$')\n", - "plt.plot(w, sigma_loc.real, '-', label=r'$Re \\Sigma(\\omega)$')\n", - "plt.plot(w, sigma_loc.imag, '-', label=r'$Im \\Sigma(\\omega)$')\n", - "plt.xlabel(r'$\\omega$')\n", - "plt.legend(loc=0)\n", + "gloc, sigma_loc = dmft_loop(gloc, w, U, beta, D, 100)\n", + "fig,ax = plt.subplots(2,1)\n", + "ax[0].plot(w, -gloc.imag, label=r'$\\pi A(\\omega)$')\n", + "ax[0].legend(loc=0)\n", + "ax[1].plot(w, sigma_loc.real, '-', label=r'$Re \\Sigma(\\omega)$')\n", + "ax[1].plot(w, sigma_loc.imag, '-', label=r'$Im \\Sigma(\\omega)$')\n", + "ax[1].set_xlabel(r'$\\omega$')\n", + "ax[1].legend(loc=0)\n", "#plt.ylim([0, 0.1])\n", - "#plt.xlim([-1, 1])\n", - "plt.show()" + "#plt.xlim([-1, 1])\n" ] }, { @@ -296,7 +309,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.3" + "version": "3.9.10" } }, "nbformat": 4,