Replace Constants with Functions from R-space

demos/burgers-hessian.py

@@ -27,10 +27,11 @@ def get_form(mesh_seq):
     def form(index, solutions):
         u, u_ = solutions["u"]
         P = mesh_seq.time_partition
-        dt = Constant(P.timesteps[index])
 
-        # Specify viscosity coefficient
-        nu = Constant(0.0001)
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(P.timesteps[index])
+        nu = Function(R).assign(0.0001)
 
         # Setup variational problem
         v = TestFunction(u.function_space())

demos/burgers.py

@@ -50,17 +50,19 @@ def get_function_spaces(mesh):
 #
 # Timestepping information associated with a given subinterval
 # can be accessed via the :attr:`TimePartition` attribute of
-# the :class:`MeshSeq`. ::
+# 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, solutions):
         u, u_ = solutions["u"]
         P = mesh_seq.time_partition
-        dt = Constant(P.timesteps[index])
 
-        # Specify viscosity coefficient
-        nu = Constant(0.0001)
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(P.timesteps[index])
+        nu = Function(R).assign(0.0001)
 
         # Setup variational problem
         v = TestFunction(u.function_space())

demos/burgers1.py

@@ -28,10 +28,11 @@ def get_form(mesh_seq):
     def form(index, solutions):
         u, u_ = solutions["u"]
         P = mesh_seq.time_partition
-        dt = Constant(P.timesteps[index])
 
-        # Specify viscosity coefficient
-        nu = Constant(0.0001)
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(P.timesteps[index])
+        nu = Function(R).assign(0.0001)
 
         # Setup variational problem
         v = TestFunction(u.function_space())

demos/burgers2.py

@@ -28,10 +28,11 @@ def get_form(mesh_seq):
     def form(index, solutions):
         u, u_ = solutions["u"]
         P = mesh_seq.time_partition
-        dt = Constant(P.timesteps[index])
 
-        # Specify viscosity coefficient
-        nu = Constant(0.0001)
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(P.timesteps[index])
+        nu = Function(R).assign(0.0001)
 
         # Setup variational problem
         v = TestFunction(u.function_space())

demos/burgers_ee.py

@@ -43,10 +43,11 @@ def get_form(mesh_seq):
     def form(index, solutions):
         u, u_ = solutions["u"]
         P = mesh_seq.time_partition
-        dt = Constant(P.timesteps[index])
 
-        # Specify viscosity coefficient
-        nu = Constant(0.0001)
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(P.timesteps[index])
+        nu = Function(R).assign(0.0001)
 
         # Setup variational problem
         v = TestFunction(u.function_space())

demos/burgers_oo.py

@@ -33,10 +33,11 @@ def get_form(self):
         def form(index, solutions):
             u, u_ = solutions["u"]
             P = self.time_partition
-            dt = Constant(P.timesteps[index])
 
-            # Specify viscosity coefficient
-            nu = Constant(0.0001)
+            # Define constants
+            R = FunctionSpace(mesh_seq[index], "R", 0)
+            dt = Function(R).assign(P.timesteps[index])
+            nu = Function(R).assign(0.0001)
 
             # Setup variational problem
             v = TestFunction(u.function_space())
@@ -84,7 +85,8 @@ def get_initial_condition(self):
 
     @annotate_qoi
     def get_qoi(self, solutions, i):
-        dt = Constant(self.time_partition[i].timestep)
+        R = FunctionSpace(self[i], "R", 0)
+        dt = Function(R).assign(self.time_partition[i].timestep)
 
         def end_time_qoi():
             u = solutions["u"]

demos/burgers_time_integrated.py

@@ -23,10 +23,11 @@ def get_form(mesh_seq):
     def form(index, solutions):
         u, u_ = solutions["u"]
         P = mesh_seq.time_partition
-        dt = Constant(P.timesteps[index])
 
-        # Specify viscosity coefficient
-        nu = Constant(0.0001)
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(P.timesteps[index])
+        nu = Function(R).assign(0.0001)
 
         # Setup variational problem
         v = TestFunction(u.function_space())
@@ -93,12 +94,13 @@ def solver(index, ic):
 #    \;\mathrm dy\;\mathrm dt.
 #
 # Note that in this case we multiply by the timestep.
-# It is wrapped in a :class:`Constant` to avoid
+# It is wrapped in a :class:`Function` from `'R'` space to avoid
 # recompilation if the value is changed. ::
 
 
 def get_qoi(mesh_seq, solutions, i):
-    dt = Constant(mesh_seq.time_partition[i].timestep)
+    R = FunctionSpace(mesh_seq[i], "R", 0)
+    dt = Function(R).assign(mesh_seq.time_partition[i].timestep)
 
     def time_integrated_qoi(t):
         u = solutions["u"]

