Use new UFL API

gusto/complex_proxy/mixed_impl.py

@@ -49,7 +49,7 @@ def FunctionSpace(V):
 
     :arg V: the real-valued FunctionSpace.
     """
-    return fd.FunctionSpace(V.mesh(), FiniteElement(V.ufl_element()))
+    return fd.FunctionSpace(V.mesh().unique(), FiniteElement(V.ufl_element()))
 
 
 def DirichletBC(W, V, bc, function_arg=None):

gusto/core/coordinates.py

@@ -51,13 +51,13 @@ def __init__(self, mesh, on_sphere=False, rotated_pole=None, radius=None):
                 self.coords_name = ['lon', 'lat']
         else:
             self.coords = SpatialCoordinate(mesh)
-            if mesh.geometric_dimension() == 1:
+            if mesh.geometric_dimension == 1:
                 self.coords_name = ['x']
-            elif mesh.geometric_dimension() == 2 and mesh.extruded:
+            elif mesh.geometric_dimension == 2 and mesh.extruded:
                 self.coords_name = ['x', 'z']
-            elif mesh.geometric_dimension() == 2:
+            elif mesh.geometric_dimension == 2:
                 self.coords_name = ['x', 'y']
-            elif mesh.geometric_dimension() == 3:
+            elif mesh.geometric_dimension == 3:
                 self.coords_name = ['x', 'y', 'z']
             else:
                 raise ValueError('Cannot work out coordinates of domain')
@@ -86,7 +86,7 @@ def register_space(self, domain, space_name):
         """
 
         comm = self.mesh.comm
-        topological_dimension = self.mesh.topological_dimension()
+        topological_dimension = self.mesh.topological_dimension
 
         if space_name in self.chi_coords.keys():
             logger.warning(f'Coords for {space_name} space have already been computed')

gusto/core/domain.py

