Hi Lawrence,

thanks, it does indeed work if I use smaller grids… It’s one of my test scripts, so I really only want to try it on small function spaces.

Eike

On 5 Jan 2015, at 13:06, Lawrence Mitchell <lawrence.mitchell@imperial.ac.uk> wrote:


On 2 Jan 2015, at 14:19, Eike Mueller <eike.h.mueller@googlemail.com> wrote:

Dear firedrakers,

I get strange error messages when I try to access the values of an assembled matrix in 3d (but it works in 2d).

The code I use is shown below. If I set dimension=2, higherorder=True (or False), it does indeed print out the values of the mass matrix. But for dimension=3 I get the following errors:

*** (1) higherorder=True ***

eikemueller@Eikes-MacBook-Pro $ python mat_assembly.py 
pyop2:WARNING No maximum assembly cache size. Install psutil or risk leaking memory!
Traceback (most recent call last):
File "mat_assembly.py", line 47, in <module>
  print mat.M.values
File "<string>", line 2, in values
File "/Users/eikemueller/PostDocBath/EllipticSolvers/PyOP2/pyop2/versioning.py", line 111, in modifies
  retval = method(self, *args, **kwargs)
File "/Users/eikemueller/PostDocBath/EllipticSolvers/PyOP2/pyop2/petsc_base.py", line 434, in values
  return self.handle[:, :]
File "Mat.pyx", line 218, in petsc4py.PETSc.Mat.__getitem__ (src/petsc4py.PETSc.c:99054)
File "petscmat.pxi", line 962, in petsc4py.PETSc.mat_getitem (src/petsc4py.PETSc.c:27761)
File "petscmat.pxi", line 877, in petsc4py.PETSc.matgetvalues (src/petsc4py.PETSc.c:26551)
ValueError: total size of new array must be unchanged


*** (2) higherorder=False ***

eikemueller@Eikes-MacBook-Pro $ python mat_assembly.py 
pyop2:WARNING No maximum assembly cache size. Install psutil or risk leaking memory!
Traceback (most recent call last):
File "mat_assembly.py", line 47, in <module>
  print mat.M.values
File "<string>", line 2, in values
File "/Users/eikemueller/PostDocBath/EllipticSolvers/PyOP2/pyop2/versioning.py", line 111, in modifies
  retval = method(self, *args, **kwargs)
File "/Users/eikemueller/PostDocBath/EllipticSolvers/PyOP2/pyop2/petsc_base.py", line 434, in values
  return self.handle[:, :]
File "Mat.pyx", line 218, in petsc4py.PETSc.Mat.__getitem__ (src/petsc4py.PETSc.c:99054)
File "petscmat.pxi", line 962, in petsc4py.PETSc.mat_getitem (src/petsc4py.PETSc.c:27761)
File "petscmat.pxi", line 876, in petsc4py.PETSc.matgetvalues (src/petsc4py.PETSc.c:26531)
File "arraynpy.pxi", line 85, in petsc4py.PETSc.empty_s (src/petsc4py.PETSc.c:7556)
ValueError: negative dimensions are not allowed

I’m using the “dmplex-1d-refinement" branch of PETSc, "columnwise_kernels" branch of PyOP2 and "multigrid-extrusion" branch of firedrake and I ran firedrake-clean.

Since it works in 2d, I suspect there is something wrong somewhere else, and I thought I let you know.

Thanks a lot (and happy New Year!),

Eike

from firedrake import *
from ffc import log
log.set_level(log.ERROR)
op2.init(log_level="WARNING")
parameters["coffee"]["O2"] = False

dimension = 3
higherorder = True

if (dimension == 2):
  ncells=16
  host_mesh = CircleManifoldMesh(ncells)
else:
  host_mesh = UnitIcosahedralSphereMesh(1)

D = 0.1
nlayers = 4
nlevel = 4
host_mesh_hierarchy = MeshHierarchy(host_mesh,nlevel)
mesh_hierarchy = ExtrudedMeshHierarchy(host_mesh_hierarchy,
                                     layers=nlayers,
                                     extrusion_type='radial',
                                     layer_height=D/nlayers)


mesh = mesh_hierarchy[-1]

if (dimension == 2):
  if (higherorder):
      U1 = FiniteElement('CG',interval,2)
  else:
      U1 = FiniteElement('CG',interval,1)
else:
  if (higherorder):
      U1 = FiniteElement('BDFM',triangle,2)
  else:
      U1 = FiniteElement('RT',triangle,1)

V1 = FiniteElement('DG',interval,0)
W2 = FunctionSpace(mesh, HDiv(OuterProductElement(U1,V1)))

u = TestFunction(W2)
v = TrialFunction(W2)


In 3D with this setup, W2 has 491520 dofs.

mat = assemble(dot(u,v)*mesh._dx)

print mat.M.values

Here you're asking for a dense array with

491520 * 491520 ~= 2^38

entries.  This is bigger than the 32bit int petsc has been built with.  So you're getting an error due to integer overflow.  You probably didn't want to look at this matrix directly on a grid this big though, no?  With 2 levels and 2 layers in 3D this dense matrix already uses 1.8 GB

Lawrence
_______________________________________________
firedrake mailing list
firedrake@imperial.ac.uk
https://mailman.ic.ac.uk/mailman/listinfo/firedrake