Mechanical Deconvolution of Elastic Moduli by Indentation of Mechanically Heterogeneous Materials

Most materials are mechanically heterogeneous on a certain length scale. In many applications, this heterogeneity is crucial for the material’s function, and exploiting mechanical heterogeneity could lead to new materials with interesting features, which require accurate understanding of the local mechanical properties. Generally used techniques to probe local mechanics in mechanically heterogeneous materials include indentation and atomic force microscopy. However, these techniques probe stresses at a region of finite size, so that experiments on a mechanically heterogeneous material lead to blurring or convolution of the measured stress signal. In this study, finite element method simulations are performed to find the length scale over which this mechanical blurring occurs. This length is shown to be a function of the probe size and indentation depth, and independent of the elastic modulus variations in the heterogeneous material, for both 1D and 2D modulus profiles. Making use of these findings, we then propose two deconvolution methods to approximate the actual modulus profile from the apparent, blurred measurements, paving the way for an accurate determination of the local mechanical properties of heterogeneous materials.


INTRODUCTION
Controlling the mechanical heterogeneity of materials could lead to new innovative materials with interesting properties [1], such as scaffolds for tissue engineering [2], stimuli-responsive materials [3], materials with soft-lithography applications [4] and materials with enhanced lifetime [5]. To design such materials, and to understand their mechanical performance, methods are needed that can accurately measure the spatial variation of mechanical properties. This is challenging, since the mechanical response of a material depends strongly on the length scale at which it is probed [6,7]. The most common method to measure local mechanical properties in heterogeneous materials is indentation, in which an indenter is pushed on a material and the local resistance force to deformation is acquired with a force sensor. Using an appropriate contact model, such as Hertzian contact theory [8], the local elastic modulus E can be obtained from the measured force-displacement relation. However, the indenting probe deforms the material in a region of finite size, so that variations in mechanical properties that occur on length scales that are similar to the size of the deformed region or smaller necessarily appear smoothed out. This is illustrated in Figure 1, which shows a material with a sharp gradient in the modulus. Probing the material using spatiallydependent indentation leads to a blurred profile, which is a convolution of the actual modulus profile and the volume of material probed by the indenter. Depending on the size of the probe, this blurring may occur on the nm scale for nano-indentation using atomic force microscopy [9][10][11][12][13], on the µm-scale for microindentation [14].
Reference [15] characterization, or on the mm-scale or larger for macroscopic probes [16][17][18]. Previous studies have shown that the extent of mechanical blurring depends on the indentation depth [6] and occurs over an area that is at least three times the contact area of the probe with the substrate [7]. However, the exact dependence of the mechanical response on the probe size and depth remain unknown. Clearly, the characterization of mechanically heterogeneous materials would benefit from a better understanding of the effects of mechanical blurring, and from ways of improving the spatial resolution of mechanical measurements.
A similar problem occurs in optical imaging, where the captured image is a convolution of the real image and the so-called point spread function, which describes the response of the imaging system to a point force and thereby characterizes the degree of blurring. If the point spread function of an optical device is known, deconvolution methods can be used to deblur the signal and increase the image quality [19,20]. In this work, we explore a similar deconvolution strategy for indentation measurements on mechanically heterogeneous materials. We first identify the mechanical analogue of the point spread function. We do this by performing finite element calculations on materials with a known modulus profile E(r), indented by a spherical probe. From the calculated normal force F(r), we obtain an apparent local modulus E a (r), which we compare to the real modulus E(r) to assess the degree of blurring and to estimate the mechanical point spread function. We find that the blurring is independent of the magnitude of E, which is a prerequisite for deconvolution to be feasible, since the real modulus is unknown. In addition, we find that the deblurring is set by a combination of the probe radius R and the indentation depth δ, which together determine the length scale over which the probe deforms the material. Based on these findings, we then describe two methods to approximate the real modulus profile E from the apparent moduli E a measured for different values of R and δ.

