Feature/new demo docs

.gitignore

@@ -28,3 +28,6 @@ cashocs.egg-info/
 /tests/temp
 /tests/out
 *cashocs_remesh_*/
+
+### Documentation
+docs/source/user/demos/*/*.md

demos/documented/cashocs_as_solver/control_solver/demo_control_solver.py

@@ -1,24 +1,62 @@
-# Copyright (C) 2020-2022 Sebastian Blauth
+# ---
+# jupyter:
+#   jupytext:
+#     text_representation:
+#       extension: .py
+#       format_name: light
+#       format_version: '1.5'
+#       jupytext_version: 1.14.4
+# ---
+
+# ```{eval-rst}
+# .. include:: ../../../global.rst
+# ```
 #
-# This file is part of cashocs.
+# (demo_control_solver)=
+# # cashocs as Solver for Optimal Control Problems
 #
-# cashocs is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
+# ## Problem Formulation
 #
-# cashocs is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-# GNU General Public License for more details.
+# As a model problem we again consider the one from {ref}`demo_poisson`, but now
+# we do not use cashocs to compute the adjoint equation and derivative of the (reduced)
+# cost functional, but supply these terms ourselves. The optimization problem is
+# given by
 #
-# You should have received a copy of the GNU General Public License
-# along with cashocs.  If not, see <https://www.gnu.org/licenses/>.
-
-"""For the documentation of this demo see https://cashocs.readthedocs.io/en/latest/demos/cashocs_as_solver/doc_control_solver.html.
-
-"""
+# $$
+# \begin{align}
+#     &\min\; J(y,u) = \frac{1}{2} \int_{\Omega} \left( y - y_d \right)^2
+#     \text{ d}x + \frac{\alpha}{2} \int_{\Omega} u^2 \text{ d}x \\
+#     &\text{ subject to } \qquad
+#     \begin{alignedat}[t]{2}
+#         -\Delta y &= u \quad &&\text{ in } \Omega,\\
+#         y &= 0 \quad &&\text{ on } \Gamma.
+#     \end{alignedat}
+# \end{align}
+# $$
+#
+# (see, e.g., [Tröltzsch - Optimal Control of Partial Differential Equations](
+# https://doi.org/10.1090/gsm/112) or [Hinze, Pinnau, Ulbrich, and Ulbrich -
+# Optimization with PDE constraints](
+# https://doi.org/10.1007/978-1-4020-8839-1)).
+#
+# For the numerical solution of this problem we consider exactly the same setting as
+# in {ref}`demo_poisson`, and most of the code will be identical to this.
+#
+# ## Implementation
+#
+# The complete python code can be found in the file {download}`demo_control_solver.py
+# </../../demos/documented/cashocs_as_solver/control_solver/demo_control_solver.py>`
+# and the corresponding config can be found in {download}`config.ini
+# </../../demos/documented/cashocs_as_solver/control_solver/config.ini>`.
+#
+# ### Recapitulation of {ref}`demo_poisson`
+#
+# For using cashocs exclusively as solver, the procedure is very similar to regularly
+# using it, with a few additions after defining the optimization problem. In particular,
+# up to the initialization of the optimization problem, our code is exactly the same as
+# in {ref}`demo_poisson`, i.e., we use
 
+# +
 from fenics import *
 
 import cashocs
@@ -42,21 +80,92 @@
 )
 
 ocp = cashocs.OptimalControlProblem(e, bcs, J, y, u, p, config=config)
+# -
+
+# ### Supplying Custom Adjoint Systems and Derivatives
+#
+# When using cashocs as a solver, the user can specify their custom weak forms of
+# the adjoint system and of the derivative of the reduced cost functional. For our
+# optimization problem, the adjoint equation is given by
+#
+# $$
+# \begin{alignedat}{2}
+#     - \Delta p &= y - y_d \quad &&\text{ in } \Omega, \\
+#     p &= 0 \quad &&\text{ on } \Gamma,
+# \end{alignedat}
+# $$
+#
+# and the derivative of the cost functional is then given by
+#
+# $$
+# dJ(u)[h] = \int_\Omega (\alpha u + p) h \text{ d}x.
+# $$
+#
+# :::{note}
+# For a detailed derivation and discussion of these objects we refer to
+# [Tröltzsch - Optimal Control of Partial Differential Equations](
+# https://doi.org/10.1090/gsm/112) or [Hinze, Pinnau, Ulbrich, and Ulbrich -
+# Optimization with PDE constraints](https://doi.org/10.1007/978-1-4020-8839-1).
+# :::
+#
+# To specify that cashocs should use these equations instead of the automatically
+# computed ones, we have the following code. First, we specify the derivative
+# of the reduced cost functional via
 
 dJ = Constant(alpha) * u * TestFunction(V) * dx + p * TestFunction(V) * dx
 
+# Afterwards, the weak form for the adjoint system is given by
+
 adjoint_form = (
     inner(grad(p), grad(TestFunction(V))) * dx - (y - y_d) * TestFunction(V) * dx
 )
 adjoint_bcs = bcs
 
+# where we can "recycle" the homogeneous Dirichlet boundary conditions used for the
+# state problem.
+#
+# For both objects, one has to define them as a single UFL form for cashocs, as with the
+# state system and cost functional. In particular, the adjoint weak form has to be in
+# the form of a nonlinear variational problem, so that
+# {python}`fenics.solve(adjoint_form == 0, p, adjoint_bcs)` could be used to solve it.
+# In particular, both forms have to include {py:class}`fenics.TestFunction` objects from
+# the control space and adjoint space, respectively, and must not contain
+# {py:class}`fenics.TrialFunction` objects.
+#
+# These objects are then supplied to the
+# {py:class}`OptimalControlProblem <cashocs.OptimalControlProblem>` via
+
 ocp.supply_custom_forms(dJ, adjoint_form, adjoint_bcs)
 
-ocp.solve()
+# :::{note}
+# One can also specify either the adjoint system or the derivative of the cost
+# functional, using the methods {py:meth}`supply_adjoint_forms
+# <cashocs.OptimalControlProblem.supply_adjoint_forms>` or {py:meth}`supply_derivatives
+# <cashocs.OptimalControlProblem.supply_derivatives>`.
+# However, this is potentially dangerous, due to the following. The adjoint system
+# is a linear system, and there is no fixed convention for the sign of the adjoint
+# state. Hence, supplying, e.g., only the adjoint system, might not be compatible with
+# the derivative of the cost functional which cashocs computes. In effect, the sign
+# is specified by the choice of adding or subtracting the PDE constraint from the
+# cost functional for the definition of a Lagrangian function, which is used to
+# determine the adjoint system and derivative. cashocs internally uses the convention
+# that the PDE constraint is added, so that, internally, it computes not the adjoint
+# state $p$ as defined by the equations given above, but $-p$ instead.
+# Hence, it is recommended to either specify all respective quantities with the
+# {py:meth}`supply_custom_forms <cashocs.OptimalControlProblem.supply_custom_forms>`
+# method.
+# :::
+#
+# Finally, we can use the {py:meth}`solve <cashocs.OptimalControlProblem.solve>` method
+# to solve the problem with the line
 
+ocp.solve()
 
-### Post Processing
+# as in {ref}`demo_poisson`.
+#
+# We visualize the results with the lines
 
+# +
 import matplotlib.pyplot as plt
 
 plt.figure(figsize=(15, 5))
@@ -78,3 +187,16 @@
 
 plt.tight_layout()
 # plt.savefig('./img_poisson.png', dpi=150, bbox_inches='tight')
+# -
+
+# The results are, of course, identical to {ref}`demo_poisson` and look as follows
+# ![](/../../demos/documented/cashocs_as_solver/control_solver/img_control_solver.png)
+
+# :::{note}
+# In case we have multiple state equations as in {ref}`demo_multiple_variables`,
+# one has to supply ordered lists of adjoint equations and boundary conditions,
+# analogously to the usual procedure for cashocs.
+#
+# In the case of multiple control variables, the derivatives of the reduced cost
+# functional w.r.t. each of these have to be specified, again using an ordered list.
+# :::