Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

backport of #837, #838, #834 plus distutils fix #835

Merged
merged 9 commits into from
Nov 29, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .github/workflows/test.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ jobs:
- name: Configure MKL
if: ${{ matrix.matrix-backend == 'mkl' }}
run: |
python -um pip install --upgrade --upgrade-strategy eager mkl
python -um pip install --upgrade --upgrade-strategy eager 'mkl<2024'
python -um devtools.gha.configure_mkl
- name: Test
run: python -um coverage run -m unittest discover -b -q -t . -s tests
Expand Down
2 changes: 1 addition & 1 deletion docs/sphinx_mods.py
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ def generate_api(app):
dst_root.mkdir(parents=True, exist_ok=True)

srcs = tuple(f for f in sorted(nutils.glob('**/*.py')) if f != nutils/'__init__.py' and (f.name == '__init__.py' or not f.name.startswith('_')))
for src in sphinx.util.status_iterator(srcs, 'generating api... ', 'purple', len(srcs), app.verbosity):
for src in sphinx.util.display.status_iterator(srcs, 'generating api... ', 'purple', len(srcs), app.verbosity):
module = '.'.join((src.parent if src.name == '__init__.py' else src.with_suffix('')).relative_to(nutils).parts)
dst = dst_root/(module+'.rst')
with dst.open('w', encoding='utf-8') as f:
Expand Down
7 changes: 0 additions & 7 deletions nutils/__init__.py
Original file line number Diff line number Diff line change
@@ -1,12 +1,5 @@
'Numerical Utilities for Finite Element Analysis'

import sys
import numpy
from distutils.version import LooseVersion

assert sys.version_info >= (3, 5)
assert LooseVersion(numpy.version.version) >= LooseVersion('1.16'), 'nutils requires numpy 1.16 or higher, got {}'.format(numpy.version.version)

__version__ = version = '8.3'
version_name = 'idiyappam'
long_version = ('{} "{}"' if version_name else '{}').format(version, version_name)
Expand Down
9 changes: 2 additions & 7 deletions nutils/function.py
Original file line number Diff line number Diff line change
Expand Up @@ -2733,12 +2733,6 @@ def get(__array: IntoArray, __axis: int, __index: IntoArray) -> Array:
return numpy.take(Array.cast(__array), Array.cast(__index, dtype=int, ndim=0), __axis)


def _range(__length: int, __offset: int) -> Array:
length = Array.cast(__length, dtype=int, ndim=0)
offset = Array.cast(__offset, dtype=int, ndim=0)
return _Wrapper(lambda l, o: evaluable.Range(l) + o, _WithoutPoints(length), _WithoutPoints(offset), shape=(__length,), dtype=int)


def _takeslice(__array: IntoArray, __s: slice, __axis: int) -> Array:
array = Array.cast(__array)
s = __s
Expand All @@ -2749,7 +2743,8 @@ def _takeslice(__array: IntoArray, __s: slice, __axis: int) -> Array:
stop = n if s.stop is None else s.stop if s.stop >= 0 else s.stop + n
if start == 0 and stop == n:
return array
index = _range(stop-start, start)
length = stop - start
index = _Wrapper(evaluable.Range, _WithoutPoints(_Constant(length)), shape=(length,), dtype=int) + start
elif isinstance(n, numbers.Integral):
index = Array.cast(numpy.arange(*s.indices(int(n))))
else:
Expand Down
4 changes: 2 additions & 2 deletions nutils/matrix/_mkl.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,8 +30,8 @@ def assemble(data, index, shape):
class Pardiso:
'''Wrapper for libmkl.pardiso.

https://www.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/
sparse-solver-routines/onemkl-pardiso-parallel-direct-sparse-solver-iface.html
https://www.intel.com/content/www/us/en/docs/onemkl/developer-reference-c/2023-2/
onemkl-pardiso-parallel-direct-sparse-solver-iface.html
'''