Finite Element Method Simulations
Finite Element simulations were performed using COMSOL Multiphysics 5.4. The 3D stationary structural mechanics module was used to model the indentation of a rigid sphere on a linear elastic substrate. The linear elastic substrate was meshed with small hexahedra swept underneath the spherical probe and coarsely meshed far away from the probe, as shown in Figure 2A, and was solved for 43 k degrees of freedom. Three probe radii were simulated; 400, 500 and 750 µm. The maximum indentation depth was set at δ max 0.4·R and the step size was 4 µm. The linear elastic substrate was modeled with a fixed constraint at the bottom. Contact between the rigid sphere and the elastic substrate was implemented by the penalty method. Further model details are provided in the Supplementary Material (SM). With F obtained from the bottom of the PDMS substrate in the simulations, the apparent modulus E a was calculated using Hertzian contact theory [8], where Y is the apparent Young's modulus, and ] is Poisson's ratio, which was here taken as 0.45 [21]. FEM simulations were performed on a total of four profiles: i) a 1-dimensional (1D) stepwise modulus profile ii) 1D sigmoidal profiles with varying widths w iii) 1D wells with various widths w, and iv) 2-dimensional (2D) wells with various widths w 2 . The modulus profiles were incorporated with a piecewise function in COMSOL's material section using the following relations: x the position, x 0 the position of the interface and H the Heaviside function defined as: | Schematic representation of a material that has a step-wise modulus profile and the resulting modulus profile measured by an indentation experiment at several positions across the modulus interface. The four insets a-d display the heterogeneous stress profiles in the material, color coded for the stress, as the material is indented with, in this example, a spherical probe.

Heterogeneity Leads to Asymmetric Stress Profile
FEM simulations are performed on materials with various elastic modulus profiles. Figure 2B shows, for the 1D sigmoidal profile with width w 10.0 mm, the elastic strain energy density profile u U V ϵ f 0 σzϵ. For mechanically heterogeneous materials, the strain energy profile is asymmetric, leading to a convolution in the locally recorded force, which becomes a weighted average of the stress profile. Hence, this calls for a method that converts the convolved profile in a deconvolved profile that approximates the actual material profile, such that mechanically heterogeneous materials can be mapped accurately.

The Mechanical Point Spread Function
To find the mechanical point spread function that sets the degree of mechanical blurring, a stepwise profile (Eq. 2a) is simulated with FEM. Results in Figures 3A,B show that the mechanical blurring is dependent on R and δ. Likely, this is because R and δ set the contact radius a Rδ √ between the probe and the substrate. The larger a, the larger the area probed by the sphere and the larger the mechanical blurring. At this point, we assume that E a is a convolution of E and the mechanical blurring function g (r −r′), with the real modulus at position r′ contributing to the modulus probed at position r, analogous to the point spread function: Next, we assume that the mechanical blurring function g (r −r′) is a Gaussian function, independent of E: where L is a characteristic length that determines the degree of blurring. Eq. 3 can be solved analytically for the stepwise profile, by inserting Eq. 2a and Eq. 4 in Eq. 3, resulting in: with E c the convolved elastic modulus, which should approximate E a . By fitting Eq. 5 to E a from the FEM simulations, the characteristic blurring length L can be determined as a function of R and δ. Indeed, the results presented in Figures 3C,D show a collapse of E a when x is rescaled with L obtained from these fits, as expected based on Eq. 5. Furthermore, to verify whether L ∼ (R, δ), Figure 4 shows that L can be fitted with L b · (Rδ) p . Collapse of L for multiple E 1 E2 ratio's shows that L is independent of E. These findings show that our assumption of a Gaussian point spread function works well for deconvolving the mechanical blurring. Now that we have obtained the mechanical point spread function and its characteristic blurring length L from the analysis of the blurred step profiles, we investigate whether this can be used to deconvolve apparent moduli measured with indentation to obtain a good approximation of the actual modulus profile of heterogeneous materials. We propose two methods to do this.

Deconvolution of 1-Dimensional Profiles
The first method to deconvolve a blurred profile is through extrapolation of E a to L 0, as at zero indentation depth the material is not deformed and there is no blurring of the mechanical response and, hence, E a E in this case. This extrapolation to obtain E e is shown in Figure 5A for a sigmoidal profile with width w 10.0 mm, for several probe positions x − x 0 . An advantage of this extrapolation method to deconvolve E a into E e that approximates E is that the precise shape of g(r − r ′ ) does not have to be known.
The second method to deconvolve a blurred profile is accomplished with an algebraic deconvolution approach. First, we write the convolution as a matrix operation: with G a matrix describing the mechanical blurring. Every row of G contains coefficients that indicate how much the modulus at position r′ contributes to the modulus probed at position r. These coefficients are determined by g (r −r′), of which L is known from the simulations on the 1D steplike profile. To obtain the unknown E from the measured E a we find a least squares solution to Eq. 6 using an iterative method, Simultaneous Iterative Reconstruction Technique (SIRT) algorithm [22,23]. Further details of this method are provided in the SM. With the two methods presented above, E a is extrapolated and deconvolved leading to E e and E d , respectively. E e and E d convert E a into an approximation of E, as shown in Figures 5B,C for a sigmoidal profile with width w 10.0 mm and a well profile with width w 2.0 mm, respectively. This shows that both presented

