preliminary DP NS code update
This commit is contained in:
NT committed
1 parent
05f99a4ed2
commit
d42fe17a6f
1 file changed
+89
-154
+89
-154
@@ -10,7 +10,7 @@
|
||||
"\n",
|
||||
"**(Note - this example is not final. TODO: update from diffphys-code-ns-v1-syncWith2.ipynb, diff display for target, adjust params so that the loss properly goes down)**\n",
|
||||
"\n",
|
||||
"Next, we'll target a more complex example with the Navier-Stokes equations as model. In line with {doc}`overview-ns-forw`, we'll target a 2D case with velocity $\\mathbf{u}$, no explicit viscosity term, and a marker density $d$ that drives a simple Boussinesq buoyancy term $\\eta d$ adding a force along the y dimension. For the velocity this gives:\n",
|
||||
"Next, we'll target a more complex example with the Navier-Stokes equations as physical model. In line with {doc}`overview-ns-forw`, we'll target a 2D case with velocity $\\mathbf{u}$, no explicit viscosity term, and a marker density $d$ that drives a simple Boussinesq buoyancy term $\\eta d$ adding a force along the y dimension. For the velocity this gives:\n",
|
||||
"\n",
|
||||
"$\\begin{aligned}\n",
|
||||
" \\frac{\\partial u_x}{\\partial{t}} + \\mathbf{u} \\cdot \\nabla u_x &= - \\frac{1}{\\rho} \\nabla p \n",
|
||||
@@ -20,10 +20,10 @@
|
||||
" \\text{s.t.} \\quad \\nabla \\cdot \\mathbf{u} &= 0,\n",
|
||||
"\\end{aligned}$\n",
|
||||
"\n",
|
||||
"With an additional transport equation for the passively advected marker density $d$:\n",
|
||||
"With an additional transport equation for the passively advected marker density $s$:\n",
|
||||
"\n",
|
||||
"$\\begin{aligned}\n",
|
||||
" \\frac{\\partial d}{\\partial{t}} + \\mathbf{u} \\cdot \\nabla d &= 0 \n",
|
||||
" \\frac{\\partial s}{\\partial{t}} + \\mathbf{u} \\cdot \\nabla s &= 0 \n",
|
||||
"\\end{aligned}$\n",
|
||||
"\n",
|
||||
"\n"
|
||||
@@ -39,7 +39,9 @@
|
||||
"\n",
|
||||
"With the notation from {doc}`overview-equations` this gives the minimization problem \n",
|
||||
"\n",
|
||||
"$\\text{arg min}_{\\mathbf{u}_{0}} \\sum_i |f(x_{t_e,i} ; \\theta)-y^*_{t_e,i}|^2$ , \n",
|
||||
"$$\n",
|
||||
"\\text{arg min}_{\\mathbf{u}_{0}} \\sum_i |f(x_{t_e,i} ; \\theta)-y^*_{t_e,i}|^2 ,\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"where $y^*_{t_e,i}$ are samples of the reference solution at a targeted time $t_e$,\n",
|
||||
"and $x_{t_e,i}$ denotes the estimate of our simulator at the same\n",
|
||||
@@ -47,8 +49,11 @@
|
||||
"The index $i$ here runs over all discretized, spatial degrees of freedom in our fluid solver (we'll have $32\\times40$ below).\n",
|
||||
"\n",
|
||||
"In contrast to before, we are not dealing with pre-computed quantities anymore, but now $x_{t_e,i}$ is a complex, non-linear function itself.\n",
|
||||
"More specifically, the simulator starts with the initial velocity $\\mathbf{u}_0$ and density $d_0$ to compute the $x_{t_e,i}$, by $n$ evaluations of the discretized pde $\\mathcal{P}$.\n",
|
||||
"This gives as simulated state and quantitiy to optimize: $x_{t_e,i} = d_{t_e} = \\mathcal{P}^n(\\mathbf{u}_0,d_0)$. We will leave $d_0$ fixed in the following, and focus on $\\mathbf{u}_0$ as our degrees of freedom.\n",
|
||||
"More specifically, the simulator starts with the initial velocity $\\mathbf{u}_0$ and density $s_0$ to compute the $x_{t_e,i}$, by $n$ evaluations of the discretized PDE $\\mathcal{P}$.\n",
|
||||
"This gives as simulated final state \n",
|
||||
"$y_{t_e,i} = s_{t_e} = \\mathcal{P}^n(\\mathbf{u}_0,s_0)$, where we will leave $s_0$ fixed in the following, and focus on $\\mathbf{u}_0$ as our degrees of freedom.\n",
|
||||
"Hence, the optimization can only change $\\mathbf{u}_0$ to align $y_{t_e,i}$ with the references $y^*_{t_e,i}$ as closely as possible.\n",
|
||||
"\n",
|
||||
"\n"
|
||||
]
|
||||
},
|
||||
@@ -58,13 +63,13 @@
|
||||
"source": [
|
||||
"## Start-up\n",
|
||||
"\n",
|
||||
"First, let's get the loading of python modules out of the way. By importing `phi.torch.flow` below, we get fluid simulation functions that work within pytorch graphs and can provide gradients (`phi.tf.flow` would be the alternative for tensorflow).\n",
|
||||
"First, let's get the loading of python modules out of the way. By importing `phi.torch.flow`, we get fluid simulation functions that work within pytorch graphs and can provide gradients (`phi.tf.flow` would be the alternative for tensorflow).\n",
|
||||
"\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 1,
|
||||
"execution_count": 13,
|
||||
"metadata": {
|
||||
"id": "da1uZcDXdVcF"
|
||||
},
|
||||
@@ -73,8 +78,8 @@
|
||||
"#!pip install --upgrade --quiet phiflow\n",
|
||||
"#!pip install --upgrade --quiet git+https://github.com/tum-pbs/PhiFlow@develop\n",
|
||||
"\n",
|
||||
"from phi.torch.flow import * # The Dash GUI is not supported on Google Colab, ignore the warning\n",
|
||||
"import pylab\n"
|
||||
"from phi.torch.flow import * \n",
|
||||
"import pylab # for visualizations later on\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -85,15 +90,15 @@
|
||||
"source": [
|
||||
"## Batched simulations\n",
|
||||
"\n",
|
||||
"Now we can set up the simulation, which will work in line with the previous \"regular\" simulation example from the {doc}`overview-ns-forw`. However, now we'll directly include an additional dimension, in line with a mini-batch used for NN training. For this, we'll include an additional named dimension called `inflow_loc`. This dimension will exist \"above\" the previous spatial dimensions `x`, `y` and the channel dimensions `vector`. As indicated by the name `inflow_loc`, the main differences for this dimension will lie in different locations of the inflow, in order to obtain different flow simulations. The named dimensions in phiflow make it very convenient to broadcast information\n",
|
||||
"Now we can set up the simulation, which will work in line with the previous \"regular\" simulation example from the {doc}`overview-ns-forw`. However, now we'll directly include an additional dimension, in line with a mini-batch used for NN training. For this, we'll include an additional named dimension called `inflow_loc`. This dimension will exist \"above\" the previous spatial dimensions `y`, `x` and the channel dimensions `vector`. As indicated by the name `inflow_loc`, the main differences for this dimension will lie in different locations of the inflow, in order to obtain different flow simulations. The named dimensions in phiflow make it very convenient to broadcast information across matching dimensions in different tensors.\n",
|
||||
"\n",
|
||||
"The `Domain` object is allocated just like before, but the `INFLOW_LOCATION` tensor now receives a string\n",
|
||||
"`'inflow_loc,vector'` that indicates the names of the two dimensions.\n"
|
||||
"`'inflow_loc,vector'` that indicates the names of the two dimensions. This leads to the creation of an `inflow_loc` dimensions in addition to the two spatial dimensions (the `vector` part).\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 12,
|
||||
"execution_count": 2,
|
||||
"metadata": {
|
||||
"id": "mphMP0sYIOz-"
|
||||
},
|
||||
@@ -118,7 +123,7 @@
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The last `print` statement verifies that our `INFLOW` grid likewise has an `inflow_loc` dimension in addition to the two spatial ones. You can test for the existence of a tensor dimension in phiflow with the `.exists` boolean, which can be evaluated for any dimension name. E.g., above `INFLOW.inflow_loc.exists` will give `True`, while `INFLOW.some_unknown_dim.exists` will give `False`.\n",
|
||||
"The last `print` statement verifies that our `INFLOW` grid likewise has an `inflow_loc` dimension in addition to the spatial dimenions `x` and `y`. You can test for the existence of a tensor dimension in phiflow with the `.exists` boolean, which can be evaluated for any dimension name. E.g., above `INFLOW.inflow_loc.exists` will give `True`, while `INFLOW.some_unknown_dim.exists` will give `False`.\n",
|
||||
"\n",
|
||||
"Phiflow tensors are automatically broadcast to new dimensions via their names, and hence typically no-reshaping operations are required. E.g., you can easily add or multiply tensors with differing dimensions. Below we'll multiply a staggered grid with a tensor of ones along the `inflow_loc` dimension to get a staggered velocity that has `x,y,inflow_loc` as dimensions via `DOMAIN.staggered_grid(0) * math.ones(inflow_loc=4)`."
|
||||
]
|
||||
@@ -132,12 +137,12 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 13,
|
||||
"execution_count": 6,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABHgAAAFtCAYAAAB80d3MAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA8bElEQVR4nO3deZhcBZX38d/p7ux7SMjGDlEWhYAxgjCsLoACLujgjPMyMzqICuqr44yir+I64os6LoiCaMBBBVkGnBdFwACyCAQIEBKWhARISAgh6aQ7S2913j/qRpvYdW53Lbfqdn8/z5Mn3fWre+/p231P3zp9q8rcXQAAAAAAAMivpnoXAAAAAAAAgMow4AEAAAAAAMg5BjwAAAAAAAA5x4AHAAAAAAAg5xjwAAAAAAAA5BwDHgAAAAAAgJxjwDOEmdk0M7vTzNrM7Ftmdr6Z/VcG273bzA6t9XZSajjYzO6pZw3AUDPEe840M1tqZiPqWQcwlAzxnsN5DoaUIX68Z3qOsfO+zmKb6D8GPIOMma00szf18+5nSVovaby7f6qGZf2ZmZ0iqc3dH04+f42Z3Wxm683M+7j/ZDO73sy2mNmzZvZ3A9jWfDPrNLP2Xv+aJcndH5XUmtQDoEyDqeeY2QgzuyzpNW1mtsjMThrAts43s66des4+kuTuL0paoOI+AFCmwdRzdlputpltH8gDUs5zMNgNxuPdzM5IhjFbzGy5mf1NP7dV1XOMZN9uS9azNuknY/u5eOb7Gv3HgGdo21PSEncvecJRA2dL+nmvz7skXS3pAyXuf5GkTknTJP29pIvN7KABbO+b7j6217+eXtmVkj40gHUBqEyj95wWSc9LOkbSBEmfl3S1me01gO1dtVPPeaZXRs8BstXoPae3iyQ9UMb2OM8Bihr+eDezN0u6QNI/SRon6WhJz/R13xKqfY5xiruPlTRH0qGSPtvP5cre12bWMtBlMDAMeAYxM/tHM7vLzC40s41mtmLHX6PNbL6kMyX9WzK5/avpuJmdamaPm1mrmd1uZgckt/+Tmf2m1/2eNrNf9/r8eTOb08f6hks6XtIdO25z9yfd/TJJj/dx/zGS3i3p/7h7u7vfJelGSf9Q5i7Z2e2STsjqckZgsMt7z3H3Le5+vruvdPeCu/+PpBWSXlf2Tnml+yTtY2Z7Vml9wJCW957Ta7kzJLVKum2g+yDF7eI8B4PEIDnevyTpy+7+p+Q8Y7W7ry5jd/Sl7HMMd18r6WYVBz2SJDM73MzuSfbXI2Z2bHL7fO20r82sycw+k1yR9LKZXW1mk5P772VmbmYfMLPnJP0huf2frXgl00YrXvW0Z69tu5mdnXwvWs3sIjOzXvm/JMu2mdkSMzssuX2mmV1rZi8lPx8fG+i+GAwY8Ax+b5D0pKQpkr4p6TIzM3f/RxUnvTv+8nNr74XM7FWSfinpE5KmSrpJ0m+SZnaHpL9JDuaZkoZLOiJZbh9JYyU92kctsyUV3H1VP2t/laRud3+q122PSBrIFTwfMbMNZvagmb27d5A01C5Jrx7A+gDE8txzXsHMpqnYh0o+MOvDKUnPedzMPtw7cPduScskHVJOPQD6lOueY2bjJX1Z0if7/RW/Euc5GEpye7xb8emTcyVNNbNlZrbKzH5gZqMG8PXX5BzDzHaTdFKyvMxslqT/J+mrkiZL+ldJ15rZ1BL7+lxJ71DxCuiZkjaqeFVib8dIOkDSW83sNEnnSXqXit+PP6r4/ent7ZJeL+lgSe+V9NaktvdIOl/S/5I0XtKpkl42syZJv1HxseIsSSdI+oSZvXWg+yPvGPAMfs+6+6XJJbuXS5qh4tOd0vytpP/n7re4e5ekCyWNkvTG5HLANhWnvEerOPF9wcz2V/Hg/aO7F/pY58Rkuf4aK2nzTrdtUvGSxv74norNd1dJ/0fSfDM7cqf7tCV1AaiOPPecPzOzYSqewFzu7k/0c7GrVTx5mSrpXyR9wczet9N96DlAdeW953xF0mVlDqI5z8FQk+fjfZqkYZJOl/Q3+svToj7fz+VrcY7x32bWpuLT09dJ+mJy+/sl3eTuNyVXGt0iaaGkk0us52xJn3P3Ve7eoeIA5nR75dOxzk+ulN6W3P8/3H1pMpj6uqQ5O1199A13b3X351R8faE5ye0fVHG49IAXLXP3Z1UcBk119y+7e2fyfb1U0hkD2B+DAs+BG/zW7vjA3bcmV7f15wW0Zkp6tteyBTN7XsWJqFScdh8rab/k41YVm+AR6nWp4k42qv/DGUlqV3Ey29t49bOZuvtDvT69ycyuVHFSfHev28epWDuA6shzz5EkJX8F+rmKr/91Tn+Xc/clvT69x8y+q+KJXO+/StFzgOrKbc9JnvbxJhUf5A0Y5zkYgnJ7vEvalvz/fXdfI0lm9m0VBzyfS1u4RucY73D3W83sGEm/UPHKqFYVX2PnPfbKF2kfpuKgpS97SrrezHoPwnr0yuHb8zvd/7v2ynfgMhW/Hzu+T2t7ZVv1l+/z7pKWl6hhppm19rqtWcWrg4YUruBBKS+oeKBIkpLnPe4uacfzRHc0wr9JPr5DxUZ4jEo3wmXJqmaVyHf2lKQWM5vd67ZDNLCnS/TmKjYPSX++/HC4ipd6AqivRug5O7Z7mYonJe9O/tJXrp17TouKJ4+PVLBOANXRCD3nWEl7SXrOzNaq+DSId5vZQ9FCAc5zgL7V/Xh3942SVql4nP755n5/BX2sUlU6x3D3OyTNV/HKJqk4jPm5u0/s9W+Mu3+jxCqel3TSTvcfudPrC/lO9//QTvcf5e739KPc5yXtW+L2FTutc5y7l7rqaNBiwINSrpb0NjM7IXmqwqckdUjaceDdIek4SaOSy4r/KOlESbtIerivFbp7p6RbVWyWkopd0cxGqngCIjMbacmLAbr7FknXSfqymY1JLjs+Tcmr1dtfXrRrr762Z2anm9nY5Dm1b1HxcsMbe93lGEl/SC4lBFBfde85iYtVvAT6lOQy4ldIes6xfW3PzE4zs0nJNuZJ+pikG3rdZZ6klcmlxADqqxF6ziUqPlCZk/z7kYqve7HjtSY4zwGqoxGOd0n6maRzzWxXM5sk6X9L+p9ey1ftHMPMjrUSb9Vewn9KerOZHSLpv1R8vZ+3mllz8nUca8XX6unLjyR9zZKnWJnZVCu+zk4pP5L0WUveGdnMJljxtXX64yeS/tXMXpfsi/2S7d4vqc3M/t3MRiV1v8bMXt/P9Q4aDHjQJ3d/UsUThe9LWi/pFBUf8HQm+VMqPoXqj8nnm1V8m7+7/ZVv0bmzH+uV74K1p4qXLO64KmebXvmXpo+o+BzZdSpegvhhd99x391VvIyv1KvPfzzJWiX9X0n/4u6398r/XsUGA6DOGqHnJCcIH1LxgdZaK747RLuZ/X2S767iU0QfK7GtM1T8i16bpCskXeDul/fK6TlAg2iEnuPuW9197Y5/yfa2u/tLyX05zwGqoBGO98RXJD2g4rMUlqo4PPqaVJNzjN31lwFWqqTvXCHpC+7+vIp/VD9P0ksqXh3zaZWeHXxXxeHy7634mj5/UvEFsUtt63oV3y7+V2a2WdJiFV/kuT91/lrFffYLFffFf0uanHyf3q7iOdwKFb/PP5E0oT/rHUzMB/729UBFzOxuSee4e58T8QGs5/OSXnL3H5ex7MGSfuzuR1RSA4DGV8We835JB7n7Z8tYdlcV/0J4qLtvr6QOAI2N8xxg6GjUcwwz+4mkX7v7zZXUhfxhwAMAAAAAAJBzPEULAAAAAAAg5xjwAAAAAAAA5BwDHgAAAAAAgJxryXJjA3yrNgCDx3p3n5r1Ruk5wJBVl54j0XeAocrdrR7bpecAQ1af5zqZDniKmrPfJIA663m2ftum5wBDTz17jkTfAYaa6J20s0DPAYaevs91eIoWAAAAAABAzjHgAQAAAAAAyDkGPAAAAAAAADnHgAcAAAAAACDnGPAAAAAAAADkHAMeAAAAAACAnGPAAwAAAAAAkHMMeAAAAAAAAHKOAQ8AAAAAAEDOMeABAAAAAADIOQY8AAAAAAAAOceABwAAAAAAIOcY8AAAAAAAAOQcAx4AAAAAAICcY8ADAAAAAACQcwx4AAAAAAAAco4BDwAAAAAAQM4x4AEAAAAAAMg5BjwAAAAAAAA5x4AHAAAAAAAg5xjwAAAAAAAA5BwDHgAAAAAAgJxjwAMAAAAAAJBzDHgAAAAAAAByjgEPAAAAAABAzjHgAQAAAAAAyDkGPAAAAAAAADnHgAcAAAAAACDnGPAAAAAAAADkHAMeAAAAAACAnGPAAwAAAAAAkHOpAx4zG2lm95vZI2b2uJl9Kbl9vpmtMLNFyb85Na8WwKBHzwGQNfoOgCzRcwDUSks/7tMh6Xh3bzezYZLuMrPfJtmn3f2a2pUHYAii5wDIGn0HQJboOQBqInXA4+4uqT35dFjyz2tZFIChi54DIGv0HQBZoucAqJV+vQaPmTWb2SJJ6yTd4u73JdHXzOxRM/uOmY0osexZZrbQzBZWp2QAgx09B0DW6DsAskTPAVALVhwg9/POZhMlXS/pXEkvS1orabikSyQtd/cvpyzvUnPZxQLIq54H3X3uQJei5wAoT3k9R6LvAChHj9zdylmSngOgPH2f6wzoXbTcvVXSAkknuvsaL+qQ9DNJ86pSJwAk6DkAskbfAZAleg6AaurPu2hNTSbLMrNRkt4s6Qkzm5HcZpLeIWlx7coEMFTQcwBkjb4DIEv0HAC10p930Zoh6XIza1ZxIHS1u/+Pmf3BzKZKMkmLJJ1duzIBDCH0HABZo+8AyBI9B0BNDOg1eCreGM8RBYao8l8PoxL0HGCoqk/Pkeg7wNBU/mvwVIqeAwxVVXgNHgAAAAAAADQeBjwAAAAAAAA5x4AHAAAAAAAg5xjwAAAAAAAA5BwDHgAAAAAAgJxjwAMAAAAAAJBzDHgAAAAAAAByjgEPAAAAAABAzjHgAQAAAAAAyDkLine truncated
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABHgAAAFtCAYAAAB80d3MAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA9V0lEQVR4nO3deZhcBZX38d/pJZ2dJGQjISQsYVXWgKAguKDAOwoqOjg6wygKuDH66jjoqKMo76jj7usr4iBBVBRZRIURUZFVwQQSIAGSAAnZmyTdSSedpNNV5/2jKtKJXae6tlt1u7+f5+kn3fWre+/p231Pbp2+VWXuLgAAAAAAAKRXU70LAAAAAAAAQGUY8AAAAAAAAKQcAx4AAAAAAICUY8ADAAAAAACQcgx4AAAAAAAAUo4BDwAAAAAAQMox4BnCzGyKmd1rZl1m9lUz+6yZ/SiB7T5gZsfVejtFajjazB6sZw3AUDPEe84UM3vSzNrqWQcwlAzxnsN5DoaUIX68J3qOsfe+TmKbGDgGPIOMmS03s9cO8O4XS9ogaay7f7SGZf2Vmb1BUpe7P5r/+iVmdqeZbTAz7+f+E8zsVjPbZmYrzOwfStjWXDPrMbOtfT6aJcndH5PUma8HQJkGU88xszYzuybfa7rMbIGZnV3Ctj5rZrv26jkHSZK7r5d0t3L7AECZBlPP2Wu52Wa2o5QHpJznYLAbjMe7mV2QH8ZsM7NnzOy0AW6rqucY+X27Pb+edfl+MnqAiye+rzFwDHiGtpmSFrt7wROOGrhU0vV9vt4l6UZJFxW4/3ck9UiaIukdkr5rZkeVsL0vu/voPh+ZPtmPJV1SwroAVKbRe06LpJWSTpe0j6RPSbrRzGaVsL2f7dVznu2T0XOAZDV6z+nrO5L+Usb2OM8Bchr+eDezMyV9SdK7JI2R9EpJz/Z33wKqfY7xBncfLelYScdJ+sQAlyt7X5tZS6nLoDQMeAYxM/tnM7vfzL5iZh1m9tzuv0ab2VxJF0r6eH5y+zfTcTN7o5ktMrNOM/ujmR2Rv/1dZvarPvdbamY/7/P1SjM7tp/1DZP0akn37L7N3Z9292skLern/qMkvUXSp919q7vfL+mXkv6xzF2ytz9Kek1SlzMCg13ae467b3P3z7r7cnfPuvuvJT0n6YSyd8qeHpJ0kJnNrNL6gCEt7T2nz3IXSOqU9PtS90ERfxTnORgkBsnx/jlJV7j7n/PnGavdfXUZu6M/ZZ9juPs6SXcqN+iRJJnZyWb2YH5/LTSzM/K3z9Ve+9rMmszs8vwVSRvN7EYzm5C//ywzczO7yMyel/SH/O3vttyVTB2Wu+ppZp9tu5ldmv9ZdJrZd8zM+uTvzS/bZWaLzez4/O3TzOxmM3sh//txWan7YjBgwDP4vUzS05ImSvqypGvMzNz9n5Wb9O7+y8/v+i5kZodKukHShyVNknSHpF/lm9k9kk7LH8zTJA2TdEp+uYMkjZb0WD+1zJaUdfdVA6z9UEm97r6kz20LJZVyBc/7zWyTmc03s7f0DfINdZekw0pYH4BYmnvOHsxsinJ9qOADs368Id9zFpnZ+/oG7t4raZmkY8qpB0C/Ut1zzGyspCsk/e8Bf8d74jwHQ0lqj3fLPX1yjqRJZrbMzFaZ2f81sxElfP81Occws/0lnZ1fXmY2XdLtkr4gaYKkj0m62cwmFdjXH5J0nnJXQE+T1KHcVYl9nS7pCEmvN7NzJX1S0puV+3ncp9zPp6+/k3SipKMlvU3S6/O1vVXSZyX9k6Sxkt4oaaOZNUn6lXKPFadLeo2kD5vZ60vdH2nHgGfwW+Hu389fsnudpP2Ue7pTMX8v6XZ3v8vdd0n6iqQRkl6evxywS7kp7yuVm/iuMbPDlTt473P3bD/rHJdfbqBGS9qy122blbukcSC+pVzznSzp05Lmmtkr9rpPV74uANWR5p7zV2bWqtwJzHXu/tQAF7tRuZOXSZLeK+kzZvb2ve5DzwGqK+095/OSrilzEM15DoaaNB/vUyS1Sjpf0ml68WlRnxrg8rU4x/iFmXUp9/T0dkn/kb/9nZLucPc78lca3SVpnqRzCqznUkn/7u6r3H2ncgOY823Pp2N9Nn+l9Pb8/f/T3Z/MD6b+j6Rj97r66Ivu3unuzyv3+kLH5m9/j3LDpb94zjJ3X6HcMGiSu1/h7j35n+v3JV1Qwv4YFHgO3OC3bvcn7t6dv7ptIC+gNU3Sij7LZs1spXITUSk37T5D0iH5zzuVa4KnqM+linvp0MCHM5K0VbnJbF9jNcBm6u6P9PnyDjP7sXKT4gf63D5GudoBVEeae44kKf9XoOuVe/2vDw50OXdf3OfLB83sm8qdyPX9qxQ9B6iu1Pac/NM+Xqvcg7yScZ6DISi1x7uk7fl/v+3uayXJzL6m3IDn34stXKNzjPPc/Xdmdrqknyh3ZVSncq+x81bb80XaW5UbtPRnpqRbzazvICyjPYdvK/e6/zdtz3fgMuV+Hrt/Tuv6ZN168ec8Q9IzBWqYZmadfW5rVu7qoCGFK3hQyBrlDhRJUv55jzMk7X6e6O5GeFr+83uUa4Snq3AjXJZf1fQC+d6WSGoxs9l9bjtGpT1doi9XrnlI+uvlh8OUu9QTQH01Qs/Zvd1rlDspeUv+L33l2rvntCh38riwgnUCqI5G6DlnSJol6XkzW6fc0yDeYmaPRAsFOM8B+lf3493dOyStUu44/evNA/4O+lmlqnSO4e73SJqr3JVNUm4Yc727j+vzMcrdv1hgFSslnb3X/Yfv9fpCvtf9L9nr/iPc/cEBlLtS0sEFbn9ur3WOcfdCVx0NWgx4UMiNkv6Xmb0m/1SFj0raKWn3gXePpFdJGpG/rPg+SWdJ2lfSo/2t0N17JP1OuWYpKdcVzWy4cicgMrPhln8xQHffJukWSVeY2aj8ZcfnKv9q9fbii3bN6m97Zna+mY3OP6f2dcpdbvjLPnc5XdIf8pcSAqivuvecvO8qdwn0G/KXEe8h33PO6G97ZnaumY3Pb+MkSZdJuq3PXU6StDx/KTGA+mqEnnO1cg9Ujs1/XKXc617sfq0JznOA6miE412SrpX0ITObbGbjJX1E0q/7LF+1cwwzO8MKvFV7Ad+QdKaZHSPpR8q93s/rzaw5/32cYbnX6unPVZKutPxTrMxskuVeZ6eQqyR9wvLvjGxm+1jutXUG4r8lfczMTsjvi0Py231YUpeZ/ZuZjcjX/RIzO3GA6x00GPCgX+7+tHInCt+WtEHSG5R7wNOTz5co9xSq+/Jfb1Hubf4e8D3fonNv39Oe74I1U7lLFndflbNde/6l6f3KPUe2XblLEN/n7rvvO0O5y/gKvfr8v+SzTkn/Jem97v7HPvk7lGswAOqsEXpO/gThEuUeaK2z3LtDbDWzd+TzGco9RfTxAtu6QLm/6HVJ+qGkL7n7dX1yeg7QIBqh57h7t7uv2/2R394Od38hf1/Oc4AqaITjPe/zkv6i3LMUnlRueHSlVJNzjBl6cYBVVL7v/FDSZ9x9pXJ/VP+kpBeUuzrmX1V4dvBN5YbLv7Xca/r8WbkXxC60rVuVe7v4n5rZFklPKPcizwOp8+fK7bOfKLcvfiFpQv7n9HfKncM9p9zP+b8l7TOQ9Q4m5qW/fT1QETN7QNIH3b3fiXgJ6/mUpBfc/XtlLHu0pO+5+ymV1ACg8VWx57xT0lHu/okylp2s3F8Ij3P3HZXUAaCxcZ4DDB2Neo5hZv8t6efufmcldSF9GPAAAAAAAACkHE/RAgAAAAAASDkGPAAAAAAAACnHgAcAAAAAACDlWpLcWIlv1QZg8Njg7pOS3ig9Bxiy6tJzJPoOMFS5u9Vju/QcYMjq91wn0QFPTnPymwRQZ5kV9ds2PQcYeurZcyT6DjDURO+knQR6DjD09H+uw1O0AAAAAAAAUo4BDwAAAAAAQMox4AEAAAAAAEg5BjwAAAAAAAApx4AHAAAAAAAg5RjwAAAAAAAApBwDHgAAAAAAgJRjwAMAAAAAAJByDHgAAAAAAABSjgEPAAAAAABAyjHgAQAAAAAASDkGPAAAAAAAACnHgAcAAAAAACDlGPAAAAAAAACkHAMeAAAAAACAlGPAAwAAAAAAkHIMeAAAAAAAAFKOAQ8AAAAAAEDKMeABAAAAAABIOQY8AAAAAAAAKceABwAAAAAAIOUY8AAAAAAAAKQcAx4AAAAAAICUY8ADAAAAAACQcgx4AAAAAAAAUo4BDwAAAAAAQMox4AEAAAAAAEg5BjwAAAAAAAApx4AHAAAAAAAg5RjwAAAAAAAApBwDHgAAAAAAgJQrOuAxs+Fm9rCZLTSzRWb2ufztc83sOTNbkP84tubVAhj06DkAkkbfAZAkeg6AWmkZwH12Snq1u281s1ZJ95vZ/+Szf3X3m2pXHoAhiJ4DIGn0HQBJoucAqImiAx53d0lb81+25j+8lkUBGLroOQCSRt8BkCR6DoBaGdBr8JhZs5ktkNQu6S53fygfXWlmj5nZ182srcCyF5vZPDObV52SAQx29BwASaPvAEgSPQdALVhugDzAO5uNk3SrpA9J2ihpnaRhkq6W9Iy7X1FkeZeayy4WQFpl5rv7nFKXoucAKE95PUei7wAoR0bubuUsSc8BUJ7+z3VKehctd++UdLeks9x9refslHStpJOqUicA5NFzACSNvgMgSfQcANU0kHfRmpSfLMvMRkg6U9JTZrZf/jaTdJ6kJ2pXJoChgp4DIGn0HQBJoucAqJWBvIvWfpKuM7Nm5QZCN7r7r83sD2Y2SZJJWiDp0tqVCWAIoecASBp9B0CS6DkAaqKk1+CpeGM8RxQYosp/PYxK0HOAoao+PUei7wBDU/mvwVMpeg4wVFXhNXgAAAAAAADQeBjwAAAAAAAApBwDHgAAAAAAgJRjwAMAAAAAAJByDHgAAAAAAABSjgEPAAAAAABAyjHgAQAAAAAASDkGPAAAAAAAACnHgAcAAAAAACDLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x360 with 4 Axes>"
|
||||
]
|
||||
@@ -156,7 +161,7 @@
|
||||
" smoke = advect.mac_cormack(smoke, velocity, dt=1) + INFLOW\n",
|
||||
" buoyancy_force = smoke * (0, 1) >> velocity\n",
|
||||
" velocity = advect.semi_lagrangian(velocity, velocity, dt=1) + buoyancy_force\n",
|
||||
" velocity, _, _, _ = fluid.make_incompressible(velocity, DOMAIN)\n",
|
||||
" velocity, _ = fluid.make_incompressible(velocity, DOMAIN)\n",
|
||||
" return smoke, velocity\n",
|
||||
"\n",
|
||||
"for _ in range(20):\n",
|
||||
@@ -166,7 +171,7 @@
|
||||
"fig, axes = pylab.subplots(1, 4, figsize=(16, 5))\n",
|
||||
"for i in range(INFLOW.shape.inflow_loc):\n",
|
||||
" axes[i].imshow(smoke_final.values.inflow_loc[i].numpy('y,x'), origin='lower', cmap='magma')\n",
|
||||
" axes[i].set_title(f\"Inflow {math.to_int(INFLOW_LOCATION).inflow_loc[i].vector.unstack_spatial('x,y')}\" + (\", Reference\" if i==3 else \"\"))\n",
|
||||
" axes[i].set_title(f\"Inflow {INFLOW_LOCATION.inflow_loc[i].vector.unstack_spatial('x,y')}\" + (\", Reference\" if i==3 else \"\"))\n",
|
||||
"pylab.tight_layout()\n"
|
||||
]
|
||||
},
|
||||
@@ -176,7 +181,7 @@
|
||||
"source": [
|
||||
"The last image shows the state of the advected smoke fields after 20 simulation steps. The final smoke shape of simulation `[3]` with an inflow at `(16,5)`, with the straight plume on the far right, will be our **reference state** below. The initial velocity of the other three will be modified in the optimization procedure below to match this reference.\n",
|
||||
"\n",
|
||||
"(As a small side note: phiflow tensors will keep track of their chain of operations using the backend they were created for. E.g. a tensor created with NumPy will keep using NumPy/SciPy operations unless a TensorFlow tensor is also passed to the same operation. Thus, it is a good idea to verify that tensors have the right `type` once in a while.)\n"
|
||||
"(As a small side note: phiflow tensors will keep track of their chain of operations using the backend they were created for. E.g. a tensor created with NumPy will keep using NumPy/SciPy operations unless a PyTorch or TensorFlow tensor is also passed to the same operation. Thus, it is a good idea to verify that tensors have the right `type` once in a while.)\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -187,17 +192,21 @@
|
||||
"\n",
|
||||
"Let's look at how to get gradients from our simulation. The first trivial step taken care of above was to include `phi.torch.flow` to import differentiable operators from which to build our simulator.\n",
|
||||
"\n",
|
||||
"Now we want to optimize the initial velocities so that all simulations arrive at a final state that is similar to the right simulation where the inflow is located at `(16, 5)`.\n",
|
||||
"To achieve this, we record the gradients during the simulation and define a simple $L^2$ based loss function. The loss function we'll use is given by $L = | u - u^* |^2$, where $s$ denotes the smoke density, and $u^*$ the reference state from the fourth simulation above. When evaluating the loss function we treat this reference state as an external constant via `field.stop_gradient()`.\n",
|
||||
"Now we want to optimize the initial velocities so that all simulations arrive at a final state that is similar to the simulation no the right, where the inflow is located at `(16, 5)`, i.e. centered along `x`.\n",
|
||||
"To achieve this, we record the gradients during the simulation and define a simple $L^2$ based loss function. The loss function we'll use is given by $L = | s_{t_e} - s_{t_e}^* |^2$, where $s_{t_e}$ denotes the smoke density, and $s_{t_e}^*$\n",
|
||||
"denotes the reference state from the fourth simulation in our batch (both evaluated at the last time step $t_e$).\n",
|
||||
"When evaluating the loss function we treat the reference state as an external constant via `field.stop_gradient()`.\n",
|
||||
"As outlined at the top, $s$ is a function of $\\mathbf{u}$ (via the advection equation), which in turn is given by the Navier-Stokes equations. Thus, via a chain of multiple time steps $s$ depends in the initial velocity state $\\mathbf{u}_0$.\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"It is important that our initial velocity has the `inflow_loc` dimension before we record the gradients, such that we have the full \"mini-batch\" of four versions of our velocity (three of which will be updated via gradients in our optimization later on). To get the approprate velocity tensor, we multiply a `staggered_grid` tensor from our domain object, with a tensor full of ones along the `inflow_loc` dimension.\n",
|
||||
"\n",
|
||||
"phiflow provides a unified API for gradients across different platforms by using functions that need to return a loss values, in addition to optional state values. Below the `simulate` function computes the $L^2$ error outlined above and returns the evolved `smoke` and `velocity` states are 20 simulation steps.\n"
|
||||
"Phiflow provides a unified API for gradients across different platforms by using functions that need to return a loss values, in addition to optional state values. Below the `simulate` function computes the $L^2$ error outlined above and returns the evolved `smoke` and `velocity` states after 20 simulation steps.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 4,
|
||||
"execution_count": 7,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
@@ -229,37 +238,39 @@
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Then, we use `field.functional_gradient()` to obtain the gradient with respect to the initial velocity. Since the velocity is the second argument of the `simulate()` function, we pass `wrt=[1]`. (Note that there's also a `field.spatal_gradient` with computes derivatives of tensors along spatial dimensions like `x,y`.)\n",
|
||||
"Phiflow's `field.functional_gradient()` function is the central function to compute gradients. Next, we'll use it to obtain the gradient with respect to the initial velocity. Since the velocity is the second argument of the `simulate()` function, we pass `wrt=[1]`. (Note that phiflow also has a `field.spatal_gradient` function which instead computes derivatives of tensors along spatial dimensions like `x,y`.)\n",
|
||||
"\n",
|
||||
"`functional_gradient` returns a gradient function, which we can evaluate with the initial states for smoke and velocity. \n"
|
||||
"`functional_gradient` in generates a gradient function, which we evaluate one term for test purposes with the initial states for smoke and velocity. The last statement prints a summary of one part of the tensor:\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 5,
|
||||
"execution_count": 12,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"Gradient info StaggeredGrid[(inflow_loc=4, x=32, y=40, vector=2), size=(32, 40) along vector, extrapolation=0]\n"
|
||||
"ename": "RuntimeError",
|
||||
"evalue": "grad can be implicitly created only for scalar outputs",
|
||||
"output_type": "error",
|
||||
"traceback": [
|
||||
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
|
||||
"\u001b[0;31mRuntimeError\u001b[0m Traceback (most recent call last)",
|
||||
"\u001b[0;32m<ipython-input-12-e357726b3491>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m 6\u001b[0m \u001b[0;31m#[velocity_grad] = sim_grad(initial_smoke, initial_velocity)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 7\u001b[0m \u001b[0;31m#(velocity_grad,) = sim_grad(initial_smoke, initial_velocity)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 8\u001b[0;31m \u001b[0msim_grad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minitial_smoke\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0minitial_velocity\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 9\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 10\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"Gradient info \"\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mformat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mvelocity_grad\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~/Dropbox/mbaDevelSelected/phiflow-v2/phi/math/_functional.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 300\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mkey\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrads\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 301\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrads\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_trace_grad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mwrt_natives\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 302\u001b[0;31m \u001b[0mnative_result\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrads\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mnatives\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 303\u001b[0m \u001b[0moutput_key\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mmatch_output_signature\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrecorded_mappings\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 304\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_output\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
|
||||
"\u001b[0;32m~/Dropbox/mbaDevelSelected/phiflow-v2/phi/torch/_torch_backend.py\u001b[0m in \u001b[0;36meval_grad\u001b[0;34m(*args)\u001b[0m\n\u001b[1;32m 635\u001b[0m \u001b[0moutput\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 636\u001b[0m \u001b[0mloss\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0maux\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0moutput\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0moutput\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0moutput\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mtuple\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlist\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32melse\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0moutput\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 637\u001b[0;31m \u001b[0mgrads\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtorch\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mautograd\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mloss\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mwrt_args\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;31m# TODO this does not work during jit for some reason\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 638\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mget_output\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 639\u001b[0m \u001b[0mloss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mloss\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdetach\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~/miniconda3/envs/tf/lib/python3.8/site-packages/torch/autograd/__init__.py\u001b[0m in \u001b[0;36mgrad\u001b[0;34m(outputs, inputs, grad_outputs, retain_graph, create_graph, only_inputs, allow_unused)\u001b[0m\n\u001b[1;32m 195\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 196\u001b[0m \u001b[0mgrad_outputs_\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_tensor_or_tensors_to_tuple\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mgrad_outputs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0moutputs\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--> 197\u001b[0;31m \u001b[0mgrad_outputs_\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_make_grads\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0moutputs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mgrad_outputs_\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 198\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 199\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mretain_graph\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
|
||||
"\u001b[0;32m~/miniconda3/envs/tf/lib/python3.8/site-packages/torch/autograd/__init__.py\u001b[0m in \u001b[0;36m_make_grads\u001b[0;34m(outputs, grads)\u001b[0m\n\u001b[1;32m 48\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mout\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrequires_grad\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 49\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mout\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnumel\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m!=\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 50\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mRuntimeError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"grad can be implicitly created only for scalar outputs\"\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 51\u001b[0m \u001b[0mnew_grads\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtorch\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mones_like\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mout\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmemory_format\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mtorch\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpreserve_format\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 52\u001b[0m \u001b[0;32melse\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: grad can be implicitly created only for scalar outputs"
|
||||
]
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"(x=33, y=40) float32 -0.9229158759117126 < ... < 2.964700937271118"
|
||||
]
|
||||
},
|
||||
"execution_count": 5,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"sim_grad = field.functional_gradient(simulate, wrt=[1], get_output=False)\n",
|
||||
"[velocity_grad] = sim_grad(initial_smoke, initial_velocity)\n",
|
||||
"\n",
|
||||
"TODO fix\n",
|
||||
"#[velocity_grad] = sim_grad(initial_smoke, initial_velocity)\n",
|
||||
"#(velocity_grad,) = sim_grad(initial_smoke, initial_velocity)\n",
|
||||
"sim_grad(initial_smoke, initial_velocity)\n",
|
||||
"\n",
|
||||
"print(\"Gradient info \" + format(velocity_grad))\n",
|
||||
"velocity_grad.inflow_loc[0].values.vector[0] # one example, location 0, x component, automatically prints size & content range\n"
|
||||
@@ -276,22 +287,9 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 11,
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABHUAAAFgCAYAAAA1sUWiAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABUEklEQVR4nO3de7hkVX3n//eXprnIRcBuEQHTqBijJGLsoJncjKJiRsVENKgxmKDEGCaT+TmTqMyoQU0wMTEmOo4dQdB4HS8DozgEvIwTE5FGUUFQkEBo5NZc5KJcGr6/P/Y+3dWn9z6n6lStql1V79fz1NOnVu3ae+06VZ9eZ9Xaa0VmIkmSJEmSpOmy06QrIEmSJEmSpMHZqSNJkiRJkjSF7NSRJEmSJEmaQnbqSJIkSZIkTSE7dSRJkiRJkqaQnTqSJEmSJElTyE4dSZIkSZKkHhFxVER8NyKuiIjXNjz+yxHx9YjYEhHHLHrsuIi4vL4d11P+pIj4dr3Pv42IGLaedupIkiRJkiTVImIV8G7g2cDjgBdHxOMWbfZvwMuBDy967n7AG4EnA0cAb4yIfeuH3wO8Eji0vh01bF3t1JEkSZIkSdrmCOCKzLwyM+8FPgoc3btBZl6Vmd8CHlj03GcB52bmLZl5K3AucFREHADsnZlfzcwEPgA8f9iK7jzsDiSN1lERuXnIfVwI52Tm0L2+kmbTsDljxkhaim0ZSaWNoC1zCXB3T9GGzNzQc/9A4Jqe+5uoRt70o+m5B9a3TQ3lQ7FTR+qYzcDGIfcRsGYUdZE0m4bNGTNG0lJsy0gqbQRtmbszc/2o6jNJdupIXbTTkFdGPrB4BKAkLTJMzpgxkpZjW0ZSaWXbMtcCB/fcP6gu68e1wFMXPfdLdflBK9xnKzt1pC6yISSpNDt1JJVkW0ZSaWXbMhcAh0bEIVQdL8cCL+lz7+cAf9YzOfIzgddl5i0RcXtEPAU4H/ht4O8GrvsiTpQsSZIkSZJUy8wtwIlUHTSXAh/PzEsi4uSIeB5ARPxcRGwCXgi8NyIuqZ97C/Bmqo6hC4CT6zKAVwPvA64Avg98bti6OlJH6pqI4b/dkqSlmDOSSjJjJJU2hpzJzLOBsxeVvaHn5wvY/nKq3u1OA05rKN8IHDbKetqpI3WRDSFJpZkzkkoyYySVZs4AdupI3eO3W5JKM2cklWTGSCrNnNnKV0GSJEmSJGkKOVJH6iJ7nSWVZs5IKsmMkVSaOQPYqSN1kwElqTRzRlJJZoyk0swZwE4dqXu8PlRSaeaMpJLMGEmlmTNb2akjdZEBJak0c0ZSSQUzJiIOBj4A7A8ksCEz31nsgJK6ybYMYKeOJEmSpOmyBXhNZn49IvYCLoyIczPzO5OumCSNm506Utc4lFBSaeaMpJIKZ0xmXgdcV/98R0RcChwI2KkjzQvbMlvZqSN1kQElqTRzRlJJw2fMmojY2HN/Q2ZuWLxRRKwDngicP+wBJU0Z2zKAnTpSNxlQkkozZySVNHzGbM7M9UttEBF7Ap8E/igzbx/2gJKmjG0ZAHwVJEmSJE2ViFhN1aHzocz81KTrI0mT4kgdqWu8PlRSaeaMpJIKZ0xEBHAqcGlm/nWxA0nqLtsyW9mpI3WRASWpNHNGUkllM+YXgJcB346Ii+qy12fm2SUPKqljbMsAdupI3WOvs6TSzBlJJZVf/eqfgCh2AEndZ1tmKzt1pC4yoCSVZs5IKsmMkVSaOQM4UbIkSZIkSdJUcqSO1EWFe50j4ijgncAq4H2ZeUrRA0rqHr/dklSSGSOpNHMGsFNH6p7yK0asAt4NPAPYBFwQEWdl5neKHVRSt3gduqSSzBhJpZkzW9mpI3VR2YA6ArgiM68EiIiPAkcDdupI88SGkKSSzBhJpZkzgHPqSLNqTURs7Lmd0PPYgcA1Pfc31WWSJEmSpCniSB2pa0YzlHBzZq4fRXUkzSCHLEsqyYyRVJo5s5WdOlIXlQ2oa4GDe+4fVJdJmic2hCSVZMZIKs2cAezUkbqpbEBdABwaEYdQdeYcC7yk5AEldZANIUklmTGSSjNnADt1pO4pPJQwM7dExInAOVRLmp+WmZcUO6Ck7nHIsqSSzBhJpZkzW9mpI82hzDwbOHvS9ZAkSZIkrZydOlIX2essqTRzRlJJZoyk0swZwE4dqXscSiipNHNGUklmjKTSzJmt5vpViIirIuLI+ufXR8T7xnTciIj3R8StEfG1cRxzFCLi9Ih4ywSP/5WIeOKkjl/X4Wci4p+LH2innYa7qRPMmMGYMRAR+0fEpRGxa/GDmTEzwZwZjDljW0aDMWMGY8bYlpmEzp5NRBwbEedHxF0RcWP986sjIkocLzP/LDNfMex+ImJdRGRELDUK6heBZwAHZeYRwx5zGvT5uiz1/OcCd2TmN+r7h0XEORGxOSKyYfv9IuLT9fvn6ojoe3WnOozvjYg7e26rADLzW8BtdX00xcyY2TLOjImIXSPi1Dpb7oiIiyLi2QMc600Rcd+ijHkkQGbeAHwROGEl56FuMWdmy7jbMj3POzQi7o6IfxjgWLZl5oAZM1smkTH1e+jS+j30/Yj4pT6PZVumQzrZqRMRrwHeCfwl8DBgf+BVwC8Au7Q8Z9XYKji8nwCuysy7Bn3iSj/kM+BVwAd77t8HfBw4vmX7dwP3Ur13Xgq8JyIeP8Dx/iIz9+y53d/z2IeA3xtgX4Pz262izJh2ZsxWS2XMzsA1wK8ADwb+K/DxiFg3wPE+tihjrux5rHzGgBlTmDnTzpzZarm2zIJ3Axes4Hi2ZWaYGdPOjNlqyYyJiGcAbwN+B9gL+GXgyqZtW9iW6YjOnU1EPBg4GXh1Zn4iM+/Iyjcy86WZeU+93ekR8Z6IODsi7gJ+NSL+fUR8IyJuj4hrIuJNi/b9svqb1Zsj4qRFj72p9xuQiHhKRPxzRNwWEd+MiKf2PPaliHhzVMPb7oiIf4yINfXDX67/va3usfz5Rcc5Hngf8PP1439al78yIq6IiFsi4qyIeHjPczIi/iAiLgcub3jNfjMi/jUi9q7vPzsiro+ItQ3bfi6q5ax7y74ZEb9R//zYiDi3rsd3I+JFO/ySRqD+/b07Ij5bv4bnR8SjWrbdBXga8H8XyjLzu5l5KrDDUtwRsQfwAuC/ZeadmflPwFnAy0ZU/S8BT49SQwoXrg81pIowY8yYhm0HypjMvCsz35SZV2XmA5n5GeBfgSeNqPrnA4+MiJ8Y0f52NGzO9HWIOKr+HV8REa8tdi4dZM6YMw3bDpQzPc87FrgN+Pxoa29bZpqZMWZMw7YryZg/BU7OzK/W7ZlrM/PaEVV/Jtoy06KLZ/PzwK7AmX1s+xLgrVQ9i/8E3AX8NrAP8O+B34+I5wNExOOA91D9Yf9w4CHAQU07jYgDgc8CbwH2A/4z8MlFH/qXUPVqPpSqN/w/1+W/XP+7T91j+S+9+64/WK8C/qV+/I0R8TTgz4EXAQcAVwMfXVSt5wNPBh63uL6Z+THgn4G/jYiHAKcCr8jMmxpO7yPAi3vO9XFUPeGfjaoz5Fzgw/V5HQv893qbEo6lCpN9gSuofpdNDgUeyMxNfe73McCWzPxeT9k3gUFG6ry6DuoLI+IFvQ/UYXcf8JMD7G8wNoRKMmPMmMUGzZjtRMT+VLnT+odZg+fWGXNJRPx+7wOZuaWu7xNWUp++FcyYqL4NfjfwbKr31IsL/p67yJwxZxYbOGfqPz5PBv6/FdbNtszsMmPMmMUGypj6/+n1wNq6o2xTRLwrInYfoG4z3ZaB5b+giuqS/I/Vj58f9ajtiHhpVJfnL9weiIjD68e+VO9z4bGHDv0yDLuDAtYAm+s3AgA9PcA/johf7tn2zMz8St2zeHdmfikzv13f/xbVB/JX6m2PAT6TmV+ue6//G/BASx1+Czg7M8+u93UusBH4tZ5t3p+Z38vMH1MNazt8iHN+KXBaZn69rtvrqHqm1/Vs8+eZeUt9vCZ/QNU7+yXgf9ffHDf5NHB4T6/pS4FP1cd9DtUwx/dn5pasrsf8JPDCIc5tKZ/OzK/Vv+sP0f4a7gPcMcB+9wRuX1T2Q6r/zPrxt1TB+FCq98npEfELi7a5o67X6PntVmlmjBmz2D4MljFbRcTqet9nZOZlfT7t48BPAWuBVwJviIgXL9qmXMbAOL7dOgK4IjOvzMx7qRreRxc7n+4xZ8yZxfZh8Jx5M3DqCjucbcvMNjPGjFlsHwbLmP2B1VS/81+q9/tEqkvK+zHzbZk+v6A6Hrg1Mx8NvIPqcjYy80OZeXhmHk7VSfqvmXlRz/NeuvB4Zt447EvRxcS8GVgTPddCZua/y8x96sd663xN7xMj4skR8cWIuCkifkjVw7swzO/hvdtndX3mzS11+AnghXUw3hYRt1FN1nVAzzbX9/z8I6qOhJV6OFVv80Ld7qzrdmDPNtcsflKvzLwN+J/AYcBfLbHdHVS96sfWRS+mCgiozvvJi877pVTX6S6p7qFdmCSrrwm26P81vJX+O2QA7gT2XlS2N30GXf2fxc11UJ9N9fr8xqLN9qIaDq3pY8aYMYsNmjELddqJ6tr1e4ETl9l8q8z8Tmb+IDPvz8x/ppoT4ZhFm017xhzI9u+pTWz/fpt15ow5s9hAOVN/o3sk1R8JA7MtM/PMGDNmsUHbMgsdX3+Xmddl5mbgr9m+U67VnLRl+vmC6mjgjPrnT1Bd1rp4ovIXs+OospHq4iRS/wLcQ/UCfXKZbRfP4v1h4F3AszPz7oj4G7aF1HVUvYkARMSDqIYUNrkG+GBmvnKwqjfWqR8/oAqIhbrtUdet95rGJfdb/+f/u1S97X8LHLXE5h8B3hgRXwZ2o5qdHKrz/r+Z+YwB609mDnJp06CuACIiDsz+rvP8HrBzRByamQvX1D6BwS6N6JXA1g9nVMNNdwG+u8L9Lc9vqEoyY8yYxQbNGOr/sE+l+qbr1zLzviGOvzhjdgYeTXXZaDnD5cyaiNjYc39Line truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x360 with 8 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"fig, axes = pylab.subplots(1, 4, figsize=(16, 5))\n",
|
||||
"for i in range(INFLOW.shape.inflow_loc):\n",
|
||||
@@ -301,57 +299,39 @@
|
||||
"pylab.tight_layout()\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Not surprisingly, the fourth gradient on the left is zero (it's already matching the reference). The other three gradients have detected variations for the intial round inflow posiitons shown as positive and negative regions around the circular shape.\n",
|
||||
"\n",
|
||||
"**TODO, double check after code update**\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {
|
||||
"id": "3mpyowRYUSS4"
|
||||
},
|
||||
"source": [
|
||||
"Not surprisingly, the fourth gradient on the left is zero (it's already matching the reference). The other three gradients have detected variations for the intial round inflow posiitons shown as positive and negative regions around the cirucular shape.\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"## Optimization\n",
|
||||
"\n",
|
||||
"This is just the linearized change that points in the direction of an increasing loss. Now we can proceed by updating the initial velocities in the opposite direction to minimize the loss.\n",
|
||||
"The gradient visualized above is just the linearized change that points in the direction of an increasing loss. Now we can proceed by updating the initial velocities in the opposite direction to minimize the loss, and iterate to find a minimizer.\n",
|
||||
"\n",
|
||||
"Note that this is a fairly tough task: we're producing different dynamics by changing the boundary conditions (the marker inflow position), and an optimizer should now find a single initial velocity state for each of the three simulations on the right, that gives the same state as simulation `0` above at $t=20$. Thus, after 20 update steps the simulation should reproduce a different set of boundary conditions from the velocity state. It would be much easier to simply change the position of the marker inflow to arrive at this goal, but -- to make things a bit more difficult and interesting here -- the inflow is _not_ a degree of freedom. The optimizer can only change the velocity $\\mathbf{u}$ at time $t=0$.\n",
|
||||
"This is a difficult task: the simulation is producing different dynamics due to the differing initial spatial density configuration, and an optimizer should now find a single initial velocity state that gives the same state as the reference simulation at $t=20$. Thus, after non-linear 20 update steps the simulation should reproduce the desired marker density state. It would be much easier to simply change the position of the marker inflow to arrive at this goal, but -- to make things more difficult and interesting here -- the inflow is _not_ a degree of freedom. The optimizer can only change the initial velocity $\\mathbf{u}_0$.\n",
|
||||
"\n",
|
||||
"The following cell re-evalutes the gradient function, and iterates a few times to optimize the initial velocities with a learning rate (step size) `LR`.\n",
|
||||
"The following cell re-evalutes the gradient function, and iterates several times to optimize $\\mathbf{u}_0$ with a learning rate (step size) `LR`.\n",
|
||||
"\n",
|
||||
"`field.functional_gradient` has a parameter `get_output` that determines whether the original results of the function (`simulate` here) are returned, or only the gradient. As it's interesting to track how the loss evolves over the course of the iterations, let's redefine the gradient function with `get_output=True`.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 7,
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"Optimization step 0, loss: 62.215240\n",
|
||||
"Optimization step 1, loss: 62.544254\n",
|
||||
"Optimization step 2, loss: 62.061436\n",
|
||||
"Optimization step 3, loss: 63.683025\n",
|
||||
"Optimization step 4, loss: 64.774864\n",
|
||||
"Optimization step 5, loss: 71.795349\n",
|
||||
"Optimization step 6, loss: 74.215714\n",
|
||||
"Optimization step 7, loss: 60.461643\n",
|
||||
"Optimization step 8, loss: 65.317680\n",
|
||||
"Optimization step 9, loss: 65.126022\n",
|
||||
"Optimization step 10, loss: 65.193939\n",
|
||||
"Optimization step 11, loss: 65.325432\n",
|
||||
"Optimization step 12, loss: 65.471474\n",
|
||||
"Optimization step 13, loss: 65.525093\n",
|
||||
"Optimization step 14, loss: 64.756439\n",
|
||||
"Optimization step 15, loss: 65.134476\n",
|
||||
"Optimization step 16, loss: 64.554466\n",
|
||||
"Optimization step 17, loss: 65.206245\n",
|
||||
"Optimization step 18, loss: 66.495720\n",
|
||||
"Optimization step 19, loss: 66.141006\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"initial_velocity = DOMAIN.staggered_grid(0) * math.ones(inflow_loc=4) # reset\n",
|
||||
"LR = 1e-02\n",
|
||||
@@ -374,29 +354,16 @@
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The loss should have gone down significantly, from above 60 to below 40, and now we can also visualize the initial velocities that were obtained in the optimization.\n",
|
||||
"The loss should have gone down significantly, from above 60 to below 40, and now we can also visualize the initial velocities that were obtained in the optimization. **TODO, double check after code update**\n",
|
||||
"\n",
|
||||
"The following images show the resulting three initial velocities (we're skipping the fourth set with `inflow_loc[0]` because it only contains zeros).\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 8,
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA1QAAAFgCAYAAAC8FTjPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABCnklEQVR4nO3dfbxldV33//ebYQYGh2GkGQfkxuEno0VcinVCvawkxBrNRCsRK4WC0IxuLi1F6dKyvMK6vMnkqo7CBZQKpCEjoQgID9ISORQaN2IjwsUgMAwwDMjdzPD5/bHWwc2Ztc7Zd9+91/ru1/Px2I85+7vXXuu7zjn7Pd/vuvkcR4QAAAAAAL3bZdwdAAAAAIC2YkIFAAAAAH1iQgUAAAAAfWJCBQAAAAB9YkIFAAAAAH3addwdALCwdXZsHnAd10qXRMS6oXQIQHbIGQCpDZozTc0YJlRAC2yWNDPgOiytHEZfAOSJnAGQ2qA509SMYUIFtMUuA16h+8QTw+kHgHyRMwBSGyRnGpox3EMFAAAAAH3iDBXQFhw5BpAaOQMgtQzPUDGhAtrAHnygAwDzIWcApJZpzjChAtoiwwAC0DDkDIDUMswZJlRAG2R6RAdAg5AzAFLLNGfy2yMAAAAAGBHOUAFtkeERHQANQ84ASC3DnGFCBbRFhgEEoGHIGQCpZZgzTKiANsj0mmMADULOAEgt05zJb48AAAAAYEQ4QwW0RYZHdAA0DDkDILUMc4YJFdAGmZ4iB9Ag5AyA1DLNGSZUQFtkGEAAGoacAZBahjmT3x4BAAAAwIhwhgpoiwyP6ABoGHIGQGoZ5gwTKqANMr3mGECDjCBnbN8q6UFJOyRtj4ippBsE0CyZjmeYUAFtkWEAAWiY0eTMz0TE5lFsCEADZTieYUIFtEGmR3QANAg5AyC1THMmvz0CAABNFZK+ZPta2yeNuzMAMAycoQLaIsMjOgAaZvCcWWl7puP5dERMdzz/yYi4w/YzJF1q+1sRcdWgGwXQIhmOZ5hQAW2RYQABaJjBc2bzfIUmIuKO8t9Nti+QdLgkJlTAJMlwPMOECmiDTK85BtAgiXPG9tMk7RIRD5Zf/6yk9yXbIIDmyXQ8w4QKAACMwmpJF9iWivHHpyLii+PtEgAMjgkV0BYZHtEB0DAJcyYibpH0/GQbANAOGY5nmFABbZDpKXIADULOAEgt05xhQgW0RYYBBKBhyBkAqWWYM/ntEQAAAICJZHud7Zttb7B9yjzL/ZLtsF1bmbRbnKEC2iLDIzoAGoacAZBa2mqiiySdLunlkjZKusb2+oi4cc5ye0r6PUlXD2O7TKiANsj0mmMADULOAEgtfc4cLmlDWQRHts+VdLSkG+cs96eSPiDpD4exUSZUQFsw0AGQGjkDILXBcmal7ZmO59MRMd3xfD9Jt3c83yjphZ0rsP1jkg6IiH+2PZQJ1UQkp+0v2D5untf/1vb/7HJdV9o+cXi9e3K9YfvgAddxg+0jhtOjnra72vZVth+0/cEhrXOV7W/ZXjqM9Q3Qj9+x/YFx9qHsSBFAgzyQFDmTVuY58wu2zxtnH8qOkDMNRsaklXnGNGMsIw2eM9LmiJjqeEwvtMmnbt67SPqQpLcPc7dam362b7V9VDfLRsQrIuLs8n3H2/7KnNffEhF/mqKfoxQRPxoRV0qS7T+2/Q8j2vRJkjZLWh4Rw/oFPUXSWRHxiCTZPsb2v9p+2PaVcxe2fZjta8vXr7V9WLcbKn+XHrH9UPn4UsfLH5f0q7afMeD+oIXImZ1Nas7Yfo7tC23fY/s+25fYfm63GyoHsI925MzNs69FxOcl/ajt5w1pv9ASZMzOJjVjOtl+Uzk57XrSy1jmSXdIOqDj+f5l26w9JR0q6Urbt0p6kaT1HrAwRWsnVGiUZ0m6MSKi1zfa3umyU9u7STpOUmeI3ifpI5JOq1h+iaQLy+WfLulsSReW7d36hYhYVj5+drYxIh6V9AVJb+phXWlw5BiTbaw5I2mFpPWSnitptaSvq8idXpzckTNzJ2OfVjGgGy9yBpNr3Bkz+76nS3q3pBt67YfaMJaRUmfMNZLW2j6oHAceqyK7JUkR8UBErIyINRGxRtLXJL06ImaqV9flLg3y5qaYPVJj+3/bvt/2d22/ouP1K22faPtHJP2tpBeXs/ct5etn2f6z8uun276oPAp5f/n1/l304ZnlkYG9O9peYHuz7cXl89+wfVO53ktsP6tmXXvZPqfsw222/6g8RTn7+m+W63nQ9o0urgV98kiX7XUqPoyvL/fzG7ZfZ/vaOdt5m+2dBgS297a90fYvlM+XuSg9udMH0fZZKgLjHeW2jrK9m+2P2P5e+fhIGSyyfUS57nfavkvS/634FrxQ0paI2DjbEBGXRcT5kr5XsfwRKu4H/EhEPBYRH5VkSUdWfX/7cKWknx/SuvrHQGesyJnJzpmI+HpEnBER90XENkkflvRc2z9U9f3tw5UiZyYaGTPZGdPhzyV9VMXZsmG6Uk3IGClpxkTEdkknS7pE0k2Szo+IG2y/z/ark+1SqhWPwQsl3SxppaS/kHSGbXcuEBE3SXqLpH8rZ+8rKtazi4oPxrMkHSjpEUkfW2jjEfE9Sf8m6Zc6mn9F0mciYpvto1UEwy9KWiXpX1Qckazy15L2kvT/SXqpiiMKvy5Jtl8n6Y/LtuWSXi3p3jl9+aKk/yXpvHI/n69idn5QGcSz3ijpnIp9uU/Sb0j6uIvTwx+WdF1EVC17vKRPSvqLcluXSTpVxSnUwyQ9X0XFlT/qeNs+kvZW8T2uOiL731T8LLv1o5K+Oeeo0jfL9m59sgz9L9l+/pzXblKxH+PDvQ1NQc78oC+TljNz/bSkuyLi3gWX/IE/LwemX/XO94jcJGmN7eUD9Gkw5EwTkDE/6MvEZYztwyVNqZgw96PZYxlpGPdQLSgiLo6I50TEsyPi/WXbeyJifcWyRwx6dkrKa0J1W0R8PCJ2qLjka18Vl2X0JCLujYjPRsTDEfGgpPerCIJufErSGySpDMBjyzapCL8/j4ibytnz/5J02NwjOy7q5x8r6V0R8WBE3CrpgyoCQ5JOVPGBvyYKGyLiti726zFJ50n6tXI7PyppjaSLapb/kqR/lHS5pFdKenNX34HCr0p6X0Rsioh7JP1JR/8l6QlJ7y3PJj1S8f4Vkh7sYXvLJD0wp+0BFdfJdtvfNSpC8QpJl9he0fH6gyr+U8ia7TNtb7J9fc3rR9h+wPZ15eM9o+5jA5Az8+9XzjnzpPJI/+mS3tbD296pYmC5n6RpSZ+3/eyO12f7sqKfPiEbZMz8+5VtxpTfs/+j4tLgJ3roZ2d/12jCxzLjktOE6q7ZLyLi4fLLZb2uxPYetv+uPD29VdJVklaUv+gL+ayKU/D7qjh6+YSKozdS8Qv+V7a3uDg9f5+Ky9L2m7OOlZIWS+oMlts6ljtA0nd63a/S2ZJ+pQzIN6o4DfrYPMtPq7hx76wej8I+Uzv3/5kdz++J4nreOver+8mQJD2k4ghXp+XqMsgi4qsR8Uj5H8+fS9oi6ac6FtlTO0/YRi/9keOzJK1bYJl/iYjDysf7Bt6n9iFnFpZrzkgqqnZJ+pKk/xMRdUfmdxIRV5cDy8eiKCzwVRUDvFmzfdnSa5+GijNU40bGLCzXjHmriqttvtbDe57UmrGMlGXGNLdn6Sx0s+HbVdx0/MKIWK4iTKQiMOZfccT9Kv6jfb2KU+TnRjx5Gdrtkt4cESs6Hksj4l/nrGazpG0qQmvWgfpBhZLbJT1bC9tpP8sP6eMqPmC/Iunv695chu60itPob3VvZVC/p53733m98EI/g29Kek4P27tB0vPmXBbxPPV3Q6dU9K9zXT8i6Rt9rms4RnApTkRcpeI/RwyOnMkvZ2ZvFv+SpPWzl5EMoCpnbo2IrQOut39c8tcmZEx+GfMySa+1fVd5X9Z/l/RB2wteqlmjeWMZaSSX/I1Dc3uWzt2S9nd9Bbg9VVxrvMXFTZnv7XH9n1JxTfAv6wenyKXieth3laenZ2/WfN3cN5en+c+X9H7be5an0d+mH1SJ+YSkP7D94y4cPPdUe8d+rnHHDaClc1RcR70tIr6y89ue9G4VH8bfkPSXks7p8siWVFxP/Ucu/v7CSknv0VOr3Czk6yqOpD15xMv2Itu7qyg+sYvt3V3eIKviRssdkn7XxU2kJ5ftXy7fe7yL0pg7sX2g7ZfYXlKu8w9VHFn7asdiL1VRHWe8mjHQebGLG4O/MPu7jErkTGY54+LepkskfTUiTpm7MheXxFYOsGyvsP1z5fp2tf2rKga4X+xYjJxBL8iYzDJG0vEqJj2HlY8ZFZcZnlq+N4+xjJRlxjS3Z+l8WcWZi7tsV1VQ+YikpSqOrnxNT/0PrxvrJa1VcbPyk0cCIuICSR+QdG55+v16Sa+oXoV+R9L3Jd0i6SsqwuzMcj3/qOJa6E+puKTtcypuipzrH8t/77X97x3tf6/i1HdtKNj+cRXB96YyFD+gIpB2GkTU+DMVQfBNSf8p6d/Ltq5ExOMqLj/7tY7mN6r4z+FvVByVekTF31WYXf41KsJ/i4rgfE3ZLhWXFnSGSqc9y3Xer+LI2TpJr5i9LKAMvlequMSg7Vbanul49Fqi+d8lPSuKG4P/WsXvHqqRM5nljKTXSvoJSb/uH/ydl4dsH1i+foCkuUfpZy0u+3aPip/576jIqG93LPMGSX/Xbf8x8ciYzDImIrZExF2zDxVn4bZGxOxleoxlGszRe7l9tJiLv9a9SdKPRcR/jbs/dVzcp/Avkl5Qc7NnL+v6kqTfi6IyUq/v/R1JB0TEOwbpw6Cm9twzZqYG+ptz8pVXXhsR867E9hpJF0XEoQuurzhSNhURwy7tipab0Jz5hKR/jIhL+njvL0h6Y0QcM0gfBjWqnAEGNaEZ0/qxjDR4zjQ1Y3b6Q2TI3m9JuqbJASRJZUWdHx7Sun524aVq3/vXw+jDwGaLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 864x360 with 6 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"fig, axes = pylab.subplots(1, 3, figsize=(12, 5))\n",
|
||||
"for i in range(INFLOW.shape.inflow_loc-1):\n",
|
||||
@@ -408,22 +375,9 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 9,
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA1QAAAFgCAYAAAC8FTjPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA/vUlEQVR4nO3de7xcdX3v//cnyU4CJBgxEQIBgwJavBDtFlF7oRZr8LRGetSCLWK9oP6Kp33IqRfweK2/4q3gaa01KgVbFSnWA6VUbsLBu+wgt4AK0lCBINlKhAi57ORz/lhrw2SzvrNnzaw167u+83o+HvPInu+sWeu79ux5Z33Xd81nzN0FAAAAAChvTtMdAAAAAIC2YkAFAAAAAH1iQAUAAAAAfWJABQAAAAB9YkAFAAAAAH2a13QHAPRmtZlPDriOddKl7r66kg4BSAoZA6Bug+ZMrBnDgApoiUlJEwOuw6SlVfQFQHrIGAB1GzRnYs0YBlRAm8wZ8CrdXbuq6QeANJExAOo2SM5EmjF8hgoAAAAA+sQMFdAmnD0GUCcyBkDdEpyhYkAFtIXZ4Ac7ABBCxgCoW6I5w4AKaJMEQwhARMgYAHVLMGfS2yMAAAAAGBJmqIC2SHSaHEAkyBgAdUs0ZxhQAW2SYAgBiAgZA6BuCeYMAyqgTRIMIQARIWMA1C3BnGFABbTFkKbJzWy1pE9Imivps+5+Ru0bBdC8RC/FARCRRHMmvT0C0Dczmyvpk5KOlXS4pBPM7PBmewUAABAvZqiANqn/rM6Rkm539zskyczOk7RG0i11bxhABBI8cwwgMgnmDAMqoC2qmSZfamYTHffXuvvajvsHSPppx/27JD1v0I0CaIFEL8UBEJFEc4YBFdAmg4fQpLuPV9EVAAlK8EAHQGQSzJn09gjAIO6WdGDH/RV5GwAAAAowQwW0Sf1nda6VdKiZHaxsIHW8pFfXvVEAkUjwzDGAyCSYMwyogLYYwnXH7j5lZqdIulRZ2fSz3X19rRsFEIdEP9sAICKJ5gwDKqBNhhBC7n6JpEtq3xCA+CR4oAMgMgnmTHp7BAAAAABDwgwV0BaJTpMDiAQZA6BuieYMAyqgTRIMIQARIWMA1C3BnGFABbRJgiEEICJkDIC6JZgz6e0RAAAAAAwJM1RAWyR63TGASJAxAOqWaM4woALaJMEQAhARMgZA3RLMGQZUQFskelYHQCTIGAB1SzRn0tsjAAAAABiSkRhQmdl/mNlJXR7/BzP7Xz2u62oze0N1vXtkvW5mhwy4jvVmdnQ1PSq13X3N7Boze9DMPl7ROpeZ2Q/NbI8q1jdAP95qZh9usg+7mTNnsBtqQ87UK/Gc+QMz+3KTfXgEGRMtMqZeiWdMOscykYq3Z7Mwsw1mdkwvy7r7se5+bv6815rZN2c8/mZ3/2Ad/Rwmd3+6u18tSWb2PjP75yFt+mRJk5L2dvdTK1rnOyWd4+4PS5KZvcrMvm1mD5nZ1TMXNrNVZrYuf3ydma3qdUP539LDZrYlv13W8fBnJP2xmT1xwP2pBgc7Q0XOPNao5oyZHWZmF5rZJjP7hZldamZP7XVD+QHs1o6c+dH0Y+7+b5KebmbPqmi/+kfGDBUZ81ijmjGdzOw1+eC050HvyBzLRCrenqFNniTpFnf3sk80s8d8js/MFkg6SVJniP5C0lmSzihYfr6kC/PlHy/pXEkX5u29+gN3X5Tffm+60d23SvoPSa8psa56TF93nGAQAT1oNGckLZF0kaSnStpX0veV5U4Zp3TkzMzB2JeUHdA1h4zBaGs6Y6af93hJp0laX7YfGoVjmUjF27MSps/UmNnHzOx+M/tPMzu24/GrzewNZvZrkv5B0vPz0fvm/PFzzOyv8p8fb2YX52ch789/XtFDH/bPzwzs09H2bDObNLOx/P7rzOzWfL2XmtmTAut6nJl9Pu/DnWb2bjOb0/H4G/P1PGhmt5jZc/L2DWZ2jJmtVvZm/KN8P28ws1ea2boZ23mbmT3mgKDksucoC4y359s6xswWmNlZZnZPfjsrDxaZ2dFmdpeZvcPM7pX0jwW/gudJ2uzud003uPsV7n6+pHsKlj9aWYGVs9x9m7v/b0km6UUFy/bjakn/raJ1oaXImdHOGXf/vrt/zt1/4e47JJ0p6alm9oSi328frhY5M9LImNHOmA5/Lel/K5stq9LVImNqk8SAKvc8ST+StFTSRyR9zsyscwF3v1XSmyV9Jx+9LylYzxxlb4wnSTpI0sOS/m62jbv7PZK+I+m/dzS/WtIF7r7DzNYoC4Y/lLRM0jeUnZEs8reSHifpyZJ+W9kZhT+VsoCQ9L68bW9JL5P08xl9+Zqk/1/Sl/P9PELZmdWD8yCedqKkzxdsv+dl3f21kr4g6SP5tq6QdLqkoyStknSEpCMlvbvjaftJ2kfZ77jojOwzlb2WvXq6pBtnnFW6MW/v1Rfy0L/MzI6Y8dityvajeZw9bho582hfRi1nZvotSfe6+89nXfJRf50fmH7LHvsZkVslrTSzvQfo0+DImKaRMY/2ZeQyxsyOlDSubMDcj/SPZSIVb8/Ku9PdP+PuO5Vd8rVc2WUZpbj7z939K+7+kLs/KOlDyoKgF1+UdIIk5QF4fN4mZeH31+5+q7tPKQuJVTPP7JjZ3Px573L3B919g6SPKwsBSXqDsjf8tZ653d3v7GG/tkn6sqQ/ybfzdEkrJV08yLIBfyzpA+5+n7tvkvT+jv5L0i5J781nkx4ueP4SSQ/2uC1JWiTplzPafilpcYn+rlQWildJutTMlnQ8/qCy/xSaxeU4MSBnuu9XyjnziPxM/yclva3E096h7MDyAElrJf2bmT2l4/Hpvizpp0+VIGNiQMZ0369kMyb/nf29skuDd/X6vBn9XanUj2UiFW/Pyrt3+gd3fyj/cVHZlZjZnmb26Xx6+gFJ10hakv+hz+Yryqbglys7e7lL2dkbKfsD/4SZbc6n53+h7LK0A2asY6mkMUmdwXJnx3IHSvpJ2f3KnSvp1XlAnijp/DxwBl12pv312P7v33F/U349b8j96n0wJElblJ3h6rS3egwyd/+Wuz+c/8fz15I2S/rNjkUW67EDtmZwsNM0cmZ2qeaMpKxql6TLJP29u4fOzD+Gu38vP7Dc5llhgW9JemnHItN92Vy2T5UiY5pGxswu1Yz5/5RdbfPdEs95xMgcy0Qq3p7VZ7YPG56q7EPHz3P3vZWFiZQFRvcVu9+v7D/aP1I2RX5ex2VoP5X0Jndf0nHbw92/PWM1k5J2KAutaQdJurtjPU/R7B6zn/mbdLuyN9irJf1Tl33pedkC9+ix/e+8Xni21+BGSYeV2N56Sc+acVnEs9TfBzqlrH+d6/o1STf0uS6MJnImvZyZ/rD4ZZIucvcPlXlugaKc2eDuDwy4XowGMia9jPldSceZ2b3557JeIOnjZjbrpZoBHMsM0SgOqH4maYWFK8AtVnat8WbLPpT53pLr/6Kya4JfoUenyKXseth35VPO0x/WfOXMJ+fT/OdL+pCZLc6n0d+mR6vEfFbS/zSzX7fMITOn2jv2c6V1fAA093ll11HvcPdvPvZpfS/b6UuS3m3Z9y8slfQe7V7lZjbfV3Ym7ZEzXmY218wWKis+McfMFlr+AVllH7TcKel/WPYh0lPy9q/nz32tmW0o2pCZHWRmLzSz+fk6/1LZmbVvdSz228qq4zSPs8dtQc4kljOWfbbpUknfcvd3zlyZZR9SLzzAMrMlZvaSfH3zzOyPlR3gfq1jsThyhoxpCzImsYyR9Fplg55V+W1C2WWGp+fP5Vgm4oyJt2f1+bqymYt7zayogspZkvZQdnblu9r9P7xeXCTpUGUfVn7kTIC7f1XShyWdl0+/3yzp2OJV6K2SfiXpDknfVBZmZ+fr+Rdl10J/Udklbf9H2YciZ/qX/N+fm9l1He3/JOkZ6i0Uyizb6a+UBcGNkm6SdF3e1hN33y7pHOXXPedOVPafw6eUnWl6WNn3Kkwv/3Jl4b9Z0uskvTxvl7JLCzpDpdPifJ33KztztlrSsZ5/0DwPvpcqu2ygWXy+oU3ImcRyRtJxkp4r6U/t0e952WJmB+WPHyhp5ln6aWN53zYpe83fqiyjftyxzAmSPt1r/2tBxrQJGZNYxrj7Zne/d/qmbGbtAXefvkyPY5mIM8a8fLl9tJhl39Z9n6TnuPttVS1bNcs+p/ANSc8OfNizzLouk/TnnlVGKvvct0o60N3fPkgfqjC+aJFPHDFYgR779rfXuft4RV0CCo1oznxW0r+4+6V9PPcPJJ3o7q8apA+DGlbGWFYO+xOS5kr6rLsHv5MHKDKiGcOxjOI9jnnMF5EheW+RdG2PoVJm2UrlFXWeVtG6fm/2pYLP/dsq+gCMmFHMmTcM8Nx/k/RvVfQjdpYVRfikpBdLukvStWZ2kbvf0mzP0DKjmDEcy0SMAdUIya+9NWWXx1W2LIZkepociBg502LDyZgjJd3u7ndkm7TzJK2RxIAKPSFjWm4IOTPbLLiZvU1Z6f4pZZdiv857KNvfDQOqEeLuK+tYFkPEgAqRI2dabvCMWWpmEx3317r72o77Byir8DbtLmVfZgv0hIxJQI3HMj3Ogv9A0ri7P2Rmb1H2Jdp/NMh2GVABbcKACkCdBs+YyRg/3wAgIvUey8w6C+7uV3Us/13tXjikLz0PqPIR34Sku939983sYEnnSXqCpHXKPlC7vds6ls6Z4yvn9vKdckD61k1NTbr7sqb7EYtKMmZszFcuXFh/Z4EWWLdlS4wZc7eyamXTVujR7yaqVSUZs3Spr1y5sva+Am2wbt26GDOm7Cz461VBOfkyM1R/LulWSXvn9z8s6Ux3P8/M/iHv0Ke6rWDl3LmaWLKkn34CybHJyXLX66b/GarBM2bhQk085zn19hJoCbvmmhgz5lpJh+aDmbslHa/sC1eHYfCMWblS11470W0RYGTMmWPlP3c0eM7Mdllxia7Line truncated
|
||||
"text/plain": [
|
||||
"<Figure size 864x360 with 6 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"fig, axes = pylab.subplots(1, 3, figsize=(12, 5))\n",
|
||||
"for i in range(INFLOW.shape.inflow_loc-1):\n",
|
||||
@@ -439,29 +393,16 @@
|
||||
"source": [
|
||||
"## Re-simulation\n",
|
||||
"\n",
|
||||
"Most importantly, we can also visualize how the full simulation over the course of 20 steps turns out given the new initial conditions for each of the inflow locations. This is what happened internally at optimization time, and was measured by our loss function, but it's good to get an intuition for which solution the optimization has found.\n",
|
||||
"We can also visualize how the full simulation over the course of 20 steps turns out, given the new initial velocity conditions for each of the inflow locations. This is what happened internally at optimization time, and what was measured by our loss function. Hence, it's good to get an intuition for which solutions the optimization has found.\n",
|
||||
"\n",
|
||||
"Below, we re-run the forward simulation with the new initial conditions from `initial_velocity`:\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 15,
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA54AAAEiCAYAAACRGA3AAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA3DklEQVR4nO3deZxcdZnv8e/TS/alsxFCEnaUTQkIUYRBwAUFF7yDKI5enMvIOCOOjssVnRm3cUHHBXdEhaCDuOIVFUFEGAZEIECAhARICJCFbCTdnbWT7nruH3UyNs1zuut0n1N9uvvzfr361d3fqjrnd6rqqVO/qlNPmbsLAAAAAICiNAz2AAAAAAAAwxsTTwAAAABAoZh4AgAAAAAKxcQTAAAAAFAoJp4AAAAAgEIx8QQAAAAAFIqJZ47M7EAzczNr6ufl/8bMfp/3uJJl/4OZrTezbWY2Ladlfs7M3pfHsgYwhtFmtszMZgzmODA0UKODw8zuNrOjBnscKD9qtP7Yj5YfdTE46rnvsqorzWyLmd1dj3UOCnfnp8ePpJMl/UlSm6TNku6QdEINlztQkktqyvO8OWxPs6Sdko7JcZkzJK2RNDb5f5Skn0t6ItmuU3uc3yR9XtIzyc/nJVmN63qHpC5J27r9nNrt9P8r6UuDfb/hp34/1GhNy8xaox+StFjSVkkrJX0ow7pOlVTpUaPndzv9XEm/GOz7DT/1+6FGa1pmphrtdrlRkpZKWp1hXexHS/BDXdS0zMx1Iek4Sbcl9+v1kt5b47py3XdJWiBpd7KczZJuknR4jZf9K0mrJY0f7PtpkT+849mDmU2S9BtJX5c0VdJsSZ+U1DGY4xqgmZLGSFqS9YLJKzDR/eQdkq53953dstslvU3SuuD8F0o6W9Ixkl4o6XWS/j7DUO509wndfm7tdtqPJJ1vZqMzLA9DFDX6bDnWqEn635KmSHq1pIvM7C0ZhrK2R41e1e206ySdZmb7Zlgehihq9NlyrNG9PiRpY9ZxiP3ooKIuni2vujCz6ZJukPQdSdMkHSopy7u7ee+7vuDuE1S9fddI+n6NlztA0hPuvj3DuiRJ/X0nfDAw8Xyu50mSu1/j7l3uvtPdf+/uD0qSmTWY2b+a2ZNmtsHMfmBmk6MFmdkTZvaKbv9/wsz+M/n3tuR3a3J4wolm9g4zu73b+V9qZveYWVvy+6XdTrvVzP7dzO4ws61m9vuk+HqO4XmSHum2rj/WuOzPmNkdknZIOjjYvNdI+q+9/7j7bne/1N1vV/VV1Z7OV/XV1NXuvkbSl1R9cBkwd18taYukl+SxPJQeNVpAjbr7F9z9PnfvdPdHJP1K0knR9ZaVu++SdK+kM/JYHkqPGi1mPyozO0jVJ+Cfi07vL/ajdUFdFFMX75d0o7tf7e4d7r7V3ZdG11tWA9l3JRPnn0qatzczs/3M7BdmttHMVprZPyX5BZK+J+nE5Db7ZJK/1swWmVmrmf3JzF7YbVlPmNmHzexBSdvNrMnMXpKcr9XMHjCzU7udv9fb1cxO7nbZVWb2jiQfbWZfNLOnrHpI9WVmNjbr9bEXE8/nelRSl5ldZWavMbMpPU5/R/JzmqoFM0HSN/qxnlOS3y3JKyx3dj/RzKZK+q2kr6n6Cs6XJf3Wnn38/Fsl/a2kfVQ9FOGDPVfi7o9K2nt8eou7n17jst+u6ruUEyU9GYz/BfrLA04tjpL0QLf/H+g2rloca2abzOxRM/s3e+6rO0tVfTcVwx81WpV3jXbfNlP1sJ8sr2Lvk+yUVprZV8xsfI/TqdGRgxqtKqJGvy7po6oe3pgV+9HBRV1U5V0XL5G0OZk0bTCzX5vZ/hkuX8i+K1nOeZKWJ/83SPq1qs9/Z0t6uaT3mdkZ7v59Se/SX45K+LiZHSvpClWPDpym6ju619mzj0o4T9JZklpUfff5t5I+reo76h+U9At79me3w9vVzA6Q9DtVH19mqDpZXpRc5hJVXzSZp+q7ybMlfSzr9bEXE88e3L1d1WPwXdJ3JW00s+vMbGZylr+R9GV3f9zdt0n6iKS3BA/gA3WWpMfc/YfJOxDXSFqm6iGqe13p7o9Gr6rksOwF7r4kOX1PsIwWVT8LVqsJqn6mYa82SROSJ7h9uU3S0aoWyl+rWmgf6nGercmYMMxRo/8j7xrt7hOq7h+urPH8y1TdtlmSTpf0IlWfcHRHjY4Q1Oj/yLVGzeyNkhrd/Ze1XqYb9qODjLr4H3nvu+aoelTdeyXtr2qPgmtqvGwR+64PmllrcrmTVZ1oS9IJkma4+6eSd3EfV/V+kPaRlgslfcfd70reIb9K1cOyux+V8DV3X5XcTm9T9RDl69294u43SVoo6cxu50+7Xd8q6Q/Ju/F73P0Zd1+UPEe/UNI/u/tmd98q6bO9jLlPTDwD7r7U3d/h7nNUfaDeT9Klycn76dmv0DwpqUnVVxry1HM9e9c1u9v/3Y9136Hq5C6vZa/qYxlbVH21qlbbJE3q9v8kSdvcq5+o7k3yILwyKaSHJH1K0jk9zjZRUmuG8WAIo0Yl5V+jkiQzu0jVz3qe5e41ffbI3de5+8NJja5UtVHJX/c4GzU6glCjknKs0eTdky9I+qcax/cs7EfLgbqQlP++a6ekX7r7PV49NPaTkl5qKYcpd1fQvuuL7t6iapOnnZKen+QHSNovOZS1NZmcflTpt+8Bkj7Q4/xzVb2O91rV4/xv6nH+k1WdVO+VdrvOlbQiGMMMSeMk3dttmTckeb8w8eyDuy9TtUvV0Um0VtUbd6/9JXWq2kWrp+2q3mB7df9wcl8Trp7r2buuNX1crha1LLuv8T2o5PMKNVqiZx+qcIz68WH0hKvaCKW7I/TsQ3kxQlCjqbLWqMzs/0i6WNLLvfqZr/5yPXf/Qo2OUNRoqiw1epiqT2T/28zWSbpW0iwzW2dmB9a4jO7Yjw4y6iJV1n3Xgz2W2ecbGr3Ibd/l7k+p+i7sV5PPRK6StNLdW7r9THT3M1MWsUrSZ3qcf1zyLnL38XY//w97nH+8u19Sw3BXSTokyDepOnk+qtsyJ3u1eVK/MPHswcwON7MPmNmc5P+5qh6S8ufkLNdI+mczO8jMJqj6lvNP3L0zWNwiVQ+TaDaz4/XsVxc3qtrCOfpgtSRdL+l5ZvZWq35g+M2SjlS1I9pA5bHs6yW9rHuQfAB5TPLvKDMb0+1Q2h9Ier+ZzTaz/SR9QNUH3L2XfcKSDzL3lHwWYmby9+GS/k3Vxid7T5+t6vHsf44uj+GFGs20jJpr1Mz+RtXr6pXJIUDqcdlbzewT0YrM7DQzO8Cq5qr6mZDuNTpG1UOYbsowfgxR1GimZdRao4tVfVdiXvLzd6pOSOYpedeD/Wi5UReZlpHl+eWVkt5oZvPMrFnV+/bt7t6WXDbXfZdVv0/11Fo2JDncda2qh6veLWmrVRsCjTWzRjM72sxOSLn4dyW9y8xenIxvvJmdZWZp7wb/p6TXmdkZybLHmNmpe+9vfbha0ivM7NzkdptmZvPcvZKM4ytmtk+y/bPNrN+NApl4PtdWSS+WdJeZbVf1AWGxqhMlqfpB3x+q+nmJlZJ2SXpPyrL+TdVXELao+tb/j/ae4O47JH1G0h1Wffv6WZ3k3P0ZSa9N1vuMqm//v9bdNw10A3Na9g8knWnP7mz1iKqvjMyWdGPy995Xvr6j6oeqH1L1+vxtksnMRqn6wem0Hd7LJT2Y3B7Xq/pK72e7nf5WSVfVelgghjxqtDZZa/TTqtbhPVbtqrfNzC7rdtm5qn7nXORYVb+bbnvy+yE9+5DA10m61d3XZhg/hi5qtDY112jyebh1e39U/Y7ASvJ/F/vRIYG6qE2mfZe7/1HVQ1Z/K2mDqg1w3trtsrntu5LJ6dbkfLX6D1WvhyZVr5t5qt6+m1TtZBseEuzuCyW9U9UGU1tUbVL0jrSVuPsqSW9Q9brYqOoLUh9SDXO95N3ZM1W93Tar+sLG3qMUP5ys+89m1i7pD/rL4cOZWQ0fsQNCZvZZSRvc/dIBLudkSe929/P6cdnRqh4CcYq7bxjIOIDhJscanSPpp+7+0j7PHF/+LkkXuPvigYwDGG7YjwLPVdZ9l5m9TdXDTj8ykHGNZEw8AQAAAACF4lBbAAAAAEChmHgCAAAAAArFxBMAAAAAUCgmngAAAACAQjXVc2VmRicjQJK79/zi7lKgRoEqahQoN2oUKLeoRus68axqrP8qgVLpGuwB9IEaxUhHjQLlRo0C5RbXKIfaAgAAAAAKxcQTAAAAAFAoJp4AAAAAgEIx8QQAAAAAFIqJJwAAAACgUEw8AQAAAACFYuIJAAAAACgUE08AAAAAQKGYeAIAAAAACsXEEwAAAABQKCaeAAAAAIBCMfEEAAAAABSKiScAAAAAoFBMPAEAAAAAhWLiCQAAAAAoFBNPAAAAAEChmHgCAAAAAArFxBMAAAAAUCgmngAAAACAQjHxBAAAAAAUioknAAAAAKBQfU48zWyMmd1tZg+Y2RIz+2SSLzCzlWa2KPmZV/hoATwHNQqUGzUKlBs1CtRHUw3n6ZB0urtvM7NmSbeb2e+S0z7k7j8vbngAakCNAuVGjQLlRo0CddDnxNPdXdK25N/m5MeLHBSA2lGjQLlRo0C5UaNAfdT0GU8zazSzRZI2SLrJ3e9KTvqMmT1oZl8xs9FFDRJA76hRoNyoUaDcqFGgeDVNPN29y93nSZojab6ZHS3pI5IOl3SCpKmSPhxd1swuNLOFZrYwnyED6IkaBcqNGgXKjRoFimfVowsyXMDsY5J2uPsXu2WnSvqgu7+2j8u61NiPYQLDSZfc3YpaOjUKDBQ1CpQbNQqUW1yjtXS1nWFmLcnfYyW9UtIyM5uVZCbpbEmL8xwugNpQo0C5UaNAuVGjQH3U0tV2lqSrzKxR1YnqT939N2b2RzObIckkLZL0ruKGCaAX1ChQbtQoUG7UKFAHmQ+1HdDKOPwAUNGHCA0ENQpI1ChLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x360 with 4 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"smoke = initial_smoke \n",
|
||||
"velocity = initial_velocity\n",
|
||||
@@ -479,7 +420,7 @@
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Naturally, the image on the right is the same (this is simulation (3), our reference), and the other three simulations now exhibit a noticeable shift towards the right. This is especially visible for the far left simulation (0), which was strongly curved in its original state, and now has a shape that tries to match the straighter reference on the right. \n",
|
||||
"Naturally, the image on the right is the same (this is the reference), and the other three simulations now exhibit a noticeable shift towards the right. This is especially visible for the far left simulation (0), which was strongly curved in its original state, and now has a shape that tries to match the straighter reference on the right. **TODO, double check after code update**\n",
|
||||
"\n",
|
||||
"Note that all three simulations need to \"work\" with a fixed inflow, hence they cannot simply \"produce\" marker density out of the blue to match the target. Also each simulation needs to take into account how the non-linear model equations change the state of the system over the course of 20 time steps. So the optimization goal is quite difficult, and it is not possible to exactly satisfy the constraints to match the reference simulation in this scenario. This is exactly what makes it an interesting test case, though.\n",
|
||||
"\n",
|
||||
@@ -488,22 +429,9 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 18,
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA3EAAAE/CAYAAADouUp5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAvt0lEQVR4nO3de3ycd3Xn8e+RNJIsW5Z8t2M7ca4kISUOmCQUwiU0NFBoQi9AtmXDQhvahe6yyxZS+mqhu4VSdhPobrvQUNKEFgiUQgksbUlDQkhbkjiJc8O54djx3fJFsmxZ97N/zLjI5jmPpNFoZn6az/v10svS78wzz08jz9FzNNI55u4CAAAAAKShqdYbAAAAAABMHUUcAAAAACSEIg4AAAAAEkIRBwAAAAAJoYgDAAAAgIRQxAEAAABAQiji5gAz+4yZ/d6Ej3/TzPaa2REzW2JmLzezZ0ofXx3cxx+Z2fuqtedgD21m9qSZLavlPgDUp7mS60r7uN/MXljrfQCoT3Ml33FtN3uMOXH1zcy2SlohaVTSmKQfSvq8pJvcfTzj9gVJhyVd6u6PlNbulHS7u/9JcI5lkjZJOsvdj5lZq6QvStog6TRJr3H3uyfc3iR9XNKvlZb+QtL1PoX/TGb2Dkmfk3RswvIbj9+/mX1A0gp3f/9k9wVg7qjTXPfbkq4txfZL+r/u/j+n+Pm8WtJ3JQ1MWH6Pu99air9F0lvd/Rencn8A5o56zHcTjmuV9IikTndfM8XP5x3i2q7qeCUuDW9y904Vn3Qfl/RBFZ8sWVZIapf0xIS10076+GTvkPRtd5/45LtX0q9K2pNx++skXS3pQkkvkvQmSe+e7JOY4F/dfcGEt7snxL4o6Voza5vG/QGYG+ot15mkfy9pkaQrJb3XzN42+afxb3adlOtunRC7XdJrzGzlNO4PwNxRb/nuuN+W1JO782xc21UZRVxC3L3P3W+X9FYVnwwXSJKZ3WJmf2hm50h6qnTzXjP7rpn9SNIZkr5Zesk96wn0eknfm3CeYXf/lLvfq+JPiE52raQb3H2Hu++UdIOKyaISn+MOSYckXVqJ+wOQnnrJde7+CXd/yN1H3f0pSd+Q9PIKfY6Dkh6U9LOVuD8AaaqXfFc65+kqFnl/VLFPUFzbzRaKuAS5+/2Sdki67KT1pyUd/xuLbne/3N3PlPS8ij/xWeDuQxl3+VP6cYKYiheq+FL7cY9MOO9UXGRm+83saTP7PTNrOSm+WcVX+QA0sDrIdf+m9Gvklyn/J98nW176G5bnzOyTZjb/pDi5DoCkusl3/0fSh3Tir0VOFdd2VUYRl65dkhZX6L66JfVP4/YLJPVN+LhP0oLSRc5k7pF0gaTlkn5R0jUqvnQ/UX9pTwBQy1w30UdU/J75l1O8/ZOS1ktaJelySS+RdONJtyHXAZioZvnOzN4sqdndv17Gubi2qwGKuHStlnSwQvd1SFLnNG5/RNLCCR8vlHRkKo1N3H2Luz/n7uPu/pik/y7pl066Waek3mnsB8DcVctcJ0kys/eq+LdxPxf8xPsnuPsed/9hKdc9J+kDKl7cTESuAzBRTfJd6bcEPiHpP5VzIq7taoMiLkFm9lIVn+j3VuguH5V0zjRu/4ROfEn8Qk3vV4wmchWbB0x0nk78dU0ADagOcp3M7J2Srpf02tLfdZTL9ZPfc8l1ACTVPN+dLWmdpO+b2R5JX5O0ysz2mNm6Ms7NtV0VUMQlxMwWmtkbJd0m6a9LP+2ohG9LetVJ52ozs/bSh61m1j7h1yU/L+m/mtlqMztF0vsl3TLh2K2ldrNZn8PrzWxF6f1zJf2eis0CjsdXq/irBD+oxCcGID31kuvM7FckfUzSFe6+JWOfd5vZR4LP4TVmdpoVrVWx+9zEXNeu4q9Y3lGBzwtAouok3z0uaa2KvwK+XsURUntL728vHcu1XZ05+Y8OUZ++aWajksZVnCVyo6TPVPD+Py9pk5nNm9CK9ikV29dK0j+W/j1d0lZJf65iV6TjieYvSmvH54ssUfxEfa2kW8xsgYoJ4q9VvEg67t9JunWqv7IEYE6pt1z3hyrmswcm/MnvX7v7b5TeXyvpn4NzXaRiflsk6YCkr0v63QnxN0m62913zeQTApCsusl37r5VE8YOmNlBSePuvqf0Mdd2dYhh35AkmdnHJO1z90/N8H5eoeJA22vKOLZNxZfaX+nu+2ayDwDIUsFct0bSV9z9p8s8/j5J73L3x2eyDwCIcG03t1HEAQAAAEBC+Js4AAAAAEgIRRwAAAAAJIQiDgAAAAASQhEHAAAAAAmp6oiBrpYOX97WVc1TAg1l31Cf+kYHTh6wiSrrKnT4CnIdMGv2DvWpb4RcV2vdhXm+qp1cB8ymJ4/s3e/uy05er2oRt7ytS5984TureUqgofyXJ26u9RYgaUVbl/70gnfUehvAnPXex2+p9RYgaVV7l2656FdqvQ1gTrv0+zduy1rn1ykBAAAAICEUcQAAAACQEIo4AAAAAEgIRRwAAAAAJIQiDgAAAAASQhEHAAAAAAmhiAMAAACAhFDEAQAAAEBCKOIAAAAAICEUcQAAAACQEIo4AAAAAEgIRRwAAAAAJIQiDgAAAAASQhEHAAAAAAmhiAMAAACAhFDEAQAAAEBCKOIAAAAAICEUcQAAAACQEIo4AAAAAEgIRRwAAAAAJIQiDgAAAAASQhEHAAAAAAmhiAMAAACAhFDEAQAAAEBCKOIAAAAAICGTFnFm1m5m95vZI2b2hJn9QWn9FjN7zsw2ld7Wz/puAWCWkOsANAJyHTA3tEzhNkOSLnf3I2ZWkHSvmf19Kfbb7v7V2dseAFQNuQ5AIyDXAXPApEWcu7ukI6UPC6U3n81NAUC1kesANAJyHTA3TOlv4sys2cw2Sdon6Q53v68U+qiZPWpmnzSzttnaJABUA7kOQCMg1wHpm1IR5+5j7r5e0hpJF5vZBZJ+R9K5kl4qabGkD2Yda2bXmdlGM9vYNzpQmV0DwCyoWK4bIdcBqF+VynW95DqgZqbVndLdeyXdJelKd9/tRUOS/lLSxcExN7n7Bnff0NXSMeMNA8Bsm3GuK5DrANS/mea6bnIdUDNT6U65zMy6S+/Pk3SFpCfNbFVpzSRdLenx2dsmAMwuch2ARkCuA+aGqXSnXCXpVjNrVrHo+4q7f8vMvmtmyySZpE2SfmP2tgkAs45cB6ARkOuAOWAq3SkflXRRxvrls7IjAKgBch2ARkCuA+aGaf1NHAAAAACgtijiAAAAACAhFHEAAAAAkBCKOAAAAABICEUcAAAAACSEIg4AAAAAEkIRBwAAAAAJoYgDAAAAgIRQxAEAAABAQijiAAAAACAhFHEAAAAAkBCKOAAAAABICEUcAAAAACSEIg4AAAAAEkIRBwAAAAAJoYgDAAAAgIRQxAEAAABAQijiAAAAACAhLbXeAACgcXmZx1lFdwEAtTPm8Wsqg2PNmeuWkz3bmsfDWLPFMaSFV+IAAAAAICEUcQAAAACQEIo4AAAAAEgIRRwAAAAAJIQiDgAAAAASQhEHAAAAAAmZdMSAmbVLukdSW+n2X3X3D5vZ6ZJuk7RE0oOS3u7uw7O52Xo14nGz6/6R7NawkjSW01u7kFNed7aMhTGz7DulWgfyketmJm9UwLaBtjD2TH+cnTpyvkOd2pHdJnvNvKHwmPbmOHcOBW28gbmGXDcz4znXfHl2HpsXxr7fUwhjh4ay89Zly+N9nNM5EMYWdcRf0pacHHlwIHv/eY9HU3BNisqYyrX9kKTL3f1CSeslXWlml0r6Y0mfdPezJB2S9K5Z2yUAzD5yHYBGQK4D5oBJizgvOlL6sFB6c0mXS/pqaf1WSVfPxgYBoBrIdQAaAbkOmBum9Ft2ZtZsZpsk7ZN0h6QfSep199HSTXZIWj0rOwSAKiHXAWgE5DogfVMq4tx9zN3XS1oj6WJJ5071BGZ2nZltNLONfaPx7+gCQK1VLNeNkOsA1K9K5bpech1QM9Pqd+HuvZLukvQySd1mdvzPztdI2hkcc5O7b3D3DV0tHTPZKwBUxYxzXYFcB6D+zTTXdZPrgJqZtIgzs2Vm1l16f56kKyRtVvFJ/0ulm10r6RuztEcAmHXkOgCNgFwHzA2TjhiQtErSrWbWrGLR9xV3/5aZ/VDSbWb2h5IelvS5Wdxn1WQ3rS46MJT9cN29N26v+uDwc2FsRHEr7EXjy8LYusKiMLZ+SfYeT+sYzVyXpPk5IwsYTYAG0lC5Lk/UFHrXsdbwmL/fGbeSfmz86TC2dHxlGLtgQVcYW96WnXcHc0YFrOw4FsbmFUbCWN9gexjrH8luDT46Xl4bcqAKyHWTyGubvy8nH3xzZ3xZ/ffH7gljy3VWGPsPp2T/aeKZC+J8lrf/4ZwcubBjMIwtXX4kc/1gz4LwmGgsgSSNM31gxiYt4tz9UUkXZaxvUfH3qAEgeeQ6AI2AXAfMDbzQAgAAAAAJoYgDAAAAgIRQxAEAAABAQijiAAAAACAhFHEAAAAAkJCpjBiYc0ZyWq9uH8huFy1Jm/uy19vjbq36mc4zw1jfcDzQYPOxA2Hs7pG7wthdu7Pv80zbEB5z9aolYezczngMAj8BAOamh3uz20J/cW88MmWdnRLG3rni3DC2oi0ecdJsebHs9YGxODM91dcZxs7qzG6fLUkvumRfGNv/VHa78c374rya1/4bQPVEz8XDwegQSfrUM3Eb/u8d/kwY27Dw2jD2gbMWh7Fzug6GsUjfYFsYe+ZwnAcP5IxPWDPYn7l+xisH4o3cG4f2H43HD2BquA4HAAAAgIRQxAEAAABAQijiAAAAACAhFHEAAAAAkBCKOAAAAABICEUcAAAAACQk6REDzeZh7JzFvWFs0aK4HeoT25aHscuWj2Suj47HtfChodYwdnA4bmH7c03dYWzXsZ8JY3+++9nM9Qf6vxQes3NPPH7g18YuDWMXdA2HsULO1wZAdeQ9C584HLeLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x360 with 3 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"fig, axes = pylab.subplots(1, 3, figsize=(16, 5))\n",
|
||||
"for i in range(3):\n",
|
||||
@@ -511,6 +439,13 @@
|
||||
" axes[i].set_title(f\"Diff {math.to_int(INFLOW_LOCATION).inflow_loc[i].vector.unstack_spatial('x,y')}\")\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Conclude... **TODO, double check after code update**"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
@@ -521,7 +456,7 @@
|
||||
"This example illustrated how the differentiable physics approach can easily be extended towards significanlty more \n",
|
||||
"complex PDEs. Above, we've optimized for a mini-batch of 20 steps of a full Navier-Stokes solver.\n",
|
||||
"\n",
|
||||
"This is a powerful basis to bring NNs into the picture! As you might have noticed, our degrees of freedom were still a regular grid, and we've jointly solved a single inverse problem. There were three cases to solve, of course, but nonetheless no question of generalization and no neural network structure. Thus, in line with the PINN example of {doc}`physicalloss-code` we've not really dealt with a _machine learning_ task here. However, that will be the topic of the next chapters.\n"
|
||||
"This is a powerful basis to bring NNs into the picture! As you might have noticed, our degrees of freedom were still a regular grid, and we've jointly solved a single inverse problem. There were three cases to solve, of course, but nonetheless the setup remains a direct optimization. Thus, in line with the PINN example of {doc}`physicalloss-code` we've not really dealt with a _machine learning_ task here. However, that will be the topic of the next chapters.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -532,7 +467,7 @@
|
||||
"\n",
|
||||
"Based on the code example above, we can recommend experimenting with the following:\n",
|
||||
"\n",
|
||||
"- Modify the setup of the simulation to differ more strongly across the four instances, run longer, or use a finer spatial discretization (i.e. larger grid size). This will make the optimization problem tougher, and might not converge directly with this simple setup.\n",
|
||||
"- Modify the setup of the simulation to differ more strongly across the four instances, run longer, or use a finer spatial discretization (i.e. larger grid size). Note that this will make the optimization problem tougher, and hence might not converge directly with this simple setup.\n",
|
||||
"\n",
|
||||
"- For more complex cases, add a multi-resolution optimization, that first solves with a coarse discretization, and then uses this solution as an initial guess for a finer discretization. (Ideally, implement it with a loop that continuously refines the discretization.)\n"
|
||||
]
|
||||
|
||||
Reference in new issue
Block a user