def FunctionSpace(mesh, family, degree=None, name=None, vfamily=None, vdegree=None): """Create a :class:`.FunctionSpace`. :arg mesh: The mesh to determine the cell from. :arg family: The finite element family. :arg degree: The degree of the finite element. :arg name: An optional name for the function space. :arg vfamily: The finite element in the vertical dimension (extruded meshes only). :arg vdegree: The degree of the element in the vertical dimension (extruded meshes only). The ``family`` argument may be an existing :class:`ufl.FiniteElementBase`, in which case all other arguments are ignored and the appropriate :class:`.FunctionSpace` is returned. """ element = make_scalar_element(mesh, family, degree, vfamily, vdegree) # Support FunctionSpace(mesh, MixedElement) if type(element) is ufl.MixedElement: return MixedFunctionSpace(element, mesh=mesh, name=name) # Check that any Vector/Tensor/Mixed modifiers are outermost. check_element(element) # Otherwise, build the FunctionSpace. topology = mesh.topology new = impl.FunctionSpace(topology, element, name=name) if mesh is not topology: return impl.WithGeometry(new, mesh) else: return new
def _coordinates_function(self): """The :class:`.Function` containing the coordinates of this mesh.""" import firedrake.functionspaceimpl as functionspaceimpl import firedrake.function as function self.init() coordinates_fs = self._coordinates.function_space() V = functionspaceimpl.WithGeometry(coordinates_fs, self) f = function.Function(V, val=self._coordinates) return f
def MixedFunctionSpace(spaces, name=None, mesh=None): """Create a :class:`.MixedFunctionSpace`. :arg spaces: An iterable of constituent spaces, or a :class:`~ufl.classes.MixedElement`. :arg name: An optional name for the mixed function space. :arg mesh: An optional mesh. Must be provided if spaces is a :class:`~ufl.classes.MixedElement`, ignored otherwise. """ if isinstance(spaces, ufl.FiniteElementBase): # Build the spaces if we got a mixed element assert type(spaces) is ufl.MixedElement and mesh is not None sub_elements = [] def rec(eles): for ele in eles: # Only want to recurse into MixedElements if type(ele) is ufl.MixedElement: rec(ele.sub_elements()) else: sub_elements.append(ele) rec(spaces.sub_elements()) spaces = [FunctionSpace(mesh, element) for element in sub_elements] # Check that function spaces are on the same mesh meshes = [space.mesh() for space in spaces] for i in range(1, len(meshes)): if meshes[i] is not meshes[0]: raise ValueError( "All function spaces must be defined on the same mesh!") # Select mesh mesh = meshes[0] # Get topological spaces spaces = tuple(s.topological for s in flatten(spaces)) # Error checking for space in spaces: if type(space) in (impl.FunctionSpace, impl.RealFunctionSpace): continue elif type(space) is impl.ProxyFunctionSpace: if space.component is not None: raise ValueError("Can't make mixed space with %s" % space) continue else: raise ValueError("Can't make mixed space with %s" % type(space)) new = impl.MixedFunctionSpace(spaces, name=name) if mesh is not mesh.topology: return impl.WithGeometry(new, mesh) return new
def setup_dump(self, t, tmax, pickup=False): """ Setup dump files Check for existence of directory so as not to overwrite output files Setup checkpoint file :arg tmax: model stop time :arg pickup: recover state from the checkpointing file if true, otherwise dump and checkpoint to disk. (default is False). """ if any([ self.output.dump_vtus, self.output.dumplist_latlon, self.output.dump_diagnostics, self.output.point_data, self.output.checkpoint and not pickup ]): # setup output directory and check that it does not already exist self.dumpdir = path.join("results", self.output.dirname) running_tests = '--running-tests' in sys.argv or "pytest" in self.output.dirname if self.mesh.comm.rank == 0: if not running_tests and path.exists( self.dumpdir) and not pickup: raise IOError("results directory '%s' already exists" % self.dumpdir) else: if not running_tests: makedirs(self.dumpdir) if self.output.dump_vtus: # setup pvd output file outfile = path.join(self.dumpdir, "field_output.pvd") self.dumpfile = File(outfile, project_output=self.output.project_fields, comm=self.mesh.comm) # make list of fields to dump self.to_dump = [field for field in self.fields if field.dump] # make dump counter self.dumpcount = itertools.count() # if there are fields to be dumped in latlon coordinates, # setup the latlon coordinate mesh and make output file if len(self.output.dumplist_latlon) > 0: mesh_ll = get_latlon_mesh(self.mesh) outfile_ll = path.join(self.dumpdir, "field_output_latlon.pvd") self.dumpfile_ll = File(outfile_ll, project_output=self.output.project_fields, comm=self.mesh.comm) # make functions on latlon mesh, as specified by dumplist_latlon self.to_dump_latlon = [] for name in self.output.dumplist_latlon: f = self.fields(name) field = Function(functionspaceimpl.WithGeometry( f.function_space(), mesh_ll), val=f.topological, name=name + '_ll') self.to_dump_latlon.append(field) # we create new netcdf files to write to, unless pickup=True, in # which case we just need the filenames if self.output.dump_diagnostics: diagnostics_filename = self.dumpdir + "/diagnostics.nc" self.diagnostic_output = DiagnosticsOutput(diagnostics_filename, self.diagnostics, self.output.dirname, self.mesh.comm, create=not pickup) if len(self.output.point_data) > 0: pointdata_filename = self.dumpdir + "/point_data.nc" ndt = int(tmax / self.timestepping.dt) self.pointdata_output = PointDataOutput(pointdata_filename, ndt, self.output.point_data, self.output.dirname, self.fields, self.mesh.comm, create=not pickup) # if we want to checkpoint and are not picking up from a previous # checkpoint file, setup the dumb checkpointing if self.output.checkpoint and not pickup: self.chkpt = DumbCheckpoint(path.join(self.dumpdir, "chkpt"), mode=FILE_CREATE) # make list of fields to pickup (this doesn't include # diagnostic fields) self.to_pickup = [field for field in self.fields if field.pickup] # if we want to checkpoint then make a checkpoint counter if self.output.checkpoint: self.chkptcount = itertools.count() # dump initial fields self.dump(t)