+
Skip to content

Extract C-contiguous arrays by default from fields #856

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.

Sign up for GitHub

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

Merged
merged 3 commits into from
Sep 23, 2024
Merged

Extract C-contiguous arrays by default from fields #856

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
95ec435
Extract field data as C-contiguous array
adtzlr Sep 23, 2024
2fc49e2
Change to default `np.einsum(..., order="C")` in `Field`
adtzlr Sep 23, 2024
78dee44
Update _container.py
adtzlr Sep 23, 2024
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
4 changes: 4 additions & 0 deletions CHANGELOG.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,10 @@ All notable changes to this project will be documented in this file. The format

### Added
- Add the hessian of the element shape functions for a quadratic quad element `QuadraticQuad.hessian()`.
- Add the `order`-argument to `FieldContainer.extract(order="C")` to return C-contiguous arrays by default.

### Changed
- Change default `np.einsum(..., order="K")` to `np.einsum(..., order="C")` in the methods of `Field`, `FieldAxisymmetric` and `FieldPlaneStrain`.

## [9.0.0] - 2024-09-06

Expand Down
33 changes: 24 additions & 9 deletions src/felupe/field/_axi.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@ def __init__(self, region, dim=2, values=0.0, dtype=None):
# in the region
self.radius = self.scalar.interpolate()

def _interpolate_2d(self, dtype=None, out=None):
def _interpolate_2d(self, dtype=None, out=None, order="C"):
"""Interpolate 2D field values at points and evaluate them at the
integration points of all cells in the region."""

Expand All @@ -101,20 +101,22 @@ def _interpolate_2d(self, dtype=None, out=None):
self.region.h,
dtype=dtype,
out=out,
order=order,
)

def interpolate(self, dtype=None, out=None):
def interpolate(self, dtype=None, out=None, order="C"):
# out-argument is not supported
# if out is not None:
# out = out[:2]

# extend dimension of in-plane 2d-gradient
# extend dimension of in-plane 2d-gradient (out-keyword can't be used here)
return np.pad(
self._interpolate_2d(dtype=dtype, out=None), ((0, 1), (0, 0), (0, 0))
self._interpolate_2d(dtype=dtype, out=None, order=order),
((0, 1), (0, 0), (0, 0)),
)

def _grad_2d(self, sym=False, dtype=None, out=None):
"""In-plane 2D gradient as partial derivative of field values at points
def _grad_2d(self, sym=False, dtype=None, out=None, order="C"):
r"""In-plane 2D gradient as partial derivative of field values at points
w.r.t. the undeformed coordinates, evaluated at the integration points
of all cells in the region. Optionally, the symmetric part of the
gradient is returned.
Expand All @@ -130,6 +132,12 @@ def _grad_2d(self, sym=False, dtype=None, out=None):
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
order : {'C', 'F', 'A', 'K'}, optional
Controls the memory layout of the output. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -148,14 +156,15 @@ def _grad_2d(self, sym=False, dtype=None, out=None):
self.region.dhdX,
dtype=dtype,
out=out,
order=order,
)

if sym:
return symmetric(g, out=g)
else:
return g

def grad(self, sym=False, dtype=None, out=None):
def grad(self, sym=False, dtype=None, out=None, order="C"):
r"""3D-gradient as partial derivative of field values at points w.r.t.
the undeformed coordinates, evaluated at the integration points of all
cells in the region. Optionally, the symmetric part of the gradient is
Expand All @@ -182,6 +191,12 @@ def grad(self, sym=False, dtype=None, out=None):
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
order : {'C', 'F', 'A', 'K'}, optional
Controls the memory layout of the output. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -197,11 +212,11 @@ def grad(self, sym=False, dtype=None, out=None):

# extend dimension of in-plane 2d-gradient
g = np.pad(
self._grad_2d(sym=sym, dtype=dtype, out=None),
self._grad_2d(sym=sym, dtype=dtype, out=None, order=order),
((0, 1), (0, 1), (0, 0), (0, 0)),
)

# set dudX_33 = u_r / R
g[-1, -1] = self.interpolate(dtype=dtype)[1] / self.radius
g[-1, -1] = self.interpolate(dtype=dtype, order=order)[1] / self.radius

return g
43 changes: 36 additions & 7 deletions src/felupe/field/_base.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -138,7 +138,7 @@ def _indices_per_cell(self, cells, dim):

return cai, ai

def grad(self, sym=False, dtype=None, out=None):
def grad(self, sym=False, dtype=None, out=None, order="C"):
r"""Gradient as partial derivative of field values w.r.t. undeformed
coordinates, evaluated at the integration points of all cells in the region.
Optionally, the symmetric part the gradient is evaluated.
Expand All @@ -159,6 +159,12 @@ def grad(self, sym=False, dtype=None, out=None):
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
order : {'C', 'F', 'A', 'K'}, optional
Controls the memory layout of the output. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -177,14 +183,15 @@ def grad(self, sym=False, dtype=None, out=None):
self.region.dhdX,
dtype=dtype,
out=out,
order=order,
)

