December 2015 - firedrake - imperial.ac.uk

Re: [firedrake] UFL errors
by Colin Cotter 21 Dec '15

21 Dec '15

Hi Anna, I'm not sure what's wrong, but it's a bit overcomplicated to write abs(grad(phi))**2, better to write inner(grad(phi),grad(phi)) all the best --cjc On 21 December 2015 at 10:51, Anna Kalogirou <A.Kalogirou(a)leeds.ac.uk> wrote: > Dear all, > > I get some UFL errors when trying to assemble the expression > > E = 0.5*rho*assemble( (Hb*abs(grad(phi0))**2 + g*eta0**2)*dx ) > > where phi0, eta0, Hb are functions and the rest are scalars. > > At first, I got this: > ----------------- > E0_w = 0.5*rho*assemble( (Hb*abs(grad(phi0))**2 + g*eta0**2)*dx ) # > Energy of the water > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble.py", > line 70, in assemble > return assemble_expressions.assemble_expression(f) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 592, in assemble_expression > evaluate_expression(Assign(result, expr), subset) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 577, in evaluate_expression > k = expression_kernel(e, args) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 524, in expression_kernel > ast_expr = _ast(expr) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 44, in _ast > return expr.ast > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 128, in ast > return self._ast(_ast(self.ufl_operands[0]), > _ast(self.ufl_operands[1])) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 48, in _ast > return f(expr) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 603, in <lambda> > ufl.algebra.Product: (lambda e: ast.Par(_to_prod(e.ufl_operands))), > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 597, in <lambda> > _to_prod = lambda o: ast.Prod(_ast(o[0]), _to_sum(o[1:])) if len(o) > > 1 else _ast(o[0]) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 596, in <lambda> > _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 > else _ast(o[0]) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 48, in _ast > return f(expr) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 602, in <lambda> > ufl.algebra.Sum: (lambda e: ast.Par(_to_sum(e.ufl_operands))), > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 596, in <lambda> > _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 > else _ast(o[0]) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 596, in <lambda> > _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 > else _ast(o[0]) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 48, in _ast > return f(expr) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 603, in <lambda> > ufl.algebra.Product: (lambda e: ast.Par(_to_prod(e.ufl_operands))), > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 597, in <lambda> > _to_prod = lambda o: ast.Prod(_ast(o[0]), _to_sum(o[1:])) if len(o) > > 1 else _ast(o[0]) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 596, in <lambda> > _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 > else _ast(o[0]) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", > line 49, in _ast > raise TypeError("No ast handler for %s" % str(type(expr))) > TypeError: No ast handler for <class 'ufl.tensoralgebra.Inner’> > ----------------- > > and then I tried to write instead > > E = 0.5*rho*assemble( (Hb*abs(grad(phi0))*abs(grad(phi0)) + > g*eta0*eta0)*dx ) > > and got a different error: > ----------------- > E0_w = 0.5*rho*assemble( (Hb*abs(grad(phi0))*abs(grad(phi0)) + > g*eta0*eta0)*dx ) # Energy of the water > File > "/Users/matak/firedrake/lib/python2.7/site-packages/ufl/exproperators.py", > line 166, in _mul > return _mult(self, o) > File > "/Users/matak/firedrake/lib/python2.7/site-packages/ufl/exproperators.py", > line 144, in _mult > error("Invalid ranks {0} and {1} in product.".format(r1, r2)) > File "/Users/matak/firedrake/lib/python2.7/site-packages/ufl/log.py", > line 158, in error > raise self._exception_type(self._format_raw(*message)) > ufl.log.UFLException: Invalid ranks 1 and 1 in product. > ----------------- > > I don’t understand what’s going wrong. It was working on Friday. > > Thank you. > > Best, Anna. > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

1 0

UFL errors
by Anna Kalogirou 21 Dec '15

21 Dec '15

