Clamped and face-clamped edges

A clamped edge of the third-order shear deformation theory (TSDT) removes the slope of \(w\) normal to the edge together with the rotations, x1wr = 0 and x1phix = x1phiy = 0 on \(x = 0\), such that the cross-section stays planar and unrotated; an edge clamped by its faces, \(u = v = 0\) at \(z = d \mp h/2\) and \(w = 0\), leaves the section free to warp and is prescribed with the keyword z of Shell.add_distr_pd_fixed_x(), see Boundary conditions. The square plates below, \(h/a = 0.1\), are clamped at \(x = 0, a\), free at \(y = 0, b\), and loaded by a uniform pressure:

def plate(model, offset=0., m=10, n=10):
    r"""Square isotropic plate, `h/a = 0.1`, all the flags ``1``, under a
    uniform pressure"""
    s = Shell(model=model, a=A, b=B, stack=[0.], plyt=H, laminaprop=(E, NU),
              offset=offset, m=m, n=n)
    for edge in ('x1', 'x2', 'y1', 'y2'):
        for field in FIELDS:
            setattr(s, edge + field, 1.)
            setattr(s, edge + field + 'r', 1.)
    s.add_pressure_load(Q)
    return s

def clamp_by_flags(s, edge):
    r"""Clamped edge by the boundary condition flags: `u, v, w` and the
    rotations, plus the slope of `w` normal to the edge for the CLPT and the
    TSDT"""
    for field in FIELDS:
        setattr(s, edge + field, 0.)
    if 'fsdt' not in s.model:
        setattr(s, edge + 'wr', 0.)

def face_clamp(s, x, h, k=1., cte=True):
    r"""Face clamp of the edge at ``x``: `u = v = 0` at the faces `z = d \mp
    h/2` and `w = 0`, prescribed with :meth:`.Shell.add_distr_pd_fixed_x`
    and the keyword ``z``"""
    zero = lambda y: 0.
    for z in (s.offset - h/2, s.offset + h/2):
        s.add_distr_pd_fixed_x(x, ku=k, kv=k, funcu=zero, funcv=zero, z=z,
                               cte=cte)
    s.add_distr_pd_fixed_x(x, kw=k, funcw=zero, cte=cte)

For the CLPT and the FSDT the face clamp is the clamp of the flags:

@pytest.mark.parametrize('model', ['plate_clpt_donnell', 'plate_fsdt_donnell'])
def test_face_clamp_equals_clamp_for_clpt_and_fsdt(model):
    r"""For the CLPT and the FSDT `u(z)` is linear in `z`: `u = v = 0` at
    both faces is the clamp of the flags"""
    s = plate(model)
    for x in (0., A):
        face_clamp(s, x, H)
    _, c = static(s)
    t = plate(model)
    clamp_by_flags(t, 'x1')
    clamp_by_flags(t, 'x2')
    _, ct = static(t)
    xs, ys = np.meshgrid(np.linspace(0, A, 7), np.linspace(0, B, 5))
    f = fields(s, c, xs.ravel(), ys.ravel())
    ft = fields(t, ct, xs.ravel(), ys.ravel())
    for k in FIELDS:
        assert np.allclose(f[k], ft[k], rtol=0, atol=1e-8*np.abs(ft['w']).max())

For the TSDT, \(u(z) = 0\) prescribed at a third coordinate gives the planar clamp of the flags, which the former clamp with the slope of \(w\) free does not:

@pytest.mark.parametrize('offset', [0., 0.3*H])
def test_planar_clamp_is_u_zero_at_three_coordinates(offset):
    r"""For the TSDT `u(z)` depends on three constants, `u`, `\phi_x` and
    `w_{,x}`: `u = v = 0` at three coordinates `z`, here the faces and `z =
    d + h/4`, gives `u(z) = v(z) = 0` for every `z`, the clamp of the flags
    ``x1u = x1v = x1w = x1phix = x1phiy = x1wr = 0``; the mid-surface `z =
    d` would not do, `u(d) = 0` following from the faces. Without ``x1wr =
    0`` the edge is softer"""
    model = 'plate_tsdt_donnell'
    s = plate(model, offset)
    zero = lambda y: 0.
    for x in (0., A):
        face_clamp(s, x, H)
        s.add_distr_pd_fixed_x(x, ku=1., kv=1., funcu=zero, funcv=zero,
                               z=offset + H/4)
    _, c = static(s)
    t = plate(model, offset)
    clamp_by_flags(t, 'x1')
    clamp_by_flags(t, 'x2')
    _, ct = static(t)
    xs, ys = np.meshgrid(np.linspace(0, A, 7), np.linspace(0, B, 5))
    f = fields(s, c, xs.ravel(), ys.ravel())
    ft = fields(t, ct, xs.ravel(), ys.ravel())
    for k in FIELDS:
        assert np.allclose(f[k], ft[k], rtol=0, atol=1e-8*np.abs(ft['w']).max())
    # the former clamp of the shear deformation theories, with w,x free
    t.x1wr = t.x2wr = 1.
    _, cf = static(t)
    assert center_w(t, cf) > 1.05*center_w(t, ct)

The face clamp lies between the planar clamp and the simple support, here with the central deflections \(w D/(q a^4)\) = 0.002901 (planar clamp), 0.002953 (face clamp) and 0.013459 (simple support): the displacements of both faces vanish along the edges, whereas the section warps, with \(\phi_x = w_{,x}/2\):