if sym:
return symmetric(g, out=g)
else:
return g

def hess(self, dtype=None, out=None):
def hess(self, dtype=None, out=None, order="C"):
r"""Hessian as second partial derivative of field values w.r.t. undeformed
coordinates, evaluated at the integration points of all cells in the region.

Expand All @@ -204,6 +211,12 @@ def hess(self, dtype=None, out=None):
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
order : {'C', 'F', 'A', 'K'}, optional
Controls the memory layout of the output. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -222,11 +235,12 @@ def hess(self, dtype=None, out=None):
self.region.d2hdXdX,
dtype=dtype,
out=out,
order=order,
)

return h

def interpolate(self, dtype=None, out=None):
def interpolate(self, dtype=None, out=None, order="C"):
r"""Interpolate field values located at mesh-points to the quadrature points
``q`` of cells ``c`` in the region.

Expand All @@ -243,6 +257,12 @@ def interpolate(self, dtype=None, out=None):
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
order : {'C', 'F', 'A', 'K'}, optional
Controls the memory layout of the output. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -259,10 +279,13 @@ def interpolate(self, dtype=None, out=None):
self.values[self.region.mesh.cells],
self.region.h,
dtype=dtype,
out=None,
out=out,
order=order,
)

def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None):
def extract(
self, grad=True, sym=False, add_identity=True, dtype=None, out=None, order="C"
):
"""Generalized extraction method which evaluates either the gradient or the
field values at the integration points of all cells in the region. Optionally,
the symmetric part of the gradient is evaluated and/or the identity matrix is
Expand All @@ -284,6 +307,12 @@ def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None)
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
order : {'C', 'F', 'A', 'K'}, optional
Controls the memory layout of the output. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -298,7 +327,7 @@ def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None)
"""

if grad:
gr = self.grad(out=out, dtype=dtype)
gr = self.grad(out=out, dtype=dtype, order=order)

if sym:
gr = symmetric(gr, out=gr)
Expand All @@ -308,7 +337,7 @@ def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None)

return gr
else:
return self.interpolate(out=out, dtype=dtype)
return self.interpolate(out=out, dtype=dtype, order=order)

def copy(self):
"Return a copy of the field."
Expand Down
19 changes: 16 additions & 3 deletions src/felupe/field/_container.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -117,7 +117,9 @@ def __repr__(self):

return "\n".join([header, size, fields_header, *fields])

def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None):
def extract(
self, grad=True, sym=False, add_identity=True, dtype=None, out=None, order="C"
):
"""Generalized extraction method which evaluates either the gradient or the
field values at the integration points of all cells in the region. Optionally,
the symmetric part of the gradient is evaluated and/or the identity matrix is
Expand All @@ -142,6 +144,12 @@ def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None)
A location into which the result is stored. If provided, it must have a
shape that the inputs broadcast to. If not provided or None, a freshly-
allocated array is returned (default is None).
orders : str or list of str, optional
Controls the memory layout of the outputs. 'C' means it should be C
contiguous. 'F' means it should be Fortran contiguous, 'A' means it should
be 'F' if the inputs are all 'F', 'C' otherwise. 'K' means it should be as
close to the layout as the inputs as is possible, including arbitrarily
permuted axes. Default is 'C'.

Returns
-------
Expand All @@ -153,13 +161,18 @@ def extract(self, grad=True, sym=False, add_identity=True, dtype=None, out=None)
if isinstance(grad, bool):
grad = (grad,)

if isinstance(order, str):
order = (order,)

if out is None:
out = [None] * len(self.fields)

grads = np.pad(grad, (0, len(self.fields) - 1))
orders = order * len(self.fields)

return tuple(
f.extract(g, sym, add_identity=add_identity, dtype=dtype, out=res)
for g, f, res in zip(grads, self.fields, out)
f.extract(g, sym, add_identity=add_identity, dtype=dtype, out=res, order=od)
for g, f, res, od in zip(grads, self.fields, out, orders)
)

def values(self):
Expand Down
Loading
Loading
点击 这是indexloc提供的php浏览器服务，不要输入任何密码和下载