demos/gray_scott.py

@@ -57,11 +57,12 @@ def form(index, sols):
         psi_a, psi_b = TestFunctions(mesh_seq.function_spaces["ab"][index])
 
         # Define constants
-        dt = Constant(mesh_seq.time_partition[index].timestep)
-        D_a = Constant(8.0e-05)
-        D_b = Constant(4.0e-05)
-        gamma = Constant(0.024)
-        kappa = Constant(0.06)
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(mesh_seq.time_partition[index].timestep)
+        D_a = Function(R).assign(8.0e-05)
+        D_b = Function(R).assign(4.0e-05)
+        gamma = Function(R).assign(0.024)
+        kappa = Function(R).assign(0.06)
 
         # Write the two equations in variational form
         F = (

demos/gray_scott_split.py

@@ -57,11 +57,12 @@ def form(index, sols):
         psi_b = TestFunction(mesh_seq.function_spaces["b"][index])
 
         # Define constants
-        dt = Constant(mesh_seq.time_partition[index].timestep)
-        D_a = Constant(8.0e-05)
-        D_b = Constant(4.0e-05)
-        gamma = Constant(0.024)
-        kappa = Constant(0.06)
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(mesh_seq.time_partition[index].timestep)
+        D_a = Function(R).assign(8.0e-05)
+        D_b = Function(R).assign(4.0e-05)
+        gamma = Function(R).assign(0.024)
+        kappa = Function(R).assign(0.06)
 
         # Write the two equations in variational form
         F_a = (

demos/point_discharge2d-goal_oriented.py

@@ -36,11 +36,16 @@ def get_form(mesh_seq):
     def form(index, sols):
         c, c_ = sols["c"]
         function_space = mesh_seq.function_spaces["c"][index]
-        D = Constant(0.1)
-        u = Constant(as_vector([1, 0]))
         h = CellSize(mesh_seq[index])
         S = source(mesh_seq[index])
 
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        D = Function(R).assign(0.1)
+        u_x = Function(R).assign(1.0)
+        u_y = Function(R).assign(0.0)
+        u = as_vector([u_x, u_y])
+
         # SUPG stabilisation parameter
         unorm = sqrt(dot(u, u))
         tau = 0.5 * h / unorm

demos/point_discharge2d-hessian.py

@@ -40,11 +40,16 @@ def get_form(mesh_seq):
     def form(index, sols):
         c, c_ = sols["c"]
         function_space = mesh_seq.function_spaces["c"][index]
-        D = Constant(0.1)
-        u = Constant(as_vector([1, 0]))
         h = CellSize(mesh_seq[index])
         S = source(mesh_seq[index])
 
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        D = Function(R).assign(0.1)
+        u_x = Function(R).assign(1.0)
+        u_y = Function(R).assign(0.0)
+        u = as_vector([u_x, u_y])
+
         # SUPG stabilisation parameter
         unorm = sqrt(dot(u, u))
         tau = 0.5 * h / unorm

demos/point_discharge2d.py

@@ -71,11 +71,16 @@ def get_form(mesh_seq):
     def form(index, sols):
         c, c_ = sols["c"]
         function_space = mesh_seq.function_spaces["c"][index]
-        D = Constant(0.1)
-        u = Constant(as_vector([1, 0]))
         h = CellSize(mesh_seq[index])
         S = source(mesh_seq[index])
 
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        D = Function(R).assign(0.1)
+        u_x = Function(R).assign(1.0)
+        u_y = Function(R).assign(0.0)
+        u = as_vector([u_x, u_y])
+
         # SUPG stabilisation parameter
         unorm = sqrt(dot(u, u))
         tau = 0.5 * h / unorm

demos/solid_body_rotation.py

@@ -151,10 +151,14 @@ def form(index, sols, field="c"):
         V = mesh_seq.function_spaces[field][index]
         mesh = mesh_seq[index]
 
+        # Define velocity field
         x, y = SpatialCoordinate(mesh)
         u = as_vector([-y, x])
-        dt = Constant(mesh_seq.time_partition[index].timestep)
-        theta = Constant(0.5)
+
+        # Define constants
+        R = FunctionSpace(mesh_seq[index], "R", 0)
+        dt = Function(R).assign(mesh_seq.time_partition[index].timestep)
+        theta = Function(R).assign(0.5)
 
         psi = TrialFunction(V)
         phi = TestFunction(V)

test_adjoint/examples/burgers.py

@@ -39,8 +39,9 @@ def form(i, sols):
         u, u_ = sols["uv_2d"]
         dt = self.time_partition[i].timestep
         fs = self.function_spaces["uv_2d"][i]
-        dtc = Constant(dt)
-        nu = Constant(0.0001)
+        R = FunctionSpace(self[i], "R", 0)
+        dtc = Function(R).assign(dt)
+        nu = Function(R).assign(0.0001)
         v = TestFunction(fs)
         F = (
             inner((u - u_) / dtc, v) * dx
@@ -104,7 +105,8 @@ def get_qoi(self, sol, i):
     norm over the right hand
     boundary.
     """
-    dtc = Constant(self.time_partition[i].timestep)
+    R = FunctionSpace(self[i], "R", 0)
+    dtc = Function(R).assign(self.time_partition[i].timestep)
 
     def time_integrated_qoi(t):
         u = sol["uv_2d"]

test_adjoint/examples/point_discharge2d.py

@@ -52,15 +52,17 @@ def source(mesh):
 
 def get_form(self):
     """
-    Advection-diffusion with SUPG
-    stabilisation.
+    Advection-diffusion with SUPG stabilisation.
     """
 
     def form(i, sols):
         c, c_ = sols["tracer_2d"]
         fs = self.function_spaces["tracer_2d"][i]
-        D = Constant(0.1)
-        u = Constant(as_vector([1.0, 0.0]))
+        R = FunctionSpace(self[i], "R", 0)
+        D = Function(R).assign(0.1)
+        u_x = Function(R).assign(1.0)
+        u_y = Function(R).assign(0.0)
+        u = as_vector([u_x, u_y])
         h = CellSize(self[i])
         S = source(self[i])
 
@@ -151,8 +153,9 @@ def analytical_solution(mesh):
     a given mesh. See [Riadh et al. 2014].
     """
     x, y = SpatialCoordinate(mesh)
-    u = Constant(1.0)
-    D = Constant(0.1)
+    R = FunctionSpace(mesh, "R", 0)
+    u = Function(R).assign(1.0)
+    D = Function(R).assign(0.1)
     Pe = 0.5 * u / D
     r = max_value(sqrt((x - src_x) ** 2 + (y - src_y) ** 2), src_r)
     return 0.5 / (pi * D) * exp(Pe * (x - src_x)) * bessk0(Pe * r)

test_adjoint/examples/point_discharge3d.py

@@ -70,8 +70,12 @@ def get_form(self):
     def form(i, sols):
         c, c_ = sols["tracer_3d"]
         fs = self.function_spaces["tracer_3d"][i]
-        D = Constant(0.1)
-        u = Constant(as_vector([1.0, 0.0, 0.0]))
+        R = FunctionSpace(self[i], "R", 0)
+        D = Function(R).assign(0.1)
+        u_x = Function(R).assign(1.0)
+        u_y = Function(R).assign(0.0)
+        u_z = Function(R).assign(0.0)
+        u = as_vector([u_x, u_y, u_z])
         h = CellSize(self[i])
         S = source(self[i])
 
@@ -164,8 +168,9 @@ def analytical_solution(mesh):
     a given mesh.
     """
     x, y, z = SpatialCoordinate(mesh)
-    u = Constant(1.0)
-    D = Constant(0.1)
+    R = FunctionSpace(mesh, "R", 0)
+    u = Function(R).assign(1.0)
+    D = Function(R).assign(0.1)
     Pe = 0.5 * u / D
     r = max_value(sqrt((x - src_x) ** 2 + (y - src_y) ** 2 + (z - src_z) ** 2), src_r)
     return 0.5 / (pi * D) * exp(Pe * (x - src_x)) * bessk0(Pe * r)

test_adjoint/examples/steady_flow_past_cyl.py

@@ -38,7 +38,8 @@ def get_form(self):
     def form(i, sols):
         up, up_ = sols["up"]
         W = self.function_spaces["up"][i]
-        nu = Constant(1.0)
+        R = FunctionSpace(self[i], "R", 0)
+        nu = Function(R).assign(1.0)
         u, p = split(up)
         v, q = TestFunctions(W)
         F = (

test_adjoint/test_fp_iteration.py

@@ -30,11 +30,10 @@ def solver(index, ic):
 
 
 def get_qoi(mesh_seq, solutions, index):
+    R = FunctionSpace(mesh_seq[index], "R", 0)
+
     def qoi():
-        if mesh_seq.fp_iteration % 2 == 0:
-            return Constant(1, domain=mesh_seq[index]) * dx
-        else:
-            return Constant(2, domain=mesh_seq[index]) * dx
+        return Function(R).assign(1 if mesh_seq.fp_iteration % 2 == 0 else 2) * dx
 
     return qoi