diff --git a/diffphys-code-ns.ipynb b/diffphys-code-ns.ipynb index afe15d5..a53d052 100644 --- a/diffphys-code-ns.ipynb +++ b/diffphys-code-ns.ipynb @@ -8,6 +8,11 @@ "source": [ "# Differentiable Fluid Simulations\n", "\n", + "\n", + "**TODO , NEXT update diffphys example, then check burgers examples**\n", + "\n", + "\n", + "\n", "Next, we'll target a more complex example with the Navier-Stokes equations as model. 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:\n", "\n", "$\\begin{aligned}\n", diff --git a/overview-ns-forw-v2.ipynb b/overview-ns-forw-v2.ipynb index 7ce9dcc..9e45afc 100644 --- a/overview-ns-forw-v2.ipynb +++ b/overview-ns-forw-v2.ipynb @@ -8,8 +8,6 @@ "source": [ "# Navier-Stokes Forward Simulation\n", "\n", - "**TODO , double check here - sync with latest local-colab-notebook version, then update diffphys example**\n", - "\n", "Now let's target a somewhat more complex example: a fluid simulations based on the Navier-Stokes equations. This is still very simple with ΦFlow, as differentiable operators for all steps exist. The Navier-Stokes equations (in their incompressible form) introduce an additional pressure field $p$, and a constraint for conservation of mass:\n", "\n", "$\\begin{align}\n", @@ -69,14 +67,14 @@ }, { "cell_type": "code", - "execution_count": 10, + "execution_count": 2, "metadata": { "id": "WrA3IXDxv31P" }, "outputs": [], "source": [ "DOMAIN = Domain(x=40, y=32, boundaries=CLOSED, bounds=Box[0:100, 0:80])\n", - "INFLOW = domain.grid(Sphere(center=(30, 15), radius=10)) * 0.2\n", + "INFLOW = DOMAIN.grid(Sphere(center=(30, 15), radius=10)) * 0.2\n", "\n", "DT = 1.5\n", "NU = 0.01" @@ -98,7 +96,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -118,7 +116,7 @@ }, { "cell_type": "code", - "execution_count": 12, + "execution_count": 4, "metadata": { "colab": { "base_uri": "https://localhost:8080/", @@ -138,10 +136,10 @@ { "data": { "text/plain": [ - "" + "" ] }, - "execution_count": 12, + "execution_count": 4, "metadata": {}, "output_type": "execute_result" }, @@ -159,15 +157,15 @@ } ], "source": [ - "def step(velocity, smoke, pressure, buoyancy_factor=1.0):\n", - " smoke = advect.semi_lagrangian(smoke, velocity, dt) + inflow\n", + "def step(velocity, smoke, pressure, dt=1.0, buoyancy_factor=1.0):\n", + " smoke = advect.semi_lagrangian(smoke, velocity, dt) + INFLOW\n", " buoyancy_force = smoke * (0, buoyancy_factor) >> velocity # resamples smoke to velocity sample points\n", " velocity = advect.semi_lagrangian(velocity, velocity, 1) + dt * buoyancy_force\n", - " velocity = diffuse.explicit(velocity, nu, dt)\n", - " velocity, pressure, iterations, divergence = fluid.make_incompressible(velocity, domain, pressure_guess=pressure)\n", + " velocity = diffuse.explicit(velocity, NU, dt)\n", + " velocity, pressure, iterations, divergence = fluid.make_incompressible(velocity, DOMAIN, pressure_guess=pressure)\n", " return velocity, smoke, pressure\n", "\n", - "velocity, smoke, pressure = step(velocity, smoke, None)\n", + "velocity, smoke, pressure = step(velocity, smoke, None, dt=DT)\n", "\n", "print(\"Max. velocity and mean density: \" + format( [ math.max(velocity.values) , math.mean(smoke.values) ] ))\n", "\n", @@ -184,36 +182,88 @@ "\n", "The pressure projection step in `make_incompressible` is typically the computationally most expensive step in the sequence above. It solves a Poisson equation for the boundary conditions of the domain, and updates the velocity field with the gradient of the computed pressure.\n", "\n", - "Just for testing, we've also printed the mean value of the velocities, and the max density after the update. As you can see in the resulting image, we have a first round region of smoke, with a slight upwards motion (which does not show here yet). " + "Just for testing, we've also printed the mean value of the velocities, and the max density after the update. As you can see in the resulting image, we have a first round region of smoke, with a slight upwards motion (which does not show here yet). \n", + "\n", + "## Datatypes and Dimensions\n", + "\n", + "The created grids are instances of the class `Grid`.\n", + "Like tensors, grids also have the `shape` attribute which lists all batch, spatial and channel dimensions.\n", + "[Shapes in ΦFlow](https://tum-pbs.github.io/PhiFlow/Math.html#shapes) store not only the sizes of the dimensions but also their names and types." ] }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 5, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Smoke: (x=40, y=32)\n", + "Velocity: (x=40, y=32, vector=2)\n", + "Inflow: (x=40, y=32)\n", + "Inflow, spatial only: (x=40, y=32)\n" + ] + } + ], "source": [ - "# ? print(\"Fluid state: \" + format(fluid.state))\n", - "# ? print(\"\\nVelocity content:\")\n", - "# ? [print(grid.data.shape) for grid in fluid.velocity.unstack()];" + "print(f\"Smoke: {smoke.shape}\")\n", + "print(f\"Velocity: {velocity.shape}\")\n", + "print(f\"Inflow: {INFLOW.shape}\")\n", + "print(f\"Inflow, spatial only: {INFLOW.shape.spatial}\")\n", + "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "# TODO\n", - "\n", - "Note that we actually created two variables, one for each velocity component. If you're interested in how this works, have a look at the [Struct documentation](https://github.com/tum-pbs/PhiFlow/blob/master/documentation/Structs.ipynb).\n", - "\n", - "If you look closely at the output from the last `print` command, you'll notice that the shapes of the variables differ. This is because the velocity is sampled in [staggered form](https://github.com/tum-pbs/PhiFlow/blob/master/documentation/Staggered_Grids.md).\n", + "The grid values can be accessed using the `values` property." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "(x=40, y=32) float32 0.0 < ... < 0.2\n", + "(x=(41, 40) along vector, y=(32, 33) along vector, vector=2) float32 -0.10927753 < ... < 0.135354\n", + "(x=40, y=32) float32 0.0 < ... < 0.2\n" + ] + } + ], + "source": [ + "print(smoke.values)\n", + "print(velocity.values)\n", + "print(INFLOW.values)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Grids have many more properties which are documented [here](https://tum-pbs.github.io/PhiFlow/phi/field/#phi.field.Grid).\n", + "Also note that the staggered grid has a [non-uniform shape](https://tum-pbs.github.io/PhiFlow/Math.html#non-uniform-tensors) because the number of faces is not equal to the number of cells. The `INFLOW` grid naturally has the same dimensions as the `smoke` grid.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Time Evolution\n", "\n", "With this setup, we can easily advance the simulation forward in time a bit more by repeatedly calling the `step` function." ] }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 7, "metadata": { "colab": { "base_uri": "https://localhost:8080/" @@ -241,7 +291,7 @@ ], "source": [ "for time_step in range(10):\n", - " velocity, smoke, pressure = step(velocity, smoke, pressure)\n", + " velocity, smoke, pressure = step(velocity, smoke, pressure, dt=DT)\n", " print('Computed frame {}, max velocity {}'.format(time_step , np.asarray(math.max(velocity.values)) ))\n" ] }, @@ -256,7 +306,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 8, "metadata": { "colab": { "base_uri": "https://localhost:8080/", @@ -269,10 +319,10 @@ { "data": { "text/plain": [ - "" + "" ] }, - "execution_count": 6, + "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, @@ -304,7 +354,7 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 9, "metadata": { "colab": { "base_uri": "https://localhost:8080/", @@ -342,7 +392,7 @@ "for time_step in range(20):\n", " if time_step<3 or time_step%10==0: \n", " print('Computing time step %d' % time_step)\n", - " velocity, smoke, pressure = step(velocity, smoke, pressure)\n", + " velocity, smoke, pressure = step(velocity, smoke, pressure, dt=DT)\n", " if time_step%5==0:\n", " steps.append(smoke.values.numpy('y,x'))\n", "\n",