Internal sites of actuation and activation in thin elastic films and membranes of finite thickness

Functionalized thin elastic films and membranes frequently feature internal sites of net forces or stresses. These are, for instance, active sites of actuation, or rigid inclusions in a strained membrane that induce counterstress upon externally imposed deformations. We theoretically analyze the geometry of isotropic, flat, thin, linearly elastic films or membranes of finite thickness, laterally extended to infinity. At the mathematical core of such characterizations are the fundamental solutions for localized force and stress singularities associated with corresponding Green's functions. We derive such solutions in three dimensions and place them into the context of previous two-dimensional calculations. To this end, we consider both no-slip and stress-free conditions at the top and/or bottom surfaces. We provide an understanding for why the fully free-standing thin elastic membrane leads to diverging solutions in most geometries and compare these situations to the truly two-dimensional case. A no-slip support of at least one of the surfaces stabilizes the solution, which illustrates that the divergences in the fully free-standing case are connected to global deformations. Within the mentioned framework, our results are important for associated theoretical characterizations of thin elastic films, whether supported or free-standing, and of membranes subject to internal or external forces or stresses.


I. INTRODUCTION
Under various circumstances, we meet thin elastic films and membranes being exposed to the action of internal forces and stresses.One obvious situation arises when such membranes are mechanically reinforced, for instance, by the inclusion of relatively firmer fibers to increase the overall stability under loading 1 .When the overall system is stretched, mechanical counterstresses are generally exerted by these fibers on the elastic environment.
The same is true for any rigid inclusion upon functionalization of thin elastic films or membranes.A specific example is given by elastic membranes containing rigid magnetic or magnetizable inclusions 2 .The action of loudspeakers can be mimicked in such cases by directly exerting magnetic forces on the membrane itself, without additional external mechanical components of stimulation.
Situations of supported elastic films on a substrate are encountered frequently during measurements of atomic force microscopy (AFM).In the field of soft matter, scales of the probes are given by lengths up to micrometers 3,4 .The tip used for this purpose directly exerts forces onto the film.Depending on the situation, one can consider such forcing as internal, with the AFM tip pressed into the membrane.
On the theoretical side, there have been some recent observations concerning two-dimensional elastic systems as representations of free-standing thin elastic membrane systems.Net forces acting in an in-plane direction within infinitely extended two-dimensional systems lead to a formal logarithmic divergence of the induced displacement field as a function of the distance from the force center 13 .The two-dimensional elastic membrane, even if infinitely extended, does not offer sufficient resistance against a net force.Instead, the whole two-dimensional system will get displaced.This is in contrast to threedimensional materials, which do not experience such diverging displacements in bulk in response to persistent point-like force centers 14 .
The divergence in two dimensions naturally cancels if the overall force on the elastic membrane vanishes 15 .Particularly, this is the case when so-called symmetric force dipoles 16 act on the material.Similarly, it vanishes in the presence of, for example, a no-slip boundary 17,18 .This means that, no matter how distant a boundary is, it always has a decisive influence on the displacement field resulting from a net force in the two-dimensional membrane system.
Here, we employ a two-dimensional Fourier transformation technique to investigate thin elastic films.Yet and notably, we here do not confine ourselves to strictly two-dimensional systems.Instead, we consider elastic films and membranes of finite thickness, that is, of finite extension in the third dimension.The twodimensional Fourier transformation technique has been widely utilized in elucidating the behavior of elastic sheets and deformable membranes.Its fundamental concept revolves around harnessing the translationally invariant symmetry along horizontal planes.The approach involves transformations of the partial differential equations governing membrane displacements or fluid flowsparticularly in fluid dynamics -into remaining ordinary differential equations with respect to the out-of-plane coordinate.Previously, this technique has been effectively used to describe in low-Reynolds-number hydrodynamics the behavior near walls and interfaces [19][20][21][22][23] .The same applies to the behavior of elastic membranes undergoing shear and/or bending deformation modes [24][25][26][27][28][29][30][31][32] .
When in the present work we consider flat elastic membranes or thin films that are infinitely extended in the lateral direction, yet of finite thickness in the normal, third direction, the situation is intermediate between genuinely two-dimensional and three-dimensional, space-filling, bulk systems.An important question from the conceptual perspective is the role of force-induced divergences in this intermediate situation.How does it compare to the results for two-and three-dimensional systems?
To address this question, the remainder of the paper is structured as follows.In Sec.II, we present the twodimensional Fourier transformation technique applied to the governing equations of elasticity, along with solutions for various boundary conditions.Specifically, these are no-slip and/or stress-free conditions at the top and/or bottom surfaces.We explore both monopole and dipole types of singularity and address the solutions in terms of resulting Green's functions.Our results are discussed in Sec.III, where we also delineate the conditions under which the solutions converge and diverge.The paper concludes in Sec.IV.Complex expressions pertaining to the solution are presented in tables within Appendix A, while additional background concerning the calculations is included in Appendix B.