Deconvolution of 2-Dimensional Profiles
After demonstrating mechanical deconvolution for 1dimensional (1D) modulus profiles, simulations were performed to showcase the potential of deconvolving a 2dimensional (2D) modulus profile. In these simulations, the same meshing was used as for the 1D profiles. Symmetry was assumed across x x 0 and y y 0 and, therefore, a quadrant was modeled to reduce computation time.
The discrepancy between the apparent modulus E a and the actual modulus profile E, shown in Figures 6A,B, once again demonstrates the necessity for the mechanical deconvolution methodology demonstrated for the 1D profiles. The convolution kernel g and L obtained from the 1D step profile are used to extrapolate and deconvolve E a .
The results from these simulations show that both the extrapolation and deconvolution methods presented in this work are able to deconvolve both 1D and 2D mechanical profiles, approximating a blurred E a into profiles close to the actual profile E. Results in Figures 3, 4 confirm that mechanical blurring is at least three times L, as observed by [7]. In addition, simulations on two E1 E2 ratio's for 2D wells show that the differences between E a , E e and E d scale approximately with the E1 E2 ratio, albeit that larger deviations are found for smaller w, see Supplementary Figure S5. Furthermore, like the 1D profile, a sharper modulus profile leads to blurring of E a , making extrapolation and deconvolution into an approximation of E more challenging. Examples are Supplementary Figure S3A  and Supplementary Figure S4A, where a sharp, (almost) discontinuous transition in E leads to the Gibbs phenomenon [24] in the deconvolved E e and E d ; oscillations of the Fourier function at sharp transitions which in imaging is known as the ringing artefact [25]. Nevertheless, the two methods presented allow deconvolution of E a into an approximate of E.

Towards Experimental Validation
Our finite element simulations show how indentation measurements on heterogeneous samples lead to mechanical blurring, and how knowledge of the mechanical point spread function can be used to deconvolve the blurred signal to obtain a more accurate modulus profile, experimental validation is required to transform this proposed method into a reliable practical procedure. To experimentally determine the mechanical point spread function and the associated characteristic length, and to compare this to simulation results shown in Figure 4, a material with a precisely known gradient in modulus is needed, ideally a x-dependent step function in modulus. Indentation measurements on this reference sample then allows the determination of the characteristic length in exactly the same way as described in Section 3.2. Unfortunately, we have not been able to produce such a material of known modulus gradient with the capabilities in our lab, without also creating gradients in the depth-direction or differences in sample height across the modulus step, which made it impossible for us to provide such validation. Our attempts are described in detail in the SM. We hope that this work will inspire others to provide the experimental validation needed to turn our approach into a robust experimental method. We further note that in our approach we have assumed the substrate to be a purely elastic material. While this is justified for elastic solids, such as crosslinked rubbers, care should be taken for viscoelastic soft materials that also have viscous (and therefore rate-dependent) contributions to the mechanical response. Furthermore, since our method is based on linear elasticity, deviations may occur for large indentation depths, where strain hardening may become relevant.

CONCLUSION AND FUTURE OUTLOOK
The work presented in this paper shows two methods to transform a blurred and measured E a towards an expected modulus E in FEM simulations in 1D and 2D. This transformation is achieved by using a characteristic length L, which is found to be a function of δ and R. The methods presented to deconvolve a blurred mechanical profile should work on smaller length as well, for example for mechanical maps obtained using spatially resolved AFM-based force Frontiers in Physics | www.frontiersin.org September 2021 | Volume 9 | Article 723768 6 spectroscopy [10,11], provided that the material can still be considered as a flat elastic continuum. For this it is necessary to obtain force-distance curves for each position, so that effective modulus data at the same penetration depth can be compared. To validate the extrapolation and deconvolution methods proposed here experimentally, indentation measurements on a material with a precisely known modulus profile are performed. Unfortunately, we did not succeed in obtaining these, but we hope that future work will demonstrate the potential of mechanical deconvolution. We anticipate that the mechanical deconvolution technique presented here can be used to approximate the actual mechanical heterogeneous material from a blurred stress response, leading to more accurate knowledge of local mechanical properties in for example biological systems [12,14], polymeric materials [6,7,15,26] and meat analogues [27].

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
JB, Conceptualization, Simulations, Data processing and analysis, Experiments, Writing-original draft, Writing-review and editing. RK, Conceptualization, Simulations, Data processing, Writing-review and editing. TK, Conceptualization, Experiments, Writing-original draft, Writing-review and editing.