Prescribed displacements#
The prescribed displacements of Shell.add_point_pd(),
Shell.add_distr_pd_fixed_x() and Shell.add_distr_pd_fixed_y()
are imposed by a MultiDomain, exactly with
pd_method='null-space' or approximately with pd_method='penalty',
see Boundary conditions and Prescribed displacements imposed exactly. The plates below
have the in-plane displacements free at the edges and are compressed by an
end shortening:
def plate(a, b, h=2e-3, m=6, n=6, model='plate_clpt_donnell', x0=0., y0=0.,
stack=(0, ), laminaprop=(E, NU), membrane_free=True):
r"""Plate with the in-plane displacements free at all the edges, ``w``
restrained at the edges"""
s = Shell(a=a, b=b, m=m, n=n, x0=x0, y0=y0, model=model,
stack=list(stack), plyt=h/len(stack), laminaprop=laminaprop)
if membrane_free:
for edge in ('x1', 'x2', 'y1', 'y2'):
setattr(s, edge + 'u', 1.)
setattr(s, edge + 'v', 1.)
#NOTE also the rotations of the shear deformation theories
setattr(s, edge + 'phix', 1.)
setattr(s, edge + 'phiy', 1.)
return s
def end_shortening(p, delta, x1=True, x2=True, cte=False):
r"""``u = 0`` at ``x = 0`` and ``u = -delta`` at ``x = a``"""
if x1:
p.add_distr_pd_fixed_x(0., ku=1., funcu=lambda y: 0., cte=True)
if x2:
p.add_distr_pd_fixed_x(p.a, ku=1., funcu=lambda y: -delta, cte=cte)
A single shell is wrapped in an assembly without connections, conn=[],
and MultiDomain.linear_static() returns the configuration with the
prescribed displacements, here the uniform strain of the patch test:
def test_single_shell_without_conn():
r"""A single panel is wrapped in a :class:`.MultiDomain` without
connections, ``conn=[]``; ``conn=None`` raises an error, as before"""
p = plate(1., 0.6)
end_shortening(p, 1e-3)
p.add_point_pd(0., 0., 0, None, 1., 0., 0, None)
md = MultiDomain([p], conn=[], pd_method='null-space')
c = md.linear_static()
assert np.isclose(c @ np.zeros_like(c), 0)
X, Y, f = fields(p, c)
assert np.allclose(f['u'], -1e-3*X, atol=1e-15)
with pytest.raises(RuntimeError, match='No connectivity'):
MultiDomain([plate(1., 1.)], pd_method='null-space').get_Z()
A cylinder of five domains under a uniform end shortening gives the same
prebuckling state as the axial load of the corresponding reaction,
MultiDomain.reaction_resultant(), and a slightly higher critical
load; with the penalty method, the springs of fkCld_xcte() are added
by the user to the stiffness matrix, and the critical load converges to the
one of the null-space method as the penalty constant increases:
def cylinder(pd_method, **kw):
md, conns = create_cylinder(height=0.5, r=0.25, plyt=0.125e-3,
stack=[0, 45, -45, 90, -45, 45], laminaprop=LAMINAPROP,
npanels=5, m=8, n=12)
md = MultiDomain(md.panels, conns, pd_method=pd_method)
for p in md.panels:
p.x2u = 1
return md
@pytest.mark.slow
def test_cylinder_end_shortening():
r"""Cylinder of 5 panels, simply supported, compressed by a uniform
end-shortening, see ``test_cylinder_compression_lb_Nxx_from_static``
- the penalty method converges to the null-space method with the error
decreasing as `1/k`, from above, the penalty springs shortening the
cylinder less than `\delta`
- the end-shortening of the axial load `N_{xx}` gives the reaction
`2 \pi r N_{xx}`, the same pre-buckling state, and a larger critical
load, the buckling mode having `u = 0` at the loaded edge under
displacement control
"""
Nxx = -100.
md = cylinder('null-space')
for p in md.panels:
p.add_distr_load_fixed_x(p.a, lambda y: Nxx, None, None)
eig_N, _ = md.linear_buckling(num_eigvalues=3)
assert np.isclose(Nxx*eig_N[0], -47162.406, rtol=2e-3)
c_N = md.linear_static()
delta = -fields(md.panels[0], c_N)[2]['u'][0, -1]
md = cylinder('null-space')
for p in md.panels:
p.add_distr_pd_fixed_x(p.a, ku=1., funcu=lambda y: -delta, cte=False)
eig_d, modes = md.linear_buckling(num_eigvalues=3)
assert modes.shape[1] == 3
c = md.linear_static()
assert np.allclose(c, c_N, rtol=0, atol=1e-6*np.abs(c_N).max())
F = md.reaction_resultant(c, [(p, 'x', p.a) for p in md.panels],
NLgeom=False)
assert np.isclose(F[0], Nxx*2*np.pi*0.25, rtol=1e-6)
assert 1.02 < eig_d[0]/eig_N[0] < 1.05
errors = []
for factor in (1e3, 1e4):
md = cylinder('penalty')
size = md.get_size()
k = factor*md.panels[0].lam.ABD[0, 0]
kCp = 0.
for p in md.panels:
p.add_distr_pd_fixed_x(p.a, ku=k, funcu=lambda y: -delta,
cte=False)
kCp = kCp + fkCld_xcte(k, 0., 0., p, p.a, size, p.row_start,
p.col_start)
K = md.calc_kC(silent=True) + finalize_symmetric_matrix(kCp)
c = md.expand(solve(md.reduce(K), md.reduce(md.calc_fext()),
silent=True))
eig, _ = lb(md.reduce(K), md.reduce(md.calc_kG(c=c, silent=True)),
silent=True, num_eigvalues=3)
errors.append(eig[0]/eig_d[0] - 1)
assert 0 < errors[1] < 2e-4
assert 9 < errors[0]/errors[1] < 11
The geometrically non-linear analyses under displacement control are shown
in Prescribed displacements imposed exactly. Without an assembly, a point of a Shell
is supported by the penalty springs of
Shell.calc_stiffness_point_constraint(), added to the stiffness
matrix:
def test_stiffness_point_constraint_fixes_the_point():
"""A w penalty at the point where the load is applied cancels w there"""
from structsolve import solve
s = _plate(a=1., b=1.)
s.add_point_load(0.5, 0.5, 0, 0, 1.)
kC = s.calc_kC()
x, y = 0.5, 0.5
c = solve(kC + s.calc_stiffness_point_constraint(x, y, u=False, v=False,
kuvw=1e12),
s.calc_fext(), silent=True)
s.uvw(c, xs=np.array([x]), ys=np.array([y]))
w_constrained = s.fields['w'][0]
c = solve(kC, s.calc_fext(), silent=True)
s.uvw(c, xs=np.array([x]), ys=np.array([y]))
w_free = s.fields['w'][0]
assert abs(w_constrained) < 1e-3*abs(w_free)