II. MATHEMATICAL FORMULATION
The displacement field u(r) of the material points in the elastic medium is governed by the Navier-Cauchy equation with µ and λ denoting the shear modulus and Lamé parameter, respectively, and f (r) the force density field acting on the elastic material.Both coefficients are related to each other via where ν ∈ (−1, 1/2] is the Poisson ratio.The components of the stress tensor read In this expression, δ ij represents the Kronecker delta and equal to the components of the unit matrix.ε ij = (∇ i u j + ∇ j u i )/2 includes the components of the linearized strain tensor.
To derive solutions for the displacement field, we utilize a two-dimensional Fourier transformation in both the xand y-direction.This transformation is defined as follows.The forward two-dimensional Fourier transformation of a function g(ρ, z) is denoted as g(k, z) and is expressed as where k represents the wavevector in the two-dimensional plane of the Fourier transformation.Here, ρ = (x, y) denotes the projection of the position vector into the xy-plane.We do not perform a Fourier transformation with respect to the z-component.Similarly, we define the inverse two-dimensional Fourier transformation as Additionally, we introduce the wavenumber k = |k|, representing the magnitude of the two-dimensional wavevector, and we define the unit vector k = k/k.
To address the elastic equations, we adopt a cylindrical coordinate system and represent the unit wavevector as k = cos ϕ x + sin ϕ ŷ.Additionally, we define the twodimensional unit vector t in the xy-plane perpendicular to k such that t = sin ϕ x − cos ϕ ŷ.Consequently, the Fourier-transformed displacement field and force density field can be projected onto the basis consisting of the unit vectors k and t.The longitudinal and transverse components of the Fourier-transformed displacement field are defined as u l = u • k and u t = u • t, respectively 33 .We refer to the z-component u z of the Fourier-transformed displacement field as the normal component.Similarly, the longitudinal, transverse, and normal components of the Fourier-transformed force density are denoted as f l , f t , and f z , respectively.
Next, in the new orthogonal basis, we derive the twodimensional Fourier-transformed elastic equations.We find that the transverse component of the displacement field is completely independent of the longitudinal and normal components, The longitudinal and normal components of the displacement field are governed by a system of secondorder differential equations, Although our focus lies in seeking solutions to the elastic equations involving force singularities concentrated at one point, we remark that Eqs. ( 6) and (7) apply to any arbitrary force density f (r), given that its Fourier transformation is known.

A. Solution in a bulk elastic medium
We begin by revisiting the scenario involving a force singularity within an infinitely extended, threedimensional, bulk elastic medium.Thus, we examine the effect of a force density represented by a Dirac delta function located at the origin of our Cartesian coordinate system, f (r) = F δ 3 (r).
The characteristic polynomial of the homogeneous equation features two distinct roots, ±k.As a consequence, the solution for the transverse component of the displacement field takes the form u t = Ae −k|z| , where the constant A is obtained from the jump in the first z-derivative due to the Dirac delta function.Specifically, This leads to The characteristic polynomials associated with the homogeneous parts of Eqs.(7) similarly yield two distinct roots ±k, each with a multiplicity of 2. Consequently, the longitudinal and normal components of the displacement field take on the forms u , where the superscripts 'plus' and 'minus' correspond to the regions above and below the singularity.Continuity of the displacements across z = 0 implies α Further relations follow by introducing these functional forms into Eqs.(7) away from the singularity position.Moreover, the Dirac delta function at z = 0 implies the jump conditions Upon solving the resulting system of linear equations, we acquire the final expressions

