{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Vo4mY-6N2yoA"
   },
   "source": [
    "# ４章 線形システムの構造"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "id": "gJYRxJfd2yoC"
   },
   "outputs": [],
   "source": [
    "from control.matlab import *\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "#plt.rcParams['font.family'] ='sans-serif' #使用するフォント\n",
    "plt.rcParams['font.family'] = 'Times New Roman' # font familyの設定\n",
    "plt.rcParams['mathtext.fontset'] = 'cm' # math fontの設定\n",
    "plt.rcParams['xtick.direction'] = 'in' #x軸の目盛線が内向き('in')か外向き('out')か双方向か('inout')\n",
    "plt.rcParams['ytick.direction'] = 'in' #y軸の目盛線が内向き('in')か外向き('out')か双方向か('inout')\n",
    "plt.rcParams['xtick.major.width'] = 1.0 #x軸主目盛り線の線幅\n",
    "plt.rcParams['ytick.major.width'] = 1.0 #y軸主目盛り線の線幅\n",
    "plt.rcParams['font.size'] = 11 #フォントの大きさ\n",
    "plt.rcParams['axes.linewidth'] = 0.5 # 軸の線幅edge linewidth。囲みの太さ\n",
    "plt.rcParams['mathtext.default'] = 'it'#'regular'\n",
    "plt.rcParams['axes.xmargin'] = '0'\n",
    "plt.rcParams['axes.ymargin'] = '0.05'\n",
    "plt.rcParams['savefig.facecolor'] = 'None'\n",
    "plt.rcParams['savefig.edgecolor'] = 'None'\n",
    "\n",
    "plt.rcParams[\"legend.fancybox\"] = True     # 丸角\n",
    "# plt.rcParams[\"legend.framealpha\"] = 1    # 透明度の指定、0で塗りつぶしなし\n",
    "# plt.rcParams[\"legend.edgecolor\"] = 'gray' # edgeの色を変更\n",
    "plt.rcParams[\"legend.handlelength\"] = 1.8  # 凡例の線の長さを調節\n",
    "plt.rcParams[\"legend.labelspacing\"] = 0.4  # 垂直方向（縦）の距離の各凡例の距離\n",
    "plt.rcParams[\"legend.handletextpad\"] = 0.7 # 凡例の線と文字の距離の長さ\n",
    "plt.rcParams[\"legend.markerscale\"] = 1.0   # 点がある場合のmarker scale"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "id": "iOA70Kcj2yoE"
   },
   "outputs": [],
   "source": [
    "def linestyle_generator():\n",
    "    linestyle = ['-', '--', '-.', ':']\n",
    "    lineID = 0\n",
    "    while True:\n",
    "        yield linestyle[lineID]\n",
    "        lineID = (lineID + 1) % len(linestyle)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "id": "EEUYxua92yoE"
   },
   "outputs": [],
   "source": [
    "def plot_set(fig_ax, *args):\n",
    "    fig_ax.set_xlabel(args[0])\n",
    "    fig_ax.set_ylabel(args[1])\n",
    "    fig_ax.grid(ls=':', lw=0.5)\n",
    "    if len(args)==3:\n",
    "        fig_ax.legend(loc=args[2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "id": "FomeQ5aF2yoE"
   },
   "outputs": [],
   "source": [
    "def bodeplot_set(fig_ax, *args):\n",
    "    fig_ax[0].grid(which=\"both\", ls=':', lw=0.5)\n",
    "    fig_ax[0].set_ylabel('Gain [dB]')\n",
    "\n",
    "    fig_ax[1].grid(which=\"both\", ls=':', lw=0.5)\n",
    "    fig_ax[1].set_xlabel('$\\omega$ [rad/s]')\n",
    "    fig_ax[1].set_ylabel('Phase [deg]')\n",
    "    \n",
    "    if len(args) > 0:\n",
    "        fig_ax[1].legend(loc=args[0])\n",
    "    if len(args) > 1:\n",
    "        fig_ax[0].legend(loc=args[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 図を保存するかどうか\n",
    "is_savefig = False\n",
    "# 図の保存パス\n",
    "figpath=\"./notebook_output/\"  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 数式処理のためにsympyをインポート\n",
    "import sympy as sp\n",
    "from sympy.matrices import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.linalg import expm\n",
    "from numpy.linalg import inv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.1〜4.4 可制御性のチェック〜可制御正準形への等価変換"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "例4.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.matrix([[0,1,2],[1,3,4],[0,5,6]])\n",
    "B = np.matrix([[1],[0],[0]])\n",
    "C = np.matrix([1,0,0])\n",
    "D = np.matrix([0])\n",
    "sys = ss(A, B, C, D)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1 0 1]\n",
      " [0 1 3]\n",
      " [0 0 5]]\n"
     ]
    }
   ],
   "source": [
    "Vc=np.block([B,A*B,A*A*B])\n",
    "print(Vc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1. 0. 1.]\n",
      " [0. 1. 3.]\n",
      " [0. 0. 5.]]\n"
     ]
    }
   ],
   "source": [
    "#controlライブラリのctrb関数を用いても同じ\n",
    "Vc=ctrb(A,B)\n",
    "print(Vc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.linalg.matrix_rank(Vc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "#特性多項式とその係数を求める\n",
    "[_,a2,a1,a0] = np.poly(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-9.0, -2.999999999999999, -3.999999999999998)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a2,a1,a0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-3. -9.  1.]\n",
      " [-9.  1.  0.]\n",
      " [ 1.  0.  0.]]\n"
     ]
    }
   ],
   "source": [
    "T0 = np.matrix([[a1,a2,1],[a2,1,0],[1,0,0]])\n",
    "print(T0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-0.  -0.   0.2]\n",
      " [-0.   1.   1.2]\n",
      " [ 1.   9.  11.2]]\n"
     ]
    }
   ],
   "source": [
    "T=np.linalg.inv(Vc*T0)\n",
    "print(T)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 8.88178420e-16  1.00000000e+00  0.00000000e+00]\n",
      " [-7.10542736e-15  0.00000000e+00  1.00000000e+00]\n",
      " [ 4.00000000e+00  3.00000000e+00  9.00000000e+00]]\n"
     ]
    }
   ],
   "source": [
    "Abar=T*A*np.linalg.inv(T)\n",
    "Bbar=T*B\n",
    "print(Abar)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.]\n",
      " [0.]\n",
      " [1.]]\n"
     ]
    }
   ],
   "source": [
    "print(Bbar)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "可制御性グラム行列と有限時間到達（式(4.4),(4.5)）については章末問題【１】を参照"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.5 Kalmanの正準構造分解"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "import scipy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-14  13  -5  25]\n",
      " [ -3   2  -3   3]\n",
      " [  3  -3   1  -6]\n",
      " [ -9   9  -2  17]]\n",
      "[[ 3]\n",
      " [-1]\n",
      " [-2]\n",
      " [ 2]]\n",
      "[[ 7 -4  7 -7]]\n"
     ]
    }
   ],
   "source": [
    "A = np.matrix([[-14,13,-5,25],[-3,2,-3,3],[3,-3,1,-6],[-9,9,-2,17]])\n",
    "B = np.matrix([[3],[-1],[-2],[2]])\n",
    "C = np.matrix([7,-4,7,-7])\n",
    "D = np.matrix([0])\n",
    "print(A)\n",
    "print(B)\n",
    "print(C)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 3.  5.  3.  5.]\n",
      " [-1.  1. -1.  1.]\n",
      " [-2. -2. -2. -2.]\n",
      " [ 2.  2.  2.  2.]]\n"
     ]
    }
   ],
   "source": [
    "#可制御性のチェック\n",
    "Vc=ctrb(A,B)\n",
    "print(Vc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "可制御性行列のランクは2で，不可制御です．"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.linalg.matrix_rank(Vc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 7. -4.  7. -7.]\n",
      " [-2. -1. -2.  2.]\n",
      " [ 7. -4.  7. -7.]\n",
      " [-2. -1. -2.  2.]]\n"
     ]
    }
   ],
   "source": [
    "#可観測性のチェック\n",
    "Vo=obsv(A,C)\n",
    "print(Vo)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "可観測性行列のランクは2で，不可観測です．"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.linalg.matrix_rank(Vo)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "部分空間の構造解析"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 準備\n",
    "# 単位行列\n",
    "I=np.eye(A.shape[0])\n",
    "# 射影行列\n",
    "def proj(X):\n",
    "    return X@inv(X.T@X)@X.T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-0.82746882 -0.25552955]\n",
      " [-0.03490737 -0.8653216 ]\n",
      " [ 0.39628072 -0.30489603]\n",
      " [-0.39628072  0.30489603]]\n"
     ]
    }
   ],
   "source": [
    "Tc = scipy.linalg.orth(Vc)  # Image Vcの基底（正規直交）=可制御部分空間\n",
    "print(Tc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-8.16496581e-01  0.00000000e+00]\n",
      " [-1.11022302e-16  7.37035603e-17]\n",
      " [ 4.08248290e-01  7.07106781e-01]\n",
      " [-4.08248290e-01  7.07106781e-01]]\n"
     ]
    }
   ],
   "source": [
    "Tno = scipy.linalg.null_space(Vo) # Ker Voの基底＝不可観測部分空間\n",
    "print(Tno)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-8.16496581e-01]\n",
      " [-5.55111512e-17]\n",
      " [ 4.08248290e-01]\n",
      " [-4.08248290e-01]]\n"
     ]
    }
   ],
   "source": [
    "#TnoのうちTcと従属な成分=可制御・不可観測部分空間\n",
    "Tc_no=np.matrix(proj(Tc)@Tno)\n",
    "# 実はこの2本の基底は線形従属なので，１本目だけ取り出します．\n",
    "Tc_no = Tc_no[:,0]\n",
    "print(Tc_no)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.        ]\n",
      " [0.        ]\n",
      " [0.70710678]\n",
      " [0.70710678]]\n"
     ]
    }
   ],
   "source": [
    "#TnoのうちTcと独立な成分=不可制御・不可観測部分空間\n",
    "Tnc_no = np.matrix((I-proj(Tc))@Tno)\n",
    "# 実はこの2本の基底は線形従属なので，2本目だけ取り出します．\n",
    "Tnc_no = Tnc_no[:,1]\n",
    "print(Tnc_no)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-0.01163579]\n",
      " [-0.03490737]\n",
      " [-0.01163579]\n",
      " [ 0.01163579]]\n"
     ]
    }
   ],
   "source": [
    "#TcのうちTnoと独立な成分=可制御・可観測部分空間\n",
    "Tc_o = np.matrix((I-proj(Tno))@Tc)\n",
    "# 実はこの2本の基底は線形従属なので，１本目だけ取り出します．\n",
    "Tc_o = Tc_o[:,0]\n",
    "print(Tc_o)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.25 -0.25  0.25 -0.25]\n",
      " [-0.25  0.25 -0.25  0.25]\n",
      " [ 0.25 -0.25  0.25 -0.25]\n",
      " [-0.25  0.25 -0.25  0.25]]\n"
     ]
    }
   ],
   "source": [
    "#TcともTnoともnoと独立な成分=不可制御・可観測部分空間\n",
    "Tnc_o = np.matrix((I-proj(Tno))@(I-proj(Tc)))\n",
    "print(Tnc_o)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.25]\n",
      " [-0.25]\n",
      " [ 0.25]\n",
      " [-0.25]]\n"
     ]
    }
   ],
   "source": [
    "# 実はこの４本の基底も線形従属なので，１本目だけ取り出します．\n",
    "Tnc_o = Tnc_o[:,0]\n",
    "print(Tnc_o)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "以上を並べて変換行列を作ります．"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-8.16496581e-01  1.08779196e-15  4.08248290e-01 -4.08248290e-01]\n",
      " [-7.16181060e+00 -2.14854318e+01 -7.16181060e+00  7.16181060e+00]\n",
      " [ 7.85046229e-17 -1.57009246e-16  7.07106781e-01  7.07106781e-01]\n",
      " [ 1.00000000e+00 -1.00000000e+00  1.00000000e+00 -1.00000000e+00]]\n"
     ]
    }
   ],
   "source": [
    "T=inv(np.block([Tc_no,Tc_o,Tnc_no,Tnc_o]))\n",
    "print(T)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "式(4.24)(4.25)のようなKalmanの正準形に整形されます．"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.00000000e+00 -5.70035028e-02 -1.73205081e+01  1.67381799e+01]\n",
      " [ 0.00000000e+00 -1.00000000e+00  5.68434189e-14  7.16181060e+01]\n",
      " [ 1.33226763e-15  7.49400542e-16  5.00000000e+00 -4.24264069e+00]\n",
      " [-5.55111512e-16 -3.57353036e-16 -4.32986980e-15  1.00000000e+00]]\n",
      "[[-4.08248290e+00]\n",
      " [ 2.86472424e+01]\n",
      " [ 4.44089210e-16]\n",
      " [ 4.21884749e-15]]\n",
      "[[-1.33226763e-15 -1.04722122e-01  8.88178420e-16  6.25000000e+00]]\n"
     ]
    }
   ],
   "source": [
    "Abar=T*A*inv(T)\n",
    "Bbar=T*B\n",
    "Cbar=C*inv(T)\n",
    "print(Abar)\n",
    "print(Bbar)\n",
    "print(Cbar)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## 章末問題"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 【１】 可制御性グラム行列を用いて目標状態に到達させる\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.matrix([[-1,-1],[1,-1]])\n",
    "B = np.matrix([[0],[1]])\n",
    "C = np.matrix([1,0])\n",
    "D = np.matrix([0])\n",
    "sys = ss(A, B, C, D)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Vc=ctrb(A,B)\n",
    "np.linalg.matrix_rank(Vc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "dt=0.001\n",
    "tf=10.00\n",
    "trange=np.arange(0,tf,dt)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "式(4.5)の近似積分"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "matrix([[ 0.125     , -0.12499992],\n",
       "        [-0.12499992,  0.37550017]])"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gc=np.zeros(A.shape)\n",
    "for tau in trange:\n",
    "    Gc = Gc + expm(A*tau) *B *B.transpose() *expm(A.transpose()*tau)*dt\n",
    "Gc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "def reach(t0,x0,tf,xf,tau):\n",
    "    u = B.transpose() * expm(A.transpose()*(tf-tau)) * inv(Gc) * (xf-expm(A*tf)*x0 )\n",
    "    return u[0,0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "初期状態と目標状態を設定（自由に変えてみて下さい）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "x0=np.array([[1],[0]])\n",
    "xf=np.array([[-0.5],[-0.5]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARoAAAEICAYAAACArTsqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAnpUlEQVR4nO2de5QcVZ34P9+efsyjQ4YQISSBJLwlJCAiYcMjRBQUFFEQweiKIrBABMETg4TfISFkc4yuiissPlijCJ5FBQxLWESQEJAguoEACmRBIm+cxDwmM9MzPXN/f9yeZh49M3W7u7q+xdzPOXPS1dVV/anq9Lfrfuve7xVjDB6PxxMmiagFPB7POx8faDweT+j4QOPxeELHBxqPxxM6PtB4PJ7Q8YHG4/GEjg80Ho8ndCILNCJykIjcJSJzhlj/WxF5RUTWiUhdrf08Hk/1SEb1xsaYZ0VkByAD14nIacAiY8xjNRfzeDxVJ+qmU+cQz88FfiUit4vIHrUU8ng81SeyK5rhMMZcKiKXA1cCNwMnDnzNzJkzzVNPPVVcnj17Nscffzw9PT10d3cDkEqlMMaQz+eLywBdXV0AJJNJRKS4XFdXRyKR6LdcV1dHZ2dnyeVEIkEymaSrqwtjzKBlESGVSpHP5+np6Rm0DJBOp+nu7i46D1xOpVL+mPwxqTqm5cuX32SM+SIOqAw0AMaYbmCpiKwrtf5DH/oQGzZsqLHV8LS0tDB+/PioNfrhnYKh0Ql0ei1fvvwV122ibjoVEZFmEUkVHicK/9YDD0Uq5vF4KiayKxoRmQIcABwtIo8B1wJrROQBYK2I3Ac8BywttX1dnb4bUel0OmqFQXinYGh0Ar1erkR512kTMLvPU/P7PD54pO19oAmGdwqGRifQ6+WKmqaTK71JMU20trZGrTAI7xQMjU6g18uV2AYaj8cTH2IbaHzTKRjeKRganUCvlys+0FSRTCYTtcIgvFMwNDqBXi9XYhtofI4mGN4pGBqdQK+XK7ENNBrRWOjdOwVDoxPo9XIltoEmkdCn3tvNXBPeKRganUCvlyv6vq0BSSb1jZ7Q2J72TsHQ6AR6vVyJbaDpHXymiba2tqgVBuGdgqHRCfR6uRLbQKOx7do7KlYT3ikYGp1Ar5crsQ00GnM0Gptz3ikYGp1Ar5cr+r6tAdH4AdTX10etMAjvFAyNTqDXy5XYBhqfowmGdwqGRifQ6+VKbAONz9EEwzsFQ6MT6PVyJbaBRmRQTfPI0Zg38k7B0OgEer1cie1RaOzI1NjYGLXCILxTMDQ6gV4vV2IbaHqLPmuio6MjaoVBeKdgaHQCvV6uxDbQaGy7agx+3ikYGp1Ar5crsQ00PkcTDO8UDI1OoNfLldgehc/RBMM7BUOjE+j1ciW2gUbjJWUul4taYRDeKRganUCvlyuxDTQaczQaOxF6p2BodAK9Xq7ENtBoRGPeyDsFQ6MT6PVyJbaBRmPR5mw2G7XCILxTMDQ6gV4vV2IbaHonHteExva0dwqGRifQ6+WKDzRVRGPBdO8UDI1OoNfLldgGGo/HEx9iG2h8jiYY3ikYGp1Ar5crkQUaETlIRO4SkTkl1s0VkSsKfzNLbe+bTsHwTsHQ6AR6vVyJrEydMeZZEdkB9Lt/JyL1wArgSCAN3Ccic8yAAjQ+0ATDOwVDoxPo9XIl6qZTqbM4G2gxlhzQCOxdWy2Px1NN9BXehQnAlj7LHcBEYFPfF61du7ZfZ6YFCxawZMkS8vl88ZZgU1MTPT09tLe3A2+PG+ktj9jQ0EAikWDnzp2AnUMnmUwWl9PpNOl0ujgtaTqdJpPJ0NraijGGVCpFJpOhra2Nnp4euru76ejoKC4nEgkaGxvp6Oggn88Xl3O5HF1dXYgI2WyWXC5X/OXKZrN0dnYWl5uamkofU1sbdRs20Pjii5DP0zFxIl1HHEFDc3O/YwJob28v+5iSyST19fVVPaa6ujp27NgR2edU6piSySRbt26t/udU4THlcjlaWloi+ZyGOqayMMZE9gesBI4f8NwJwB19lp8E9h247RVXXGG0sX379tq80aOPGnPoocbst58x8+YZ88//bMysWcbssosx555rzLPP1t7JAe8UHI1ewGLj+F2PuulURESaRSQFPAxMKjyXAbqMMS8MfL3GHE1NOlf993/DqafClVfC88/Dz34GP/kJrFsHGzfCXnvBscfCZZdBn6sGTXin4Gj1ciXKu05TgAOAo0WkAbgWOM3YvMxiEVkIfAU4LypHdWzcCJ//PNx1F5x5JgwcB7P77nD11fDnP8PWrTBzJsn16yNR9Xj6EuVdp03YxG8v8/usWw2sHm57jfVompqawn2DL30JrrgCZs0a/nXjx8OPfwy3387YefNg2TI4T0+8Dv08lYFGJ9Dr5YrGZHAgzGibbuXxx+Evf7FXM0H5xCdo33dfGs88E15+GZYsGXwVFAEaS3xodAK9Xq6oydG4orHwVe8dhlC4/np7ReN4Jdc2aRKsXQurV8O//AsoyG2Fep7KRKMT6PVyJbaBZlTR2QmrVsHZZ5e3/e67w+9+B889BxdeCAqvBj3vbGIbaDTmaEKr7/rww7D//jBpkvOmRacxY2yza8MGuPzySIONxjq4Gp1Ar5crsQ00o4pHHoE5g4aEuTNmDNxzj726ueaayvfn8QQktoFGYy3V0CZkX7cO/umfytp0kNOuu8K998LKlXDLLZW7VcNJARqdQK+XK7G96zSq2LABDjusevvbYw/bjHr/+2HaNJg9e+RtPJ4KiO0VTTKpL0Y2NDRUf6ft7fD3v8Pe5Y0rHdLpkEPsVc3pp8NLL5WtVw6hnKcK0egEer1ciW2g0VgdPpRZBV98EaZOhbq6sjYf1unkk2HBAjjjDKjhHM8aZ1/U6AR6vVyJ7VFozNH0HTFdNV55xY5fKpMRnS67zAayyy8v+z1cCeU8VYhGJ9Dr5UpsA82ooaXFDikICxG46Sa47z649dbw3sczqoltoKkrsykRJplMpvo73by5okATyGnsWPjFL+DSS+HZZ8t+r6o61RiNTqDXy5XYBhqNbddQEtSbN8O4cWVvHtjpsMNg6VKYN8/2RA4RjYl8jU6g18sVfd/WgIyaHE1HB1TQO9TJ6YILYMIEG3BCRGPeQaMT6PVyJbaBZtTQ2Qm1mlqmN1/zwx/aToIeT5WIbaDRmKMJZa6pCgONs9OECXDDDfDZz0KhBm+10Tgnl0Yn0Ovlig80VSS0QFPBANKynD7xCTj6aFi4sOz3HQ6NXx6NTqDXy5XYBhqN8920hnEFkEhABcWPynb69rfhzjvtgM4qE8p5qhCNTqDXy5XYBppRQzoNUSS+d90VrrvOlgB9hxTI9kRHbAPNqGk6pdMV3W6uyOn00+GAA2D58vL3UQKNzQGNTqDXyxUfaKpIKJ2rUqmKAk1FTiK2hOj118Mzz5S/n2o6hYRGJ9Dr5UpsA82oydGk0xU1XSp2mjTJFsk6//yKckVVdQoBjU6g18uV2AYajYQyM8Muu1R0m7kqThdcAPk83Hxz5ftC5wwWGp1Ar5crsQ00GocghFLHuLnZTgZXJlVxSiTg3/8dvvY12Lat4t1prPes0Qn0ermi79saEI1jQEJpT1cYaKrmdOSR8OEPV6XWsMa8g0Yn0OvlSmwDjcaxTqHUd60w0FTVafly+OlP7ZS7FaCxDq5GJ9Dr5UpsA43GtmsoswqOHVtRoKmq0+67w1VXwSWXVDRdi8bZFzU6gV4vV2IbaDTmaEJpzu26K2zZUvbmVXe6+GJ4803ba7hMNDZ7NTqBXi9X9H1bA6LxA6ivr6/+TidMgDfeKPsKoupOySR885twxRVl91gO5TxViEYn0OvlSmSBRkQuFpHzRORqEWkusf63IvKKiKwTkUG980ZNjiabtYXJt28va/NQnE480c7K8KMflbW5xryDRifQ6+VKJIFGRA4HZhljfgjcD1w9YP1pwCJjzGRjzFHGmEEz04+aHA28fVVTBqE4icCKFfYO1I4dOpwqRKMT6PVyJar2x2lAb5/2J4Cbgcv6rJ8LnC4ifwAuNMa8OXAHa9as6TflyoIFC1iyZAn5fJ5coSdtU1MTPT09tLe3A2/PY9z7K9HQ0EAikShWMctkMiSTyeJyOp0mnU4Xe2em02kymQytra0YY0ilUmQyGdra2ujp6WHnzp1ks9niciKRoLGxkY6ODvL5fHE5l8vR1dWFiJDNZsnlcsWeztlsls7OzuJyU1MTqd13Z+ezz9K1227Ox9TR0UF7e3vZx5RMJqmvrx98TNOmUX/ssZilS5FrrnE6pq6uLnbs2BHZ51TqmPL5PFu3bq3ocwrj/972wpVs2Z9Thf/3Sh1TWRhjav4H/AA4r/A4CeRKvKYO+H/Ab0rtY9GiRUYb7e3t4ez4zDONufXWsjYNzckYY156yZhx44x57TWnzUJ1KhONTsbo9AIWG8fvfFQ5mjeB3kK4Y4C3Br7AGNNtjFkK7FJqB/l8Pjy7MukIaxK2CRPg9dfL2jQ0J4ApU+ALX4DFi502C9WpTDQ6gV4vV6IKNLcDhxUezwBWiUiziKQARCRR+LceeKjUDjS2XUMLfnvtBS+/XNamoQfkr30NfvUrO6NmQDT+SGh0Ar1erkQSaIwx64E/i8gXgQ8AVwHXAqeJyG7A0yJyHfAFoGRJ/lEzJS7YmSTLnB879P5G48bZvjXXXht4E419oDQ6gV4vVyLrjGKM+caAp+b3eXzwSNtrHGzWWMG0KMMybRr89a9lbRqaU18uuwz22w82boT999fh5IhGJ9Dr5Upsw6XGS8pcWCUvK7iiCc2pL83NdpbLgPNB1cTJEY1OoNfLldgGGo05mtA6EY4bZ4tO/eMfzpvWrGPjpZfC//wPPPfciC/V2NlSoxPo9XIltoFGI6HljUTKvqqpWS5rl11sEypAGQmN+TWNTqDXy5XYBhqNRZuz2Wx4Oy8z0ITqNJD58+G+++D554d9WU2dAqLRCfR6uRLbQNPdPWhUQuSE2p6eNs3pFnIvNW3jjxlj70B9Y2Cevz8a8w4anUCvlys+0FSRUAumH3hgoPzHQGpexH3+fNuv5tVXh3yJxsLyGp1Ar5crsQ00o44yA03N2W03+Nzn7EyXtSKXg5YW+Nvf7Dl67jnbzHzjDTvoU+EA3NGGvqIuARl1OZoyA00kbfzLL4fDDoMrr7R3zKrhZIzN/axfDxs2wNNP28Dyyiu2hMaYMdDYCA0NNnne0WH/duywyxMm2L/99oODDoJ3vxsOOQT23RcKAww1otXLldgGGq1Np9AKFU2aBDt32rKezc06nIZir73gYx+DG26wpT/Lddq+HVatsrfNH3jATqZ3xBFw6KHw+c/bsVaTJ8P48XamhqFobbVjxV5/3XYqfPZZW0tnwwZob4ejjrKB8UMfglmzbHEvJUTy+YWAnjPqiNZAExoidnra556zXwYNTsPx1a/CnDnwla/Yq4ygTsbAvffaQHDffXYfH/kILFkC++xjz4Mr2aztsbz//nDccf3XvfoqPPoo5v774aKL7BXSBz8IJ58MH/2oLaUaIT5H46k9Bx0UjzwNWNcjj4Rbbgn2+nwebroJpk+3ZUJPOskOu1i1ys6SWWjiVJ1Jk+CMM9i5dCk8+aS9yvngB+GOO2yXgtNOg9tug3dIpbuoiG2g0TjWqampKdw3KCNPE7rTcFxyCXz3u4OSsYOc7roLZsyAn/3MzvO9fj2cd17J/E5YFJ0mTYJzz7WB5m9/s4HmRz+yTbQvfQn+8peaOfXzijmxDTQahyCEPv7qwANtfsGBSMeEfeAD9krlwQf7PV10+vvf4cwzYcECW/D8gQdg7txwrlxGoOR5GjsWzjkHfvMbe7XT3Azvf7/9u/POqs1F7uwVQ2IbaDTmaELvXDV9OjzzzMiv60OkHb5E3r6q6UMul4M1a2DmTJvQfeIJOOWUSAJMP6fh2GsvO2h00ybblLv2Wut/yy02mEblFRNiG2hGJQceaP+jF+rQxoLPfhbWru1X5iJz8832SuanP7W9iON0VyWdhrPOgscfh3/7N/j+920+6uaba3KFE1diG2hGZY4mlbJ3nhympI28jd/UZG9FX3+9Xb7uOpq++114+GGbdFWC83kSsQnrhx6ySewbboD3vhd++9tovZQS29vbRmFvz5rkjWbMgKeesv+pA6Ail3XhhfaW/L77wnXX0bF6NY0BCmTVkorO05w58Pvfw+2322MtHCcHHhitlyJie0WjMUnWXosmzYwZ9hZsQGriNBL77GNvD190EdxzD23velfURoOo+DyJwOmn26vNk06Co4+Gq6+2vZOj9FJCbAPNqGXmTHtFEyd27ny7H0oVfuVVk0rZujxPPGET9zNmVL05FUdiG2g05mhqUt/V8YpGRc3Za66BM86wQwVefFGH0wCq7jR5Mvzyl/Cd79gpaebPtwE3aq+IiG2gGbVMmgRdXfDmoMk7dfLmm/CDH9hb3J/5DPznf0ZtVFtOOcX+MGzbBu95Dzz2WNRGkRDbQKOxlmpNJmQXeTshHIDIJ4m/8Ub41Kdgzz1tj9uVK2krTF2riVDPU3Ozvf39r/8Kp55q++METPJG/vlVidgGmlHNzJm2p2ocuO0227sWbFmGceNIrlsXqVJknHGGHV7xm9/YQZstLVEb1YzYBpqkoqH8vTQMGKUcGocfDn/6U6CX1sypFK+/bptORx759nNnn03TqlXROQ1Bzc7TxIl2qMXMmbaLwghNqUg/vyoS20CjsTp8zWYVfO97AweaSGc6fOop+4Xq63DWWSR//WubZ1JETc9TKgUrVti81Uc/CitX6vAKkdgehcYczc4y7iqUxcEHv11ZbgRq5lSKLVtgYJ+ZadPIT5tma80oIpLz9LGP2TFfS5faOcxL5G0i/fyqSGwDzagmmbRXCuvXR20yPMkklCjc1HnKKbbOjMeWFH3sMTsk45OffMfWvRkx0IhIVkTOFJHZheVpInJ82GIjUVdXF7XCIDKZTO3eLGDzqaZOA9lnn9JzPH3kI3D33aqKhkd6nsaPt536sllbgmLLFh1eVSTIFc0VQDswXUQuBzYBV1b6xiJysYicJyJXi0jzgHVzReSKwt/MkuIK2641TVAHDDSRJs0PPdTWnNm4sd/TiXe/246CVtTDOfKbC5mMzdUcdxwcdxy3f+81pk6FqVOTTJ0avFChVoJ8W/8EPG2M+SHwH8CpQEWpcBE5HJhV2Of9wNV91tUDK4CvA98GviclMr+jOkcDgQNNpG38ujo7yPCqq/pdvexsa7OFwBXlaVTkQkRgxQrWH/JZ3nPJMdRteoGJE3cWS+DEOdgECTTPAmcDGGPajTF3AjdV+L6nAb0VnJ4oLPcyG2gxlhzQCOxd4fu98zj4YHj5ZTudiGa++lVbfnThwv4FombPhkcfjc5LMR9ft5Dl5goe4jia85sBm7pZtChisQoIcr34BeAOABE5FlhnjFlZ4ftOAB4vPO4AJg5Yt6XPcu/6TX138NBDD/W7xb1gwQKWLFlCPp8vViVramqip6enOAK2d9xIb2/LhoYGEolE8dcsk8mQTCaLy+l0mnQ6TWuhJ2s6nSaTydDa2ooxhlQqRSaToa2tjZ6eHnK5HB0dHcXlRCJBY2MjHR0d5PP54nIul6OrqwspzCeUy+WK1e6z2SydnZ3F5aampiGPaezBB9P96KP0HHPMkMeUz+dpb28v+5iSyST19fWVHdOqVSTPPRfZf3/aL7gAPvUpWqdPp2HhQv7R0lLzz6nUMQFs3bo1lM/J9ZjGjIF7Jp/F01sPp37Peg5JtrBjR5pt2zJs3hzi5xTwmMrCGDPsH3AOMK7wOAV8ZqRtAuxzKXBp4fGuwMt91p0A3NFn+Ulg34H7WLRokdFGe3t7bd/w4ouN+da3hn1JzZ2G4+GHjfnkJ03PuHHG3HmnMamUMZ2dUVsZY3SdpylTjLFtTWPGjm0vPp4yJWozC7DYOH7ngzSdpgA/EZEvA+8BmssLaf24HTis8HgGsEpEmkUkBTwMTAIQkQzQZYx5YeAONM5301rrMTzve58tKTkMNXcajqOPhttuY9vPf26TDl1dkdYJ7oum87RsmZ10E2CvvaxXY6N9Pq4ECTT/B5wLvIJtRlV8v80Ysx74s4h8EfgAcBVwLXCasXmZxSKyEPgKcF6l7/eO5aijYpnnyB9+uO0DtG6dqlkhtTBvnh3wPmWKXZ4yxS7PmxetVyUE+ZTbgW5jzC9F5I/AydV4Y2PMNwY8Nb/PutXA6uG219iPpubzge+/vy0/8MYbdl5pDU4BSKfTtu/IxIkjv7hGaDtP8+bZv+3b0+yyS9Q2lTPiFY0x5nZgz8LiWEBFtWSNgabmnasSCVuLd5iBeRo7fHmn4Gj1ciVQrzdjzNOFf58scSUSCT5HU+Coo2wTZAg05R568U7B0erlir7utTHGRNGlfoRAE4nTCHin4Gj1ciW2gUbjEIRI6hjPmgV//OOQsyVqrK3snYKj1csVfd/WgEQ+NqUEkbSnm5vtdK1DTJWrsY3vnYKj1cuV2AYajWOdIqvvOkzzSWPNWe8UHK1ersQ20Ghsu0Y2q+Aw/Wk0znTonYKj1cuV2AYajTmayJpzw1zRaGxieqfgaPVyRd+3NSAaP4D6+vpo3nj6dHjttX4Fk3qJzGkYvFNwtHq5EttA43M0fairs3effv/7Qas0tvG9U3C0erkS20DjczQDOPZYWLt20NMa2/jeKThavVyJbaAZ1dOtlOKYY0oGGo25LO8UHK1ersT2KDR2ZIp0QvajjrKzVxYKLfWicZJ47xQcrV6uxDbQ5IfoCRslHR0d0b15Y6Odk/sPf+j3dKROQ+CdgqPVy5XYBhqNbdfIg1+J5lPkTiXwTsHR6uVKbAONz9GUoERCOHKnEnin4Gj1ciW2R+FzNCU45hjbca/Pr2DkTiXwTsHR6uVKbAONxkvK3mrxkbHbbjB5MmzYUHwqcqcSeKfgaPVyJbaBRmOORkUnwgF5GhVOA/BOwdHq5UpsA41GVOSNBuRpVDgNwDsFR6uXK7ENNNqKSYOdgCtyegNNoee0CqcBeKfgaPVyJbaBpru7O2qFQahoT++9t50wfuNGQInTALxTcLR6ueIDTRVRUTBdBObMgQcfBJQ4DcA7BUerlyuxDTSeYZg7txhoPB4NxDbQ+BzNMBx/vA00xuhx6oN3Co5WL1diG2h802kYpk2zU80+/7wepz54p+Bo9XLFB5oqouY/hUix+aTGqQ/eKThavVyJbaDxjEBv88njUUDNA42IzBWRKwp/M4d4zbUi8oqI/E1E9i31Go1jnZqaVExLDsCdE7by97tv48AbD2Dqd6Zyy1O3RK1URNN56kWjE+j1cqWmFb5FpB5YARwJpIH7RGSO6VOXU0T2BN4yxkwebl8ahyBoGX91y1O3cP76q3gm0cP0LUnWdm/i/LvOB2DejHkR2+k5T33R6AR6vVyp9RXNbKDFWHJAI7D3gNecBFwkIv8rIkcMtSONORotnasW3b+Itq42HpwKc161d+fautpYdP+iaMUKaDlPfdHoBHq9XAntikZEFgP7DXj608DP+yx3ABOBTb1PGGNWAitF5GTgDhGZbozZPnD/a9as6TcOZMGCBSxZsoR8Pl/8cJqamujp6aG9UN6yd8h9b2X5hoYGEokEO3fuBOz0o8lksricTqdJp9O0trYWlzOZDK2trRhjSKVSZDIZ2tra6OnpobW1lWw2W1xOJBI0NjbS0dFBPp8vLudyObq6uhARstksuVyumPTLZrN0dnYWl5uampyPaWfbTg7JHsL3zuphTMO72LW1jUn1k6Abtm/f7nRMyWSS+vr6qh5TLpdjx44dkX1OpY6ps7OTrVu31vRzCnJM27Zti+xzGuqYykFqOZuAiJwAzDfGfLyw/CTwCWPMC0O8/kbgRmPMEwPXLVq0yCxbtixMXWfa29tpaGiIWoOp35nKpm02do9PjaelqwWAKWOn8NKXX4rQzKLlPPVFoxPo9BKRJcaYxS7b1Lrp9DAwCUBEMkCXMeYFEUmJSHPh+b5OO4GSs9f76VaGZtkJy2hM2V/QpNiL1sZUI8tO0BGYtZynvmh0Ar1ertQ0GWyMyYnIYhFZCNQB5xVWfRw4DpgP3C0i/wAew17NlCzIoTFJ1t7eruIuQW/Cd9H9ixjTPYZMQ4ZlJyxTkQgGPeepLxqdQK+XKzWfV9YYsxpYPeC524DbCo8/XGundyLzZsxj3ox5tLS0MH78+Kh1PKOc2HbY09iPRmN9V+8UDI1OoNfLldgGGo/HEx9iG2g01lLVOCG7dwqGRifQ6+VKbAONx+OJD7ENNMlkzfPYI6KtvwN4p6BodAK9Xq7ENtBorA6vcVZB7xQMjU6g18uV2B6FxhxNb/dxTXinYGh0Ar1ersQ20Hg8nvgQ20BTV1cXtcIgMplM1AqD8E7B0OgEer1ciW2g0dh21Zig9k7B0OgEer1c0fdtDYjP0QTDOwVDoxPo9XIltoHG4/HEh9gGGo05Go1zTXmnYGh0Ar1ervhAU0U0/qfwTsHQ6AR6vVyJbaDRON9NbylJTXinYGh0Ar1ersQ20Hg8nvgQ20Djm07B8E7B0OgEer1c8YGmimjsXOWdgqHRCfR6uRLbQONzNMHwTsHQ6AR6vVyJbaDRiMaZGbxTMDQ6gV4vV2IbaDQOQdBYx9g7BUOjE+j1ckXftzUgGseAaGxPe6dgaHQCvV6uxDbQaBzrpLG+q3cKhkYn0OvlSmwDjca2q8ZZBb1TMDQ6gV4vV2IbaDTmaDQ257xTMDQ6gV4vV/R9WwOi8QOor6+PWmEQ3ikYGp1Ar5crsQ00PkcTDO8UDI1OoNfLldgGGp+jCYZ3CoZGJ9Dr5UrN2x8isjvwVWCzMWZ5ifWTgYuAl4Btxpj/GmI/YWqWhca8kXcKhkYn0OvlSs2PwhjzFvA8MFQHgeuBm4wxPwA+LSJ7lnqRxo5MGidk907B0OgEer1ciSpclhyoJCIZYK4x5oXCUxuBE0u9Np/Ph6RWPh0dHVErDMI7BUOjE+j1ciW0ppOILAb2G/i8MeYzw2y2G7Cjz3IHMLHUC3/3u9/1az4tWLCAJUuWkM/nyeVyADQ1NdHT00N7ezvw9q9Db4KtoaGBRCJRLACdyWRIJpPF5XQ6TTqdLg5sS6fTZDIZWltbMcaQSqXIZDK0tbXR09NDa2sr9fX1xeVEIkFjYyMdHR3k8/nici6Xo6urCxEhm82Sy+WKg0Sz2SydnZ3F5aampoqOqa2tjUwmU/YxJZPJqh9TLpejrq4uss+p1DH1fqGj+pyGOqbNmzeTz+cj+ZyGOqZykCiSqiJyDjDVGLN4wPNp4DVjzPjC8nXABmPMTQP3sXDhQvP1r3+9BrbB2bJlC+PGjYtaox/eKRganUCnl4gsGfjdHQkVmSYRSYlIszGmE3hERKYWVh0A3FtqG5+jCYZ3CoZGJ9Dr5UrNA42INAPvA2aIyB6Fpz8OXFt4/GXgIhG5ALjZGPNKqf1ozNH0XmJqwjsFQ6MT6PVypea3t40xW4GLBzx3G3Bb4fFfsbe/h0Vj/wKNnQi9UzA0OoFeL1dUNJ3eKWjs2+OdgqHRCfR6uRLbQKOxaHM2m41aYRDeKRganUCvlyuxDTTd3d1RKwxCY3vaOwVDoxPo9XLFB5oqorFguncKhkYn0OvlSmwDjcfjiQ+xDTQ+RxMM7xQMjU6g18uV2AYa33QKhncKhkYn0Ovlig80VUTjfwrvFAyNTqDXy5XYBhqPxxMfYhtoNI51ampqilphEN4pGBqdQK+XK7ENNBqHIGgcf+WdgqHRCfR6uRLbQKMxR6Oxc5V3CoZGJ9Dr5UpsA43H44kPsQ00PkcTDO8UDI1OoNfLldgGGj/dSjC8UzA0OoFeL1diG2g0Jsl668NqwjsFQ6MT6PVyJbaBxuPxxIfYBhqNORqN9V29UzA0OoFeL1diG2g8Hk98iG2g0VhLVeOE7N4pGBqdQK+XK7ENNB6PJz7ENtAkkzWfwGFEGhoaolYYhHcKhkYn0OvlSmwDjcbq8ImEvtPpnYKh0Qn0erkS26PQmKPpnTdZE94pGBqdQK+XK7ENNB6PJz7ENtDU1dVFrTCITCYTtcIgvFMwNDqBXi9XYhtoNLZdNSaovVMwNDqBXi9X9H1bA+JzNMHwTsHQ6AR6vVyJbaB55JFHolYYxIoVK6JWGIR3CoZGJ1DrdbzrBlLrcgsisjvwVWCzMWb5EK+5FjgH6AHmGmNeKPEao61UhIioK1/hnYKh0Ql0ehWcnPqX1LwBaIx5S0SeByaWWi8iewJvGWMm19bM4/GERVSZpuEmqzkJuEhEzgHON8b8cYjXvSEiE/osrwEerI5e2RwvIlE7DMQ7BUOjE+j0OtV1g9CaTiKyGNhv4PPGmM8UgshUY8ziYbY/Gfg+MN0Ysz0USY/HUxNCu6IZLogE3H61iNwN7AM8UQ0nj8cTDSruOolISkSaC4/7Ou0EnolEyuPxVI2a52gKAeV9wAQR2cMY8ybwceA4YD5wt4j8A3gMuNEYo6/DjMfjcaLmt7fLpVq3xWvpJCKTgYuAl4Btxpj/CtOn8J5zgVmFxdXGmA0lXlOz8yQiF2OT/xOB64wxW11ca+1UWP9b4CDgFeBoY0zosxWKyEHAN4BvGmPWDFgX1Xka0qmwPvh5MsbE5g84H1g8xLo9gUuUOf0a2LfP4z1DdqkHHgcEyAAPUfgxieI8AYcDPy08Pgb4totrrZ0Kz50GzKr1/6PCe98KHO/6mdbaqZzzpCJH40CQ2+L/KyJH1EqIIZxEJEP/q4WNwIkhu8wGWowlBzQCew94TS3P02m8nWN7orDcSxDXWjsBzAV+JSK3i8geNfDpS6n/S1Gdp+GcwPE8qRuxNdxt8eG2M8asBFYWbovfISJVuy1eptNuwI4+yx0M0Umxik6fBn5e4j039T4R5nkqwQTsr3Ffl77rtgznGoETxphLReRy4ErgZsL/cRiJqM7TsLieJ3WBxii8LV6mUwv2UreXMcBfq+EDpZ1E5MfYhHrf93xriO1r0X3gTewvcCmXvutKrQ+L4ZwAMDbXsFRE1tXAZySiOk8j4nKe4tZ06ofG2+K9TsaYTuAREZlaWHUAcG/Ib/8wMKngkQG6jDEvRHiebgcOKzyeAawSkWYRSQ3lGqJLEKfi+RGRemw+JBIUnKfhnJzPU5zuOjUDy7CXkhcZY94UkTOB44wx80XkHqD3tvhqY8xGBU7TgAuBF4Adxphba+B0MvYLVAfcY4xZH+V5EpEFhffbG/g2sBRYY4z5RSnXMF1GcgIeANYC9wHPATcbY3YMtZ8qO03BNnvvBr6FvdsT9Xkq6UQZ5yk2gcbj8cSXWDedPB5PPPCBxuPxhI4PNB6PJ3R8oPF4PKHjA43H4wkdH2g8Hk/o+EDj8XhCR90QBM87DxE5HWjAjs3aBEzHdpT7HPCEMWZ1hHqeGuADjSdURGQGtot6Ajgb+B6w0Rjzqoi8hr+qHhX4QOMJFWPMUwAi8kngfmNMC3b8jmcU4X9NPKEiIocVxnx9APi9iKRF5LiovTy1xV/ReMLmRKAN+D/gSOwo9tBLmnp04QONJ1SMMSUnjxaRScD7ge0i8gdjjIoaK55w8KO3PR5P6PgcjcfjCR0faDweT+j4QOPxeELHBxqPxxM6PtB4PJ7Q8YHG4/GEjg80Ho8ndHyg8Xg8oeMDjcfjCZ3/D1ktpYo+0jXAAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 288x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(figsize=(4, 4))\n",
    "\n",
    "ax.set_xlim(-1.5,1.5)\n",
    "ax.set_ylim(-1.5,1.5)\n",
    "ax.set_xticks(np.arange(-1.5, 1.51, step=0.5))\n",
    "ax.set_yticks(np.arange(-1.5, 1.51, step=0.5))\n",
    "ax.grid(ls = ':')\n",
    "plot_set(ax, '$x_1$', '$x_2$')\n",
    "\n",
    "u=np.array(list(map(lambda tau:reach(0,x0,tf,xf,tau), trange)))  # 理論式\n",
    "_, _, x = lsim(sys, U=u, X0=x0, T=trange)    \n",
    "ax.plot(x[:,0], x[:,1], lw=1, color='r',)     #　シミュレーション結果\n",
    "ax.scatter(x0[0],x0[1],c='b')\n",
    "ax.scatter(xf[0],xf[1],c='g')\n",
    "\n",
    "if (is_savefig):\n",
    "    fig.savefig(figpath+\"ans/ch4_1_phase_portrait.pdf\", transparent=True, bbox_inches=\"tight\", pad_inches=0.0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAIeCAYAAAClYUjYAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAB87ElEQVR4nO2deXgT1frHv2+apW1SKFBK2XdEZRNRuAjIovJTAQXEBVwAhSvgRVELIijUFXevK4JeF7iobIKC1wUQERBBERRBNgHZadlK2yxNc35/nCZNS5tm0iSTzLyf58nTTDLLO59M8+acM+ccEkKAYRiGYSrCoHYADMMwTGzDiYJhGIYJCCcKhmEYJiCcKBiGYZiAcKJgGIZhAsKJgmEYhgmIUe0AwsXIkSNFo0aN1A5DdYqKipCQkKB2GDEBu5CwBwl7KMHfRVZW1ntCiHsCra+ZRJGWlobp06erHYbq5OTkIC0tTe0wYgJ2IWEPEvZQgr+LrKysQ5Wtr5mqJ5PJpHYIMUFycrLaIcQM7ELCHiTsoQSlLjSTKBiGYZjIoJlEUVhYqHYIMUFBQYHaIcQM7ELCHiTsoQSlLjSTKBiGYZjIoJlEYTRqpl2+SiQlJakdQszALiTsQcIeSlDqQjOJgojUDiEmMBg085FWGXYhYQ8S9lCCUheqmCOi1kT0BRFdWc57vYjokeJHu+LXTET0FBHdRUSTytsnt1FI8vPz1Q4hZmAXEvYgYQ8lKHWhSqIQQvwJ4ByAUsUAIkoE8DyA5wC8AuANkkWF+wDsE0J8CMBERDeU3ee3334Lt9sd8dgZhmH0hpplMVc5r3UFkCMkTgDJABoBGAhgW/E6WwDcWHbDLVu24KeffopMpHGExWJRO4SYgV1I2IMkljwUFhaq+sNWqYtYawHOAHDKb9kBoF6Z172vnUe3bt18zzMzM5GVlQW32w2n0wkAsFqt8Hg8sNvtAEo6nXhvFUtKSoLBYPAVyywWC4xGo2/ZbDbDbDYjLy/Pt2yxWJCXlwchBEwmEywWCwoKCuDxeGA0GpGYmOhbNhgMSE5OhsPhgNvt9i07nU4UFhaCiGCz2eB0OuFyyTxqs9ngcrl8y1arNeA5GQwGFBQUaOqcQv2ccnNz4XQ6NXVOoXxOdrsdhYWFmjqnUD4np9MJk8mk+jl9+umn+Oc//wkiQps2bTBixAjcdtttMBqNUbv2CgsLkZOTg6ARQqjyAPABgJ5lXusD4DO/5a0AmgP4AUD74tcGAphTzv7E5ZdfLvROdna22iHEDOxCwh4kanrIzc31Pd+4caMwGAyCiAQAAUB07NhRHDhwIGrx+LsAMF1U8n0dE7cBEFEqEZkArAVQv/g1C4BCIcReAIsBXFK8elsAn5fdh8FgwKZNm3Dy5MkoRc0wDBMYIQSefPJJ9OjRw/eLv2PHjrDb7cjPz8eCBQvQpEkTbN68Gd27d8eBAwdUjrh81LrrqTGAVgCuIKIkAE8BuFHIdonpxXc2PQRgVPEmbwG4mIhGQLZtLCy7zwYNGkAIgRUrVkTlHGIVs9msdggxA7uQsAeJGh6mTJmCxx9/HL/99hs2btwIAEhISIDZbEZSUhJuuukm/Prrr+jSpQv+/vtvvPfee1GJS7GLyooc8fLo1auXACBGjBhRpSJZvGO329UOIWZgFxL2IIm2hzfeeEMAEEajUSxatCjguqdOnRIzZ84UHo8nKrH5u0C8VD2FA+9cFF9//bW3zUKXeBvRGHbhhT1Ioulh9erVGD9+PABg9uzZGDRoUMD1a9So4WvgjgZKXWgmUdSuXRsZGRk4cuQI/vjjD7XDYRhGp2RnZ2Po0KHweDyYPHkyhg8frmj7HTt2YPr06TH1g1czicJoNOKaa64BIEsVeoXro0tgFxL2IImWhxkzZuDo0aPo1q0bnnjiCUXbOp1O9OjRA1lZWVi5cmWEIlTuQjOJIiEhAX379gWg70QRS52K1IZdSNiDJFoennvuObz88suYN2+e4sFKLRYLHnzwQQBAVlZWJMLzHUcJmkkULpcLV199NYgIa9as0e3Y81wfXQK7kLAHSbQ8GI1GTJgwAQ0bNgxp+/vuuw/Vq1fH2rVr8euvv4Y5Oolu2ygA2U7RsWNHOJ1OrFmzRu1wVCGW6jXVhl1I2IMk0h7eeecdHD16tMr7SUlJ8bVrvPnmm1XeX3kodaGZROEdNlfv1U88d3gJ7ELCHiSR9LBx40bce++9uOSSS8JSmzF27FgAwH//+1+cOXOmyvsri1IXmkkU3rpAvScKro8ugV1I2IMkUh6EEHjooYcAAMOHD/eND1UVWrVqhV69esHhcGDx4sVV3l9ZdNtG4Z2P4h//+AdSUlKwY8cOHDx4UOWooo9e22bKg11I2IMkUh4+++wzrF27FmlpaZg8eXLY9nvvvffivvvuw2WXXRa2fXrR7ZzZ3jo3k8mE3r17A9BnqcLj8agdQszALiTsQRIJDy6XC5MmybnUpk+fjurVq4dt3zfffDNef/11tG3bNmz79KLUhWYShf/UfnqufuK5w0tgFxL2IImEh9mzZ2PPnj244IILMHr06LDvP1IodaGZROF/4t5EsWLFCt3NepeYmKh2CDEDu5CwB0m4PRQWFuLZZ58FADz99NMRaSw/d+4cZs2ahccffzys+1XqQjOJwn/O7GbNmqFFixY4c+YMNm3apGJU0Yfro0tgFxL2IAm3B5PJhEWLFmHs2LEYOHBgWPftxeVyYcyYMZgxYwZyc3PDtl/dt1F40Wv1E9dHl8AuJOxBEgkPnTt3xptvvlmq6juc1KpVC127dkVhYWFYv8t020ZRdtRFvSaKSF2w8Qi7kLAHSTg9RHOCtP79+wMAvvjii7DtU6kLzVxBZesHe/XqBZPJhI0bN+L06dMqRRV9wnEPt1ZgFxL2IAmXB4fDgXbt2uG6666LynfLddddBwD49ttvw9a7XKkLzSSKso3WNpsNV1xxBTwej65mvXM4HGqHEDOwCwl7kITLw5w5c3DkyBEcOXIkrLfDVsTFF1+M9PR0HDt2DDt37gzLPpW60EyiKK/OTY/VT3q7yysQ7ELCHiTh8CCEwKuvvgoAyMzMjEq1HhH5+oatWrUqLPtU6kIziaK8maH8E4VeBkbj+ugS2IWEPUjC4eHbb7/F9u3bUa9ePQwZMiQMUQXHddddh759+6J+/fph2Z9SF5rpiVPePczt27dHeno6Dh06hB07duCiiy5SIbLowvXRJbALCXuQhMPDK6+8AgAYN25cVCeEuuOOO3DHHXeEbX/cRuGHwWDQ3ax3TqdT7RBiBnYhYQ+SqnrYsWMHvvrqKyQmJuKf//xnmKJSB6UuNJMoKrov2Fv99NVXX0UzHNXw73iod9iFhD1Iqurhl19+gcViwZ133olatWqFKargsdvtWLt2LbZu3VrlfSl1oUqiIKJxRDSKiKYRUarf60RE24hoLxHtIaJcImpe/N71RHSo+DE02GP17dsXRITVq1frYqav8tpq9Aq7kLAHSVU93H777Th48CCmT58enoAU8u6776J79+6+6q+qoNRF1BMFEXUE0FkIMRvASgDT/N5uBmCAEKK5EKIFgEVCiL0kz6orgIZCiAZCiHll91tRfWHt2rXRuXNnuFyuiE5WHivYbDa1Q4gZ2IWEPUjC4aF27dqoW7duGKJRTpcuXQAAGzZsqPK+lLpQo0RxI4A/ip9vKV4GAAgh9goh/gIAIroEgHfC2PYA/g/AASIaVt5Oi4qKKjzg9ddfDwBYtmxZVeKOC7g+ugR2IWEPklA9uFwuLFy4UPUqvPbt2yMxMRE7d+7EqVOnqrQvpS7UuOspA4B3pD4HgHoVrHczgFcAQAixBcClRNQSwHIi2i+EWOe/8qpVq0oVpzIzM5GVlQW3242uXbsCAJYvX468vDxfZxNvy793gKykpCQYDAbk5+cDkLNAGY1G37LZbIbZbPZVYZnNZlgsFuTl5UEIAZPJBIvFgoKCAng8HhiNRiQmJvqWDQYDkpOT4XA44Ha7fctOpxOFhYUgIthsNjidTrhcLgAy87tcLt+y1WqF2+32fdBWqxUejwd2u10KdThgNBo1dU6hfk7Z2dlwuVyaOqdQPqfTp0/D7XZr6pxC+ZzOnj0Ls9ms+Jw+/fRTjBw5Ej179sSyZctUPacOHTpgw4YNWLVqFXr27Bny55Sbm+tbPyiEEFF9AHgSwP3Fz2sAOFjOOgRgdgXb3wrggbKvZ2ZmiorweDyiXr16AoDYvHlzhetpgezsbLVDiBnYhYQ9SELx4PF4RKdOnQQAMWvWrAhEpYzMzEwBQDz22GNV2o+/CwDTRSXf22pUPS0G0KH4eVsAnxNRKhH5d4S4EsD33gUi8o+zOoDzxuQIdE8zEfmqn5YvXx5q3HEB10eXwC4k7EESiof169fj559/Rq1atXD77bdHICpleNspfvrppyrtJ+bbKIQQvwLYTkT3ALgKwFQAT8GvraL4+RK/5ceIaBURPQxglxBiW9n9BmqjAKCbRKGoOKlx2IWEPUhC8eAdruPee+9FUlJSmCNSTseOHQEABw8erNJ+lLpQpWe2EOKFMi/dV+b9B8osZ1W2z8oSRZ8+fWA2m/HTTz8hOzsbtWvXDjLa+IK/FEpgFxL2IFHqYf/+/Vi8eDFMJhPGjh0boaiU0bhxYxw8eLDKQ3kodaGZDneVYbPZ0LNnTwghdNP5jmGY0HnjjTfg8Xhwyy23oF69iu65iS5EhAYNGkS9b4xmEkUw89Xq4TZZq9WqdggxA7uQsAeJUg8XXXQRmjdvjgceeCAyAVURUYWBTpW60EyiCGZqP2+i+Prrr1W/JzpS8JDSJbALCXuQKPUwcuRI7Nq1C5deemmEIgqN9evXo23btrj55ptD3oduhxmvrI0CAJo3b44LLrgAZ8+exfr166MQVfThzlUlsAsJe5CE4iEWh2hPTU3Ftm3bsGnTpspXrgDdDgoYLP369QOg/bufGIYJja+//hqPPvooDh8+rHYo5XLBBRcgKSkJBw4cqHIP7WDRTKIIpo0C0H47BddHl8AuJOxBEqyHZ599Fs8++ywWLFgQ4YhCIyEhAe3atQOAkEeS1W0bRbANO926dUP16tWxY8cO7N69O8JRRZ9g2mr0AruQsAdJMB5+/fVXfP/990hJScHIkSOjEFVotGnTBgDwxx9/VLJm+Si9JjSTKIJtnDGZTL5SxdKlSyMZkip4x3Nh2IUX9iAJxoO3g93dd9+NatWqRTii0Ln44osBhJ4olF4TmkkUShg4cCAAYMmSJeoGwjBMzHD06FF8/PHHMBgMGD9+vNrhBKSqiUIpmp4zuyL69u0Li8WC9evX4/jx46hTp04EI4suPD9yCexCwh4klXl4++23UVhYiIEDB6Jp06ZRiio0OnTogKlTp6JTp04hba/bObOVkJKSgquuugpCCHzxxRdqh8MwjMo4HA68/fbbAIAJEyaoHE3lpKen48knn8QNN9wQleNpJlEo7UB34403AgA+++yzCESjHt65ABh24YU9SAJ5sFgsWLhwIR588EF069YtilGpg9JrQjNVT0rp378/iAgrVqzAuXPnkJKSonZIDMOoBBHhyiuvxJVXXql2KEGzc+dOrFy5EhdddJFvEqNIoZkShdGoLOfVqVMHXbt2hcvl0tQggbEwFHKswC4k7EFSkYd4Hc7nq6++wrhx4/DJJ58o3lbpNaGZRBHKaIre6ict3f0Ui0MOqAW7kLAHSUUeBg4ciH79+mHv3r1RjqhqXHTRRQBCu/NJ6TWhmSsolF8F3kSxfPlyzYzZ7527l2EXXtiDpDwPf/75J5YvX46VK1eievXqKkQVOhdeeCEAWQWlFKXXhGYSRSi0aNECbdq0wdmzZ/H9999XvgHDMJri3//+NwDgzjvvRFpamsrRKKNevXpITk5GdnY2zpw5E9FjaSZRJCQkhLSdt1SxePHiMEajHhaLRe0QYgZ2IWEPkrIeTp48iQ8//BAAcP/996sRUpUwGAxo0aIFACgejkjpNaGZRBFqPeygQYMAyEQRzFDlsY7SRn0twy4k7EFS1sPs2bNht9vRt29fX31/vNGqVSsAwK5duxRtp/Sa0EyiCPXOhQ4dOqBFixY4ceIE1qxZE+aoog/XR5fALiTsQeLvweVy4fXXXweAmJ3BLhhatmyJtLQ0xf0iuI1CIUTkmylq/vz5KkfDMEw0+Pnnn5GTk4OLLroIffv2VTuckHniiSeQnZ2NUaNGRfQ4mkkUobZRAPAlikWLFsX9tJFms1ntEGIGdiFhDxJ/D127dsWBAwfw4YcfhnRrfawQarWi0mtClURBROOIaBQRTSOi1HLeX0FEh4hoAxElEJGJiJ4ioruIaFJ5+6xKomjXrh1atWqF7OzsuL/7ib8USmAXEvYgKeshIyMj5EH1Yg232x30nDxAHCQKIuoIoLMQYjaAlQCmlXn/RgBThBANhBBdhBBFAO4DsE8I8SEAExGdNxJWVfpBaKn6KS8vT+0QYgZ2IWEPEq+HtWvXauLGFS9dunRBUlISTpw4EfQ2Sq8JNUoUNwLwdiXcUrzsTy8Ai4hoMRF5x/8eCGBbgG2qjJaqnxiGKZ9du3ahe/fuaNu2rWb+z4UQcLvdiu98UoIa981lANhU/NwBoJ7/m0KI+4noQQCPApgD4JribU5VtA0ArFmzplRdY2ZmJrKysuB2u+F0OgHIeWI9Ho9vdifvmOwFBQXIyMhAq1atsGvXLixZsgQ9e/aExWKB0Wj03SFgNpthNpt92dhsNsNisSAvLw9CCJhMJlgsFhQUFMDj8cBoNCIxMdG3bDAYkJycDIfDAbfb7Vt2Op0oLCwEEcFms8HpdPpKSDabDS6Xy7dstVoDnlNRUREKCgp8d0EkJSXBYDD4ziEez8n/c1JyTv7LWjmnUD6n/Px8zZ1TKJ9Tfn4+XnrpJQDApZdeijNnzsT9OQFA8+bNsXHjRvzyyy+48MILgzqnoqIi5OTkIGiEEFF9AHgSwP3Fz2sAOBhg3Q3Ff38A0L74+UAAc8quO2XKFFFVHn/8cQFAjBo1qsr7UguHw6F2CDEDu5CwB8nhw4dFUlKSACD++OMPtcMJG1lZWQKAeOSRR4Lexv+aADBdVPK9rUbV02IAHYqftwXwORGlEpEJAIjIUPw3EcAav20u8d+m7E7DMVaTf/VTvI79xPXRJbALCXuQvPXWW7Db7bj22mvjtoNdeXhn49u3b1/Q28R8G4UQ4lcA24noHgBXAZgK4CkANxJRLQDbiOjfAEZClj4A4C0AFxPRCAAuAAsjEdvFF1+MNm3a4NSpU3E79LhQcOeD1mEXEvYAOJ1OvPvuuwCABx98UOVowksoiULpNaFK334hxAtlXrrP7/l5qV4I4QSQGWif4RpK+Y477sCkSZMwZ84cDBgwICz7jCZK5g7XOuxCwh6AuXPn4vjx42jbti369OmjdjhhJZREofSa0EyHu3CNZzN06FAQEb744ouIj8gYCXgAuBLYhYQ9yAbchg0bYtKkSXHdwa486tati6effhqvvfZa0CUF3Q4KGK5Zqho0aIBevXrB6XRi4cKI1HBFFJ4fuQR2IWEPwG233YZNmzbhlltuUTuUsGMwGPDoo4/i1ltvDToJKr0mNJMowlkPe8cddwAA5syZE7Z9RguPx6N2CDEDu5CwB0lCQgKPpFuM0mtCM4kinNM9Dho0CElJSVizZg32798ftv1GA/5HKIFdSPTsYfXq1Xj88ceRnZ2taQ+//fYbXnvttaBHwNbtMOPhvAiqVauGG26Qo4TMmzcvbPuNBomJiWqHEDOwC4mePTz55JN48skn8e6772raw1dffYX7778/6AnYlLrQTKIIVxuFF//qp3i6vZDro0tgFxK9eti4cSNWrVqFlJQUjBkzRtMelN75xG0UYeLqq69Geno6/vzzT2zcuDGs+44kXB9dAruQ6NXDjBkzAABjx45Famqqpj00a9YMQPCJQrdtFOG+5c1kMuHOO+8EAF9HnXggnG018Q67kOjRw7Zt27BkyRJYLBbfDHZa9uBfogjmR7NSF5oxF4lORXfffTcA4OOPP8a5c+fCvv9I4B3EjGEXXvToISsrC0IIjBo1ChkZGQC07aFGjRqoVq0a8vLyghrsT6kLzSSKSAwZ3Lp1a3Tv3h35+fn49NNPw77/SOBwONQOIWZgFxK9edi9ezcWLlwIi8WCRx55xPe6lj0QkaJ2CqUuNJMoIlX/eM899wAAZs+eHZH9hxutjLEfDtiFRG8eWrRogVWrVuGll15C/fr1fa9r3UPjxo2Rnp4eVO2HUheaSRSR6pZ/0003oXr16ti4cSN+++23iBwjnGi5HlYp7EKiNw9EhF69emHcuHGlXte6h88++wzHjx8PaiwrbqMIM8nJyRg2bBgA4L333ovIMcKJluthlcIuJHry8Pfff1f4ntY9KPny5zaKCOCtfpozZ45vlqhYxTurFcMuvOjFw9atW9G0aVPceeed5d75oxcPwdz1pNSFZhJFJO+RvuSSS9CpUyecPn0aH3/8ccSOEw7C3fEwnmEXEr14eOyxx+DxeFCzZs1yq6K17uG3335D06ZNceWVV1a6rlIXmkkUkaSoqMg36uSrr74a0z21tTaEclVgFxI9ePjhhx/wxRdfwGq1YvLkyeWuo3UP1atXx/79+7F3795K11XqQjOJwmw2h3V/p06dwttvv43evXujevXqyMyU8yb9/vvvWLt2LQDg8ccfx8CBA/HOO+/g6NGjYT1+qNhsNrVDiBnYhUTrHoQQmDRpEgDg4YcfRp06dcpdT+se6tWrByLC0aNHKy0xKHWhmURRVFQUlv3k5ubi0UcfRePGjTF27Fh89913yM/PR/369VG3bl0AwOuvvw5A/opZsmQJ7r33XjRs2BADBgzAsmXLVB0qQC/1sMHALiRa97B06VL8+OOPqF27Nh566KEK19O6B5PJhHr16kEIgcOHDwdcV7dtFOFKFKdOncJrr72GvLw8XHXVVZgzZw6OHz+OQ4cOYdOmTTAajVi8eDEOHjyI//znP3jnnXfQv39/36x4/fv3R4cOHfDtt9+GJR6luFwuVY4bi7ALiZY9uN1uX1XT448/jpSUlArX1bIHLw0bNgQAHDx4MOB6Sl1oJlGEiyZNmuCjjz7Cjz/+iG+//Ra333470tPTAQD169fH4MGDUVRUhJkzZ6Jp06YYPXo0Pv/8cxw+fBgvvPACGjRogN9//13TI1UyTKxgt9vRo0cPtGrVCqNHj1Y7HNUJNlEoRTOJItQ2CiEExo8fj//85z++1wYNGoQuXbqUu/6//vUvAMCsWbNKJYP09HQ8/PDD2LNnD+bNm4f+/fv73vvwww8D3t8dTrReD6sEdiHRsoeUlBS888472Lp1a6XfAVr24KVRo0YAKk8U3EahACEEMjMz8frrr+P+++/H6dOnK92ma9euuOyyy5CTk4P333//vPctFgtuu+02X+eXHTt24O6778bFF1+Mt99+O+LtF3ooXgcLu5Bo1YP/3YfBTMSjVQ/+9OvXDzNmzEDv3r0DrqfbqqdQEsXrr7+Ol156CUajER9//DFq1KhR6TZE5Bto7IUXXqj07oLq1avjhhtuQF5eHsaOHYtevXph9+7dimMNFj38MwQLu5Bo0cPu3bvRrl07fPnll0Fvo0UPZenZsycmTZqEyy67LOB6cZEoiGgcEY0iomlElFrmveeJ6Fci2kxEF/q9fj0RHSp+DK1qDN9//z0efPBBAMBHH32Efv36Bb3tDTfcgFatWuHAgQOYP39+wHXr1auHRYsWYeHChahTpw7WrFmDDh064N13343p/hgME8s8+OCD2LZtGxYuXKh2KPpACBHVB4COAD4qft4NwCt+77UH0Kb4+RgAC4ufE4CnAVBF+3300UdFsGRnZ4v09HQBQEyaNCno7fx59913BQDRtm1b4fF4gtrm5MmTYtiwYQKAACCysrJCOnYgCgoKwr7PeIVdSLTm4X//+58AIFJSUsTRo0eD3k5rHsrD7XaLxYsXi7feeivgev4uAEwXlX1vV7ZCuB8AngAwqfi5DcC+CtZrD+C14ucdAPwC4G8Aw8pb/5FHHgla5q233ioAiCuvvFK43e6gt/PH4XCIevXqCQBi+fLliradO3euaNy4sThw4EBIxw5Ebm5u2PcZr7ALiZY8OBwO0apVKwFAvPDCC4q21ZKHiigqKhImk0kAEHa7vcL1/F0EkyiMkSqpBCADwKbi5w4A9SpY72oALwOAEGILgEuJqCWA5US0Xwixzn/lVatWleqWnpmZiaysLLjdbl/nEqvVCo/Hg3vuuQd79uzBW2+9BafT6bt7KSkpCQaDAfn5+QBkw7TRaPQtm81mmM1m5OXlAQDGjBmDxx57DFlZWbjssstgNpthsVhQUFAAj8cDo9GIxMRE37LBYEBycjKuv/569OrVC4mJiXA4HLDb7ZgzZw5uueUWpKamwul0+uoQbTYbXC6Xb9lqtZZ7Tt7BCh0OBxISEkI+J+855OXlQQgBk8kU1Dk5HA643W7fstPpRGFhIYgINputSufkHelS6TmdOHECTqdTU+cUyud0+vRpFBYWauKcXnzxRezatQstW7bE7bffDofDEfQ5nT17FiaTKebOKdzXXkZGBg4ePIg//vgDjRs3Lveczp49q6zTXWWZJNwPAE8CuL/4eQ0AB8tZpz2AnhVsfyuAB8q+npmZGXzaFSLo6qJA5Obmilq1agkA4n//+1/I+3nhhRcEANG1a1exf//+KsWUnZ1dpe21BLuQaMXD9u3bfb+Wf/jhB8Xba8VDZXTp0kUAEGvXrq1wHX8XCKJEoUZj9mLIqiQAaAvgcyJKJSITABBRGwDVhRCriagaEdmIyD/O6gBWlN1pMPNR+N9bHI4BwlJSUnxjzDz22GPeRKaYyy67DPXr18f69evRoUMHLF68OOSYrFZryNtqDXYh0YqHY8eOIS0tDf/85z/RrVs3xdtrxUNl1KsnK2mOHDlS4TpKXUQ9UQghfgWwnYjuAXAVgKkAngJwIxH1ALAMwMtE9DOAlQDyATxGRKuI6GEAu4QQ28rZb8DjHjhwAK1atcKQIUPCOtzwuHHjUKdOHfz888/4/PPPQ9rHlVdeiS1btqB///44c+YMBg8ejDFjxoQ094Wa40zFGuxCohUPvXr1wo4dO/Dcc8+FtL1WPFRGMIlCqQtVbo8VQrwghHhXCPG4EOK0EOI+IcQCIcQaIUQTIUSn4sdlxaWjLCFEbyHEi0KI78rbZ2UTFz355JNwOBwwmUxhnQ0vOTkZjz76KAA51kyoF2NaWhqWLl2Kf//73zCbzZg5cyYuv/zygB92ecT6xErRhF1I4t2Dfx+p6tWro3r16iHtJ949BEswiUKpC810uAvE/v378eGHH8JgMCArKyvs+x89ejQaNGiA3377rdJ+FYEgIowfPx4bNmxAq1atYLPZULt27TBGyjDxhcfjQd++ffHoo4/qosNcOPAmitzc3LDtUzOJIlAp4ZlnnoHb7cbQoUPRsmXLsB87MTERjz/+OABg8uTJcDgcVdrfJZdcgl9++QWLFi3yndfJkydx9uzZSrfV+rzASmAXknj28Oqrr2LlypV47733grr+AxHPHpRw6623wul04u23365wHd3OmV0Rhw4dwgcffACDwYCpU6dG7DgjR45E27ZtsX//frz66qtV3p/NZvP9MvB4PLjzzjtxySWX4KeffqryvhkmHvjjjz981brvvfcel66DxGKxhH0iN80kiooaqN966y0UFhbipptuwgUXXBCx4yckJODll18GIEswx48fD9u+z5w5g2PHjmHfvn3o1q0bnn/++QrbQnh48xLYhSQePdjtdgwdOhROpxP33HOPoiF2KiIePUQKpS40kygqwuVywWKx4P7774/4sa666ipcf/31OHfuHB577LGw7bdmzZpYv349JkyYALfbjUmTJuHaa68NazJimFhi/Pjx+O2339CyZUvfDzAmeK666io0bdrU1/GvylTW0aK8B4DaoWwXycfkyZMr7FySk5MTlg52wbBjxw5hNBoFEYmNGzeGff/Lli3zdfKrU6eO+Oabb0q9n5eXF/ZjxivsQhJvHj7++GMBQCQmJootW7aEbb/x5qEqNGvWTAAQO3fuLPd9fxcIV4c7IrrR73kqgA/Ck6bCR6AOdLVq1QpLB7tgaN26NR544AEIITB69OhKb9tVyvXXX4+tW7eiZ8+eOH78+HlTrnrnwWDYhZd489CnTx/07dsXb7zxBtq3bx+2/cabh6pQ2S2ySl0Eu3ZvIrqFiK4DsB5AmMoz4aNsG8Wff/6Jzz77TJVb6qZPn47GjRtjy5YteO2118K+//r162PFihWYNWsWnnrqKd/rRUVFvjFnGLCLYuLNQ+3atfHll19i5MiRYd1vvHmoCpUlCqUugk0UnwK4A8BMALcIIW5RdBQVePPNNzFo0KCI3ulUEVarFW+99RYAObTHgQMHwn6MhIQEjBo1ynd3w8mTJ9GuXTt89NFHuumBymgHl8uFV1991feDz2AwRK0WQIsE0+lOCcEmim8gh9ZoBuAKIno8LEcPIwkJCb7nTqcT8+bNAwAMHVrlOY5C4rrrrsOQIUNQUFCA0aNHhzwOVLD85z//wfbt2/HQQw+hZ8+e2LlzZ0SPFw9YLBa1Q4gJYt2DEAJjx47FhAkTwl6K8CfWPYSTyhKFUhfBJoqhQoiZQgi3EGIm5LwQMYV/nduyZctw6tQptG/fHh06dFAtptdeew21atXCN998gzfffDOix3r44YfxySefID09HT/88APatWuHp59+Wte9WY1GNUbRjz1i3cOzzz6L9957D4mJiRG9OzHWPYSTyhKFUhdBJQohxNIyyx8oOkoU8G+j+OCDDwAAw4cPVyeYYjIyMjBr1iwAcn6MHTt2ROxYRIRbbrkFa9euxYgRI+ByuTB16lR06NAB69atq3wHGkRPddKBiGUPb7/9NqZMmQIiwpw5c9CpU6eIHSuWPYSb9u3bY/z48RgwYEC57yt2UdltUfHy8M5HcerUKWE0GkVCQoI4fvx45feRRYG77rpLABAdO3YUTqczosfyjjO/YsUK0aJFi5DH7tcCepl/oDJi1cPHH38siEgAEO+8807EjxerHtQgHuajiAjeNoqlS5fC7XajZ8+eSE9PVzkqyWuvvYYmTZpg8+bNePDBByN6LG/jdp8+fbBt2zYsWbKk1Nj9c+fOxalTpyIaQ6wQ7mEM4pVY9PD999/jjjvugBACzz77LEaPHh3xY8aiB7VQ6kJziSIxMRFt2rTBTTfdpHJEJVSrVg2ffvopzGYz3nzzTcyZMydix/K/ACwWC2644Qbf8g8//IA77rgDTZs2xRNPPBHW0SVjEf5ikMSih86dO6Nnz5545JFHfJN/RZpY9BBJNm7ciCVLlpTbTqnYRWVFjnh5lJ0KtaioKKQiWSR55513BACRlJQU1h6n/gQqXm/dulVcc801AoAAIGrWrClmzJih2UnnuapBEkseXC6X77nT6YzaiAlCxJaHaNCgQQMBQOzbt++893Rb9VSWWOyFOWrUKIwcORJ2ux39+/cP2z3OwdKuXTt8/fXX+P7779GtWzecOnUKjzzyCBo2bIgnn3wyqrEw+sLj8eCRRx7B9ddf77vxxGw2c1+JCFKnTh0ACMuYcLH3bRoiCQkJWLRoEfbs2aN2KBVCRHjzzTfRtWtXHDx4ENddd13Yq3+CKVL26NEDa9aswVdffYXu3bvj7NmzpdottNJhT29VDRWhtoeCggLccssteO6557Bq1Sr8+OOPqsShtodoEyhR6LaNQgiBO+64Ay1btozpUVUTExOxdOlStGrVClu3bsXgwYPhdDrDtv9gO9IQEfr27Ys1a9Zgw4YNePjhh33vvfPOO7j00ksxc+bMuG7H0FMHq0Co6WHXrl3o2rUrFi5ciGrVquGrr75Cjx49VIlFb9dDRkYGgPITRaQ63MU8e/fuhd1uxyWXXOLLpLFKWloa/ve//yE9PR0rVqzAkCFDwpYsQhlWuHPnzqhfv75vecGCBdi8eTPGjBmDunXrYtiwYfj888/DmtCiQdiGWI5z1PLw6aef4tJLL8XWrVvRokULrF+/HldddZUqsQD6ux6834PHjh077z2lLjTTVXHv3r0A5NAZ8UCzZs3w9ddfo0+fPvjiiy8wZMgQLFiwoMq/ekQYhgr58ssvsXjxYsyePRurV6/GvHnzMG/ePFSvXh1ZWVlRmdsjHITDRSzjcDhw5swZnDlzBqdPnz7vuffviRMnUFRUBLvdjoKCAhQUFPie2+122O32Ug2XQIk7i8WC5ORkJCUlITk5GcnJyUhNTUV6enqpR5MmTdCiRQs0bNgQCQkJ+Prrr3HrrbcCAG6++WbMnj0b1apVU80VoP3roSyBqp6UutBMoti3bx8AOQx3vNChQwesXLnSlywGDBiABQsWVOkfKtDc4cGSmJiIoUOHYujQofjrr78wf/58fPrpp9iyZQtq1KjhW2/9+vX48ccfcc011+Diiy+OuRsIwuEi0rjdbpw5cwYnT57EqVOncOrUqVLPT506hdOnT5/35X/mzJkqz80eDIWFhYp+fZpMJjRt2hQXXXQRWrRogWuuuQaPPfaY6knCG5ueCJQolLogNbIsEY0D4AJQD8C/hRBn/N7rBaBz8eKXQojfiMgEYBqA3QAyhBDPlbNPUbNmTZw4caLUAIHxwJYtW3DNNdcgOzsb7du3x/Lly0tVBSnB4XAgMTExzBFKdu7cifr168NmswEAxowZg5kzZwIAUlNTccUVV6Bbt27o1q0bLrnkElit1ojEESyRdFGWoqIinD59utwv+0DPz549G/IxTSYTatSogdTUVN/f8p4nJSUhNTW1VMnA/29iYqLvf4aISt2J5HK5fKWQgoIC5Ofn+0op3sfRo0fx+++/Y/v27RVWTzZv3hx9+vRB79690atXL1U6w0bzeogFzp07h9OnTyM9Pf288/Z3QURZQojpgfYV9URBRB0BPCCEuJOIugEYLISYUPxeIoAfAFwOwAzgWwBXAngAQK4Q4j0imgrgd1Fm/CkiEkOHDsV///vfKJ5N+Ni7dy+uvfZa7N69Gw0aNMCCBQvQpUsXxfs5deoUatasGYEIz2f58uWYP38+vvvuOxw8eLDUe927d8eaNWsAyF/NX3/9NVq1aoUmTZpE7ZddqC6Kiopw8uRJZGdnIzs7GydOnPA9P3nyZLlf/GfOnAkpRiJCamoqatWqhZo1a6JmzZqlnnsf/gnA+zcpKSmo20sjeU3s3LkTzzzzDObMmQMhBIYNG4aJEyfit99+w+bNm/HLL79g8+bNpUolRISuXbti8ODBGDx4MBo1ahSR2MoSzf+NWMffRawmiicA5AshniMiG+SXftPi93oDyBRCXFu8/DOAwQDmFL/+ExH1g0wuI/z3azKZxAcffIBhw4ZF9XzCSU5ODm688UasW7cORqMRzz//PB544AFF95rn5OQgLS0tglGWz99//421a9di7dq1WL9+Pfr06YOXXnoJALBt2za0bdsWgLyNuUmTJmjZsqWvTvuuu+7yFZOLiorCViL0d+F2u3H8+HEcOXLE9zh69Oh5iSA7OxunTp0KqT67Ro0a5X7ZB3pevXr1iJeAI3FN/Pzzz3j++eexcOFCCCFgMplw33334bHHHitVPQlI95s3b8aqVauwatUqrFmzplTJo1u3bhg1ahSGDBmCpKSksMbpj1r/G7GIv4tYTRSzAGwSQswmIiNk0rAUvzcUwPVCiGHFy2sBZAL4sPj13UR0FWTS6Ou/36ZNm4r9+/f7lrt27YqePXvC4/GgqKgIgCyqCyF805N6f9l6OwAZjUYQkW85ISEBBoOh1HJCQoKvS3zZZYPBAKPRiMLCQtmbscwyEcFkMsHtdsPj8Zy37N3nihUr8PPPPwOQjd59+/b1VfeYTKaA5+TxeGCxWGLqnE6ePIkVK1bg9OnTOHfu3HnXxL333ovU1FQUFRVh+fLl2LNnj69KJDEx0VdFkpaWhvbt28NoNMLj8WDnzp0wmUwwm82+uvSCggI4HA5fQ6+3ukTpaJn+x01OTobNZvPFk5SUhJSUFJjNZpjNZt+yEKLU5xIr157T6YTFYqn0czKbzSgqKvLFXHbZe07r16/3lRYNBgPatGmDyy+/3FedVNk55efnY+/evdi9ezf27t3rc2KxWNCuXTt07NgR1apVC8u1538OLpcLNput3HOKhc8plHMK9DkVFRXhiy++QF5eHm677TbftW0ymeB0On1tii+88EKliUKNxuzjAJKLn6cAOFHBe/7vB9oGADBs2LBS04LGM88++yw+++wzjBw5En/99RfmzJmDZ555BmPGjKl0HPlYr4e12+3466+/sGfPHuzZsweHDx/Gk08+6WvP+O2337B9+3a4XK7z6u87d+6MFi1aYP/+/di1axeWLFmi6NhGo9F3F4/VakVKSgrS09MxevRo1K5dG7Vr18a3336L5ORkpKSkwGq1+u70sVqtqF+/PmrXrg2g5K6ReOhZXJVrIjs7G1988QVSU1MxaNAgAMD27dvRvXt3DB8+HBMmTECDBg1Cju3cuXP45JNPMHv2bGzatAmbNm3C5s2bMXToUEyaNAkXX3xxyPsuS6z/b0SCefPm4fDhwxg3bhwaN27se93fxQsvvFD5jiob4yPcDwCXAHi/+HkPAG8CSAVgAmABsLH4PQuAn4ufTwAwvPj54wCGlN3vxIkTzxvPJN45evSoGDx4sG9sptatW4vPPvss4Pg4J0+ejGKE4aOoqEgcPHhQrFixQrz00kti5MiR4oorrhANGjQQCQkJPgcVPUwmk7BYLMJisQij0Sj+9a9/iTlz5oiNGzeKsWPHVrhd48aNS8VRq1atCtd98sknfet99tlnIiEhQVSvXl00adJEdOjQQfTq1UsMGjRI3H333eLMmTO+dX/44Qfxv//9T2zYsEHs3LlTnDhxotSYR5Em2GuiqKhI7N27VyxYsEDcd999ok2bNr5z79KlS6l1IzFc/ubNm8WwYcN8nzcRibvuukv8/fffYdl/vP5vVIWOHTsKAOKnn34q9bq/CwQx1lPUSxRCiF+JaDsR3QOgEYCpAJ4C8L0QYgERTSeiSQASAIwq3uwtAE8R0QjIu6UWlrPf6JxAFMnIyMDChQuxZMkSPPzww/jzzz8xcOBAdOzYEQ899BCGDBlyXsNwrA+/cfr0aezcuRN//vkndu3ahV27dmH37t3YvXs37HZ7hds1aNAALVq0QJMmTXyPxo0bo0mTJqhfv365DeTeetimTZti7Nixvmoo/zt4yg5lMHz4cJw+fbrUut6/3lnDADksRVFREc6ePVvunUuvvPKK7/ljjz2G1atXn7eOzWbD0KFD8c477wAATpw4gUmTJqFGjRq+h7fxukaNGmjbti1SUlIqdVwW7zVRWFiIs2fP4syZMzhy5AgOHDiA7t27o0mTJgDk5Fovv/xyqW0tFgt69uyJwYMH+6pGgMgMh3HJJZdg7ty5eOKJJ/Diiy/i3XffxYcffohPP/0UDzzwAKZMmeKrgg2FWP/fiAQV3SKr2EVlmSReHlosUfjjdDrF66+/LtLT032/8urXry8mT54s/vjjD996sfCrye12iz179ohly5aJl156SYwaNUp07969VOzlPWrXri2uuOIKMXz4cPHMM8+IBQsWiK1bt4q8vLyQ4oi0C5fLJU6ePCn27t0rfv75Z/Htt9+K+fPni3feeadUqW/SpEni6quvFp06dRLNmzcXNWvWFAaDQQAQd999t2+9X3/9NaAf/wmoxo8fL2w2m6hVq5aoV6+eaNq0qWjdurVo166duPXWW33rud1ukZGRIWw2W7n7/PDDD33rvvvuuyIjI0Ncc8014oknnhDff/+9sNvtEXUYiD179oibb77ZF2ujRo3EsmXLQt5fLPxvRJvhw4cLAGL27NmlXo/5EkWk0HpnGrPZjPvuuw9333035s6di5dffhl//vknnn32WTz77LNo06YN/u///g89e/ZEr169kJycXPlOq0BRUREOHz6Mv/76C3/99Rf27t3rKyns3r27wrm6k5OTccEFF6B169Zo1aoVWrVqhZYtW6Jly5ZITU0Na4yRdmAymXx3MAVixowZ570mhMC5c+dK/bKrX78+3nvvPV8Hu7IP/6Fp8vLyfI+y+N9FZTAYfEM4GAwGVK9eHTVq1ECdOnXQuHHjUv11RowYgbvvvjt4ARGmefPm+PTTT/Hggw9i7Nix2Lx5M/r164ebb74Zb7zxhq+9KFgifT3EIt7xnsoO46HUhSod7iLBI488Isr7h9QqHo8Ha9euxdy5c7FgwYJS9/EnJCTgoosuQseOHXHxxRf7qmrq169f6f33hYWFOHfuHHJzc3H8+HEcO3YMR48e9T0OHDiAv/76CwcOHCg1T3lZGjRo4EsIrVu39j2vX79+1HpwnzlzJuzJJ1YoLCyE3W6H0+ks9XC5XLBarWjZsqVv3W3btqFRo0aw2Wwx13s+WNxuN15//XVMnToVBQUFyMjIwEcffYSrr7466H1o+XqoiFdffRUTJkzAfffdh9dff933ur+LYG6P1UyJQm/1jwaDAT169ECPHj3w+uuvY926dfjmm2/w5Zdf4o8//sDvv/+O33//vdxtTSYTqlWrVuqXZ1FREfLy8hQN/Fe3bl00bdoUzZo1Q7NmzUqVFKpSlxwuvLc4ahGTyRR0KTojIyMmhtCoCkajERMmTMCgQYNw11134fvvv8c111yDhx56CM8880xQbSZavh4qom3bthg8eDDat29f6nWlLjRTopg0aZJ47rnzRvbQHadOnUJiYqKvZ+zu3buxf/9+7Nu3D8eOHcPZs2cDjhGUkJCAlJQUpKSkoHbt2qhbty7q1q2LjIwM1K1bF40aNUKzZs3QpEmTmC/Kc09cidY8FBUVYcaMGZg2bRqKiorQtWtXLFq0yFfNUhFa81AVlPbM1kyJQuttFMGSnJyMxMREdOnSpcIhQJxOJ3Jzc0uVwogIKSkpSExMjIu+AcEQ64ksWmjNQ0JCAqZMmYI+ffpgyJAhWL9+PTp16oQlS5agU6dOFW6nNQ9VQamL+KywLAc9FivLI5iqI4vFgtq1a6NOnTq+R3p6etBjB8UL8TZ/RqTQqocuXbpg06ZNuOKKK3D48GF0794d8+fPr3B9rXoIhMfjwYEDB7Bly5ZSryt1oZlEobc2iooI1MCsN9iFRMseMjIysGrVKowePRoOhwO33nprqUZbf7TsoSIKCwvRpEkTXHbZZaX6mil1oZlEwUi0VCKoKuxConUPZrMZM2fOxIwZMyCEwPjx4zFlypTzOuFq3UN5WCwWVKtWzTfviRelLjSTKPQ2cXpFxMLdRrECu5DowQMRYdKkSXj//feRkJCAZ555BuPGjSuVLPTgoTy8/U2ys7N9ryl1oZlE4R05Ue/osR62ItiFRE8ehg8fjiVLlsBiseDtt98ulSz05MEf73DiOTk5vtd020bBiUJSUY9oPcIuJHrz0K9fPyxduvS8ZKE3D17KK1EodaGZRMEwDOOlb9++pZLFAw88oMmBQ4OhvEShFM0kCm6jkOi1HrY82IVErx68ycJsNuO1117zze+uN7iNwg+uepLotXhdHuxComcPffv2xdy5c0FEePzxx/Hee++pHVLUGT16NNavX19qwEel14RmemZzopDo+UuhLOxConcPQ4YMQXZ2NsaNG+ebzXDAgAFqhxU1mjdvjubNm5d6jdsoGIZhyjB27Fg89NBD8Hg8uO22287rqcwERjOJgsd6knjnnmbYhRf2IHniiSdw1113oaCgAAMGDDhvjgatcuTIEdx///2YOnWq7zWl14RmEgUP4SHhMa9KYBcS9iApKirCO++8g3/84x84ePAgBg4cGHAkZa3gcDjw2muvYe7cub7XlF4TmkkU3EYh0WunovJgFxL2IHE6nbBYLPjss8/QqFEjbNiwAaNHj9b8bbPl3fWk2w53DMMwwVCnTh18/vnnSE5Oxpw5c/DWW2+pHVJEsdlssFgsKCgoQEFBQUj70Eyi4DYKCddHl8AuJOxB4u+hffv2eP/99wEAEyZMwKZNm9QKK+IQ0XnDeOi2jULrxcdg4baaEtiFhD1Iynq4+eabcd9996GwsBBDhgzBqVOnVIos8pStflJ6TUQ1URBRLyJ6pPjRrpz3hxLRj0S0l4gG+b1eh4j+IqJDRFTuYPPcYCex2+1qhxAzsAsJe5CU5+HFF1/E5ZdfjgMHDuDOO+/UbFItmyiUXhNRSxRElAjgeQDPAXgFwBvkNyg6EZkAuIUQ/wBwIwD/isNRANoKIRoIIf4VrZgZhtE2FosF8+fPR40aNbB8+XK89NJLaocUEdq1a4cuXbrAYrGEtH00SxRdAeQIiRNAMoBG3jeFEIVCiPnFz38HsAsAiCgJQE8Ae4loBhGV25uc2ygkPC9wCexCwh4kFXlo3Lgx5syZAwCYMmUKNm/eHM2wosKLL76IH3/8Eb169QKg/JqIyBAeRDQdQIsyLw8F8LHfsgNAPQAHytm+F4B/A4AQwg7gKiKqDuAjABMBPFN2m1WrVpWatSkzMxNZWVlwu92+W8GsVis8Ho+v2OWV5b0TICkpCQaDAfn5+QDkrw2j0ehbNpvNMJvNyMvL8y1bLBbk5eVBCAGTyeS7u8Dj8cBoNCIxMdG3bDAYkJycDIfDAbfb7Vt2Op0oLCwEEcFms8HpdPq62NtsNrhcLt+y1WoNeE5eB1o6p1A/p1OnTqGgoEBT5xTK51RQUACHw6Gpcwrlc3I4HDAYDOWeU+fOnXH33Xfjvffew6233ooVK1YgOTk55s8p1M/JbrcruwNKCBGVB4A+AD7zW94KoHk569UEcEcF+8gAsKS89zIzMwUjRHZ2ttohxAzsQsIeJJV5KCgoEBdddJEAIMaMGROlqKKHy+US586dE0KUdgFguqjk+zuaVU9rAdQHACKyACgUQuwlIhMRpRa/XhPANUKIOcXL9Yr/eosK6QCWRzFmhmF0QlJSEv773//CZDLh7bffxhdffKF2SGFj8eLFMJvNGD58eEjbRy1RCNkuMZ2IJgF4CLKBGgAGAniKiOoC+A7Aw0T0MxHtA2AhoqsBbCGiaZAN2rPL27/RqJmBcKtEUlKS2iHEDOxCwh4kwXjo0KEDnn32WQDA3XffjRMnTkQ6rKhQvXp1AMDJkycBKL8movrtKoT4EsCXZV6bD2B+8WL7cjbbV8HrpfBvn9AzBoNmusZUGXYhYQ+SYD1MmDABX375JVatWoVx48ZhwYIFEY4s8tSqVQsAfH1FlF4TmrmCCgsL1Q4hJvA2qjHswgt7kATrwWAw4L333oPVasXChQuxcOHCCEcWeWrWrAmgpESh9JrQTKJgGIYJF02aNMELL7wAQM5l4R36Il7xlii8iUIpmkkUCQkJaocQE4TaoUaLsAsJe5Ao9fDPf/4TPXv2RHZ2Nv71r/ju55ucnAyz2QyHwwG73a7YhWYSBdfDSrhRvwR2IWEPEqUevFVQycnJ+OSTT7BkyZLIBBYFiKhUqUKxi0gEpQbcRiHh+ugS2IWEPUhC8dCsWTPMmDEDAHDvvffG9cCBr7zyCubPn4/U1FRuo2AYhgkn48aNQ/fu3XH8+HE88sgjaocTMrfccguGDBkCm82meFvNJApuo5CYzWa1Q4gZ2IWEPUhC9WAwGDBr1iyYTCbMnj0bP/zwQ5gjiz5KXXCi0Bj8pVACu5CwB0lVPLRu3RqTJ08GIBu5veMqxRM//PADnn/+efz000/6TRTx+MFFAu9gZAy78MIeJFX1MHnyZLRs2RI7duzw3TobTyxbtgyTJk3CqlWrFLvQTKJgGIaJJImJiZg5cyYA4Mknn8Tu3btVjkgZ3k53oTTIayZRcNWThKsZSmAXEvYgCYeH3r17484774TT6cSYMWPiagpm/9tjdVv1xIlCwp2rSmAXEvYgCZeHl156CbVq1cLKlSvx3//+Nyz7jAb+4z3ptsMdt1FIuD66BHYhYQ+ScHlIS0vDiy++CEAOIBgvfSv8x3viNgqdE09F4UjDLiTsQRJOD3fddReuvPJK5OTkxE3fCv+qJ6UuNJMoeAgPCc8dXgK7kLAHSTg9EBHefvttX9+K9evXh23fkaJWrVpISkpCYmKiYhea+Xbl8WwkXB9dAruQsAdJuD1ceOGFmDhxIgDZtyLWhxHKyMhAQUEBNm/erN82ilj/kKKFognTNQ67kLAHSSQ8TJkyBc2aNcO2bdvw6quvhn3/4cR/cjelLjSTKLgeVuLxeNQOIWZgFxL2IImEh6SkJLz55psAgOnTp+PAgQNhP0YkUOpCM/U1FbVRFBYW4tChQ3A4HFGOqOokJiaiQYMGiuoTuQquBHYhYQ+SSHn4v//7PwwZMgQLFizA+PHjsXTp0ogcJxzcdNNNWLt2LebPn48ePXoEvZ1mrqCKLoJDhw4hJSUFTZo0iat5tYUQOHnyJA4dOoSmTZsGvV1iYmIEo4ov2IWEPUgi6eHVV1/FV199hc8//xxLly7FDTfcELFjVYXTp0/j+PHjyM3NVbSdZqqeKmqjcDgcqFWrVlwlCaBkohGlJSGujy6BXUjYgySSHurVq4enn34aAPCvf/0rZvuueG+RPXbsmKLtNJMoArVRxFuS8BJK3FwfXQK7kLAHSaQ9jB07FpdeeikOHjyI6dOnR/RYoRLq3NmaSRTxmgzCDfcnKYFdSNiDJNIeEhISMHPmTBARXn31VWzdujWixwsFb+/ss2fPKtouqlcQEfUiokeKH+0qWOcpIjpERH8TUfPi18YR0SgimkZEqeVtx52KJMnJyWqHEDOwCwl7kETDQ6dOnTBu3DgUFRXh3nvvjbnSnLdEEbNtFESUCOB5AM8BeAXAG1SmGEBEdQGcEEI0EEI0EkLsJaKOADoLIWYDWAlgWnn7d7vdkT2BOCEe7+6KFOxCwh4k0fLw1FNPISMjAxs2bMC7774blWMGS/Xq1VG7dm3FCSyaJYquAHKExAkgGUCjMuv0BTCWiDYTUafi124E8Efx8y3Fy+cRa5m7LIsWLcLcuXMxffp0vP/++3j44Ydx+PBhPPPMM/jyyy/DdhxOmCWwCwl7kETLQ/Xq1X2d7yZNmoQTJ05E5bjBUK1aNWRnZyuuFovI7bFENB1AizIvDwXwsd+yA0A9AL4eKkKIDwB8QETXAfiMiC4GkAFgU5ltzuP7778v1U6RmZmJrKwseDwe3x1RkRqXv7CwEEIIEBEMBgM8Hk+p5a1bt+If//gHhBCYN28e7r33XjRt2hR16tRBRkYGCgsLUVhY6FtfCOFLfEII5Obm+kbHtVqtcLvdcDqdvmWPxwO73Q5AjqJbUFDgu8MjKSkJBoMB+fn5AOQwBkaj0bdsNpthNpt9d2mYzWZYLBbk5eVBCAGTyQSLxYKCggJ4PB4YjUYkJib6lg0GA5KTk+FwOOB2u33LTqfTd042mw1Op9N3DjabDS6XK+hz8lYZKD0nb/FaS+cUyueUm5uruXMK5XPKzc2N2jn17t0bvXr1wnfffYcHHngAr7/+ekxce4sXLwYA9OrVCzk5OQgaIURUHgD6APjMb3krgOYB1p8JoAOAJwHcX/xaDQAHy1t/ypQpojy2b9/uew4gIo9gmT9/vnjppZdKvfb++++LL774osJt/OMPBrvdrmh9LcMuJOxBEm0Pe/bsERaLRQAQK1eujOqxy6OoqEjUqVNHABCbNm3yvQ5guqjk+zuaVU9rAdQHACKyACgUsg3C5G2gJiL/ePIhq5wWQyYMAGgL4PPydh5MsbIyGaE+KmPLli3Yt28fVqxYga5du8LlcmHNmjWVbhcK3l8RDLvwwh4k0fbQvHlzTJ06FQAwZswY1T+HLVu24Pjx42jQoAFatChb4ROYqCUKIdslphPRJAAPARhV/NZAAE8VP19ORPOI6H4AM4UQhUKIXwFsJ6J7AFwFYGp5+4/lNopvvvkGy5cvR4sWLbBx40Z88skn6Ny5c0SOxYMjlsAuJOxBooaHzMxMXHDBBdi1axeef/75qB/fH+/QItddd53i9pqoDuEhhPgSwJdlXpsPYH7x82sr2O6FyEcXObxDEZfl8OHDWLVqFapVq4bLL78c6enpVT4W9ycpgV1I2INEDQ8WiwVvv/02evfujaeffhq33Xab4l/z4WLhwoUAgEGDBil2oZmeOPE4gXz9+vXx0Ucf4Y033ghLkgBkwxYjYRcS9iBRy0OvXr1wxx13wOl0Yty4caqMdP3nn39i+/btSE1NRa9evRS70EyiKCoqUjuEmEDtetBYgl1I2INETQ8vvvgiatSogW+++Qaffvpp1I+/aNEiAMCAAQNgNpsVu+BEoTG8t8gx7MILe5Co6SE9PR3PPfccAGDChAk4c+ZMVI+/YMECAMDgwYMBKHehmUTBMAwTy9x99934xz/+gWPHjmHKlClRO+7WrVuxdetWpKam4pprrglpH5pJFPHYRhEJuD66BHYhYQ8StT0YDAbMnDkTCQkJeOutt7By5cqoHPeDDz4AAAwdOtQ3Jwe3UegcrmYogV1I2IMkFjy0a9cOjz32GADgrrvuwunTpyN6PJfLhblz5wIAhg8fXup1JXCi0Bix8M8QK7ALCXuQxIqHKVOmoEuXLjh8+DDuvffeiN4FtXTpUuTk5ODiiy9Gp06dfK/rNlEwDMPEA0ajEXPnzoXVasX8+fPx3//+N2LHeuWVVwDInuFV6UeimUTB81FIrFar2iHEDOxCwh4kseShefPmeO211wDImfH+/PPPsB/jxx9/xI8//ogaNWqUqnYClLvQTKKI5SE8ogkPKV0Cu5CwB0mseRgxYgRuueUWnDt3DoMGDcK5c+fCun/vkCH//Oc/z0sMSl1oJlFwG4WEO1eVwC4k7EESax6ICO+++y4uuugi7NixAyNHjgxbe8XGjRuxZMkSJCUlYfz48ee9r9sOd8FCRBU+Zs2a5Vtv1qxZAddVSnkTF23fvh2LFy/Go48+qnhqQoZh4h+bzYbFixcjJSUFCxcuxFNPPVX5RpUghMDkyZMBAOPHj0fdunWrvE/NJIpYbqP4/fff0aNHD1x99dXYtGkT+vfvjxtvvBGLFy9G+/bt0aVLl7DNchdL9bBqwy4k7EESqx4uuOACzJ07FwaDAY8//niVp0/99NNPsWrVKqSmpmLSpEnlrqPbNopgi2yB5pUYPXq0b73Ro0dXaQ4Kf9q2bYvatWtjzZo16NOnD9LS0tCtWzdMnToVzZs3x7Fjx9CgQQNF+6wIbqspgV1I2IMklj0MGDAAb7zxBgDZpvDZZ5+FtJ+TJ0/i/vvvByDbKGrUqFHuerE8Z3ZEibWGKn8CTVzknaqyW7duYTmWdwpHhl14YQ+SWPcwZswYTJ06FR6PB0OGDMG8efMUbV9UVIRhw4bhxIkT6NGjB+6+++4K11XqIqrzUeiVb775BsnJyb6Ji3bt2oVbbrkFAPD1119jwoQJ+Ouvv9CsWTOVI2UYRk2eeOIJFBUV4dlnn8Xtt9+OgwcPYuLEiZW2iwoh8MADD+Drr79GWlqaryorbERqetBoPx599FFRHkrnnI4m//nPf0TXrl1Fnz59KpxTV2n8+fn54QhNE7ALCXuQxJOHGTNmCAACgBgwYIA4dOhQhes6HA4xatQoAUCYzeag5uf2d4EYmzObKcOIESOwbt06rFixAr1791Y7HIZhYoRJkyZh6dKlqF69Oj7//HNccMEFmDx5Mvbs2eNrI3W5XFiyZAk6deqE2bNnw2w2Y9GiRRH5LtFMouB5gSUFBQVqhxAzsAsJe5DEm4cBAwZg69atGDhwIPLz8zFjxgy0bNkSGRkZaNmyJVJTUzFw4EBs27YNTZs2xbp169CvX7+g9q3UhWYSBcMwjNZo3LgxFi9ejHXr1mHYsGGoWbMmTpw4gT179sBut6NNmzZ46aWXsH379lKD/oUbzTRmG42aOZUqkZSUpHYIMQO7kLAHSTx76Nq1K7p27QqPx4PDhw8jPz8fdevWRfXq1UPan1IXmvl2rcrIiFoirHc6xDnsQsIeJFrwYDAY0LBhw7DsR9H6VT6iAoioFxE9UvxoV8777xHRfiLaQ0THiWhk8et1iOgvIjpERK+Xt+9AbRTexp94I5S48/PzIxBJfMIuJOxBwh5KUOoiaiUKIkoE8DyAywGYAXxLRFcW354FIjICmCuEuLt4+RUAi4s3HwWgrRBC8SedmJiIkydPolatWnFV6hBC4OTJk76pCxmGYdQimlVPXQHkFCcGJxElA2gE4AAACCHcAL4DfEkjWQhxhoiSAPQEcB8RfQBgavG6pUhISCj3oA0aNMChQ4eQnZ0d/jOKMImJiYqH9rBYLBGKJv5gFxL2IGEPJSh1EZFEQUTTAbQo8/JQAB/7LTsA1ENxoihDXwBfA4AQwg7gKiKqDuAjABMBPFN2g9WrV5cqMWRmZiIrKwtutxspKSlISUmB1WqFx+PxdV9PTk4GUHKrWFJSEgwGg69YZrFYYDQafctmsxlmsxl5eXm+ZYvFgry8PAghYDKZYLFYUFBQAI/HA6PRiMTERN+ywWBAcnIyHA4H3G63b9npdKKwsBBEBJvNBqfT6ZuqsKioCHa73bdstVrhdrt9wwSXPSeDwYCCgoKYPiebzQaXyxX0OYX6OeXm5sLpdGrqnEL5nOx2OwoLCzV1TqF8Tk6nEyaTSVPnFOrnVFhYiJycHAQLRav+noj6ALhPCDGweHkrgEFCiL3lrPsWgAlCCGeZ1zMAzBRC3Fh2m4kTJwrvRB16JicnB2lpaWqHEROwCwl7kLCHEvxdEFGWEGJ6oPWj2Zi9FkB9ACAiC4BCIcReIjIRUap3JSKyFr/n9HvNW1RIB7C8vJ2vW7cuUnHHFZwsS2AXEvYgYQ8lKHURtRIFABDRdQDaAkgA8D8hxK9EdDOAHkKI+4rXGQbgiBDC215xNYAXIRu29wghyp2JnIhENM8lViGiuL3LK9ywCwl7kLCHEvxdBFOiiGqiiCScKCT8z1ACu5CwBwl7KEFpotBMhzsAvxLR52oHEQP0JKLVagcRI7ALCXuQsIcS/F1UemulZkoUDMMwTGSI/z7tDMMwTEThRMEwDMMERBOJgojGEdEoIprmf6utniCiakS0iIj+JKLPi28z1i1ElExEW4ioidqxqA0RXU1EQ8sbX00PEJGBiKYS0WAimkJEbdWOKZoQUWsi+oKIrixebkBEzxDRaCK6JZh9xH2iIKKOADoLIWYDWAlgmsohqUUPACMBXAggH8AIdcNRnTGQvf91DRGNANBECDFPCPGb2vGoRDsAtYUQiwAsBXC3yvFEFSHEnwDOAfD2R3sTwHtCiFkAhhJR3cr2EfeJAsCNAP4ofr6leFl3CCGWCSHOFt8jvAHAEbVjUgsiugPAfOg8URR/AUwDUERE7xf/qNIjOwB0I6LOADoDeFvleNTABfg6O/fyGxFjN4BrKttYC4kiA8Cp4ufe8aP0ThMAurxVmIh6ANgrhDiodiwxwADIjq3/gRwn7XMi0sL/vCKKR3kYA+BhAN0A7FM3IlWpBVm68BLUd6YWLprjAJKLn6cAOKFiLKpDRAMBvFLeCLs64T4AzxTfI94BwCdEVFvViNSjBoBcACge6SAZ8otCVxSXrIYLIYYA+BNAuXPa6IQcAP5Dxwb1namFRLEY8gsBkMOD6PKXNAAQ0TUAfhdC/E1E9dWORw2EEDcLIXoKIXpCVkXeKoSIvzHmw8NqAJcBABElANivUxedAXjHjnseQCsVY1EVIYQLwDq/mzxaoXik7kBoosMdEWUCOA05v8UrQojTKocUdYjoAQATAGRD/gDYKIS4V9WgVKa4VDFcCLFf5VBUg4gmQtZPOwD8JIT4VeWQog4RmQG8BOB7yEnTjgkhVqkbVfQgosaQUzwsB/AyZHX9GAB7AZwTQsyrdB9aSBQMwzBM5NBC1RPDMAwTQThRMAzDMAHhRMEwDMMEhBMFwzAMExBOFAzDMExAOFEwDMMwAeFEwTARhoiaqR0Dw1QFThQME0GIqDmAp9WOg2GqAicKhoksVwK4JNhx/xkmFuFEwTCRZRXkkBGfqh0Iw4QKJwqGYRgmIJwoGCayCAAGPc4DwWgHvngZJrIch5wXYrTagTBMqPDosQzDMExAuETBMAzDBIQTBcMwDBMQThQMwzBMQDhRMAzDMAHhRMEwDMMEhBMFwzAMExBOFAzDMExAOFEwDMMwAeFEwTAMwwSEEwXDMAwTEE4UDMMwTEA4UTAMwzAB4UTBMAzDBIQTBcMwDBMQThQMwzBMQDhRMAzDMAHhRMEwDMMEhBMFwzAMExBOFAzDMExAjGoHEC5uv/120aJFC7XDUB273Y6kpCS1w4gJ2IWEPUjYQwn+LrKyst4TQtwTaP24SBREVA3A+wAuBrALwG1CiHz/dRo3bozp06erEF1sUVBQgOTkZLXDiAnYhYQ9SNhDCf4usrKyDlW2frxUPfUAMBLAhQDyAYxQNxyGYRj9EBeJQgixTAhxVgghAGwAcKTsOoWFhdEPLAYpKChQO4SYgV1I2IOEPZSg1EVcVD2VoQmAN8u++P3334OIfMuZmZnIysqC2+2G0+kEAFitVng8HtjtdgDwFb280pKSkmAwGJCfL2u1LBYLjEajb9lsNsNsNiMvL8+3bLFYkJeXByEETCYTLBYLCgoK4PF4YDQakZiY6Fs2GAxITk6Gw+GA2+32LTudThQWFoKIYLPZ4HQ64XK5AAA2mw0ul8u3bLVaA56Tw+FAQUGBps4p1M/p7NmzmjunUD6ns2fPau6cQvmczp49q7lzCvVzcjgcyMnJQbCQ/JEeHxDRQAC/CCH+Lvveo48+Kp555hkVooot8vPzYbVa1Q4jJmAXEvYgYQ8l+LsgoiwhxPRA68dF1RMAENE1AH4XQvxNRPXLeV+FqGIPgyFuPtKIwy4k7EHCHkpQ6iIuzBHRAwBmA/iEiDYDeKzsOtxGIfEWgRl24YU9SNhDCUpdxEUbhRDiVQCvqhwGwzCMLomLRBEMCQkJaocQE1gsFrVDiBnYhYQ9SCLpobCwEAcOHIDT6UR6ejrS0tJiujpcqQvNJAquf5QYjZr5SKsMu5CwB0m4PZw9exYfffQRFixYgE2bNiEjIwOJiYk4duwYrFYrbrrpJjzwwANo0qRJWI8bDpS60My3K7dRSLgetgR2IWEPknB5sNvteOKJJ9C8eXOsW7cOEydOxMmTJ7Fv3z7s2LEDp06dwooVK5CUlIRLL70U06ZNg9vtDsuxw4VSF5pJFAzDMJHmu+++Q7t27fDbb7/hp59+wieffIJ+/fqVGhqEiNC6dWs8++yz2LZtGzZs2ICrr74aubm5KkZeNTSTKLiNQmI2m9UOIWZgFxL2IKmKB4/HgyeeeALDhg3DK6+8goULF6J58+aVble3bl18+eWXuPDCC9G7d++YSRZKXXCi0Bj8pVACu5CwB0moHvLy8tC/f3+sXLkSv/zyC/r166do+4SEBLz55pvo1KkTbr755pioJtdtovB2X9c73qEDGHbhhT1IQvFw4sQJ9OrVC/Xq1cPKlStRt27dkI5NRHjjjTcAANOmTQtpH+FEqQvNJAqGYZhwcuDAAVxxxRW49tprMWvWrCrfNWU0GvHRRx/hgw8+wOrVq8MTZJTQTKLgqicJVzOUwC4k7EGixMOhQ4fQu3dv3HfffXjiiSfC1iciPT0ds2bNwj333AOHwxGWfYaCbqueOFFIuHNVCexCwh4kwXo4cuQIevfujbFjx+L+++8Pexz9+vVDmzZt8PLLL4d938Gi9JrQTKLgNgoJ10eXwC4k7EESjIecnBz06dMHI0aMwEMPPRSxWF555RW8/PLLOHr0aMSOEQhuo9A58TRsfKRhFxL2IKnMQ0FBAfr3748bb7wRkydPjmgsTZs2xZ133okZM2ZE9DgVofSa0Eyi4CE8JCaTSe0QYgZ2IWEPkkAeioqKMGzYMDRv3hzRmtdm4sSJmDNnDg4fPhyV4/mj9JrQzLcrj2cj4froEtiFhD1IKvIghMCECRNw9uxZ/Oc//4naYH4ZGRkYMWIEXnzxxagczx/dtlHEQieWWIDnBS6BXUjYg6QiD6+88gpWrVqFxYsXR/0Osfvvvx8ffvghzp07F9XjKr0mNJMouB5W4vF41A4hZmAXEvYgKc/Dl19+iRdffBFffvklUlNTox5To0aN0KdPH3z44YdRPa7Sa0IziYLbKCRcBVcCu5CwB0lZD7t27cLw4cMxf/58NGrUSKWogPHjx+P111+PakLX7TDj/M8gSUxMVDuEmIFdSNiDxN9Dbm4ubrjhBjz11FPo1q2bilEB3bp1g9lsxpo1a6J2TKXXhGYSBbdRSLg+ugR2IWEPEq8Hj8eD22+/HT179sTo0aNVjkqOA3XXXXdFtfqJ2yh0DtdHl8AuJOxB4vUwbdo0nD59Gv/+979VjqiEYcOG4bPPPovaJFO6baOI5flpowm31ZTALiTsQWIwGLBw4UJ89NFHWLhwYUyNgVW3bl1cccUVWLx4cVSOp/Sa0MwVxJ2KJP4zbekddiFhD5K9e/dizJgxWLx4MerUqaN2OOdxxx13YN68eVE5ltJrQjOJItbmpFULNUekjDXYhYQ9yDGchgwZgn//+9+49NJL1Q6nXK6//nqsW7cOZ86cifixlF4TmkkUXA8r4YRZAruQ6N2D2+3GLbfcgv79+2Po0KFqh1MhKSkp6NmzJ5YtWxbxYym9JjSTKLiNQsL10SWwC4nePTz88MMwmUx4/PHH1Q6lUgYPHhyVdgpuo9A5XB9dAruQ6NnDhx9+iGXLluHjjz9GSkqK2uFUindu7kjf/cRtFDrH6XSqHULMwC4kevWwadMmPPzww1i6dClq1KgRFx5q1qyJyy67DN9++21Ej6PURVwlCiJqTURfENGVZd/jNgoJdzwsgV1I9Ojh+PHjGDx4MGbNmoWLL74YQPx4uO666/DVV19F9BhKXcRVohBC/AngHABukKgAbqspgV1I9ObB5XJh8ODBGDFiBAYOHOh7PV48/N///R+++uqriHYiVuoirhJFMeXOeRpLnWfUxGazqR1CzMAuJHrzcP/996NWrVqYNm1aqdfjxcOFF14Ij8eDnTt3RuwYSl1oZiS9FStW4Omnn/YtZ2ZmIisrC26321cfZ7Va4fF4YLfbAZQ06HjHPUlKSoLBYPA1JFksFhiNRt+y2WyG2Wz2zTdrNpthsViQl5cHIQRMJhMsFgsKCgrg8XhgNBqRmJjoWzYYDEhOTobD4YDb7fYtO51OFBYWgohgs9ngdDp9c4DbbDa4XC7fstVqDXhORUVFSElJ0dQ5hfo5HT9+HFarVVPnFMrndO7cOaSkpGjqnCr6nD7++GOsWrUKX331FU6fPl3qnPLz81GnTp24OKfevXtj4cKFGDt2bEQ+p3PnziEhIQFBI4SIqweADwD0LPt6ZmamYITIzs5WO4SYgV1I9OJh3bp1onbt2mLnzp3lvh9PHhYtWiSuueaaiO3f3wWA6aKS7914rHpiGIYpxcGDBzFkyBB88MEHaNWqldrhVJk+ffpg3bp1MXOnVlwlCiJqDKAVgCuIKMn/PW6jkMRLPWw0YBcSrXs4d+4c+vXrhwcffBDXXXddhevFk4fq1aujdevW2LhxY0T2r9RFXCUKIcQBIURXIcTTQgi7/3tFRUVqhRVTeOspGXbhRcseioqKcNttt+Hyyy/Hgw8+GHDdePNw5ZVX4vvvv4/IvpW6iKtEEQhOFJJ4+2eIJOxComUPDz/8MOx2O956661Kb/mMNw89e/bE6tWrI7Jv3SYKhmH0xcyZM/G///0PCxcu1OQQPt27d8dPP/0UEwlOM4lCixdKKFitVrVDiBnYhUSLHr755htMnz4dy5YtQ40aNYLaJt48pKamokWLFvj555/Dvm+lLjSTKHgIDwmPeVUCu5BozcPPP/+MYcOGYcGCBWjRokXQ28Wjh0i1U+h2mHFuo5DEyu10sQC7kGjJw+7du9G/f3/Mnj0b3bt3V7RtPHro0aMHfvjhh7DvV9ODAjIMo1+OHj2Kvn374oknnsCNN96odjhRoUuXLvjpp58iOu5TMGgmUXAbhSTe6mEjCbuQaMHD2bNnce2112LkyJEYNWpUSPuIRw/16tWDzWbDnj17wrpf3bZRqJ1xYwVuqymBXUji3UN+fj769++Pbt26YcqUKSHvJ149dOnSBRs2bAjrPpW60EyiiMeGqkjgHfiLYRde4tlDQUEB+vXrh5YtW+K1116r0lDh8eqhc+fOYU8USl1oJlEwDKMtHA4HbrjhBjRs2BCzZs3S7dzfkShRKEUz5rmNQqLn+ZHLwi4k8ejB6XRi0KBBSEtLw/vvv69sSOwKiEcPAHDJJZdgx44dvqHOw4Fu58xmGEYb5OXloV+/frDZbJgzZ05YkkQ8k5SUhDZt2uCXX35RLQbNJIp4mQ830oTzV0e8wy4k8eTh9OnTuPrqq9GoUSN8/PHHMBrDN7daPHkoS+fOncM6kqxSF5pJFAzDxDfHjh1Dz5490aVLF8yePVv3JQl/OnbsiF9//VW142smUYTzl0c8k5SUVPlKOoFdSOLBw86dO9GtWzcMGjQIL7/8ckQaruPBQ0VccsklYU0USl1oJlFU5bY5LaHXO0PKg11IYt3D6tWr0aNHDzzyyCOYNm1axP6XY91DIC666CLs27fPNzd3VVHqIn7NlYHbKCThupC0ALuQxLKHjz76CLfccgvmzZuHe+65J6LHimUPlWE2m3HRRRfht99+C8v+lLrg+hqGYaJOYWEhJk2ahKVLl2L16tW48MIL1Q4p5vFWP/3jH/+I+rE1kyi44UtisVjUDiFmYBeSWPNw5MgR3HzzzahevTo2bdqEmjVrRuW4seZBKR07dgzbLbJKXWim6ime6x/DCTfql8AuJLHkYeXKlejUqRP+7//+D1988UXUkgQQWx5CIZx3Pil1Ed/m/OA2Ckl+fn5c390RTtiFJBY82O12TJ48GQsXLsSHH36Iq6++OuoxxIKHqtC2bVvs2LEDLpcLZrO5SvtS6oJ/hjMME1F++eUXXHrppTh69Ci2bt2qSpLQAsnJyWjWrBm2b98e9WNrJlFwG4Wkqr80tAS7kKjl4dy5c3jwwQdx7bXXYurUqfjkk09Qq1YtVWIBtHE9tGvXDr///nuV96PUBScKjaGFf4ZwwS4k0fYghMCCBQtw4YUX4vTp0/jjjz8wdOhQ1fs6aeF6aNOmDbZt21bl/Sh1oZk2CpfLpXYIMUFeXh4SExPVDiMmYBeSaHr48ccfMWnSJJw6dQoff/yx4nmtI4kWroc2bdrgnXfeqfJ+lLrQTImCYRj1+OOPPzBw4EDcfPPNGDFiBLZu3RpTSUIrhKtEoRTNJAquepJooXgdLtiFJJIe1q5diwEDBqB3797o2rUrdu3ahREjRsTk/6MWrocmTZrg5MmTyM3NrdJ+uI1C58R7p6Jwwi4k4fZgt9sxd+5cdO3aFXfddReuvfZa7N+/H5mZmTF9+6kWrgeDwYCLLroIf/zxR5X2o+kOd0Q0johGEdE0Ikr1f4/bKCR5eXlqhxAzsAtJODwIIfDzzz/jX//6Fxo0aIC5c+fioYcewq5duzBmzJiYThBetHI9tGnTpsp3Pil1ETeN2UTUEUBnIcSdRNQNwDQAE1QOK+YQQqgdQszALiShevB4PNi0aRMWLlyIhQsXwmg0YujQodi8eTMaN24c5igjj1auh3C0Uyh1ETeJAsCNALzlrS0A5sAvUaxZswaPPPJIuRsGui2vovdC2SYW9udwOAL+ulM7vkDvGQwGmEwmmM1mmM1mxc8tFguSkpKQnJzse48Jfj55IQT++usvrFy5EitWrMCqVatQp04dDB48GEuXLkXbtm1Vv8W1KmjlemjTpg2WLVtWpX0odRFPiSIDwKbi5w4A9fzf3L17N3766SffcqNGjdCsWTMIIeDxeADIdgz/Ze/4UEVFRb73/ZcNBgOIyLdMRDAYDHC73b73DQYDioqKIITwve/xeEote7f3HsP7vncf/stlYyxv/UDn5L9P/3P07qPsOVV0jv7LRHTeOZWNuaIYvesHOifvORQVFfmOW1RUhMLCQt96RUVFcLvd8Hg88Hg8vnWLiop8rxUWFsLtdvvWMxqNMJlMvkdCQgKMRqMviXiXTSYTEhMTkZCQ4Es6ycnJMBqNMBqNsFgsvuWEhATfev7xll02mUy+OL3LQgjfteP9R/UOPWM0GkFEvuWEhAQYDIZSywkJCb4q1rLLBoMBRqMRhYWFvs/Zu+x2u33n6XUjhEB+fj6OHz+Oo0eP4tixYzh27BgSEhLQqFEjNG7cGEOHDoXVaoXD4cDcuXNj6py816L/OZVdLvu5FBUVISkpKWY/p2DPyel0YsOGDZg4cWLI157b7VZUqqB4KY4R0ZMAcoQQ/yaiGgB+E0I09L4/adIk8dxzz6kXYIxw6tSpqA60Fqu43W4cPnwYiYmJsNvtKCgo8P31f26325Gfn4+8vDycO3eu0kdubi5MJhNSUlJKPapVq4aUlBTYbLZSj/JeK/tITk4O26CWQgg4HA7k5ubi1KlTOHr0KHbt2oXc3FwcPXoU+/btw86dO7Fv3z5kZGSgdevW6NSpEy677DJ06tQJ9evXD0scsYhW/jeEEKhZsyZ27tyJ9PT0kPbh74KIsoQQ0wOtH08lisUAxhc/bwvgc/834yXhRRr/koWeMRqNsFqtSEtLC+t+hRCw2+0VJpK8vDzf4+zZszh8+HCp18quk5eXB7vdDqvVCqvV6ivtlH2YzWYYjUZ4PB5fqcntdvtKWvn5+cjNzUVubi4SEhJQrVo1pKamol69eqhVqxaaNm2Khg0b4oorrsAFF1yAFi1axEUDdDjRyv8GEaF169ZVShRKXcRNohBC/EpE24noHgCNAEz1f5+HGZfE+1DK4SQSLogIycnJSE5ORp06dcKyz6KiIhQUFCAvLw8ulwuFhYXnPVwuF9xut6+6zFsN5v1rtVp9pZqytz6eOXMGqampYYk1ntHS/4Y3UYTaqVHTw4wLIV6o6D0tXQRVId6HKAgn8eIiISHBV4UVCeLFQ6TRkocLLrgAf/75Z8jbK3WhmZ/hPB+FpKCgQO0QYgZ2IWEPEi15uOCCC7Bz586Qt1fqQjOJgtsoJFqphw0H7ELCHiRa8tC6desqlSiUutBMoojn+7vDCbfVlMAuJOxBoiUPzZs3x8GDB+F0OkPaXqkLzZjTSmeaqpKcnKx2CDEDu5CwB4mWPJjNZjRq1Ah79+4NaXulLjSTKLydSfSOw+FQO4SYgV1I2INEax6q0k6h1IVmEoWW6h+rAifMEtiFhD1ItOahKolCqQvNJApuo5BoqR62qrALCXuQaM1DVRq0uY1C52ipHraqsAsJe5BozUNVShTcRqFzQr0LQouwCwl7kGjNg7dEEUrXAKUuNJMouI1Cwh0PS2AXEvYg0ZqHtLQ0CCFw8uRJxdsqdaGZRMFIuK2mBHYhYQ8SrXkgIjRv3jykW2SVutBMotDCxOnhwGazqR1CzMAuJOxBokUPoSYKpS40kyj8JwfSM1qrh60K7ELCHiRa9NCiRYuQEoVu2yg4UUi8M2kx7MILe5Bo0UPz5s2xZ88exdspdaGZRMEwDKM3Qq16UopmEgW3UUi0WA8bKuxCwh4kWvTAbRQK4aoniRaL16HCLiTsQaJFD/Xr18eZM2eQn5+vaDvdVj1xopBo8Z8hVNiFhD1ItOjBYDCgadOm+OuvvxRtp9tEwTAMo0dCbdBWgmYSBY/1JLFarWqHEDOwCwl7kGjVQyi3yCp1oZlEwUN4SHjMqxLYhYQ9SLTqIZQGbd0OM85tFBItdioKFXYhYQ8SrXoIJVHotsMdwzCMHuE2CgVwG4VEq/WwocAuJOxBolUPTZo0weHDhxWNCKvbNopQxmTXItxWUwK7kLAHiVY9mM1mpKen4/Dhw0Fvo9SFZhKFVhuqlGK329UOIWZgFxL2INGyhyZNmuDAgQNBr6/UhWYSBcMwjF5p0qQJ9u/fH7H9ayZRcBuFRGvzAlcFdiFhDxIte2jcuLGiRKHJObOJKJ2IXiSiyWrHwjAME2sorXpSSlwkCiHECQC7AFgqWkdr8+GGSkFBgdohxAzsQsIeJFr2oLTqSamLuEgUxWhvRC+GYZgwEOk2CmPE9lwORJQhhDgW4P3pAFqUfV0IcXtl+16zZk2pCcMzMzORlZUFt9vt64VotVrh8Xh8Lf7eejpvdk1KSoLBYPAN2WuxWGA0Gn3LZrMZZrMZeXl5vmWLxYK8vDwIIWAymWCxWFBQUACPxwOj0YjExETfssFgQHJyMhwOB9xut2/Z6XSisLAQRASbzQan0+kb3dFms8HlcvmWrVZrwHPyno+WzinUz8lutyMnJ0dT5xTK52S323Hq1ClNnVMon5PdbofD4dDUOXk/p+rVq+Pw4cM4fvw4qlWrFtQ55eTkIFgokv0PiOh9AN4DGACkCyGuC3FfwwE0EUJML+/9KVOmiKeffjqUXWsKu92OpKQktcOICdiFhD1ItO6hfv362LBhAxo2bFjpuv4uiCirou9VL5GuetoN4MPixxcAvo/UgbiNQqJ0AhMtwy4k7EGidQ9Kqp+UuohoohBCPCOE+L74sQjApaHsh4hSAVwGoC0R1QlnjAzDMFqgcePGEbvzKaJtFET0HUqqnqwAKmyfCIQQ4gyAcYHWSUhICGXXmsNiqfDGMN3BLiTsQaJ1D0pKFEpdRLoxey6AFcXPXUKIo5E6kMEQTzdwRQ6jMar3J8Q07ELCHiRa99CkSRNs3LgxqHWVuoh01dN7QogDxY+IJQmA2yi8aL0eVgnsQsIeJFr3oKTTXUy1UTAMwzDRQekwHkrQTKLgNgqJ2WxWO4SYgV1I2INE6x4aNWqEgwcPBjWEuFIXnCg0htb/GZTALiTsQaJ1D0lJSahWrRpOnDhR6bq6TRTeXop6x9sjlGEXXtiDRA8eGjZsiIMHD1a6nlIXmkkUDMMweqdBgwY4dOhQ2PermUTBVU8SrRevlcAuJOxBogcPwZYodFv1xIlCovVORUpgFxL2INGDh2BLFEpdaCZRcBuFRA/1sMHCLiTsQaIHD9xGwQRFJEcDjjfYhYQ9SPTgIdgShVIXmkkUPISHhOcOL4FdSNiDRA8egi1RKHWhmW9XrY/jEix6qIcNFnYhYQ8SPXioX78+jh49WmmnO922UfBYTxItzwusFHYhYQ8SPXiwWCxITU3F8ePHA66n5TmzA6KH+sdgCKb7vl5gFxL2INGLh2DaKZS60Eyi4DYKCVfBlcAuJOxBohcPwbRTxNQw49FELxdBZSQmJqodQszALiTsQaIXDw0aNKg0USh1oZlEwW0UEj3UwwYLu5CwB4lePDRs2LDSqiduo9A5eqmHDQZ2IWEPEr14CKZEods2CiJSO4SYgNtqSmAXEvYg0YuHYEoUSl1oxpweOtMEQ3JystohxAzsQsIeJHrxEEyJQqkLzSQKt9utdggxgcPhUDuEmIFdSNiDRC8evJ3uioqKKlxHqQvNJAq91D9WBifMEtiFhD1I9OLBYrGgRo0aATvdKXWhmUTBbRQSvdTDBgO7kLAHiZ481KtXD0ePHq3wfW6j0Dl6qYcNBnYhYQ8SPXmoV68ejhw5UuH73Eahc5xOp9ohxAzsQsIeJHryUFmiUOpCM4mC2ygk3PGwBHYhYQ8SPXmorOpJqYu4SBRE9DwR/UpEm4noQrXjiWW4raYEdiFhDxI9eahbt27AEoVSFzGfKIioPYCPhBCXAJgN4Mny1tPDxOnBYLPZ1A4hZmAXEvYg0ZOHyqqelLqI+UQhhNgqhNhWvLgeQLlnH+ieYT2hp3rYymAXEvYg0ZOHcLdRxNSQq0Q0HUCLsq8LIW4vfno1gJfL23bVqlWlilOZmZnIysqC2+32SbFarfB4PLDb7QBKWv69A2QlJSXBYDAgPz8fgLwf2Wg0+pbNZjPMZrNvYnKz2QyLxYK8vDwIIWAymWCxWFBQUACPxwOj0YjExETfssFgQHJyMhwOB9xut2/Z6XSisLAQRASbzQan0wmXywVAZn6Xy+VbtlqtAc/J4XDAaDRq6pxC/Zyys7Phcrk0dU6hfE6nT5+G2+3W1DmF8jmdPXsWZrNZU+dU0edUrVo1HD58GDk5OeWeU25urm/9YKB4GUyvuAqqhhBidXnvT5w4UTz//PPRDSoGycnJQVpamtphxATsQsIeJHry4P1hUFBQUO4UDP4uiChLCDE90P5ivuoJAIioDYDqQojVRFSNiM6rYOM2Come6mErg11I2INETx6MRiNq1apVYe9szbVREFEPAMsAvExEPwNYCSC/7HrcRiFRUpzUOuxCwh4kevMQqJ1CqYuYaqMoDyHEGgBNKluPE4VEb/8MgWAXEvYg0ZuHQH0plLqI+RIFwzAMo5zK7nxSgmYSBY/1JLFarWqHEDOwCwl7kOjNQ6BEodSFZhIFD+Eh4TGvSmAXEvYg0ZuHunXrVlj1pNthxrmNQqKnTkWVwS4k7EGiNw+BShS6HRSQYRiGKYHbKMqB2ygkequHDQS7kLAHid48cBtFOcRLD/NIw201JbALCXuQ6M1D7dq1cfr06XKHFFfqQjOJQm8NVRXhHc+FYRde2INEbx4SEhJQu3btcntnK3WhmUTBMAzDlCZc7RSaSRTcRiHR07zAlcEuJOxBokcPGRkZOHbs2Hmv63bObIZhGKY0FSUKpWgmUehpPtxAeMfNZ9iFF/Yg0aOHOnXqlNtGodSFZhIFwzAMU5qKEoVSNJMoypucQ48kJSWpHULMwC4k7EGiRw8VVT0pdaGZROE/DaqeMRg085FWGXYhYQ8SPXqoqESh1IVmzHEbhcQ7dy/DLrywB4kePVSUKJS60EyiYBiGYUrDdz2VISEhQe0QYgKLxaJ2CDEDu5CwB4kePVSrVg0ul+u8u5yUutBMotBj/WN5cKN+CexCwh4kevRAROVWPyl1oZlvV26jkOixHrYi2IWEPUj06iEjI+O8RMFtFAzDMIyPcPSl0Eyi4DYKidlsVjuEmIFdSNiDRK8e6tSpc16DtlIXnCg0hl7/GcqDXUjYg0SvHsqretJtonC5XGqHEBPk5eWpHULMwC4k7EGiVw/lVT0pdaGZRMEwDMOcT3lVT0rRTKLgqieJXovX5cEuJOxBolcPXPXkBycKiR47FVUEu5CwB4lePZRXotBthztuo5DotR62PNiFhD1I9OpBN20URHQPEX1JRFuJqJXa8cQyQgi1Q4gZ2IWEPUj06qFatWpwu92lOtkpdRHziYKIzADWCCGuA7AAwDXlrcdDeEh47vAS2IWEPUj06qG8YTyUuoj5b1chhEsIsat40Q1gXnnr6XEcl/LQaz1sebALCXuQ6NlD2USh1EVMfbsS0XQALcp56w4A9wAYD+AYgA/KrvDtt9/i6aef9i1nZmYiKysLbrcbTqcTAGC1WuHxeGC32wEAycnJAErmj01KSoLBYPAV0SwWC4xGo2/ZbDbDbDb76vfMZjMsFgvy8vIghIDJZILFYkFBQQE8Hg+MRiMSExN9ywaDAcnJyXA4HHC73b5lp9OJwsJCEBFsNhucTqevzcVms8HlcvmWrVZrwHNyuVxITU3V1DmF+jkdOXIE1apV09Q5hfI5nTlzBqmpqZo6p1A+p9zcXNSrV09T5xTs51S7dm3s3r0bLVu2RHJyMs6cOaPozieKp3o7ImoE4AshRPuy702cOFE8//zzKkQVW+Tk5CAtLU3tMGICdiFhDxI9exg9ejQ6duyIe++9F0BpF0SUJYSYHmj7mK96KsNRAJvKe4PbKCRcBVcCu5CwB4mePaSnp+PEiRO+ZaUuYt4cEaUDWATgMwC5AKaWt56eLwJ/EhMT1Q4hZmAXEvYg0bOH2rVrY8+ePb5lpS5i/ttVCHECQPfK1uP5KCQFBQW6/ofwh11I2INEzx7S09Px448/+paVutBMfU08tbVEEo/Ho3YIMQO7kLAHiZ491K5dG9nZ2b5lpS40kyiISO0QYgJuqymBXUjYg0TPHmrXrl2qjUKpC82Y02tnmrJ4b+dj2IUX9iDRs4eyJQqlLjSTKNxut9ohxAQOh0PtEGIGdiFhDxI9e0hLS8PJkyd9VU5KXWgmUei5/tEfTpglsAsJe5Do2YPZbIbNZsOZM2cAKHehmUTBbRQSPdfDloVdSNiDRO8e/KufuI1C5+i5HrYs7ELCHiR69+DfoM1tFDrHO74Lwy68sAeJ3j34lyiUutBMouA2Cgl3PCyBXUjYg0TvHvwThVIXmkkUjITbakpgFxL2ING7h/T0dF+iUOpCM4lCrxOnl8Vms6kdQszALiTsQaJ3D/4lCqUuNJMoioqK1A4hJtB7Paw/7ELCHiR69+DfmK3bNgpOFBLv5CUMu/DCHiR69+BfolDqQjOJgmEYhqkY/zYKpWgmUXAbhUTv9bD+sAsJe5Do3QO3UYCrnrzovXjtD7uQsAeJ3j2kpaUhOzsbQgjFLmJ+4qJg4UQh0fs/gz/sQsIeJHr3YLFYkJycjDNnzij+vtRMiYJhGIYJTNnhxoNFM4mCx3qSWK1WtUOIGdiFhD1I2ENJg7ZSF5pJFDyEh4THvCqBXUjYg4Q9lJQodDvMOLdRSPTeqcgfdiFhDxL2UNLpTrcd7hiGYZjAcBsFt1EA4HpYf9iFhD1I2ENJotBtG4UQQu0QYgJuqymBXUjYg4Q9lDRmK3WhmUTBDVUSu92udggxA7uQsAcJeygpUSh1oZlEwTAMwwQmLS0NOTk5irfTTKLgNgqJ3ucF9oddSNiDhD2UJArdzpnNMAzDBEYXJQoimkVEw8t7T+/z4XopKChQO4SYgV1I2IOEPchSlRBCcbKIm0RBRLcA4PvbGIZhQoSIUKtWLZw6dUrRdnGRKIioM4AcALsrWsdo1MxAuFUiKSlJ7RBiBnYhYQ8S9iCZNm0aatSooWibmPp2JaLpAFqUebkugI+EEB8SUfeKtv3uu+9ARL7lzMxMZGVlwe12+7qrW61WeDwe361h3gYdb5E0KSkJBoMB+fn5AOSwvEaj0bdsNpthNpuRl5fnW7ZYLMjLy4MQAiaTCRaLBQUFBfB4PDAajUhMTPQtGwwGJCcnw+FwwO12+5adTicKCwtBRLDZbHA6nb4hkW02G1wul2/ZarUGPCeDwQAi0tQ5hfo5nTlzBna7XVPnFMrnZLfb4XQ6NXVOoXxOTqcTCQkJmjqnUD6nYcOGweFwKKp+oljvqEZEkwH0LV5sAsAB4F9CiG/915s4caJ4/vnnoxxd7JGTk4O0tDS1w4gJ2IWEPUjYQwn+LogoSwgxPdD6MVWiKA8hxLMAngV8JY79ZZMEwzAMEznioo0iGBISEtQOISawWCxqhxAzsAsJe5CwhxKUuoj5EoU/gYpHBoNmcl6V4Eb9EtiFhD1I2EMJSl1o5tuV+1FIvI1qDLvwwh4k7KEEpS40kyjWrVundggxATfol8AuJOxBwh5KUOoi5u96ChYiElo5l6pARDzkejHsQsIeJOyhBH8Xwdz1pJkSBcMwDBMZtNS68ysRfa52EDFATyJarXYQMQK7kLAHCXsowd9Fg8pW1kzVE8MwDBMZuOqJYRiGCQgnCoZhGCYgnCgYhmGYgGgiURDROCIaRUTTiChV7XjUgIiqEdEiIvqTiD4nIl3P3UFEyUS0hYiaqB2L2hDR1UQ0lIjaqR2LGhCRgYimEtFgIppCRG3VjimaEFFrIvqCiK4sXm5ARM8Q0ejieX4qJe4TBRF1BNBZCDEbwEoA01QOSS16ABgJ4EIA+QBGqBuO6oyBHGlY1xDRCABNhBDzhBC/qR2PSrQDUFsIsQjAUgB3qxxPVBFC/AngHADvPAxvAnhPCDELwFAiqlvZPuI+UQC4EcAfxc+3FC/rDiHEMiHE2eJehxsAHFE7JrUgojsAzIfOE0XxF8A0AEVE9H7xjyo9sgNAt+IJ0DoDeFvleNTABQBEZAHQSwixt/j13QCuqWxjLSSKDADeef0cAOqpGEus0ASALvuUEFEPAHuFEAfVjiUGGADgf0KI/wD4CMDnRKSF/3lFCCGckCXMhwF0A7BP3YhUpRZk6cJLUN+ZWrhojgNILn6eAuCEirGoDhENBPCKEMKtdiwqcR+AZ4o7E3UA8AkR1VY1IvWoASAXAIQQ30H+n9RSNSIVKC5ZDRdCDAHwJ4DXVQ5JTXIA+I8xHtR3phYSxWLILwQAaAud/pIGACK6BsDvQoi/iai+2vGogRDiZiFETyFET8iqyFuFENnqRqUaqwFcBgBElAA56ZceXXQG4Cx+/jyAVirGoipCCBeAdX43ebQC8HVl22miZzYRZQI4DaAR5K/p0yqHFHWI6AEAEwBkQ/4A2CiEuFfVoFSmuFQxXAixX+VQVIOIJkLWTzsA/CSE+FXlkKIOEZkBvATgewBmAMeEEKvUjSp6EFFjAB8DWA7gZcjq+jEA9gI4J4SYV+k+tJAoGIZhmMihhaonhmEYJoJwomAYhmECwomCYRiGCQgnCoZhGCYgnCgYhmGYgHCiYBiGYQLCiYJhIgwRNVM7BoapCpwoGCaCEFFzAE+rHQfDVAVOFAwTWa4EcEmw4/4zTCzCiYJhIssqyCEjPlU7EIYJFU4UDMMwTEA4UTBMZBEADHqcB4LRDnzxMkxkOQ45L8RotQNhmFDh0WMZhmGYgHCJgmEYhgkIJwqGYRgmIJwoGIZhmIBwomAYhmECwomCYRiGCQgnCoZhGCYgnCgYhmGYgHCiYBiGYQLy/ynNnOpxC3UJAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x662.4 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(2, 1, figsize=(6, 9.2))\n",
    "\n",
    "u=np.array(list(map(lambda tau:reach(0,x0,tf,xf,tau), trange)))  # 理論式\n",
    "_, _, x = lsim(sys, U=u, X0=x0, T=trange)     # シミュレーション\n",
    "ax[0].plot(trange, x[:,0], lw=2, ls = '-',  color='k', label='$x_1$')                    \n",
    "ax[0].plot(trange, x[:,1], lw=2, ls = '--', color='k', label='$x_2$')           \n",
    "ax[0].axhline(0, color=\"k\", linewidth=0.5)  \n",
    "ax[0].grid(ls = ':')\n",
    "plot_set(ax[0], 't', 'x', 'best')\n",
    "ax[0].set_xlim(0,tf+0.1)\n",
    "\n",
    "ax[1].plot(trange, u, lw=1, color='k')\n",
    "ax[1].axhline(0, color=\"k\", lw=0.5)  \n",
    "ax[1].grid(ls = ':')\n",
    "plot_set(ax[1], 't', 'u')\n",
    "ax[1].set_xlim(0,tf+0.1)\n",
    "\n",
    "if (is_savefig):\n",
    "    fig.savefig(figpath+\"ans/ch4_1_time_response.pdf\", transparent=True, bbox_inches=\"tight\", pad_inches=0.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 【３】 Kalmanの正準分解"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [],
   "source": [
    "import scipy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0 -1  0]\n",
      " [ 1 -1  0]\n",
      " [ 1  0 -2]]\n"
     ]
    }
   ],
   "source": [
    "A = np.matrix([[0,-1,0],[1,-1,0],[1,0,-2]])\n",
    "B = np.matrix([[0],[0],[1]])\n",
    "C = np.matrix([1,0,0])\n",
    "D = np.matrix([0])\n",
    "sys = ss(A, B, C, D)\n",
    "print(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.  0.  0.]\n",
      " [ 0.  0.  0.]\n",
      " [ 1. -2.  4.]]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "1"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Vc=ctrb(A,B)\n",
    "print(Vc)\n",
    "np.linalg.matrix_rank(Vc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.]\n",
      " [ 0.]\n",
      " [-1.]]\n"
     ]
    }
   ],
   "source": [
    "Tc = scipy.linalg.orth(Vc)  # Image Vcの基底\n",
    "print(Tc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.  0.  0.]\n",
      " [ 0. -1.  0.]\n",
      " [-1.  1.  0.]]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "2"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Vo=obsv(A,C)\n",
    "print(Vo)\n",
    "np.linalg.matrix_rank(Vo)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.]\n",
      " [0.]\n",
      " [1.]]\n"
     ]
    }
   ],
   "source": [
    "Tno = scipy.linalg.null_space(Vo) # Ker Voの基底\n",
    "print(Tno)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "この例では$T_c$と$T_{no}$は（たまたま）平行である．これを補う基底として以下を選ぶ．"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1 0]\n",
      " [0 1]\n",
      " [0 0]]\n"
     ]
    }
   ],
   "source": [
    "T_comp = np.matrix([[1,0,0],[0,1,0]]).transpose()\n",
    "print(T_comp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0. 0. 1.]\n",
      " [1. 0. 0.]\n",
      " [0. 1. 0.]]\n"
     ]
    }
   ],
   "source": [
    "T=inv(np.block([Tno,T_comp]))\n",
    "print(T)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$T$を用いて等価変換"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-2.  1.  0.]\n",
      " [ 0.  0. -1.]\n",
      " [ 0.  1. -1.]]\n"
     ]
    }
   ],
   "source": [
    "Abar = T*A*inv(T)\n",
    "print(Abar)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1.]\n",
      " [0.]\n",
      " [0.]]\n"
     ]
    }
   ],
   "source": [
    "Bbar=T*B\n",
    "print(Bbar)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0. 1. 0.]]\n"
     ]
    }
   ],
   "source": [
    "Cbar=C*inv(T)\n",
    "print(Cbar)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 【５】伝達関数への変換，等価変換"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### (1) 状態方程式から伝達関数への変換"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[-2,1],[2,-3]])\n",
    "B = np.array([[0],[1]])\n",
    "C = np.array([1,1])\n",
    "D = np.array([0])\n",
    "sys = ss(A, B, C, D)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\frac{s + 3}{s^2 + 5 s + 4}$$"
      ],
      "text/plain": [
       "TransferFunction(array([1., 3.]), array([1., 5., 4.]))"
      ]
     },
     "execution_count": 54,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P=ss2tf(sys)\n",
    "P"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### (2) 行列$T$による座標変換"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [],
   "source": [
    "T=np.array([[1,0],[-2,1]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.,  1.],\n",
       "       [-4., -5.]])"
      ]
     },
     "execution_count": 56,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Abar=T@A@inv(T)\n",
    "Abar"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0],\n",
       "       [1]])"
      ]
     },
     "execution_count": 57,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Bbar=T@B\n",
    "Bbar"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([3., 1.])"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Cbar=C@inv(T)\n",
    "Cbar"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### (3)変換後の伝達関数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [],
   "source": [
    "sys_bar = ss(Abar,Bbar,Cbar,D)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\frac{s + 3}{s^2 + 5 s + 4}$$"
      ],
      "text/plain": [
       "TransferFunction(array([1., 3.]), array([1., 5., 4.]))"
      ]
     },
     "execution_count": 60,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Pbar=ss2tf(sys_bar)\n",
    "Pbar"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "name": "chap07.ipynb",
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}