Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Do not require get_form in user code #228

Merged
merged 13 commits into from
Nov 17, 2024
16 changes: 2 additions & 14 deletions demos/burgers-hessian.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,8 +23,8 @@ def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}


def get_form(mesh_seq):
def form(index):
def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define constants
Expand All @@ -39,17 +39,6 @@ def form(index):
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form


def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define form
F = mesh_seq.form(index)["u"]

# Time integrate from t_start to t_end
tp = mesh_seq.time_partition
Expand Down Expand Up @@ -91,7 +80,6 @@ def get_initial_condition(mesh_seq):
meshes,
get_function_spaces=get_function_spaces,
get_initial_condition=get_initial_condition,
get_form=get_form,
get_solver=get_solver,
)

Expand Down
62 changes: 17 additions & 45 deletions demos/burgers.py
Original file line number Diff line number Diff line change
Expand Up @@ -40,25 +40,29 @@ def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}


# In order to solve PDEs using the finite element method, we
# require a weak form. For this, Goalie requires a function
# that maps the :class:`MeshSeq` index and a dictionary of
# solution data to the form. The form should be
# returned inside its own dictionary, with an entry for each equation
# being solved for.
#
# The solution :class:`Function`\s are automatically built on the function spaces given
# by the :func:`get_function_spaces` function and are accessed via the :attr:`fields`
# attribute of the :class:`MeshSeq`. This attribute provides a dictionary of tuples
# containing the current and lagged solutions for each field.
# Similarly, timestepping information associated with a given subinterval
# can be accessed via the :attr:`TimePartition` attribute of
# the :class:`MeshSeq`. For technical reasons, we need to create a :class:`Function`
# in the `'R'` space (of real numbers) to hold constants. ::
#
# In order to solve the PDE, we need to choose a time integration routine and solver
# parameters for the underlying linear and nonlinear systems. This is achieved below by
# using a function :func:`solver` whose input is the :class:`MeshSeq` index. The
# function should return a generator that yields the solution at each timestep, so
# that Goalie can efficiently track the solution history. This is done by using the
# `yield` statement before progressing to the next timestep.
#
# The lagged solution is assigned the initial condition for the current subinterval
# index. For the :math:`0^{th}` index, this will be provided by the initial conditions,
# otherwise it will be transferred from the previous mesh in the sequence.
# Timestepping information associated with a given subinterval can be accessed via the
# :attr:`TimePartition` attribute of the :class:`MeshSeq`. For technical reasons, we
# need to create a :class:`Function` in the `'R'` space (of real numbers) to hold
# constants.::


def get_form(mesh_seq):
def form(index):
def get_solver(mesh_seq):
def solver(index):
# Get the current and lagged solutions
u, u_ = mesh_seq.fields["u"]

Expand All @@ -74,37 +78,6 @@ def form(index):
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form


# We have a weak form. The dictionary usage may seem cumbersome when applied to such a
# simple problem, but it comes in handy when solving adjoint problems associated with
# coupled systems of equations.

# In order to solve the PDE, we need to choose
# a time integration routine and solver parameters for the underlying
# linear and nonlinear systems. This is achieved below by using a function
# :func:`solver` whose input is the :class:`MeshSeq` index. As noted above, the solution
# :class:`Function`\s are automatically initialised and accessed via the
# :attr:`fields` attribute of the :class:`MeshSeq`. The lagged solution is assigned
# the initial condition for the current subinterval index. For the :math:`0^{th}` index,
# this will be provided by the initial conditions, otherwise it will be transferred
# from the previous mesh in the sequence.
#
# The function should return a generator that yields the solution at each timestep, so
# that Goalie can efficiently track the solution history. This is done by using the
# `yield` statement before progressing to the next timestep. ::


def get_solver(mesh_seq):
def solver(index):
# Get the current and lagged solutions
u, u_ = mesh_seq.fields["u"]

# Define form
F = mesh_seq.form(index)["u"]

# Time integrate from t_start to t_end
tp = mesh_seq.time_partition
Expand Down Expand Up @@ -164,7 +137,6 @@ def get_initial_condition(mesh_seq):
meshes,
get_function_spaces=get_function_spaces,
get_initial_condition=get_initial_condition,
get_form=get_form,
get_solver=get_solver,
)
solutions = mesh_seq.solve_forward()
Expand Down
18 changes: 3 additions & 15 deletions demos/burgers1.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@
from goalie_adjoint import *

# For ease, the list of field names and functions for obtaining the
# function spaces, forms, solvers, and initial conditions
# function spaces, solvers, and initial conditions
# are redefined as in the previous demo. The only difference
# is that now we are solving the adjoint problem, which
# requires that the PDE solve is labelled with an
Expand All @@ -29,8 +29,8 @@ def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}


def get_form(mesh_seq):
def form(index):
def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define constants
Expand All @@ -45,17 +45,6 @@ def form(index):
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form


def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define form
F = mesh_seq.form(index)["u"]