B. Image solution technique
Our goal is to quantify the static deformational response of an elastic membrane of finite thickness H to an applied static force acting at position h ∈ (−H/2, H/2), see Fig. 1(a).Initially, in the undeformed state, the membrane features two flat, parallel, infinitely extended (top and bottom) surfaces at z = ±H/2.We examine three scenarios concerning the possible confinement of these surfaces when a force singularity is applied between these surfaces.These are (i) no-slip conditions both at the top and bottom boundaries, see Fig. 1(b), (ii) one noslip condition (top boundary) and one stress-free condition (bottom boundary), see Fig. 1(c), and (iii) stress-free conditions both at the top and bottom boundaries, see Fig. 1(d).
Utilizing the bulk solution derived in Fourier space in Sec.II A, we express the result as a combination of the bulk term and an additional contribution necessary to meet the prescribed boundary conditions.This technique is commonly known as the image solution technique.Specifically, with u * representing image solution.
We remark that the position of the force center in the bulk elastic medium considered in Sec.II A was located at the origin.Obtaining the corresponding expressions for the force acting at z = h is straightforward in bulk through simple translation.To this end, z is replaced by z − h in the expressions for u t , u l , and u z in Eqs. ( 9) and (11).
The image solution u * satisfies the homogeneous parts of Eqs. ( 6) and (7).It can be expressed in the form for the transverse component and for the longitudinal and normal components, where The wavenumber-dependent coefficients c ± 1 , c ± 2 , and c ± 3 will subsequently be determined based on the boundary conditions specified at z = ±H/2, see Fig. 1.
At each no-slip surface, we enforce the vanishing components of the displacement field, namely u = 0, ensuring u t = u l = u z = 0 at these boundaries.For a stress-free surface, we prescribe σ • ẑ = 0, implying the boundary conditions Particularly, we note that, if ν = 1/2, Eq. (16c) becomes the incompressibility condition for the elastic medium.Subsequently, we express lengths in dimensionless units by normalizing them with respect to the thickness of the membrane H.Effectively we thus set a unit membrane thickness.In this way, the resulting expressions simplify.
Imposing the relevant boundary conditions on the top and bottom surfaces, we obtain a system of six linear equations for each combination of boundary conditions.
Their solutions yield the expressions for the remaining unknown coefficients c ± j , j ∈ {1, 2, 3}.Due to their complexity and lengthiness, and since it is straightforward to formulate the system of equations using computer algebra systems, we do not include the details of evaluation here.
First, as previously indicated, the transverse component of the displacement field, see Eq. ( 13), remains unaffected by normal and longitudinal forces, see Eq. ( 8).This results in a set of equations for c ± 1 when imposing the relevant boundary conditions.For both types of boundary conditions, the expressions for the coefficients c ± 1 follow in the form Here, the specific expressions for w 0 , w 1 , and d are determined by the specified set of boundary conditions.As depicted in Fig. 1, we introduce the abbreviations 2NOS (two no-slip walls, NOS-STF (bottom no-slip and top stress-free), and 2STF (two stress-free surfaces).We find that d = µ e 2k − 1 for 2NOS and 2STF, and d = µ e 2k + 1 for NOS-STF.In addition, (w 0 , w 1 ) = (1, −1), (−1, ±1), and (1, 1) for 2NOS, NOS-STF, and 2STF, respectively.
The coefficients c ± 2 and c ± 3 entail more complex expressions, as they specify both the longitudinal and normal components of the displacement field, see Eqs. ( 14) and (15).This results from their combination in the coupled differential equations, see Eqs. (7).For ease of reference, we introduce the abbreviations together with δ ± = 1 ± h and η ± = 1 ± 2h.Imposing the different sets of boundary conditions, the expressions for the coefficients c ± 2 and c ± 3 determining the expressions for the longitudinal and normal components of the image displacement field, see Eqs. ( 14) and (15), for all sets of boundary conditions can be formulated as Here, the specific expressions for the coefficients a ± n , b ± n , g ± n , and q ± n together with the denominator D depend on the specific boundary conditions.They are provided in Tab.A1 of Appendix A.
We now define the displacement Green's function as It is composed of the Green's function associated with the bulk displacement field, G ∞ (r), and the Green's function associated with the image displacement field, G * (r), Based on Eqs. ( 9) and ( 11), we deduce that in Fourier space the components of the Green's function associated with the bulk displacement field for a force acting at a distance h above the center plane are These expressions can simply be read off from Eqs. ( 9) and ( 11), recalling the vertical shift due to the position of the force center at z = h instead of z = 0 and regarding Eq. ( 18).In addition, it follows that Next, in a similar way, we infer the the corresponding expressions for the Fourier-transformed image solution of the Green's function from the previous results for the image displacement field.Utilizing Eqs. ( 13) and ( 17), we may express the tt-component of the image Green's function for all sets of boundary conditions as where Moreover, it follows from Eqs. ( 14), (15), and ( 19) that for i, j ∈ {l, z} for all sets of boundary conditions.We here defined the abbreviations The expressions for w 0 , w 1 , d, a ± n , b ± n , g ± n , q ± n , and D again depend on the specific set of boundary conditions and are listed below Eq. ( 17) and in Tab.A1 of Appendix A.
It is beneficial to have available the expressions of the Green's functions in the standard Cartesian coordinate system as well.To this end, we invert the transformation outlined below Eq. ( 5).This implies and G zy = G zl sin ϕ.Moreover, G yx = G xy , which is dictated by the geometric symmetry in the xy-plane.
We next proceed to the inverse Fourier transformation back to real space 34 .Given that the expressions in Fourier space are formulated utilizing polar coordinates (k, ϕ), we anticipate that the expressions in real space, after inverse Fourier transformation, can likewise be expressed in polar coordinates (ρ, θ).Here, ρ := x 2 + y 2 represents the radial distance from the origin, while θ := arctan(y/x) denotes the angle formed between the positive x-direction and the radial direction of the considered position in space.A comprehensive collection of associated formulas can be found in Ref. 35.In Appendix B, we include an overview of the pertinent formulas utilized in the present computation.For later reference and for the ease of computing the inverse Fourier transformation, we present the expressions of the Green's function itself in Cartesian coordinates.Conversely, its dependence on the position for convenience is maintained in polar coordinates ρ and θ.
First, we list the components of the Green's function associated with the z-component of the applied force.A force only featuring a z-component is oriented perpendicular to the surfaces, implying an axisymmetric geometry.The inverse Fourier transformation yields In our notation, J n refers to the Bessel function of the first kind of order n 36 .It follows that G xz = G ρz cos θ and G yz = G ρz sin θ in Cartesian coordinates.
The remaining components of the Green's function are calculated in analogy to that.We introduce the notation Then, we obtain Moreover, Again, G yx = G xy , which agrees with the geometric symmetry within the xy-plane.Furthermore, G zx = G zρ cos θ and G zy = G zρ sin θ, where In a bulk elastic medium, which is infinitely extended in all space directions, the Green's function reads 37 where r := |r| = ρ 2 + (z − h) 2 denotes the distance from the position of the applied force singularity.The result in Eq. ( 34) is obtained by inverse Fourier transformation of Eqs.(22).For ν = 1/2, implying σ = 1 in Eq. ( 18), representing an incompressible elastic medium that does not allow any divergence of the displacement field at any position in space, we recover the classical solution recognized in fluid mechanics as the Oseen tensor 38,39 .This tensor characterizes the behavior of an incompressible fluid under the influence of a force singularity under low-Reynolds-number conditions.There, µ assumes the role of the (dynamic) shear viscosity and u(r) takes the role of the flow field of the fluid.

C. Convergence of the Green's functions for force monopoles
We here consider the displacement fields arising in a laterally infinitely extended elastic film or membrane of finite thickness.The expressions for the Green's functions quantifying these displacements in real space in response to the localized application of a force were derived in Sec.II B, see Eqs. ( 28)- (33).They are provided by integral expressions that result from the inverse two-dimensional Fourier transformation back to real space, see Eq.( 5).We discuss the divergence of these expressions for the different types of boundary conditions.
Considering that the Green's functions exhibit a screened exponential decay as the wavenumber k approaches infinity, which results from the Bessel functions in Eqs. ( 28)- (33), problems of convergence of the integrals do not arise when k approaches infinity.Our primary concern remains the convergence of the integrals determining the Green's functions in Eqs. ( 28)- (33) as k tends towards zero.Near k = 0, the Bessel functions J n (ρk) can be approximated as (ρk/2) Thus, G must exhibit O k −1 or higher-order behavior in k near k = 0 to warrant convergence of the integrals.The cause lies with the factor k that the Green's functions G are multiplied with in polar coordinates when performing the two-dimensional inverse Fourier transformation in Eq. ( 5).
Overall, this requires an exponent of k of at least −1 in G to hinder divergence of the integrals at k = 0.
For boundary conditions involving at least one no-slip surface, that is, for the categories 2NOS and NOS-STF, we observe that the diagonal components G tt , G ll , and G zz all exhibit behavior of the order O k −1 , see Eqs. ( 34)-( 26) together with the expressions of the parameters listed in Tab.A1 of Appendix A. Moreover, the off-diagonal components G lz and G zl show behavior on the order of O (1), as equally inferred from this set of equations.Overall, this indicates convergence of the Green's function for all cases involving at least one no-slip boundary.
The situation becomes markedly different in the case of two stress-free boundaries (2STF).In this situation, from Eqs. ( 23)-( 26) together with the coefficients listed in Tab.A1 in Appendix A, we observe that G tt and G ll exhibit behavior on the order of O k −2 .Even more so, G lz and G zl show behavior of the order of O k −3 , and G zz is of the order of O k −4 .These relations indicate a non-convergent Green's function in real space when the integrals in Eqs. ( 28)-( 33) are carried out.

D. Convergence of the Green's functions for force dipoles
Having analyzed the behavior in response to singularities consisting of force monopoles, our focus shifts to singularities of force dipoles.Again, we investigate how the behavior varies for the different boundary conditions enforced at the surfaces.
First, we introduce the â-directed Green's function as where a is a unit vector and r 0 is the position of applying the force singularity.Next, we calculate the gradient ∇ 0 of the â-directed Green's function with respect to the position r 0 .In this way, we define the Green's function for the force dipole composed of â-directed forces as This force dipole extends along the direction given by the unit vector b.
In this expression, the differentiation ∇ 0 (•) with respect to the position r 0 of the force singularity stems from casting the discrete setup of two nearby force centers into a differential form.During this process, as a limit, the distance between the two centers of the point forces tends to zero.This is how we turn from the two infinitesimally separated force centers to the force dipole.In contrast to ∇ 0 (•), the regular ∇(•) is the differentiation with respect to the regular space coordinate.Due to the translational symmetry along the xy plane, it follows from Eq. ( 36) in Cartesian coordinates that ∂/∂x 0 = −∂/∂x and ∂/∂y 0 = −∂/∂y.However, in the presence of a boundary, ∂/∂h does not equate to −∂/∂z.In fact, G(r, r 0 ) in Eq. ( 35) depends on x − x 0 and y − y 0 , but not on z − z 0 , where r 0 = (x 0 , y 0 , z 0 ).Next, we define the stresslet as the symmetric part of the force dipole, expressed as Likewise, we define the rotlet as the antisymmetric part, where ĉ = â × b.
It is more convenient to derive the expressions of the force dipole in Fourier space according to Eq. ( 4).Then, we obtain recalling that the Fourier transformation is performed only in two dimensions so that the arguments are given by the two-dimensional wavevector k together with the real-space z coordinate.The minus signs in Eqs.(39a) and (39b) result from the relations ∂/∂x 0 = −∂/∂x and ∂/∂y 0 = −∂/∂y noted above.When turning to the inverse two-dimensional Fourier transformation according to Eq. ( 5), the solution is expressed in integral form Integration with respect to the variable ϕ has already been performed in the function Γ(r, k; â, b).This calculation can be supported by conventional computer algebra systems.In Cartesian coordinates, the expression for Γ(r, k; â, b) needs to be evaluated with respect to all â, b ∈ { x, ŷ, ẑ}.Each such combination of â and b implies three Cartesian components of G D (r; â, b) in Eq. (40).We present the corresponding components of G D (r; â, b) in Tab.A2 of Appendix A. Next, we investigate the convergence of the integral in Eq. ( 40) in analogy to the procedure outlined in Sec.II C. To this end, we examine whether Γ is of order k −1 or higher as k approaches 0. From the relations in Eqs. ( 35) and (39), we infer convergence for G D whenever the integrals for G converge.We found convergence for the pure force singularity in all cases that involve at least one no-slip boundary (2NOS and NOS-STF).Thus, we may assume converge of the solution under these boundary conditions also for G D .
Therefore, the remaining task is to consider the situation of two stress-free boundaries (2STF).It turns out that in this case the solution for the force dipole generally can still show divergence.For certain geometries, it converges while the solution for the pure force singularity diverges.40).Here, we refer to the scenario of a free-standing elastic film or membrane of finite thickness featuring two parallel surfaces of stress-free boundary conditions (2STF).The combinations Γ(r, k; x, ŷ) and Γ(r, k; ŷ, x) are omitted because both are of O (1) to leading order in k at k = 0 and thus imply convergent integrals in Eq. ( 40).
The leading-order behavior of Γ with respect to k near k = 0 can be determined from the expressions in Tab.A2 in Appendix A systematically using computer algebra.Upon examination, we observe that only the configurations Γ(r, k; x, ŷ) and Γ(r, k; ŷ, x) satisfy the criteria necessary for the integral in Eq. ( 40) to generally converge.Thus, only the Green's functions G D (r; x, ŷ) and G D (r; ŷ, x) for the force dipole are well defined in general.We present in Tab.I the leading-order terms in the series expansion of the remaining combinations of Γ(r, k; â, b) around k = 0 with â, b ∈ { x, ŷ, ẑ}.From there, it is obvious that no other combination of â and b provides the required leadingorder behavior of k −1 or higher so that convergence of the integral in Eq. ( 40) does not generally result.
Concerning the Green's function for the stresslet G D , see Eq. ( 37), it is evident that convergence is achieved, nevertheless, if the stresslet is positioned precisely at the center plane of the elastic film or membrane.This implies h = 0.Then, we find from Tab.I that the leading orders in k −2 cancel for Γ(r, k; x, x), Γ(r, k; ŷ, ŷ), and Γ(r, k; ẑ, ẑ).Moreover, upon symmetrization in Eq. ( 37), the leading orders in k −2 cancel for Γ(r, k; x, ẑ) + Γ(r, k; ẑ, x) and Γ(r, k; ŷ, ẑ) + Γ(r, k; ẑ, ŷ).Therefore, the leading order is shifted to at least k −1 around k = 0, thus satisfying the required criterion of convergence for the integral in Eq. (40).Besides, this implies that for the Green's function of the rotlet G R , see Eq. (38), only the configuration of ĉ = ẑ implies general convergence from the integrals in Eq. (40).

