Skip to content

Commit

Permalink
[mod] misc
Browse files Browse the repository at this point in the history
  • Loading branch information
NaokiHori committed Jul 17, 2023
1 parent 06e9313 commit 27c9af0
Show file tree
Hide file tree
Showing 17 changed files with 266 additions and 146 deletions.
7 changes: 3 additions & 4 deletions README.rst
Original file line number Diff line number Diff line change
Expand Up @@ -82,7 +82,8 @@ Quick start
#. Set initial condition

A solenoidal velocity field should be given.
Although the scalar field can be arbitrary, the velocity field should be solenoidal.
``main.py`` offers several examples.

.. code-block:: console
Expand All @@ -106,7 +107,7 @@ Quick start
bash exec.sh
This may take a few minutes, depending on your machine spec.
This may take a few minutes, depending on your machine spec.

The flow fields are saved under ``output/save/`` in `NPY <https://numpy.org/devdocs/reference/generated/numpy.lib.format.html>`_ format.
Note that these velocities are in the spectral domain; you need to perform the inverse Fourier transform (and the normalisation) to recover the velocities in the physical domain.
Expand All @@ -117,8 +118,6 @@ If proper Python libraries are installed, you can visualise the flow fields by
python3 visualise/2d.py
This should give you a similar movie to the one you can see above.

*************
3D simulation
*************
Expand Down
9 changes: 5 additions & 4 deletions docs/source/equation/main.rst
Original file line number Diff line number Diff line change
Expand Up @@ -15,10 +15,10 @@ and the momentum balance:
.. math::
\der{u_i}{t}
+
u_j \der{u_i}{x_j}
=
-
u_j \der{u_i}{x_j}
-
\der{p}{x_i}
+
\frac{1}{Re}
Expand All @@ -32,9 +32,10 @@ Also a passive scalar field :math:`T` is transported, which is governed by the a
.. math::
\der{T}{t}
+
u_j \der{T}{x_j}
=
-
u_j \der{T}{x_j}
+
\frac{1}{Re Sc}
\der{}{x_j} \der{T}{x_j},
Expand Down
61 changes: 33 additions & 28 deletions docs/source/equation/spectral.rst
Original file line number Diff line number Diff line change
@@ -1,3 +1,6 @@

.. _spectral:

################################
Equations in the spectral domain
################################
Expand Down Expand Up @@ -44,46 +47,49 @@ The external forcing terms :math:`\wav{a_x}{\kx \ky}` and :math:`\wav{a_y}{\kx \
\lx \ly
\left[
-
I \kx \frac{2 \pi}{\lx}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\ux}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\ux}{\kx_1^{\prime} \ky_1^{\prime}}
+
-
I \ky \frac{2 \pi}{\ly}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\uy}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\ux}{\kx_1^{\prime} \ky_1^{\prime}}
-
+
\wav{a_x}{\kx \ky}
\right],
\lx \ly
\left[
-
I \kx \frac{2 \pi}{\lx}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\ux}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\uy}{\kx_1^{\prime} \ky_1^{\prime}}
+
-
I \ky \frac{2 \pi}{\ly}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\uy}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\uy}{\kx_1^{\prime} \ky_1^{\prime}}
-
+
\wav{a_y}{\kx \ky}
\right],
\lx \ly
\left[
-
I \kx \frac{2 \pi}{\lx}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\ux}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{ T}{\kx_1^{\prime} \ky_1^{\prime}}
+
-
I \ky \frac{2 \pi}{\ly}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
Expand Down Expand Up @@ -134,27 +140,24 @@ Summary
.. math::
\der{\wav{\ux}{\kx \ky}}{t}
+
\wav{h_x}{\kx \ky}
=
\wav{h_x}{\kx \ky}
-
I \kx \frac{2 \pi}{\lx} \wav{p}{\kx \ky}
-
\frac{1}{Re} \left[ \left( \kx \frac{2 \pi}{\lx} \right)^2 + \left( \ky \frac{2 \pi}{\ly} \right)^2 \right] \wav{\ux}{\kx \ky},
\der{\wav{\uy}{\kx \ky}}{t}
+
\wav{h_y}{\kx \ky}
=
\wav{h_y}{\kx \ky}
-
I \ky \frac{2 \pi}{\ly} \wav{p}{\kx \ky}
-
\frac{1}{Re} \left[ \left( \kx \frac{2 \pi}{\lx} \right)^2 + \left( \ky \frac{2 \pi}{\ly} \right)^2 \right] \wav{\uy}{\kx \ky},
\der{\wav{T}{\kx \ky}}{t}
+
\wav{g}{\kx \ky}
=
\wav{g}{\kx \ky}
-
\frac{1}{Re Sc} \left[ \left( \kx \frac{2 \pi}{\lx} \right)^2 + \left( \ky \frac{2 \pi}{\ly} \right)^2 \right] \wav{T}{\kx \ky},
Expand All @@ -164,60 +167,63 @@ where
\wav{h_x}{\kx \ky}
\equiv
-
I \kx \frac{2 \pi}{\lx}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\ux}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\ux}{\kx_1^{\prime} \ky_1^{\prime}}
+
-
I \ky \frac{2 \pi}{\ly}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\uy}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\ux}{\kx_1^{\prime} \ky_1^{\prime}}
-
+
\wav{a_x}{\kx \ky},
\wav{h_y}{\kx \ky}
\equiv
-
I \kx \frac{2 \pi}{\lx}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\ux}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\uy}{\kx_1^{\prime} \ky_1^{\prime}}
+
-
I \ky \frac{2 \pi}{\ly}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\uy}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{\uy}{\kx_1^{\prime} \ky_1^{\prime}}
-
+
\wav{a_y}{\kx \ky},
\wav{g}{\kx \ky}
\equiv
-
I \kx \frac{2 \pi}{\lx}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\ux}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{ T}{\kx_1^{\prime} \ky_1^{\prime}}
+
-
I \ky \frac{2 \pi}{\ly}
\sum_{\kx_0^{\prime} + \kx_1^{\prime} = \kx}
\sum_{\ky_0^{\prime} + \ky_1^{\prime} = \ky}
\wav{\uy}{\kx_0^{\prime} \ky_0^{\prime}}
\wav{ T}{\kx_1^{\prime} \ky_1^{\prime}}.
To eliminate the pressure from the momentum equation, I consider to compute the inner product of the wave vector and the momentum balance, namely the sum of
To eliminate the pressure from the momentum equation, I consider the inner product of the wave vector and the momentum balance, namely the sum of

.. math::
I \kx \frac{2 \pi}{\lx}
\der{\wav{\ux}{\kx \ky}}{t}
+
=
I \kx \frac{2 \pi}{\lx}
\wav{h_x}{\kx \ky}
=
+
\left( \kx \frac{2 \pi}{\lx} \right)^2 \wav{p}{\kx \ky}
-
\frac{1}{Re} \left[
Expand All @@ -233,10 +239,10 @@ and
I \ky \frac{2 \pi}{\ly}
\der{\wav{\uy}{\kx \ky}}{t}
+
=
I \ky \frac{2 \pi}{\ly}
\wav{h_y}{\kx \ky}
=
+
\left( \ky \frac{2 \pi}{\ly} \right)^2 \wav{p}{\kx \ky}
-
\frac{1}{Re} \left[
Expand All @@ -250,9 +256,10 @@ Since the temporal derivative terms and the diffusive terms are zero because of

.. math::
-
I \kx \frac{2 \pi}{\lx}
\wav{h_x}{\kx \ky}
+
-
I \ky \frac{2 \pi}{\ly}
\wav{h_y}{\kx \ky}
=
Expand All @@ -269,9 +276,8 @@ which is used to eliminate the pressure, giving
\der{\wav{\ux}{\mkx \mky}}{t}
=
-
\wav{h_x}{\mkx \mky}
+
-
\frac{\mkx}{\mkx^2 + \mky^2}
\left(
\mkx \wav{h_x}{\mkx \mky}
Expand All @@ -283,19 +289,18 @@ which is used to eliminate the pressure, giving
\der{\wav{\uy}{\mkx \mky}}{t}
=
-
\wav{h_y}{\mkx \mky}
+
-
\frac{\mky}{\mkx^2 + \mky^2}
\left(
\mkx \wav{h_x}{\mkx \mky}
+
\mky \wav{h_y}{\mkx \mky}
\right)
-
\frac{1}{Re} \left( \mkx^2 + \mky^2 \right) \wav{\uy}{\mkx \mky},
\frac{1}{Re} \left( \mkx^2 + \mky^2 \right) \wav{\uy}{\mkx \mky}.
where
Recall that

.. math::
Expand Down
35 changes: 34 additions & 1 deletion docs/source/numerical_method/adv.rst
Original file line number Diff line number Diff line change
Expand Up @@ -2,5 +2,38 @@
Advective terms
###############

TBA
********
Overview
********

To evaluate the convolution sums efficiently, instead of computing

.. math::
\sum_{\ix_0^{\prime} + \ix_1^{\prime} = \ix}
\sum_{\iy_0^{\prime} + \iy_1^{\prime} = \iy}
\wav{p}{\ix_0^{\prime} \iy_0^{\prime}}
\wav{q}{\ix_1^{\prime} \iy_1^{\prime}},
where :math:`p` and :math:`q` are the arguments and whose computational cost is :math:`\mathcal{O} \left( N_x^2 \times N_y^2 \right)`, I consider their representations in the physical domain, i.e. computing the product:

.. math::
p q
directly, whose cost is :math:`\mathcal{O} \left( N_x \log N_x \times N_y \log N_y \right)`, which is known as the transform method.

**************
Implementation
**************

#. ``src/fluid/physical.c``

This file contains a function to compute the description of each flow field in the physical place by performing the multi-dimensional inverse Fourier transform.

#. ``src/fluid/slope.c``

This file contains several functions which evaluate the right-hand-side terms of the Runge-Kutta scheme.
In particular, ``convolute`` computes the convolution sum (product of the given two arrays), while ``compute_adv`` computes the advective terms by multiplying pre-factors corresponding to the spatial derivative in each direction and the pre-factors.
Finally ``project_velocity`` is necessary for the velocity field to satisfy the divergence-free condition (see the :ref:`equations <spectral>`, where two terms are involved in each direction to describe the advection).

Loading

0 comments on commit 27c9af0

Please sign in to comment.