Wettability of semispherical droplets on layered elastic gradient soft substrates

Research on the wettability of soft matter is one of the most urgently needed studies in the frontier domains, of which the wetting phenomenon of droplets on soft substrates is a hot subject. Scholars have done considerable studies on the wetting phenomenon of single-layer structure, but it is noted that the wetting phenomenon of stratified structure is ubiquitous in nature, such as oil exploitation from geological structural layers and shale gas recovery from shale formations. Therefore, the wettability of droplets on layered elastic gradient soft substrate is studied in this paper. Firstly, considering capillary force, elastic force and surface tension, the constitutive equation of the substrate in the vector function system is derived by using the vector function system in cylindrical coordinates, and the transfer relation of layered structure is obtained. Further, the integral expressions of displacement and stress of double Bessel function are given. Secondly, the numerical results of displacement and stress are obtained by using the numerical formula of double Bessel function integral. The results show that the deformation of the substrate weakens with the increase of the elastic modulus, also the displacement and stress change dramatically near the contact line, while the variation is flat when the contact radius is far away from the droplet radius.


Scientific Reports
| (2021) 11:2236 | https://doi.org/10.1038/s41598-020-80575-9 www.nature.com/scientificreports/ low adhesion of lotus leaves 31 , and the geological structure of oil exploitation is also layered. The water and ice repellent properties of layered flexible superhydrophobic substrate composed of poly (dimethylsiloxane) and zinc oxide were studied, and this study offered insights into the design of multi-structured materials 32 . So it had potential application value to study the wettability of layered structure. It is noted that most anisotropic materials can be synthesized by lamination technology, so layered isotropy can be used to approximate the wetting condition of anisotropy. For layered elastic structures, Pan et al. 33,34 expressed the displacement and stress by vector function system, and gave the transfer matrix of layered elastic structure by using the transfer matrix method, then solved the displacement and stress of each layer. In this paper, we follow Pan's work and develop their method into the wetting phenomenon, and study the deformation of elastic layered isotropic thin substrate caused by hemispherical droplets without considering the gravity of droplet under micro-scale. The partial differential equation can be transformed into ordinary differential matrix equation by algebraic operation when the displacement and stress are expressed by vector function system, and then the inter-layer displacement and stress are given by using transfer relation. Finally, the expressions of displacement and stress at any field point are given by using layered technology. The results have potential application value in designing and fabricating soft devices.

Model establishment and basic equations
The semispherical droplet is placed on a layered elastic soft substrate, and the bottom is fixed with a rigid surface or a semi-infinite space. It is assumed that the layered structure is not affected by body force, and the layers are in complete contact, that is the displacement and stress at the interface are continuous, and the thickness, elastic modulus and Poisson's ratio of each layer are known. The contact model of elastic substrate deformation caused by equilibrium droplet is established, where the deformation of substrate is axisymmetric, as shown in Fig. 1.
In the cylindrical coordinate system (r, θ , z) , the equilibrium equation of layered elastic medium is where σ ij is the stress components. The geometric equation is  www.nature.com/scientificreports/ where ε ij , u i are strain components and displacement components, respectively. The constitutive relation of layered isotropic materials in cylindrical coordinate system is with where E is the elastic modulus, ν is the Poisson's ratio, and the following relation is given by The inter-layer contact conditions are as follows where u − i (r, θ, h k ), u + i (r, θ, h k ) denote the displacements of the upper surface of the k-layer and the lower surface of the (k − 1)-layer, respectively, and the stress has similar expression.
In this paper, the deformation of the substrate caused by only the surface normal load is considered, and the stress boundary conditions can be expressed as where δ(x) and H(x) are Dirac delta function and Heaviside step function, respectively, P is the Laplace pressure inside the droplet. For semispherical droplets,P = 2γ lv R , where R is the droplet radius, and γ lv is the liquid-air surface tension. Since the bottom surface is fixed on a rigid surface or half space, the displacement boundary conditions are