def test_face_clamp_tsdt_thick_plate():
    r"""TSDT plate, `h/a = 0.1`, clamped at `x = 0, a` and free at `y = 0,
    b`: the face clamp is softer than the planar clamp and stiffer than the
    simple support; `u` and `v` vanish at both faces along the clamped
    edges, whereas the section warps, `\phi_x = w_{,x}/2`, i.e. `\gamma_{xz}
    = 3 \phi_x`, with `u(z) \neq 0` between the faces"""
    model = 'plate_tsdt_donnell'
    s = plate(model)
    for x in (0., A):
        face_clamp(s, x, H)
    _, c = static(s)
    planar = plate(model)
    clamp_by_flags(planar, 'x1')
    clamp_by_flags(planar, 'x2')
    _, cpl = static(planar)
    ss = plate(model)
    for edge in ('x1', 'x2'):
        for field in ('u', 'v', 'w', 'phiy'):
            setattr(ss, edge + field, 0.)
    _, css = static(ss)
    w_face = center_w(s, c)
    w_planar = center_w(planar, cpl)
    w_ss = center_w(ss, css)
    assert w_planar < 0.99*w_face
    assert w_face < 0.99*w_ss
    # central deflections w D/(q a^4) of the documentation
    D = E*H**3/(12*(1 - NU**2))
    for w, ref in ((w_planar, 0.002901), (w_face, 0.002953), (w_ss, 0.013459)):
        assert np.isclose(w*D/(Q*A**4), ref, rtol=2e-4)
    ys = np.linspace(0, B, 11)
    for x in (0., A):
        f = fields(s, c, np.full_like(ys, x), ys)
        # the face displacement inside the plate
        scale = np.abs(uv_at(s, fields(s, c, np.full_like(ys, A/4), ys),
                             H/2, H)[0]).max()
        for z in (-H/2, H/2):
            u, v = uv_at(s, f, z, H)
            assert np.abs(u).max() < 1e-8*scale
            assert np.abs(v).max() < 1e-8*scale
        assert np.abs(f['w']).max() < 1e-8*w_face
        # warping: phi = w,x/2 and u(z) = (w,x/2) zeta (1 - 4 zeta^2/h^2)
        assert np.allclose(f['gxz'], 3*f['phix'], rtol=1e-6, atol=0)
        u, _ = uv_at(s, f, H/4, H)
        assert np.abs(u).max() > 1e3*np.abs(uv_at(s, f, H/2, H)[0]).max()
        assert np.allclose(u, f['gxz']/3*H/4*(1 - 4*(H/4)**2/H**2),
                           rtol=1e-6, atol=0)

Both clamps are defined by the physical section, and the results do not depend on the reference surface of the laminate:

def test_face_clamp_does_not_depend_on_the_reference_surface():
    r"""The faces are physical surfaces: the face-clamped homogeneous plate
    gives the same deflection for the reference surface at the
    mid-surface, `d = 0`, and at `d = 0.3 h` or `d = h/2`, as the planar
    clamp of the flags"""
    model = 'plate_tsdt_donnell'
    xs, ys = np.meshgrid(np.linspace(0, A, 5), np.linspace(0, B, 5))
    res = {}
    for offset in (0., 0.3*H, 0.5*H):
        s = plate(model, offset)
        for x in (0., A):
            face_clamp(s, x, H)
        _, c = static(s)
        t = plate(model, offset)
        clamp_by_flags(t, 'x1')
        clamp_by_flags(t, 'x2')
        _, ct = static(t)
        res[offset] = (fields(s, c, xs.ravel(), ys.ravel())['w'],
                       fields(t, ct, xs.ravel(), ys.ravel())['w'])
    for offset in (0.3*H, 0.5*H):
        for i in (0, 1):
            assert np.allclose(res[offset][i], res[0.][i], rtol=1e-7, atol=0)
    assert res[0.][1].max() < 0.99*res[0.][0].max()

With the penalty method, the kernels of the prescribed displacements are given the same z:

@pytest.mark.parametrize('offset', [0., 0.3*H])
def test_face_clamp_penalty_and_null_space_agree(offset):
    r"""Face clamp at `x = 0` and the faces of `x = a` pushed by `\delta`,
    `u(d \mp h/2) = -\delta`, with the pressure: the penalty method, with the
    kernels of the prescribed displacements at the same `z`, converges to
    the null-space method"""
    model = 'plate_tsdt_donnell'
    delta = 1e-5

    def build(k):
        s = plate(model, offset)
        face_clamp(s, 0., H, k=k)
        for z in (offset - H/2, offset + H/2):
            s.add_distr_pd_fixed_x(A, ku=k, kv=k, funcu=lambda y: -delta,
                                   funcv=lambda y: 0., z=z)
        s.add_distr_pd_fixed_x(A, kw=k, funcw=lambda y: 0.)
        s._rebuild()
        return s

    s = build(1.)
    _, c = static(s)
    k = 1e4*E
    p = build(k)
    md = MultiDomain([p], conn=[])
    size = md.get_size()
    kCp = 0
    for x in (0., A):
        for z in (offset - H/2, offset + H/2):
            kCp = kCp + fkCld_xcte(k, k, 0., p, x, size, 0, 0, z=z)
        kCp = kCp + fkCld_xcte(0., 0., k, p, x, size, 0, 0)
    K = md.calc_kC() + finalize_symmetric_matrix(kCp)
    cp = solve(K, md.calc_fext(), silent=True)
    xs, ys = np.meshgrid(np.linspace(0, A, 6), np.linspace(0, B, 5))
    f = fields(s, c, xs.ravel(), ys.ravel())
    fp = fields(p, cp, xs.ravel(), ys.ravel())
    for key in ('u', 'w'):
        assert np.allclose(fp[key], f[key], rtol=0,
                           atol=1e-3*np.abs(f[key]).max())
    # the faces of x = a are at -delta
    f = fields(s, c, np.full(5, A), np.linspace(0, B, 5))
    for z in (offset - H/2, offset + H/2):
        assert np.allclose(uv_at(s, f, z, H)[0], -delta, rtol=1e-8)