Tria3DSG - Triangular element with discrete shear gap (pyfe3d.tria3dsg)

Three-node triangular shell element with 6 degrees-of-freedom per node, \(u\), \(v\), \(w\), \(r_x\), \(r_y\), \(r_z\), first-order shear deformation theory, and a transverse shear field obtained with the discrete shear gap (DSG) method:

Bletzinger, K.-U., Bischoff, M., and Ramm, E., 2000, “A unified approach for shear-locking-free triangular and rectangular shell finite elements”, Computers & Structures, 75(3), pp. 321-334. https://doi.org/10.1016/S0045-7949(99)00140-6

Why this element exists

pyfe3d.tria3r.Tria3R takes its transverse shear from a single point at the centroid, which still locks, and cures the locking by dividing the transverse shear stiffnesses by \(1 + \alpha \ell^2/h^2\) with \(\alpha = 0.7\), following Lyly, Stenberg and Vihinen (1993). That works for plate bending, where the error it introduces is \(O(1/n^2)\) in the number of elements per span, but it leaves the element unable to reproduce a state of constant transverse shear: measured on an 8 by 8 patch, its stiffness against a transverse displacement field with all rotations at zero is \(3 \times 10^{-1}\), \(2 \times 10^{-2}\) and \(1 \times 10^{-3}\) of the pyfe3d.quad4.Quad4 value at \(\ell/h\) of 1.25, 6.25 and 25. A Donnell geometric stiffness matrix works on \(\partial w/\partial x\), which is exactly that subspace, so thin-shell buckling suffers badly. See tests/test_transverse_shear_stiffness.py.

The DSG removes the locking in the operator instead of in the constitutive matrix, so no stabilisation parameter is needed and the transverse shear stiffness is used as the laminate gives it.

The algebraic stabilisation of Bischoff and Bletzinger (2004) and Castro et al. (2019) is nevertheless available through the alpha_shear_locking attribute, which is 0. by default and then has no effect whatsoever. It is there because in those papers the factor is applied on top of a discrete shear gap field exactly like this one, so setting it reproduces the element they describe, and because on a coarse mesh of a thin plate it does improve the answer: \(\alpha = 0.07\), the lower bound of the range recommended by Castro et al., takes the first natural frequency of a \(5 \times 7\) mesh from \(12.06\) to \(7.93\) per cent above the analytical value and the buckling load of a thin plate from \(+1.71\) to \(+0.11\) per cent. That is the opposite of its role in pyfe3d.tria3r.Tria3R, where the same factor is applied to a centroid-sampled shear field and has to do the unlocking itself, which is why that element defaults to \(\alpha = 0.7\).

The default here is zero because the factor is unbounded in \(\ell/h\), so any non-zero \(\alpha\) removes nearly all of the transverse shear stiffness on a practical shell mesh and with it the three properties listed above. The alpha_shear_locking attribute documentation carries the measurements and the recommended value.

The discrete shear gap

With \(\pmb{\phi} = (-r_y, r_x)\) the transverse shear strains of this element’s convention are

\[\gamma_{xz} = \frac{\partial w}{\partial x} + r_y \qquad \gamma_{yz} = \frac{\partial w}{\partial y} - r_x\]

that is \(\pmb{\gamma} = \nabla w - \pmb{\phi}\). The shear gap \(\Delta w\) is the function whose gradient is the shear strain, and its nodal values are obtained by integrating along the straight edges from node 1, with \(\pmb{\phi}\) interpolated linearly,

\[\Delta w_1 = 0 \qquad \Delta w_i = (w_i - w_1) - \frac{1}{2}(\pmb{\phi}_1 + \pmb{\phi}_i) \cdot (\pmb{r}_i - \pmb{r}_1)\]

The assumed shear strain is then the gradient of the interpolated gap,

\[\pmb{\gamma} = \nabla \left( \sum_i N_i \Delta w_i \right) = \sum_i \Delta w_i \nabla N_i\]

which is constant over the element because the shape functions are linear. Four properties follow, all of them verified symbolically and then numerically in tests/test_tria3dsg.py:

• the three rigid-body modes that involve \(w\) give \(\pmb{\gamma} = 0\) exactly;

• a state of constant transverse shear is reproduced exactly, which is what Tria3R cannot do;

