Adjoint

.github/workflows/build.yml

@@ -11,13 +11,6 @@ on:
     # Scheduled build at 0330 UTC on Monday mornings to detect bitrot.
     - cron:  '30 3 * * 1'
 
-concurrency:
-  # Cancels jobs running if new commits are pushed
-  group: >
-    ${{ github.workflow }}-
-    ${{ github.event.pull_request.number || github.ref }}
-  cancel-in-progress: true
-
 jobs:
   build:
     name: "Build Gusto"

examples/compressible_euler/unsaturated_bubble.py

@@ -216,11 +216,16 @@ def unsaturated_bubble(
     R_v = eqns.parameters.R_v
     epsilon = R_d / R_v
 
+    # make expressions for determining water_v0
+    exner = thermodynamics.exner_pressure(eqns.parameters, rho_averaged, theta0)
+    p = thermodynamics.p(eqns.parameters, exner)
+    T = thermodynamics.T(eqns.parameters, theta0, exner, water_v0)
+    r_v_expr = thermodynamics.r_v(eqns.parameters, rel_hum, T, p)
+
     # make expressions to evaluate residual
     exner_expr = thermodynamics.exner_pressure(eqns.parameters, rho_averaged, theta0)
     p_expr = thermodynamics.p(eqns.parameters, exner_expr)
     T_expr = thermodynamics.T(eqns.parameters, theta0, exner_expr, water_v0)
-    water_v_expr = thermodynamics.r_v(eqns.parameters, rel_hum, T_expr, p_expr)
     rel_hum_expr = thermodynamics.RH(eqns.parameters, water_v0, T_expr, p_expr)
     rel_hum_eval = Function(Vt)
 
@@ -242,7 +247,7 @@ def unsaturated_bubble(
 
         # first solve for r_v
         for _ in range(max_inner_solve_count):
-            water_v_eval.interpolate(water_v_expr)
+            water_v_eval.interpolate(r_v_expr)
             water_v0.assign(water_v0 * (1 - delta) + delta * water_v_eval)
 
             # compute theta_vd

gusto/core/configuration.py

@@ -1,7 +1,7 @@
 """Some simple tools for configuring the model."""
-from abc import ABCMeta, abstractproperty
+from abc import ABCMeta, abstractmethod
 from enum import Enum
-from firedrake import sqrt, Constant
+from firedrake import sqrt, Constant, Function, FunctionSpace
 
 
 __all__ = [
@@ -12,7 +12,7 @@
     "EmbeddedDGOptions", "ConservativeEmbeddedDGOptions", "RecoveryOptions",
     "ConservativeRecoveryOptions", "SUPGOptions", "MixedFSOptions",
     "SpongeLayerParameters", "DiffusionParameters", "BoundaryLayerParameters",
-    "SubcyclingOptions"
+    "SubcyclingOptions", "convert_parameters_to_real_space"
 ]
 
 
@@ -79,11 +79,7 @@ def __setattr__(self, name, value):
 
         # Almost all parameters should be Constants -- but there are some
         # specific exceptions which should be kept as integers
-        non_constants = [
-            'dumpfreq', 'pddumpfreq', 'chkptfreq',
-            'fixed_subcycles', 'max_subcycles', 'subcycle_by_courant'
-        ]
-        if type(value) in [float, int] and name not in non_constants:
+        if type(value) in [float, int] and name not in ['dumpfreq', 'pddumpfreq', 'chkptfreq']:
             object.__setattr__(self, name, Constant(value))
         else:
             object.__setattr__(self, name, value)
@@ -167,7 +163,7 @@ class ShallowWaterParameters(Configuration):
 class WrapperOptions(Configuration, metaclass=ABCMeta):
     """Base class for specifying options for a transport scheme."""
 
-    @abstractproperty
+    @abstractmethod
     def name(self):
         pass
 
@@ -269,19 +265,6 @@ class BoundaryLayerParameters(Configuration):
     mu = 100.                   # Interior penalty coefficient for vertical diffusion
 
 
-class HeldSuarezParameters(Configuration):
-    """
-    Parameters used in the default configuration for the Held Suarez test case.
-    """
-    T0stra = 200               # Stratosphere temp
-    T0surf = 315               # Surface temperature at equator
-    T0horiz = 60               # Equator to pole temperature difference
-    T0vert = 10                # Stability parameter
-    sigmab = 0.7               # Height of the boundary layer
-    tau_d = 40 * 24 * 60 * 60  # 40 day time scale
-    tau_u = 4 * 24 * 60 * 60   # 4 day timescale
-
-
 class SubcyclingOptions(Configuration):
     """
     Describes the process of subcycling a time discretisation, by dividing the
@@ -308,3 +291,17 @@ def check_options(self):
             raise ValueError(
                 "Cannot provide both fixed_subcycles and subcycle_by_courant"
                 + "parameters.")
+
+
+def convert_parameters_to_real_space(parameters, mesh):
+    """Convert the float type parameters attributes to real space.
+
+    Args:
+        parameters (:class:`Configuration`): the configuration object containing
+            the parameters to convert
+        mesh (:class:`firedrake.Mesh`): the mesh object to use for the real space.
+    """
+    R = FunctionSpace(mesh, 'R', 0)
+    for name, value in vars(parameters).items():
+        if isinstance(value, (float, Constant)):
+            setattr(parameters, name, Function(R, val=float(value)))

gusto/core/domain.py

@@ -6,10 +6,9 @@
 
 from gusto.core.coordinates import Coordinates
 from gusto.core.function_spaces import Spaces, check_degree_args
-from firedrake import (
-    Constant, SpatialCoordinate, sqrt, CellNormal, cross, inner, grad,
-    VectorFunctionSpace, Function, FunctionSpace, perp, curl
-)
+from firedrake import (Constant, SpatialCoordinate, sqrt, CellNormal, cross,
+                       inner, grad, VectorFunctionSpace, Function, FunctionSpace,
+                       perp)
 import numpy as np
 
 
@@ -64,15 +63,16 @@ def __init__(self, mesh, dt, family, degree=None,
         # -------------------------------------------------------------------- #
 
         # Store central dt for use in the rest of the model
+        R = FunctionSpace(mesh, "R", 0)
         if type(dt) is Constant:
-            self.dt = dt
+            self.dt = Function(R, val=float(dt))
         elif type(dt) in (float, int):
-            self.dt = Constant(dt)
+            self.dt = Function(R, val=dt)
         else:
             raise TypeError(f'dt must be a Constant, float or int, not {type(dt)}')
 
         # Make a placeholder for the time
-        self.t = Constant(0.0)
+        self.t = Function(R, val=0.0)
 
         # -------------------------------------------------------------------- #
         # Build compatible function spaces
@@ -125,14 +125,12 @@ def __init__(self, mesh, dt, family, degree=None,
                 V = VectorFunctionSpace(mesh, "DG", sphere_degree)
                 self.outward_normals = Function(V).interpolate(CellNormal(mesh))
                 self.perp = lambda u: cross(self.outward_normals, u)
-                self.divperp = lambda u: inner(self.outward_normals, curl(u))
         else:
             kvec = [0.0]*dim
             kvec[dim-1] = 1.0
             self.k = Constant(kvec)
             if dim == 2:
                 self.perp = perp
-                self.divperp = lambda u: -u[0].dx(1) + u[1].dx(0)
 
         # -------------------------------------------------------------------- #
         # Construct information relating to height/radius