Dear all, I get some UFL errors when trying to assemble the expression E = 0.5*rho*assemble( (Hb*abs(grad(phi0))**2 + g*eta0**2)*dx ) where phi0, eta0, Hb are functions and the rest are scalars. At first, I got this: ----------------- E0_w = 0.5*rho*assemble( (Hb*abs(grad(phi0))**2 + g*eta0**2)*dx ) # Energy of the water File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble.py", line 70, in assemble return assemble_expressions.assemble_expression(f) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 592, in assemble_expression evaluate_expression(Assign(result, expr), subset) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 577, in evaluate_expression k = expression_kernel(e, args) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 524, in expression_kernel ast_expr = _ast(expr) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 44, in _ast return expr.ast File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 128, in ast return self._ast(_ast(self.ufl_operands[0]), _ast(self.ufl_operands[1])) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 48, in _ast return f(expr) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 603, in <lambda> ufl.algebra.Product: (lambda e: ast.Par(_to_prod(e.ufl_operands))), File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 597, in <lambda> _to_prod = lambda o: ast.Prod(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 else _ast(o[0]) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 596, in <lambda> _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 else _ast(o[0]) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 48, in _ast return f(expr) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 602, in <lambda> ufl.algebra.Sum: (lambda e: ast.Par(_to_sum(e.ufl_operands))), File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 596, in <lambda> _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 else _ast(o[0]) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 596, in <lambda> _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 else _ast(o[0]) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 48, in _ast return f(expr) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 603, in <lambda> ufl.algebra.Product: (lambda e: ast.Par(_to_prod(e.ufl_operands))), File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 597, in <lambda> _to_prod = lambda o: ast.Prod(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 else _ast(o[0]) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 596, in <lambda> _to_sum = lambda o: ast.Sum(_ast(o[0]), _to_sum(o[1:])) if len(o) > 1 else _ast(o[0]) File "/Users/matak/firedrake/lib/python2.7/site-packages/firedrake/assemble_expressions.py", line 49, in _ast raise TypeError("No ast handler for %s" % str(type(expr))) TypeError: No ast handler for <class 'ufl.tensoralgebra.Inner’> ----------------- and then I tried to write instead E = 0.5*rho*assemble( (Hb*abs(grad(phi0))*abs(grad(phi0)) + g*eta0*eta0)*dx ) and got a different error: ----------------- E0_w = 0.5*rho*assemble( (Hb*abs(grad(phi0))*abs(grad(phi0)) + g*eta0*eta0)*dx ) # Energy of the water File "/Users/matak/firedrake/lib/python2.7/site-packages/ufl/exproperators.py", line 166, in _mul return _mult(self, o) File "/Users/matak/firedrake/lib/python2.7/site-packages/ufl/exproperators.py", line 144, in _mult error("Invalid ranks {0} and {1} in product.".format(r1, r2)) File "/Users/matak/firedrake/lib/python2.7/site-packages/ufl/log.py", line 158, in error raise self._exception_type(self._format_raw(*message)) ufl.log.UFLException: Invalid ranks 1 and 1 in product. ----------------- I don’t understand what’s going wrong. It was working on Friday. Thank you. Best, Anna.

1 0

gmsh: Accessing PhysicalNames for boundary conditions
by Buesing, Henrik 18 Dec '15

18 Dec '15

Dear all, I create a very simple gmsh file: just a square with "left", "bottom", "right" and "top" boundary. I have attached two versions of the mesh. 1) neu3.msh: This has the "PhysicalNames" declared at the beginning. I get a "File is not a valid Gmsh file" error. 2) neu4.msh: This has the "PhysicalNames" declared at the end. But if I now try to set Dirichlet bc for the "left" and "right" boundary, e.g. with bc0 = DirichletBC(W.sub(0), Expression("0.0"), "right") nothing seems to happen. I guess I'm missing sth. obvious here... Thank you! Henrik Error message for neu3.msh Traceback (most recent call last): File "charge.py", line 7, in <module> mymesh = Mesh("neu3.msh") File "<decorator-gen-274>", line 2, in __init__ File "/home/hb111949/.local/lib/python2.7/site-packages/PyOP2-0.11.0_425_g64db689-py2.7-linux-x86_64.egg/pyop2/profiling.py", line 203, in wrapper return f(*args, **kwargs) File "/home/hb111949/Code/firedrake/firedrake/mesh.py", line 351, in __init__ plex = _from_gmsh(meshfile) File "/home/hb111949/Code/firedrake/firedrake/mesh.py", line 152, in _from_gmsh gmsh_plex = PETSc.DMPlex().createGmsh(gmsh_viewer) File "PETSc/DMPlex.pyx", line 115, in petsc4py.PETSc.DMPlex.createGmsh (src/petsc4py.PETSc.c:214898) petsc4py.PETSc.Error: error code 62 [0] DMPlexCreateGmsh() line 140 in /rwthfs/rz/cluster/home/hb111949/Code/petsc/src/dm/impls/plex/plexgmsh.c [0] Invalid argument [0] File is not a valid Gmsh file -- Dipl.-Math. Henrik Büsing Institute for Applied Geophysics and Geothermal Energy E.ON Energy Research Center RWTH Aachen University ------------------------------------------------------ Mathieustr. 10 | Tel +49 (0)241 80 49907 52074 Aachen, Germany | Fax +49 (0)241 80 49889 ------------------------------------------------------ http://www.eonerc.rwth-aachen.de/GGE hbuesing(a)eonerc.rwth-aachen.de ------------------------------------------------------

