Issue 
4open
Volume 6, 2023
Special issue on GDR Cosm’actifs



Article Number  7  
Number of page(s)  12  
Section  Physics  Applied Physics  
DOI  https://doi.org/10.1051/fopen/2023006  
Published online  13 June 2023 
Research Article
Dwell time in contactfree creep tests plays an agedependent role in the viscoelastic behavior of in vivo human skin
Université de Lyon, ECLENISE – Laboratoire de Tribologie et Dynamique des Systèmes (LTDS) – UMR 5513 CNRS, 58 rue Jean Parot, 42023 SaintEtienne, France
^{*} Corresponding author: marieangele.abellan@enise.fr
Received:
20
June
2022
Accepted:
10
May
2023
The experimental characterization of human skin in vivo may be highly dependent on the experimental protocol and on the procedure employed for analyzing the experimental data. One way of overcoming these drawbacks is the contactfree indenter at the Laboratoire de Tribologie et Dynamique des Systèmes, Lyon, France. However, it is still necessary to determine the possible influence of the experimental protocol in terms of the dwell time chosen. This paper describes a coupled theoretical/experimental/numerical study on the viscoelastic responses of the human skin at two different ages during contactfree creep tests with dwell times of 0.1 s, 0.5 s, 1 s, 2 s and 100 s, respectively. The 3D finite element simulations are conducted with the SYSTUS^{®} software. The numerical results point to different behaviors for the young and elderly subjects. The young subject does not appear influenced by dwell times during the numerical creep tests. On the contrary, the response of the elderly subject is highly dependent on dwell time.
Key words: Human skin in vivo / Contactfree indenter / Numerical simulation / Viscoelastic parameters / Large displacements / Nonlinear analysis
© M. Abellan et al., Published by EDP Sciences, 2023
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Introduction
The human skin is a living organ constantly renewed by new cells. Its physiology and its evolutions during aging influence its responses to coupled chemicobiologicalthermohydromechanical influences. Moreover the response can be different depending whether the load is applied on the outer surface, in the volume through the surface or directly in the tissues. Understanding these mechanisms is of prime importance for surgical and medical purposes and in dermatology for both curative and care needs. Nowadays, biorobots and devices are being developed and used to improve patientdesigned therapies. They must rely on robust, objective, qualitative and quantitative parameters.
In vivo human skin has a multilayered structure well adapted to its barrier functions. It is a stratified medium with four main layers: the stratum corneum (SC), the viable epidermis (VE), the dermis (DE) and the hypodermis (HY) [1–3]. The stratum corneum is made of hard dead corneocytes embedded in a lipid mortar [4, 5]. On the contrary, the viable epidermis contains evolving soft cells bathed in a physiological fluid. Its cohesion stems from the large amount of desmosomes that bridge the cells together. The dermis has dense fiber networks of collagen and elastin surrounded by an amorphous ground substance [6]. The adsorbed part of this ground substance provides cohesion between the fibers [7]. The hypodermis is composed of lipid cells separated by fibrous walls [8]. Together the fibers and the adsorbed ground substance initiate a structural process which gives skin its mechanical resistance and structural integrity. The strong desmosome bridges, the adsorbed part of the lipid mortar, and the fibrous walls also contribute to this structural process. These particularities make human skin a highly complex material (nonhomogeneous, anisotropic, nonlinear viscoelastic material subjected to prestress in vivo [9–12]) difficult to model.
In this context, attempts have been made in the past to simulate the mechanical properties of skin. Models, numerical simulations and/or innovative experiments have been developed. In vivo experimental tests include tensile testing [13], suction methods [14–16], torsion tests [7, 17], and indentation experiments [18, 19]. These works were confronted by the fact that the device had to be attached to the skin throughout the tests and thus disturbed its natural state. To overcome this drawback, the present study is based on contactfree creep tests performed with a contactfree indenter: the Tonoderm^{®} device [20, 21].
Constitutive models are based on simplifying assumptions that consider skin as either homogeneous isotropic linear elastic [15, 22] or nonlinear isotropic hyperelastic [14, 19]. Bischoff et al. [23] coupled a hyperelastic model for fibrous tissue with a viscoelastic reptation model for soft materials. Azzez et al. [24] proposed using a homogeneous linear viscoelastic model to simulate indentation tests at two ages. To capture the in vivo behavior of the skin of the arm, Flynn et al. [25] used an Ogden strain energy function coupled with viscoelasticity and anisotropic prestress. Abellan et al. [26] proposed comparing four viscoelastic models to simulate indentation tests and to describe the coupled transient water flow and ion transport through skin. But these studies were still based on experimental contact tests.
This paper goes further by proposing a coupled theoreticalexperimentalnumerical simulation approach to study the influence of the dwell time of numerical contactfree creep tests on the viscoelastic responses of human skin in vivo at two ages. This study is intended to improve understanding of the state (healthy or diseased) of the skin before and after treatment in order to provide practitioners with clues for adapting the protocol of contactfree creep tests to the special needs of their patients, based on the knowledge of their state to design the appropriate curative/care therapy.
The experimental contactfree creep tests were performed on the volar forearm of two volunteers: a young adult and an elderly subject. An inverse approach was taken to determine the viscoelastic characteristics of the skin of both subjects. To ensure admissible physical consistency between the experimental data and the numerical simulations, the numerical tests were based on these experimental data. Five numerical contactfree creep tests were formulated with dwell times of 0.1 s, 0.5 s, 1 s, 2 s and 100 s, respectively.
The 3D numerical simulations of the five numerical contactfree creep tests were conducted using the SYSTUS^{®} software. The skin is seen as a stratified medium described by a viscoelastic law with a large displacement formulation. This model allows a highly nonlinear description of the skin’s behavior. The main advantage is to gain insights into the mechanical properties of each layer. However, this is not the sole goal. It appears that a simple observation of the response of the skin surface to an air flow provides knowledge of what occurs in the volume. The mechanical theoretical model considered in conjunction with the constitutive viscoelastic law presented in this paper makes it possible to analyze the evolutions of each layer on its own and the relative influences of each layer on the others. The aim is to deduce the contribution of each layer, or more generally of each material point of the human skin, to the overall response of the skin to the different numerical contactfree creep tests, and to find out if certain people respond differently as a function of dwell time and/or age.
In the following, the materials and methods are presented. In particular the contactfree indenter is described and the numerical contactfree creep tests used to study the possible influence of the dwell time chosen are defined. The theoretical model, the numerical simulations and the viscoelastic parameters estimation procedure used in this work are presented. Then the results of the numerical simulations are given. A discussion analyses the main features of the results. The conclusion will end this paper.
Materials and methods
Contactfree creep tests in vivo
Experimental device
Air flow as an external load applied has been used for a long time by ophtalmologists in order to measure eye pressure. Based on this principle, the contactfree indenters (Tonoderm^{®} and Waveskin^{®}) at the Laboratoire de Tribologie et Dynamique des Systèmes (LTDS, Lyon, France), has been developed by Prof. H. Zahouani’s team [12, 20]. This device allows performing very difficult tests even on damaged or painful skins in vivo. Indeed, in such cases, it is compulsory to avoid any contact with the skin tissue in order to prevent contamination or destruction of the default studied. Moreover, in the case of contactfree indentations, the load applied can be removed instantaneously. It is then possible to observe and measure the undisturbed and free evolutions of the skin in vivo coming back to its natural stress state.
The experimental setup of the Tonoderm^{®} and the sketch of its network are given in Figures 1a and 1b, respectively. A compressor (A) provides the source of compressed air necessary to the system. The regulation of this compressed air is firstly done by two expansion valves (not shown on Fig. 1) connected in series. They convert the pressure of 6 bars provided by the compressor at the working pressure of the regulator which is 4.5 bars. The choice of two instead of one regulator is based on the fact that such assembly can effectively minimize the change in the output pressure due to a variation of the inlet pressure. The air flow is then controlled by a mass flow controller (B) with a flow range from 2 to 100 NL/min (Normalized litters per minute) expressed in normal conditions of pressure P_{N} (Pa) and temperature T_{N} (K) for a standard density ρ_{N} (kg/Nm^{3}) which are for the air: T_{N} = 273.15 °K, P_{N} = 101,325 Pa and ρ_{N} = 1.293 kg/Nm^{3}. At this regulator’s exit are connected two solenoid valves, the first one (C) aims to empty the air network before a test. The second one (D) is placed right before the exit of the air flow. Therefore this air flow can be regulated (diameter 2 mm, distance exitskin 10 mm). The air flow reaches the skin surface with an invertedbell shape (Fig. 1c). It is equivalent to a spherical indenter with a tip surface of 6.28E6 m^{2} (Fig. 1d). The deflection is measured by a laser beam (E) (LKG82, Keyence). The experimental setup is controlled with a computer (F) and the experimental data are recorded with an analogical/numerical card (16 bits, 6221 M, National Instruments).
Figure 1 (a) Principle of the Tonoderm^{®} [20]. (b) Forearm setup [31]. (c) Photo of the pressured area of the skin surface under the air flow [20]. (d) Invertedbell shape of the air flow and typical contact profile which may develop. 
Protocol and experimental data
The experimental tests of this study are based on the following protocol:
Before the tests, the subjects were resting seated in a medical chair during 10–15 min in a room with controlled temperature (22 °C) and humidity (RH 50%).
The imposed force tests were performed on the inner forearm of the volunteers.
The tested area was 6 cm from the elbow.
For the tests, the system applies an air flow onto the outer surface of the forearm until the pressure equivalent to a flux of 20 NL/min is obtained (loading branch of 9.1 s). Then the flux is kept constant during 1 s. Finally the load is removed (the flux is stopped) and the system discharges. The freereturn is monitored for 14 s. The mean of 3 measurements has been done.
The experimental data recorded during the indentation tests give profiles of the load (Fig. 2) and of the deflection with respect to time.
Figure 2 Experimental loading curve as a function of time (Exp). 
Numerical contactfree creep tests
The possible influence of the dwell time chosen when performing the in vivo contactfree creep tests described before is analyzed based on five numerical contactfree creep tests numerically designed and presented here. This is done in order to ensure admissible physical consistency between the experimental records and the numerical findings.
The protocols of the five numerical contactfree creep tests are: a loading branch of 9.1 s; dwell times of respectively 0.1 s, 0.5 s, 1 s, 2 s and 100 s; then the applied pressure simulating the air flow at the upper surface of the REV is stopped and the freereturn is monitored for 3600 s for each test.
The choice of monitoring the freereturn for 3600 s lies on the observation of the numerical results obtained. After 3600 s, the human skin has completed its freereturn. There is no more evolution in the responses of the skin. It will be shown, discussed later in the paper, that the numerical results exhibit close hysteresis cycles after 3600 s both for the young subject and for the elderly subject.
The loading curves of the numerical tests are given in Figure 3.
Figure 3 Loading curves of the five numerical contactfree creep tests with dwell times of 0.1 s, 0.5 s, 1 s, 2 s and 100 s. 
Theoretical model
Within the framework of the thermodynamics of irreversible processes [27] and following Bathe [28], in an updatedLagrangian formulation, when neglecting inertia forces, the total balance of momentum reads:
(1)where S(x, t) is the second stress tensor of Piola–Kirchhoff of the skin in vivo in Pa, Ω is the configuration of the REV gathering all the material points M(t) defining the REV for the updatedLagrangian description and x are the spatial coordinates of the material point M(t). Following Latorre and Montans [29], the constitutive equations are derived here in terms of the strain measure: the Eulerian rate of deformation tensor D and the corresponding constitutive stress measure: the (objective) Green–Naghdi stress rate σ^{△} defined by:
(2)R being the rotation tensor coming from the polar decomposition of the deformation gradient.
Skin is supposed to follow an isotropic behavior and a Zener rheological model is chosen for modelling the skin behavior (Fig. 4). The Zener model is constituted by an elastic branch in parallel with a viscoelastic Maxwell’s branch. Therefore, the Cauchy stress is defined as the sum of an elastic part σ^{E} and a viscoelastic Maxwell’s part σ^{M}. It reads: σ = σ^{E} + σ^{M}.
Figure 4 Rheological Zener model. 
In the elastic branch, the generalized Hooke’s law (Eq. (2)) is written in a rate form:
(3)where ν is the Poisson’s ratio, E_{∞} is the long term modulus in Pa, δ is the Kronecker delta.
The long term modulus E_{∞} is defined by:
(4)where E_{i} is the Young modulus of layer i, i standing for layers SC, VE, DE and HY, respectively.
In the viscoelastic Maxwell’s branch, for an isotropic material, we suppose an additive decomposition of the Eulerian rate of deformation tensor:
The stressstrain rate relations then separate volumetric and deviatoric responses of the skin in vivo. Toward this end, the Maxwell part σ^{M} of the Cauchy stress tensor and the Eulerian rate of deformation tensor D are decomposed into a spheric and a deviatoric part:
(7)where σ′^{M} and are respectively the deviatoric part of the stress tensor and the mean stress in the Maxwell’s branch, D′ is the deviatoric part of the eulerian rate of deformation tensor and D_{m}, the volumetric part.
The volumetric and deviatoric responses of the material are then expressed by:
(9)with η_{v} the volumetric viscosity in Pa s, , the compressibility modulus in Pa, η_{d} the deviatoric viscosity in Pa s, G^{M} = E_{1}/2(1 + ν) the shear modulus in Pa and E_{1} the modulus in the Maxwell’s branch in Pa.
The balance equation (1) and the behavior law (2), (8) and (9) are complemented by physically admissible boundary conditions set on additional parts of the boundary ∂Ω of the REV:
(10)with ∂Ω = ∂Ω_{t} ∪ ∂Ω_{u} and ∂Ω_{t} ∩ ∂Ω_{u} = ∅, σ(x, t) the Cauchy stress tensor, n the outward unit normal of the boundary ∂Ω_{t}, t_{p} the prescribed effort on the boundary ∂Ω_{t} in Pa, u(x, t) the displacement vector of the material point M(t) in meter and u_{p} the prescribed displacement on the boundary ∂Ω_{u} in meter.
Numerical simulations
The theoretical model is applied to the numerical simulations of the experimental contactfree creep tests in vivo and also to the numerical simulations of the numerical contactfree creep tests. The boundary value problem is solved, 6 times for the young subject and 6 times for the elderly subject (one time for each protocol for each subject), using the finite element method with the SYSTUS^{®} software.
When defining the representative elementary volume (REV), the structural architecture spread in all layers is of importance for defining a userdefined REV and work at the soft cutaneous tissue scale. The thicknesses of the forearm were extracted from echograph images available at the Laboratory for a panel of 28 volunteers (14 women 23.2 ± 1.6 yearsold and 14 women 60.4 ± 2.4 yearsold) [20] from which were extracted the records of the young and the elderly volunteers of this study. It has been found for the thicknesses of the three layers SC + VE + DE, a mean of 1.094 mm ± 0.128 mm for the thickness of the young volunteers and 1.033 mm ± 0.108 mm for the thickness of the elderly volunteers. An analysis of variance conducted on these results has led to a statistically nonsignificant difference between the two groups. These results are similar to the measurements reported by Agache [7]. Hence, the userdesigned REV considered here (Fig. 5) is a stratified volume composed of four layers: SC, VE, DE and HY with dimensions 30 mm × 30 mm × 11.11 mm where the thicknesses of each layer is taken equal to: 0.01 mm for SC, 0.1 mm for VE, 1 mm for DE and 10 mm for HY. The same REV appears to be admissible for both subjects.
Figure 5 Userdesigned REV deduced from the echograph images and 3D mesh. 
A large displacement analysis with an updated Lagrangian formulation allows handling the geometrical nonlinearities as well as the nonlinearities resulting from the behavior law. For the time increment, the choice was made of 0.1 s for the loading branch, 0.001 s during the dwell time and 1 s for monitoring the freereturn. These choices for the formulations and for the time increments ensured the convergence of the computations.
For symmetry reasons, all the simulations were carried out on ¼ of the REV (Fig. 5) with dimensions 15 mm × 15 mm × 11.11 mm (Fig. 6). Furthermore, it is important to work with a mesh adapted to the mechanical settings of our model. In order to do so, SC, VE and DE are meshed with prisms while HY is meshed with tetrahedrons. The numerical results have not shown any mesh dependence on the elements chosen. The mesh is refined in and close to the loading zone (blue front corner in Fig. 6) for reasons of convergence and precision. The mesh consists of 9689 nodes spread into 25,181 elements.
Figure 6 3D mesh with the boundary conditions. 
The boundary conditions prescribed are: a Gaussian function P_{0} is used to model the pressure distribution given in Figures 2 and 3. It is applied on the top corner (blue zone) of the mesh (Fig. 6); a zero displacement vector for all the nodes of the bottom surface; zero horizontal components in both directions for the displacement vector for all the nodes of the vertical surfaces for reasons of symmetry (Fig. 6).
The mechanical characteristics should be experimentally determined. For the Poisson ratio, Zahouani et al. [30] have reported of contact indentation and static friction tests performed on volunteers aged 55–70 years. They found a Poisson ratio between 0.44 and 0.48. Delalleau et al. [18] performed numerical simulations of contact indentation tests which lead to a Poisson ratio of 0.48. In a first theoretical model [21], we have worked with a Poisson ratio of 0.3 for obtaining convergence of the numerical process. In 2017, Ayadh et al. [31] proposed an extension of the first model with a Poisson ratio of 0.45 which presented no problem of convergence. Based on that information, the choice is done here to take a Poisson ratio of 0.45 for all the layers both for the young subject and for the elderly subject. The Poisson ratio should be complemented with viscoelastic parameters experimentally determined using the inverse procedure presented here after.
Inverse procedure
The input file supplies the finite element program with initial values for all the parameters. The relative error e_{rel} between the experimental and the numerical deflection is calculated using Equation (11):
(11)where u_{exp} and u_{num} are respectively the experimental and the numerical displacements in meter. If this error e_{rel} is less than a given tolerance, the procedure stops: the parameters used in the calculation are the parameters sought. If the error e_{rel} is greater than the tolerance, the procedure continues with new parameters. The procedure stops when the calculations converge at an estimate of the parameters that minimize the error e_{rel}. Validation of this approach with a sensitivity analysis has already been documented [32].
Results
Experimental contactfree creep records in vivo
One volunteer is chosen in each group, i.e. a young healthy adult 23 yearsold and an elderly subject 62 yearsold. Figure 7 presents the deflection of the outer skin surface of the forearms in the test zone recorded for the young and the elderly subjects.
Figure 7 Profiles of the experimental deflection as a function of time for the young subject (blue line) and for the elderly subject (red line) for a flux of 20 NL/min. 
Determination of the mechanical characteristics for the young and the elderly subject: inverse procedure
The estimation of the mechanical parameters for the young and the elderly subjects was carried out using the inverse procedure presented in section “Materials and methods”. Figure 8 shows an example of comparisons between the experimental and the numerical displacements of the central point under the air flow for three different sets of parameters. The relative errors are also calculated and drawn in Figure 8 (Bottom) for the different sets of Figure 8 (Top). It can be noticed the quite good fit between the experimental and the numerical curve 5 for the young subject and between the experimental and the numerical curve 111 for the elderly subject. Table 1 gives the sets of mechanical parameters associated with the numerical results of Figure 8. The set 5 for the young subject and the set 111 for the elderly subject are showing the best fit with the minimum error. They will provide the viscoelastic parameters used in the numerical simulations of the numerical contactfree creep tests for the young and the elderly subjects presented here after.
Figure 8 (Top) Comparison of the experimental and numerical curves of the displacement for three sets of parameters (Bottom) Curves of the relative error for the same sets of parameters. 
Examples of sets of parameters checked during the inverse approach showing the best fit for the young subject (Set 5) and for the elderly subject (Set 111). The instantaneous Young modulus is defined by E_{0} = E_{1} + E_{∞} in MPa.
Numerical results of the numerical simulations of the numerical contactfree creep tests
Taking into account the mechanical parameters estimated in Table 1 for each volunteer, ten boundary value problems (five for the young subject and five for the elderly subject) are now solved for the boundary conditions given in Figures 2 and 3. A key point appeared when studying the numerical results. For a given layer, each material point of the layer exhibits the same kind of stress state as the point of the layer belonging to the central vertical under the air flux but with lower amplitudes. For this reason, the obtained numerical results are presented and analyzed in terms of evolutions of layers exemplified with the numerical results of seven nodes of the central vertical under the air flux: top node of SC (6966), middle node of SC (5832), boundary between SC and VE (4688), boundary between VE and DE (4116), middle node of DE (3515), boundary between DE and HY (2920) and middle node of HY (4251) (Fig. 9).
Figure 9 Sketch giving the position of the seven nodes of the central vertical under the air flux used in the analysis. 
Figure 10 presents the numerical isovalues of the vertical component of the displacement of the central point under the air flux at the end of the dwell times 0.1 s, 0.5 s and 1 s for the young subject. It gives an example of what happens in the volume and more specifically inside the layers for the young subject. The top surface sinks under the applied air flow. The large blue zone faces the maximum applied load and its response gives the highest values for the displacements. Then for a given layer, the amplitudes of the displacements decrease slowly for material points becoming increasingly distant from the central vertical. At the end of the different dwell times 0.1 s, 0.5 s, 1 s, 2 s and 100 s, the isovalues have the same shape as the ones proposed in Figures 10a–10c for the three dwell times 0.1 s, 0.5 s and 1 s. Moreover the isovalues exhibit equal amplitudes with a maximum deflection of the central point under the air flow of 1.038 mm.
Figure 10 Isovalues of the vertical displacements at the end of the dwell times: (a) 0.1 s, (b) 0.5 s and (c) 1 s for the young subject. (d) Contour lines of the mean stress for the dwell time of 0.5 s for the elderly subject. 
Equivalently for the elderly subject, the most highly stressed area is under the air flow. The influence of the external mechanical applied load decreases with increasing distance (Fig. 10d). However depending on the dwell time the maximum deflection of the central point under the air flow is different. It equals 1.246 mm for 0.1 s, 1.251 mm for 0.5 s, 1.255 mm for 1 s, 1.264 mm for 2 s and 1.470 mm for 100 s.
Figure 11 presents the numerical normal component of the stress tensor σ_{zz} drawn as a function of the numerical deflection for the central point of the top surface of SC under the air flux. It should be noticed the peculiar shape of the graphs of these numerical results available, in Civil Engineering, they are called “hysteresis cycle”. These hysteresis cycles are given for the young subject and for the elderly subject for the dwelltimes 0.1 s, 0.5 s, 1 s, 2 s and 100 s These hysteresis cycles give information on the exchanges of energy taking place between the REV and the outside when the studied medium is subjected to an external mechanical load, i.e. when the material points of the REV are provided with external mechanical energy. They are available here because the theoretical model considers that the responses of the cutaneous soft tissues of the different layers of the REV can be modelled reasonably well with the rheological Zener model. As said before, this viscoelastic model treated with a large displacement analysis with an updated Lagrangian formulation allows to follow the nonlinear evolutions of the responses and in particular the hysteresis cycles.
Figure 11 Numerical normal stress σ_{zz} as a function of the numerical deflection of the central point under the air flux for the young subject and for the elderly subject for the dwell times of 0.1 s, 0.5 s, 1 s, 2 s and 100 s. 
Figure 12 gives the comparison between the numerical normal component of the stress tensor σ_{zz} as a function of the numerical deflection for dwell times 1 s and 100 s for the seven nodes of the central vertical under the air flux. It displays the results obtained for the young subject (Fig. 12a) and for the elderly subject (Fig. 12b).
Figure 12 Numerical normal stress σ_{zz} as a function of the numerical deflection for the dwell times 1 s and 100 s for seven nodes of the central vertical under the air flux: (a) for the young subject, (b) for the elderly subject. 
During the numerical contactfree creep tests, the pressure distribution (simulating an air flow) applied at the top corner of the mesh provides mechanical energy to the material points of the REV for both subjects. For the dwell times 1 s and 100 s, the amounts of energy exchanged are extracted from the hysteresis curves of Figure 12 for the young subject and for the elderly subject. These amounts are given in Table 2 for the seven nodes of the central vertical under the air flow for the two subjects. The order of amplitude of the values given in Table 2 should be analyzed keeping in mind the dimensions of the REV (15 mm × 15 mm × 11.11 mm) and the amplitude of the external mechanical load (from 0 Pa to 7277 Pa). It should be noticed that for the elderly subject, the amount of mechanical energy of the dwell time of 100 s is 24.17 times higher than the amount for the dwell time of 1 s for the top of the SC, 21.06 times higher for the middle of SC, 20.51 times higher for the boundary SCVE, 21.51 times higher for the boundary VEDE, 20.14 times higher for the middle of DE, 17.62 times higher for the boundary DEHY and 18.5 times higher for the middle of HY.
Surface of the hysteresis cycles for the dwell times 1 s and 100 s for seven nodes of the central vertical under the air flux for the young subject and for the elderly subject.
Discussion
As explained when introducing the numerical results, they are analyzed in terms of evolutions of layers exemplified with the numerical results of seven nodes of the central vertical under the air flux.
In civil engineering, hysteresis cycles (Figs. 11 and 12) are studied to quantify the amount of mechanical energy supplied to the material during loading and the amount of mechanical energy given back or dissipated by the material during unloading.
Here the hysteresis cycles of the young subject for the central point of the top surface of SC under the air flux (Fig. 11) overlap with each other whatever the dwell time during which the load was applied. Each test supplied the same amount of mechanical energy (equal to 1.01 mm^{2}) to the top surface of SC. Further there is no storage of energy during the creep tests i.e. the starting and ending points of the cycles are the same. These mechanical findings mean that the evolutions of the material points of the top surface of SC are reversible and unique. In other words, let us consider a material point of the top surface of SC. Before one numerical creep test, it is in a given unknown initial state. During the test, the mechanical energy supplied by the load changes this state. After the test, the free evolution of this material point causes it to return to what appears to be its initial state (Fig. 11). This observation can be done for all material points of the top surface of SC and for all dwelltimes and whatever the amount of energy provided.
At the opposite, the hysteresis cycles (Fig. 11) of the elderly subject appear to be dwell time dependent. The starting and ending points are still the same for the five tests. However the surface of the cycles of the top node of SC increases when the dwell time increases. As the amount of mechanical energy provided by the test is equal to the surface of the cycle, increasing the surfaces is equivalent to increasing these amounts: 0.77 mm^{2} for 0.1 s, 0.86 mm^{2} for 0.5 s, 0.97 mm^{2} for 1 s, 1.20 mm^{2} for 2 s and 23.45 mm^{2} for 100 s. But the material points of the top surface of SC return to their initial state, i.e. the dissipated mechanical energy is equal to the mechanical energy supplied. Hence, each material point of the top surface of SC of the elderly subject uses all the supplied energy to return to its initial state after the load has been removed.
These findings are coherent with the qualitative observations done during the tests for both subjects i.e. the outer surface of the skin seems to come back to its initial position more or less quickly after the tests. The plus of the numerical simulations is to quantify the amount of energy exchanged during the tests and to be able to confirm that the skin outer surface comes really back to its mechanical initial state after the tests whatever can be the amount of energy provided. The numerical results have also shown that the evolutions of the outer surface of the skin are not influenced by the dwell time for the young subject while they are dwell time dependent for the elderly subject.
These questions remain open for the evolutions in the volume and will be discussed here after based on the other available numerical results.
It can be observed (Fig. 12) that the starting and ending points of a given hysteresis cycle are the same for both subjects. However each sublayer of the skin appears to have its own particular evolution as emphasized by the different hysteresis cycles drawn in Figure 12. They all contribute differently to the overall response of the skin. In a sense, mechanically speaking, each sublayer of the skin starts from an unknown initial state and returns to it.
For the young subject (Fig. 12a) each layer follows overlapping paths whatever the dwell time. This observation is strengthened for the two dwell times (1 s and 100 s) by the amounts of energy exchanged in each case by the seven nodes of the central vertical (Table 2). It should be noticed that either the amounts are equal or when they differ, the difference is not significant. Again the numerical results indicate that in the volume also the responses are not influenced by dwell times for the young subject. In other words, for the chosen REV, the response of the skin of the young subject is not influenced by the dwell time in all its volume when subjected to contactfree creep tests.
Observations are different for the elderly subject (Fig. 12b). It appears that each sublayer has a dwell time dependent evolution. When the dwell time is increased, the surfaces of the hysteresis cycles of the different sublayers increase (Table 2). These surfaces give a way to quantify the amount of mechanical energy supplied by the load. For each sublayer, the amount of the mechanical energy of the dwell time of 100 s is 20.50 ± 2.13 times higher than the one of the dwell time of 1 s. Recalling the main features of the physiology of the different sublayers it appears that the comparison of the different curves shows a link between the dwell time of the creep test and the viscous behavior of the skin. When the dwell time is increased, the viscous contribution of a given sublayer has time to develop and drastically change the response and the surfaces of its hysteresis cycles (Table 2). Again the surfaces give a way to quantify the amount of dissipated mechanical energy used by the sublayer to return to its unknown initial state after the load has been removed. The top surface of the SC faces the external load first and it is provided with the highest amount of supplied external mechanical energy. It is also, the sublayer which needs the highest amount of energy for coming back to its unknown initial state after the load is removed because it dissipates all the energy supplied. Contrary, in a Civil engineering sense, HY is the last sublayers reached by the charge transfer process. This sublayer receives the smallest amount of supplied external mechanical energy. When the load is removed, HY dissipates all this energy and comes back to its unknown initial state without needing extra energy for completing this evolution. In between these two responses, the other sublayers receive from the charge transfer process different amounts of external mechanical energy. Each sublayer uses entirely this amount for coming back when the load is removed. More, each sublayer satisfies itself with the amount available and complete its evolution taking it back to its unknown initial state without needing extraenergy. Hence, for the chosen REV, every layer of the elderly subject exhibits its own evolutions that are strongly dependent on the dwell time. These evolutions are triggered by the viscosity coefficient.
This work is a first contribution for establishing the causal relationship existing between mechanical fields and physiology associating the age of a healthy subject, his physiology, his mechanical viscoelastic characteristics and his numerical stress fields in the volume obtained through numerical simulations of in vivo contactfree creep tests.
Conclusion
In this paper we presented 3D finite element simulations of five numerical contactfree creep tests with dwell times of respectively 0.1 s, 0.5 s, 1 s, 2 s and 100 s . These numerical tests were based on experimental data recorded during contactfree creep tests performed on the volar forearm of a young adult subject and an elderly adult subject.
The numerical results showed different behaviors for the young subject and for the elderly subject.
For the young subject, the upper surface of the skin faces the external load with the same viscoelastic response whatever the dwell time of the numerical creep test. This evolution was also retrieved in the volume of the skin. The choice of dwell time appeared irrelevant and the response obtained was equivalent to that of an indentation test (dwell time of 0.1 s).
On the contrary, the skin of the elderly subject exhibited responses highly dependent on the amount of energy supplied during the tests. The sensorial feelings associated with the tests can be very different also. Therefore the optimal choice may be to work with a patient’s specific experimental protocol with a given dwell time and keep to it throughout the study.
Therefore, for a practitioner, the contactfree indenter is an easy way of characterizing the state of their patient’s skin before and after treatments.
In the context of medical survey, a practitioner needs an easy and robust way to characterize the state of his patient. He can perform an in vivo contactfree creep test with a short dwell time. A 1 s dwell time looks to be a good compromise between the needs of clinical tests and to be sure to provide sufficient amount of information to the coupled experimentaltheoreticalnumerical approach.
Conflict of interest
The authors declare that they have no conflicts of interest in relation to this article.
Ethical statement
Each volunteer has provided his informed consent before participating.
References
 Sheppard C, Shotton D (1997), Confocal laser scanning microscopy. Royal Microsc. Soc. Microsc. Handbooks, Bios Scientific Publishers, Oxford. [Google Scholar]
 Ginefri J, Darasse I, Crozat P (2001), Hightemperature superconducting surface coil for in vivo microimaging of the human skin. Magnetic Res Med 45, 376–382. [CrossRef] [PubMed] [Google Scholar]
 Tran HV (2007), Caractérisation des propriétés mécaniques de la peau humaine in vivo via l’IRM. Mechanics, Université de Technologie de Compiègne, Compiègne, France. [Google Scholar]
 Norlén L, AlAmoudi A (2004), Stratum corneum keratin structure, function, and formation: The cubic rodpacking and membrane templating model. J Invest Dermatol 123, 715–732. [CrossRef] [PubMed] [Google Scholar]
 Myer K, Maibach H (2013), Stratum corneum evaluation methods: Overview. Skin Res Technol 19, 213–219. [CrossRef] [PubMed] [Google Scholar]
 Daly CH, Odland GF (1979), Agerelated change in mechanical properties of human skin. J Invest Dermatol 73, 84–87. [CrossRef] [PubMed] [Google Scholar]
 Agache P, Monneur C, Lévêque J, De Rigal J (1980), Mechanical properties and young’s modulus of human skin in vivo. Arch Dermatol Res 269, 221–232. [CrossRef] [PubMed] [Google Scholar]
 Geerligs M (2009), Skin layer mechanics, Dissertation, Eindhoven University of Technology, Eindhoven, Netherlands. [Google Scholar]
 Leveque JL, de Rigal J, Agache PG, Monneur C (1980), Influence of ageing on the in vivo extensibility of human skin at a low stress. Arch Dermatol Res 269, 2, 127–135. [CrossRef] [PubMed] [Google Scholar]
 Rochefort A (1986), Propriétés biomécaniques du stratum corneum : – Modélisation rhéologique – Application à la cosmetology, Thèse de Doctorat, Université de FrancheComté, Besançon. [Google Scholar]
 Benitez JM, Montans FJ (2017), The mechanical behavior of skin: Structures and models for the finite element analysis. Comput Struct 190, 75–107. [CrossRef] [Google Scholar]
 Ayadh M, Abellan MA, Didier C, Bigouret A, Zahouani H (2020), Methods for characterizing the anisotropic behavior of the human skin’s relief and its mechanical properties in vivo linked to age effects. Surf Topogr: Metrol Prop 8, 1, 014002. [CrossRef] [Google Scholar]
 Moerman KM, Holt CA, Evans SL, Simms CK (2009), Digital image correlation and finite element modeling as a method to determine mechanical properties of human soft tissue in vivo. J Biomech 42, 1150–1153. [CrossRef] [PubMed] [Google Scholar]
 Hendriks FM, Brokken D, Oomens CWJ, Bader DL, Baaijens FPT (2006), The relative contributions of different skin layers to the mechanical behavior of human skin in vivo using suction experiments. Med Eng Phys 28, 259–266. [CrossRef] [PubMed] [Google Scholar]
 Khatyr KF, Imberdis C, Varchon D, Lagarde JM, Josse G (2006), Measurement of the mechanical properties of the skin using the suction test. Comparison between three methods: Geometric, Timoshenko and finite elements. Skin Res Technol 12, 24–31. [CrossRef] [PubMed] [Google Scholar]
 Krueger N, Luebberding S, Oltme M, Streker M, Kerscher M (2011), Agerelated changes in skin mechanical properties: A quantitative evaluation of 120 female subjects. Skin Res Technol 17, 141–148. [CrossRef] [PubMed] [Google Scholar]
 Hara Y, Masuda Y, Hirao T, Yoshikawa N (2013), The relationship between the Young’s modulus of the stratum corneum and age: A pilot study. Skin Res Technol 19, 339–345. [CrossRef] [PubMed] [Google Scholar]
 Delalleau A, Josse G, Lagarde JM, Zahouani H, Bergheau JM (2006), Characterization of the mechanical properties of skin by inverse analysis combined with the indentation test. J Biomech 39, 1603–1610. [CrossRef] [PubMed] [Google Scholar]
 Abellan MA, Zahouani Z, Bergheau JM (2013), Contribution to the determination of in vivo mechanical characteristics of human skin by indentation test. Comput Math Meth Med 2013, 11, Article ID 814025. [Google Scholar]
 Boyer G (2010), Modélisation du comportement mécanique de la peau humaine in vivo : Application au vieillissement et aux gestes du clinicien, Thèse de Doctorat, Ecole Nationale Supérieure des Mines, SaintEtienne. [Google Scholar]
 Abellan MA, Ayadh M, Feulvarch E, Zahouani Z, Bergheau JM (2015), Caractérisation des paramètres mécaniques de la peau humaine jeune in vivo par essais d’indentation sans contact, in: Congrès Français de Mécanique, 24–28 août, Lyon. [Google Scholar]
 Diridollou S, Patat F, Gens F, Vaillant L, Black D, Lagarde JM, Gall Y, Berson M (2000), In vivo model of the mechanical properties of the human skin under suction. Skin Res Technol 6, 214–221. [CrossRef] [PubMed] [Google Scholar]
 Bischoff JE, Arruda EM, Grosh K (2004), A rheological network model for the continuum anisotropic and viscoelastic behavior of soft tissue. Biomech Mod Mechanobiol 3, 56–65. [CrossRef] [PubMed] [Google Scholar]
 Azzez K, Bergheau JM, Dogui A, Zahouani H, Abellan MA, Chaabane M (2015), Contribution à l’étude du vieillissement de la peau humaine in vivo par simulations numériques d’essais d’indentation, in: Colloque National en Calculs des Structures, 18–22 mai, Giens. [Google Scholar]
 Flynn C, Taberner A, Nielsen P (2011), Modeling the mechanical response of in vivo human skin under a rich set of deformations. Ann Biomed Eng 39, 7, 1935–1946. [CrossRef] [PubMed] [Google Scholar]
 Abellan MA, Bergheau JM, Zahouani Z (2014), Comparaison de différents modèles viscoélastiques pour la caractérisation des propriétés mécaniques de la peau humaine in vivo par essai d’indentation. Comp Meth Biomech Biomed Eng 17, Suppl. 1, 22–23. [CrossRef] [PubMed] [Google Scholar]
 Jouanna P, Abellan MA (1995), Chapter 1. Generalized approach to heterogeneous media, Modern issues in nonsaturated soils, in: A Gens, P Jouanna, BA Schrefler (Eds.), CISM courses and lectures N°357 (International Center for Mechanical Sciences), SpringerVerlag, Wien NewYork, ISBN 5211827838. [Google Scholar]
 Bathe KJ (1996), Finite element procedures, Prentice Hall, Upper Saddle River, New Jersey. [Google Scholar]
 Latorre M, Montans FJ (2016), Stress and strain mapping tensors and general workconjugacy in large strain continuum mechanics. Appl Math Model 40, 3938–3950. [CrossRef] [Google Scholar]
 Zahouani H, PaillerMattei C, Sohm B, Vargiolu R, Cenizo V, Debret R (2009), Characterization of the mechanical properties of a dermal equivalent compared with human skin in vivo by indentation and static friction tests. Skin Res Technol 15, 68–76. [CrossRef] [PubMed] [Google Scholar]
 Ayadh M, Abellan MA, Chatelin R, Bergheau JM, Zahouani H (2017), Contribution à la simulation numérique du comportement viscoélastique de la peau humaine jeune et âgée in vivo, in: Colloque National en Calculs des Structures, 15–19 mai, Giens. [Google Scholar]
 Azzez K, Chaabane M, Abellan MA, Bergheau JM, Zahouani H, Dogui A (2015), In vivo identification of mechanical properties of human skin: Effect of initial viscoelastic parameters on indentation response, in: VI International Congress on Design and Modeling of Mechanical Systems (CMSM), 23–25 march, Hammamet, Tunisia, paper 296. [Google Scholar]
