guanjihuan/GitHub_collection_pykan
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity
GitHub_collection_pykan/tutorials/Example_9_singularity.ipynb
kindxiaoming 994d482d49 update tutorials
2024-07-13 22:17:48 -04:00

388 lines
47 KiB
Plaintext
Raw Blame History

{
"cells": [
{
"cell_type": "markdown",
"id": "134e7f9d",
"metadata": {},
"source": [
"# Example 9: Singularity"
]
},
{
"cell_type": "markdown",
"id": "2571d531",
"metadata": {},
"source": [
"Let's construct a dataset which contains singularity $f(x,y)=sin(log(x)+log(y))\n",
" (x>0,y>0)$"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "2075ef56",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"train loss: 5.01e-03 | test loss: 5.28e-03 | reg: 7.54e+00 : 100%|██| 20/20 [00:06<00:00, 3.16it/s]\n"
]
}
],
"source": [
"from kan import *\n",
"import torch\n",
"\n",
"# create a KAN: 2D inputs, 1D output, and 5 hidden neurons. cubic spline (k=3), 5 grid intervals (grid=5).\n",
"model = KAN(width=[2,1,1], grid=20, k=3, seed=0)\n",
"f = lambda x: torch.sin(2*(torch.log(x[:,[0]])+torch.log(x[:,[1]])))\n",
"dataset = create_dataset(f, n_var=2, ranges=[0.2,5])\n",
"\n",
"# train the model\n",
"model.fit(dataset, opt=\"LBFGS\", steps=20);"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "3f95fcdd",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAFICAYAAACcDrP3AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAtiklEQVR4nO3deVBUZ74+8OftbvYdRJBFpaGDuC9BXEDUqLhEb6Le0sp1bhnNxGtpMuNMJuY6cXQck2jGuUETU5PSqRtNcsMkQccYvOINyaACccFdWQVUQESBZmtolj6/P2L3D1xRDvT2fKpSU8U5ffg249tPv8t5j5AkSQIREZGMFOYugIiIbA/DhYiIZMdwISIi2TFciIhIdgwXIiKSHcOFiIhkx3AhIiLZMVyIiEh2DBciIpIdw4WIiGTHcCEiItkxXIiISHYMFyIikh3DhYiIZMdwISIi2anMXQCRNZAkCVVVVWhoaIC7uzv8/PwghDB3WUQWiz0XokfQarXYvn07NBoN/P39ERYWBn9/f2g0Gmzfvh1ardbcJRJZJMEnURI9WGpqKhYsWACdTgfg596LkbHX4urqiuTkZCQkJJilRiJLxXAheoDU1FTMmTMHkiTBYDA89DyFQgEhBFJSUhgwRB0wXIjuodVqERISgqampkcGi5FCoYCLiwtKS0vh7e3d8wUSWQHOuRDdY8+ePdDpdF0KFgAwGAzQ6XTYu3dvD1dGZD3YcyHqQJIkaDQaFBUV4UmahhACarUaBQUFXEVGBIYLUSd37tyBv79/t17v5+cnY0VE1onDYkQdNDQ0dOv19fX1MlVCZN0YLkQduLu7d+v1Hh4eMlVCZN0YLkQd+Pn5ITw8/InnTYQQCA8Ph6+vbw9VRmRdGC5EHQgh8Nprrz3Va19//XVO5hPdxQl9onvwPhei7mPPhege3t7eSE5OhhACCsWjm4jxDv19+/YxWIg6YLgQPUBCQgJSUlLg4uICIcR9w13Gn7m4uODQoUOYMWOGmSolskwMF6KHSEhIQGlpKRITE6FWqzsdU6vVSExMRFlZGYOF6AE450LUBZIk4ccff8Rzzz2HtLQ0TJkyhZP3RI/AngtRFwghTHMq3t7eDBaix2C4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQPUZrayvKysqQk5MDALh69Sqqq6thMBjMXBmR5eJjjokeQqvVIjk5GV988QUuX76M+vp6tLS0wNnZGf7+/oiLi8Py5csxceJEqFQqc5dLZFEYLkQPkJWVhTVr1uDChQuIjo7GnDlzMHz4cLi7u0Or1SI7OxsHDx5EYWEhFi1ahM2bN8Pf39/cZRNZDIYL0T2OHDmCpUuXwt3dHe+99x5mz56NlpYWJCUlQa/Xw9PTE4sXL0ZrayuSkpKwceNGDBkyBJ999hkCAgLMXT6RRWC4EHWQn5+PmTNnws3NDUlJSRg8eDCEECgqKsLo0aNRW1uLsLAwZGdnw8fHB5Ik4fjx43jppZcwefJk7N69G05OTuZ+G0Rmxwl9orva29vx7rvvoqamBh999JEpWB5FCIHY2Fi8//77OHDgAA4fPtxL1RJZNoYL0V2FhYU4ePAg5s+fj9jY2McGi5EQAi+88ALGjRuHXbt2oa2trYcrJbJ8XOJCdFdmZiYaGhqwYMEClJSUoLGx0XSstLQU7e3tAICWlhZcvnwZnp6epuNBQUGYP38+Nm7ciIqKCoSEhPR6/USWhOFCdFdubi5cXV2hVquxYsUKZGRkmI5JkgS9Xg8AKC8vx/Tp003HhBD4y1/+gmHDhkGn06G8vJzhQnaP4UJ0V1NTE1QqFZycnKDX69Hc3PzA8yRJuu9YW1sbXFxcOoUQkT1juBDd1bdvXzQ1NUGr1SImJgZubm6mY01NTcjMzDSFyIQJE0w3Tgoh0L9/f1RWVkKhUMDHx8dcb4HIYjBciO4aM2YMWltbcfLkSWzdurXTsaKiIkRHR6O2thYBAQH4+9//Dm9vb9NxIQTWrVuHwMBADokRgavFiEzGjh0LtVqNPXv2oLGxEUqlstN/RkIIKBQK088VCgVu3ryJb775BnPmzIGXl5cZ3wWRZWC4EN3l5+eH1atX48yZM9ixY0eXlxTr9Xr86U9/QlNTE1asWNHlJcxEtozDYkQdLF26FEePHsXWrVvh6uqKlStXwtnZGQCgUqmgUqlMvRhJklBfX4933nkHSUlJ+OCDDxAZGWnO8oksBrd/IbrH7du3sWrVKnz33XdISEjAmjVrEBUVhby8PBgMBjg6OiIiIgInT57Etm3bcO7cOWzatAkrV67sNHxGZM8YLkQP0NjYiF27dmHHjh24desW1Go1NBoNPDw8UFNTg7y8PJSXl2PMmDHYsGED4uPjoVBwlJnIiOFC9AgVFRVIS0tDeno6zp8/j5MnTyIuLg4TJ07EjBkzEBMTA1dXV3OXSWRxGC5EXXTq1CmMHTsWp06dwrPPPmvucogsGvvxRF2kVCpNy5CJ6NHYSoiISHYMFyIikh3DhYiIZMdwISIi2TFciIhIdgwXIiKSHcOFiIhkx3AhIiLZMVyIiEh2DBciIpIdw4WIiGTHcCEiItkxXIiISHYMFyIikh2f50LURZIkwWAwQKFQQAhh7nKILBp7LkRPgM9yIeoalbkLIJKLJEkoKChAVVWVuUvpFoVCgaFDh8LNzc3cpRA9NQ6Lkc0wGAxYtWoVQkND4e7ubu5yHqu9vR1KpfK+nx87dgzr16/H8OHDzVAVkTzYcyGb4uTkhOXLlyMgIKBb15EkCefPn0dqairGjBmDKVOmPDAInobBYMBXX32Fzz//HH/4wx8QHR1tmsORJAkNDQ3gdz6ydgwXontIkoQzZ85g/vz5KC0thbu7O3bv3o2FCxfKMpFfWlqKtWvXorS0FA4ODkhKSoKTk5MMlRNZDs5OEt2jtbUVW7ZsQWlpKZydndHQ0IDNmzejurq629eWJAnfffcdysrKAABpaWk4ffo0eypkcxguRPfIycnB999/D0dHR7zzzjvo27ev6WfdDYG2tjYcPnwYkiRBqVSisbERn3zyCdrb22WqnsgyMFyIOpAkCYcOHUJdXR2GDh2Kl19+GTNmzEB7ezuSk5O7HQKVlZU4c+YMlEolfvGLX0ClUiElJQWXLl1i74VsCsOFqIOWlhYcOXIEADBr1ix4enpi3rx5UCqV+Omnn1BZWfnU15YkCZcuXcLt27fh7++PN954A6NGjYJWq8Xf/vY3hgvZFIYLUQelpaW4fPkyHB0d8dxzz0EIgbFjx8Lf3x+3bt3CuXPnunX948ePo62tDYMHD0Z4eDheffVVKBQKJCcn4+rVq/K8CSILwHAhukuSJJw7dw41NTUICgrC4MGDAQCBgYEYNmwY2tracPz48afuYbS0tCAjIwMAEBsbCwcHB8ybNw+DBg3CrVu3sHfvXvZeyGYwXIg6yMrKgsFgwLBhw+Dj4wMAUKlUGD9+PADg5MmTaG1tfapr37x5E1euXIGDgwNiY2MhhICfnx+WLVsGIQS+/PLLbg27EVkShgvRXa2trThz5gwAICYmxnTTpBACMTExUKlUyMvLe6rtZSRJwtmzZ1FdXY3AwEAMGTLEdO2FCxciKCgI165dQ1pamnxviMiMGC5Ed1VVVaGgoAAqlQqjR4/udMPkoEGD4O3tjTt37jz13MixY8fQ3t6OESNGoE+fPqafBwUFYfr06TAYDPj666/R1tbW7fdCZG4MF6K7iouLUVVVBS8vL0RGRnY6FhAQgIEDB6KlpQUXLlx44rmR5uZmZGVlAQDi4uI6bSWjUCiwaNEijBgxAtHR0bznhWwCw4UI/3+ZsF6vR2hoKPr27dvpuJOTk2mC/+zZs098/dLSUuTn58PJyQkTJkzo1CsSQmDKlClIT0/HunXruBUM2QSGC9Fd58+fBwBERUXB2dm50zEhBEaNGgUAuHLlCvR6fZeva5xvqaurQ3BwMJ555pn7znFwcICHh0c3qieyLAwXIvw8mX/lyhUAwPDhw+/boFIIgSFDhkClUuHatWtPvM/Y8ePHYTAYMHLkSNMqNCJbxnAhAlBbW4vi4mIolUoMGTLkgbsfq9VqeHl5obq6GtevX+/ytfV6vWkV2vjx4/k0S7IL/FdOBKC8vBxVVVVwc3NDeHj4A8/x9/dHv3790NLSgvz8/C5P6ldXV6O4uBgqlQojRoyQZdt+IkvHcCECUFBQgKamJgQEBCAwMPCB57i4uECj0Zgm/7uquLgY1dXV8Pb2RkREhFwlE1k0hgvZPUmSkJOTA0mSEBYW9tBHJCsUCtOKsZycnC4tGTZeu6WlBaGhofD395e1diJLxXAhuydJkmkyPzIy8pGPMx48eDCEECgqKoJOp+vS9Y3XjoiIuG8VGpGtYriQzZEkCU1NTaiqqoLBYHjs+Xq93nTXvbFn8iBCCERERMDJyQm3bt3C7du3H3ttg8GAwsJCAD8HF+dbyF4wXMjmnDx5EvHx8Zg7dy7q6uoee35NTQ3Kysrg4OCAZ5555pEBEBwcDB8fH9TX16OkpOSx125qakJJSQmEEI+9NpEtYbiQzXFxcUFubi6uXLmCmzdvPvb80tJS1NTUwMPDAwMGDHjkuT4+PggJCUFbW1uXVoxVV1ejsrISDg4OCAsLe6L3QWTNGC5kc4KDg+Hv74/GxkYUFBQ88lxJklBQUICWlhYEBgZ22lDyQRwdHU0rvi5fvvzYWsrLy1FXVwcPDw8EBwd3/U0QWTmGC9kcDw8PDBw4EAaDwbQK7FGuXLkCSZIQHh4OV1fXR54rhDDNy+Tn5z92B+Pi4mK0tLQgICAAvr6+T/ZGiKwYw4VsjoODAwYNGgTg597Fo8LFYDAgNzcXwM8T7o+7e14IgUGDBkGhUKCkpAQNDQ0PPVeSJOTl5UGSJPTv3x8uLi5P8W6IrBPDhWyOcR8wAKYhr4dpamrC1atXO73mccLDw+Hs7Izbt2/j1q1bDz1PkiTk5+cD+HkZ8qOWOBPZGoYL2aTIyEioVCqUlpZCq9U+9LyqqircvHkTDg4O0Gg0XVrNFRQUBD8/PzQ2NqKoqOih57W0tJiOG3tSRPaC4UI2aeDAgfDw8EB1dTVKS0sfet7169dRW1sLb29vhIaGdunanp6e6N+/P9rb203DXg9SW1uL8vJyqFQqREREcBky2RWGC9kkf39/BAQEQK/Xm25ivJckScjNzUVrayuCg4O7POHu6OgIjUYD4NErxioqKlBTUwM3N7fHLnEmsjUMF7JJrq6uUKvVpq1dHta7MG5A+aRbswwdOhQAkJeX99A5neLiYjQ3N6NPnz7cU4zsDsOFbJJSqTTNc+Tk5DxwG5i2tjbk5OQA+DksujpsZVyOrFQqce3atQfO6RjvnzEYDAgNDX3oZphEtorhQjbL2Lswbqd/r7q6OhQVFUGhUDxRuAA/PzjMw8MDVVVVuHHjxgPPMS5xjoiIgEqleop3QGS9GC5kk4QQiIyMhJOTE8rKyh64yWRZWRkqKyvh5ub2wOfaP0pgYCCCgoKg1+uRm5t737Cb8YFiABAVFfX0b4TISjFcyGYNGDAAvr6+qKuru29S3ziZr9PpTEHxJFxdXREZGQlJknDu3Ln7jldXV6OkpAQqleqhj00msmUMF7JZvr6+CAsLQ1tbGy5cuHBf7+LMmTOQJAnPPPPME8+JKBQKjBw5EgBw/vz5+7aBuXbtGqqqquDp6fnQxyYT2TKGC9ksR0dHDB8+HACQnZ3dKVza2tpMPY5Ro0Y98d3zQgjT6/Lz81FdXW06JkkSzp8/D71ej/79+yMgIKD7b4bIyjBcyKZFR0dDCIGLFy922gespqYGubm5UCqVGDNmzFMNW0VFRcHLywu3b982PWwM+DlcMjIyAPwcXNxTjOwRw4VslhACI0eOhKurK65fv47r16+bjhUUFKCyshJeXl5d3lPsXoGBgVCr1dDr9Z16RvX19Thz5gyEEIiNjeV8C9klhgvZtLCwMISGhqKxsRGnT5+GJEmQJAknTpyAXq9HRETEE0/mG7m4uODZZ58FABw7dsx0L01JSQmuX78OV1dXjB49muFCdonhQjbNw8MD0dHRkCQJaWlpMBgMaG9vxz//+U8AwPjx45/ozvx7TZo0CQqFAtnZ2aiuroYkScjKyoJOp8PAgQMxcOBAed4IkZVhuJBNE0Jg+vTpUCgUyMzMRFVVFcrLy5GdnQ2VSoWpU6d269rR0dHw9fVFWVkZzpw5g7a2Nhw6dAiSJCE2NpZ35pPdYriQTTPOe/j7++PGjRtIT09HWloaKisrERwcbJrwf1qhoaGIjo5Ga2sr9u/fj/z8fGRlZcHBwQFz587lkBjZLe5JQTYvJCQE06ZNwxdffIEdO3ZAr9fDYDBg1qxZ3d5QUqVSYeHChUhNTcW3336L6upqVFdXY/jw4Rg/fjzDhewWey5k8xQKBVasWAEPDw9kZmYiOzsbPj4+WLZsWbc//IUQmDNnDgYPHoxbt24hOTkZSqUSr776Kry8vGR6B0TWh+FCNk8IgZiYGKxduxaenp7w8vLCunXrMHLkSFl6Fn369MGWLVsQHBwMZ2dnLF68GEuWLGGvhewah8XIpkiShJqaGjg4ONx3bPny5YiLi4MQAhqNBrW1tbL93rFjx+If//gHqqqqEBUVhdbW1k537T+J5uZm2eoiMheGC9kMIQQGDBiADz/88Im3c5HTwYMHu/X6pqYmDqmR1RPSwx7RR2RljDdI2gIhBIfVyKoxXIiISHac0CciItlxzoWoizp28jlkRfRo7LkQddHZs2ehVCpx9uxZc5dCZPEYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EHWBJEmoqakBANTU1IBPByd6NIYL0SNotVps374dGo0G06ZNgyRJmDZtGjQaDbZv3w6tVmvuEokskpD4FYzogVJTU7FgwQLodDoAD37MsaurK5KTk5GQkGCWGoksFcOF6AFSU1MxZ84cSJIEg8Hw0PMUCgWEEEhJSWHAEHXAcCG6h1arRUhICJqamh4ZLEYKhQIuLi4oLS2Ft7d3zxdIZAU450J0jz179kCn03UpWADAYDBAp9Nh7969PVwZkfVgz4WoA0mSoNFoUFRU9EQrwoQQUKvVKCgoMM3HENkzhgtRB3fu3IG/v3+3Xu/n5ydjRUTWicNiRB00NDR06/X19fUyVUJk3RguRB24u7t36/UeHh4yVUJk3RguRB34+fkhPDz8iedNhBAIDw+Hr69vD1VGZF0YLkQdCCHw2muvPdVrX3/9dU7mE93FCX2ie/A+F6LuY8+F6B7e3t5ITk6GEAIKxaObiPEO/X379jFYiDpguBA9QEJCAlJSUuDi4gIhxH3DXcafubi44NChQ5gxY4aZKiWyTAwXoodISEhAaWkpEhMToVarOx1Tq9VITExEWVkZg4XoATjnQtQFkiThxx9/xHPPPYe0tDRMmTKFk/dEj8CeC1EXCCFMcyre3t4MFqLHYLgQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBA9RmtrK8rKypCTkwMAuHr1Kqqrq2EwGMxcGZHl4mOOiR5Cq9UiOTkZX3zxBS5fvoz6+nq0tLTA2dkZ/v7+iIuLw/LlyzFx4kSoVCpzl0tkURguRA+QlZWFNWvW4MKFC4iOjsacOXMwfPhwuLu7Q6vVIjs7GwcPHkRhYSEWLVqEzZs3w9/f39xlE1kMhgvRPY4cOYKlS5fC3d0d7733HmbPno2WlhYkJSVBr9fD09MTixcvRmtrK5KSkrBx40YMGTIEn332GQICAsxdPpFFYLgQdZCfn4+ZM2fCzc0NSUlJGDx4MIQQKCoqwujRo1FbW4uwsDBkZ2fDx8cHkiTh+PHjeOmllzB58mTs3r0bTk5O5n4bRGbHCX2iu9rb2/Huu++ipqYGH330kSlYHkUIgdjYWLz//vs4cOAADh8+3EvVElk2hgvRXYWFhTh48CDmz5+P2NjYxwaLkRACL7zwAsaNG4ddu3ahra2thyslsnxc4kJ0V2ZmJhoaGrBgwQKUlJSgsbHRdKy0tBTt7e0AgJaWFly+fBmenp6m40FBQZg/fz42btyIiooKhISE9Hr9RJaE4UJ0V25uLlxdXaFWq7FixQpkZGSYjkmSBL1eDwAoLy/H9OnTTceEEPjLX/6CYcOGQafToby8nOFCdo/hQnRXU1MTVCoVnJycoNfr0dzc/MDzJEm671hbWxtcXFw6hRCRPWO4kN0rKSlBeno6MjIyoNPpoNVqERMTAzc3N9M5TU1NyMzMNIXIhAkTTDdOCiHQv39/VFZWoq2tDYWFhYiOjoazs7O53hKR2XEpMtmdGzdu4OjRo0hPT0d6ejquX79uCojCwkLs3LkTr7zySqfXFBUVITo6GrW1tRg4cCBOnz4Nb29v03EhBNatW4dt27ZBpVLB2dkZMTExmDRpEuLj4xEdHc0lymRXGC5k88rLy5Genm4KlOLiYgDA8OHDMWnSJEyaNAlxcXFob29HbGwsfHx8cPjw4U4T9g+7zwX4eZisvLwc8fHxmDt3LpYuXYpjx47h6NGjOHr0KGpra+Hi4oJx48aZwmbMmDFwdHQ0y9+DqDcwXMjm3Lp1q1OYFBYWAgCGDBli+nCPi4uDr6/vfa/duXMnfvvb3+Ltt9/GW2+9ZRr6elS4NDc349e//jUOHjyIH374AZGRkabrtbe34+LFi6Zajh8/jrq6Ori6umL8+PGIj49HfHw8Ro0aBQcHh1746xD1DoYLWb3bt2+begnp6enIy8sDAAwaNKhTmHRl76/GxkYsW7YMhw4dwh//+EesXLkSzs7OKC4uxtixY03DYidPnoS3tzfq6+vxzjvv4JNPPsEHH3yAl19++ZHXb2trw/nz503hl5GRgYaGBri7u2PChAmmsBkxYgQ3wySrxnAhq1NVVWUadkpPT8eVK1cAABqNxhQmkyZNeup9vm7fvo1Vq1bhu+++Q0JCAtasWYOoqCjk5eXBYDDA0dEREREROHnyJLZt24Zz585h06ZNWLlyJZRK5RP9rtbWVpw9e9YUNpmZmdDpdPD09MTEiRNNYTNs2LAnvjaROTFcyOJptVocO3bM9AF88eJFAIBare4UJkFBQbL9zsbGRuzatQs7duzArVu3oFarodFo4OHhgZqaGuTl5aG8vBxjxozBhg0bEB8fD4Wi+xtetLS0IDs72xScWVlZaG5uhpeXF+Li4kzvdejQobL8PqKewnAhi1NbW4vjx4+bwuTChQuQJAkDBgwwfbjGx8f3yo2KFRUVSEtLQ3p6OoqKitDc3AwfHx8MHToUM2bMQExMDFxdXXvs9+v1epw6dcoUNidOnIBer4ePj48pbOLj47u0DxpRb2K4kNnV19cjIyPD9AF67tw5GAwGBAcHmz484+PjMWDAALPW2d7eDkmSoFAozNZraG5uxokTJ0x/q5MnT6K1tRV9+vTpFDaRkZEMGzIrhgv1uoaGBmRlZZl6JmfOnEF7ezv69evXqWcSFhbGD8jH0Ol0+Omnn0xhc/r0abS1taFv376mv2N8fDwiIiL4t6RexXChHmf8ADSGSccPwI5hwg/A7mtoaDD9rdPT0zsFd8ewYXBTT2O4kOyMQzfGMOk4dGO8aZFDN72jrq7O1EvsOOQYEhLSKWzMPeRItofhQt1mnHQ2hknHSeeOYcJJZ/MzLpYwDqN1XCzRMWy4qzN1F8OFnphxuazx2/BPP/1kWi5rDBMul7UONTU1ncLGuMw7LCys02KKfv36mblSsjYMF3qsjjf6Ge+9MN7oFxsba5o34Y1+1q+qqsq0DLzjDaoRERGmoOnODapkPxgudJ+OW5Skp6ebntDo7u6OiRMnmoZPuEWJ7bt9+7bpBtaOW+tERkZ2Cps+ffqYuVKyNAwXMm2uaPwAycjIMG2uaNzvatKkSdxckVBRUdEpbIybgg4ePNgUNg/bFJTsC8PFDhkMBly6dMk0zn7s2DHU1tbC2dnZtFPvpEmTuC08PVZ5eXmnZ+MYH2cwbNgw07+juLi4Ts++IfvAcLEDkiThypUrpg+AY8eOoaamBk5OToiJiTF9CPCBVtRdD3sQ2/Dhw009m4kTJ8LLy8vcpVIPY7jYIEmSkJeX1ylM7ty5AwcHB4wdO9YUJjExMXwUL/WokpKSTmFTVlYGhUKBkSNHmsJmwoQJ8PDwMHepJDOGiw2QJAmFhYWmBnz06FFUVlZCpVLh2WefNYXJuHHjenSTRaJHkSQJxcXFpn+n6enpqKiogFKpxOjRo01hM378eLi5uZm7XOomhosVureRHj16FDdv3ryvkY4bNw7u7u7mLpfogbrypSg+Pr7Hd56mnsFwsRLXrl3r1AhLS0vvG14YP358p+e+E1mTe4dzjx49iqqqKjg6OiI6OprDuVaG4WIlRowYgYKCAk6Mkt0wGAzIycm5byHKZ599hoULF5q7PHoMhouVMBgMEEJwby6yW5IkQZIktgMrwXAhIiLZce8OmRgnJ6uqqsxdSrcoFAoMGTKEq3XoibENUEcMF5lIkoQdO3YgODgYLi4uaG1thZOTk9V13zMyMrBu3ToMGzbM3KWQlWEboI4YLjJycnLC6NGjsWPHDgDAF198YVUT7pIkobGxERwppafl5OSEMWPGIDExEUIIfP7551a19QvbgHwYLj0gMzMTjo6OuHPnjlWFC5EcJElCVlYWnJycUF1dbVXhQvLhk5xkFhAQAG9vbzQ2NqK0tNTc5RD1uj59+sDDwwM6nQ63bt0ydzlkJgwXmXl7e6Nfv35oa2tDQUEBu9dkd7y9vdGnTx+0trbi2rVrbAN2iuEiM0dHR4SHhwMAcnJyzFwNUe9zdnbGgAEDTHfck31iuMjMuIwRAK5cuYK2tjYzV0TUuxQKBTQaDQAgLy+PPRc7xXDpAUOHDoVSqcTVq1eh1WrNXQ5RrxJCYMiQIRBCoKCgADqdztwlkRkwXHpAZGQkPDw8cPv2bVy7ds3c5RD1uqioKDg5OaG0tBQVFRXmLofMgOHSAwIDA9G/f380Nzfj3LlzHBYguxMWFoaAgADU1dXh0qVLbAN2iOHSA1xdXTFy5EgAP9/zwoZF9sbb2xtDhw6FwWBARkaGucshM2C49AAhBGJjYyGEQHZ2Nurq6sxdElGvUiqVmDRpEgDg+PHjaGxsNHNF1NsYLj1ACIHo6Gh4enri+vXryM3NNXdJRL1KCIFJkybBzc0N+fn5bAN2iOHSQwYMGIBBgwahubkZ33//PYfGyO4888wziIqKgk6nw//93/+xDdgZhksPcXZ2xowZMwAAqampXI5JdsfFxQUJCQkAgP/93/9lG7AzDJceIoTAzJkz4ebmhsuXL+PixYv85kZ2RQiBWbNmwc3NDZcuXWIbsDMMlx4UFRWF0aNHo6mpCV999RUbFtmdwYMHY9SoUWhqakJycjLbgB1huPQgZ2dnLF68GAqFAt9++y3KysrMXRJRr3J2dsa//uu/QgiBgwcP8oZKO8Jw6UFCCMyePRsDBw5EeXk5ey9kd4xtIDg4GDdu3MC3337LNmAnGC49LCAgAP/2b/8GSZLw3//937h586a5SyLqVUFBQVi4cCEMBgP+9re/obq62twlUS9guPQwIQSWLFmC0NBQFBUVYc+ePTAYDOYui6jXCCHw8ssvo2/fvrhy5QqSkpLYe7EDDJdeEBoaildeeQUAsGvXLj5EjOyKEAIRERFYsmQJDAYDPvzwQz5EzA4wXHqBEAJLly7FkCFDcPPmTWzZsgUtLS3mLouo1ygUCqxcuRIRERG4du0atm7ditbWVnOXRT2I4dJL/P398dZbb8HJyQn79+/HgQMH+M2N7EpwcDDeeustODo64u9//ztSUlLYBmwYw6WXCCHw/PPPY+HChdDr9diwYQPy8/PZuMhuCCGwYMEC/Mu//Auampqwfv16FBUVsQ3YKIZLL3J0dMQf/vAHDBo0CCUlJXjjjTdQU1PDxkV2w8nJCRs3bkRERASuXr2KN998E3V1dWwDNojh0ouEEAgNDcX7778Pb29vpKWlYcOGDWhubjZ3aUS9QgiBsLAwbNmyBR4eHjh8+DA2b94MvV5v7tJIZgyXXiaEwNSpU7F+/Xo4Ojri008/xdatWznBT3bDuO/ef/7nf0KpVOKTTz5BYmIiJ/htDMPFDJRKJV555RWsXr0aAPDBBx8gMTGRAUN2Q6lUYuXKlVixYgUMBgO2bt2Kjz/+mAFjQxguZuLo6Ih169Zh+fLlaG9vx7vvvott27ahubmZ489kF5ydnbFhwwYsWbIELS0t2LhxIxITE6HX69kGbADDxYxcXFywefNmLFu2DG1tbdiyZQt+//vfc4KT7Iabmxvef/99/OIXv0BLSws2b96M9evXo76+nm3AyjFczEgIATc3N7z33ntYtWoVhBD461//imXLlqGkpISNi2yeEAIeHh7485//jP/4j/+AJEnYuXMnXn31Vdy4cYNtwIoxXMxMCAFXV1ds2rQJmzZtgru7Ow4dOoQXX3wR33//Pdrb281dIlGPEkLA3d0d77zzDjZs2ABXV1ccOHAACxYswNGjR7kXn5ViuFgAIQScnJywevVqfPrpp1Cr1cjLy8OSJUvw7rvvQqvV8hsc2TQhBJydnfHrX/8au3fvRv/+/XHp0iUsXrwY27ZtQ21tLduAlWG4WBClUomZM2di//79mDlzJhobG7FlyxYsXLgQmZmZaG9vZwMjm6ZUKjFv3jzs378fU6dORV1dHf74xz9i0aJFOHHiBHvyVoThYmGEENBoNNi7dy82btwIHx8fZGRkYP78+Vi7di3HocnmCSEQFRWF//mf/8Hvf/97eHl5IT09HS+++CLefvttlJWVsQ1YAYaLBTJOcv7mN7/B/v37MW3aNOh0OuzcuRPTp0/Hzp07UV1dzQZGNksIAS8vL6xduxb79u1DfHw86uvrkZiYiBkzZuCvf/0rqqqq2AYsGMPFgikUCjz77LP48ssv8dFHH0Gj0eDGjRtYu3YtZs+ejT179rCBkU1TKpWIiYnB119/jf/6r/9CWFgYioqK8MYbb2DWrFnYs2cP5yQtFMPFwhlX0vz7v/87UlNTsXbtWvj7++P8+fNYtWoVpk+fjg8//BAVFRVsYGSTjD35X/7ylzhy5Ah+97vfoW/fvrh48SJWr16NWbNm4dNPP+UXLQvDcLESQgj069cP69evx5EjR/Daa6+hb9++yMnJwVtvvYVp06bhz3/+M4qKijjxTzZJCIGQkBBs3LjR1AZ8fX1x7tw5rF69GtOmTcMHH3yAa9euwWAwsA2YGcPFyigUCmg0GmzZsgU//PAD1q1bh5CQEBQVFWHDhg2YPHkyfvnLX+KHH35AY2MjGxjZHIVCgYiICGzZsqXTF63c3Fy8/fbbmDJlCn71q18hMzMTOp2ObcBMGC5WSAgBhUKBsLAwvP322/jhhx+wadMmREVFoaamBl9++SXmz5+PWbNmYefOnSgsLERbWxsbGdkMYxuIjIzEli1b8OOPP2Ljxo3QaDSorKzE7t278fzzz2P27Nn4+OOPUVxczB59L2O4WDnjUMFvf/tbpKWl4fPPP8fcuXPh4uKCU6dO4c0338TkyZPx0ksv4auvvsLNmzfZyMhmGENm4MCBePPNN5GWlobdu3dj2rRpcHBwwIkTJ/C73/0OkydPxtKlS7Fv3z5UVFSwDfQClbkLIHkIIeDj44N58+Zh9uzZyM3Nxddff43vvvsOhYWFOHjwIFJSUhAYGIgJEybg+eefx8SJE9GvXz8oFAoIIcz9Foi6RQiBPn36YNGiRXjxxReRk5OD5ORkHDx4EFevXsU333yD/fv3o1+/foiNjcXzzz+PcePGITAwkG2gBzBcbIwQAg4ODhg2bBiGDh2K3/zmNzh9+jSSk5ORlpaG8vJyfPPNN9i3bx/69u2LmJgYTJ06FZMnT4ZarTZ3+UTdZtxOaeTIkRgxYgTWrFmDrKws7Nu3D+np6bh58yaSkpLw1VdfISAgADExMZg1axbi4uIQGhpq7vJtBsPFhgkh4O3tjeeeew5Tp05FZWUlTpw4gUOHDiE9PR1lZWU4cOAADhw4AF9fXyxYsAA+Pj7mLptINkII+Pr6Yvbs2Zg5cyYqKyuRlZWFlJQUHD9+HGVlZfjHP/6BAwcOwM/PD/PmzYO/v7+5y7YJDBc7IISAEAKBgYGYN28e5s6di6qqKpw6dQrff/89jh49ivz8fNTV1SEgIMDc5RLJTggBpVKJfv364cUXX8QLL7yAO3fu4NSpUzh8+DDS09NRXFyMxsZGBAcHm7tcm8BwkZEkSdBqtXBwcDB3KY+lVCoxbtw4xMTE4Fe/+hXS09OhVquRlpZm7tLIillTG1CpVBg/fjzGjRuH2tpanDhxAkFBQfjnP/9p7tJsAsNFJkII9O/fHx9//DGUSqW5y3kq586dQ1NTE7y8vMxdClkhW2gDANgGZCIkrseThSRJNrO00TiMRvQk2AaoI4YLERHJjjdREhGR7BguVkKSJG7GR3aP7cB6MFysxPnz5+Hm5obz58+buxQis2E7sB4MFyIikh3DhYiIZMdwISIi2TFciIhIdgwXIiKSHcOFiIhkx3AhIiLZMVyIiEh2DBciIpIdw4WIiGTHcCEiItkxXIiISHYMFyIikh3DhYiIZMdwISIi2TFcrIAkSaipqen0v0T2hu3AujBcLJhWq8X27duh0WgwdepU6PV6TJ06FRqNBtu3b4dWqzV3iUQ9ju3AOgmJ8W+RUlNTsWDBAuh0OgDo9C1NCAEAcHV1RXJyMhISEsxSI1FPYzuwXgwXC5Samoo5c+aYnhf+MAqFAkIIpKSksGGRzWE7sG4MFwuj1WoREhKCpqamRzYoI4VCARcXF5SWlsLb27vnCyTqBWwH1o9zLhZmz5490Ol0XWpQAGAwGKDT6bB3794eroyo97AdWD/2XCyIJEnQaDQoKip6opUwQgio1WoUFBSYxqGJrBXbgW1guFiQO3fuwN/fv1uv9/Pzk7Eiot7HdmAbOCxmQRoaGrr1+vr6epkqITIftgPbwHCxIO7u7t16vYeHh0yVEJkP24FtYLhYED8/P4SHhz/xeLEQAuHh4fD19e2hyoh6D9uBbWC4WBAhBF577bWneu3rr7/OSUyyCWwHtoET+haG6/uJ2A5sAXsuFsbb2xvJyckQQkChePT/PcY7k/ft28cGRTaF7cD6MVwsUEJCAlJSUuDi4gIhxH3dfOPPXFxccOjQIcyYMcNMlRL1HLYD68ZwsVAJCQkoLS1FYmIi1Gp1p2NqtRqJiYkoKytjgyKbxnZgvTjnYgUkSUJ1dTXq6+vh4eEBX19fTlqS3WE7sC4MFyIikh2HxYiISHYMFyIikh3DhYiIZMdwISIi2TFciIhIdgwXIiKSHcOFiIhkx3AhIiLZMVyIiEh2DBciIpIdw4WIiGTHcCEiItkxXIiISHYMFyIikt3/A4U2ubFgiRu7AAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 500x400 with 6 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"model.plot()"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "ccb7ec43",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"r2 is 0.9999974370002747\n",
"r2 is 0.9999890923500061\n",
"r2 is 0.999965488910675\n"
]
},
{
"data": {
"text/plain": [
"tensor(1.0000)"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"model.fix_symbolic(0,0,0,'log')\n",
"model.fix_symbolic(0,1,0,'log')\n",
"model.fix_symbolic(1,0,0,'sin')"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "0937db67",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"train loss: 2.85e-07 | test loss: 2.82e-07 | reg: 7.54e+00 : 100%|██| 20/20 [00:02<00:00, 9.03it/s]\n"
]
}
],
"source": [
"model.fit(dataset, opt=\"LBFGS\", steps=20);"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "e959cda3",
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"$\\displaystyle - 1.0 \\sin{\\left(2.0 \\log{\\left(2.205 x_{1} \\right)} + 2.0 \\log{\\left(2.018 x_{2} \\right)} + 0.156 \\right)}$"
],
"text/plain": [
"-1.0*sin(2.0*log(2.205*x_1) + 2.0*log(2.018*x_2) + 0.156)"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"ex_round(model.symbolic_formula()[0][0], 3)"
]
},
{
"cell_type": "markdown",
"id": "16e4da06",
"metadata": {},
"source": [
"We were lucky -- singularity does not seem to be a problem in this case. But let's instead consider $f(x,y)=\\sqrt{x^2+y^2}$. $x=y=0$ is a singularity point."
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "1ce52cec",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"train loss: 4.94e-03 | test loss: 5.23e-03 | reg: 5.98e+00 : 100%|██| 20/20 [00:03<00:00, 5.04it/s]\n"
]
}
],
"source": [
"from kan import *\n",
"import torch\n",
"\n",
"# create a KAN: 2D inputs, 1D output, and 5 hidden neurons. cubic spline (k=3), 5 grid intervals (grid=5).\n",
"model = KAN(width=[2,1,1], grid=5, k=3, seed=0)\n",
"f = lambda x: torch.sqrt(x[:,[0]]**2+x[:,[1]]**2)\n",
"dataset = create_dataset(f, n_var=2)\n",
"\n",
"# train the model\n",
"model.fit(dataset, opt=\"LBFGS\", steps=20);"
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "3a69ec41",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAFICAYAAACcDrP3AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAshUlEQVR4nO3deVRUV54H8O+tKpZiLUVEDdEGJLgSN0RlNSIYTNJGTexxMt1O7E7iMWZccrqzTVqTNplkNKLGyURjMhqT0J24i4rRCIgbrtGgIrhEAXEBiq2KYqk3f7RVB1xRXtUriu/nHE/O8Vm8H4Rffeve+959QpIkCURERDJSKV0AERE5H4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESy0yhdAFFbIEkSSktLUV1dDS8vL/j5+UEIoXRZRA6LIxeie9Dr9Vi8eDFCQ0Ph7++PoKAg+Pv7IzQ0FIsXL4Zer1e6RCKHJPgkSqI7S09Px4QJE2AwGAD8c/RiYRm1eHh4YO3atUhKSlKkRiJHxXAhuoP09HSMHTsWkiTBbDbf9d+pVCoIIZCWlsaAIWqC4UJ0C71ej8DAQBiNxnsGi4VKpYJWq0VhYSF0Op3tCyRqA7jmQnSLVatWwWAwtChYAMBsNsNgMGD16tU2royo7eDIhagJSZIQGhqK8+fP40FaQwiB4OBg5Ofn8yoyIjBciJq5ceMG/P39W/V6Pz8/GSsiaps4LUbURHV1dateX1VVJVMlRG0bw4WoCS8vr1a93tvbW6ZKiNo2hgtRE35+fggJCXngdRMhBEJCQtCxY0cbVUbUtjBciJoQQmDGjBkP9drXXnuNi/lEN3FBn+gWvM+FqPU4ciG6hU6nw9q1ayGEgEp17xax3KG/bt06BgtREwwXojtISkpCWloatFothBC3TXdZ/k6r1WLr1q1ITExUqFIix8RwIbqLpKQkFBYWIiUlBcHBwc2OBQcHIyUlBUVFRQwWojvgmgtRC0iShN27d2PUqFHYtWsXRo4cycV7onvgyIWoBYQQ1jUVnU7HYCG6D4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgvRfdTX16OoqAinT58GAJw7dw5lZWUwm80KV0bkuPiYY6K70Ov1WLt2Lb755hvk5uaiqqoKdXV1cHd3h7+/P2JiYjB16lRERUVBo9EoXS6RQ2G4EN3B/v37MWvWLJw4cQIREREYO3YswsPD4eXlBb1ejyNHjmDz5s0oKCjApEmT8Le//Q3+/v5Kl03kMBguRLfYsWMHpkyZAi8vL3z44YdITk5GXV0dUlNTYTKZ4OPjg9/97neor69Hamoq5s6di759++Lrr79GQECA0uUTOQSGC1ETZ8+exZgxY+Dp6YnU1FT06dMHQgicP38egwYNQkVFBYKCgnDkyBF06NABkiQhOzsbkydPRnx8PL744gu4ubkp/W0QKY4L+kQ3NTY24oMPPkB5eTk+/fRTa7DcixAC0dHR+Pjjj7Fx40Zs377dTtUSOTaGC9FNBQUF2Lx5M8aPH4/o6Oj7BouFEALjxo3DsGHDsGLFCjQ0NNi4UiLHx0tciG7at28fqqurMWHCBFy8eBE1NTXWY4WFhWhsbAQA1NXVITc3Fz4+Ptbj3bp1w/jx4zF37lyUlJQgMDDQ7vUTORKGC9FNZ86cgYeHB4KDg/Hyyy9j79691mOSJMFkMgEAiouLMXr0aOsxIQQWLlyI/v37w2AwoLi4mOFC7R7Dhegmo9EIjUYDNzc3mEwm1NbW3vHfSZJ027GGhgZotdpmIUTUnjFciG7q3LkzjEYj9Ho9IiMj4enpaT1mNBqxb98+a4iMGDHCeuOkEALdu3fHtWvXoFKp0KFDB6W+BSKHwXAhumnw4MGor69HTk4OPvroo2bHzp8/j4iICFRUVCAgIAB///vfodPprMeFEHjrrbfQpUsXTokRgVeLEVkNHToUwcHBWLVqFWpqaqBWq5v9sRBCQKVSWf9epVLhypUr+OGHHzB27Fj4+voq+F0QOQaGC9FNfn5+ePXVV3H06FEsWbKkxZcUm0wmvP/++zAajXj55ZdbfAkzkTPjtBhRE1OmTEFWVhY++ugjeHh4YNq0aXB3dwcAaDQaaDQa6yhGkiRUVVVh/vz5SE1NxaJFixAWFqZk+UQOg9u/EN3i+vXrmD59OrZs2YKkpCTMmjULvXv3Rl5eHsxmM1xdXdGzZ0/k5ORgwYIFOH78ON577z1Mmzat2fQZUXvGcCG6g5qaGqxYsQJLlizB1atXERwcjNDQUHh7e6O8vBx5eXkoLi7G4MGD8de//hVxcXFQqTjLTGTBcCG6h5KSEuzatQuZmZn4+eefkZOTg5iYGERFRSExMRGRkZHw8PBQukwih8NwIWqhQ4cOYejQoTh06BCGDBmidDlEDo3jeKIWUqvV1suQieje2CVERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmOz3MhaiFJkmA2m6FSqSCEULocIofGkQvRA+CzXIhaRqN0AURykSQJ+fn5KC0tVbqUVlGpVOjXrx88PT2VLoXooXFajJyG2WzG9OnT8eijj8LT0xONjY3Wp0e2JXv27MF//ud/Ijw8XOlSiB4aRy7kVNzc3DB06FB8+eWXKCkpwcCBA/Hss89i4MCB8PDwcPigkSQJ1dXV4Gc+aus4gUxOx2QyYdOmTcjMzERKSgoSExORlJSE9PR0NDQ08I2byA4YLuR0QkJC8MYbb+CDDz5AQkIC3NzcsH//fkyaNAlvvfUWSktLGTBENsZwIafToUMHvP322/jzn/+MzZs3Y9euXXjmmWdQW1uLTz75BE899RR27tzJUQyRDTFcyCkJISCEgKurKwYOHIg1a9ZgwYIF6Ny5Mw4dOoQJEyZgxowZuHDhAgOGyAYYLuT0hBDw9PTE9OnTsW3bNjzzzDOoq6vD8uXLMWrUKCxZsgR6vZ4hQyQjhgu1GyqVCuHh4fj222+xcuVK9OrVC5cvX8brr7+O5ORkbN26FXV1dQwZIhkwXKhdEUJAq9Vi8uTJ2LVrF95880107NgRBw8exPPPP4+pU6fi7NmzDBiiVmK4ULskhECXLl0wb9487NixA8899xwkScK3336LxMRE/O///i9qamoYMkQPieFC7ZplqmzVqlX45ptvEB4ejqKiIsycOROTJ09Gbm4uA4boITBcqN0TQsDNzQ3jxo3D9u3bMWvWLLi7u2PLli1ITk7Gl19+idraWoYM0QNguBDdJIRA586d8eGHH+L777/HgAEDUFxcjFdffRUvv/wyLl++zIAhaiGGC1ETQghoNBqMHj0aaWlpeOWVV6BSqbBmzRqMHTsW6enpaGxsVLpMIofHcCG6AyEEAgICsGjRInz55Zf4zW9+g1OnTuFf/uVf8MEHH6CqqoqjGKJ7YLgQ3YUQAi4uLnj++eeRlpaG5ORk1NTU4P3338cLL7yAgoICBgzRXTBciO5DCIGwsDB8++23ePfdd+Hp6YktW7bgqaeewtatWzlNRnQHDBeiFhBCwNvbG2+++Sa+/fZbPPbYYygoKMALL7yAlJQUGI1GjmKImmC4ED0AtVqNMWPGYPPmzUhOTkZ1dTXefvttzJw5E2VlZQwYopsYLkQPSAiBkJAQrFmzBjNnzoRKpcLKlSvxb//2b7h06RIDhggMF6KHIoSAr68v5s+fj08++QS+vr5IT0/HpEmTcPr0aQYMtXsMF6JWcHV1xUsvvYSvvvoK3bp1w6FDhzBp0iQcO3aMAUPtGsOFqJVUKhWefvppfPfddwgODsapU6cwefJkHD58mAFD7RbDhUgGQghERUXhu+++Q1hYGPLz8/H73/8ex48fZ8BQu8RwIZKJEAKDBw/GmjVr8Nhjj+Hs2bP4wx/+gDNnzjBgqN1huBDJSAiBgQMHYvXq1QgKCkJubi5efPFFbnpJ7Q7DhUhmQggMGTIEK1euRJcuXZCTk4Np06bxPhhqVxguRDYghEBsbCyWLVsGnU6H9PR0vPHGG6itrVW6NCK7YLgQ2YgQAk8//TTmz58PNzc3rF69GikpKWhoaFC6NCKbY7gQ2ZBarcaLL76IV199FWazGf/1X/+FjRs3cnqMnB7DhcjGXF1d8c477+C3v/0tqqurMWfOHJw4cYIBQ06N4UJkB15eXli4cCHCw8Nx+fJlzJ49G3q9XumyiGyG4UJkB0IIdO/eHSkpKejYsSOysrKwZMkSmM1mpUsjsgmGC5GdCCEQExODOXPmAACWLVvGPcjIaTFciOxIpVLhlVdeQWRkJEpLS/Hxxx+jvr5e6bKIZMdwIbIzX19f/OUvf4Gbmxu2bduGAwcOcPRCTofhQmRnQggkJCQgNjYWNTU1WL58ORobG5Uui0hWDBciBbi7u+NPf/oTNBoN0tPTkZ+fr3RJRLJiuBApQAiBUaNGoXfv3igrK8OGDRs4NUZOheFCpBBfX188++yzAIANGzbAYDAoXBGRfBguRAqx7D3m5eWF3Nxc5ObmKl0SkWwYLkQK6tWrF/r16wej0YgdO3YoXQ6RbBguRArSarVISEgAAOzatQt1dXUKV0QkD4YLkYKEEHjiiSfg6uqKX375BZcvX1a6JCJZMFyIFNa3b18EBgZCr9fjyJEjSpdDJAuGC5HCOnTogIEDB8JsNmPPnj1Kl0MkC4YLkcJUKhWio6OhVquRl5cHk8mkdElEraZRugAiOUmShPLycri4uChdygOJjY3FypUr0b9/f/zwww9Kl0PUagwXchpCCPTo0QNLly6FWq1WupyHcujQIRiNRvj6+ipdClGrCIl7TpCTkCTJabZQEUJACKF0GUQPjeFCRESy44I+ERHJjmsuRC3UdJDPKSuie+PIhaiFjh07BrVajWPHjildCpHDY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4ULUApIkoby8HABQXl4OPh2c6N4YLkT3oNfrsXjxYoSGhiIhIQGSJCEhIQGhoaFYvHgx9Hq90iUSOSQh8SMY0R2lp6djwoQJMBgMAO78mGMPDw+sXbsWSUlJitRI5KgYLkR3kJ6ejrFjx0KSJJjN5rv+O5VKBSEE0tLSGDBETTBciG6h1+sRGBgIo9F4z2CxUKlU0Gq1KCwshE6ns32BRG0A11yIbrFq1SoYDIYWBQsAmM1mGAwGrF692saVEbUdHLkQNSFJEkJDQ3H+/PkHuiJMCIHg4GDk5+db12OI2jOGC1ETN27cgL+/f6te7+fnJ2NFRG0Tp8WImqiurm7V66uqqmSqhKhtY7gQNeHl5dWq13t7e8tUCVHbxnAhasLPzw8hISEPvG4ihEBISAg6duxoo8qI2haGC1ETQgjMmDHjoV772muvcTGf6CYu6BPdgve5ELUeRy5Et9DpdFi7di2EEFCp7t0iljv0161bx2AhaoLhQnQHSUlJSEtLg1arhRDitukuy99ptVps3boViYmJClVK5JgYLkR3kZSUhMLCQqSkpCA4OLjZseDgYKSkpKCoqIjBQnQHXHMhagFJkrB7926MGjUKu3btwsiRI7l4T3QPHLkQtYAQwrqmotPpGCxE98FwISIi2TFciIhIdgwXIiKSHcOFiIhkx3AhIiLZMVyIiEh2DBciIpIdw4WIiGTHcCEiItkxXIiISHYMFyIikh3DhYiIZMdwISIi2TFciIhIdgwXIiKSHcOFiIhkx3Ahuo/6+noUFRXh9OnTAIBz586hrKwMZrNZ4cqIHBcfc0x0F3q9HmvXrsU333yD3NxcVFVVoa6uDu7u7vD390dMTAymTp2KqKgoaDQapcslcigMF6I72L9/P2bNmoUTJ04gIiICY8eORXh4OLy8vKDX63HkyBFs3rwZBQUFmDRpEv72t7/B399f6bKJHAbDhegWO3bswJQpU+Dl5YUPP/wQycnJqKurQ2pqKkwmE3x8fPC73/0O9fX1SE1Nxdy5c9G3b198/fXXCAgIULp8IofAcCFq4uzZsxgzZgw8PT2RmpqKPn36QAiB8+fPY9CgQaioqEBQUBCOHDmCDh06QJIkZGdnY/LkyYiPj8cXX3wBNzc3pb8NIsVxQZ/opsbGRnzwwQcoLy/Hp59+ag2WexFCIDo6Gh9//DE2btyI7du326laIsfGcCG6qaCgAJs3b8b48eMRHR1932CxEEJg3LhxGDZsGFasWIGGhgYbV0rk+HiJC9FN+/btQ3V1NSZMmICLFy+ipqbGeqywsBCNjY0AgLq6OuTm5sLHx8d6vFu3bhg/fjzmzp2LkpISBAYG2r1+IkfCcCG66cyZM/Dw8EBwcDBefvll7N2713pMkiSYTCYAQHFxMUaPHm09JoTAwoUL0b9/fxgMBhQXFzNcqN1juBDdZDQaodFo4ObmBpPJhNra2jv+O0mSbjvW0NAArVbbLISI2jOGC7V7Fy9eRGZmJrKzs2EwGKDX6xEZGQlPT0/rvzEajdi3b581REaMGGG9cVIIge7du+PatWtoaGhAQUEBIiIi4O7urtS3RKQ4XopM7c7ly5eRlZWFjIwMZGZm4tKlSxBCoEePHigoKMCyZcvwxz/+sdlrzp8/j4iICFRUVOA3v/kNDh8+DJ1OZz0uhMBbb72FBQsWQK1Ww93dHZGRkYiNjUV8fDwiIiJ4iTK1KwwXcnrFxcXIzMy0/rlw4QIAIDw8HLGxsYiLi0N0dDTMZjOio6PRoUMHbN++vdmC/d3ucwH+OU1WXFyMuLg4PP3005gyZQqysrKQlZWFPXv2QK/XQ6vVYtiwYdawGTx4MFxdXRX5eRDZA8OFnM7Vq1eRkZGBrKwsZGZmoqCgAADQt29fa5jExsaiY8eOt7122bJlmDNnDt555x288cYb1qmve4VLbW0tZs6cic2bN+Onn35CWFiY9es1NjbixIkT1lqys7NRWVkJDw8PDB8+HHFxcYiPj8fAgQPh4uJih58OkX0wXKjNu379OjIzM61v4Hl5eQCAXr16NQuTluz9VVNTgxdffBFbt27FvHnzMG3aNLi7u+PChQsYOnSodVosJycHOp0OVVVVmD9/Pj7//HMsWrQI//7v/37Pr9/Q0IDjx49b6927dy+qq6vh5eWFESNGIC4uDnFxcRgwYAA3w6Q2jeFCbU5paal12ikzMxOnTp0CAISGhlrDJC4u7qH3+bp+/TqmT5+OLVu2ICkpCbNmzULv3r2Rl5cHs9kMV1dX9OzZEzk5OViwYAGOHz+O9957D9OmTYNarX6gc9XX1+PYsWPWkda+fftgMBjg4+ODqKgo6/cSHh7+wF+bSEkMF3J45eXl2LNnj/XT/smTJwEAwcHBzcKkW7dusp2zpqYGK1aswJIlS3D16lUEBwcjNDQU3t7eKC8vR15eHoqLizF48GD89a9/RVxcHFSq1m94UVdXhyNHjljXh/bv34/a2lrodDpER0dbv9d+/frJcj4iW2G4kMOpqKhAdna2NUx+/vlnSJKEHj16WN9c4+Li7HKjYklJCXbt2oXMzEycP38etbW16NChA/r164fExERERkbCw8PDZuc3mUw4dOiQNWwOHjwIk8mEjh07IiYmxvqzaMk+aET2xHAhxVVVVWHv3r3Waa5jx47BbDbjkUcesS54x8XFoUePHorW2djYCEmSoFKpFBs11NbW4uDBg9ZLqXNyclBfX49OnTo1G8WFhYUxbEhRDBeyu+rqauzfv9/6afzo0aNobGxE165dER8fb32TDAoK4hvkfRgMBhw4cMAaNocPH0ZDQwM6d+7cbJTXs2dP/izJrhguZHOWN0BLmDR9A2waJnwDbL3q6mrrzzojI6NZcDcNGwY32RrDhWRnmbqxhEnTqRvLZcGcurGPysrKZqNEy5RjYGCg9YZOR5hyJOfDcKFWsyw6Wy6nbbroHBsbaw0TLjorT6/XN1vfanqxRNOw4a7O1FoMF3pglstlLXtzHThwwHq5bExMjDVMeLms4ysrK8PevXutIxvLZd5BQUHWKbT4+Hh07dpV4UqprWG40H01vdHPcu+F5UY/y70XsbGxvNHPCZSWllrvKWp6g2rPnj2brdk87A2q1H4wXOg2TbcoyczMtD6h0cvLC1FRUdaRCbcocX7Xr1+3TqE13VonLCysWdh06tRJ4UrJ0TBcyLq5ouUNZO/evdbNFZvud8XNFamkpKRZ2Fg2Be3Tp4/19+Rum4JS+8JwaYfMZjN++eUX6xtEdnY29Ho93N3drTv1xsXFcVt4uq+7Pc6gf//+1t8jy2MMqH1huLQDkiTh1KlT1jeAPXv2oKysDG5uboiMjLS+CfCBVtRad3sQ2+OPP279PYuKioKvr6/SpZKNMVyckCRJyMvLs4ZJVlYWbty4ARcXFwwdOtR642JkZCQfxUs2dfHixWZhU1RUBJVKhYEDB1rDZsSIEfD29la6VJIZw8UJSJKEgoKCZtMT165dg0ajwZAhQ6xhMmzYMJtuskh0L5Ik4cKFC9bf0YyMDJSUlECtVmPQoEHWe2yGDx8OT09PpculVmK4tEG3NmlmZiauXLlyW5MOGzYMXl5eSpdLdEct+VAUFxdn852nyTYYLm3Er7/+av20l5WVhcLCwtumF4YPH97sue9Ebcmt07mZmZkoLS2Fq6srIiIiOJ3bxjBc2ojw8HDk5+c3WxgdMWIEdDqd0qUR2YTZbMbp06dvuxDl66+/xnPPPad0eXQfDJc2wmw2QwjBvbmo3ZIkCZIksQ/aCIYLERHJjnt3yMSyOFlaWqp0Ka2iUqnQt29fXq1DD4w9QE0xXGQiSRKWLFmCwMBAm/1SWh6zq1arbTYtsHfvXrz11lvo37+/Tb4+OS979oAtHzXNHpAHw0VGbm5umDJlimw7xkqSBJPJhJycHGzZsgW5ubkwGo3o0qULYmJiMHbsWDz66KOyBY0kSaipqQFnSulh2aIH6urqcPjwYWzZsgUnT56EwWBAQEAAoqKi8NRTT6F79+6yBQ17QD4MFwclSRJyc3Mxb9487Ny5E7W1tc2Ob9iwAQsXLsSMGTMwdepUeHp6cpGTnIokSTh9+jTee+897NixA0ajsdnxDRs2YNGiRZg+fTr+9Kc/wcvLiz3gQBguDshsNmPLli2YPXs2ioqK4OHhgaSkJMTFxcHHxwd5eXnYvn07CgoK8Pbbb+PgwYNYuHAhunTpwuYip2A2m7F161bMnj0bly9fhlarRUJCAuLj46HT6ZCfn4/t27fj7NmzePfdd3Ho0CF88skn7AEHwnBxMGazGevXr8eMGTNQXl6Oxx9/HPPnz0dMTEyz7e7nzJmDzz77DMuWLcP69euh1+uxcuVKNhe1eWazGevWrcNrr72G8vJy9O/fH++//z7i4+Ob7dJt6YGlS5diw4YNqKysxMqVKxEQEMAecAB8Bq0DkSQJmZmZmDlzJsrLyzFq1Ch8//33eOKJJ+Dq6mq9vl8IgYCAALzzzjtYtmwZOnbsiN27d2PmzJmorKzkfDG1WZIk4aeffmrWAz/88AMSExPh5ubWrAf8/f3x9ttvW3vgp59+wuzZs1FVVaX0t0FguDgMy35hM2fOxI0bNxAVFYXly5cjMDDwrp/CNBoNJk6ciJSUFHh7e2PLli346KOP0NDQYOfqiVpPkiTk5+dj5syZKC0tRUxMDJYvX37Pi1bUajUmTpyIhQsXwsvLCxs3bsTChQvZAw6A4eIgjEYj3nnnHZw9exbBwcH49NNP0bVr1/sO71UqFcaPH4833ngDKpUKn3/+OdLS0jh6oTanpqYGb775Js6dO4fQ0FAsW7asxT0wceJEvP766xBC4H/+53/w448/sgcUxnBxAJIk4bvvvsPmzZvh4eGBDz/8EGFhYS2eN1ar1XjllVfwzDPPwGAwYO7cuSgqKmJzUZthNpvx1VdfIT09Hd7e3vj444/Rs2fPFveARqPBq6++iuTkZFRXV+Pdd99FSUmJjaume2G4KEySJFy8eBH//d//jYaGBvz+97/Hk08++cALklqtFvPmzUP37t2Rl5eHRYsWobGx0UZVE8lHkiScPXsWKSkpMJvN+OMf/4iEhIQH7gEPDw/MmzcP3bp1Q25uLpYuXcoeUBDDRWGNjY1YvHgxfv31V4SGhmLOnDnQaB78Ij4hBEJCQvDnP/8ZGo0Ga9asweHDhzl6IYfX0NCABQsWoLi4GP369cN//Md/PHQP9OrVC7Nnz4ZKpcL//d//4eeff2YPKIThoiBJknDs2DGkpqZCo9Fgzpw5eOSRRx76MkohBCZNmoSoqChUVlZi4cKFMJlMMldNJB9JkrB3715s2LABrq6u+Mtf/oLOnTs/9NcTQuCFF15AREQEysvL8cknn6C+vl7GiqmlGC4Kqq+vx5IlS1BRUYHIyEiMHz++1dfne3p64vXXX4dWq8XOnTuRmZnJT27ksGpra5GSkoKamhqMHDkSycnJre4BHx8fvP7663Bzc8O2bduwb98+9oACGC4KkSQJhw8fxrZt2+Dq6oqZM2fK8khiIQRiYmKQlJSE2tpaLF269LatY4gcgSRJyMjIQEZGBjw8PDB79mxZnjAphMCoUaMwcuRIGAwGLF26FHV1dTJUTA+C4aKQhoYGfPbZZ6ipqUFUVBSeeOIJ2e4qdnV1xfTp0+Hh4YHs7GxkZ2fzkxs5HJPJhM8++wwmkwmJiYkYNmyYbD3g5uaGGTNmQKvVYvfu3di/fz97wM4YLgqQJAknT57Ejh074OLigldeeQVarVa2ry+EwNChQ/HEE0/AZDJhxYoVnHcmhyJJEg4cOIDs7GxotVpMmzat2fZGrSWEwIgRIxAbGwuj0Yjly5fzxko7Y7goQJIkrFq1CpWVlRg0aJCsoxYLFxcXvPTSS3Bzc8Pu3btx/PhxfnIjh9HQ0IAvvvgCRqMRsbGxiIyMlL0H3Nzc8NJLL8HV1RU7d+7EiRMn2AN2xHBRwKVLl7Bp0yaoVCq8+OKLNnmwkhACUVFRGDJkCKqrq7FmzRo2FjkEy1b6O3fuhIuLC6ZOndpsQ0q5CCEQGxuLQYMGoaqqij1gZwwXO5MkCevXr8fVq1cREhLyUDdMtpRWq8Uf/vAHqFQqbNmyBZcvX7bJeYgehGVHioqKCoSHhyMuLs5mPeDp6WntgU2bNqGwsNAm56HbMVzsrLKyEqmpqZAkCRMnTkSnTp1sdi4hBJKSkhAUFISSkhJs2rSJn9xIcdeuXcOGDRsghMC//uu/wtvb22bnEkLgySefRFBQEK5cuYINGzawB+yE4WJHlhvGTp8+jQ4dOmDixIk2f+6Ev78/xo0bB0mS8P3336Ompsam5yO6F0mS8OOPP+LSpUvo2rUrnn76aZv3QOfOnfHss89CkiT84x//QHV1tU3PR//EcLGjxsZG/P3vf0d9fT1iY2MRGhpq83MKITBx4kR4e3vjxIkTOHLkCD+5kWLq6uqQmpoKs9mM5ORkdOvWzebnFELgueeeg4+PD3755Rfk5OSwB+yA4WJHly5dQkZGBtRqNSZNmvRQ+yc9jN69e2Po0KEwmUz44Ycf2FikmFOnTuHQoUNwd3fH888/b7cnRvbq1QsjRoyAyWTCP/7xD5jNZructz1juNiJJEnYvn07rl+/jqCgIMTExNitsVxdXfHcc89BpVJhx44duHr1ql3OS9SUJEnYuHEjqqqq0K9fPwwaNMhuPeDi4oLnn38eKpUKO3fuxJUrV+xy3vaM4WInJpMJ69evhyRJePLJJ+Hn52e3c1u2w+jSpQsKCwuxZ88ejl7I7iorK5GWlgYAGDduHDw8POx2biEE4uPjERgYiJKSEuzatYs9YGMMFzs5c+YMjh8/Dnd3d4wbN87u5+/atSvi4uJgNpuxfv16PueC7EqSJBw5cgRnz56FTqeTZYPKBxUQEICEhARrD3DXCttiuNiBJEnYtm0bqqur0adPH4SHh9u9sVQqFcaNGweNRoN9+/ahqKjIruen9k2SJGzatAl1dXUYMmQIQkJC7F6DEALjxo2Di4sLcnJycPHiRbvX0J4wXOzAYDBg69atAIDk5GSb3JF/P0IIDB8+HI8++iiuX7+OjIwMTguQ3ZSXl2PXrl0QQuC3v/2trPuItZQQwhpser0eO3fuZA/YEMPFDs6cOYNTp07B09MTY8aMsfuoxcLPzw8jR46EJElIS0vj1BjZhWVK7Ndff7X+DirVAzqdDqNHj7b2ALfitx2Gi41JkoQdO3bAYDCgT58+6N27t2K1CCEwduxYaDQaHDx4kFNjZBeWaeH6+noMGTIE3bt3V7Se5ORkuLq64tixY5wasyGGi43V1tYiPT0dAJCYmCjr1voPyjItEBgYiBs3bvA5L2QXFRUV2L17N4QQSE5Ottv9XXcihMCAAQMQHBwMvV7P6WEbYrjYWEFBAU6dOgWtVovRo0crNh1g0alTJ0RHR8NsNmPbtm28mYxs7uTJk7hw4QJ0Op1NN6lsKV9fX+v08Pbt2/mcFxthuNiQ5TGuVVVVCA0NRZ8+fZQuCUIIjBkzBiqVCgcOHMC1a9eULomcmGUvsbq6OvTv3x89evRQuiQA/5xF0Gg0OHr0KIqLi5UuxykxXGyooaEBP/74IwAgLi4OXl5eClf0z3CJjIxE586dUVJSwr3GyKaMRiMyMjIAAKNHj7bJc1selBACgwYNQrdu3XDjxg0cPHiQPWADDBcbKi4uxvHjx6HRaDB69Gily7Hq0qULBg8ejMbGRmv4EdlCQUEBzpw5A61Wi/j4eMWnxCz8/PwQGRkJs9mMHTt2MFxsgOFiI5IkIScnB6WlpejWrRsef/xxh2kstVptDbvs7GxUVlYqXBE5I0mSkJ2djerqaoSEhCAsLEzpkqxUKpV1DfTAgQPQ6/VKl+R0GC42IkkSdu7cCbPZjIiICLvuJXY/QghER0fDx8cHFy5cwJkzZ5QuiZxQY2Mjdu/eDQCIjo52iGlhCyEEhg0bBp1Oh8uXLyM3N1fpkpwOw8VGKisrceDAAQghkJCQAJXKsX7UQUFBCAsLg9Fo5EaWZBPXr1/H0aNHoVarMXLkSKXLuc2jjz6K3r17o66uDllZWewBmTnWO54TOXPmDC5dugQfHx8MGzbMYabELLRaLaKjowEAGRkZ3MSPZCVJEk6cOIFr167B39/frtvrt5SbmxtiY2MBAFlZWbxbX2YMFxuwPM64trYWYWFhDnP55a1GjhwJjUaDkydP8hkvJLvMzEw0NDSgf//+CAgIULqc2wghEBsbCxcXF5w6dQolJSVKl+RUGC420NDQgKysLADAiBEj4O7urnBFtxNCWJu+tLQUR48eVbokciK1tbXYu3cvACA2NlbRu/LvpU+fPujSpQvKysrYAzJjuNjAtWvXcPLkSWg0GsTFxSldzl116tQJgwYNQmNjozUMieRQVFSE/Px8uLm5ISoqyuGmxCz8/PwwYMAAmM1mZGdnK12OU2G42MAvv/yC69evw9/fH/369XPYxlKr1YiJiQEA7N+/HwaDQeGKyFkcPXoUFRUVeOSRR/DYY48pXc5dqdVq69rjwYMHYTQaFa7IeTBcZGa5tr+hoQF9+/Z1yLlmC8szXrRaLc6dO4dff/1V6ZLICUiShD179sBsNmPAgAHQ6XRKl3RXlh0r3N3dce7cORQWFipdktNguMisvr4e+/fvBwBERUU57FyzRc+ePdG9e3dUVVXh8OHDSpdDTsBoNFp/l6Kjox3uMvxbhYaGolu3bqioqMDx48eVLsdpOPb/9TboypUryMvLg4uLC4YPH+6wU2IWPj4+GDx4sPXTJq/1p9a6dOkSzp8/D61Wi4iICIfvAV9fXzz++OOQJAn79u1TuhynwXCR2blz51BbW4uAgAD06tVL6XLuy3I5po+PD6qrq3mtP7XajRs30LlzZ4SEhCAkJETpcu5LpVJhxIgRUKvVyMvLg8lkUrokp+DYczZtjCRJ6N27N1JTU3H16lWoVCqUlZUpXdZ9xcTEYOPGjQgMDMQXX3yhdDnUhkmShL59+2LTpk0oKytDY2Njm+iB6OhofP755+jXrx/WrVundDlOgeEiEyEEunfvjq+++gpqtRoA2uQahtFohK+vr9JlUBtk6YHPPvvM2gNtzZEjR9gDMhESJ9llIUmS06xXCCEcfp6cHA97gJpiuBARkey4oE9ERLJjuLQRkiTBbDY7zbQD0cNgH7QdDJc24vjx4/Dw8OBNXtSusQ/aDoYLERHJjuFCRESyY7gQEZHsGC5ERCQ7hgsREcmO4UJERLJjuBARkewYLkREJDuGCxERyY7hQkREsmO4EBGR7BguREQkO4YLERHJjuFCRESyY7gQEZHsGC5tgCRJKC8vb/ZfovaGfdC2MFwcmF6vx+LFixEaGopRo0ahrq4Oo0aNQmhoKBYvXgy9Xq90iUQ2xz5om4TE+HdI6enpmDBhAgwGAwA0+5QmhAAAeHh4YO3atUhKSlKkRiJbYx+0XQwXB5Seno6xY8danxd+NyqVCkIIpKWlsbHI6bAP2jaGi4PR6/UIDAyE0Wi8Z0NZqFQqaLVaFBYWQqfT2b5AIjtgH7R9XHNxMKtWrYLBYGhRQwGA2WyGwWDA6tWrbVwZkf2wD9o+jlwciCRJCA0Nxfnz5x/oShghBIKDg5Gfn2+dhyZqq9gHzoHh4kBu3LgBf3//Vr3ez89PxoqI7I994Bw4LeZAqqurW/X6qqoqmSohUg77wDkwXByIl5dXq17v7e0tUyVEymEfOAeGiwPx8/NDSEjIA88XCyEQEhKCjh072qgyIvthHzgHhosDEUJgxowZD/Xa1157jYuY5BTYB86BC/oOhtf3E7EPnAFHLg5Gp9Nh7dq1EEJApbr3/x7Lncnr1q1jQ5FTYR+0fQwXB5SUlIS0tDRotVoIIW4b5lv+TqvVYuvWrUhMTFSoUiLbYR+0bQwXB5WUlITCwkKkpKQgODi42bHg4GCkpKSgqKiIDUVOjX3QdnHNpQ2QJAllZWWoqqqCt7c3OnbsyEVLanfYB20Lw4WIiGTHaTEiIpIdw4WIiGTHcCEiItkxXIiISHYMFyIikh3DhYiIZMdwISIi2TFciIhIdgwXIiKSHcOFiIhkx3AhIiLZMVyIiEh2DBciIpIdw4WIiGT3/7T7NL8LxJUNAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 500x400 with 6 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"model.plot()"
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "abef7aa9",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"r2 is 0.9999871253967285\n",
"r2 is 0.9999728798866272\n",
"r2 is 0.9998090863227844\n"
]
},
{
"data": {
"text/plain": [
"tensor(0.9998)"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"model.fix_symbolic(0,0,0,'x^2')\n",
"model.fix_symbolic(0,1,0,'x^2')\n",
"model.fix_symbolic(1,0,0,'sqrt')"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "e14000d8",
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"$\\displaystyle 1.01262303277627 \\sqrt{\\left(8.59418125232242 \\cdot 10^{-5} - x_{2}\\right)^{2} + 0.999965395886852 \\left(- x_{1} - 2.26198704007758 \\cdot 10^{-5}\\right)^{2} + 0.00768977463773129} - 0.0159889459609985$"
],
"text/plain": [
"1.01262303277627*sqrt((8.59418125232242e-5 - x_2)**2 + 0.999965395886852*(-x_1 - 2.26198704007758e-5)**2 + 0.00768977463773129) - 0.0159889459609985"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"formula = model.symbolic_formula()[0][0]\n",
"formula"
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "c56ee3d5",
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"$\\displaystyle 1.01 \\sqrt{1.0 x_{1}^{2} + x_{2}^{2} + 0.01} - 0.02$"
],
"text/plain": [
"1.01*sqrt(1.0*x_1**2 + x_2**2 + 0.01) - 0.02"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"ex_round(formula, 2)"
]
},
{
"cell_type": "markdown",
"id": "1fd57d41",
"metadata": {},
"source": [
"w/ singularity avoiding"
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "de708f21",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"train loss: 4.85e-08 | test loss: 4.84e-08 | reg: 5.95e+00 : 100%|██| 20/20 [00:01<00:00, 14.88it/s]\n"
]
}
],
"source": [
"model.fit(dataset, opt=\"LBFGS\", steps=20, update_grid=False, singularity_avoiding=True);"
]
},
{
"cell_type": "markdown",
"id": "6fd34c4c",
"metadata": {},
"source": [
"w/o singularity avoiding, nan may appear"
]
},
{
"cell_type": "code",
"execution_count": 18,
"id": "031fabd6",
"metadata": {
"scrolled": true
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"train loss: nan | test loss: nan | reg: nan : 25%|████▌ | 5/20 [00:01<00:03, 3.90it/s]"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"Intel MKL ERROR: Parameter 6 was incorrect on entry to SGELSY.\n",
"\n",
"Intel MKL ERROR: Parameter 6 was incorrect on entry to SGELSY.\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\n"
]
},
{
"ename": "RuntimeError",
"evalue": "false INTERNAL ASSERT FAILED at \"/Users/runner/work/pytorch/pytorch/pytorch/aten/src/ATen/native/BatchLinearAlgebra.cpp\":1540, please report a bug to PyTorch. torch.linalg.lstsq: (Batch element 0): Argument 6 has illegal value. Most certainly there is a bug in the implementation calling the backend library.",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mRuntimeError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m/var/folders/6j/b6y80djd4nb5hl73rv3sv8y80000gn/T/ipykernel_33275/1949812002.py\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mmodel\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfit\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdataset\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mopt\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"LBFGS\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msteps\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m20\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m;\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
"\u001b[0;32m~/Desktop/2022/research/code/pykan/kan/MultKAN.py\u001b[0m in \u001b[0;36mfit\u001b[0;34m(self, dataset, opt, steps, log, lamb, lamb_l1, lamb_entropy, lamb_coef, lamb_coefdiff, update_grid, grid_update_num, loss_fn, lr, start_grid_update_step, stop_grid_update_step, batch, small_mag_threshold, small_reg_factor, metrics, save_fig, in_vars, out_vars, beta, save_fig_freq, img_folder, device, singularity_avoiding, y_th, reg_metric)\u001b[0m\n\u001b[1;32m 804\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 805\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0m_\u001b[0m \u001b[0;34m%\u001b[0m \u001b[0mgrid_update_freq\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m0\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0m_\u001b[0m \u001b[0;34m<\u001b[0m \u001b[0mstop_grid_update_step\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0mupdate_grid\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0m_\u001b[0m \u001b[0;34m>=\u001b[0m \u001b[0mstart_grid_update_step\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 806\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mupdate_grid_from_samples\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdataset\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'train_input'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mtrain_id\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mto\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdevice\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 807\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 808\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mopt\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;34m\"LBFGS\"\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/Desktop/2022/research/code/pykan/kan/MultKAN.py\u001b[0m in \u001b[0;36mupdate_grid_from_samples\u001b[0;34m(self, x)\u001b[0m\n\u001b[1;32m 286\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0ml\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdepth\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 287\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mforward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 288\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mact_fun\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0ml\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mupdate_grid_from_samples\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0macts\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0ml\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 289\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 290\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0minitialize_grid_from_another_model\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmodel\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/Desktop/2022/research/code/pykan/kan/KANLayer.py\u001b[0m in \u001b[0;36mupdate_grid_from_samples\u001b[0;34m(self, x)\u001b[0m\n\u001b[1;32m 228\u001b[0m \u001b[0mgrid\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrid_eps\u001b[0m \u001b[0;34m*\u001b[0m \u001b[0mgrid_uniform\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m \u001b[0;34m-\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrid_eps\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m*\u001b[0m \u001b[0mgrid_adaptive\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 229\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrid\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdata\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mextend_grid\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mgrid\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mk_extend\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mk\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 230\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcoef\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdata\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mcurve2coef\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx_pos\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0my_eval\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrid\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mk\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdevice\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdevice\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 231\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 232\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0minitialize_grid_from_parent\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mparent\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/Desktop/2022/research/code/pykan/kan/spline.py\u001b[0m in \u001b[0;36mcurve2coef\u001b[0;34m(x_eval, y_eval, grid, k, device)\u001b[0m\n\u001b[1;32m 164\u001b[0m \u001b[0;31m# coef shape: (in_dim, outdim, G+k)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 165\u001b[0m \u001b[0my_eval\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0my_eval\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpermute\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0munsqueeze\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdim\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;31m# y_eval: (in_dim, out_dim, batch, 1)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 166\u001b[0;31m coef = torch.linalg.lstsq(mat.to(device), y_eval.to(device),\n\u001b[0m\u001b[1;32m 167\u001b[0m driver='gelsy' if device == 'cpu' else 'gels').solution[:,:,:,0]\n\u001b[1;32m 168\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mcoef\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mto\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdevice\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;31mRuntimeError\u001b[0m: false INTERNAL ASSERT FAILED at \"/Users/runner/work/pytorch/pytorch/pytorch/aten/src/ATen/native/BatchLinearAlgebra.cpp\":1540, please report a bug to PyTorch. torch.linalg.lstsq: (Batch element 0): Argument 6 has illegal value. Most certainly there is a bug in the implementation calling the backend library."
]
}
],
"source": [
"model.fit(dataset, opt=\"LBFGS\", steps=20);"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.7"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Reference in New Issue View Git Blame Copy Permalink