5 11

PyOP2 error
by Justin Chang 17 Dec '15

17 Dec '15

Hi all, Just made an update to firedrake, and my code breaks down with this PyOP2 error: Traceback (most recent call last): File "2D_OS_boundary_ex3.py", line 286, in <module> op2.par_loop(reactions3_kernel,Q.dof_dset,u0_A.dat(op2.READ),u0_B.dat(op2.READ),c_A.dat(op2.WRITE),c_B.dat(op2.WRITE),c_C.dat(op2.WRITE),k1.dat(op2.READ),ctx(op2.READ)) File "<decorator-gen-272>", line 2, in par_loop File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/versioning.py", line 154, in modifies_arguments retval = func(*args, **kwargs) File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/op2.py", line 278, in par_loop return backends._BackendSelector._backend.par_loop(kernel, iterset, *args, **kwargs) File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/base.py", line 4422, in par_loop return _make_object('ParLoop', kernel, it_space, *args, **kwargs).enqueue() File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/backends.py", line 72, in _make_object return _BackendSelector._backend.__dict__[obj](*args, **kwargs) File "<decorator-gen-29>", line 2, in __init__ File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/utils.py", line 129, in wrapper self.check_args(args, kwargs) File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/utils.py", line 155, in check_args self.check_arg(arg, argcond, exception) File "/Users/justin/Software/firedrake/src/PyOP2/pyop2/utils.py", line 170, in check_arg % (self.file, self.line, arg, argtype)) pyop2.exceptions.SetTypeError: /Users/justin/Software/firedrake/src/PyOP2/pyop2/base.py:4008 Parameter OP2 DataSet: None_nodes_dset on set OP2 Set: None_nodes with size 20301, with dim (1,) must be of type <class 'pyop2.base.Set'> It worked previously, I am guessing the latest change to PyOP2 changed a few things here and there. Attached is the code and the PyOP2 kernel. Thanks, Justin

2 5

Gmsh help
by Anna Kalogirou 17 Dec '15

17 Dec '15

Dear all, I need to import my own mesh which was not created using Gmsh, so I am trying to convert it to Gmsh format. I have some questions and I didn't find any documentation on this on the Gmsh website. For example, in the sample msh file attached there should be 40 elements and not 72. I don't understand what the first 32 lines (i.e. lines 64-95 in the file) under "Elements" actually are. Also, in the actual elements (lines 96-135), what is the first quadruplet {3 2 0 5} after the global element number? Thanks, Anna.

3 9

Icosahedral mesh
by Angwenyi David 17 Dec '15

17 Dec '15