III. DISCUSSION
We touched here the subtle question of possible divergences in the displacements of elastic films and membranes when exposed to mechanical force density fields.Conversely to many previous works, we here focused on systems of finite, non-vanishing thickness.In deriving associated Green's functions, we concentrated on the response to imposed forces, see Secs.II A-II C, stresslets (symmetric contributions to force dipoles), see Sec.II D and rotlets implying torques (antisymmetric contributions to force dipoles), see Sec.II D as well.
During these investigations, we encountered possible divergences in the resulting displacement fields.Generally, we found that the presence of at least one no-slip boundary stabilizes the system so that no divergence is observed, see Sec.II C.However, divergences in the displacement field emerge for free-standing elastic membranes that are infinitely extended to the lateral sides and of finite thickness in the normal direction.
That is, stressfree boundary conditions apply on both the top and bottom surface.The divergences are connected to small wavenumbers and thus large length scales, see Sec.II C. In fully three-dimensional, bulk systems such divergences are not observed 14 .
We first compare to the situation of two-dimensional elastic systems, that is, infinitely thin, flat elastic membranes.In that case, logarithmic divergences are found for displacement fields in response to net forces 40 , but not in response to force dipoles 15 , all applied along in-plane directions.
In our considerations, the two-dimensional situation is recovered when applying forces and force dipoles that are oriented in the xy-plane and positioned in the center plane at z = 0 (h = 0), see Fig. 1.In that case, in strictly two dimensions, net forces lead to diverging displacement fields, but force dipoles do not 15 .Indeed, we recover this convergence for force dipoles, if G D (r; x, x) and G D (r; ŷ, ŷ) are evaluated for positioning the force dipoles in the center plane.This is observed from Eq. ( 40) together with entries 1 and 3 in Tab.I.For h = 0, the leading order in k around k = 0 for Γ(r, k; x, x) and Γ(r, k; ŷ, ŷ) shifts to at least k −1 , thus ensuring convergence of the integral in Eq. (40), see also Sec.II D. Along the same lines, from the orders of Γ mentioned in the caption of Tab.I, we infer that the expressions of G D (r; ŷ, x) and G D (r; ŷ, x) are always well defined.Therefore, the symmetric positioning of force dipoles at h = 0 recovers the situation of strictly two-dimensional considerations 15 , where all resulting displacement fields remain finite.Conversely, off-center positioning of force dipoles, implying h ̸ = 0, can lead to diverging displacement fields, signaled by the nonconverging integral in Eq. ( 40) when k approaches 0 as highlighted by the expressions in Tab.I.
Additionally, we note that displacements according to stresslets corresponding to the Green's functions G S (r; x, ẑ) and G S (r; ŷ, ẑ) are always well defined, also for off-center positioning at h ̸ = 0.The cause lies with the symmetrization in Eq. (37).In combination, as is obvious from entries 2 and 5 in Tab.I, the leading order of the symmetrized Γ shifts to at least order k −1 , thus ensuring convergence of the integral in Eq. ( 40) around k = 0.The analogous behavior is observed when considering entries 4 and 6 in Tab.I. Therefore, the displacements associated with G S (r; x, ẑ) and G S (r; ŷ, ẑ) remain finite.The remaining question is, how the diverging displacement fields can arise even in the case of symmetric force dipoles.
In fact, considering elastic films or membranes in three spatial dimensions, a bending contribution can generally be involved, as soon as the force dipole is not placed to the center plane at h = 0. We recall absence of anchoring at infinity.Consequently, we expect bending modes associated with infinite magnitudes of normal displacements at infinite distance from the center of the force dipole.Moreover, antisymmetric parts of force dipoles that involve the normal direction lead to overall rotations of the membrane, and therefore again to diverging displacement fields at infinity.Once more, the cause lies with the absence of anchoring of the lateral ends at infinity.We refer to the remark at the end of Sec.II D. Only for in-plane rotations, we found converging integrals in Eq. ( 40), implying finite displacements.Thus, again, for the scenario of in-plane rotations, we recover the twodimensional situation of well-defined displacement fields 15 .
Finally, we remark that in the context of linear elasticity, small strains do not necessarily equate to small displacements.Rather, they signify a small symmetrized gradient in displacement, that is, strain.
In other words, we may reside well in the regime of linear elasticity (small strains) but still encounter divergences in the displacement field (for instance, at infinity in an infinitely extended system).Nonlinear elastic response may suppress such divergences.To our knowledge, there is limited analytical exploration of nonlinear elasticity theory concerning singularities within the framework of free-standing elastic films and membranes.Investigating the challenging nonlinear problem in future research could provide further insight.