Constitutive equation under vector function system
The vector functions in cylindrical coordinate system can be expressed as 33,34 where e r , e θ , e z are unit vectors in r, θ , z directions, and S satisfies Helmholtz function, i.e. and where J m ( r) is the m-th order Bessel function of the first kind, j = √ −1 , and , m correspond to the transformation variables of horizontal physical variables r, θ . m = 0 corresponds to axisymmetric deformation. For brevity, let (2) ε rr = ∂u r ∂r , ε θθ = ∂u θ r∂θ + u r r , ε zz = ∂u z ∂z , σ rz (r, θ, 0) = σ θ z (r, θ , 0) = 0, σ zz (r, θ , 0) = −γ lv δ(r − R) + PH(R − r), u i (r, θ , z n ) = 0, (i = r, θ , z).  www.nature.com/scientificreports/ Due to the orthogonality and completeness of the vector function system, the displacement and traction vectors of any field point on every layer can be expanded as follows 33 : where U K , T K (K = L, M, N) are the expansion coefficients of displacement and stress vector functions in cylindrical coordinates, respectively. These values can be expressed as 35 where the horizontal line above vectors represents the complex conjugate vector. From Eqs. (7), (12) and (13), the expressions of the transformation domain of displacement and stress are as follows L(r, θ ; , m) = e z S(r, θ ) L(r, θ ), M(r, θ ; , m) = e r ∂ ∂r + e θ ∂ r∂θ S(r, θ ) M(r, θ ),  www.nature.com/scientificreports/  In this paper, the coefficient formula with N subscript will not appear in the following calculation when only the normal load is considered. Rewriting Eqs. (26), (27), (34) and (35) in matrix form as follows www.nature.com/scientificreports/ Rewriting the former formula gives In this case, with the help of the vector function system in cylindrical coordinates, the partial differential equations are transformed into the following ordinary differential equations as The vector function of the stress boundary condition at the structural surface layer z = 0 in the cylindrical coordinate system is expressed as 33 The expansion coefficient is expressed as follows    .  Since [E(z 0 )] has been solved, Eq. (42) becomes a set of quaternion linear equations about the expansion coefficients U L (z k ), U M (z k ), T L (z k ), T M (z k ) of the displacement and stress in the k-layer. Thus, the normal displacement of each inter-layer deformation can be obtained as follows The radial displacement of the inter-layer is According to the expansion formula of stress, the normal stress of the inter-layer is The radial shear stress of the inter-layer is (2) If the bottom boundary is a semi-infinite space, the displacement at infinity is zero. The overall transfer relation is where [a k ] is the transfer matrix corresponding to k-layer.
In a continuous homogeneous semi-infinite space at z = z 0 , the displacement at infinity is zero, so the two coefficients corresponding to the displacement in matrix [K] are zero, and the coefficients corresponding to the stress are K 31 , K 41 . As the traction force exerted by the droplet on the substrate at the surface is solved, the coefficient containing T L , T M in [E(z 0 )] will be determined, then Eq. (45) becomes where Equation (46) is also a system of quaternion linear equations, which can be easily solved as At this point, the coefficients of the displacement in [E(z 0 )] and the stress in matrix [K] are obtained, the normal displacement of substrate surface deformation is C σ z ( ) = −a 44 γ lv a 33 a 44 − a 43 a 34 . 34 .
· · · [a n ][B n ] e ν * n h n = b l,m 4×4 , (l, m = 1, 2, 3, 4). 34 . The normal stress on the bottom of the substrate is The radial shear stress of the bottom surface of the substrate is The solution of displacement and stress between layers is similar to that of rigid layer, and the expressions of the displacement and stress are as follows.
The normal displacement of the inter-layer surface is The radial displacement of the inter-layer is The normal stress of the inter-layer is The radial shear stress of the inter-layer is As a result, we have solved the problems of deformation and stress on the surface, the bottom and the interlayer, respectively. For the displacement and stress of any point in the field, the layered technology is implemented to discuss. Taking the k-layer as an example, the local coordinate system O 1 − r − ξ is established. From the point to be solved, the layer is divided into two layers with thickness h k1 , h k2 , respectively, as shown in Fig. 1. In this case, the boundary conditions of the substrate surface and the bottom are given, so it is easy to establish the transfer relation from the point to the last layer.
(1) The bottom boundary is rigid layer, and the transfer relation is www.nature.com/scientificreports/ Using the bottom boundary conditions, we can easily obtain the expansion coefficients U L (z hk2 ), U M (z hk2 ), T L (z hk2 ), and T M (z hk2 ) of the displacement and stress.
The normal displacement is The radial displacement is expressed as The normal stress has the form as follow The radial shear stress is (2) The bottom boundary is semi-infinite space, and the transfer relation is It is similar to the process of case (1), the stress and displacement can be easily obtained by using the bottom boundary condition, so we will not repeat it here.