Dear all, I have had an error while trying to solve the Helmholtz eqn on Icosahedral grids. The code used is: >>> from firedrake import * >>> m = IcosahedralSphereMesh(radius=3, refinement_level=3) >>> mesh = ExtrudedMesh(m, 5, layer_height=0.1, extrusion_type='radial') >>> horiz_elt = FiniteElement("CG", triangle, 1) >>> vert_elt = FiniteElement("CG", interval, 1) >>> elt = OuterProductElement(horiz_elt, vert_elt) >>> V = FunctionSpace(mesh, elt) >>> u = TrialFunction(V) >>> v = TestFunction(V) >>> f = Function(V) >>> f.interpolate(Expression("(2*x[0])+(sin(3*x[1]))+(cos(0.5*x[2]))")) Coefficient(FunctionSpace(Mesh(OuterProductVectorElement(OuterProductElement(FiniteElement('Lagrange', triangle, 1), FiniteElement('Lagrange', interval, 1), OuterProductCell(triangle, interval)), dim=3), 4), OuterProductElement(FiniteElement('Lagrange', triangle, 1), FiniteElement('Lagrange', interval, 1), OuterProductCell(triangle, interval))), 3) >>> a = (dot(grad(v), grad(u)) + v * u) * dx >>> L = f * v * dx >>> solve(a == L, u, solver_parameters={'ksp_type': 'cg'}) Traceback (most recent call last): File "<stdin>", line 1, in <module> File "/Users/master/firedrake/lib/python2.7/site-packages/firedrake/solving.py", line 120, in solve _solve_varproblem(*args, **kwargs) File "/Users/master/firedrake/lib/python2.7/site-packages/firedrake/solving.py", line 132, in _solve_varproblem nest = _extract_args(*args, **kwargs) File "/Users/master/firedrake/lib/python2.7/site-packages/firedrake/solving.py", line 259, in _extract_args u = _extract_u(args[1]) File "/Users/master/firedrake/lib/python2.7/site-packages/firedrake/solving.py", line 299, in _extract_u raise RuntimeError("Expecting second argument to be a Coefficient") RuntimeError: Expecting second argument to be a Coefficient I will be glad to get help on how to fix this, and any other advice on implementation of icosahedral grids. Receive thanks in advance. David

2 4

mesh from gmsh
by Buesing, Henrik 14 Dec '15

14 Dec '15

Dear all, I tried to create a simple mesh in gmsh with a hole in it and would like to import this into firedrake. When I try to use this mesh with mymesh = Mesh("test1.msh") elt = FiniteElement("CG", "triangle", 1) CG1 = FunctionSpace(mymesh, elt) I get the error below. Any ideas why? I attached the mesh. Thank you! Henrik Traceback (most recent call last): File "charge.py", line 25, in <module> CG1 = FunctionSpace(mymesh, elt) File "/home/hb111949/.local/lib/python2.7/site-packages/PyOP2-0.11.0_425_g64db689-py2.7-linux-x86_64.egg/pyop2/caching.py", line 161, in __new__ obj = make_obj() File "/home/hb111949/.local/lib/python2.7/site-packages/PyOP2-0.11.0_425_g64db689-py2.7-linux-x86_64.egg/pyop2/caching.py", line 142, in make_obj obj.__init__(*args, **kwargs) File "/home/hb111949/Code/firedrake/firedrake/functionspace.py", line 578, in __init__ mesh.init() File "/home/hb111949/Code/firedrake/firedrake/mesh.py", line 480, in init self._callback(self) File "/home/hb111949/Code/firedrake/firedrake/mesh.py", line 404, in callback old_to_new = self._plex.getOrdering(PETSc.Mat.OrderingType.RCM).indices File "PETSc/DMPlex.pyx", line 601, in petsc4py.PETSc.DMPlex.getOrdering (src/petsc4py.PETSc.c:224945) petsc4py.PETSc.Error: error code 77 [0] DMPlexGetOrdering() line 136 in /rwthfs/rz/cluster/home/hb111949/Code/petsc/src/dm/impls/plex/plexreorder.c [0] DMPlexCreateOrderingClosure_Static() line 43 in /rwthfs/rz/cluster/home/hb111949/Code/petsc/src/dm/impls/plex/plexreorder.c [0] Petsc has generated inconsistent data [0] Number of depth 1 faces 40 does not match permuted nubmer 38 -- Dipl.-Math. Henrik Büsing Institute for Applied Geophysics and Geothermal Energy E.ON Energy Research Center RWTH Aachen University ------------------------------------------------------ Mathieustr. 10 | Tel +49 (0)241 80 49907 52074 Aachen, Germany | Fax +49 (0)241 80 49889 ------------------------------------------------------ http://www.eonerc.rwth-aachen.de/GGE hbuesing(a)eonerc.rwth-aachen.de ------------------------------------------------------

2 4

Periodic boundary conditions
by Buesing, Henrik 14 Dec '15

14 Dec '15