• a state of constant curvature gives \(\pmb{\gamma} = 0\) exactly, so the element develops no parasitic shear in pure bending. This is the locking-free property, and it is the reason no stabilisation parameter appears anywhere in this element;

• the operator is a linear functional of the nodal values only and is constant over the element, so it is formed once per element.

The gap is integrated away from one node, which makes the plain operator depend on which node that is: measured on an irregular triangle, the element matrix moves by 14 and 17 per cent under the two cyclic relabellings, against \(3 \times 10^{-16}\) for pyfe3d.tria3r.Tria3R, whose shear comes from a symmetric centroid evaluation. A mesh generator orders the nodes of a triangle as it pleases, so this element averages the three operators, one per starting node. The averaging costs nothing: each of the three is separately exact for all of the properties listed above, and those are linear in the operator, so the mean inherits them, which the tests confirm. What it buys is invariance, measured at \(6 \times 10^{-16}\) in tests/test_tria3dsg.py::test_node_numbering_invariance.

Membrane, bending and drilling

The in-plane and curvature fields are the constant-strain ones of the linear triangle, so no locking arises there and one point integrates them exactly. The drilling rotation \(r_z\) is given stiffness with the regularisation of Hughes and Brezzi (1989), tying \(r_z\) to the in-plane rotation \(\theta_z = (v_{,x} - u_{,y})/2\), exactly as in pyfe3d.quad4.Quad4 and pyfe3d.tria3r.Tria3R:

Hughes, T. J. R., and Brezzi, F., 1989, “On drilling degrees of freedom,” Computer Methods in Applied Mechanics and Engineering, 72(1), pp. 105-121. https://doi.org/10.1016/0045-7825(89)90124-2

drilling_model = 0, the default, uses the physics-based coefficient \(\gamma_{r_z} = A_{66}\), a modulus of the laminate and not an adjustable parameter. Any other value uses the fictitious penalty \(\gamma_{r_z} = K6ROT \cdot 10^{-6} A_{66}\).

Note

This element does not carry the hierarchical edge-mode enrichment of Allman (1984) that pyfe3d.tria3r.Tria3R and pyfe3d.quad4.Quad4 apply to the in-plane displacement field, so its membrane response in in-plane bending is the stiff one of the constant-strain triangle. The element is aimed at the transverse shear behaviour; for in-plane dominated problems prefer the enriched elements.

Geometrically nonlinear analysis

The element carries the von Karman membrane strains,

\[\epsilon_{xx} = u_{,x} + \frac{1}{2} w_{,x}^2, \quad \epsilon_{yy} = v_{,y} + \frac{1}{2} w_{,y}^2, \quad \gamma_{xy} = u_{,y} + v_{,x} + w_{,x} w_{,y}\]

with the same split of the tangent stiffness as pyfe3d.quad4.Quad4 and pyfe3d.tria3r.Tria3R,

\[[K_T] = [K_{C_0}] + [K_{CNL}(u)] + [K_G(u)]\]

which is the exact Jacobian of the internal forces of Tria3DSG.update_fint() called with nonlinear=1. \([K_G]\) is left homogeneous of degree one in the displacements so that a linear buckling analysis can use it directly, the stress of the nonlinear membrane strain being collected in \([K_{CNL}]\) instead. See Tria3DSG.update_KCNL().

Two things are simpler here than in the quadrilateral. First, the transverse shear strains of first-order shear deformation theory have no von Karman terms, so the discrete shear gap operator appears in \([K_{C_0}]\) alone and the nonlinear methods never touch it. Second, every operator of this element is constant over the triangle, so \([K_{CNL}]\), its internal forces and \([K_G]\) are all integrated exactly by a single multiplication by the area, with no quadrature loop and no quadrature to keep consistent between the matrix and the residual. The drilling term is the only one that needs more than one point, and it is linear.

class pyfe3d.tria3dsg.Tria3DSG

Three-node triangular shell element with a discrete shear gap

Attributes:

eid,int

Element identification number.

pid,int

Property identification number.

area,double

Element area.

n1, n2, n3int

Node identification numbers.

c1, c2, c3int

Position of each node’s first degree-of-freedom in the global displacement vector.

init_k_KC0, init_k_KCNL, init_k_KG, init_k_Mint

Position of this element’s first entry in the sparse vectors.

drilling_model,int