_errorcodes = {
Expand Down
40 changes: 30 additions & 10 deletions nutils/topology.py
Original file line number Diff line number Diff line change
Expand Up @@ -1866,6 +1866,19 @@
def __rsub__(self, other):
return other

@property
def connectivity(self):
return ()

Check warning on line 1871 in nutils/topology.py

View check run for this annotation

Codecov / codecov/patch

nutils/topology.py#L1871

Added line #L1871 was not covered by tests

@property
def boundary(self):
if self.ndims:
return EmptyTopology(self.space, self.transforms.todims, self.ndims - 1)
else:
raise ValueError('A 0D topology has no boundary.')

Check warning on line 1878 in nutils/topology.py

View check run for this annotation

Codecov / codecov/patch

nutils/topology.py#L1878

Added line #L1878 was not covered by tests

interfaces = boundary


def StructuredLine(space, root: transform.TransformItem, i: int, j: int, periodic: bool = False, bnames: Optional[Tuple[str, str]] = None):
assert isinstance(i, int), f'i={i!r}'
Expand Down Expand Up @@ -2350,22 +2363,28 @@
class SimplexTopology(TransformChainsTopology):
'simpex topology'

def _renumber(simplices):
simplices = numpy.asarray(simplices)
keep = numpy.zeros(simplices.max()+1, dtype=bool)
keep[simplices.flat] = True
return types.arraydata(simplices if keep.all() else (numpy.cumsum(keep)-1)[simplices])

def __init__(self, space: str, simplices: numpy.ndarray, transforms: transformseq.Transforms, opposites: transformseq.Transforms):
assert isinstance(space, str), f'space={space!r}'
assert isinstance(simplices, numpy.ndarray), f'simplices={simplices!r}'
assert len(simplices), f'simplices is empty'
assert simplices.shape == (len(transforms), transforms.fromdims+1)
self.simplices = numpy.asarray(simplices)
assert numpy.greater(self.simplices[:, 1:], self.simplices[:, :-1]).all(), 'nodes should be sorted'
assert not numpy.equal(self.simplices[:, 1:], self.simplices[:, :-1]).all(), 'duplicate nodes'
references = References.uniform(element.getsimplex(transforms.fromdims), len(transforms))
super().__init__(space, references, transforms, opposites)

@cached_property
def contiguous_simplices(self):
simplices = numpy.asarray(self.simplices)
keep = numpy.zeros(simplices.max()+1, dtype=bool)
keep[simplices.flat] = True
return simplices if keep.all() else (numpy.cumsum(keep)-1)[simplices]

@cached_property
def nverts(self):
return self.contiguous_simplices.max() + 1

def take_unchecked(self, indices):
space, = self.spaces
return SimplexTopology(space, self.simplices[indices], self.transforms[indices], self.opposites[indices])
Expand All @@ -2374,6 +2393,8 @@
def boundary(self):
space, = self.spaces
ielem, iedge = (self.connectivity == -1).nonzero()
if len(ielem) == 0:
return self.empty_like().boundary
nd = self.ndims
edges = numpy.arange(nd+1).repeat(nd).reshape(nd,nd+1).T[::-1]
simplices = self.simplices[ielem, edges[iedge].T].T
Expand All @@ -2398,7 +2419,7 @@
def basis_std(self, degree):
if degree == 1:
coeffs = element.getsimplex(self.ndims).get_poly_coeffs('bernstein', degree=1)
return function.PlainBasis([coeffs] * len(self), self.simplices, self.simplices.max()+1, self.f_index, self.f_coords)
return function.PlainBasis([coeffs] * len(self), self.contiguous_simplices, self.nverts, self.f_index, self.f_coords)
return super().basis_std(degree)

def basis_bubble(self):
Expand All @@ -2410,9 +2431,8 @@
coeffs[:-1] = poly.change_degree(bernstein, self.ndims, 1 + self.ndims) - bubble[None] / (self.ndims+1)
coeffs[-1] = bubble
coeffs = types.frozenarray(coeffs, copy=False)
nverts = self.simplices.max() + 1
ndofs = nverts + len(self)
nmap = [types.frozenarray(numpy.hstack([idofs, nverts+ielem]), copy=False) for ielem, idofs in enumerate(self.simplices)]
ndofs = self.nverts + len(self)
nmap = [types.frozenarray(numpy.hstack([idofs, self.nverts+ielem]), copy=False) for ielem, idofs in enumerate(self.contiguous_simplices)]
return function.PlainBasis([coeffs] * len(self), nmap, ndofs, self.f_index, self.f_coords)


Expand Down
2 changes: 1 addition & 1 deletion pyproject.toml
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,7 @@ classifiers = [
[project.optional-dependencies]
docs = ["Sphinx>=1.8"]
export_mpl = ["matplotlib>=1.3", "pillow>2.6"]
matrix_mkl = ["mkl"]
matrix_mkl = ["mkl<2024"]
matrix_scipy = ["scipy>=0.13"]
import_gmsh = ["meshio"]

Expand Down
7 changes: 7 additions & 0 deletions tests/test_function.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
import warnings as _builtin_warnings
import functools
import fractions
import pickle


class Array(TestCase):
Expand Down Expand Up @@ -171,6 +172,12 @@ def test_lower_eval(self):
desired = self.n_op(*self.args)
self.assertArrayAlmostEqual(actual, desired, decimal=15)

def test_pickle(self):
f = self.op(*self.args)
s = pickle.dumps(f)
f_ = pickle.loads(s)
self.assertEqual(f.as_evaluable_array, f_.as_evaluable_array)


def generate(*shape, real, imag, zero, negative):
'generate array values that cover certain numerical classes'
Expand Down
21 changes: 21 additions & 0 deletions tests/test_topology.py
Original file line number Diff line number Diff line change
Expand Up @@ -1352,9 +1352,30 @@ def setUp(self):
self.desired_references = [element.TriangleReference()]*4
self.desired_vertices = coords[simplices].tolist()

def test_contiguous_simplices(self):
# Same structure as the simplices in `setUp`, but renumbered.
simplices = numpy.array([[10,13,25],[10,18,25],[13,21,25],[18,21,25]])
transforms = transformseq.IndexTransforms(2, len(simplices))
topo = topology.SimplexTopology('X', simplices, transforms, transforms)
self.assertEqual(self.topo.contiguous_simplices.tolist(), self.topo.simplices.tolist())
self.assertEqual(self.topo.nverts, 5)

def test_boundary(self):
self.assertIsInstance(self.topo.boundary, topology.SimplexTopology)

def test_empty_boundary(self):
simplices = numpy.array([[0, 1, 2, 3]])
transforms = transformseq.IndexTransforms(3, len(simplices))
topo = topology.SimplexTopology('X', simplices, transforms, transforms)
self.assertEqual(len(topo.boundary), 4)
self.assertEqual(len(topo.boundary.boundary), 0)

def test_empty_interfaces(self):
simplices = numpy.array([[0, 1, 2, 3]])
transforms = transformseq.IndexTransforms(3, len(simplices))
topo = topology.SimplexTopology('X', simplices, transforms, transforms)
self.assertEqual(len(topo.interfaces), 0)

def test_getitem(self):
self.assertIsInstance(self.topo[numpy.arange(4) < 2], topology.SimplexTopology)
self.assertIsInstance(self.topo[numpy.arange(2)], topology.SimplexTopology)
Expand Down
Loading