Dear all, I tried to build a mixed space with periodic bc on two of the subspaces: meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", quadrilateral, 0) DG1 = FunctionSpace(mesh, horiz_elt) meshbase = PeriodicRectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", quadrilateral, 0) DG2 = FunctionSpace(mesh, horiz_elt) DG3 = FunctionSpace(mesh, horiz_elt) W = DG1 * DG2 * DG3 but got an error message (see below). Is there a way to set periodic bc on two of the three subspaces? Henrik Traceback (most recent call last): File "twophase.py", line 43, in <module> W = DG1 * DG2 * DG3 File "/home/hb111949/Code/firedrake/firedrake/functionspace.py", line 536, in __mul__ return MixedFunctionSpace((self, other)) File "/home/hb111949/.local/lib/python2.7/site-packages/PyOP2-0.11.0_425_g64db689-py2.7-linux-x86_64.egg/pyop2/caching.py", line 161, in __new__ obj = make_obj() File "/home/hb111949/.local/lib/python2.7/site-packages/PyOP2-0.11.0_425_g64db689-py2.7-linux-x86_64.egg/pyop2/caching.py", line 142, in make_obj obj.__init__(*args, **kwargs) File "/home/hb111949/Code/firedrake/firedrake/functionspace.py", line 742, in __init__ self._ufl_element = ufl.MixedElement(*[fs.ufl_element() for fs in self._spaces]) File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/finiteelement/mixedelement.py", line 57, in __init__ domains = tuple(sorted(set(chain(*[element.domains() for element in elements])))) File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/domain.py", line 216, in __lt__ return self.hash_data() < other.hash_data() File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/exproperators.py", line 55, in _lt return LT(left, right) File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/conditional.py", line 144, in __init__ BinaryCondition.__init__(self, "<", left, right) File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/conditional.py", line 72, in __init__ "Expecting scalar arguments.") File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/assertions.py", line 37, in ufl_assert if not condition: error(*message) File "/home/hb111949/.local/lib/python2.7/site-packages/ufl/log.py", line 151, in error raise self._exception_type(self._format_raw(*message)) ufl.log.UFLException: Expecting scalar arguments. -- Dipl.-Math. Henrik Büsing Institute for Applied Geophysics and Geothermal Energy E.ON Energy Research Center RWTH Aachen University ------------------------------------------------------ Mathieustr. 10 | Tel +49 (0)241 80 49907 52074 Aachen, Germany | Fax +49 (0)241 80 49889 ------------------------------------------------------ http://www.eonerc.rwth-aachen.de/GGE hbuesing(a)eonerc.rwth-aachen.de ------------------------------------------------------

2 4

Re: [firedrake] matrix/function/mesh entries numbering
by David Ham 13 Dec '15

13 Dec '15

> > Hi Tomasz, > > > On Mon, 14 Dec 2015 at 04:46 Tomasz Salwa [RPG] <mmtjs(a)leeds.ac.uk> wrote: > > > > > > > Dear Firedrakers, > > > > > I use utility meshes, including the extruded ones. I need to access the > entries of system matrices, but I am not sure about the node/entries > numbering. Say, I have the code: > > > m = UnitIntervalMesh(3) > mesh = ExtrudedMesh(m, 3) > V = FunctionSpace(mesh, "CG", 1) > > > > f = Function(V) > trial = TrialFunction(V) > test = TestFunction(V) > A = assemble( trial * test * dx ).M.values > > > > > > Does the element A[i,j] correspond to nodes mesh.coordinates.dat.data[i] & > [j], same > > > as f.dat.data[i] & j? > Yes. This is exactly what happens. Note, however that coordinates is vector-valued so mesh.coordinates.dat.data[i] is itself a vector. What is the case for VectorFunctionSpace? > For a matrix assembled over a VectorFunctionSpace, the entries are interleaved. For example, if the mesh were 2D then the row corresponding to f.dat.data[i,0] would be row 2*i and that corresponding to f.dat.data[i,1] would be 2*i + 1. Regards, David > > > > > Thanks, > > Tomasz >

1 0

matrix/function/mesh entries numbering
by Tomasz Salwa [RPG] 13 Dec '15

13 Dec '15

Dear Firedrakers, I use utility meshes, including the extruded ones. I need to access the entries of system matrices, but I am not sure about the node/entries numbering. Say, I have the code: m = UnitIntervalMesh(3) mesh = ExtrudedMesh(m, 3) V = FunctionSpace(mesh, "CG", 1) f = Function(V) trial = TrialFunction(V) test = TestFunction(V) A = assemble( trial * test * dx ).M.values Does the element A[i,j] correspond to nodes mesh.coordinates.dat.data[i] & [j], same as f.dat.data[i] & j? What is the case for VectorFunctionSpace? Thanks, Tomasz

1 0