0, the default, uses \(\gamma_{r_z} = A_{66}\) of Hughes and Brezzi. Any other value uses \(K6ROT \cdot 10^{-6} A_{66}\).

gamma_rz,double

When non-negative and drilling_model = 0, overrides \(A_{66}\) as the regularisation coefficient, which is what makes the sensitivity study recommended in the literature possible. Negative means use \(A_{66}\).

K6ROT,double

Coefficient of the fictitious penalty, used when drilling_model != 0.

alpha_shear_locking,double

Optional algebraic stabilisation of the transverse shear stiffness, \(\hat{A}_{ts} = A_{ts}/(1 + \alpha \ell^2/h^2)\) with \(\ell\) the longest edge and \(h\) the total thickness, applied after the rotation into the element coordinate system and after the shear correction. The default 0. disables it entirely and is the element as described above: the discrete shear gap is already locking-free, so nothing has to be rescaled.

A non-zero value reproduces the element of Bischoff and Bletzinger (2004) and Castro et al. (2019) as those papers actually define it, that is the stabilisation applied on top of a discrete shear gap field, which is not how pyfe3d.tria3r.Tria3R uses the same symbol. There the factor is applied to a centroid-sampled shear field and has to do the unlocking itself, which is why its default is 0.7, seven times the literature value, and why it is load bearing there and merely helpful here.

Recommended value. 0.07, the lower bound of the range recommended by Castro et al. (2019), on a thin plate discretised coarsely, and only there. Measured with this element, with the remaining columns from tests/test_tria3dsg.py and the study described in the module documentation:

Methods