Cite this article as: Abellan MA, Ayadh M, Bergheau JM & Zahouani H 2023. Dwell time in contactfree creep tests plays an agedependent role in the viscoelastic behavior of in vivo human skin. 4open, 6, 7.
All Tables
Examples of sets of parameters checked during the inverse approach showing the best fit for the young subject (Set 5) and for the elderly subject (Set 111). The instantaneous Young modulus is defined by E_{0} = E_{1} + E_{∞} in MPa.
Surface of the hysteresis cycles for the dwell times 1 s and 100 s for seven nodes of the central vertical under the air flux for the young subject and for the elderly subject.
All Figures
Figure 1 (a) Principle of the Tonoderm^{®} [20]. (b) Forearm setup [31]. (c) Photo of the pressured area of the skin surface under the air flow [20]. (d) Invertedbell shape of the air flow and typical contact profile which may develop. 

In the text 
Figure 2 Experimental loading curve as a function of time (Exp). 

In the text 
Figure 3 Loading curves of the five numerical contactfree creep tests with dwell times of 0.1 s, 0.5 s, 1 s, 2 s and 100 s. 

In the text 
Figure 4 Rheological Zener model. 

In the text 
Figure 5 Userdesigned REV deduced from the echograph images and 3D mesh. 

In the text 
Figure 6 3D mesh with the boundary conditions. 

In the text 
Figure 7 Profiles of the experimental deflection as a function of time for the young subject (blue line) and for the elderly subject (red line) for a flux of 20 NL/min. 

In the text 
Figure 8 (Top) Comparison of the experimental and numerical curves of the displacement for three sets of parameters (Bottom) Curves of the relative error for the same sets of parameters. 

In the text 
Figure 9 Sketch giving the position of the seven nodes of the central vertical under the air flux used in the analysis. 

In the text 
Figure 10 Isovalues of the vertical displacements at the end of the dwell times: (a) 0.1 s, (b) 0.5 s and (c) 1 s for the young subject. (d) Contour lines of the mean stress for the dwell time of 0.5 s for the elderly subject. 

In the text 
Figure 11 Numerical normal stress σ_{zz} as a function of the numerical deflection of the central point under the air flux for the young subject and for the elderly subject for the dwell times of 0.1 s, 0.5 s, 1 s, 2 s and 100 s. 

In the text 
Figure 12 Numerical normal stress σ_{zz} as a function of the numerical deflection for the dwell times 1 s and 100 s for seven nodes of the central vertical under the air flux: (a) for the young subject, (b) for the elderly subject. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.