Numerical results and discussion
Taking the radius of droplet R = 200 µm , thickness of substrate h = 50 µm , elastic modulus of substrate E = 3 kPa , Poisson's ratio ν = 0.48 and gradient index a = 0 , Fig. 2 shows the comparisons between the normal displacement of our results and those in Ref. 20 . It can be seen that the deformation displacement of a single-layer structure in this work is in good agreement with that in Ref. 20 , and the error range is within 10 −7 , which shows that the model established in this paper is correct.  Figure 3 shows the variations of substrate surface deformation with modulus. It can be seen that the normal and radial deformations of the substrate decrease with the increase of the elastic modulus of the substrate, which is consistent with the experimental results. Figure 4 shows the variations of the normal and radial displacements of the single-layer structure surface deformation with the contact radius and substrate thickness. The results indicate that the normal displacement increases with the increase of substrate thickness, and the change range of wetting ridge also increases. Similarly, with the increase of substrate thickness, the radial deformation is becoming ever more larger. However, compared with the normal deformation, the radial deformation is very weak, so the normal displacement can well describe the surface profile of the substrate. Figure 5 demonstrates the variations of the normal stress and tangential normal stress of the bottom surface with contact radius and substrate thickness under single layer structure. It shows that the substrate thickness has little effect on normal stress, and it can be seen from local enlargement that the whole normal stress increases with the decrease of the substrate thickness. On the contrary, the shear stress decreases with the decrease of substrate thickness. Comparing Fig. 4 with Fig. 5, it is found that the changes of radial displacement and shear stress show opposite trends. www.nature.com/scientificreports/ (II) The droplet wets the bilayer structure, and let R = 200 µm , a = 20 , h 1 = 20 µm , E 1 = 2 kPa , ν 1 = 0.42 , h 2 = 60 µm , E 2 = 4 kPa , ν 2 = 0.44 , and γ lv = 0.05 Nm −1 . The following numerical results are obtained. Figure 6 shows the variations of the normal and radial displacement of the surface and inter-layer surface of the double-layer structure. It can be seen that the displacement of inter-layer surface is larger than that of surface. In addition, it also indicated that the closer the contact radius r is droplet radius R, the greater the difference in the displacement of the two layers will be. This is because the surface tension near the three contact lines has a greater influence on the deformation. As r is far away from R, the deformation of the substrate becomes weaker, finally, the displacements of the two layers almost coincide. Figure 7 illustrates the stress changes on the surface and inter-layer surfaces of the bilayer structure. It can be seen that the closer r is to R, the greater the difference in the stress changes between the two layers will be. On the contrary, the stress changes of the two layers almost coincide. As described in Fig. 6, the reason is that the surface tension near the three contact lines leads to the phenomenon. Figures 6 and 7 also show that the change trends of tangential displacement and tangential stress are opposite. www.nature.com/scientificreports/

Conclusion
In this paper, the wettability of a semispherical droplet on gradient layered structure was researched. Both theoretical formulas and numerical results showed that surface was producing weak radial shear stress when the droplet was placed on a thin substrate, which resulted in the very small radial deformation of the substrate. Due to unequal contact angles during the deformation process, the substrate led to tangential deformation. The numerical results also illustrated that the deformation of the substrate decreased with the increase of elastic modulus, also the variation amplitude of displacement and shear stress increased with the increase of substrate thickness, while the change of normal stress was opposite. In addition, the changes of displacement and stress between inter-layer and surface were different with the relationship of the size between contact radius and droplet radius. The contact radius was closer to droplet radius, producing the larger difference of displacement and stress between the surface and inter-layer. On the other hand, because of surface tension near the three contact lines, the changes of the two layers almost overlapped. The results of this paper have potential application value for the wetting of thin films and soft materials.