# Time integrate from t_start to t_end
tp = mesh_seq.time_partition
Expand Down Expand Up @@ -125,7 +114,6 @@ def end_time_qoi():
mesh,
get_function_spaces=get_function_spaces,
get_initial_condition=get_initial_condition,
get_form=get_form,
get_solver=get_solver,
get_qoi=get_qoi,
qoi_type="end_time",
Expand Down
16 changes: 2 additions & 14 deletions demos/burgers2.py
Original file line number Diff line number Diff line change
Expand Up @@ -24,8 +24,8 @@ def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}


def get_form(mesh_seq):
def form(index):
def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define constants
Expand All @@ -40,17 +40,6 @@ def form(index):
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form


def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define form
F = mesh_seq.form(index)["u"]

# Time integrate from t_start to t_end
tp = mesh_seq.time_partition
Expand Down Expand Up @@ -105,7 +94,6 @@ def end_time_qoi():
meshes,
get_function_spaces=get_function_spaces,
get_initial_condition=get_initial_condition,
get_form=get_form,
get_solver=get_solver,
get_qoi=get_qoi,
qoi_type="end_time",
Expand Down
22 changes: 8 additions & 14 deletions demos/burgers_ee.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,10 @@
set_log_level(DEBUG)

# Redefine the ``field_names`` variable and the getter functions as in the first
# adjoint Burgers demo. ::
# adjoint Burgers demo. The only difference is the inclusion of the
# :meth:`GoalOrientedMeshSeq.read_forms()` method in the ``get_solver`` function. The
# method is used to communicate the variational form to the mesh sequence object so that
# Goalie can utilise it in the error estimation process described above. ::

field_names = ["u"]

Expand All @@ -40,8 +43,8 @@ def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}


def get_form(mesh_seq):
def form(index):
def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define constants
Expand All @@ -56,17 +59,9 @@ def form(index):
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form


def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define form
F = mesh_seq.form(index)["u"]
# Communicate variational form to mesh_seq
mesh_seq.read_forms({"u": F})

# Time integrate from t_start to t_end
P = mesh_seq.time_partition
Expand Down Expand Up @@ -122,7 +117,6 @@ def end_time_qoi():
meshes,
get_function_spaces=get_function_spaces,
get_initial_condition=get_initial_condition,
get_form=get_form,
get_solver=get_solver,
get_qoi=get_qoi,
qoi_type="end_time",
Expand Down
15 changes: 4 additions & 11 deletions demos/burgers_oo.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,8 +30,8 @@ class BurgersMeshSeq(GoalOrientedMeshSeq):
def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}

def get_form(self):
def form(index):
def get_solver(self):
def solver(index):
u, u_ = self.fields["u"]

# Define constants
Expand All @@ -46,16 +46,9 @@ def form(index):
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form

def get_solver(self):
def solver(index):
u, u_ = self.fields["u"]

# Define form
F = self.form(index)["u"]
# Communicate variational form to mesh_seq
self.read_forms({"u": F})

# Time integrate from t_start to t_end
t_start, t_end = self.subintervals[index]
Expand Down
40 changes: 14 additions & 26 deletions demos/burgers_time_integrated.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,35 +12,14 @@

from goalie_adjoint import *

# Redefine the ``get_initial_condition``, ``get_function_spaces``,
# and ``get_form`` functions as in the first Burgers demo. ::
# Redefine the ``get_initial_condition`` and ``get_function_spaces``,
# functions as in the first Burgers demo. ::


def get_function_spaces(mesh):
return {"u": VectorFunctionSpace(mesh, "CG", 2)}


def get_form(mesh_seq):
def form(index):
u, u_ = mesh_seq.fields["u"]

# Define constants
R = FunctionSpace(mesh_seq[index], "R", 0)
dt = Function(R).assign(mesh_seq.time_partition.timesteps[index])
nu = Function(R).assign(0.0001)

# Setup variational problem
v = TestFunction(u.function_space())
F = (
inner((u - u_) / dt, v) * dx
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)
return {"u": F}

return form


def get_initial_condition(mesh_seq):
fs = mesh_seq.function_spaces["u"][0]
x, y = SpatialCoordinate(mesh_seq[0])
Expand All @@ -61,8 +40,18 @@ def get_solver(mesh_seq):
def solver(index):
u, u_ = mesh_seq.fields["u"]

# Define form
F = mesh_seq.form(index)["u"]
# Define constants
R = FunctionSpace(mesh_seq[index], "R", 0)
dt = Function(R).assign(mesh_seq.time_partition.timesteps[index])
nu = Function(R).assign(0.0001)

# Setup variational problem
v = TestFunction(u.function_space())
F = (
inner((u - u_) / dt, v) * dx
+ inner(dot(u, nabla_grad(u)), v) * dx
+ nu * inner(grad(u), grad(v)) * dx
)

# Time integrate from t_start to t_end
t_start, t_end = mesh_seq.subintervals[index]
Expand Down Expand Up @@ -126,7 +115,6 @@ def time_integrated_qoi(t):
meshes,
get_function_spaces=get_function_spaces,
get_initial_condition=get_initial_condition,
get_form=get_form,
get_solver=get_solver,
get_qoi=get_qoi,
qoi_type="time_integrated",
Expand Down
Loading
Loading