@@ -107,12 +107,12 @@ def __init__(self, mesh, dt, family, degree=None,
         # WARNING: if we ever wanted to run on other domains (e.g. circle, disk
         # or torus) then the identification of domains would no longer be unique
         if hasattr(mesh, "_base_mesh") and hasattr(mesh._base_mesh, 'geometric_dimension'):
-            self.on_sphere = (mesh._base_mesh.geometric_dimension() == 3 and mesh._base_mesh.topological_dimension() == 2)
+            self.on_sphere = (mesh._base_mesh.geometric_dimension == 3 and mesh._base_mesh.topological_dimension == 2)
         else:
-            self.on_sphere = (mesh.geometric_dimension() == 3 and mesh.topological_dimension() == 2)
+            self.on_sphere = (mesh.geometric_dimension == 3 and mesh.topological_dimension == 2)
 
         #  build the vertical normal and define perp for 2d geometries
-        dim = mesh.topological_dimension()
+        dim = mesh.topological_dimension
         if self.on_sphere:
             x = SpatialCoordinate(mesh)
             R = sqrt(inner(x, x))
@@ -223,13 +223,13 @@ def construct_domain_metadata(mesh, coords, on_sphere):
         metadata['domain_type'] = 'extruded_spherical_shell'
     elif on_sphere:
         metadata['domain_type'] = 'spherical_shell'
-    elif mesh.geometric_dimension() == 1 and mesh.topological_dimension() == 1:
+    elif mesh.geometric_dimension == 1 and mesh.topological_dimension == 1:
         metadata['domain_type'] = 'interval'
-    elif mesh.geometric_dimension() == 2 and mesh.topological_dimension() == 2 and mesh.extruded:
+    elif mesh.geometric_dimension == 2 and mesh.topological_dimension == 2 and mesh.extruded:
         metadata['domain_type'] = 'vertical_slice'
-    elif mesh.geometric_dimension() == 2 and mesh.topological_dimension() == 2:
+    elif mesh.geometric_dimension == 2 and mesh.topological_dimension == 2:
         metadata['domain_type'] = 'plane'
-    elif mesh.geometric_dimension() == 3 and mesh.topological_dimension() == 3 and mesh.extruded:
+    elif mesh.geometric_dimension == 3 and mesh.topological_dimension == 3 and mesh.extruded:
         metadata['domain_type'] = 'extruded_plane'
     else:
         raise ValueError('Unable to determine domain type')

gusto/core/function_spaces.py

@@ -224,13 +224,13 @@ def build_dg1_equispaced(self):
         """
 
         if self.extruded_mesh:
-            cell = self.mesh._base_mesh.ufl_cell().cellname()
+            cell = self.mesh._base_mesh.ufl_cell().cellname
             hori_elt = FiniteElement('DG', cell, 1, variant='equispaced')
             vert_elt = FiniteElement('DG', interval, 1, variant='equispaced')
             V_elt = TensorProductElement(hori_elt, vert_elt)
             continuity = {'horizontal': False, 'vertical': False}
         else:
-            cell = self.mesh.ufl_cell().cellname()
+            cell = self.mesh.ufl_cell().cellname
             V_elt = FiniteElement('DG', cell, 1, variant='equispaced')
             continuity = False
 
@@ -285,9 +285,9 @@ def build_base_spaces(self, family, horizontal_degree, vertical_degree=None):
         """
 
         if self.extruded_mesh:
-            cell = self.mesh._base_mesh.ufl_cell().cellname()
+            cell = self.mesh._base_mesh.ufl_cell().cellname
         else:
-            cell = self.mesh.ufl_cell().cellname()
+            cell = self.mesh.ufl_cell().cellname
 
         hdiv_family = hcurl_hdiv_dict[family]
         hcurl_family = hdiv_hcurl_dict[family]
@@ -359,7 +359,7 @@ def build_compatible_spaces(self):
 
             return Vcg, Vcurl, Vu, Vdg, Vth
 
-        elif self.mesh.topological_dimension() > 1:
+        elif self.mesh.topological_dimension > 1:
             # 2D: two de Rham complexes (hcurl or hdiv) with 3 spaces
             # 3D: one de Rham complexes with 4 spaces
             # either way, build all spaces

gusto/core/meshes.py

@@ -618,12 +618,12 @@ def get_flat_latlon_mesh(mesh):
     coords_fs = coords_orig.function_space()
 
     if coords_fs.extruded:
-        cell = mesh._base_mesh.ufl_cell().cellname()
+        cell = mesh._base_mesh.ufl_cell().cellname
         DG1_hori_elt = FiniteElement("DG", cell, 1, variant="equispaced")
         DG1_vert_elt = FiniteElement("DG", interval, 1, variant="equispaced")
         DG1_elt = TensorProductElement(DG1_hori_elt, DG1_vert_elt)
     else:
-        cell = mesh.ufl_cell().cellname()
+        cell = mesh.ufl_cell().cellname
         DG1_elt = FiniteElement("DG", cell, 1, variant="equispaced")
     vec_DG1 = VectorFunctionSpace(mesh, DG1_elt)
     coords_dg = Function(vec_DG1).interpolate(coords_orig)

gusto/diagnostics/compressible_euler_diagnostics.py

@@ -55,7 +55,7 @@ def setup(self, domain, state_fields):
         gradu = state_fields("u_gradient")
 
         denom = 0.
-        z_dim = domain.mesh.geometric_dimension() - 1
+        z_dim = domain.mesh.geometric_dimension - 1
         u_dim = state_fields("u").ufl_shape[0]
         for i in range(u_dim-1):
             denom += gradu[i, z_dim]**2

gusto/diagnostics/diagnostics.py

@@ -492,7 +492,7 @@ def setup(self, domain, state_fields):
         """
         f = state_fields(self.fname)
 
-        mesh_dim = domain.mesh.geometric_dimension()
+        mesh_dim = domain.mesh.geometric_dimension
         try:
             field_dim = state_fields(self.fname).ufl_shape[0]
         except IndexError:
@@ -602,7 +602,7 @@ def setup(self, domain, state_fields):
             domain (:class:`Domain`): the model's domain object.
             state_fields (:class:`StateFields`): the model's field container.
         """
-        dim = domain.mesh.topological_dimension()
+        dim = domain.mesh.topological_dimension
         e_x = as_vector([Constant(1.0)]+[Constant(0.0)]*(dim-1))
         super().setup(domain, state_fields, e_x)
 
@@ -624,7 +624,7 @@ def setup(self, domain, state_fields):
         """
         assert domain.metadata['domain_type'] not in ['interval', 'vertical_slice'], \
             f'Y-component diagnostic cannot be used with domain {domain.metadata["domain_type"]}'
-        dim = domain.mesh.topological_dimension()
+        dim = domain.mesh.topological_dimension
         e_y = as_vector([Constant(0.0), Constant(1.0)]+[Constant(0.0)]*(dim-2))
         super().setup(domain, state_fields, e_y)
 
@@ -646,7 +646,7 @@ def setup(self, domain, state_fields):
         """
         assert domain.metadata['domain_type'] not in ['interval', 'plane'], \
             f'Z-component diagnostic cannot be used with domain {domain.metadata["domain_type"]}'
-        dim = domain.mesh.topological_dimension()
+        dim = domain.mesh.topological_dimension
         e_x = as_vector([Constant(0.0)]*(dim-1)+[Constant(1.0)])
         super().setup(domain, state_fields, e_x)
 
@@ -682,7 +682,7 @@ def _check_args(self, domain, field):
         """
 
         # check geometric dimension is 3D
-        if domain.mesh.geometric_dimension() != 3:
+        if domain.mesh.geometric_dimension != 3:
             raise ValueError('Spherical components only work when the geometric dimension is 3!')
 
         if np.prod(field.ufl_shape) != 3:
@@ -1043,7 +1043,7 @@ def setup(self, domain, state_fields):
                 horiz_degree = m_X_horiz.degree() + rho_d_horiz.degree()
                 vert_degree = m_X_vert.degree() + rho_d_vert.degree()
 
-                cell = domain.mesh._base_mesh.ufl_cell().cellname()
+                cell = domain.mesh._base_mesh.ufl_cell().cellname
                 horiz_elt = FiniteElement('DG', cell, horiz_degree)
                 vert_elt = FiniteElement('DG', cell, vert_degree)
                 elt = TensorProductElement(horiz_elt, vert_elt)
@@ -1052,7 +1052,7 @@ def setup(self, domain, state_fields):
                 rho_d_degree = rho_d_space.ufl_element().degree()
                 degree = m_X_degree + rho_d_degree
 
-                cell = domain.mesh.ufl_cell().cellname()
+                cell = domain.mesh.ufl_cell().cellname
                 elt = FiniteElement('DG', cell, degree)
 
             tracer_density_space = FunctionSpace(domain.mesh, elt, name='tracer_density_space')

gusto/equations/common_forms.py

@@ -257,7 +257,7 @@ def advection_equation_circulation_form(domain, test, q, ubar):
         class:`LabelledForm`: a labelled transport form.
     """
 
-    if domain.mesh.topological_dimension() == 3:
+    if domain.mesh.topological_dimension == 3:
         L = inner(test, cross(curl(q), ubar))*dx
 
     else:

gusto/recovery/recovery.py

@@ -124,12 +124,12 @@ def __init__(self, x_inout, method=BoundaryMethod.extruded, eff_coords=None):
         elif self.method == BoundaryMethod.taylor:
             # Create DG1 space ----------------------------------------------- #
             if V_inout.extruded:
-                cell = mesh._base_mesh.ufl_cell().cellname()
+                cell = mesh._base_mesh.ufl_cell().cellname
                 DG1_hori_elt = FiniteElement("DG", cell, 1, variant="equispaced")
                 DG1_vert_elt = FiniteElement("DG", interval, 1, variant="equispaced")
                 DG1_element = TensorProductElement(DG1_hori_elt, DG1_vert_elt)
             else:
-                cell = mesh.ufl_cell().cellname()
+                cell = mesh.ufl_cell().cellname
                 DG1_element = FiniteElement("DG", cell, 1, variant="equispaced")
 
             vec_DG1 = VectorFunctionSpace(mesh, DG1_element)
@@ -323,12 +323,12 @@ def find_eff_coords(V0):
 
     mesh = V0.mesh()
     if V0.extruded:
-        cell = mesh._base_mesh.ufl_cell().cellname()
+        cell = mesh._base_mesh.ufl_cell().cellname
         DG1_hori_elt = FiniteElement("DG", cell, 1, variant="equispaced")
         DG1_vert_elt = FiniteElement("DG", interval, 1, variant="equispaced")
         DG1_element = TensorProductElement(DG1_hori_elt, DG1_vert_elt)
     else:
-        cell = mesh.ufl_cell().cellname()
+        cell = mesh.ufl_cell().cellname
         DG1_element = FiniteElement("DG", cell, 1, variant="equispaced")
 
     vec_CG1 = VectorFunctionSpace(mesh, "CG", 1)
@@ -397,12 +397,12 @@ def correct_eff_coords(eff_coords):
     vec_CG1 = VectorFunctionSpace(mesh, "CG", 1)
 
     if vec_CG1.extruded:
-        cell = mesh._base_mesh.ufl_cell().cellname()
+        cell = mesh._base_mesh.ufl_cell().cellname
         DG1_hori_elt = FiniteElement("DG", cell, 1, variant="equispaced")
         DG1_vert_elt = FiniteElement("DG", interval, 1, variant="equispaced")
         DG1_element = TensorProductElement(DG1_hori_elt, DG1_vert_elt)
     else:
-        cell = mesh.ufl_cell().cellname()
+        cell = mesh.ufl_cell().cellname
         DG1_element = FiniteElement("DG", cell, 1, variant="equispaced")
 
     vec_DG1 = VectorFunctionSpace(mesh, DG1_element)

gusto/recovery/recovery_kernels.py

@@ -268,7 +268,7 @@ def __init__(self, DG1):
             DG1 (:class:`FunctionSpace`): The equispaced DG1 function space.
         """
         shapes = {"nDOFs": DG1.finat_element.space_dimension(),
-                  "dim": DG1.mesh().topological_dimension()}
+                  "dim": DG1.mesh().topological_dimension}
 
         # EFF_COORDS are the effective coordinates
         # ACT_COORDS are the actual coordinates

gusto/recovery/recovery_options.py

@@ -51,7 +51,7 @@ def __init__(self, domain, boundary_method=None, use_vector_spaces=False):
                 theta_boundary_method = boundary_method['theta']
             else:
                 theta_boundary_method = None
-            cell = mesh._base_mesh.ufl_cell().cellname()
+            cell = mesh._base_mesh.ufl_cell().cellname
             DG_hori_ele = FiniteElement('DG', cell, 1, variant='equispaced')
             DG_vert_ele = FiniteElement('DG', interval, (domain.vertical_degree + 1), variant='equispaced')
             CG_hori_ele = FiniteElement('CG', cell, 1)
@@ -66,7 +66,7 @@ def __init__(self, domain, boundary_method=None, use_vector_spaces=False):
                                                  recovered_space=VCG_theta,
                                                  boundary_method=theta_boundary_method)
         else:
-            cell = self.mesh.ufl_cell().cellname()
+            cell = self.mesh.ufl_cell().cellname
 
         # ----------------------------------------------------------------------
         # Building the DG options

gusto/spatial_methods/augmentation.py

@@ -106,7 +106,7 @@ def __init__(
         sign_u = 0.5*(sign(dot(u, n)) + 1)
         upw = lambda f: (sign_u('+')*f('+') + sign_u('-')*f('-'))
 
-        if domain.mesh.topological_dimension() == 2:
+        if domain.mesh.topological_dimension == 2:
             mix_test = test_F - domain.perp(grad(test_Z))
             F_cross_u = Z*domain.perp(u)
         elif domain.mesh.topological_dimension == 3:
@@ -151,14 +151,14 @@ def __init__(
             DG0 = FunctionSpace(domain.mesh, 'DG', 0)
             ones = Function(DG0).interpolate(Constant(1.0))
             area = assemble(ones*dx)
-            mean_dx = (area/DG0.dof_count)**(1/domain.mesh.geometric_dimension())
+            mean_dx = (area/DG0.dof_count)**(1/domain.mesh.geometric_dimension)
 
             # Divide by approximately (1 + c)
             tau /= (1.0 + sqrt(dot(u, u))*domain.dt/Constant(mean_dx))
 
             dxqp = dx(degree=3)
 
-            if domain.mesh.topological_dimension() == 2:
+            if domain.mesh.topological_dimension == 2:
                 time_deriv_form -= inner(mix_test, tau*Z*domain.perp(u)/domain.dt)*dxqp
                 transport_form -= inner(
                     mix_test, tau*domain.perp(u)*domain.divperp(Z*domain.perp(u))
@@ -168,7 +168,7 @@ def __init__(
                         mix_test,
                         tau*domain.perp(u)*domain.divperp(u_dot_nabla_F)
                     )*dxqp
-            elif domain.mesh.topological_dimension() == 3:
+            elif domain.mesh.topological_dimension == 3:
                 time_deriv_form -= inner(mix_test, tau*cross(Z, u)/domain.dt)*dxqp
                 transport_form -= inner(
                     mix_test, tau*cross(curl(Z*u), u)

gusto/spatial_methods/diffusion_methods.py

@@ -143,7 +143,7 @@ def __init__(self, equation, variable, diffusion_parameters):
 
         super().__init__(equation, variable)
 
-        if equation.domain.mesh.topological_dimension() == 1:
+        if equation.domain.mesh.topological_dimension == 1:
             self.form = interior_penalty_diffusion_form_1d(
                 equation.domain, self.test, self.field, diffusion_parameters)
         else:
@@ -157,7 +157,7 @@ def __init__(self, equation, variable, diffusion_parameters):
 
         super().__init__(equation, variable)
 
-        if equation.domain.mesh.topological_dimension() == 1:
+        if equation.domain.mesh.topological_dimension == 1:
             kappa = diffusion_parameters.kappa
             self.form = diffusion(kappa * self.test.dx(0) * self.field.dx(0) * dx)
         else:

gusto/spatial_methods/limiters.py

@@ -52,12 +52,12 @@ def __init__(self, space, subspace=None):
 
         # Create equispaced DG1 space needed for limiting
         if space.extruded:
-            cell = mesh._base_mesh.ufl_cell().cellname()
+            cell = mesh._base_mesh.ufl_cell().cellname
             DG1_hori_elt = FiniteElement("DG", cell, 1, variant="equispaced")
             DG1_vert_elt = FiniteElement("DG", interval, 1, variant="equispaced")
             DG1_element = TensorProductElement(DG1_hori_elt, DG1_vert_elt)
         else:
-            cell = mesh.ufl_cell().cellname()
+            cell = mesh.ufl_cell().cellname
             DG1_element = FiniteElement("DG", cell, 1, variant="equispaced")
 
         DG1_equispaced = FunctionSpace(mesh, DG1_element)
@@ -123,7 +123,7 @@ def __init__(self, space):
         self.Vt_brok = FunctionSpace(mesh, BrokenElement(space.ufl_element()))
 
         # Create equispaced DG1 space needed for limiting
-        cell = mesh._base_mesh.ufl_cell().cellname()
+        cell = mesh._base_mesh.ufl_cell().cellname
         DG1_hori_elt = FiniteElement("DG", cell, 1, variant="equispaced")
         DG1_vert_elt = FiniteElement("DG", interval, 1, variant="equispaced")
         CG2_vert_elt = FiniteElement("CG", interval, 2)

gusto/spatial_methods/transport_methods.py

@@ -239,7 +239,7 @@ def __init__(self, equation, variable, ibp=IntegrateByParts.ONCE,
         # Determine appropriate form to use
         # -------------------------------------------------------------------- #
         # first check for 1d mesh and scalar velocity space
-        if equation.domain.mesh.topological_dimension() == 1 and len(equation.domain.spaces("HDiv").shape) == 0:
+        if equation.domain.mesh.topological_dimension == 1 and len(equation.domain.spaces("HDiv").shape) == 0:
             assert not vector_manifold_correction
             if self.transport_equation_type == TransportEquationType.advective:
                 form = upwind_advection_form_1d(
@@ -350,7 +350,7 @@ def __init__(self, equation, variable, ibp=IntegrateByParts.ONCE,
         # Determine appropriate form to use
         # -------------------------------------------------------------------- #
         # first check for 1d mesh and scalar velocity space
-        if equation.domain.mesh.topological_dimension() == 1 and len(equation.domain.spaces("HDiv").shape) == 0:
+        if equation.domain.mesh.topological_dimension == 1 and len(equation.domain.spaces("HDiv").shape) == 0:
             assert not vector_manifold_correction
             raise ValueError('You cannot do horizontal and vertical splitting in 1D')
         else:
@@ -831,7 +831,7 @@ def upwind_circulation_form(domain, test, q, ibp=IntegrateByParts.ONCE):
     n = FacetNormal(domain.mesh)
     Upwind = 0.5*(sign(dot(ubar, n))+1)
 
-    if domain.mesh.topological_dimension() == 3:
+    if domain.mesh.topological_dimension == 3:
 
         if ibp != IntegrateByParts.ONCE:
             raise NotImplementedError

gusto/time_discretisation/wrappers.py

@@ -337,7 +337,7 @@ def setup(self, field_name):
         # -------------------------------------------------------------------- #
 
         # construct tau, if it is not specified
-        dim = domain.mesh.topological_dimension()
+        dim = domain.mesh.topological_dimension
         if self.options.tau is not None:
             # if tau is provided, check that is has the right size
             self.tau = self.options.tau

integration-tests/transport/test_dg_transport.py

@@ -53,7 +53,7 @@ def test_dg_transport_scalar(tmpdir, geometry, equation_form, tracer_setup):
 def test_dg_transport_vector(tmpdir, geometry, equation_form, tracer_setup):
     setup = tracer_setup(tmpdir, geometry)
     domain = setup.domain
-    gdim = domain.mesh.geometric_dimension()
+    gdim = domain.mesh.geometric_dimension
     f_init = as_vector([setup.f_init]*gdim)
     V = VectorFunctionSpace(domain.mesh, "DG", 1)
     if equation_form == "advective":