photonAI/MZM_MoE_PINN_Model.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "E7c5HVjM1J5Q"
   },
   "source": [
    "# MZM MoE PINN Model\n",
    "\n",
    "## Overview\n",
    "This notebook implements a Physics-Informed Neural Network (PINN) architecture for modeling a Mach-Zehnder Modulator (MZM). The model is built using PyTorch and includes:\n",
    "\n",
    "- Data preprocessing\n",
    "- Custom activation functions\n",
    "- Flexible neural network design\n",
    "- Training loop with learning rate scheduling\n",
    "- Performance visualization\n",
    "\n",
    "The goal is to approximate nonlinear relationships between device parameters and output performance metrics.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 418
    },
    "id": "vazB3wCHBm1c",
    "outputId": "19953917-22ca-48ed-fbb3-b3b976ce6819"
   },
   "outputs": [],
   "source": [
    "# ==============================\n",
    "# Import Required Libraries\n",
    "# ==============================\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from io import StringIO\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "import re\n",
    "import torch\n",
    "import torch.nn as nn\n",
    "from torch.optim.lr_scheduler import StepLR\n",
    "from torch.utils.data import DataLoader, TensorDataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "id": "b459hcyW1J5S"
   },
   "outputs": [],
   "source": [
    "# random seed for selecting test set\n",
    "random_state = 123\n",
    "\n",
    "# Custum functions are defined here\n",
    "class GaussianActivation(nn.Module):\n",
    "    def forward(self, x):\n",
    "        return torch.exp(-x**2)\n",
    "\n",
    "def weights_init(layer_in):\n",
    "  # Initialize the parameters with He initialization\n",
    "  if isinstance(layer_in, nn.Linear):\n",
    "    nn.init.kaiming_uniform_(layer_in.weight)\n",
    "    # nn.init.xavier_uniform_(layer_in.weight)\n",
    "    layer_in.bias.data.fill_(0.0)\n",
    "\n",
    "def create_flexible_nn(input_dim, output_dim, hidden_dims, activation_fn=nn.ReLU(),positive_output=False):\n",
    "    \"\"\"\n",
    "    Create a flexible neural network with variable number of hidden layers\n",
    "\n",
    "    Args:\n",
    "        input_dim (int): Number of input features\n",
    "        output_dim (int): Number of output targets\n",
    "        hidden_dims (list): List of integers specifying hidden layer dimensions\n",
    "        activation_fn: Activation function to use between layers\n",
    "    \"\"\"\n",
    "    layers = []\n",
    "\n",
    "    # Input layer\n",
    "    layers.append(nn.Linear(input_dim, hidden_dims[0]))\n",
    "    layers.append(activation_fn)\n",
    "\n",
    "    # Hidden layers\n",
    "    for i in range(len(hidden_dims) - 1):\n",
    "        layers.append(nn.Linear(hidden_dims[i], hidden_dims[i+1]))\n",
    "        layers.append(activation_fn)\n",
    "\n",
    "    # Output layer\n",
    "    layers.append(nn.Linear(hidden_dims[-1], output_dim))\n",
    "\n",
    "    # If positivity is requested, append a Softplus\n",
    "    if positive_output:\n",
    "        layers.append(nn.Softplus())\n",
    "\n",
    "    return nn.Sequential(*layers)\n",
    "\n",
    "class ExpertNN(nn.Module):\n",
    "    def __init__(self, input_dim, output_dim, hidden_dims, activation_fn=nn.ReLU(),\n",
    "             dropout_rate=0.0, use_bn=False):\n",
    "        super().__init__()\n",
    "        layers = []\n",
    "\n",
    "        prev_dim = input_dim\n",
    "        for h in hidden_dims:\n",
    "            layers.append(nn.Linear(prev_dim, h))\n",
    "            if use_bn:\n",
    "                layers.append(nn.BatchNorm1d(h))\n",
    "            layers.append(activation_fn)\n",
    "            if dropout_rate > 0:\n",
    "                layers.append(nn.Dropout(p=dropout_rate))\n",
    "            prev_dim = h\n",
    "        layers.append(nn.Linear(prev_dim, output_dim))\n",
    "\n",
    "        self.net = nn.Sequential(*layers)\n",
    "\n",
    "    def forward(self, x):\n",
    "        return self.net(x)\n",
    "\n",
    "\n",
    "class MixtureOfExperts(nn.Module):\n",
    "    def __init__(self, input_dim, output_dim, hidden_dims, n_experts=3,\n",
    "             activation_fn=nn.ReLU(), gating_hidden=32,\n",
    "             dropout_rate=0.0, use_bn=False):\n",
    "        super().__init__()\n",
    "        self.experts = nn.ModuleList([\n",
    "            ExpertNN(input_dim, output_dim, hidden_dims,\n",
    "                    activation_fn, dropout_rate, use_bn)\n",
    "            for _ in range(n_experts)\n",
    "        ])\n",
    "        self.gating = nn.Sequential(\n",
    "            nn.Linear(input_dim, gating_hidden),\n",
    "            nn.ReLU(),\n",
    "            nn.Linear(gating_hidden, n_experts),\n",
    "            nn.Softmax(dim=1)\n",
    "        )\n",
    "\n",
    "    def forward(self, x):\n",
    "        # Gating weights\n",
    "        gate_weights = self.gating(x)  # [batch, n_experts]\n",
    "        # Experts’ outputs\n",
    "        expert_outputs = torch.stack([expert(x) for expert in self.experts], dim=2)  # [batch, output_dim, n_experts]\n",
    "        # Weighted sum\n",
    "        out = torch.bmm(expert_outputs, gate_weights.unsqueeze(2)).squeeze(2)  # [batch, output_dim]\n",
    "        return out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "SbevfhHm1J5T"
   },
   "source": [
    "## Neural Network Architecture Design\n",
    "\n",
    "In this section, we define:\n",
    "\n",
    "- Input dimension (`D_i`)\n",
    "- Output dimension (`D_o`)\n",
    "- Hidden layer structure (`D_h`)\n",
    "- Activation function\n",
    "\n",
    "The architecture is intentionally flexible so it can be tuned for bias-variance tradeoff optimization.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "id": "d0Hkww28B5Da"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MixtureOfExperts(\n",
      "  (experts): ModuleList(\n",
      "    (0-49): 50 x ExpertNN(\n",
      "      (net): Sequential(\n",
      "        (0): Linear(in_features=8, out_features=8, bias=True)\n",
      "        (1): BatchNorm1d(8, eps=1e-05, momentum=0.1, affine=True, track_running_stats=True)\n",
      "        (2): ReLU()\n",
      "        (3): Linear(in_features=8, out_features=6, bias=True)\n",
      "        (4): BatchNorm1d(6, eps=1e-05, momentum=0.1, affine=True, track_running_stats=True)\n",
      "        (5): ReLU()\n",
      "        (6): Linear(in_features=6, out_features=6, bias=True)\n",
      "        (7): BatchNorm1d(6, eps=1e-05, momentum=0.1, affine=True, track_running_stats=True)\n",
      "        (8): ReLU()\n",
      "        (9): Linear(in_features=6, out_features=3, bias=True)\n",
      "      )\n",
      "    )\n",
      "  )\n",
      "  (gating): Sequential(\n",
      "    (0): Linear(in_features=8, out_features=8, bias=True)\n",
      "    (1): ReLU()\n",
      "    (2): Linear(in_features=8, out_features=50, bias=True)\n",
      "    (3): Softmax(dim=1)\n",
      "  )\n",
      ")\n",
      "\n",
      "Number of parameters: 11972\n"
     ]
    }
   ],
   "source": [
    "###### Design the NN architecture\n",
    "\n",
    "acivation_fun = nn.ReLU()\n",
    "# acivation_fun = nn.Tanh()\n",
    "# acivation_fun = GaussianActivation()\n",
    "\n",
    "D_i = 8        # Input dimensions\n",
    "D_o = 3        # Output dimensions\n",
    "# Hidden layer dimensions\n",
    "# D_h = [64, 128, 256, 128, 64]\n",
    "# D_h = [200, 300, 350, 300, 200]\n",
    "# D_h = [8,7,6,5,4]\n",
    "# D_h = [8,6,6]\n",
    "D_h = [64, 128, 64]\n",
    "\n",
    "#positive_output = False # True if we want a constraint to ensure +ve \"BW_3dB\", \"IL\", \"V_pi\" DOESN'T WORK WITH NORMALIZATION\n",
    "#model = create_flexible_nn(D_i, D_o, D_h, acivation_fun,positive_output)\n",
    "n_experts = 60  # number of experts in MoE model\n",
    "gating_hidden_layer_width = 8 #width of the gating hidden layer\n",
    "dropout_rate = 0.0 # dropout rate (regularization parameter)\n",
    "use_bn = True # Batch normalization flag\n",
    "\n",
    "model_data = MixtureOfExperts(\n",
    "    D_i, D_o, D_h,\n",
    "    n_experts=n_experts,\n",
    "    activation_fn=acivation_fun,\n",
    "    gating_hidden=gating_hidden_layer_width,\n",
    "    dropout_rate=dropout_rate,\n",
    "    use_bn=use_bn\n",
    ")\n",
    "\n",
    "print(model_data)\n",
    "print(f\"\\nNumber of parameters: {sum(p.numel() for p in model_data.parameters())}\")\n",
    "\n",
    "# Design the NN training module\n",
    "\n",
    "# SGD parameters\n",
    "batch_size = 128\n",
    "learning_rate = 1e-3 #1e-3\n",
    "weight_decay = 0.05 # Regularization parameter 0.05\n",
    "momentum = 0.9 # used in momentum SGD optimizer\n",
    "betas = (0.9, 0.999) # used in Adam optimizer beta1=0.9, beta2=0.999\n",
    "\n",
    "# Schedular parameters\n",
    "LR_scheduler_gamma = 0.5 # decrease learning by 0.5 every N steps\n",
    "LR_scheduler_step = 100 # step size parameter N of the LR scheduler\n",
    "n_epoch = 100 # loop over the dataset n_epoch, e.g., 200 times\n",
    "\n",
    "# define MSE as loss function (regression problem)\n",
    "loss_function = nn.MSELoss()\n",
    "\n",
    "# construct SGD optimizer\n",
    "# optimizer = torch.optim.SGD(model.parameters(), lr = learning_rate, weight_decay=weight_decay, momentum=momentum)\n",
    "optimizer = torch.optim.AdamW(model_data.parameters(), lr=learning_rate, weight_decay=weight_decay, betas=betas)\n",
    "\n",
    "# learning rate schedular by half every 10 epochs\n",
    "# scheduler = StepLR(optimizer, step_size=LR_scheduler_step, gamma=LR_scheduler_gamma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "id": "zYO1hW9L1J5U"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "      PN_offset  Bias_V    Core_width      P+_width      N+_width  \\\n",
      "0 -2.150000e-07   -10.0  4.500000e-07  2.062650e-07  1.083980e-07   \n",
      "1 -2.150000e-07   -10.0  4.500000e-07  3.059360e-07  3.005950e-07   \n",
      "2 -2.150000e-07   -10.0  5.081920e-07  2.085450e-07  1.936980e-07   \n",
      "3 -2.150000e-07   -10.0  5.146500e-07  2.985000e-07  3.407150e-07   \n",
      "4 -2.150000e-07   -10.0  5.180810e-07  1.589440e-07  2.348480e-07   \n",
      "\n",
      "        P_width       N_width  Phase_length     BW_3dB       IL     V_pi  \n",
      "0  1.000000e-06  7.551270e-07      0.002402  51.200000  2.58569  30.6728  \n",
      "1  1.000000e-06  6.000000e-07      0.003196  36.179487  3.00849  39.7924  \n",
      "2  8.262140e-07  7.634800e-07      0.002991  42.461539  2.82829  33.4881  \n",
      "3  6.387450e-07  8.128130e-07      0.002953  38.692308  2.67160  40.2870  \n",
      "4  6.000000e-07  6.000000e-07      0.001487  60.800000  1.39794  75.1647  \n"
     ]
    }
   ],
   "source": [
    "# Mount Drive\n",
    "# from google.colab import drive\n",
    "# drive.mount('/content/drive')\n",
    "\n",
    "# Load and scale data\n",
    "# file_path = \"/content/drive/MyDrive/MZM Data/Sim_generated_dataset.txt\"\n",
    "file_path = \"Sim_generated_dataset.txt\"\n",
    "\n",
    "with open(file_path) as f:\n",
    "    cleaned = [re.sub(r'[\\[\\]]', '', line.strip()) for line in f]\n",
    "\n",
    "df = pd.read_csv(StringIO(\"\\n\".join(cleaned)), header=None)\n",
    "df.columns = [\n",
    "    \"PN_offset\", \"Bias_V\", \"Core_width\", \"P+_width\", \"N+_width\",\n",
    "    \"P_width\", \"N_width\", \"Phase_length\", \"BW_3dB\", \"IL\", \"V_pi\"\n",
    "]\n",
    "print(df.head()) # for debugging only to check the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "id": "JjcbkxaR1J5U"
   },
   "outputs": [],
   "source": [
    "# Remove the rows where Vpi is above a threshold (data cleaning)\n",
    "df_cleaned = df[df['V_pi'] < 500].copy()\n",
    "\n",
    "# Split data into features/targets\n",
    "feature_cols = df_cleaned.columns[:8]  # first 8 columns\n",
    "target_cols = df_cleaned.columns[8:]   # last 3 columns\n",
    "\n",
    "x = df_cleaned[feature_cols].values\n",
    "y = df_cleaned[target_cols].values\n",
    "# print(y[:3]) # for debugging only to check the data\n",
    "\n",
    "# Train-test split\n",
    "x_train_raw, x_test_raw, y_train_raw, y_test_raw = train_test_split(\n",
    "    x, y, test_size=0.1, random_state=random_state)\n",
    "\n",
    "# Fit scalers on TRAINING data only\n",
    "scaler_x = StandardScaler()\n",
    "\n",
    "# Scale each target variable SEPARATELY\n",
    "scaler_y1 = StandardScaler()  # For BW_3dB\n",
    "scaler_y2 = StandardScaler()  # For IL\n",
    "scaler_y3 = StandardScaler()  # For V_pi\n",
    "\n",
    "# Scale features (x) - one scaler is fine\n",
    "x_train_scaled = scaler_x.fit_transform(x_train_raw)\n",
    "\n",
    "# Scale each target dimension INDEPENDENTLY\n",
    "y_train_scaled_1 = scaler_y1.fit_transform(y_train_raw[:, 0:1])  # BW_3dB\n",
    "y_train_scaled_2 = scaler_y2.fit_transform(y_train_raw[:, 1:2])  # IL\n",
    "y_train_scaled_3 = scaler_y3.fit_transform(y_train_raw[:, 2:3])  # V_pi\n",
    "\n",
    "# Combine targets back together\n",
    "y_train_scaled = np.hstack([y_train_scaled_1, y_train_scaled_2, y_train_scaled_3])\n",
    "\n",
    "# Transform TEST data using training parameters\n",
    "x_test_scaled = scaler_x.transform(x_test_raw)\n",
    "\n",
    "# Transform each test target separately\n",
    "y_test_scaled_1 = scaler_y1.transform(y_test_raw[:, 0:1])  # BW_3dB\n",
    "y_test_scaled_2 = scaler_y2.transform(y_test_raw[:, 1:2])  # IL\n",
    "y_test_scaled_3 = scaler_y3.transform(y_test_raw[:, 2:3])  # V_pi\n",
    "y_test_scaled = np.hstack([y_test_scaled_1, y_test_scaled_2, y_test_scaled_3])\n",
    "\n",
    "# Convert to tensors\n",
    "x_train = torch.from_numpy(x_train_scaled).float()\n",
    "y_train = torch.from_numpy(y_train_scaled).float()\n",
    "x_test = torch.from_numpy(x_test_scaled).float()\n",
    "y_test = torch.from_numpy(y_test_scaled).float()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "vZvTrgGptLKj"
   },
   "source": [
    "The outputs are:\n",
    "\n",
    "- BW_3dB\n",
    "\n",
    "- IL\n",
    "\n",
    "- V_pi\n",
    "\n",
    "These are macroscopic performance metrics, not modeling the field profile directly, so solving full Maxwell is overkill. The most reasonable PINN constraint here is Phase Accumulation Relation enforced by gradients.\n",
    "\n",
    "## ✔ Constraint 1: Monotonic Bandwidth vs Length\n",
    "\n",
    "Physically:\n",
    "\n",
    "$$\n",
    "BW \\propto \\frac{1}{L}\n",
    "$$\n",
    "\n",
    "So:\n",
    "\n",
    "$$\n",
    "\\frac{\\partial BW}{\\partial L} \\le 0\n",
    "$$\n",
    "\n",
    "This must be enforced via autograd, not sorting.\n",
    "\n",
    "---\n",
    "\n",
    "## ✔ Constraint 2: Insertion Loss increases with Length\n",
    "\n",
    "$$\n",
    "\\frac{\\partial IL}{\\partial L} \\ge 0\n",
    "$$\n",
    "\n",
    "Longer waveguide → more absorption.\n",
    "\n",
    "---\n",
    "\n",
    "## ✔ Constraint 3: Vπ scales inversely with L\n",
    "\n",
    "$$\n",
    "V_\\pi \\cdot L = constant\n",
    "$$\n",
    "\n",
    "Instead of using dataset mean, enforce:\n",
    "\n",
    "$$\n",
    "\\frac{\\partial (V_\\pi L)}{\\partial L} \\approx 0\n",
    "$$\n",
    "\n",
    "---\n",
    "\n",
    "## ✔ Constraint 4: Smoothness Constraint\n",
    "\n",
    "Device metrics should be smooth functions of geometry.\n",
    "\n",
    "So penalize second derivatives:\n",
    "\n",
    "$$\n",
    "\\left|\\frac{\\partial^2 y}{\\partial L^2}\\right|\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "id": "z03ciqUp1J5V"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using device: cuda\n"
     ]
    }
   ],
   "source": [
    "# ==============================\n",
    "# GPU Setup\n",
    "# ==============================\n",
    "device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')\n",
    "print(f\"Using device: {device}\")\n",
    "\n",
    "# Move model and all data to GPU\n",
    "model_data.to(device)\n",
    "x_train = x_train.to(device)\n",
    "y_train = y_train.to(device)\n",
    "x_test  = x_test.to(device)\n",
    "y_test  = y_test.to(device)\n",
    "\n",
    "# Precompute constant for Vpi*L on GPU\n",
    "with torch.no_grad():\n",
    "    VpiL_const = (y_train[:, 2] * x_train[:, 7]).mean()\n",
    "\n",
    "# Store polynomial fit coefficients on GPU\n",
    "L_train_np = x_train[:,7].cpu().numpy()  # polyfit uses numpy, move to CPU\n",
    "BW_train_np = y_train[:,0].cpu().numpy()\n",
    "coeffs = np.polyfit(L_train_np, BW_train_np, deg=2)  # [a2, a1, a0]\n",
    "BW_coeffs = torch.tensor(coeffs, dtype=torch.float32, device=device)\n",
    "\n",
    "# Training errors storage\n",
    "errors_train = np.zeros((n_epoch))\n",
    "errors_test  = np.zeros((n_epoch))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "id": "Z5Fa_mZttHEc"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "MixtureOfExperts(\n",
       "  (experts): ModuleList(\n",
       "    (0-49): 50 x ExpertNN(\n",
       "      (net): Sequential(\n",
       "        (0): Linear(in_features=8, out_features=8, bias=True)\n",
       "        (1): BatchNorm1d(8, eps=1e-05, momentum=0.1, affine=True, track_running_stats=True)\n",
       "        (2): ReLU()\n",
       "        (3): Linear(in_features=8, out_features=6, bias=True)\n",
       "        (4): BatchNorm1d(6, eps=1e-05, momentum=0.1, affine=True, track_running_stats=True)\n",
       "        (5): ReLU()\n",
       "        (6): Linear(in_features=6, out_features=6, bias=True)\n",
       "        (7): BatchNorm1d(6, eps=1e-05, momentum=0.1, affine=True, track_running_stats=True)\n",
       "        (8): ReLU()\n",
       "        (9): Linear(in_features=6, out_features=3, bias=True)\n",
       "      )\n",
       "    )\n",
       "  )\n",
       "  (gating): Sequential(\n",
       "    (0): Linear(in_features=8, out_features=8, bias=True)\n",
       "    (1): ReLU()\n",
       "    (2): Linear(in_features=8, out_features=50, bias=True)\n",
       "    (3): Softmax(dim=1)\n",
       "  )\n",
       ")"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# -- Precompute constant for VpiL\n",
    "with torch.no_grad():\n",
    "    VpiL_const = (y_train[:, 2] * x_train[:, 7]).mean()\n",
    "\n",
    "# Define weights for each physics term\n",
    "lambda_vpiL    = 0.0\n",
    "lambda_bw_mon  = 0.0\n",
    "lambda_IL_mon  = 0.0\n",
    "lambda_pos     = 0.0\n",
    "lambda_bw_poly = 0.2\n",
    "\n",
    "data_loader = DataLoader(\n",
    "    TensorDataset(x_train, y_train),\n",
    "    batch_size=batch_size,\n",
    "    shuffle=True,\n",
    "    worker_init_fn=np.random.seed(1)\n",
    ")\n",
    "\n",
    "model_data.apply(weights_init)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "id": "zhkJ08bl1J5W"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoch     0, train MSE 0.5561, test MSE 0.5682\n",
      "Epoch    10, train MSE 0.0699, test MSE 0.0889\n",
      "Epoch    20, train MSE 0.0549, test MSE 0.0779\n",
      "Epoch    30, train MSE 0.0506, test MSE 0.0749\n",
      "Epoch    40, train MSE 0.0478, test MSE 0.0732\n",
      "Epoch    50, train MSE 0.0461, test MSE 0.0711\n",
      "Epoch    60, train MSE 0.0444, test MSE 0.0705\n",
      "Epoch    70, train MSE 0.0439, test MSE 0.0708\n",
      "Epoch    80, train MSE 0.0425, test MSE 0.0699\n",
      "Epoch    90, train MSE 0.0406, test MSE 0.0715\n"
     ]
    }
   ],
   "source": [
    "# ==============================\n",
    "# Training Loop\n",
    "# ==============================\n",
    "for epoch in range(n_epoch):\n",
    "    model_data.train()\n",
    "    for x_batch, y_batch in data_loader:\n",
    "        # Move batch to GPU\n",
    "        x_batch = x_batch.to(device)\n",
    "        y_batch = y_batch.to(device)\n",
    "\n",
    "        optimizer.zero_grad()\n",
    "        pred = model_data(x_batch)\n",
    "\n",
    "        # Data loss\n",
    "        data_loss = loss_function(pred, y_batch)\n",
    "\n",
    "        # Physics-informed losses\n",
    "        physics_loss_total = 0.0\n",
    "\n",
    "        # 1) Vpi*L constant\n",
    "        if lambda_vpiL != 0:\n",
    "            Vpi_pred = pred[:,2]\n",
    "            L_batch  = x_batch[:,7]\n",
    "            residual_vpi = Vpi_pred * L_batch - VpiL_const\n",
    "            physics_loss_vpiL = torch.mean(residual_vpi**2)\n",
    "            physics_loss_total += lambda_vpiL * physics_loss_vpiL\n",
    "\n",
    "        # 2) BW monotonic decreasing w.r.t L\n",
    "        if lambda_bw_mon != 0:\n",
    "            L_sorted, idx = torch.sort(x_batch[:,7])\n",
    "            BW_sorted = pred[idx][:,0]\n",
    "            violation_bw = torch.relu(BW_sorted[1:] - BW_sorted[:-1])\n",
    "            physics_loss_bw = torch.mean(violation_bw**2)\n",
    "            physics_loss_total += lambda_bw_mon * physics_loss_bw\n",
    "\n",
    "        # 3) IL monotonic increasing w.r.t L\n",
    "        if lambda_IL_mon != 0:\n",
    "            IL_sorted = pred[idx][:,1]\n",
    "            violation_il = torch.relu(IL_sorted[:-1] - IL_sorted[1:])\n",
    "            physics_loss_il = torch.mean(violation_il**2)\n",
    "            physics_loss_total += lambda_IL_mon * physics_loss_il\n",
    "\n",
    "        # 4) Positivity penalty\n",
    "        if lambda_pos != 0:\n",
    "            neg_penalty = torch.relu(-pred)\n",
    "            physics_loss_pos = torch.mean(neg_penalty**2)\n",
    "            physics_loss_total += lambda_pos * physics_loss_pos\n",
    "\n",
    "        # 5) BW polynomial fit consistency\n",
    "        if lambda_bw_poly != 0:\n",
    "            L_batch = x_batch[:,7]\n",
    "            BW_pred = pred[:,0]\n",
    "            BW_fit = BW_coeffs[0] * L_batch**2 + BW_coeffs[1] * L_batch + BW_coeffs[2]\n",
    "            residual_bw_poly = BW_pred - BW_fit\n",
    "            physics_loss_bw_poly = torch.mean(residual_bw_poly**2)\n",
    "            physics_loss_total += lambda_bw_poly * physics_loss_bw_poly\n",
    "\n",
    "        # Total loss\n",
    "        loss = data_loss + physics_loss_total\n",
    "        loss.backward()\n",
    "        optimizer.step()\n",
    "\n",
    "    # ==============================\n",
    "    # Evaluation (no_grad, GPU)\n",
    "    # ==============================\n",
    "    model_data.eval()\n",
    "        with torch.no_grad():\n",
    "        pred_train = model_data(x_train)\n",
    "        pred_test  = model_data(x_test)\n",
    "\n",
    "        train_loss = loss_function(pred_train, y_train).item()\n",
    "        test_loss  = loss_function(pred_test, y_test).item()\n",
    "\n",
    "        errors_train[epoch] = train_loss\n",
    "        errors_test[epoch]  = test_loss\n",
    "\n",
    "        # Print every 10 epochs + first epoch\n",
    "        if epoch % 10 == 0 or epoch == 0:\n",
    "            print(f\"Epoch {epoch:5d}, train MSE {train_loss:.4f}, test MSE {test_loss:.4f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "P3fpdbKU1J5X"
   },
   "source": [
    "## Training Performance Visualization\n",
    "\n",
    "This section plots the Mean Squared Error (MSE) for both training and testing sets across epochs.\n",
    "\n",
    "Monitoring both curves allows us to detect:\n",
    "- Overfitting (train ↓, test ↑)\n",
    "- Underfitting (both high)\n",
    "- Proper convergence (both ↓ and stable)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "id": "Z35AHOJlCFRj"
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAHHCAYAAABDUnkqAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAXPhJREFUeJzt3Xl4E9X+BvA3SZuk6b5AS6HQslOQIhTKDkoVgauCIMtFgYqoP0HlVryKXEFxAZeLuIKoiDsCgriCWNnkIvu+ylKKQDeW7k3b5Pz+OCRp6EJLk0ybvp/nmSfNZDI5mQby9nvOnFEJIQSIiIiI3IRa6QYQERERORLDDREREbkVhhsiIiJyKww3RERE5FYYboiIiMitMNwQERGRW2G4ISIiIrfCcENERERuheGGiIiI3ArDDbmFCRMmIDIyUulmEBFRLcBwQ06lUqmqtGzYsEHpptrZsGEDVCoVVqxYoXRTrstoNOLpp59GeHg4vLy8EBcXh3Xr1lX5+efOncPIkSMREBAAPz8/3H333Th16pTdNgUFBZg4cSI6dOgAf39/+Pj4ICYmBm+99RaKi4vL7HPdunXo3bs3DAYDAgMDMWLECCQnJ9ttc/HiRbz++uvo27cvGjRogICAAHTv3h3ffPNNmf0dOnQI9957L5o3bw6DwYCQkBD07dsXP/zwQ5XfZ2n9+/ev0ufy+eefv6H9X+v999/HkiVLqry9SqXClClTHPLazvbxxx+jXbt20Ov1aNWqFd55550qP7cqn93k5ORKf0eTJk2q9j4BoLi4GC+88AKaN28OnU6H5s2b46WXXkJJSUmlbX755ZehUqnQoUOHKr9Pcj0PpRtA7u3zzz+3u//ZZ59h3bp1Zda3a9euRq/z4Ycfwmw212gfddWECROwYsUKTJ06Fa1atcKSJUswePBgrF+/Hr179670ubm5ubjllluQlZWFZ599Fp6ennjzzTfRr18/7N27F8HBwQBkuDl06BAGDx6MyMhIqNVq/O9//8O//vUvbNu2DV999ZV1nz/++CPuvvtudO7cGXPnzkV2djbeeust9O7dG3v27EGDBg0AAFu3bsWMGTMwePBg/Oc//4GHhwe+/fZbjB49GocPH8YLL7xg3eeZM2eQk5OD8ePHIzw8HPn5+fj2229x11134YMPPsBDDz1UrWM2Y8YMPPjgg9b7O3bswNtvv41nn33W7rPYsWPHau23Iu+//z5CQkIwYcIEh+yvtvjggw/wyCOPYPjw4UhMTMTmzZvx+OOPIz8/H08//fR1n1+Vz26DBg3K/H8BAGvWrMGXX36J22+/vdr7BID77rsPy5cvxwMPPIDY2Fj8+eefeO6555CSkoJFixaV296///4br7zyCry9vatzmEgJgsiFJk+eLKryscvLy3NBayq2fv16AUAsX75c0XZcz7Zt2wQA8frrr1vXFRQUiBYtWogePXpc9/mvvvqqACC2b99uXXfkyBGh0WjE9OnTr/v8KVOmCADiwoUL1nXR0dGiZcuWwmg0Wtft3btXqNVqkZiYaF136tQpkZycbLc/s9ksbr31VqHT6URubm6lr11SUiJiYmJEmzZtrtvO61m+fLkAINavX1/jfZWnffv2ol+/flXeHoCYPHmyU9riKPn5+SI4OFgMGTLEbv3YsWOFt7e3uHTpUqXPr+lnd8CAAcLPz08UFBRUe5/bt28XAMRzzz1nt88nn3xSqFQqsW/fvnJfc9SoUeLWW28V/fr1E+3bt79uG0k57JYixfXv3x8dOnTArl270LdvXxgMBjz77LMAgNWrV2PIkCEIDw+HTqdDixYt8OKLL8JkMtnt49oxN5ZS9htvvIFFixahRYsW0Ol06Nq1K3bs2OGwtp86dQr33nsvgoKCYDAY0L17d/z0009ltnvnnXfQvn17azdNbGysXbUjJycHU6dORWRkJHQ6HRo2bIjbbrsNu3fvrvT1V6xYAY1GY1e50Ov1mDhxIrZu3YqzZ89e9/ldu3ZF165drevatm2LAQMGYNmyZdd9/5ZjfuXKFQDApUuXcPjwYQwbNgxarda6XUxMDNq1a4elS5da10VFRaFZs2Z2+1OpVBg6dCiMRmOZrrFraTQaREREWF/bGX755Rf06dMH3t7e8PX1xZAhQ3Do0CG7bVJTU5GQkIAmTZpAp9OhUaNGuPvuu63dcJGRkTh06BA2btxo7Urp379/jduWl5eHJ598EhEREdDpdGjTpg3eeOMNCCHstrN0EQYEBMDHxwdt2rSx/vuyuN7nszzr16/HxYsX8eijj9qtnzx5MvLy8sr9d1BaTT67Fy5cwPr163HPPfdAr9dXe5+bN28GAIwePdpuv6NHj4YQotyu0U2bNmHFihWYP39+pe+Lagd2S1GtcPHiRQwaNAijR4/Gfffdh9DQUADAkiVL4OPjg8TERPj4+OD333/HzJkzkZ2djddff/26+/3qq6+Qk5ODhx9+GCqVCq+99hruuecenDp1Cp6enjVqc1paGnr27In8/Hw8/vjjCA4Oxqeffoq77roLK1aswLBhwwDILrPHH38cI0aMwBNPPIHCwkLs378f27Ztwz//+U8AwCOPPIIVK1ZgypQpiI6OxsWLF/HHH3/gyJEj6Ny5c4Vt2LNnD1q3bg0/Pz+79d26dQMA7N27FxEREeU+12w2Y//+/XjggQfKPNatWzf8+uuvyMnJga+vr3V9UVERsrOzUVBQgJ07d+KNN95As2bN0LJlSwByvAMAeHl5ldmnwWDAoUOHkJqairCwsArfU2pqKgAgJCSkzGN5eXkoKChAVlYWvv/+e/zyyy8YNWpUhfuqic8//xzjx4/HwIED8eqrryI/Px8LFiywdq9Zgt3w4cNx6NAhPPbYY4iMjER6ejrWrVuHlJQUREZGYv78+Xjsscfg4+ODGTNmAID1832jhBC46667sH79ekycOBGdOnXC2rVr8dRTT+HcuXN48803AcixSv/4xz/QsWNHzJ49GzqdDidOnMCWLVus+6rK57M8e/bsAQDExsbare/SpQvUajX27NmD++67r9Ln3+hnd+nSpTCbzRg7duwN7bOiz6nBYAAA7Nq1y269yWTCY489hgcffBA33XRThe+JahGFK0dUz5TXLdWvXz8BQCxcuLDM9vn5+WXWPfzww8JgMIjCwkLruvHjx4tmzZpZ758+fVoAEMHBwXbl8dWrVwsA4ocffqi0nVXplpo6daoAIDZv3mxdl5OTI6KiokRkZKQwmUxCCCHuvvvu65aw/f39b6gbon379uLWW28ts/7QoUMVHlOLjIwMAUDMnj27zGPvvfeeACCOHj1qt/7rr78WAKxLbGys2L9/v/Vxk8kkAgICxIABA+yel5mZKby9vQUAsXPnzgrbdPHiRdGwYUPRp0+fch9/+OGHra+tVqvFiBEjrtv9URXXdkvl5OSIgIAAMWnSJLvtUlNThb+/v3X95cuXy3SDlMfR3VLfffedACBeeuklu/UjRowQKpVKnDhxQgghxJtvvikAiIyMjAr3VZXPZ3kmT54sNBpNuY81aNBAjB49utLn1+Sz26VLF9GoUSPrv7Hq7vPbb78VAMTnn39ut93ChQsFANGhQwe79e+++67w9/cX6enpQgjBbqk6gN1SVCvodDokJCSUWV/6L6ucnBxkZmaiT58+yM/Px9GjR6+731GjRiEwMNB6v0+fPgBw3S6Pqvj555/RrVs3u0GKPj4+eOihh5CcnIzDhw8DAAICAvD3339X2h0WEBCAbdu24fz589VqQ0FBAXQ6XZn1llJ9QUFBpc8FUK3n33LLLVi3bh2WL1+ORx55BJ6ensjLy7M+rlar8fDDDyMpKQnTp0/HX3/9hV27dmHkyJEoKiqqtE2Wv8SvXLlS4Rk3U6dOxbp16/Dpp59i0KBBMJlM1v060rp163DlyhWMGTMGmZmZ1kWj0SAuLg7r168HID+fWq0WGzZswOXLlx3ejor8/PPP0Gg0ePzxx+3WP/nkkxBC4JdffgEgP1eA7N6taMB9VT6f5SkoKLDreixNr9dX+tmzPP9GPrvHjx/Hrl27MHr0aKjV9l9hVd3n4MGD0axZM0ybNg0rV67EmTNnsGzZMsyYMQMeHh52r33x4kXMnDkTzz33nHUwPNV+DDdUKzRu3Ljc/ygPHTqEYcOGwd/fH35+fmjQoIG11J2VlXXd/TZt2tTuviXoOOKL6MyZM2jTpk2Z9Zazbc6cOQMAePrpp+Hj44Nu3bqhVatWmDx5sl23AAC89tprOHjwICIiItCtWzc8//zzVQpgXl5e1hJ7aYWFhdbHK3sugGo9PzQ0FPHx8RgxYgQWLFiAf/zjH7jtttusXUkAMHv2bEycOBGvvfYaWrdujdjYWHh4eGDixIkAZAAsz2OPPYY1a9bgo48+QkxMTLnbtG3bFvHx8Rg3bhx+/PFH5Obm4s477ywzzqSm/vrrLwDArbfeigYNGtgtv/76K9LT0wHIYPjqq6/il19+QWhoKPr27YvXXnvN7ng4w5kzZxAeHm7XZQiU/eyNGjUKvXr1woMPPojQ0FCMHj0ay5Ytsws6Vfl8lsfLy6vCYFlYWFjpZ8/y/Bv57H755ZcAUKZLqjr71Ov1+OmnnxAcHIzhw4cjMjIS48aNw8yZMxEUFGT3Gf3Pf/6DoKAgPPbYY5W+H6pdGG6oVijvP7IrV66gX79+2LdvH2bPno0ffvgB69atw6uvvgoAVTr1W6PRlLve0V+GlWnXrh2OHTuGpUuXonfv3vj222/Ru3dvzJo1y7rNyJEjcerUKbzzzjsIDw/H66+/jvbt21v/Aq9Io0aNcOHChTLrLevCw8MrfG5QUBB0Ot0NPx8ARowYgdzcXKxevdq6TqvV4qOPPsL58+exadMmHDt2DGvXrkVWVhbUarV1fE5pL7zwAt5//33MnTsX999/f6Wvee3r79ixA8ePH6/yc6rC8tn6/PPPsW7dujJL6fc7depUHD9+HHPmzIFer8dzzz2Hdu3aWcekKMnLywubNm3Cb7/9hvvvvx/79+/HqFGjcNttt1kH5Vfl81meRo0awWQyWYOeRVFRES5evHjdz86Nfna/+uortGnTBl26dKnRPtu3b4+DBw/i4MGD2Lx5M86fP49JkyYhMzMTrVu3BiBD7qJFi/D444/j/PnzSE5ORnJyMgoLC1FcXIzk5GRcunSp0vdJClG4W4zqmYrG3JTXf71q1SoBQGzcuNFu/aJFi8qctlvRmJvyxkIAELNmzaq0nVUZc9O6dWvRrVu3Muvnzp0rAIgDBw6U+zyj0SiGDBkiNBqN3WmspaWlpYnGjRuLXr16VdrOadOmCY1GI7KysuzWv/zyywKASElJqfT5sbGxomvXrmXW33bbbaJ58+aVPlcIeYo3APHqq69Wul1JSYlo1KhRuaf4vvvuuwKAmDp16nVf71rz588XAMS2bduq/dzSrh1zs2zZMgFArF27ttr7On78uDAYDGLs2LHWdR06dHDomJuHHnpIaDQakZ2dbbf+zz//FADEO++8U+FzLZ+NdevWlft4VT6fQgjx448/CgDip59+slu/ZcsWAUB89tlnFT5XiBv77FreX3njxG50n6X99NNPAoD44IMPhBC2/wcqW5544olK90nKYOWGai1L1UWUqrIUFRXh/fffV6pJdgYPHozt27dj69at1nV5eXlYtGgRIiMjER0dDUD22Zem1WoRHR0NIQSKi4thMpnKdLE1bNgQ4eHh5ZbYSxsxYgRMJpPdpGNGoxGffPIJ4uLi7M42SUlJKTNOyVL52Llzp3XdsWPH8Pvvv+Pee++1rsvMzCy32vXRRx8BKHvGzLXeeOMNXLhwAU8++aTd+m+++QaPP/44xo4di3nz5lX4/GurA4CcYfazzz6Dl5eX9Vg7ysCBA+Hn54dXXnml3BmYMzIyAAD5+fnWLg+LFi1awNfX1+535+3t7dBT1gcPHgyTyYR3333Xbv2bb74JlUqFQYMGAUC5VYVOnToBsHVHXu/zWZFbb70VQUFBWLBggd36BQsWwGAwYMiQIdZ1mZmZOHr0KPLz863rqvPZtbCcnl7RWVw3sk+LgoICPPfcc2jUqBHGjBkDAOjQoQNWrVpVZmnfvj2aNm2KVatWWbtbqXbhqeBUa/Xs2ROBgYEYP348Hn/8cahUKnz++ecu7VL69ttvyx24PH78eDzzzDP4+uuvMWjQIDz++OMICgrCp59+itOnT+Pbb7+1Dna8/fbbERYWhl69eiE0NBRHjhzBu+++iyFDhsDX1xdXrlxBkyZNMGLECMTExMDHxwe//fYbduzYgf/+97+Vti8uLg733nsvpk+fjvT0dLRs2RKffvopkpOT8fHHH9ttO27cOGzcuNHu+D366KP48MMPMWTIEEybNg2enp6YN28eQkND7YLIF198gYULF2Lo0KFo3rw5cnJysHbtWqxbtw533nknbr31Vrttv/32W/Tt29f6XpYtW4YHH3wQw4cPt263fft2jBs3DsHBwRgwYIB1LIVFz5490bx5cwDAww8/jOzsbPTt2xeNGzdGamoqvvzySxw9ehT//e9/7cZILFmyBAkJCfjkk09ueEZgPz8/LFiwAPfffz86d+6M0aNHo0GDBkhJScFPP/2EXr164d1338Xx48cxYMAAjBw5EtHR0fDw8MCqVauQlpZmN4dKly5dsGDBArz00kto2bIlGjZsaHfMyrNz50689NJLZdb3798fd955J2655RbMmDEDycnJiImJwa+//orVq1dj6tSpaNGiBQA5/mnTpk0YMmQImjVrhvT0dLz//vto0qSJdSD89T6fFfHy8sKLL76IyZMn495778XAgQOxefNmfPHFF3j55ZcRFBRk3fbdd9/FCy+8gPXr11vn+KnOZxeQp2N/88036N69u/X9Xas6+xw5ciTCw8MRHR2N7OxsLF68GKdOncJPP/1kfd8hISEYOnRomdexzHVT3mNUSyhZNqL6pzrdUkLIEnf37t2Fl5eXCA8PF//+97/F2rVrXdYtVdFiOf375MmTYsSIESIgIEDo9XrRrVs38eOPP9rt64MPPhB9+/YVwcHBQqfTiRYtWoinnnrKWjo3Go3iqaeeEjExMcLX11d4e3uLmJgY8f7771faRouCggIxbdo0ERYWJnQ6nejatatYs2ZNme0sp9xf6+zZs2LEiBHCz89P+Pj4iH/84x/ir7/+sttmx44d4t577xVNmzYVOp1OeHt7i86dO4t58+aJ4uJiu223bdsm+vbtKwIDA4VerxcxMTFi4cKFwmw22233ySefVHqMP/nkE+u2X3/9tYiPjxehoaHCw8NDBAYGivj4eLF69eoy7+edd94RAMo9BhWpaIbi9evXi4EDBwp/f3+h1+tFixYtxIQJE6yns2dmZorJkyeLtm3bCm9vb+Hv7y/i4uLEsmXL7PaTmpoqhgwZInx9fQWA63ZRVXZcXnzxRSGEPF39X//6lwgPDxeenp6iVatW4vXXX7c7zklJSeLuu+8W4eHhQqvVivDwcDFmzBhx/Phx6zbX+3xez6JFi0SbNm2EVqsVLVq0EG+++WaZ3/WsWbPKPb5V/ewKIcSaNWsEAPH2229X2p6q7vPVV18Vbdu2FXq9XgQGBoq77rpL7Nmzp0rvmaeC134qIVz4ZzARkZONHDkSycnJ2L59u9JNISKFsFuKiNyGEAIbNmzAF198oXRTiEhBrNwQERGRW+HZUkRERORWGG6IiIjIrTDcEBERkVthuCEiIiK3Uu/OljKbzTh//jx8fX2hUqmUbg4RERFVgRACOTk5CA8PL3NF+GvVu3Bz/vz5SqfgJiIiotrr7NmzaNKkSaXb1LtwY5lW++zZs/Dz81O4NURERFQV2dnZiIiIqPSyIBb1LtxYuqL8/PwYboiIiOqYqgwp4YBiIiIicisMN0RERORWGG6IiIjIrdS7MTdERETOZDKZUFxcrHQz6iStVnvd07yrguGGiIjIAYQQSE1NxZUrV5RuSp2lVqsRFRUFrVZbo/0w3BARETmAJdg0bNgQBoOBE8VWk2WS3QsXLqBp06Y1On4MN0RERDVkMpmswSY4OFjp5tRZDRo0wPnz51FSUgJPT88b3g8HFBMREdWQZYyNwWBQuCV1m6U7ymQy1Wg/DDdEREQOwq6omnHU8WO4ISIiIrfCcENEREQOERkZifnz5yvdDA4oJiIiqs/69++PTp06OSSU7NixA97e3jVvVA0x3DhIURGQlgYIATRtqnRriIiIHEMIAZPJBA+P60eGBg0auKBF18duKQfZvl2Gmvh4pVtCRERUNRMmTMDGjRvx1ltvQaVSQaVSYcmSJVCpVPjll1/QpUsX6HQ6/PHHHzh58iTuvvtuhIaGwsfHB127dsVvv/1mt79ru6VUKhU++ugjDBs2DAaDAa1atcL333/v9PfFcOMgOp28NRqVbQcREdUCQgB5ecosQlS5mW+99RZ69OiBSZMm4cKFC7hw4QIiIiIAAM888wzmzp2LI0eOoGPHjsjNzcXgwYORlJSEPXv24I477sCdd96JlJSUSl/jhRdewMiRI7F//34MHjwYY8eOxaVLl2p0eK+H3VIOotfL28JCZdtBRES1QH4+4OOjzGvn5gJVHPfi7+8PrVYLg8GAsLAwAMDRo0cBALNnz8Ztt91m3TYoKAgxMTHW+y+++CJWrVqF77//HlOmTKnwNSZMmIAxY8YAAF555RW8/fbb2L59O+64445qv7WqYuXGQVi5ISIidxIbG2t3Pzc3F9OmTUO7du0QEBAAHx8fHDly5LqVm44dO1p/9vb2hp+fH9LT053SZgtWbhyElRsiIrIyGGQFRanXdoBrz3qaNm0a1q1bhzfeeAMtW7aEl5cXRowYgaKiokr3c+1lFFQqFcxms0PaWBGGGwcpXbkRAuAklURE9ZhKVeWuIaVptdoqXe5gy5YtmDBhAoYNGwZAVnKSk5Od3Lobw24pB7FUbgB5WjgREVFdEBkZiW3btiE5ORmZmZkVVlVatWqFlStXYu/evdi3bx/++c9/Or0Cc6MYbhxEl3nO+jPH3RARUV0xbdo0aDQaREdHo0GDBhWOoZk3bx4CAwPRs2dP3HnnnRg4cCA6d+7s4tZWjUqIapwz5gays7Ph7++PrKws+Pn5OWy/Ysv/oO7dE4CczK9hQ4ftmoiIarnCwkKcPn0aUVFR0Jcu5VO1VHYcq/P9zcqNg6j0OmghSzas3BARESmH4cZR9HroIU+V4hlTREREymG4cRSdDjpWboiIiBTHcOMorNwQERHVCgw3jlK6cpN//fkCiIiIyDkYbhyldOUmp1jhxhAREdVfDDeOUrpyk8dwQ0REpBSGG0fx9CxVuSlRuDFERET1F8ONo6hU0KllxcaYx3BDRESkFIYbB9JfDTeFuQw3RERESmG4cSCd5mrlhmdLERFRHdG/f39MnTrVYfubMGEChg4d6rD93QiGGwfSa2TFpjCP4YaIiEgpDDcOpPOQ4cZYwHBDRES134QJE7Bx40a89dZbUKlUUKlUSE5OxsGDBzFo0CD4+PggNDQU999/PzIzM63PW7FiBW666SZ4eXkhODgY8fHxyMvLw/PPP49PP/0Uq1evtu5vw4YNLn9fHi5/RTem95ChpjC/Xl1onYiIriEEkJ+vzGsbDIBKVbVt33rrLRw/fhwdOnTA7NmzAQCenp7o1q0bHnzwQbz55psoKCjA008/jZEjR+L333/HhQsXMGbMGLz22msYNmwYcnJysHnzZgghMG3aNBw5cgTZ2dn45JNPAABBQUHOeqsVYrhxIJ2nDDfGArPCLSEiIiXl5wM+Psq8dm4u4O1dtW39/f2h1WphMBgQFhYGAHjppZdw880345VXXrFut3jxYkREROD48ePIzc1FSUkJ7rnnHjRr1gwAcNNNN1m39fLygtFotO5PCQw3DqS/Gm4KC1i5ISKiumnfvn1Yv349fMpJZydPnsTtt9+OAQMG4KabbsLAgQNx++23Y8SIEQgMDFSgteVjuHEgnacMNcZChhsiovrMYJAVFKVeuyZyc3Nx55134tVXXy3zWKNGjaDRaLBu3Tr873//w6+//op33nkHM2bMwLZt2xAVFVWzF3cQhhsH0mtld1ShUeGGEBGRolSqqncNKU2r1cJksp0I07lzZ3z77beIjIyEh0f5MUGlUqFXr17o1asXZs6ciWbNmmHVqlVITEwssz8l8GwpB9Jpr1ZuGG6IiKiOiIyMxLZt25CcnIzMzExMnjwZly5dwpgxY7Bjxw6cPHkSa9euRUJCAkwmE7Zt24ZXXnkFO3fuREpKClauXImMjAy0a9fOur/9+/fj2LFjyMzMRHGx66+3yHDjQHqdDDeFxioOUyciIlLYtGnToNFoEB0djQYNGqCoqAhbtmyByWTC7bffjptuuglTp05FQEAA1Go1/Pz8sGnTJgwePBitW7fGf/7zH/z3v//FoEGDAACTJk1CmzZtEBsbiwYNGmDLli0uf0/slnIgnV7eGosYboiIqG5o3bo1tm7dWmb9ypUry92+Xbt2WLNmTYX7a9CgAX799VeHte9GsHLjQHqdvC008rASEREphd/CDqTTy4qNsZiHlYiISCn8FnYgvZcMN4UMN0RERIrht7ADWSs3JRqFW0JERFR/Mdw4kN4gD2dhMcMNEVF9JAQnca0JRx0/hhsH0nnJw2k08SQ0IqL6xNPTEwCQr9TVMt1EUVERAECjqVmRgN/CDqT3lr+MwhIeViKi+kSj0SAgIADp6ekAAIPBAFVVL81NAACz2YyMjAwYDIYKZ0auKn4LO5DuargxmjwVbgkREbma5SrYloBD1adWq9G0adMaB0OGGwfSe8vDWWhmuCEiqm9UKhUaNWqEhg0bKnLJAXeg1WqhVtd8xAzDjQPproYbo9kTQsgLpxERUf2i0WhqPGaEaqZWDCh+7733EBkZCb1ej7i4OGzfvr3CbZcsWQKVSmW36PV6F7a2YnofGW4E1GBoJyIiUobi4eabb75BYmIiZs2ahd27dyMmJgYDBw6stM/Sz88PFy5csC5nzpxxYYsrpvPVWn/mlcGJiIiUoXi4mTdvHiZNmoSEhARER0dj4cKFMBgMWLx4cYXPUalUCAsLsy6hoaEubHHFdD62sTaFhQo2hIiIqB5TNNwUFRVh165diI+Pt65Tq9WIj48v9wqlFrm5uWjWrBkiIiJw991349ChQxVuazQakZ2dbbc4i9qghyeKrr6u016GiIiIKqFouMnMzITJZCpTeQkNDUVqamq5z2nTpg0WL16M1atX44svvoDZbEbPnj3x999/l7v9nDlz4O/vb10iIiIc/j6sdDroIUs2rNwQEREpQ/Fuqerq0aMHxo0bh06dOqFfv35YuXIlGjRogA8++KDc7adPn46srCzrcvbsWec1Tq+HDrJkw8oNERGRMhQ9FTwkJAQajQZpaWl269PS0qyTIV2Pp6cnbr75Zpw4caLcx3U6HXQ6XY3bWiWs3BARESlO0cqNVqtFly5dkJSUZF1nNpuRlJSEHj16VGkfJpMJBw4cQKNGjZzVzKorXbkp5MXTiIiIlKD4JH6JiYkYP348YmNj0a1bN8yfPx95eXlISEgAAIwbNw6NGzfGnDlzAACzZ89G9+7d0bJlS1y5cgWvv/46zpw5gwcffFDJtyHpdNAjAwBQmFMMQFv59kRERORwioebUaNGISMjAzNnzkRqaio6deqENWvWWAcZp6Sk2E3FfPnyZUyaNAmpqakIDAxEly5d8L///Q/R0dFKvQWb0pWbvBIw3BAREbmeSghRr/pPsrOz4e/vj6ysLPj5+Tl252Yzemm24n/ohZWfXMGwCQGO3T8REVE9VZ3v7zp3tlStplZDp7o6z01eicKNISIiqp8YbhxMr5YXlSrMZbghIiJSAsONg+k0MtwY80wKt4SIiKh+YrhxML1GVmwKGW6IiIgUwXDjYDoPGW6MBWaFW0JERFQ/Mdw4mN5DVmwK8xluiIiIlMBw42C6q+GGlRsiIiJlMNw4mN7zauWmoF5NH0RERFRrMNw4mM5TVmx4bSkiIiJlMNw4mF4rww2vCk5ERKQMhhsH02llxabQqHBDiIiI6imGGwfT62S4MRpVCreEiIiofmK4cTCdTt4WMtwQEREpguHGwfR6eWssZrghIiJSAsONg+n0MtQUFvHQEhERKYHfwA5mq9zw0BIRESmB38AOpvOSh7SwWKNwS4iIiOonhhsH03vJbiljCcMNERGREhhuHExnkKGmsMRD4ZYQERHVTww3DqY3yENqNDHcEBERKYHhxsGslRuTp8ItISIiqp8YbhxM7yMrNqzcEBERKYPhxsEs4abQrIXghcGJiIhcjuHGwXQ+sjtKQI2SEoUbQ0REVA8x3DiYpXIDAIWFCjaEiIionmK4cTCdr9b6s9GoYEOIiIjqKYYbB1N76eCJIgCs3BARESmB4cbR9HroIEs2rNwQERG5HsONo+l00EOWbFi5ISIicj2GG0dj5YaIiEhRDDeOxsoNERGRohhuHE2ns1Vu8k0KN4aIiKj+YbhxNL3eVrnJLlK4MURERPUPw42jla7c5BUr3BgiIqL6h+HG0Tw8oL8abgpzef0FIiIiV2O4cQKdWlZsjAw3RERELsdw4wR6jQw3hXkcUExERORqDDdOoNPIig3PliIiInI9hhsn0Htcrdww3BAREbkcw40T6DxkqDHmmxVuCRERUf3DcOME+qvhprCA4YaIiMjVGG6cQOcpQ42xQCjcEiIiovqH4cYJ9J5XKzeFDDdERESuxnDjBDqtDDVGXjiTiIjI5RhunECvld1SvCo4ERGR6zHcOIFOd7Vyw+tmEhERuRzDjRPor3ZLFRpVCreEiIio/mG4cQKdXoYaYxHDDRERkasx3DiBXi9vC4t4eImIiFyN375OYK3cFPPwEhERuRq/fZ1A7yXDTWGxRuGWEBER1T8MN06g85KH1VjCcENERORqtSLcvPfee4iMjIRer0dcXBy2b99epectXboUKpUKQ4cOdW4Dq0lvkIe1sMRD4ZYQERHVP4qHm2+++QaJiYmYNWsWdu/ejZiYGAwcOBDp6emVPi85ORnTpk1Dnz59XNTSqmPlhoiISDmKh5t58+Zh0qRJSEhIQHR0NBYuXAiDwYDFixdX+ByTyYSxY8fihRdeQPPmzV3Y2qrR+8iKTaHJU+GWEBER1T+KhpuioiLs2rUL8fHx1nVqtRrx8fHYunVrhc+bPXs2GjZsiIkTJ7qimdWmM8iKjdHMcENERORqig4KyczMhMlkQmhoqN360NBQHD16tNzn/PHHH/j444+xd+/eKr2G0WiE0Wi03s/Ozr7h9lYVKzdERETKUbxbqjpycnJw//3348MPP0RISEiVnjNnzhz4+/tbl4iICCe3EtB5y3BjhgYlJU5/OSIiIipF0cpNSEgINBoN0tLS7NanpaUhLCyszPYnT55EcnIy7rzzTus6s1legdvDwwPHjh1DixYt7J4zffp0JCYmWu9nZ2c7PeDofW0Vm8JCwMfHqS9HREREpShaudFqtejSpQuSkpKs68xmM5KSktCjR48y27dt2xYHDhzA3r17rctdd92FW265BXv37i03tOh0Ovj5+dktzqbzsYWbUj1iRERE5AKKT8SSmJiI8ePHIzY2Ft26dcP8+fORl5eHhIQEAMC4cePQuHFjzJkzB3q9Hh06dLB7fkBAAACUWa8kjbceHihGCTxRWKh0a4iIiOoXxcPNqFGjkJGRgZkzZyI1NRWdOnXCmjVrrIOMU1JSoFbXqaFBgE4HHYwogScrN0RERC6mEkIIpRvhStnZ2fD390dWVpbzuqgOHkTITWG4iBAcOgRERzvnZYiIiOqL6nx/17GSSB2h10MHWbJh5YaIiMi1GG6cQaeDHnKwDcfcEBERuRbDjTOUrtwUmBVuDBERUf3CcOMMpSs3OcUKN4aIiKh+YbhxhtKVmzyGGyIiIldiuHEGT89SlRtef4GIiMiVGG6cQaWCTi0rNsY8hhsiIiJXYrhxEv3VcFOYy3BDRETkSgw3TqLTyFBjzDcp3BIiIqL6heHGSfQeVys3+TwVnIiIyJUYbpxEp5EVG1ZuiIiIXIvhxkn0nrJbqrCgXl26i4iISHEMN06i85DdUZyhmIiIyLUYbpxEr5XdUYWFrNwQERG5EsONk+g8ZagxMtwQERG5FMONk+i1sjuqsFClcEuIiIjqF4YbJ9HprlZujAo3hIiIqJ5huHESvU7eFhpZuSEiInIlhhsn0V0NN8YihhsiIiJXYrhxEv3VbqnCIh5iIiIiV+I3r5PovOShNRbzEBMREbkSv3mdRK+Xt4UMN0RERC7Fb14n0Rk0AABjiUbhlhAREdUvDDdOoveSA4kLiz0UbgkREVH9wnDjJJbKTaGJ4YaIiMiVGG6cRO99tVuK4YaIiMilGG6cROctQ02hyVPhlhAREdUvDDdOYq3cmBluiIiIXInhxkkslRuT0KCkROHGEBER1SMMN06i97VVbHjxTCIiItdhuHESnY8t3BQWKtgQIiKieobhxkk8vHXQQPZHsXJDRETkOgw3zqLTQQeZali5ISIich2GG2fR66GHTDWs3BAREbkOw42zsHJDRESkCIYbZ2HlhoiISBEMN85SunJTIBRuDBERUf3BcOMspSs3eZzFj4iIyFUYbpyldOUml+GGiIjIVRhunEWns1VucosVbgwREVH9wXDjLGo1dKoiAKzcEBERuRLDjRPp1bJiwzE3RERErsNw40QGjRxzk5djUrglRERE9QfDjRMFeuYAAC5f4qngRERErlKtcPPaa6+hoKDAen/Lli0wlpqhLicnB48++qjjWlfHBWtluLl4UaVwS4iIiOqPaoWb6dOnIycnx3p/0KBBOHfunPV+fn4+PvjgA8e1ro4L0uYBAC5dYYGMiIjIVar1rSuEqPQ+2QvS5wMALmUx3BAREbkKv3WdyBJuLmZ5KNwSIiKi+oPhxomCveUkfpeyPRVuCRERUf1R7ZLCRx99BB8fHwBASUkJlixZgpCQEACwG49DQJDharjJZbghIiJylWqFm6ZNm+LDDz+03g8LC8Pnn39eZhuSgnzkDMVX8nUwmQCNRuEGERER1QPVCjfJyclOaoZ7Cgq0Dbi+fBm4WuAiIiIiJ+KYGyfyCAmAH7IAAJcuKdwYIiKieqJa4Wbr1q348ccf7dZ99tlniIqKQsOGDfHQQw/ZTepXVe+99x4iIyOh1+sRFxeH7du3V7jtypUrERsbi4CAAHh7e6NTp05lusZqjaAgBEGmGoYbIiIi16hWuJk9ezYOHTpkvX/gwAFMnDgR8fHxeOaZZ/DDDz9gzpw51WrAN998g8TERMyaNQu7d+9GTEwMBg4ciPT09HK3DwoKwowZM7B161bs378fCQkJSEhIwNq1a6v1ui5RKtxcvKhwW4iIiOqJaoWbvXv3YsCAAdb7S5cuRVxcHD788EMkJibi7bffxrJly6rVgHnz5mHSpElISEhAdHQ0Fi5cCIPBgMWLF5e7ff/+/TFs2DC0a9cOLVq0wBNPPIGOHTvijz/+qNbrukRQEIIhUw0rN0RERK5RrXBz+fJlhIaGWu9v3LgRgwYNst7v2rUrzp49W+X9FRUVYdeuXYiPj7c1SK1GfHw8tm7det3nCyGQlJSEY8eOoW/fvuVuYzQakZ2dbbe4DLuliIiIXK5a4SY0NBSnT58GIIPJ7t270b17d+vjOTk58PSs+pwumZmZMJlMdoHJ8jqpqakVPi8rKws+Pj7QarUYMmQI3nnnHdx2223lbjtnzhz4+/tbl4iIiCq3r8aCgxluiIiIXKxa4Wbw4MF45plnsHnzZkyfPh0GgwF9+vSxPr5//360aNHC4Y28lq+vL/bu3YsdO3bg5ZdfRmJiIjZs2FDuttOnT0dWVpZ1qU5lqcZKdUtxzA0REZFrVGuemxdffBH33HMP+vXrBx8fHyxZsgRardb6+OLFi3H77bdXeX8hISHQaDRIS0uzW5+WloawsLAKn6dWq9GyZUsAQKdOnXDkyBHMmTMH/fv3L7OtTqeDTqercpscqnS3VKYZPPOeiIjI+aoVbkJCQrBp0yZrt5Dmmil3ly9fDl9f3yrvT6vVokuXLkhKSsLQoUMBAGazGUlJSZgyZUqV92M2m2/oFHSnCwiwhZv0YgAKhSwiIqJ6pFrh5oEHHqjSdhWd6VSexMREjB8/HrGxsejWrRvmz5+PvLw8JCQkAADGjRuHxo0bW08xnzNnDmJjY9GiRQsYjUb8/PPP+Pzzz7FgwYLqvBXX0GgQ5F0E5AEXM8xKt4aIiKheqFa4WbJkCZo1a4abb74ZQojrP6EKRo0ahYyMDMycOROpqano1KkT1qxZYx1knJKSArXa1p2Tl5eHRx99FH///Te8vLzQtm1bfPHFFxg1apRD2uNowf4lQB5w6ZJK6aYQERHVCypRjZQyefJkfP3112jWrBkSEhJw3333ISgoyJntc7js7Gz4+/sjKysLfn5+Tn+9ox1GoN2hFQjwLsLlXO31n0BERERlVOf7u1ojXN977z1cuHAB//73v/HDDz8gIiICI0eOxNq1ax1WyXE3QSHyEF/J08JkUrgxRERE9UC1T9/R6XQYM2YM1q1bh8OHD6N9+/Z49NFHERkZidzcXGe0sU4LCrNVay5fVrAhRERE9USNzk1Wq9VQqVQQQsDEskS5eGVwIiIi16p2uDEajfj6669x2223oXXr1jhw4ADeffddpKSkwMfHxxltrNt4CQYiIiKXqtbZUo8++iiWLl2KiIgIPPDAA/j6668REhLirLa5h6vhJhlRnKWYiIjIBaoVbhYuXIimTZuiefPm2LhxIzZu3FjuditXrnRI49wCrwxORETkUtUKN+PGjYNKxflaqoXdUkRERC5V7Un8qJqCghCEFAAMN0RERK7AKzk6G68MTkRE5FIMN84WHGzrlrrI60sRERE5G8ONswUG2sJNWonCjSEiInJ/DDfO5uGBIK9CALwyOBERkSsw3LhAsL+s2Fzi5ReIiIicjuHGBSwXTr90RaNsQ4iIiOoBhhsXsF0Z3JNXBiciInIyhhsXCGpom06IVwYnIiJyLoYbF+CVwYmIiFyH4cYVSs91w3BDRETkVAw3rlDq+lKcpZiIiMi5GG5cgVcGJyIichmGG1fglcGJiIhchuHGFRhuiIiIXIbhxhV4ZXAiIiKXYbhxhdKVm4tC4cYQERG5N4YbVygdbtJ5ZXAiIiJnYrhxBa0WQfoCAMDFDF5/gYiIyJkYblzEemVwDigmIiJyKoYbFwkKlGNtLmXxyuBERETOxHDjItYrg+fyyuBERETOxHDjIrwyOBERkWsw3LgIrwxORETkGgw3rsJZiomIiFyC4cZVeGVwIiIil2C4cZXgYF4ZnIiIyAUYblyF3VJEREQuwXDjKgw3RERELsFw4yq8MjgREZFLMNy4Cq8MTkRE5BIMN64SGGgLNxm8MjgREZGzMNy4ipcXgrR5AICL6WaFG0NEROS+GG5cKNivGAAHFBMRETkTw40LBQXJ24tXeNiJiIichd+yLhTaQHZHZeV5Ii9P4cYQERG5KYYbFwoM1cIfVwAAycmKNoWIiMhtMdy4UlAQmuMUAOD0aYXbQkRE5KYYblwpKAhRkKnm1CmF20JEROSmGG5ciZUbIiIip2O4cSVWboiIiJyO4caVWLkhIiJyOoYbV7qmciN4iSkiIiKHY7hxpeBgNMMZqGBGXh6Qmal0g4iIiNwPw40rBQdDDyPCcR4Ax90QERE5A8ONK4WFAZ6eHHdDRETkRLUi3Lz33nuIjIyEXq9HXFwctm/fXuG2H374Ifr06YPAwEAEBgYiPj6+0u1rFY0GaNaMZ0wRERE5keLh5ptvvkFiYiJmzZqF3bt3IyYmBgMHDkR6enq522/YsAFjxozB+vXrsXXrVkREROD222/HuXPnXNzyGxQZycoNERGREykebubNm4dJkyYhISEB0dHRWLhwIQwGAxYvXlzu9l9++SUeffRRdOrUCW3btsVHH30Es9mMpKQkF7f8BkVFsXJDRETkRIqGm6KiIuzatQvx8fHWdWq1GvHx8di6dWuV9pGfn4/i4mIEBQWV+7jRaER2drbdoqioKFZuiIiInEjRcJOZmQmTyYTQ0FC79aGhoUhNTa3SPp5++mmEh4fbBaTS5syZA39/f+sSERFR43bXSKnKTUoKUFKibHOIiIjcjeLdUjUxd+5cLF26FKtWrYJery93m+nTpyMrK8u6nD171sWtvEZUFBrhAnQohMkEKN0cIiIid+Oh5IuHhIRAo9EgLS3Nbn1aWhrCwsIqfe4bb7yBuXPn4rfffkPHjh0r3E6n00Gn0zmkvQ4RFQU1BCKRjGNoi1OngKgopRtFRETkPhSt3Gi1WnTp0sVuMLBlcHCPHj0qfN5rr72GF198EWvWrEFsbKwrmuo4DRoABgPH3RARETmJopUbAEhMTMT48eMRGxuLbt26Yf78+cjLy0NCQgIAYNy4cWjcuDHmzJkDAHj11Vcxc+ZMfPXVV4iMjLSOzfHx8YGPj49i76PKVCogMhJRh3nGFBERkTMoHm5GjRqFjIwMzJw5E6mpqejUqRPWrFljHWSckpICtdpWYFqwYAGKioowYsQIu/3MmjULzz//vCubfuOiotD8MCs3REREzqB4uAGAKVOmYMqUKeU+tmHDBrv7ycnJzm+Qs5U6Y4rhhoiIyLHq9NlSdRYn8iMiInIahhsllJrILyMDyM1VuD1ERERuhOFGCVFR8Ec2AlWXAbBrioiIyJEYbpQQGQkAaC5OAmC4ISIiciSGGyUEBAABARx3Q0RE5AQMN0rhBTSJiIicguFGKTxjioiIyCkYbpTCyg0REZFTMNwo5ZqJ/IRQuD1ERERuguFGKVFRaIYzUMGM/HwgPV3pBhEREbkHhhulREVBi2I0UZ0DwHE3REREjsJwo5RmzQBwrhsiIiJHY7hRisEAhIZax92cPKlwe4iIiNwEw42SoqIQjcMAgL17lW0KERGRu2C4UVJUFOKwDQCwbZvCbSEiInITDDdKiopCF+yCWmXGuXPAuXNKN4iIiKjuY7hRUlQUvJGPDj7JAIDt25VtDhERkTtguFFSVBQAIE69EwC7poiIiByB4UZJkZEAgLi83wEw3BARETkCw42SmjYF1Gp0K9kCANi5EzCZFG4TERFRHcdwoyRPT6BJE0TjMHy8TMjNBY4cUbpRREREdRvDjdKioqCBGbGRGQDYNUVERFRTDDdKuzqouFuwvLgUz5giIiKqGYYbpbVrBwCIM8lxN6zcEBER1QzDjdK6dAEAxP29EgBw4ACQl6dkg4iIiOo2hhulde4MAGh89k+Eh5lhNgO7dyvcJiIiojqM4UZpgYFAixYAgLgWHFRMRERUUww3tUFsLAAgzldeIZzhhoiI6MYx3NQGV8fddMvbAIBnTBEREdUEw01tcDXcxJ75FioVkJICpKYq3CYiIqI6iuGmNrg6qNg35RDaty0BwK4pIiKiG8VwUxsEBAAtWwIAujWVJRt2TREREd0YhpvawjKo2OsAAFZuiIiIbhTDTW1hmcwvex0AYMcOwGxWskFERER1E8NNbXG1ctP+xGr4+gLZ2cDWrQq3iYiIqA5iuKktbr4ZAOCRcgpDBxUCAJYuVbJBREREdRPDTW3h7w+0bg0AGN3hEABg2TKgpETJRhEREdU9DDe1ydVxN/HmXxEUBKSnAxs3KtwmIiKiOobhpja5Gm60+3Zg+HC5il1TRERE1cNwU5tcHVSMXbswerT8ceVKoKhIuSYRERHVNQw3tcnVQcVISUG/6AyEhgKXLgG//aZss4iIiOoShpvaxM8PaNMGAKDZuwsjR8rV7JoiIiKqOoab2ubquBvs3GntmvruO6CgQLEWERER1SkMN7VNqXE33bsDTZsCOTnAL78o2ywiIqK6guGmtilVuVGrgVGj5F12TREREVUNw01tc/PNgEoF/P03kJxs7Zr68UdZwSEiIqLKMdzUNr6+QP/+8uelS3HzzUCrVnLMzfffK9oyIiKiOoHhpjYaO1befvklVCpYqzf//jdw+rRyzSIiIqoLGG5qo+HDAa0WOHgQ2L8f//oX0L49cP48EB8PXLigdAOJiIhqL4ab2iggABgyRP781VcIDAR+/RVo3hw4dQq4/XY5uR8RERGVxXBTW/3zn/L2q68Asxnh4XKm4kaNZEFn0CAOMCYiIioPw01tNWSInLH47FlgyxYAQFQUsG4dEBQEbN8ODB0KFBYq20wiIqLaRvFw89577yEyMhJ6vR5xcXHYvn17hdseOnQIw4cPR2RkJFQqFebPn++6hrqalxdwzz3y5y+/tK5u3x5Yswbw8QF+/x0YMwYoKVGojURERLWQouHmm2++QWJiImbNmoXdu3cjJiYGAwcORHp6ernb5+fno3nz5pg7dy7CwsJc3FoFWM6aWr7c7tLgXbvK08J1OnlphkmTALNZmSYSERHVNoqGm3nz5mHSpElISEhAdHQ0Fi5cCIPBgMWLF5e7fdeuXfH6669j9OjR0Ol0Lm6tAm65BQgLk6OH164t89DSpYBaDSxZAjz1FCCEMs0kIiKqTRQLN0VFRdi1axfi4+NtjVGrER8fj61btzrsdYxGI7Kzs+2WOkOjsU1y89VXZR4eOhT4+GP587x5wJw5rmsaERFRbaVYuMnMzITJZEJoaKjd+tDQUKSmpjrsdebMmQN/f3/rEhER4bB9u4TlrKnVq8s9PWrCBOC//5U/z5ghf2YFh4iI6jPFBxQ72/Tp05GVlWVdzp49q3STqic21nb9he++K3eTxETg2Wflz9OmyWIPTxMnIqL6SrFwExISAo1Gg7S0NLv1aWlpDh0srNPp4OfnZ7fUKSqVbWDxv/8NJCeXu9lLLwFvvgl4eADLlslBxwcPuq6ZREREtYVi4Uar1aJLly5ISkqyrjObzUhKSkKPHj2Ualbt9K9/ATfdBKSmAnfcAVy8WGYTlQqYOhXYtAlo0gQ4dgzo1g348EMgLY1dVUREVH8o2i2VmJiIDz/8EJ9++imOHDmC//u//0NeXh4SEhIAAOPGjcP06dOt2xcVFWHv3r3Yu3cvioqKcO7cOezduxcnTpxQ6i24hp8f8MsvQESETC133SW7qcrRowewe7e8RENBAfDQQ/KEq5AQoE8fef+aM8uJiIjcikoIZf+mf/fdd/H6668jNTUVnTp1wttvv424uDgAQP/+/REZGYklS5YAAJKTkxEVFVVmH/369cOGDRuq9HrZ2dnw9/dHVlZW3euiOnQI6N0buHIFGDZMphSNptxNTSbgtdfk2VSnTpWt3ISGyvlxHnpIZiYiIqLarDrf34qHG1er0+EGADZulGWZoiJg8mTgnXdkn1QlCgpkwefIEVnV+eIL2cMFyHly7rwTePhhudsKshIREZGiGG4qUefDDSBHDI8eLcsxr74qBxpXQ3GxPPFqwQJg/Xrb+qZNgYkTgQcekON2iIiIaguGm0q4RbgBgPnz5UBjQF57yjIfTjUdOSIHHX/6qZwIGZDVnNtuk5e2uvtu2YVlYTIBe/cCGzYAWi0wahTQsGFN3ggREdH1MdxUwm3CDQA8+aScmtjTUw44HjDghndVWAisXAksWiR7vixUKjnMp39/YP9++diVK7bHPTxkt9bEicDAgfI+ERGRozHcVMKtwo3ZLC8LvmyZPKNq82agY8ca7/avv4AVK2TY2bmz7ON+fkDfvkB6OlD6Iu7h4fLlw8JsS7NmQJcuspvrOkODiIiIKsRwUwm3CjeALLkMHCgnuAkPB7ZulYNnHCQlBVi1CtixQ061c+utwM032yo0Bw/KM7I++8zWrVWehg3lZMtdusihQunpcv6d9HSgpASIiZGPdekiA1J2NrBnjxwAvXu3nNpn7FjgvvtkdxgREdUvDDeVcLtwAwCXL8tJbA4dkgNkliyRk/25kNEo89Xff8szsSzL8ePAgQNyrE5VaTQVb9+kieyNmzQJ8PZ2TNuJiKj2Y7iphFuGGwA4exYYNEgGHAB44glg7lxAr1e2XZCnou/bJ7u49u2TlZfQUFnNCQ2VvWu7dwO7dsltLBMwt2oFdO4sK0Vmszzr/cIF+VhwsJywMDtbjgG6cgXIzweCguQ+LfuPiADatJFLy5bytS9elAWu//1PLpcuAVFRQPPmQIsW8jY4GPDxsS3e3oBOx641IiKlMNxUwm3DDSBTxNNPyxQAAB06AF99JfuT6gghgHPn5Liea389hYWy++u114CTJ6u/b7Vahh5LQLqR53t7AwaDvG3USAawli3lEhUFNGggA5afX8VBSAjg/HmZQw8flvsbPJin3xMRVYbhphJuHW4sfv4ZSEiQA1o0Gnle96hRwNChQECA0q2rMZMJ+OknOWYnIADw95eLwSCrMJaxPGlpwOnTcgLDY8fsr5Teti3Qs6dcwsPl9UhPnpSzOZ8+LStBublyKSysfhs9PGTI8fGRFR/LYjbLrrqsrLLP6dxZnnk2aJB8/qVLssp08aKsOMXGypzKM9KIqD5iuKlEvQg3gPx2nzQJ+P572zqtVg4+HjdOTmDj6alc+1xMCFmxSUmRVZaQkKo/t6QEyMuT3V6W29xc2RN44oQ8u+zECRmQLl2q8LJfdjQa2Y7oaBnCtm6t2sVNvbzkoOuuXeXPBQW2xWwGfH3l4ucnb4uKZNddVpa8BWTWHTKkbGUMkKFu715ZWUpLs42dysqS+zIa5WI2y4HfffrI5UbGsJeUyHmWvLyAyEiGNiKqHMNNJepNuLE4fhz45hu5WMbjALJc8cgjMgCFhSnXPjdUUGCruOTn2wKBJRS0bCm7s3Q623PS02U16ocf5MBsvV6O+wkKkrdZWfKMtfIqPjdCq7VN1OjhAWzZIpfSH5HqaNoU6N5dvrcWLeRt8+YyPxcXy2BUXCwHnG/ZAvzxhwx0ubny+R4ecvtWreRtUBAQGCgrcwEB8hhYxmhd2+VnMsljnpMjw1lWlry1VOo0GtmlqNHIY24JgL6+stp3/rwMpydOyOqd0SinOrjtNvvQZpnAcv16OXbs0iU5lv/yZfl6Pj7247aaNZOv4eVlW8xmW0UwN9c2TswydULDhtf/m8Nslm0sKbEtJpP8nRoM9mPDhJDHwdLGwEA5Dq0qY8fy8+XvKyNDhvDAwPK3KymR/82YTLYKpVYr72dkyM92ejqQmSl/dxER8rhGRFR+UkBJifxjxPIHxPHjtj8kAgPl9YOHDQPatSv7fgoK5DHy96/5OLmSEnksSv8hodHYxuJ5e1ctmJtMwJkz8n0cOyZv8/Jsn+vQUNmt7eEhf2+WpbDQvop7+bL8XIWHy6VRI/k+L1+Wj2dmyu1btpR/xHh5lW3LyZPyK6GwEIiLk0tV/uArLJR/9Pj6yt/Btce2uFgOK0hJke+jZ8+qHeOqYripRL0LN6UdPCjH4Hz8sfzfBpD/k44YAUyfXqfG5tRHli6tbdvkafJms/0Xp1otvzBzcmSVJjtbftH4+cn//Pz85Bfcd9/J/1wrYvliDguT/+GGhcmAUbp7raREznG0aZMcDF6ds+FK8/WV/yFWp+tPq5X/sRYV2cKjs7RqJac/uHBBvtfSE1g6i6+vfI+enrZbo1G+1/z86x8ry9gwDw8Z9Mxm+8d9fGRYiY6WQSw/3xbSLl+WlbqzZ+XPFiqVrNT17Qv06yfbaAnE27bZQmp1BQTIL1VLiA8Kkq974oTsHi4uvv4+WreWJ4fm5Mgv7ZMn5RcsII9dgwa2xdtb/uGg09nfWn7WauV/jcnJtiUz8/pt0OnkZzIw0BbMhZDvxRKCL12q2vtxJB8fORrhn/8EunWT//aXLJF/XFyrRQtZFTYY5GfGsmRny5D799/2x8LT0xbKPTzkZ+b8edvnrX9/+8v7OALDTSXqdbixMBqBb78F3n1X/vkMyP+9Ro4Enn9eDkght3b4sJyk8ccf5V+hvXrZxiBV93IaubnyY7R/v/xisVRAzpyxVRQsX9QBAfIst9695Wu2by8/eufOyb/K//pLPu/KFVu14fJl+Z9qerr9uKlrqdUyxFnGYVn+eZvNsh0mkwwGOTm2xWiUX6qlK05mM5CUJMPbtaHNMoFlz54y+Fm+0Pz9ZZA4dcq2pKTIv8xL/8WvVstg4OMjb3U62zixtDQZGm+EWl02xJSm08njUt0vWB8f+by//77+dgaDfZVSrZaBwlKZCAmRxyglRf6OK/tdlm63paLXqpUMMi1byuO7ahXw228y5LqKl5cMQiaT/N1WN9TrdLL9ljM4fX1ldcsyTjAjQ+5TpbItWq38jFqWwEAZOM6fty3Z2baAGBJim9P1zJny26FWA/HxsuqzbRtw9GjV34NWW/kx12plVa5HD+Dzz6t3fK6H4aYSDDfX2L1bnjK+fLm8r1bLmfKefVb+6yO6QWaz7T9oRykokF8Cly/LLxmDwbbo9fLjWx0lJRV3KWRlyWuobd4sv6BvucV+AktHM5tl+Lhyxb4rr7hYfil6edm/V09P2RaNRh5jS/dJXp5ciott3XqWroniYhk+Dx2Sy5kzti4GSzegZQqFiAhbF2BqqjwOGzfKClZ+vuyG7NXLFlI1Gtt7sXyrVPa7z8qSofbiRVu3y6VLsj2WMNO4ceW/0+xsYM0aWYlo0EAGVMvi4yPDgqVrLCPDvpu4sFAupe8bjTIcREbKpVkz2QZLxaf0+xFC/o7y8mxdf5bl0iW5belqTlCQ7EYqfZycSQj5R8fXX8tJ7NPTZRfe+PHyv/jGjW3bXr4su70tc5Kp1XJRqeR7b9JEfh6aNJGfEaNR7s8yJq+oSHY1Nm0qPz/V/XdYVQw3lWC4qcC+fcCsWcDq1bZ13brJaYFHjbK/eiYREdUZJSUyjDRqVLfn6mK4qQTDzXXs2AG8+KI8ndxSc9VoZAdqVJTtvGs/P/mnTZ8+8k8SIiIiJ2K4qQTDTRWlpcnh9F9+aX91zGupVPLCUP37y6Vr17r/5wEREdU6DDeVYLi5AX/9BaxbJzuSs7Jsy4EDcqKSa4WEyMATEyPPwLJM49uwIUMPERHdEIabSjDcOFhamhx1aRl5efRoxacQ+PnZJkCJjJTdXJGRcuBy8+YMPkREVCGGm0ow3DhZQYE8DWPfPjnj2ZEjsvJz9mzlU/C2bAn84x/y+gN9+sj9bNokz8lNSpLPv+ce4LHHgE6dXPVuiIiolmC4qQTDjUIKC+XkFH/9JWfGOn3aNkvW4cP2k2/4+MhwU1EFqHdvGXL69pXbWWY3KymR525yzA8RkdthuKkEw00tlJMjx/T8+KO8BoFl9uSWLeX0sAMGyEksFi0CVqy4/kxnvr5ytq/WreV4n9JdYI0by1nnSl9ds6BABqrSS1iYfE2GJCKiWoHhphIMN7Wc2SwHKgcEyCrMtc6fBz74QAadtDTbrGZeXjKI/P33jV8L4Fo6nW32qtBQ24xpHh5yadFCzmAWG2t/oSiLrCzb/O5ERFQjDDeVYLhxExVNgVpUJOf+t1yd7uRJWzfYmTO2ecP9/GRgadhQhiPLtK65uXLa08zMql2mG5ABJjYW6NBBTtdp6W7LypJTdTZvLqdwjY6W1aO0NPn4mTPyVquVU9927iyXmBg5/WxBgW0aVb1eBi1nTf1JRFTLMdxUguGmHjOb5Rzvvr7Xr6YUFckq0dmzcrFc9MVkkt1iRqO8mNKWLbZuNGczGOSZZe3ayet/abX2c8ir1bLK1KyZbSnv0r01IYR8v2fPym7DgADH7ZuIqBIMN5VguCGHEkJerGfLFnnbuLHtojTNmsnxRIcPyzPIDh+WVw0MC7PfJj9fXuPLsiQn2/Zv6dbKz7/xSwpbLkyk18vxRJarGTZsKH82mWwX4cnIkJWr0hccCgyU7+PYMblkZdn23aYNEBcnLzQUHW27gqJlsVxG3NPTVnWyXHUwL0++L61Wvp6vr+suvENEdQ7DTSUYbqjWy8+XQUCrtQWC4mJ5ttmRI3IuoePHZbDS620BqKjIdsnlM2dkUHEGlUpO1Fjd/Xt6yudWdklhLy95eQ/LpaQti4eHLQzl5clqlUolw5Bl0WptV5f08rIter3tVq22VbosVa9r/wv08pID0Js3l0vpq05mZMiux/R02XWZkSFvMzPlMeneXYa98i6tbrlsN7sWiW5Idb6/nXR9WyK6YQZD2XWenrJKUp0rtRcUyIqLZexOQYEcU2S5TLJl8fCQFZyQEHnr71/2Msc6ne31W7aUX/gZGfLSHNu2yeX0afkalhBiNNq359rKk0ol92O59LWlzQUFMkDUFo0ayfZlZlb9OVFRQMeO8jhajrNlHFdgIBAcbFtKXzo6MFD+Pq5csS1ZWfL3V1RkWwDbJb8tz7dUzby9bZUzHx9ZEbOcBajTyXCl0dhuPTzsK2slJfJ3m5Yml8uX5eWsW7Ysf5qFkhLZTj8/GTBLEwK4cEF24R4+LKuWPXrIyiXPRCQnYuWGiJzDZLJ9GRuN8tZstn356vW2LzijUQaBnBz5RWn5YrUsludZFr1efnFaxkFZXssy75Hl1lKdsYQ7s1k+11LxsnzZl5aTI6tkJ0/KLrrS1Gr7Lj1LKAwOllWzbdvkl3hdpNHIkGM0VjyY3mCQZwn6+dkqV5cv27YPCpIBJixMrjtwoPxQGBoqQ05srPy5dNjTam1j20wm+Xu0/D5OnJA/GwxAz57ybMXu3WWAA+TvNyNDjpfLz7eFOstnzmSyH6em0chxauWNwRPC9h4tAdByK4T8vFmCuRDyfQQFVR7ahJCfj99+k8vly/J4tmxpu23dmmPZKsBuqUow3BBRlQghr6d2+rQMQZYv4euNC8rKAnbskN2HgYG2MGS5ttrFi/ZL6QrZ5cvySz0wUFbQAgLkrZeXbINWKxezWYbA0s+zdNlZJrW0nP2Xk2O7reo0CWq1DG6hobINf/8tx4JZutaqQ62WFb/oaDkQfc+eGx8/Vtn+c3Nlleh682CVxzIOrmlTGWgtXbv5+dXbj14vq1zh4fJ3qNXKMGQZ/L9pU9WqkmFh8qSBtm1l4GncWC7h4bJ6ZjbL32d2tq066+lpey2tVn7eAgOrfyxqMYabSjDcEFG9JYT8YjSb7c/+K12F8PIqP8QVFckv/BMnZHAqXbkKCpJftBcuyC/v1FS53w4dZKixjFsC5Bfx7t3A1q3AwYP2QS8zUz7PMqeURiNDXWSkrbLRvLkMnVu2yKX0AHxABsjQUFmtsQS8vDz76SMsYdFovH6ACQyUx6u4WLatuFgGqtJhwmyWAbMq9Hp5iZnbbpNVo1On5DG1LBcuVG0/VREQII9XixYyMF2+bKuGpqfLwGWZt8vDQ76XBg1kgLIs3t62wGxZLJUvy9g1T0/ZFduihVyaN5f3S//eHYDhphIMN0REbuT8eTmmJzBQVjdCQ+WXbWlC2H+RW7qOhJChylKpSUmRXVzNmslAFRFR/gSd5TEaZVvOnwfOnZNhr/Q4KSGArl1ld1plU1FkZ8uzEo8elcupU/b7LSiQ26nVtrMMvbzsQ2pRkazsKaltW3kChAMx3FSC4YaIiOokIWQ3lEYjxxBVNr4nL09WtU6elAEpNVVW5EqfhWgJRZalsFBWdC5csC0FBfaD1C1nIVrO0rRMVXHqlG1s1MmT8hqAP/3k0LfPs6WIiIjcjUolB3NXhbe3nBm9fXvntqk8lkqZgjjhAhERETmOZVyTghhuiIiIyK0w3BAREZFbYbghIiIit8JwQ0RERG6F4YaIiIjcCsMNERERuRWGGyIiInIrDDdERETkVhhuiIiIyK0w3BAREZFbYbghIiIit8JwQ0RERG6F4YaIiIjciofSDXA1IQQAIDs7W+GWEBERUVVZvrct3+OVqXfhJicnBwAQERGhcEuIiIiounJycuDv71/pNipRlQjkRsxmM86fPw9fX1+oVCqH7js7OxsRERE4e/Ys/Pz8HLpvssdj7To81q7DY+06PNau46hjLYRATk4OwsPDoVZXPqqm3lVu1Go1mjRp4tTX8PPz4z8WF+Gxdh0ea9fhsXYdHmvXccSxvl7FxoIDiomIiMitMNwQERGRW2G4cSCdTodZs2ZBp9Mp3RS3x2PtOjzWrsNj7To81q6jxLGudwOKiYiIyL2xckNERERuheGGiIiI3ArDDREREbkVhhsiIiJyKww3DvLee+8hMjISer0ecXFx2L59u9JNqvPmzJmDrl27wtfXFw0bNsTQoUNx7Ngxu20KCwsxefJkBAcHw8fHB8OHD0daWppCLXYfc+fOhUqlwtSpU63reKwd59y5c7jvvvsQHBwMLy8v3HTTTdi5c6f1cSEEZs6ciUaNGsHLywvx8fH466+/FGxx3WQymfDcc88hKioKXl5eaNGiBV588UW7axPxWN+4TZs24c4770R4eDhUKhW+++47u8ercmwvXbqEsWPHws/PDwEBAZg4cSJyc3Nr3jhBNbZ06VKh1WrF4sWLxaFDh8SkSZNEQECASEtLU7ppddrAgQPFJ598Ig4ePCj27t0rBg8eLJo2bSpyc3Ot2zzyyCMiIiJCJCUliZ07d4ru3buLnj17Ktjqum/79u0iMjJSdOzYUTzxxBPW9TzWjnHp0iXRrFkzMWHCBLFt2zZx6tQpsXbtWnHixAnrNnPnzhX+/v7iu+++E/v27RN33XWXiIqKEgUFBQq2vO55+eWXRXBwsPjxxx/F6dOnxfLly4WPj4946623rNvwWN+4n3/+WcyYMUOsXLlSABCrVq2ye7wqx/aOO+4QMTEx4s8//xSbN28WLVu2FGPGjKlx2xhuHKBbt25i8uTJ1vsmk0mEh4eLOXPmKNgq95Oeni4AiI0bNwohhLhy5Yrw9PQUy5cvt25z5MgRAUBs3bpVqWbWaTk5OaJVq1Zi3bp1ol+/ftZww2PtOE8//bTo3bt3hY+bzWYRFhYmXn/9deu6K1euCJ1OJ77++mtXNNFtDBkyRDzwwAN26+655x4xduxYIQSPtSNdG26qcmwPHz4sAIgdO3ZYt/nll1+ESqUS586dq1F72C1VQ0VFRdi1axfi4+Ot69RqNeLj47F161YFW+Z+srKyAABBQUEAgF27dqG4uNju2Ldt2xZNmzblsb9BkydPxpAhQ+yOKcBj7Ujff/89YmNjce+996Jhw4a4+eab8eGHH1ofP336NFJTU+2Otb+/P+Li4nisq6lnz55ISkrC8ePHAQD79u3DH3/8gUGDBgHgsXamqhzbrVu3IiAgALGxsdZt4uPjoVarsW3bthq9fr27cKajZWZmwmQyITQ01G59aGgojh49qlCr3I/ZbMbUqVPRq1cvdOjQAQCQmpoKrVaLgIAAu21DQ0ORmpqqQCvrtqVLl2L37t3YsWNHmcd4rB3n1KlTWLBgARITE/Hss89ix44dePzxx6HVajF+/Hjr8Szv/xQe6+p55plnkJ2djbZt20Kj0cBkMuHll1/G2LFjAYDH2omqcmxTU1PRsGFDu8c9PDwQFBRU4+PPcEN1wuTJk3Hw4EH88ccfSjfFLZ09exZPPPEE1q1bB71er3Rz3JrZbEZsbCxeeeUVAMDNN9+MgwcPYuHChRg/frzCrXMvy5Ytw5dffomvvvoK7du3x969ezF16lSEh4fzWLs5dkvVUEhICDQaTZmzRtLS0hAWFqZQq9zLlClT8OOPP2L9+vVo0qSJdX1YWBiKiopw5coVu+157Ktv165dSE9PR+fOneHh4QEPDw9s3LgRb7/9Njw8PBAaGspj7SCNGjVCdHS03bp27dohJSUFAKzHk/+n1NxTTz2FZ555BqNHj8ZNN92E+++/H//6178wZ84cADzWzlSVYxsWFob09HS7x0tKSnDp0qUaH3+GmxrSarXo0qULkpKSrOvMZjOSkpLQo0cPBVtW9wkhMGXKFKxatQq///47oqKi7B7v0qULPD097Y79sWPHkJKSwmNfTQMGDMCBAwewd+9e6xIbG4uxY8daf+axdoxevXqVmdLg+PHjaNasGQAgKioKYWFhdsc6Ozsb27Zt47Gupvz8fKjV9l9zGo0GZrMZAI+1M1Xl2Pbo0QNXrlzBrl27rNv8/vvvMJvNiIuLq1kDajQcmYQQ8lRwnU4nlixZIg4fPiweeughERAQIFJTU5VuWp32f//3f8Lf319s2LBBXLhwwbrk5+dbt3nkkUdE06ZNxe+//y527twpevToIXr06KFgq91H6bOlhOCxdpTt27cLDw8P8fLLL4u//vpLfPnll8JgMIgvvvjCus3cuXNFQECAWL16tdi/f7+4++67eXryDRg/frxo3Lix9VTwlStXipCQEPHvf//bug2P9Y3LyckRe/bsEXv27BEAxLx588SePXvEmTNnhBBVO7Z33HGHuPnmm8W2bdvEH3/8IVq1asVTwWuTd955RzRt2lRotVrRrVs38eeffyrdpDoPQLnLJ598Yt2moKBAPProoyIwMFAYDAYxbNgwceHCBeUa7UauDTc81o7zww8/iA4dOgidTifatm0rFi1aZPe42WwWzz33nAgNDRU6nU4MGDBAHDt2TKHW1l3Z2dniiSeeEE2bNhV6vV40b95czJgxQxiNRus2PNY3bv369eX+Hz1+/HghRNWO7cWLF8WYMWOEj4+P8PPzEwkJCSInJ6fGbVMJUWqqRiIiIqI6jmNuiIiIyK0w3BAREZFbYbghIiIit8JwQ0RERG6F4YaIiIjcCsMNERERuRWGGyIiInIrDDdEVO+pVCp89913SjeDiByE4YaIFDVhwgSoVKoyyx133KF004iojvJQugFERHfccQc++eQTu3U6nU6h1hBRXcfKDREpTqfTISwszG4JDAwEILuMFixYgEGDBsHLywvNmzfHihUr7J5/4MAB3HrrrfDy8kJwcDAeeugh5Obm2m2zePFitG/fHjqdDo0aNcKUKVPsHs/MzMSwYcNgMBjQqlUrfP/9985900TkNAw3RFTrPffccxg+fDj27duHsWPHYvTo0Thy5AgAIC8vDwMHDkRgYCB27NiB5cuX47fffrMLLwsWLMDkyZPx0EMP4cCBA/j+++/RsmVLu9d44YUXMHLkSOzfvx+DBw/G2LFjcenSJZe+TyJykBpfepOIqAbGjx8vNBqN8Pb2tltefvllIYS8Ovwjjzxi95y4uDjxf//3f0IIIRYtWiQCAwNFbm6u9fGffvpJqNVqkZqaKoQQIjw8XMyYMaPCNgAQ//nPf6z3c3NzBQDxyy+/OOx9EpHrcMwNESnulltuwYIFC+zWBQUFWX/u0aOH3WM9evTA3r17AQBHjhxBTEwMvL29rY/36tULZrMZx44dg0qlwvnz5zFgwIBK29CxY0frz97e3vDz80N6evqNviUiUhDDDREpztvbu0w3kaN4eXlVaTtPT0+7+yqVCmaz2RlNIiIn45gbIqr1/vzzzzL327VrBwBo164d9u3bh7y8POvjW7ZsgVqtRps2beDr64vIyEgkJSW5tM1EpBxWbohIcUajEampqXbrPDw8EBISAgBYvnw5YmNj0bt3b3z55ZfYvn07Pv74YwDA2LFjMWvWLIwfPx7PP/88MjIy8Nhjj+H+++9HaGgoAOD555/HI488goYNG2LQoEHIycnBli1b8Nhjj7n2jRKRSzDcEJHi1qxZg0aNGtmta9OmDY4ePQpAnsm0dOlSPProo2jUqBG+/vprREdHAwAMBgPWrl2LJ554Al27doXBYMDw4cMxb948677Gjx+PwsJCvPnmm5g2bRpCQkIwYsQI171BInIplRBCKN0IIqKKqFQqrFq1CkOHDlW6KURUR3DMDREREbkVhhsiIiJyKxxzQ0S1GnvOiai6WLkhIiIit8JwQ0RERG6F4YaIiIjcCsMNERERuRWGGyIiInIrDDdERETkVhhuiIiIyK0w3BAREZFbYbghIiIit/L/sLHtRbQ6Q/4AAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the results\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(errors_train,'r-',label='train')\n",
    "ax.plot(errors_test,'b-',label='test')\n",
    "#ax.set_ylim(0,100); ax.set_xlim(0,n_epoch)\n",
    "ax.set_xlabel('Epoch'); ax.set_ylabel('MSE')\n",
    "ax.set_title('Train Loss %3.5f, Test Loss %3.5f'%(errors_train[-1],errors_test[-1]))\n",
    "ax.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Physics Informed Model\n",
    "model_pinn = MixtureOfExperts(\n",
    "    D_i, D_o, D_h,\n",
    "    n_experts=n_experts,\n",
    "    activation_fn=acivation_fun,\n",
    "    gating_hidden=gating_hidden_layer_width,\n",
    "    dropout_rate=dropout_rate,\n",
    "    use_bn=use_bn\n",
    ")\n",
    "model_pinn.apply(weights_init)\n",
    "model_pinn.to(device)\n",
    "\n",
    "lambda_bw_mon = 0.1\n",
    "lambda_IL_mon = 0.1\n",
    "lambda_vpiL   = 0.05\n",
    "lambda_smooth = 0.01\n",
    "\n",
    "optimizer = torch.optim.AdamW(model_pinn.parameters(), lr=learning_rate, weight_decay=weight_decay, betas=betas)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoch     0 | Train 0.577587 | Test 0.602078\n",
      "Epoch    10 | Train 0.078204 | Test 0.101025\n",
      "Epoch    20 | Train 0.053335 | Test 0.077052\n",
      "Epoch    30 | Train 0.046076 | Test 0.070806\n",
      "Epoch    40 | Train 0.043260 | Test 0.068207\n",
      "Epoch    50 | Train 0.040504 | Test 0.065854\n",
      "Epoch    60 | Train 0.038483 | Test 0.064187\n",
      "Epoch    70 | Train 0.037348 | Test 0.064485\n",
      "Epoch    80 | Train 0.037371 | Test 0.066420\n",
      "Epoch    90 | Train 0.034940 | Test 0.064381\n"
     ]
    }
   ],
   "source": [
    "# ==============================\n",
    "# Training Loop\n",
    "# ==============================\n",
    "\n",
    "for epoch in range(n_epoch):\n",
    "\n",
    "    model_pinn.train()\n",
    "\n",
    "    for x_batch, y_batch in data_loader:\n",
    "\n",
    "        x_batch = x_batch.to(device)\n",
    "        y_batch = y_batch.to(device)\n",
    "\n",
    "        x_batch.requires_grad_(True)\n",
    "\n",
    "        optimizer.zero_grad()\n",
    "\n",
    "        pred = model_pinn(x_batch)\n",
    "\n",
    "        data_loss = loss_function(pred, y_batch)\n",
    "\n",
    "        BW_pred  = pred[:, 0]\n",
    "        IL_pred  = pred[:, 1]\n",
    "        Vpi_pred = pred[:, 2]\n",
    "        L_batch  = x_batch[:, 7]\n",
    "\n",
    "        physics_loss_total = 0.0\n",
    "\n",
    "        # =========================================\n",
    "        # 1) dBW/dL <= 0\n",
    "        # =========================================\n",
    "        if lambda_bw_mon != 0:\n",
    "\n",
    "            grads = torch.autograd.grad(\n",
    "                BW_pred,\n",
    "                x_batch,\n",
    "                grad_outputs=torch.ones_like(BW_pred),\n",
    "                create_graph=True\n",
    "            )[0]\n",
    "\n",
    "            dBW_dL = grads[:, 7]\n",
    "\n",
    "            physics_bw = torch.mean(torch.relu(dBW_dL) ** 2)\n",
    "            physics_loss_total += lambda_bw_mon * physics_bw\n",
    "\n",
    "\n",
    "        # =========================================\n",
    "        # 2) dIL/dL >= 0\n",
    "        # =========================================\n",
    "        if lambda_IL_mon != 0:\n",
    "\n",
    "            grads = torch.autograd.grad(\n",
    "                IL_pred,\n",
    "                x_batch,\n",
    "                grad_outputs=torch.ones_like(IL_pred),\n",
    "                create_graph=True\n",
    "            )[0]\n",
    "\n",
    "            dIL_dL = grads[:, 7]\n",
    "\n",
    "            physics_il = torch.mean(torch.relu(-dIL_dL) ** 2)\n",
    "            physics_loss_total += lambda_IL_mon * physics_il\n",
    "\n",
    "\n",
    "        # =========================================\n",
    "        # 3) d(Vpi*L)/dL ≈ 0\n",
    "        # =========================================\n",
    "        if lambda_vpiL != 0:\n",
    "\n",
    "            VpiL = Vpi_pred * L_batch\n",
    "\n",
    "            grads = torch.autograd.grad(\n",
    "                VpiL,\n",
    "                x_batch,\n",
    "                grad_outputs=torch.ones_like(VpiL),\n",
    "                create_graph=True\n",
    "            )[0]\n",
    "\n",
    "            dVpiL_dL = grads[:, 7]\n",
    "\n",
    "            physics_vpiL = torch.mean(dVpiL_dL ** 2)\n",
    "            physics_loss_total += lambda_vpiL * physics_vpiL\n",
    "\n",
    "\n",
    "        # =========================================\n",
    "        # 4) Smoothness (second derivative of BW)\n",
    "        # =========================================\n",
    "        if lambda_smooth != 0:\n",
    "\n",
    "            grads1 = torch.autograd.grad(\n",
    "                BW_pred,\n",
    "                x_batch,\n",
    "                grad_outputs=torch.ones_like(BW_pred),\n",
    "                create_graph=True\n",
    "            )[0]\n",
    "\n",
    "            dBW_dL = grads1[:, 7]\n",
    "\n",
    "            grads2 = torch.autograd.grad(\n",
    "                dBW_dL,\n",
    "                x_batch,\n",
    "                grad_outputs=torch.ones_like(dBW_dL),\n",
    "                create_graph=True\n",
    "            )[0]\n",
    "\n",
    "            d2BW_dL2 = grads2[:, 7]\n",
    "\n",
    "            physics_smooth = torch.mean(d2BW_dL2 ** 2)\n",
    "            physics_loss_total += lambda_smooth * physics_smooth\n",
    "\n",
    "\n",
    "        loss = data_loss + physics_loss_total\n",
    "\n",
    "        loss.backward()\n",
    "        optimizer.step()\n",
    "\n",
    "\n",
    "    # =========================\n",
    "    # Evaluation\n",
    "    # =========================\n",
    "    model_pinn.eval()\n",
    "    with torch.no_grad():\n",
    "\n",
    "        pred_train = model_pinn(x_train)\n",
    "        pred_test  = model_pinn(x_test)\n",
    "\n",
    "        train_loss = loss_function(pred_train, y_train).item()\n",
    "        test_loss  = loss_function(pred_test, y_test).item()\n",
    "\n",
    "        errors_train[epoch] = train_loss\n",
    "        errors_test[epoch]  = test_loss\n",
    "\n",
    "        if epoch % 10 == 0 or epoch == 0:\n",
    "            print(f\"Epoch {epoch:5d} | Train {train_loss:.6f} | Test {test_loss:.6f}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAHHCAYAAABDUnkqAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAXVtJREFUeJzt3XlcVOXiP/DPDDAzLLLDIAqCW+6YomiaWuJuqWVaWS7tqaU/rrf0etOyb2nZ7Wrq1ZtldkuzNG1PM1xKM3fcd0VIZXNhh4GZ5/fH4wyMLILMzIHh8369zmvgzJkzzxwG5sOzqoQQAkREREROQq10AYiIiIhsieGGiIiInArDDRERETkVhhsiIiJyKgw3RERE5FQYboiIiMipMNwQERGRU2G4ISIiIqfCcENEREROheGGnML48eMRERGhdDGIiKgWYLghu1KpVFXatm3bpnRRrWzbtg0qlQrr1q1Tuii3VVhYiFdffRWhoaFwd3dHTEwMNm/eXOXHX7p0CaNGjYKvry+8vb0xbNgwnD9/3uqY/Px8PP3002jXrh18fHzg5eWFqKgoLFy4EEVFRZWe/9lnn4VKpcLQoUPLvf+7775Dp06doNPpEB4ejtmzZ6O4uLjMcZs3b0bPnj3h4eEBPz8/jBw5EomJiVV+naU58n2Zl5eH119/vcrn4nvvfLnHpqam4vnnn0ejRo2g0+kQERGBp59+2uqY119/vdyfo06nszquOu/nPn36VPj+cHNzq/JrJcdyVboA5Nw+++wzq+//97//YfPmzWX2t27dukbPs3z5cphMphqdo64aP3481q1bh6lTp6JFixZYuXIlBg8ejK1bt6Jnz56VPjYnJwf33XcfMjMz8Y9//ANubm7497//jd69eyMhIQEBAQEA5IfBsWPHMHjwYERERECtVuOPP/7A//t//w+7d+/G6tWryz3/vn37sHLlyjIfLmY///wzhg8fjj59+mDRokU4cuQI/u///g9paWlYunSp5bgffvgBw4YNQ6dOnTBv3jxkZWVh4cKF6NmzJw4ePIigoKBqXTNHvS8BGW7eeOMNAPKD0pk44r0HAMnJyejRowcA4IUXXkCjRo1w+fJl7Nmzp9xzL126FF5eXpbvXVxcrO6vzvt55syZeOaZZ6wen5ubixdeeAH9+/ev2oUixxNEDjRp0iRRlbddbm6uA0pTsa1btwoAYu3atYqW43Z2794tAIj58+db9uXn54tmzZqJ7t273/bx77zzjgAg9uzZY9l34sQJ4eLiImbMmHHbx0+ePFkAEFeuXClzn8lkEt27dxdPPfWUaNKkiRgyZEiZY9q0aSOioqJEUVGRZd/MmTOFSqUSJ06csDquefPmorCw0LIvISFBqNVqERcXd9ty3k5V35d3Ij09XQAQs2fPrtLxfO+Vfe8NGjRIREZGioyMjErPOXv2bAFApKenV/PVSJW9n0v77LPPBACxatWqO3oesj82S5Hi+vTpg3bt2mH//v3o1asXPDw88I9//AMA8O2332LIkCEIDQ2FVqtFs2bN8Oabb8JoNFqd49Y+N4mJiVCpVHjvvffw4YcfolmzZtBqtejSpQv27t1rs7KfP38ejzzyCPz9/eHh4YFu3brhxx9/LHPcokWL0LZtW0uTSnR0tNV/h9nZ2Zg6dSoiIiKg1WoRHByMfv364cCBA5U+/7p16+Di4oLnnnvOsk+n0+Hpp5/Grl27kJycfNvHd+nSBV26dLHsa9WqFfr27Yuvvvrqtq/ffM1v3LhR5r7PPvsMR48exVtvvVXuY48fP47jx4/jueeeg6trSSXyxIkTIYSwNMtcu3YNx48fx4gRI6DRaCzHRUVFoXXr1lizZs1ty3knTCYTFixYgLZt20Kn00Gv1+P555/H9evXrY7bt28fBgwYgMDAQLi7uyMyMhJPPfUUAPk+NNcqvfHGG5bmjNdff73G5asv772TJ0/i559/xt///ncEBASgoKDgtk2hQghkZWVBCFHpcbeq7P1c2urVq+Hp6Ylhw4ZV6/zkOGyWolrh6tWrGDRoEB599FE88cQT0Ov1AICVK1fCy8sLcXFx8PLywpYtWzBr1ixkZWVh/vz5tz3v6tWrkZ2djeeffx4qlQrvvvsuHnroIZw/f77G7eWpqam45557kJeXh5dffhkBAQH49NNP8eCDD2LdunUYMWIEANlk9vLLL2PkyJGYMmUKCgoKcPjwYezevRuPP/44AFnVvm7dOkyePBlt2rTB1atXsWPHDpw4cQKdOnWqsAwHDx5Ey5Yt4e3tbbW/a9euAICEhASEhYWV+1iTyYTDhw9bPohvffwvv/yC7OxsNGjQwLLfYDAgKysL+fn52LdvH9577z00adIEzZs3t3p8dnY2Xn31VfzjH/9ASEhIhWUHgOjoaKv9oaGhaNy4seX+wsJCAIC7u3uZc3h4eODYsWNISUmp8Hnu1PPPP4+VK1diwoQJePnll3HhwgUsXrwYBw8exM6dO+Hm5oa0tDT0798fQUFBmD59Onx9fZGYmIj169cDAIKCgrB06VK8+OKLGDFiBB566CEAQIcOHWpUtvr03vv1118BAHq9Hn379sWWLVvg4uKCfv36YenSpeUOJGjatClycnLg6emJ4cOH41//+pflb0ppVX0/l5aeno7Nmzdj9OjR8PT0rPA4UpjCNUdUz5RX/d+7d28BQCxbtqzM8Xl5eWX2Pf/888LDw0MUFBRY9o0bN040adLE8v2FCxcEABEQECCuXbtm2f/tt98KAOL777+vtJxVaRqYOnWqACB+//13y77s7GwRGRkpIiIihNFoFEIIMWzYMNG2bdtKn8/Hx0dMmjSp0mPK07ZtW3H//feX2X/s2LEKr6mZublkzpw5Ze5bsmSJACBOnjxptf+LL74QACxbdHS0OHz4cJnHT5s2TURGRlp+RuU1S82fP18AEElJSWUe36VLF9GtWzchhBBGo1H4+vqKvn37Wh2TkZEhPD09BQCxb9++Cl9nVdz6vvz999/LbXbYuHGj1f4NGzYIAGLv3r0VntsezVL16b338ssvW36XBw4cKL788ksxf/584eXlJZo1a2bVhL1gwQIxefJksWrVKrFu3ToxZcoU4erqKlq0aCEyMzPLPFdV38+lLVq0SAAQP/30U6XHkbLYLEW1glarxYQJE8rsL/3fenZ2NjIyMnDvvfciLy8PJ0+evO15R48eDT8/P8v39957LwBUOCKjOn766Sd07drVquOkl5cXnnvuOSQmJuL48eMAAF9fX/z111+VNof5+vpi9+7duHz5crXKkJ+fD61WW2a/uQNvfn5+pY8FUK3H33fffdi8eTPWrl2LF154AW5ubsjNzbU65vTp01i4cCHmz59f7rmr+vzm+9VqNZ5//nnEx8djxowZOHPmDPbv349Ro0bBYDDc9nXeibVr18LHxwf9+vVDRkaGZevcuTO8vLywdetWAPLnBsgOz7drKrGl+vTey8nJAQCEhITgxx9/xKhRozBt2jQsX74c586ds2pimzJlChYtWoTHH38cDz/8MBYsWIBPP/0UZ86cwX/+858yz1WV9/OtVq9ejaCgIPTr16/S40hZDDdUKzRq1MiqP4XZsWPHMGLECPj4+MDb2xtBQUF44oknAACZmZm3PW94eLjV9+agc2u/iTtx8eJF3HXXXWX2m0fYXLx4EQDw6quvwsvLC127dkWLFi0wadIk7Ny50+ox7777Lo4ePYqwsDB07doVr7/+epUCmLu7u6XZprSCggLL/ZU9FkC1Hq/X6xEbG4uRI0di6dKlGDp0KPr164eUlBTLMVOmTME999yDhx9++LZlr+z5Sz/3nDlz8PTTT+Pdd99Fy5YtER0dDVdXV8tQ4NIjY2zhzJkzyMzMRHBwMIKCgqy2nJwcpKWlAQB69+6Nhx9+GG+88QYCAwMxbNgwfPLJJ+W+JluqT+898+2oUaOgVpd8ZD3yyCNwdXXFH3/8UWk5H3/8cYSEhFiat0qryvu5tPPnz2PXrl0YPXq0VT8xqn0YbqhWKO8P4Y0bN9C7d28cOnQIc+bMwffff4/NmzfjnXfeAYAqDf2+dQiomahmR8OaaN26NU6dOoU1a9agZ8+e+Prrr9GzZ0/Mnj3bcsyoUaNw/vx5LFq0CKGhoZg/fz7atm2Ln3/+udJzN2zYEFeuXCmz37wvNDS0wsf6+/tDq9Xe8eMBYOTIkcjJycG3334LANiyZQs2btyIKVOmIDEx0bIVFxcjPz8fiYmJyMrKspS99HPd+vyln1uj0eCjjz7C5cuX8dtvv+HUqVPYtGkTMjMzoVarK+0jcSdMJhOCg4OxefPmcrc5c+YAgGU+ml27dmHy5Mm4dOkSnnrqKXTu3NlS46AkZ3jvmW9v7TPj4uKCgICAKv2jEhYWhmvXrt32uFvfz7cy1xKNGTPmtuciZTHcUK21bds2XL16FStXrsSUKVMwdOhQxMbGWjUzKalJkyY4depUmf3m5rImTZpY9nl6emL06NH45JNPkJSUhCFDhuCtt96y/JcKyA+LiRMn4ptvvsGFCxcQEBBQ4Ugjs44dO+L06dOWwGC2e/duy/0VUavVaN++Pfbt21fmvt27d6Np06ZWnYnLY246MNeiJSUlAQAeeughREZGWrZLly5hy5YtiIyMxIoVK6zKduvzX758GX/99Ve5Zdfr9bj33nvRsmVLGI1GbNu2DTExMTavuWnWrBmuXr2KHj16IDY2tswWFRVldXy3bt3w1ltvYd++fVi1ahWOHTtmGcWlUqlsWjagfr33OnfuDEBO+FeawWBARkbGbec4EkJYjVqrzK3v51utXr0azZo1Q7du3W57LlIWww3VWuZal9K1LAaDody2cyUMHjwYe/bswa5duyz7cnNz8eGHHyIiIgJt2rQBIEeClabRaNCmTRsIIVBUVASj0Vjmj2lwcDBCQ0Nv27wxcuRIGI1GfPjhh5Z9hYWF+OSTTxATE2M1WiUpKalMP6WRI0di7969Vh8yp06dwpYtW/DII49Y9mVkZJRb2/XRRx8BKBnxdP/992PDhg1ltqCgIERHR2PDhg144IEHAABt27ZFq1at8OGHH1oN7V+6dClUKhVGjhxZ6Wt/7733cOXKFfztb3+r9Lg7MWrUKBiNRrz55ptl7isuLrYMFb5+/XqZ62L+UDf/7Dw8PADcfnhxddSn916fPn0QHByMVatWWQWylStXwmg0WvV9SU9PL1POpUuXIj09HQMHDrTsq+r7ubSDBw/ixIkTllFmVLux0ZBqrXvuuQd+fn4YN24cXn75ZahUKnz22WcObVL6+uuvy+24PG7cOEyfPh1ffPEFBg0ahJdffhn+/v749NNPceHCBXz99deW/gH9+/dHSEgIevToAb1ejxMnTmDx4sUYMmQIGjRogBs3bqBx48YYOXIkoqKi4OXlhV9//RV79+7Fv/71r0rLFxMTg0ceeQQzZsxAWloamjdvjk8//RSJiYn4+OOPrY4dO3Ystm/fbnX9Jk6ciOXLl2PIkCGYNm0a3Nzc8P7770Ov11uFhs8//xzLli3D8OHD0bRpU2RnZ2PTpk3YvHkzHnjgAdx///0AZB+nW/s5AcDUqVOh1+sxfPhwq/3z58/Hgw8+iP79++PRRx/F0aNHsXjxYjzzzDNWswN//vnn+Prrr9GrVy/L9fnqq6/wzDPPlOnbM378eMvP4U7XG+vduzeef/55zJ07FwkJCejfvz/c3Nxw5swZrF27FgsXLsTIkSPx6aef4j//+Q9GjBiBZs2aITs7G8uXL4e3tzcGDx4MQDa5tmnTBl9++SVatmwJf39/tGvXDu3atau0DHzvSVqtFvPnz8e4cePQq1cvPPnkk0hKSsLChQtx7733WobXA7LGavTo0Wjfvj10Oh127NiBNWvWoGPHjnj++ectx1X1/VzaqlWrALBJqs5QaJQW1VMVDQWvaLjqzp07Rbdu3YS7u7sIDQ0Vr7zyiti0aZMAILZu3Wo5rqKh4KVnTzVDFYblmofjVrSZh+CeO3dOjBw5Uvj6+gqdTie6du0qfvjhB6tz/fe//xW9evUSAQEBQqvVimbNmom///3vlqGphYWF4u9//7uIiooSDRo0EJ6eniIqKkr85z//qbSMZvn5+WLatGkiJCREaLVa0aVLF7Fx48Yyx5mH3N8qOTlZjBw5Unh7ewsvLy8xdOhQcebMGatj9u7dKx555BERHh4utFqt8PT0FJ06dRLvv/++1ezCFalohmIh5HDqjh07Cq1WKxo3biz++c9/CoPBYHXM7t27Ra9evYSfn5/Q6XQiKipKLFu2TJhMpjLne/jhh4W7u7u4fv36bctlVtEMxR9++KHo3LmzcHd3Fw0aNBDt27cXr7zyirh8+bIQQogDBw6Ixx57zHJdgoODxdChQ8sMTf/jjz9E586dhUajue37j++9M2WOE0IO246KihJarVbo9XoxefJkkZWVZXXMM888I9q0aSMaNGgg3NzcRPPmzcWrr75a5rjqvp+NRqNo1KiR6NSpU5WuCylPJYQD/w0mIrIzvV6PsWPHVmmSRyJyTgw3ROQ0jh07hu7du+P8+fMIDAxUujhEpBCGGyIiInIqHC1FREREToXhhoiIiJwKww0RERE5FYYbIiIicir1bhI/k8mEy5cvo0GDBnaZFp2IiIhsTwiB7OxshIaGWi2iWp56F24uX75sNS04ERER1R3Jyclo3LhxpcfUu3BjXowtOTkZ3t7eCpeGiIiIqiIrKwthYWG3XdAXqIfhxtwU5e3tzXBDRERUx1SlSwk7FBMREZFTYbghIiIip8JwQ0RERE6l3vW5ISIisiej0YiioiKli1EnaTSa2w7zrgqGGyIiIhsQQiAlJQU3btxQuih1llqtRmRkJDQaTY3Ow3BDRERkA+ZgExwcDA8PD04UW03mSXavXLmC8PDwGl0/hhsiIqIaMhqNlmATEBCgdHHqrKCgIFy+fBnFxcVwc3O74/OwQzEREVENmfvYeHh4KFySus3cHGU0Gmt0HoYbIiIiG2FTVM3Y6vox3BAREZFTYbghIiIim4iIiMCCBQuULoby4WbJkiWIiIiATqdDTEwM9uzZU+nxN27cwKRJk9CwYUNotVq0bNkSP/30k4NKS0RE5Fz69OmDqVOn2uRce/fuxXPPPWeTc9WEoqOlvvzyS8TFxWHZsmWIiYnBggULMGDAAJw6dQrBwcFljjcYDOjXrx+Cg4Oxbt06NGrUCBcvXoSvr6/jC1+mbEBqKiAEEB6udGmIiIhsQwgBo9EIV9fbR4agoCAHlOj2FK25ef/99/Hss89iwoQJaNOmDZYtWwYPDw+sWLGi3ONXrFiBa9eu4ZtvvkGPHj0QERGB3r17IyoqysElL2vPHhlqYmOVLgkREVHVjB8/Htu3b8fChQuhUqmgUqmwcuVKqFQq/Pzzz+jcuTO0Wi127NiBc+fOYdiwYdDr9fDy8kKXLl3w66+/Wp3v1mYplUqFjz76CCNGjICHhwdatGiB7777zu6vS7FwYzAYsH//fsSWSgNqtRqxsbHYtWtXuY/57rvv0L17d0yaNAl6vR7t2rXD22+/XemQscLCQmRlZVlt9qDVmp/PLqcnIqK6RAggN1eZTYgqF3PhwoXo3r07nn32WVy5cgVXrlxBWFgYAGD69OmYN28eTpw4gQ4dOiAnJweDBw9GfHw8Dh48iIEDB+KBBx5AUlJSpc/xxhtvYNSoUTh8+DAGDx6MMWPG4Nq1azW6vLejWLNURkYGjEYj9Hq91X69Xo+TJ0+W+5jz589jy5YtGDNmDH766SecPXsWEydORFFREWbPnl3uY+bOnYs33njD5uW/lU4nbwsK7P5URERU2+XlAV5eyjx3Tg7g6VmlQ318fKDRaODh4YGQkBAAsHwGz5kzB/369bMc6+/vb9VS8uabb2LDhg347rvvMHny5AqfY/z48XjssccAAG+//TY++OAD7NmzBwMHDqz2S6sqxTsUV4fJZEJwcDA+/PBDdO7cGaNHj8bMmTOxbNmyCh8zY8YMZGZmWrbk5GS7lI01N0RE5Eyio6Otvs/JycG0adPQunVr+Pr6wsvLCydOnLhtzU2HDh0sX3t6esLb2xtpaWl2KbOZYjU3gYGBcHFxQWpqqtX+1NRUS3q8VcOGDeHm5gYXFxfLvtatWyMlJQUGg6Hchba0Wi205uRhR6y5ISIiCw8PWYOi1HPbgOcttT/Tpk3D5s2b8d5776F58+Zwd3fHyJEjYTAYKj3PrcsoqFQqmEwmm5SxIoqFG41Gg86dOyM+Ph7Dhw8HIGtm4uPjK6ze6tGjB1avXg2TyWRZEv306dNo2LBhjVcQranSNTdCAJykkoioHlOpqtw0pDSNRlOl5Q527tyJ8ePHY8SIEQBkTU5iYqKdS3dnFG2WiouLw/Lly/Hpp5/ixIkTePHFF5Gbm4sJEyYAAMaOHYsZM2ZYjn/xxRdx7do1TJkyBadPn8aPP/6It99+G5MmTVLqJViYa24AOSyciIioLoiIiMDu3buRmJiIjIyMCmtVWrRogfXr1yMhIQGHDh3C448/bvcamDulaLgZPXo03nvvPcyaNQsdO3ZEQkICNm7caOlknJSUhCtXrliODwsLw6ZNm7B371506NABL7/8MqZMmYLp06cr9RIsSocbNk0REVFdMW3aNLi4uKBNmzYICgqqsA/N+++/Dz8/P9xzzz144IEHMGDAAHTq1MnBpa0alRDVGDPmBLKysuDj44PMzEx4e3vb7LxCADdbypCaCpQzByERETmpgoICXLhwAZGRkdCV/m+XqqWy61idz+86NVqqNlOpSvrdsOaGiIhIOQw3NsTh4ERERMpjuLEhDgcnIiJSHsONDbHmhoiISHkMNzbEmhsiIiLlMdzYyqlT0KbLpR1Yc0NERKQchhtbuXYNuhspAFhzQ0REpCSGG1vR6aCDTDUMN0RERMphuLEVrRZayPYoNksREREph+HGVlhzQ0REVCsw3NgKa26IiKgO6tOnD6ZOnWqz840fPx7Dhw+32fnuBMONrWi1JTU3ubdfOp6IiIjsg+HGVnS6kpqbPIYbIiKq/caPH4/t27dj4cKFUKlUUKlUSExMxNGjRzFo0CB4eXlBr9fjySefREZGhuVx69atQ/v27eHu7o6AgADExsYiNzcXr7/+Oj799FN8++23lvNt27bN4a/L1eHP6KxYc0NERDcJAeTlKfPcHh5yMeeqWLhwIU6fPo127dphzpw5AAA3Nzd07doVzzzzDP79738jPz8fr776KkaNGoUtW7bgypUreOyxx/Duu+9ixIgRyM7Oxu+//w4hBKZNm4YTJ04gKysLn3zyCQDA39/fXi+1Qgw3tuLqCi0MAIDCvGKFC0NERErKywO8vJR57pwcwNOzasf6+PhAo9HAw8MDISEhAID/+7//w9133423337bctyKFSsQFhaG06dPIycnB8XFxXjooYfQpEkTAED79u0tx7q7u6OwsNByPiUw3NiKSgWdaxFQDBTkmpQuDRER0R05dOgQtm7dCq9y0tm5c+fQv39/9O3bF+3bt8eAAQPQv39/jBw5En5+fgqUtnwMNzakcy2W4SaP4YaIqD7z8JA1KEo9d03k5OTggQcewDvvvFPmvoYNG8LFxQWbN2/GH3/8gV9++QWLFi3CzJkzsXv3bkRGRtbsyW2E4caGtK4y1BQWMNwQEdVnKlXVm4aUptFoYDSW9BXt1KkTvv76a0RERMDVtfyYoFKp0KNHD/To0QOzZs1CkyZNsGHDBsTFxZU5nxI4WsqGdG7yh1mQLxQuCRERUdVERERg9+7dSExMREZGBiZNmoRr167hsccew969e3Hu3Dls2rQJEyZMgNFoxO7du/H2229j3759SEpKwvr165Geno7WrVtbznf48GGcOnUKGRkZKCoqcvhrYrixIa2bueaG4YaIiOqGadOmwcXFBW3atEFQUBAMBgN27twJo9GI/v37o3379pg6dSp8fX2hVqvh7e2N3377DYMHD0bLli3xz3/+E//6178waNAgAMCzzz6Lu+66C9HR0QgKCsLOnTsd/prYLGVDOo0MN1x+gYiI6oqWLVti165dZfavX7++3ONbt26NjRs3Vni+oKAg/PLLLzYr351gzY0NaTWyxobLLxARESmH4caGdFoZbgoKqjh7EhEREdkcw40NaW+Gm0IDww0REZFSGG5sSKeTtwUMN0RERIphuLEhnVbeFhh4WYmI6iMhOFq2Jmx1/fgpbENanayxKSziZSUiqk/c3NwAAHlKrZbpJAwGuUaji4tLjc7DoeA2pPOQoaagqGY/FCIiqltcXFzg6+uLtLQ0AICHhwdUVV2amwAAJpMJ6enp8PDwqHBm5KpiuLEhS81NMcMNEVF9Y14F2xxwqPrUajXCw8NrHAwZbmxI5ylDTUExLysRUX2jUqnQsGFDBAcHK7LkgDPQaDRQq2vetYOfwjakdZc/EIPRFULIhdOIiKh+cXFxqXGfEaoZ9ny1IXPNDcBZiomIiJTCcGNDWs+SijCGGyIiImUw3NiQxqMk3HDxTCIiImUw3NiQyl0HHfIBMNwQEREpheHGlrRaaCHbo9gsRUREpAyGG1vS6aCDrLJhzQ0REZEyGG5siTU3REREimO4sSWtljU3RERECmO4sSWdjjU3RERECmO4sSXW3BARESmO4caWWHNDRESkOIYbW2LNDRERkeIYbmyJ4YaIiEhxDDe2xGYpIiIixTHc2FLpmpt8oXBhiIiI6ieGG1sqXXOTb1S4MERERPUTw40tla65yWG4ISIiUkKtCDdLlixBREQEdDodYmJisGfPngqPXblyJVQqldWm0+kcWNpKlF5+IY/hhoiISAmKh5svv/wScXFxmD17Ng4cOICoqCgMGDAAaWlpFT7G29sbV65csWwXL150YIkroVZDpy4CABTkMtwQEREpQfFw8/777+PZZ5/FhAkT0KZNGyxbtgweHh5YsWJFhY9RqVQICQmxbHq93oElrpzWtRgAUJhvUrgkRERE9ZOi4cZgMGD//v2IjY217FOr1YiNjcWuXbsqfFxOTg6aNGmCsLAwDBs2DMeOHavw2MLCQmRlZVlt9qRzlTU2BQw3REREilA03GRkZMBoNJapedHr9UhJSSn3MXfddRdWrFiBb7/9Fp9//jlMJhPuuece/PXXX+UeP3fuXPj4+Fi2sLAwm7+O0rRuN8NNHoeCExERKUHxZqnq6t69O8aOHYuOHTuid+/eWL9+PYKCgvDf//633ONnzJiBzMxMy5acnGzX8uluhpvCAoYbIiIiJbgq+eSBgYFwcXFBamqq1f7U1FSEhIRU6Rxubm64++67cfbs2XLv12q10Gq1NS5rVek0MtQUMNwQEREpQtGaG41Gg86dOyM+Pt6yz2QyIT4+Ht27d6/SOYxGI44cOYKGDRvaq5jVonWTfW24/AIREZEyFK25AYC4uDiMGzcO0dHR6Nq1KxYsWIDc3FxMmDABADB27Fg0atQIc+fOBQDMmTMH3bp1Q/PmzXHjxg3Mnz8fFy9exDPPPKPky7DQaW/W3BSqFC4JERFR/aR4uBk9ejTS09Mxa9YspKSkoGPHjti4caOlk3FSUhLU6pIKpuvXr+PZZ59FSkoK/Pz80LlzZ/zxxx9o06aNUi/BirkFrJDhhoiISBEqIUS96hySlZUFHx8fZGZmwtvb2+bn39LlVfTd9w7aNrqOo3/52fz8RERE9VF1Pr/r3Gip2k6rkzU2hQZeWiIiIiXwE9jGzMtcFRTx0hIRESmBn8A2pnWXl7SgyEXhkhAREdVPDDc2pnO/2SxVzHBDRESkBIYbG9N5sOaGiIhISQw3Nqb1kKGmyOQKE9fOJCIicjiGGxsz19wAnKWYiIhICQw3Nqb1LJkXkeGGiIjI8RhubMzNww0qyPaoggKFC0NERFQPMdzYmEqnhRayyoY1N0RERI7HcGNrOh10kFU2rLkhIiJyPIYbW9OW1Nww3BARETkew42tlaq5YbMUERGR4zHc2JpWy2YpIiIiBTHc2JqWHYqJiIiUxHBja+xQTEREpCiGG1tjzQ0REZGiGG5sjTU3REREimK4sTXW3BARESmK4cbWWHNDRESkKIYbW+MkfkRERIpiuLG1UvPcFBaYFC4MERFR/cNwY2ulm6VyjAoXhoiIqP5huLG10h2K8xhuiIiIHI3hxtbc3KAz97nJZbghIiJyNIYbW1OpoHUtBgAU5rPPDRERkaMx3NiB7ma4KchjuCEiInI0hhs70LrJUFNYIBQuCRERUf3DcGMHOlfZ16aAzVJEREQOx3BjB1qNrLHhJH5ERESOx3BjBzqNuVlK4YIQERHVQww3dmAONwVcOJOIiMjhGG7sQKuVt4UGlbIFISIiqocYbuxAp73Z56aQ4YaIiMjRGG7sQKuToaawiJeXiIjI0fjpawc6nbwtMPDyEhERORo/fe1A686aGyIiIqXw09cOdO7yshYUuShcEiIiovqH4cYOtOZwU+yqcEmIiIjqH4YbO9B5yhqbYpMLjEaFC0NERFTPMNzYgc6j5LIWciI/IiIih2K4sQOtR0lfG4YbIiIix2K4sQNXDw3UuLkyONeXIiIiciiGGztQ6bTQQlbZsOaGiIjIsRhu7EGngw6yyoY1N0RERI7FcGMPWtbcEBERKYXhxh60WtbcEBERKYThxh50OkvNDcMNERGRY9WKcLNkyRJERERAp9MhJiYGe/bsqdLj1qxZA5VKheHDh9u3gNVVquaGzVJERESOpXi4+fLLLxEXF4fZs2fjwIEDiIqKwoABA5CWllbp4xITEzFt2jTce++9DippNbDmhoiISDGKh5v3338fzz77LCZMmIA2bdpg2bJl8PDwwIoVKyp8jNFoxJgxY/DGG2+gadOmDixtFbHmhoiISDGKhhuDwYD9+/cjNjbWsk+tViM2Nha7du2q8HFz5sxBcHAwnn766ds+R2FhIbKysqw2u2OHYiIiIsUoGm4yMjJgNBqh1+ut9uv1eqSkpJT7mB07duDjjz/G8uXLq/Qcc+fOhY+Pj2ULCwurcblvq1SzFGtuiIiIHEvxZqnqyM7OxpNPPonly5cjMDCwSo+ZMWMGMjMzLVtycrKdSwnW3BARESnIVcknDwwMhIuLC1JTU632p6amIiQkpMzx586dQ2JiIh544AHLPpPJBABwdXXFqVOn0KxZM6vHaLVaaLVaO5S+Eqy5ISIiUoyiNTcajQadO3dGfHy8ZZ/JZEJ8fDy6d+9e5vhWrVrhyJEjSEhIsGwPPvgg7rvvPiQkJDimyakqStfc5AuFC0NERFS/KFpzAwBxcXEYN24coqOj0bVrVyxYsAC5ubmYMGECAGDs2LFo1KgR5s6dC51Oh3bt2lk93tfXFwDK7FdUqeUXCvJMAFyULQ8REVE9oni4GT16NNLT0zFr1iykpKSgY8eO2Lhxo6WTcVJSEtTqOtU1yGrhzMI8IxhuiIiIHEfxcAMAkydPxuTJk8u9b9u2bZU+duXKlbYvUE2VrrnJNSpcGCIiovqljlWJ1BEuLtCpDQCAwnyTwoUhIiKqXxhu7ETnUgzA3OeGiIiIHIXhxk60bjLUFBYw3BARETkSw42d6NxkXxsOBSciInIshhs7Kam5YbghIiJyJIYbO9FpZbjh8gtERESOxXBjJ1o3WWNTUKBSuCRERET1C8ONnei0MtwUGhQuCBERUT3DcGMn5rU6CwpZc0NERORIDDd2otPJ20IDLzEREZEj8ZPXTnTu8ragiJeYiIjIkfjJaydanWyOKmS4ISIicih+8tqJzl1e2oIirghORETkSAw3dqK9GW4KixluiIiIHInhxk50HvLSFptcUFyscGGIiIjqEYYbO9F6lNTYFBYqWBAiIqJ6huHGTnSeDDdERERKYLixE1cPDdS4uTI415ciIiJyGIYbe9FqoYNMNay5ISIichyGG3spFW5Yc0NEROQ4DDf2otNBC1llw5obIiIix2G4sRfW3BARESmC4cZeWHNDRESkCIYbe2HNDRERkSIYbuxFq7XU3DDcEBEROQ7Djb3odBwKTkREpACGG3thzQ0REZEiGG7shTU3REREimC4sRetFh7IAwDk5ChcFiIionqE4cZetFoEIgMAkJGhcFmIiIjqEYYbe9HpGG6IiIgUUK1w8+677yI/P9/y/c6dO1FYqkNJdnY2Jk6caLvS1WWla27ShcKFISIiqj+qFW5mzJiB7Oxsy/eDBg3CpUuXLN/n5eXhv//9r+1KV5fpdAhCOgCGGyIiIkeqVrgRQlT6PZVSquYmneGGiIjIYdjnxl40mlJ9blQKF4aIiKj+YLixF5UKgRrZhHf1ugomk8LlISIiqidcq/uAjz76CF5eXgCA4uJirFy5EoGBgQBg1R+HgEBtNmAATCYVbtwA/P2VLhEREZHzU4lqdJyJiIiASnX7JpYLFy7UqFD2lJWVBR8fH2RmZsLb29u+T6bXwyftNLLgg5Mngbvusu/TEREROavqfH5Xq+YmMTGxJuWqf252Ks6CDzIyGG6IiIgcgX1u7EmrLRkOzon8iIiIHKJa4WbXrl344YcfrPb973//Q2RkJIKDg/Hcc89ZTepX73GWYiIiIoerVriZM2cOjh07Zvn+yJEjePrppxEbG4vp06fj+++/x9y5c21eyDrLaq4bhctCRERUT1Qr3CQkJKBv376W79esWYOYmBgsX74ccXFx+OCDD/DVV1/ZvJB1FmtuiIiIHK5a4eb69evQ6/WW77dv345BgwZZvu/SpQuSk5NtV7q6jn1uiIiIHK5a4Uav11uGeRsMBhw4cADdunWz3J+dnQ03NzfblrAuK714JsMNERGRQ1Qr3AwePBjTp0/H77//jhkzZsDDwwP33nuv5f7Dhw+jWbNmNi9knVWqWYp9boiIiByjWvPcvPnmm3jooYfQu3dveHl5YeXKldBoNJb7V6xYgf79+9u8kHWWtzcCcQYAa26IiIgcpVrhJjAwEL/99hsyMzPh5eUFFxcXq/vXrl2LBg0a2LSAdVpQEILwBwCGGyIiIkepVrh56qmnqnTcihUrqlWIJUuWYP78+UhJSUFUVBQWLVqErl27lnvs+vXr8fbbb+Ps2bMoKipCixYt8Le//Q1PPvlktZ7TIQIDLc1SWVmAwQCUqugiIiIiO6hWuFm5ciWaNGmCu+++G9VYkqpSX375JeLi4rBs2TLExMRgwYIFGDBgAE6dOoXg4OAyx/v7+2PmzJlo1aoVNBoNfvjhB0yYMAHBwcEYMGCATcpkM0FB8MUNqGGECS7IyABCQ5UuFBERkXOr1sKZkyZNwhdffIEmTZpgwoQJeOKJJ+Bfw6WuY2Ji0KVLFyxevBgAYDKZEBYWhpdeegnTp0+v0jk6deqEIUOG4M0337ztsQ5dOPO774BhwxDsehXpxf44dAjo0MG+T0lEROSMqvP5Xa3RUkuWLMGVK1fwyiuv4Pvvv0dYWBhGjRqFTZs23VFNjsFgwP79+xEbG1tSILUasbGx2LVr120fL4RAfHw8Tp06hV69epV7TGFhIbKysqw2hwkKkjcqDgcnIiJylGovnKnVavHYY49h8+bNOH78ONq2bYuJEyciIiICOTk51TpXRkYGjEaj1cSAgJxPJyUlpcLHmTs0azQaDBkyBIsWLUK/fv3KPXbu3Lnw8fGxbGFhYdUqY40EBsobYxoAhhsiIiJHqNGq4Gq1GiqVCkIIGI1GW5Xptho0aICEhATs3bsXb731FuLi4rBt27Zyj50xYwYyMzMtm0NnUL5ZcxNoSgXAuW6IiIgcoVodigHZzLN+/XqsWLECO3bswNChQ7F48WIMHDgQanX1slJgYCBcXFyQmppqtT81NRUhISEVPk6tVqN58+YAgI4dO+LEiROYO3cu+vTpU+ZYrVYLrVZbrXLZjI8P4OqKwGI2SxERETlKtdLIxIkT0bBhQ8ybNw9Dhw5FcnIy1q5di8GDB1c72ACARqNB586dER8fb9lnMpkQHx+P7t27V/k8JpMJhYWF1X5+u1OpgMBAri9FRETkQNWquVm2bBnCw8PRtGlTbN++Hdu3by/3uPXr11f5nHFxcRg3bhyio6PRtWtXLFiwALm5uZgwYQIAYOzYsWjUqBHmzp0LQPahiY6ORrNmzVBYWIiffvoJn332GZYuXVqdl+I4gYEITGHNDRERkaNUK9yMHTsWKpXKpgUYPXo00tPTMWvWLKSkpKBjx47YuHGjpZNxUlKSVa1Qbm4uJk6ciL/++gvu7u5o1aoVPv/8c4wePdqm5bKZUhP5sc8NERGR/VVrnhtn4NB5bgBg1ChsXJuFQdiIqCggIcH+T0lERORs7DbPDd0B9rkhIiJyKIYbewsKsmqWql/1ZERERI7HcGNvpfrcGAxANec5JCIiompiuLG3oCB4IA86tRyqzqYpIiIi+2K4sbfAQKgABLlcA8BwQ0REZG8MN/ZmXoJByE7FHA5ORERkXww39mZePLNYLjHBmhsiIiL7YrixN3O44XBwIiIih2C4sTetFmjQgHPdEBEROQjDjSNwCQYiIiKHYbhxhFIT+bHmhoiIyL4YbhyhVM0Nww0REZF9Mdw4QlAQ+9wQERE5CMONI7DPDRERkcMw3DhCqT43164BRqPC5SEiInJiDDeOEBiIAFwFIFcFv35d4fIQERE5MYYbRwgKghuK4euSBYD9boiIiOyJ4cYRzLMUq+Timex3Q0REZD8MN45gXjzTlAaANTdERET2xHDjCOaaGxMXzyQiIrI3hhtH8PUFXFw41w0REZEDMNw4gkrFuW6IiIgchOHGUbgEAxERkUMw3DgKl2AgIiJyCIYbR2HNDRERkUMw3DhKqSUY2OeGiIjIfhhuHOWWDsVCKFweIiIiJ8Vw4yhBQWiESwCA3FyuL0VERGQvDDeOEhgID+QjRCMX0LxwQeHyEBEROSmGG0e5uQRDU9ckAMD580oWhoiIyHkx3DjKzSUYmprOAWC4ISIisheGG0cx19wUngDAcENERGQvDDeOEhAAAIgUrLkhIiKyJ4YbR9HpAC8vNIVMNQw3RERE9sFw40hBQZZwc/EiUFyscHmIiIicEMONIwUGIhSXoXE1wmgE/vpL6QIRERE5H4YbRwoKghoCkYHZANg0RUREZA8MN45kHg7uIyfyY7ghIiKyPYYbRzIPB3dPAcBwQ0REZA8MN450s+Ym0jUZAMMNERGRPTDcOJK55oZz3RAREdkNw40jmfvccJZiIiIiu2G4caSbNTeROUcAAFevAllZShaIiIjI+TDcONLNmhvvqxfMX+LCBQXLQ0RE5IQYbhypUSN5m52NpuFyemI2TREREdkWw40jeXoCDRsCAJoGZgJguCEiIrI1hhtHa94cAOe6ISIisheGG0dr0QIAEMnVwYmIiOyiVoSbJUuWICIiAjqdDjExMdizZ0+Fxy5fvhz33nsv/Pz84Ofnh9jY2EqPr3XMNTe5RwEw3BAREdma4uHmyy+/RFxcHGbPno0DBw4gKioKAwYMQFpaWrnHb9u2DY899hi2bt2KXbt2ISwsDP3798elS5ccXPI7ZA43V/cCABITAZNJwfIQERE5GZUQQihZgJiYGHTp0gWLFy8GAJhMJoSFheGll17C9OnTb/t4o9EIPz8/LF68GGPHjr3t8VlZWfDx8UFmZia8vb1rXP5qS0gA7r4bxYEhcL9xBcXFQHIy0Lix44tCRERUV1Tn81vRmhuDwYD9+/cjNjbWsk+tViM2Nha7du2q0jny8vJQVFQEf39/exXTtpo1AwC4ZqSgSZgRAJumiIiIbEnRcJORkQGj0Qi9Xm+1X6/XIyUlpUrnePXVVxEaGmoVkEorLCxEVlaW1aaoBg2AkBAAQNPgXAAMN0RERLakeJ+bmpg3bx7WrFmDDRs2QKfTlXvM3Llz4ePjY9nCwsIcXMpymPvdNJD9ihhuiIiIbEfRcBMYGAgXFxekpqZa7U9NTUXIzdqNirz33nuYN28efvnlF3To0KHC42bMmIHMzEzLlpycbJOy18jNcBOpugiA4YaIiMiWFA03Go0GnTt3Rnx8vGWfyWRCfHw8unfvXuHj3n33Xbz55pvYuHEjoqOjK30OrVYLb29vq01xN+e6aWo4CYDhhoiIyJZclS5AXFwcxo0bh+joaHTt2hULFixAbm4uJkyYAAAYO3YsGjVqhLlz5wIA3nnnHcyaNQurV69GRESEpW+Ol5cXvLy8FHsd1WJulrpxAAAXzyQiIrIlxcPN6NGjkZ6ejlmzZiElJQUdO3bExo0bLZ2Mk5KSoFaXVDAtXboUBoMBI0eOtDrP7Nmz8frrrzuy6HfOHG4u/Q4ASEkB8vIADw8lC0VEROQcFJ/nxtEUn+dGFgLw8QEA+PmYcCNThaNHgbZtlSkOERFRbVdn5rmpt7y9geBgAEDThvkA2O+GiIjIVhhulGJumvK5CoDhhoiIyFYYbpRyc8RUS10SAODwYSULQ0RE5DwYbpRys+amh5tc0fz335UsDBERkfNguFGKOdzkbIJKBZw5A1y5onCZiIiInADDjVJuNkv5JB5Cx45yF2tviIiIao7hRik3a26QkoJe3QwAgN9+U7A8REREToLhRik+PkBQEACgV7NLABhuiIiIbIHhRkk3a2/u9T0CADhyBLh2TckCERER1X0MN0q6GW6C0o6hdWu5a8cOBctDRETkBBhulHSzUzHOnkWvXvJLNk0RERHVDMONksydihluiIiIbIbhRknmcHPmDO69V3554ACQna1ckYiIiOo6hhslmcPNlSsI889FZCRgNAK7dilbLCIiorqM4UZJfn5AQID8+tw59O4tv9y+XbkiERER1XUMN0or1TTFfjdEREQ1x3CjNPOIqdOnLeFmzx4gP1+5IhEREdVlDDdKi4qStzt3omlTIDQUMBhkwCEiIqLqY7hRWr9+8nbrVqgMhWyaIiIiqiGGG6V16ACEhAB5ecDOnQw3RERENcRwozSVCujfX369aZMl3Pzxh2yeIiIiouphuKkNBgyQt7/8gtatSypyNm5UtlhERER1EcNNbWDud5OQAHV6Kh5/XH77v/8pVyQiIqK6iuGmNggKAjp1kl//8gvGjZNffv89cO2acsUiIiKqixhuagtz09SmTejQQY4QNxiAr75StlhERER1DcNNbWEON5s3AyYTxo6V37JpioiIqHoYbmqL7t0BLy8gLQ04dAiPPw6o1XIRzTNnlC4cERFR3cFwU1toNMB998mvN21CSEhJZc5nnylXLCIiorqG4aY2KdXvBoClaeqzzwCTSaEyERER1TEMN7WJOdzs3Ank5GDYMMDbG0hMBH7/XdGSERER1RkMN7VJ8+ZA06ZAURGwbRvc3YFHHpF3sWMxERFR1TDc1DYVNE2tXStnLSYiIqLKMdzUNqXWmQKAnj2BiAggOxv49lvlikVERFRXMNzUNvffD7i6yvHfR49CrS6pvXn7bdliRURERBVjuKltvL2BBx+UX3/0EQDg5ZeBwEDg6FHgX/9SsGxERER1AMNNbfTMM/L2s8+AggIEBADvvy93vfEGcP68ckUjIiKq7RhuaqP+/YHwcLlq5oYNAIAnngD69gUKCoCJEwEhFC4jERFRLcVwUxu5uABPPSW/Xr4cAKBSAUuXAlqt7Gu8Zo2C5SMiIqrFGG5qqwkTZKLZuhU4exYA0KIF8M9/yrunTgWuX1eueERERLUVw01tFR4ODBwov/74Y8vuV14BWreW62tOn65Q2YiIiGoxhpva7Nln5e0nn1jGgGs0wH//K3d/+CEX1SQiIroVw01tNnQooNcDqanADz9Ydt97LzBtmvx6/Hjg66+VKR4REVFtxHBTm7m5yfQCWDoWm73zjuyWYzIBjz0G/Pyz44tHRERUGzHc1HbmOW82bgSSky271WqZd0aPli1WDz0EbN+uUBmJiIhqEYab2q55c+C+++TENubONje5uMg+N0OHyvlvhg4Fdu9WqJxERES1BMNNXTBpkrx9/30gKcnqLjc3uWL4/fcDOTly/r8//1SgjERERLUEw01d8NBDQK9eQH4+EBdX5m6dTq4Y3qsXkJUlA87OnQqUk4iIqBZguKkLVCpg8WLZDvX118DmzWUO8fICfvpJtmBlZwMDBgC//aZAWYmIiBSmeLhZsmQJIiIioNPpEBMTgz179lR47LFjx/Dwww8jIiICKpUKCxYscFxBlda+PfDSS/Lrl14CDIYyh3h6yhHj/foBubnAoEFygmMiIqL6RNFw8+WXXyIuLg6zZ8/GgQMHEBUVhQEDBiAtLa3c4/Py8tC0aVPMmzcPISEhDi5tLfD663Lem1OngH//u9xDPDxkE9XAgUBenlxsMyxMzo3z5JPAa6/JGp5yshEREZFTUAmh3PrSMTEx6NKlCxYvXgwAMJlMCAsLw0svvYTpt1lbICIiAlOnTsXUqVOr9ZxZWVnw8fFBZmYmvL2977ToyvnsM2DsWJliTp6UyaUcBQUyzKxbV/5p/PyAESPkUPL77wdcXe1YZiIiohqqzue3YjU3BoMB+/fvR2xsbElh1GrExsZi165dShWr9nviCaBnT1kt87e/VXiYTidHUaWmytFTa9YA8+YBTz8tK3+uXwdWrJB9cxo1Aj79VI42JyIiqusU+389IyMDRqMRer3ear9er8fJkydt9jyFhYUoLCy0fJ+VlWWzcyvC3Lm4UyeZXj78EHjuuQoPDw6WW0xMyT6jUXY2/uorWbOTllayjMN//ws0bGj/l0FERGQvincotre5c+fCx8fHsoVV0IxTp0RFAbNny69ffBH4/vtqPdzFRY6qWroUuHwZmDtXLsj5/fdA27bAF1+wFoeIiOouxcJNYGAgXFxckJqaarU/NTXVpp2FZ8yYgczMTMuWXGoJgzrttddkG5PJJDvO3OHMfW5uwPTpwP79sjLo+nXg8ceBhx+WwYeIiKiuUSzcaDQadO7cGfHx8ZZ9JpMJ8fHx6N69u82eR6vVwtvb22pzCioVsGwZMHiwnNxv6FA5iuoOtWsn89Ebb8jOxRs2AK1by6cwmWxYbiIiIjtTtFkqLi4Oy5cvx6effooTJ07gxRdfRG5uLiZMmAAAGDt2LGbMmGE53mAwICEhAQkJCTAYDLh06RISEhJw9uxZpV6CslxdZceZrl2Bq1fl+O8rV+74dG5uwKxZshana1c52/GLL8qZj48ft2G5iYiI7EjRoeAAsHjxYsyfPx8pKSno2LEjPvjgA8Tc7P3ap08fREREYOXKlQCAxMREREZGljlH7969sW3btio9X50fCl6e9HSgRw/gzBk52d+2bYC/f41OaTQCS5YA//iHnBDQzQ145BHgqadkfx210/fWIiKi2qQ6n9+KhxtHc8pwAwDnz8sh4leuyGqXX38FGjSo8WmTkuS6nT/8ULKvSRNgwgTZN6d5c9lCRkREZE8MN5Vw2nADAMeOAb17yyaqXr2An3+Wk/3VkBDA3r3AJ58Aq1fL5iqzRo3kU5q3li0ZdoiIyPYYbirh1OEGkB1m7r9fJpCBA4FvvgG0WpudPj9fdjb+5BNg+3agqMj6/tBQueSDeWvcWC71kJ4u59PJzATuvhvw8bFZkYiIqB5guKmE04cbANixQ049nJcHPPggMGcO0KGDzatU8vLkCKvt2+X2559AqfkSAciWsexs633+/nKZrBdekH15iIiIbofhphL1ItwAss/NkCElK2Tq9XK58P795fDxgACbP2V+PvDHH0B8vHz6/ftLhpG7uABBQTJfmQd03XUX8N57sphsyiIiosow3FSi3oQbQI6amj9f3ubllezXaOQsfc89JzvK2ClZXL8u17YKCpILdarVQHEx8PHHcg7C9HR5XLdusmKpcWPZh6dxY9kn2tfXLsUiIqI6iOGmEvUq3JgVFgK7dgG//AL8+CNw+HDJfS1bAs8+K9uIvLwcVqTMTLnsw7//XVK5VJq3NzBtGjB1qk0GfRERUR3HcFOJehlubrV/P7B8ObBqFZCTI/eFhsplw8eMcegkNhcvAps2AZcuye2vv+REy4mJ8v7AQLk8xMSJchDYwYPAgQMyn0VEAHFxsraHiIicG8NNJRhuSsnJkatkzpsn58kBZHvQwoWyrUghJpOceHn2bOD0abnPza3syCxADgR74QUZgMxLkhUWArt3A1u3yhqixo2BsDC5NWnCVc+JiOoihptKMNyUo6AAWLAAeOutkpqcsWPlPj8/xYpVXAz8739yvaukJNkpuU0bOZS8XTvg22+BnTvlse7ustLp4kU5WCw/v+LzxsTIGp+HHpIrWBARUe3HcFMJhptKXLkCzJwJrFwpZ+5r1Ej2/h0wQNFiGQyyYikiAtDpSvYLIUdlvfaarKkpTa+X0/2EhsqmruRkuV26VDKCKzwceOklueZocjJw7pzcLl6Ug8latJBdklq0kOe7ehXIyCjZcnJkiCookLdFRTIsubmVbN7e8rHBwfI2JIR9iIiI7gTDTSUYbqpg1y5g3Di5VhUAPP+8HLPtwA7H1SEE8NNPsq9069Zy8sDWrcsfBJaaCvznP8DSpSWjtRwtKMg6OLVoATRtCjRrducjxISQA+IKCuRWWCg3tVoGQq1W3np6cm4hIqqbGG4qwXBTRXl5siPLokXy+6ZNZa3OqFG1NuRUR0GB7E/9wQeyX09kpAwXzZrJGqKMDLn/9GmZ8fLy5EoWQUGyk3NAgKyV0elkk5hOJ0OD0ShrcMzbjRtyZubUVLmZW/0q4u8vyxIQIGdx9vaWm7+/7CsUGiq3wEDg5ElZY2XeUlNv/7rd3OTCpyNGAMOGVd7/KDMTSEiQnbjz8uSxpTd/f2Wa9QwGOc2Avz+DGtV+QsgBEL/+Kmfh8POT/8T4+cm1+fR6pUtYdzDcVILhppq2bJGrZCYlye+9vIBHHwWeflp2XnGC2feEqPxlCCHDkLt7zZ8rOxs4e7YkNJ0+Lb8/f75q4aQqXFxk2NJoZBNcYaEs/61UKtlvvF07eZzJJMNZTg5w6JBsorudBg3kH2nz5u0t95lvi4tlk11eXkk/qODgkma64GC5kH3TpmV/BkajnBRywwYZ5C5fllvpGrfAwJKw1aiRbGps0kTehofL83t7lz13QYFshU1JkWHJ/FdQCPm85povg0FuXl4y2Jo3X1/52sz3Gwxyny3eI7aUlwf89huwebO8DQwERo+W4fbWJVCuXpV92PLySq5jw4ZlB0+ar1V1f/Vv93t2q8REOU1XfDzQqhVwzz1y69y54uuckSHD+MGD8vliY2UfvVtfw/XrcvqvpCT5HgkJkVtwsAz1ycnyvuRkuZJNeLj8pyMyUv7zk5cnf3/PnpW3aWnyH6O2beXWuLFsDl+9Gvj8c+Do0fLL6+YmZ+KYOVP+01JdQsjRpVeuAJ063dmyNkLIf44+/VQGME9P+R4PDpa3LVvKpvvw8Oqf29YYbirBcHMHsrJkO87HH5c0VQHyt7l3b+Dee+VW3icUVVlODnDhgvyjfv26vOxZWfKPbUaG/AN2+bK8TU+XHz4xMXKAW0yMDAkeHuXXpggha5IuXJDLjW3YULafUnnCw+UfTV9f+bzmzdZNemFhQJ8+cmvUCPj+e+Drr2X4qClXV1nLExgo356XL8vra2sqlbxeLVvK2bebNpXPbTLJ628yyQDYtq3sGH9rBajJJH/OV6/K90JubsmtORyag6KLi/yQbdpUbkFB8tiTJ4Hjx+W2d68MK+XNI6XVyg+sfv1krcJvv5X/AezmJj90TSbrkOruLl+n+bU2aybLfvZsyZaaKgOg0VjSzy0yEujRQ4aUHj3ktXBxsX7O48flAM7Vq+VjyytTRIT8EDZvarVcNzg5uezxer3sNnjfffKfiVtnT7eHW5ed0Wjk5PDu7vK9d+OGvF7mKS90OmDyZODVV+X7tDzFxfK9e/68/N3duVOG/6tX5f1qtfxd7dNH/llu0ED2MTT/U3DjRsk/AyEh8rqYQ82pU7d/TR07ytreBx6QNcvm8F9YCFy7VtJn0by1by9rx22J4aYSDDc1IATw++8y5KxdW3ZIUsOG8reqTx/5l6RFC4YdOzGZaj4d0aVLsp9Sero8l3nTaktGpVW0Soe5ye369ZLtxg35B90cynJy5Ie7u3vJJoT8L9e8Xb4sa4nKG+YPyP9Ehw+X2blRo5JmOT8/+Qe1dOD66y/53/bFiyX/defmVvz6tVr5R16nk29T8+biIu8zb25u8nWVXvz1Vmp19T8smzSRvyJZWfI6pKTID7A74e5e8QjB8HAZYvr2lR+Mq1YBJ06Uf2zr1vID8OJF+f4oL1zYkkYjf8YNGsjN1VUGD7N+/eRE6hcvyg/ynTtvX8PZvLn8kC8okLU+Fb0H7rpL1lpmZMhzpqTI97BOVzJ1RHi4LFdSkvzH4Pz5kqblRo3kc7VoIcPlmTMyYJ0+XXLdevcGnnhCTghf3sDTbdtkrc0ff8jvvbxkYNRq5bXRamV4uHhRvp/L+3nodLKWxVy5fifc3WUZH31U/gzM7/XUVNkFc+fO6r+/27UDjhy58zKVh+GmEgw3NpKVJcdc//abDDx795b9hGrYEIiOLmmv8PGRVQADB8p3PhFkbcAff8g/9OamgthY4JFH5AeyRnPn5y4okP/Zmke6mUwlfZd8fe8sexsM8u2v0cjNzU2Gm4wM+R+webt4UR5vDo0qlfwAPXas4hoplUqWy8tLbuaaCQ8PuZlDosFQ8mF76VJJU1FIiAymbdrIX7H775cfwKVfp7kPyOrVwJ49cumTXr1kgAwOLjmuuLjkP383N+uQmp1d8jpPn5blCAyUz9W8uazJadRIPs7FRW7m5zXXOPz5Z/l90FQq2Ww2Y4b881GaELK249Il65otg0GGgo4d5Z+a0j+rHTuAn3+Wz9msmXxP9e0rm47K+9m6uVX8vhBChmpz5/zyGAwy6Pj6Vm2CUSFk+f75T9mcVhlXVxm6OnaUNV89e8p/QjQaGe63b5e/Qzt2yPOa/xkwv9/NNcDmrVEjGb5GjrS+brdKT5f/CH37reypUFxc8v7XauVjzYMizFuLFvK9YEsMN5VguLGT/PySJcK3bZNxv7z6cED+5Rg1Si4N3qqVI0tJVCtcvSqbX86dk9k/NFSGLr2++p2kCwvlf/X+/nKrK4qLZblzcmRYysmRW/v28oOxvjGZ5P+I165Z9/lydbXuA3VrM159wnBTCYYbB8nPlw26p06VdBzJypJ/zX/6SR6jVsuZ92bOlHXEREREFWC4qQTDTS1w6JBcW+Hbb0v26fWyfrVjR3nbuzfHSBIRkQXDTSUYbmqR/ftlyPnpp5JOA2Zqtex4MWaMbIDntL5ERPUaw00lGG5qodxc2a3ePGPc3r3WPevc3WXQ0emsZ8iLiACGDJG9Ayvq3UdERE6B4aYSDDd1xNmzcjjHqlUlS4NXRKuVQ8/79ZMBqHRvPJNJ9sBzdS2Z3e7++zlai4iojmG4qQTDTR0jhGy++uMP2VRlXpHSxUXu/+EHOSa2ujp3BsaPBx57zHoyF/OMVEFBXDKciKgWYbipBMONkxFCTsv6ww9y+Ll5FjrzLFhqtZz5yjxValqanIvePFuamxsQFSVnoUtPlyO6AFnD06GDnA2sUyc5iYZKVTLVrMkkm9Nu3CjZCgpkJ2jzbHPmWy6ARERUYww3lWC4IaSnA198AaxceftZs2rKPE++eXazyEg545V5ZjZPTzk5SePGsraoptMOExE5KYabSjDckJWjR2X/nsDAkpXivL1lU9eBAyVbYmLJNLPmWy8vOe2nedNo5NSzly+XTO1a0USG5dFoZG1P48Yy8Pj4lMzq7Osrm8/Mm7+/LEdxsdyKimTNUem1EMxzyZuPDwiQNUtNm9bvmcCIqE5iuKkEww05jMkkw07p1QQTE2VzVuktI0Me56hfRQ8POZ9Qp05yTiHzKoTm+f7VajnZonluffNEjKUXXNLpSlajNG8NG8raqdBQrilGRDbHcFMJhhuqlYqKZE1PcrKs9blxQ87qnJlZUgtjXiTp6lXZ6RmQnZ7Nm1Yr5/L39ZW3Pj6yg3Tpx12+XPEKi7bi7l6yoqC/v/zevDiSRmM9nL+4WAYj8/LWTZvK8ufmlqzel54ur4l5RcyLF+U+Hx9ZG1U6YAUFlWyBgfIYb29ZhjsNXEVF8vpW9njzn1GGOiK7qc7nN4eDENUGbm5y8ZgmTez7PEajrI0xN7cdPCgDhHkFQvMqhBERckkM8xYQUDLEvrBQBqTr12Wtk3n76y/ZnJefL+ctutMlgd3cKl4m/E6p1TLkeHjI85tXvDSvegmUBJPiYrnYkXmJc4NB1lSV7iTu4yOXTL5ypWRJb51O9qkyb2FhJc1/5vDj6yvXU2vVqvxlos3Pf+WKDLrJyfLc5tUSIyOrHqAyMuTP1/xzPnRIvv5OneRowU6d5EJO5mXRiZwIa26IyJoQNavlSEyUyyKfOyfDQV5eyWYwWAcLFxcZEs6fl8enppacS6cr6QfVqJFcPdC86fXy3Oblvs1benrJlpEhj6mtf+KCg2UNl8lUsmpkTo4MjUZj+Y/x9paj+5o1k6+rdJ+rzEz5ms3XpKCg6mW5tS9ZcHDJFhAgn8tgkFtRkfzZBQfLn4NeL39GJpMMtnl58tZgkOczn1utLlkl1LyVN/O4EDIo794tly0/fVqGydK1cuHhMiA2alT2vWoyycCuVlsvc16a0ShHWd64IUdF1mQGdINBBtAGDWQ5Sw8KMK9qmpQkA6uLi6xF1OnkbViYDKxUJWyWqgTDDVEtZu6DFBAg+//UtEbBZJIftllZsiYmL09+OJs/pIuK5IepeQPkB1CDBnIzj2zLzJQ1NObO4pmZ8kPdvJx3w4ay7BculGyXLpUERfPrSEuTH6p//VV5uV1dSwKdXi/Pd+RI9TqoA7Jp0Ny3qmNHGZ4OHJBzRO3fL6+1ktzd5TVu0KCkz9f58zIIVIWnpww5EREy1CUlyWtrvk7BwTK8dOggjzt3Tgamffvk+wGQP5tWrYAuXYDoaFkOcxg312RqNNbb5cvAiRPyZ3nuXEkYValkwPHzkwEvJeX2r+Guu4AHHpDbPfdUfX6t/Hx5na5ckWVt2VKGpTsZcSmEfL+ePSu/d3WVAdbVFQgJke/FWoDhphIMN0SkuOxsWSNx/rz8sPTyKtn8/GSguXVEW1GR7Nx98KD8ADfPvG3ezH2QzP2QgoMrX5ZECFlLVFxcMn+T0SgDUFqa3FJTZWhwcSmpcdNoZI1EamrJlp4uy2DuX+XuLo8zP495bihzv6/Ll0vmlCqPi4tsMouJkbc5OSU1cmlpMuydPVtxDZc5TFb28ebpKa/17YJmVbi6lsyddSt3d9ncHBoqy5OfL2vV8vNlMCr9OF9fOVrSPE2Ep6c8d+mwlZsrr8GNG+W/ptat5abTlfTby8yUz+fnZ91H7do1GZqPHpXvhYqEhwM9ewI9egDduskAdf26fPz167Js5jBkrpXV6+WyOTbEcFMJhhsiolrAHKLM/ZtycuRtw4aytsnDo/LHGwwyHJ48KTuZBwaWNFuGhsoweOwYcPiw7G906pQMGTExQNeuQJs2Jc2i+/aVrGlnNJZ0gC/dCd68pEthoXwuc4ho1ark+Up/4Gu18vkCAiqugczMBH75Bfj+e7mA8NWr1buG5r5gWq0MezXpq+biIjv0u7qWdPYvKpK1TxWFyMp06yYnVrUhhptKMNwQEVGtYzTKEHb9uqydMdfUFBeXTPpp3oKCZAj08SkJTkVFMuwdPy6bzIxG67mydLqSQQDmPmmenrJmrH17GdJ0urLlysmR/Z927JDb/v0y8Pn7y5ogf39ZJnMYMt+2bg188IFNLxHDTSUYboiIiOqe6nx+c653IiIicioMN0RERORUGG6IiIjIqTDcEBERkVNhuCEiIiKnwnBDREREToXhhoiIiJwKww0RERE5FYYbIiIicioMN0RERORUGG6IiIjIqTDcEBERkVNhuCEiIiKnwnBDRERETsVV6QI4mhACgFw6nYiIiOoG8+e2+XO8MvUu3GRnZwMAwsLCFC4JERERVVd2djZ8fHwqPUYlqhKBnIjJZMLly5fRoEEDqFQqm547KysLYWFhSE5Ohre3t03PTdZ4rR2H19pxeK0dh9facWx1rYUQyM7ORmhoKNTqynvV1LuaG7VajcaNG9v1Oby9vfnL4iC81o7Da+04vNaOw2vtOLa41rersTFjh2IiIiJyKgw3RERE5FQYbmxIq9Vi9uzZ0Gq1ShfF6fFaOw6vtePwWjsOr7XjKHGt612HYiIiInJurLkhIiIip8JwQ0RERE6F4YaIiIicCsMNERERORWGGxtZsmQJIiIioNPpEBMTgz179ihdpDpv7ty56NKlCxo0aIDg4GAMHz4cp06dsjqmoKAAkyZNQkBAALy8vPDwww8jNTVVoRI7j3nz5kGlUmHq1KmWfbzWtnPp0iU88cQTCAgIgLu7O9q3b499+/ZZ7hdCYNasWWjYsCHc3d0RGxuLM2fOKFjiusloNOK1115DZGQk3N3d0axZM7z55ptWaxPxWt+53377DQ888ABCQ0OhUqnwzTffWN1flWt77do1jBkzBt7e3vD19cXTTz+NnJycmhdOUI2tWbNGaDQasWLFCnHs2DHx7LPPCl9fX5Gamqp00eq0AQMGiE8++UQcPXpUJCQkiMGDB4vw8HCRk5NjOeaFF14QYWFhIj4+Xuzbt09069ZN3HPPPQqWuu7bs2ePiIiIEB06dBBTpkyx7Oe1to1r166JJk2aiPHjx4vdu3eL8+fPi02bNomzZ89ajpk3b57w8fER33zzjTh06JB48MEHRWRkpMjPz1ew5HXPW2+9JQICAsQPP/wgLly4INauXSu8vLzEwoULLcfwWt+5n376ScycOVOsX79eABAbNmywur8q13bgwIEiKipK/Pnnn+L3338XzZs3F4899liNy8ZwYwNdu3YVkyZNsnxvNBpFaGiomDt3roKlcj5paWkCgNi+fbsQQogbN24INzc3sXbtWssxJ06cEADErl27lCpmnZadnS1atGghNm/eLHr37m0JN7zWtvPqq6+Knj17Vni/yWQSISEhYv78+ZZ9N27cEFqtVnzxxReOKKLTGDJkiHjqqaes9j300ENizJgxQghea1u6NdxU5doeP35cABB79+61HPPzzz8LlUolLl26VKPysFmqhgwGA/bv34/Y2FjLPrVajdjYWOzatUvBkjmfzMxMAIC/vz8AYP/+/SgqKrK69q1atUJ4eDiv/R2aNGkShgwZYnVNAV5rW/ruu+8QHR2NRx55BMHBwbj77ruxfPlyy/0XLlxASkqK1bX28fFBTEwMr3U13XPPPYiPj8fp06cBAIcOHcKOHTswaNAgALzW9lSVa7tr1y74+voiOjrackxsbCzUajV2795do+evdwtn2lpGRgaMRiP0er3Vfr1ej5MnTypUKudjMpkwdepU9OjRA+3atQMApKSkQKPRwNfX1+pYvV6PlJQUBUpZt61ZswYHDhzA3r17y9zHa20758+fx9KlSxEXF4d//OMf2Lt3L15++WVoNBqMGzfOcj3L+5vCa10906dPR1ZWFlq1agUXFxcYjUa89dZbGDNmDADwWttRVa5tSkoKgoODre53dXWFv79/ja8/ww3VCZMmTcLRo0exY8cOpYvilJKTkzFlyhRs3rwZOp1O6eI4NZPJhOjoaLz99tsAgLvvvhtHjx7FsmXLMG7cOIVL51y++uorrFq1CqtXr0bbtm2RkJCAqVOnIjQ0lNfaybFZqoYCAwPh4uJSZtRIamoqQkJCFCqVc5k8eTJ++OEHbN26FY0bN7bsDwkJgcFgwI0bN6yO57Wvvv379yMtLQ2dOnWCq6srXF1dsX37dnzwwQdwdXWFXq/ntbaRhg0bok2bNlb7WrdujaSkJACwXE/+Tam5v//975g+fToeffRRtG/fHk8++ST+3//7f5g7dy4AXmt7qsq1DQkJQVpamtX9xcXFuHbtWo2vP8NNDWk0GnTu3Bnx8fGWfSaTCfHx8ejevbuCJav7hBCYPHkyNmzYgC1btiAyMtLq/s6dO8PNzc3q2p86dQpJSUm89tXUt29fHDlyBAkJCZYtOjoaY8aMsXzNa20bPXr0KDOlwenTp9GkSRMAQGRkJEJCQqyudVZWFnbv3s1rXU15eXlQq60/5lxcXGAymQDwWttTVa5t9+7dcePGDezfv99yzJYtW2AymRATE1OzAtSoOzIJIeRQcK1WK1auXCmOHz8unnvuOeHr6ytSUlKULlqd9uKLLwofHx+xbds2ceXKFcuWl5dnOeaFF14Q4eHhYsuWLWLfvn2ie/fuonv37gqW2nmUHi0lBK+1rezZs0e4urqKt956S5w5c0asWrVKeHh4iM8//9xyzLx584Svr6/49ttvxeHDh8WwYcM4PPkOjBs3TjRq1MgyFHz9+vUiMDBQvPLKK5ZjeK3vXHZ2tjh48KA4ePCgACDef/99cfDgQXHx4kUhRNWu7cCBA8Xdd98tdu/eLXbs2CFatGjBoeC1yaJFi0R4eLjQaDSia9eu4s8//1S6SHUegHK3Tz75xHJMfn6+mDhxovDz8xMeHh5ixIgR4sqVK8oV2oncGm54rW3n+++/F+3atRNarVa0atVKfPjhh1b3m0wm8dprrwm9Xi+0Wq3o27evOHXqlEKlrbuysrLElClTRHh4uNDpdKJp06Zi5syZorCw0HIMr/Wd27p1a7l/o8eNGyeEqNq1vXr1qnjssceEl5eX8Pb2FhMmTBDZ2dk1LptKiFJTNRIRERHVcexzQ0RERE6F4YaIiIicCsMNERERORWGGyIiInIqDDdERETkVBhuiIiIyKkw3BAREZFTYbghonpPpVLhm2++UboYRGQjDDdEpKjx48dDpVKV2QYOHKh00YiojnJVugBERAMHDsQnn3xitU+r1SpUGiKq61hzQ0SK02q1CAkJsdr8/PwAyCajpUuXYtCgQXB3d0fTpk2xbt06q8cfOXIE999/P9zd3REQEIDnnnsOOTk5VsesWLECbdu2hVarRcOGDTF58mSr+zMyMjBixAh4eHigRYsW+O677+z7oonIbhhuiKjWe+211/Dwww/j0KFDGDNmDB599FGcOHECAJCbm4sBAwbAz88Pe/fuxdq1a/Hrr79ahZelS5di0qRJeO6553DkyBF89913aN68udVzvPHGGxg1ahQOHz6MwYMHY8yYMbh27ZpDXycR2UiNl94kIqqBcePGCRcXF+Hp6Wm1vfXWW0IIuTr8Cy+8YPWYmJgY8eKLLwohhPjwww+Fn5+fyMnJsdz/448/CrVaLVJSUoQQQoSGhoqZM2dWWAYA4p///Kfl+5ycHAFA/PzzzzZ7nUTkOOxzQ0SKu++++7B06VKrff7+/pavu3fvbnVf9+7dkZCQAAA4ceIEoqKi4Onpabm/R48eMJlMOHXqFFQqFS5fvoy+fftWWoYOHTpYvvb09IS3tzfS0tLu9CURkYIYbohIcZ6enmWaiWzF3d29Sse5ublZfa9SqWAymexRJCKyM/a5IaJa788//yzzfevWrQEArVu3xqFDh5Cbm2u5f+fOnVCr1bjrrrvQoEEDREREID4+3qFlJiLlsOaGiBRXWFiIlJQUq32urq4IDAwEAKxduxbR0dHo2bMnVq1ahT179uDjjz8GAIwZMwazZ8/GuHHj8PrrryM9PR0vvfQSnnzySej1egDA66+/jhdeeAHBwcEYNGgQsrOzsXPnTrz00kuOfaFE5BAMN0SkuI0bN6Jhw4ZW++666y6cPHkSgBzJtGbNGkycOBENGzbEF198gTZt2gAAPDw8sGnTJkyZMgVdunSBh4cHHn74Ybz//vuWc40bNw4FBQX497//jWnTpiEwMBAjR4503AskIodSCSGE0oUgIqqISqXChg0bMHz4cKWLQkR1BPvcEBERkVNhuCEiIiKnwj43RFSrseWciKqLNTdERETkVBhuiIiIyKkw3BAREZFTYbghIiIip8JwQ0RERE6F4YaIiIicCsMNERERORWGGyIiInIqDDdERETkVP4/g9UU6raNFKgAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the results\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(errors_train,'r-',label='train')\n",
    "ax.plot(errors_test,'b-',label='test')\n",
    "#ax.set_ylim(0,100); ax.set_xlim(0,n_epoch)\n",
    "ax.set_xlabel('Epoch'); ax.set_ylabel('MSE')\n",
    "ax.set_title('Train Loss %3.5f, Test Loss %3.5f'%(errors_train[-1],errors_test[-1]))\n",
    "ax.legend()\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "accelerator": "GPU",
  "colab": {
   "gpuType": "T4",
   "provenance": []
  },
  "jupytext": {
   "main_language": "python"
  },
  "kernelspec": {
   "display_name": "Python (mzmvenv)",
   "language": "python",
   "name": "mzmvenv"
  },
  "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.14.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}