Dolfinx.cpp.mesh.h not working

dolfinx
1

Using fenicsx 2:0.5.0.1

I took the function test_cell_h_prism() from https://github.com/FEniCS/dolfinx/blob/main/python/test/unit/mesh/test_mesh.py
to create the MWE below

import numpy as np
import ufl
from dolfinx.mesh import (CellType, create_unit_cube,)
from dolfinx import cpp as _cpp
from mpi4py import MPI


def test_cell_h_prism():
    N = 3
    mesh = create_unit_cube(MPI.COMM_WORLD, N, N, N, CellType.tetrahedron)
    tdim = mesh.topology.dim
    num_cells = mesh.topology.index_map(tdim).size_local
    cells = np.arange(num_cells, dtype=np.int32)
    h = _cpp.mesh.h(mesh, tdim, cells)
    assert np.allclose(h, np.sqrt(2 / (N**2)))

test_cell_h_prism()

but it yields the following error message

Traceback (most recent call last):
  File "/home/bruno/Documents/Fenics/Heat_Spreading/test_h.py", line 17, in <module>
    test_cell_h_prism()
  File "/home/bruno/Documents/Fenics/Heat_Spreading/test_h.py", line 14, in test_cell_h_prism
    h = _cpp.mesh.h(mesh, tdim, cells)
        ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
TypeError: h(): incompatible function arguments. The following argument types are supported:
    1. (mesh: dolfinx::mesh::Mesh, dim: int, entities: numpy.ndarray[numpy.int32]) -> numpy.ndarray[numpy.float64]

Invoked with: <dolfinx.mesh.Mesh object at 0x7f79f4786d50>, 3, array([  0,   1,   2,   3,   4,   5,   6,   7,   8,   9,  10,  11,  12,
        13,  14,  15,  16,  17,  18,  19,  20,  21,  22,  23,  24,  25,
        26,  27,  28,  29,  30,  31,  32,  33,  34,  35,  36,  37,  38,
        39,  40,  41,  42,  43,  44,  45,  46,  47,  48,  49,  50,  51,
        52,  53,  54,  55,  56,  57,  58,  59,  60,  61,  62,  63,  64,
        65,  66,  67,  68,  69,  70,  71,  72,  73,  74,  75,  76,  77,
        78,  79,  80,  81,  82,  83,  84,  85,  86,  87,  88,  89,  90,
        91,  92,  93,  94,  95,  96,  97,  98,  99, 100, 101, 102, 103,
       104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116,
       117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129,
       130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142,
       143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155,
       156, 157, 158, 159, 160, 161], dtype=int32)

Has something changed in the h()function?

2

I cannot reproduce your error using docker ( docker run --rm -it -v $(pwd):/root/shared -w /root/shared dolfinx/dolfinx:v0.5.0).

How did you install DOLFINx?

import numpy as np
import ufl
from dolfinx.mesh import (CellType, create_unit_cube,)
from dolfinx import cpp as _cpp
from mpi4py import MPI


def test_cell_h_prism():
    N = 3
    mesh = create_unit_cube(MPI.COMM_WORLD, N, N, N, CellType.tetrahedron)
    tdim = mesh.topology.dim
    num_cells = mesh.topology.index_map(tdim).size_local
    cells = np.arange(num_cells, dtype=np.int32)
    h = _cpp.mesh.h(mesh, tdim, cells)
    print(h)
    assert np.allclose(h, np.sqrt(3 / (N**2)))

test_cell_h_prism()
3

This might be the irritating mesh vs mesh._cpp_object issue that came up somewhere around version 0.6.0.

4

I installed it with Debian 12 in a virtual python 11 environment. Legacy Fenics is also installed:

fenics/stable,now 2:0.5.0.1 amd64 [installed]
fenics/stable 2:0.5.0.1 i386
fenicsx-performance-tests-source/stable,stable 0.5.0~git20220731.821823b-1 all
fenicsx-performance-tests/stable 0.5.0~git20220731.821823b-1+b1 amd64
fenicsx-performance-tests/stable 0.5.0~git20220731.821823b-1+b1 i386
fenicsx/stable,stable,now 2:0.5.0.1 all [installed]
5

I can reproduce it with dolfinx 0.6.0 (debian package). It might be a v0.6.0 issue, if not a python version issue (python 3.11)

6

I couldn’t reproduce it with the docker image docker run --rm -it -v $(pwd):/root/shared -w /root/shared dolfinx/dolfinx:v0.6.0-r1 which uses Python 3.10.

I know we have some issues with Python 3.11 for older versions of DOLFINx, see; https://github.com/FEniCS/dolfinx/pull/2500
and
Test on Python 3.11 · Issue #2423 · FEniCS/dolfinx · GitHub

1 Like
7

Using the links you provided I could make the workaround below adding mesh = _cpp.mesh.Mesh(mesh.comm, mesh.topology, mesh.geometry):

import numpy as np
import ufl
from dolfinx.mesh import (CellType, create_unit_cube,)
from dolfinx import cpp as _cpp
from mpi4py import MPI


def test_cell_h_prism():
    N = 3
    mesh = create_unit_cube(MPI.COMM_WORLD, N, N, N, CellType.tetrahedron)
    tdim = mesh.topology.dim
    num_cells = mesh.topology.index_map(tdim).size_local
    cells = np.arange(num_cells, dtype=np.int32)
    mesh = _cpp.mesh.Mesh(mesh.comm, mesh.topology, mesh.geometry)
    h = _cpp.mesh.h(mesh, tdim, cells)
    print("hmin = ", h.min())
    assert np.allclose(h, np.sqrt(3 / (N**2)))

test_cell_h_prism()

yielding hmin = 0.5773502691896257

Thank’s a lot four your help!