IV. CONCLUSIONS
In summary, we here consider the mechanical displacements and deformations of linearly elastic, initially flat films and membranes exposed to forces and force dipoles.The systems are infinitely extended to the sides, yet of finite, non-vanishing thickness in their normal direction.We address supported films and membranes subject to no-slip boundary conditions on at least one of the two surfaces (top and/or bottom).In these situations, the solutions for resulting displacements in response to both net applied forces and/or stresslets are well defined.For free-standing elastic films and membranes, that is, for stress-free boundary conditions applying to both the top and bottom surface, we observe diverging solutions for most geometries.These divergences can be related to global deformations and large-scale displacements of the whole (infinitely extended) system.When positioning the force dipoles at the center plane of the elastic sheet, our expressions are in line with corresponding considerations for genuinely two-dimensional systems 15,18 .
Our results are significant for any theoretical approach that addresses the deformation of thin, isotropic elastic films and membranes of infinite lateral extension within the linear regime.Particularly, we find that, unless for very symmetric situations of stresslets included symmetrically at the center plane of the membrane, infinitely extended free-standing membranes cannot serve as immediate theoretical model systems without further precaution.In contrast to genuinely two-dimensional treatments, even off-center stresslets imply divergence of the solutions.Therefore, some stabilization of the system, presumably by introduction of lateral boundaries and/or counterforces, is necessary in this case.
Overall, the analysis in this study relies on linear elasticity, primarily due to the employed formalism of Green's functions and the principle of superposition.Yet, linear elasticity does not always accurately capture the observed response of various materials, particularly in the fields of soft matter and biology.Defining parameters like Young's moduli for biological tissues poses significant challenges, especially when their stiffness is very small in the absence of loading or pre-stretch.We opt for a linear framework as an initial approach to understand the behavior within the regime of small deformations and to facilitate analytical advancements.Addressing nonlinearities represents an important avenue for future investigations, potentially through perturbative methods.Still, we note that the observed divergences in displacement do not necessarily imply that the regime of linear elasticity is violated.Since strains follow as symmetrized gradients in displacements, strains can be small in magnitude even if displacements diverge in infinitely extended systems.
Finally, it will be interesting to analyze experimental data in view of our results.Particularly, this concerns the analysis and quantification of observed displacement fields.Naturally, real-space experimental systems are always of finite extent, in contrast to our mathematical picture.
It is crucial to acknowledge that the role of thermal fluctuations may become significant with increasing extension.Potentially, they obscure certain aspects of the described effects, especially in large systems.Introducing lateral boundaries into our investigations, such as clamping of elastic films and membranes, comparison with experiments is facilitated.Mathematically speaking, such additional boundary conditions play a pivotal role.
Our study offers valuable insight in this regard, aiming to draw attention to this aspect.We infer that in experiments the role of lateral clamping is dominant for free-standing systems of stress-free top and bottom surfaces, as they are the ones that maintain displacements finite.Eventually, if the system is stabilized by a no-slip boundary condition at the bottom surface, while stress-free conditions apply at the top surface, we recover the situation of atomic force microscopy of thin films 3,4 .It is our plan to provide further theoretical tools in this framework.Table A1.Expressions for a ± n , b ± n , g ± n , and q ± n , where n ∈ {0, 1, 2, 3}, defining the coefficients c ± 2 and c ± 3 introduced in Eqs.(19).These coefficients specify the longitudinal and normal components of the two-dimensionally Fourier-transformed image solutions for the displacement field, u * l and u * z , respectively, see Eqs. ( 14) and (15).The expressions vary depending on the boundary conditions imposed at the top and bottom surface.We address geometries of two no-slip boundaries (2NOS), two stress-free boundaries (2STF), and the combination of one bottom no-slip and one top stress-free boundary (NOS-STF).For convenience, the definitions of σ, δ±, and η±, see Eq. ( 18) and thereafter, are repeated on the bottom right.Table A2.Expressions for the various combinations of the function Γ(r, k; â, b), where â, b ∈ { x, ŷ, ẑ}, are listed componentwise by projecting them for each combination on the Cartesian unit vectors.The function Γ upon integration defines the Green's function of the force dipole in real space, see Eq. (40).We here express Γ in the system of polar coordinates (ρ, θ).In these expressions, Gij represents the two-dimensionally Fourier-transformed Green's function quantifying the displacements induced in the elastic system in response to a force acting at position (0, 0, h) on the elastic film or membrane of finite thickness, see Eqs. ( 21)- (26).When dealing with a radially symmetric function f (k) independent of the azimuthal angle ϕ, only the term n = 0 persists in the series described by Eq. (B1).As a result, the inverse Fourier transformation simplifies to

2NOS
This specific formulation applies to the component G zz of the Green's function, as indicated by Eq. ( 29).For the other components, we observe that only terms of n ∈ {0, 1, 2} remain for the force monopole and terms of n ∈ {0, 1, 2, 3} for the force dipole.The entries of Tab.A2 were derived using Eqs.(B1)-(B3).

Figure 1 .
Figure 1.(a) Illustration of the system setup.A force F acts onto the thin elastic film or membrane of finite thickness H at the position (0, 0, h).The bottom surface of the film or membrane is located at z = −H/2, the top surface at z = H/2.For this situation, three different sets of boundary conditions are addressed.These are (b) no-slip boundary conditions at both the top and bottom surface (2NOS), (c) the bottom boundary exhibiting a no-slip condition while the top boundary remains stress-free (NOS-STF), and (d) stress-free boundary conditions at both the top and bottom surfaces (2STF).Noslip boundary conditions imply vanishing displacements u = 0, while stress-free conditions imply σ • ẑ = 0.
TableI.Leading-order expressions in k of the function Γ(r, k; â, b) around k = 0 with â, b ∈ { x, ŷ, ẑ} up to O k −3 .The function Γ determines the Green's function GD for the force dipole, see Eq. (