update_KC0(self, long[, long[, double[, ...)	Update the sparse vectors of the constitutive stiffness matrix
update_KCNL(self, long[, long[, double[, ...)	Update the sparse vectors of the nonlinear constitutive stiffness matrix
update_KG(self, long[, long[, double[, ...)	Update the geometric stiffness from the current displacements
update_KG_given_stress(self, double Nxx, ...)	Update the sparse vectors of the geometric stiffness matrix
update_M(self, long[, long[, double[, ...)	Update the sparse vectors of the mass matrix
update_area(self)	Update the element area from the probe coordinates
update_fint(self, double[, ShellProp prop, ...)	Update the internal force vector
update_probe_finte(self, ShellProp prop, ...)	Update the internal force vector of the probe
update_probe_ue(self, double[)	Update the element-frame nodal displacements in the probe
update_probe_xe(self, double[)	Update the element-frame nodal coordinates in the probe
update_rotation_matrix(self, double[, ...)	Update the rotation matrix of the element

K6ROT

K6ROT: ‘double’

alpha_shear_locking

alpha_shear_locking: ‘double’

area

area: ‘double’

c1

c1: ‘int’

c2

c2: ‘int’

c3

c3: ‘int’

drilling_model

drilling_model: ‘int’

eid

eid: ‘int’

gamma_rz

gamma_rz: ‘double’

init_k_KC0

init_k_KC0: ‘int’

init_k_KCNL

init_k_KCNL: ‘int’

init_k_KG

init_k_KG: ‘int’

init_k_M

init_k_M: ‘int’

m11

m11: ‘double’

m12

m12: ‘double’

m21

m21: ‘double’

m22

m22: ‘double’

n1

n1: ‘int’

n2

n2: ‘int’

n3

n3: ‘int’

pid

pid: ‘int’

probe

probe: pyfe3d.tria3dsg.Tria3DSGProbe

r11

r11: ‘double’

r12

r12: ‘double’

r13

r13: ‘double’

r21

r21: ‘double’

r22

r22: ‘double’

r23

r23: ‘double’

r31

r31: ‘double’

r32

r32: ‘double’

r33

r33: ‘double’

update_KC0(self, long[: :1] KC0r, long[: :1] KC0c, double[: :1] KC0v, ShellProp prop, int update_KC0v_only=0) → void

Update the sparse vectors of the constitutive stiffness matrix

The element-frame matrix is built by _update_probe_KC0ve() into the KC0ve attribute of the probe and rotated to global coordinates here.

Before this method is called, update_probe_xe() must have been used to set the element-frame coordinates of the probe.

Parameters:

KC0r, KC0carray-like

Row and column positions of the sparse values.

KC0varray-like

Sparse values.

proppyfe3d.shellprop.ShellProp

Shell property.

update_KC0v_onlyint, optional

0, the default, also fills KC0r and KC0c.

update_KCNL(self, long[: :1] KCNLr, long[: :1] KCNLc, double[: :1] KCNLv, ShellProp prop, int update_KCNLv_only=0) → void

Update the sparse vectors of the nonlinear constitutive stiffness matrix

Assuming that KCNL = KC0L + KCL0 + KCLL + KGNL, built from the von Karman membrane strains

\[\epsilon_{xx} = u_{,x} + \frac{1}{2} w_{,x}^2, \quad \epsilon_{yy} = v_{,y} + \frac{1}{2} w_{,y}^2, \quad \gamma_{xy} = u_{,y} + v_{,x} + w_{,x} w_{,y}\]

whose nonlinear part is \(\{\epsilon_{NL}\} = \frac{1}{2} [B_{mL}] \{u_e\}\), with \([B_{mL}]\) its variation. With \([B_m]\) and \([B_b]\) the linear membrane and bending strain-displacement matrices, \([G]\) the gradient of \(w\), and \([A]\), \([B]\) the laminate matrices:

• KC0L = \([B_m]^T [A] [B_{mL}] + [B_b]^T [B] [B_{mL}]\)

• KCL0 = KC0L`^T`

• KCLL = \([B_{mL}]^T [A] [B_{mL}]\)

• KGNL = \([G]^T [N_{NL}] [G]\), with \(\{N_{NL}\} = [A] \{\epsilon_{NL}\}\)

The first three groups are the constitutive terms coupling the linear and the nonlinear parts of the membrane strain. KGNL is geometric, carrying the stress of the nonlinear membrane strain. It is collected here so that update_KG() stays homogeneous of degree one in the displacements, which is what a linear buckling analysis needs. With it here,

\[K_T = K_{C0} + K_{CNL}(u) + K_G(u)\]

is the exact Jacobian of the internal forces of update_fint() with nonlinear=1, and a Newton-Raphson iteration built on them converges quadratically.

Note

The transverse shear strains of first-order shear deformation theory carry no von Karman terms, so the discrete shear gap operator of this element enters \(K_{C0}\) alone and needs nothing here. The membrane rows have no \(r_z\) entries, this element having no Allman enrichment, so the drilling regularisation is likewise untouched. Both are unlike pyfe3d.quad4.Quad4, whose enriched membrane rows have to be carried into its nonlinear tangent.

Before this method is called, update_probe_xe() must have been used with the node coordinates and update_probe_ue() with the current displacements.

Parameters:

KCNLr, KCNLcarray-like

Row and column positions of the sparse values.

KCNLvarray-like

Sparse values.

proppyfe3d.shellprop.ShellProp

Shell property.

update_KCNLv_onlyint, optional

0, the default, also fills KCNLr and KCNLc.

update_KG(self, long[: :1] KGr, long[: :1] KGc, double[: :1] KGv, ShellProp prop, int update_KGv_only=0) → void

Update the geometric stiffness from the current displacements

The membrane stress resultants are recovered from the probe displacements, which update_probe_ue() must have updated, using the same strain operators and the same membrane constitutive relation that update_fint() uses, and then handed to update_KG_given_stress(). Both the strains and the curvatures are constant over the triangle, so the resulting stress state is constant too, which is what that method assumes.

Only the linear strains contribute, which leaves this matrix homogeneous of degree one in the displacements, as a linear buckling analysis needs. The stress of the von Karman terms of the membrane strain is carried by update_KCNL() instead.

Parameters:

KGr, KGcarray-like

Row and column positions of the sparse values.

KGvarray-like

Sparse values.

proppyfe3d.shellprop.ShellProp

Shell property.

update_KGv_onlyint, optional

0, the default, also fills KGr and KGc.

update_KG_given_stress(self, double Nxx, double Nyy, double Nxy, long[: :1] KGr, long[: :1] KGc, double[: :1] KGv, int update_KGv_only=0) → void

Update the sparse vectors of the geometric stiffness matrix

A constant membrane stress state \(N_{xx}, N_{yy}, N_{xy}\) is assumed within the element. The operator is the Donnell one, involving only the gradients of the element-frame \(w\), which are taken from the Gwx, Gwy rows of the probe and turned into global translations with the third row of the rotation matrix, so only the translational degrees-of-freedom appear.

Warning

The stress resultants are in the element coordinate system, and for a triangle that is easy to get wrong. The element \(x\) axis runs from node 1 to node 2, so in a structured mesh whose cells are split along a diagonal the first triangle of each cell has its frame along the mesh while the second has it along the diagonal, some 45 degrees away. Handing the same \((N_{xx}, 0, 0)\) to both then loads half the elements along the diagonal, and the answer comes to depend on how the mesh generator happened to order the nodes of each triangle.

A uniaxial state \(N\) along a global direction \(\{d\}\) has to be rotated. With \([R]\) the element-to-global rotation of this object, \(\{a\} = [R]^T \{d\}\) is that direction in element coordinates and the state \(N \{a\} \{a\}^T\) has components

\[N^e_{xx} = N a_x^2, \quad N^e_{yy} = N a_y^2, \quad N^e_{xy} = N a_x a_y\]

which reduces to \((N, 0, 0)\) when the element \(x\) axis is along \(\{d\}\), so a quadrilateral mesh aligned with the load needs nothing. Measured on the cylinder of tests/test_quad4_linear_buckling_cylinder_displ.py, whose cells are nearly square so the diagonal sits at 44.3 degrees, the three cyclic numberings of one triangle give buckling loads spanning a factor of fifteen without this rotation and agree to 2e-3 with it. This is checked by test_buckling_load_does_not_depend_on_triangle_node_numbering in tests/test_tria3dsg.py.

update_KG() is free of this, recovering the stress from the element’s own strains in its own frame.

Parameters:

Nxx, Nyy, Nxydouble

Membrane stress resultants in the element coordinate system.

KGr, KGcarray-like

Row and column positions of the sparse values.

KGvarray-like

Sparse values.

update_KGv_onlyint, optional

0, the default, also fills KGr and KGc.

update_M(self, long[: :1] Mr, long[: :1] Mc, double[: :1] Mv, ShellProp prop, int mtype=0) → void

Update the sparse vectors of the mass matrix

Parameters:

Mr, Mcarray-like

Row and column positions of the sparse values.

Mvarray-like

Sparse values.

proppyfe3d.shellprop.ShellProp

Shell property.

mtypeint, optional

0 for the consistent mass matrix, for which \(\int N_i N_j dA\) is \(A/6\) when \(i = j\) and \(A/12\) otherwise, and 2 for the lumped one, which puts \(A/3\) on each node. Any other value is treated as 0.

update_area(self) → void

Update the element area from the probe coordinates

update_fint(self, double[: :1] fint, ShellProp prop, int nonlinear=0) → void

Update the internal force vector

Parameters:

fintarray-like

Array updated in place with the internal forces, in global coordinates. update_probe_finte() is called to obtain the element-frame internal forces, available afterwards as .probe.finte.

proppyfe3d.shellprop.ShellProp

Shell property.

nonlinearint, optional

The default 0 gives the linear internal forces, KC0*u. Any other value adds the geometrically nonlinear terms of the von Karman strains, for which the exact Jacobian of the internal forces is KC0 + KCNL + KG, see update_KCNL().

update_probe_finte(self, ShellProp prop, int nonlinear=0) → void

Update the internal force vector of the probe

The attribute finte of the Tria3DSGProbe is updated with the internal forces in element coordinates. Mind that the probe can be shared amongst more than one finite element, so it always holds the values of the last update.

Before this method is called, update_probe_xe() and update_probe_ue() must have been used.

Parameters:

proppyfe3d.shellprop.ShellProp

Shell property.

nonlinearint, optional

The default 0 gives the linear internal forces, KC0*u. Any other value adds the geometrically nonlinear terms of the von Karman strains, for which the exact Jacobian of the internal forces is KC0 + KCNL + KG, see update_KCNL().

update_probe_ue(self, double[: :1] u) → void

Update the element-frame nodal displacements in the probe

Parameters:

uarray-like

Global displacement vector.

update_probe_xe(self, double[: :1] x) → void

Update the element-frame nodal coordinates in the probe

Parameters:

xarray-like

Global nodal coordinates.

update_rotation_matrix(self, double[: :1] x, double xmati=0., double xmatj=0., double xmatk=0.) → void

Update the rotation matrix of the element

The element \(x\) axis runs from node 1 to node 2 and the element \(z\) axis is the triangle normal, which is the same construction used by pyfe3d.tria3r.Tria3R.

Parameters:

xarray-like

Global nodal coordinates, x_1, y_1, z_1, ..., x_M, y_M, z_M.

xmati, xmatj, xmatkdouble, optional

Material direction in global coordinates, projected onto the element. Leaving it at zero means the material and element directions coincide.

class pyfe3d.tria3dsg.Tria3DSGData

Sizes needed to allocate the sparse matrices

Attributes:

KC0_SPARSE_SIZE,int

KC0_SPARSE_SIZE = 324, the full 18 by 18 element matrix stored row by row.

KCNL_SPARSE_SIZE,int

KCNL_SPARSE_SIZE = 324, also the full element matrix, since the nonlinear constitutive terms couple every degree-of-freedom.

KG_SPARSE_SIZE,int

KG_SPARSE_SIZE = 81, the 9 by 9 translational block, the Donnell geometric stiffness involving only the gradients of \(w\).

M_SPARSE_SIZE,int

M_SPARSE_SIZE = 324

KC0_SPARSE_SIZE

KC0_SPARSE_SIZE: ‘int’

KCNL_SPARSE_SIZE

KCNL_SPARSE_SIZE: ‘int’

KG_SPARSE_SIZE

KG_SPARSE_SIZE: ‘int’

M_SPARSE_SIZE

M_SPARSE_SIZE: ‘int’

class pyfe3d.tria3dsg.Tria3DSGProbe

Probe carrying the element-frame coordinates, displacements and operators

The idea behind using a probe is to avoid allocating one set of buffers per finite element. The buffers belong to the probe, and one probe can be shared amongst many elements, the contents being updated and read on demand. Mind that the probe always holds the values of the last update.

Attributes:

xe,array-like

Element-frame nodal coordinates, x1, y1, z1, x2, ..., z3.

ue,array-like

Element-frame nodal displacements, 18 values, in the order \(u_1, v_1, w_1, {r_x}_1, {r_y}_1, {r_z}_1, u_2, \ldots, {r_z}_3\).

finte,array-like

Element-frame internal forces corresponding to ue, 18 values, updated by Tria3DSG.update_probe_finte().

KC0ve,array-like

The 18 by 18 element-frame constitutive stiffness matrix stored row by row, 324 values.

KCNLve,array-like

The 18 by 18 element-frame nonlinear constitutive stiffness matrix stored row by row, see Tria3DSG.update_KCNL().

BLexx, BLeyy, BLgxyarray-like

Rows of 18 values giving the membrane strains \(\epsilon_{xx}, \epsilon_{yy}, \gamma_{xy}\).

BLkxx, BLkyy, BLkxyarray-like

Rows of 18 values giving the curvatures \(\kappa_{xx}, \kappa_{yy}, \kappa_{xy}\).

BLgxz, BLgyzarray-like

Rows of 18 values giving the transverse shear strains \(\gamma_{xz}, \gamma_{yz}\) of the discrete shear gap, symmetrised over the three starting nodes.

BLdrillingarray-like

BLdrilling: ‘double[::1]’

Gwx, Gwyarray-like

Rows of 18 values giving \(w_{,x}\) and \(w_{,y}\), used by the geometrically nonlinear terms.

.. note:: Every one of these operators is constant over the triangle,

because all three shape functions are linear, so each is formed once per element and no quadrature loop appears anywhere in this module except for the drilling term.

BLdrilling

BLdrilling: ‘double[::1]’

BLexx

BLexx: ‘double[::1]’

BLeyy

BLeyy: ‘double[::1]’

BLgxy

BLgxy: ‘double[::1]’

BLgxz

BLgxz: ‘double[::1]’

BLgyz

BLgyz: ‘double[::1]’

BLkxx

BLkxx: ‘double[::1]’

BLkxy

BLkxy: ‘double[::1]’

BLkyy

BLkyy: ‘double[::1]’

Gwx

Gwx: ‘double[::1]’

Gwy

Gwy: ‘double[::1]’

KC0ve

KC0ve: ‘double[::1]’

KCNLve

KCNLve: ‘double[::1]’

finte

finte: ‘double[::1]’

ue

ue: ‘double[::1]’

xe

xe: ‘double[::1]’