Biomechanical Modeling of Wounded Skin

Skin injury is the most common type of injury, which manifests itself in the form of wounds and cuts. A minor wound repairs itself within a short span of time. However, deep wounds require adequate care and sometime clinical interventions such as surgical suturing for their timely closure and healing. In literature, mechanical properties of skin and other tissues are well known. However, the anisotropic behavior of wounded skin has not been studied yet, specifically with respect to localized overstraining and possibilities of rupture. In the current work, the biomechanics of common skin wound geometries were studied with a biofidelic skin phantom, using uniaxial mechanical testing and Digital Image Correlation (DIC). Global and local mechanical properties were investigated, and possibilities of rupture due to localized overstraining were studied across different wound geometries and locations. Based on the experiments, a finite element (FE) model was developed for a common elliptical skin wound geometry. The fidelity of this FE model was evaluated with simulation of uniaxial tension tests. The induced strain distributions and stress-stretch responses of the FE model correlated very well with the experiments (R2 > 0.95). This model would be useful for prediction of the mechanical response of common wound geometries, especially with respect to their chances of rupture due to localized overstraining. This knowledge would be indispensable for pre-surgical planning, and also in robotic surgeries, for selection of appropriate wound closure techniques, which do not overstrain the skin tissue or initiate tearing.


Introduction
Skin wounds are caused due to rupture of one or more layers of the skin [1]. Figure 1 shows the three layers of the skin, namely the epidermis (outermost layer), dermis (middle layer) and hypodermis (bottom layer) and a dermal skin wound. Bleeding occurs mainly due to the rupture of blood vessels which extend from the hypodermis through some parts of the dermis. Wounds may be classified broadly as shallow (epidermal and dermal) and deep (dermal and hypodermal). Typically, deep wounds require surgical intervention and suturing for timely healing [2]. To date, suturing is based on a surgeon's experience and there is a wide variation in suturing techniques from one surgeon to another [3]. Better understanding and characterization of the biomechanics of wounds may help in valuable improvements in wound healing and also provide benchmarks and guidelines for robotic surgeries involving wound closure with sutures.
In literature, the mechanics of wound closure has been studied using computational modeling [4]. One of the earliest numerical wound closure models was developed by Danielson and Natarajan et al. [5], which involved the estimation of maximum and minimum stress zones using the tension field theory. The first 2D finite element (FE) model of wound closure was developed by time. Skin pre-stress which occurs due to the natural tension of the skin has been observed to affect skin wrinkling and wound closure significantly [8]. Cavicchi et al. [7] considered the skin as a 2D surface with isotropic material property and investigated the effect of equi-biaxial pre-stress on two cases of circular wounds (small and large). From the analyses, it was observed that maximum stress around the large circular wound was approximately 200 kPa higher than that of the small circular wound. It was concluded that the wound size also affects the stresses around the wounds. Flynn et al. [9] in 2010 developed an FE model to quantify the wound closure force requirements for symmetric wound geometries. A skin pre-stress of 10 kPa (x-axis) and 5 kPa (y-axis) was applied for four wound geometries namely circular, fusiform, elliptical and lazy S-plasty. From the analyses, the elliptical geometry was found to have the lowest wound closure force requirement. Capek et al. [10] in 2012 further extended the study on wound closure force estimation of elliptical wounds with varying sizes, with and without skin pre-stress effect. It was observed that the wound closure force requirements had a strong dependence on both the size of the elliptical wound and skin pre-stress. Feng et al. [11] investigated the impact of neck incision location on total laryngectomy procedure outcome, using FE analysis. Two separate neck incision models namely low-neck apron (involving the tracheostoma) and mid-neck apron (not involving the tracheostoma), were investigated under different neck loading conditions. From the study, the mid-neck apron model was observed to result in up to 11 times lower von Mises stresses compared to the low-neck apron model, which complied with clinical observations. Chanda et al. [12,13] in 2017 developed a novel 3D computational model of a multilayer skin. A symmetric dermal wound was created and sutures were placed using a novel suture pulling technique, and in the presence of biaxial skin pre-stress. A seven-step interrupted suturing process was simulated and force requirement for each suture were estimated. A mathematical relation was also developed relating the wound flap displacement and suture force requirement, which is applicable to a range of symmetric wound geometry and sizes. A few experimental studies have been conducted in literature to model wound closure and suture. Melis et al. [14] in 2001 conducted suturing experiments on 9 cm × 9 cm rectangular wounds created on skin samples of piglets, and manually estimated suture force requirements. From the study, a mean tensile force of 22.66 N was recorded for sutures with the maximum at the center suture. Freck et al. [15] in 2001 investigated suturing on sheep skin. The average axial force estimated for suture was 1.9 N. Subcuticular suturing was performed by Cacou et al. [16,17] on human skin specimens with different wound geometries. It was concluded from the study that the wound closure force has strong dependence on wound geometry, with elliptical geometry leading to lowest suture force requirements. Corr et al. [18] reviewed studies on the similarities and differences in mechanical properties of scar tissue and uninjured skin. Stress-relaxation experiments and tensile failure testing were found to indicate the changes in the directional properties of the tissue in healthy state versus in injured state. Specifically, compared with uninjured skin, the scar tissue was observed to exhibit similar high-load stiffness, significantly less resistance to failure, lower low-load compliance, and changes in the material directionality. Chanda et al. [19] in 2017 also developed a biofidelic human skin tissue phantom material and created different wound geometries to study the effect of wound geometry and size on the average skin response in a tensile test. Actual sutures were placed to close the wounds, and the average stress-strain responses were compared with the mechanical behavior of wounded specimens without sutures. It was observed from the experiments that, with sutures, the A few experimental studies have been conducted in literature to model wound closure and suture. Melis et al. [14] in 2001 conducted suturing experiments on 9 cm × 9 cm rectangular wounds created on skin samples of piglets, and manually estimated suture force requirements. From the study, a mean tensile force of 22.66 N was recorded for sutures with the maximum at the center suture. Freck et al. [15] in 2001 investigated suturing on sheep skin. The average axial force estimated for suture was 1.9 N. Subcuticular suturing was performed by Cacou et al. [16,17] on human skin specimens with different wound geometries. It was concluded from the study that the wound closure force has strong dependence on wound geometry, with elliptical geometry leading to lowest suture force requirements. Corr et al. [18] reviewed studies on the similarities and differences in mechanical properties of scar tissue and uninjured skin. Stress-relaxation experiments and tensile failure testing were found to indicate the changes in the directional properties of the tissue in healthy state versus in injured state. Specifically, compared with uninjured skin, the scar tissue was observed to exhibit similar high-load stiffness, significantly less resistance to failure, lower low-load compliance, and changes in the material directionality. Chanda et al. [19] in 2017 also developed a biofidelic human skin tissue phantom material and created different wound geometries to study the effect of wound geometry and size on the average skin response in a tensile test. Actual sutures were placed to close the wounds, and the average stress-strain responses were compared with the mechanical behavior of wounded specimens without sutures. It was observed from the experiments that, with sutures, the average skin response is similar to unwounded skin. Local changes around the wounds were not quantified and the influence of geometrical parameters of the wound on the surrounding skin straining was not characterized.
To date, there is a lack of knowledge on localized straining of skin wounds. This inhibits the understanding of localized overstraining and possibilities of rupture in wounds during skin surgeries. To address this key gap in literature, a novel experimental framework was developed to study skin straining due to common wound geometries, using uniaxial mechanical testing and Digital Image Correlation (DIC). A biofidelic skin phantom [20][21][22][23] was employed, and common wound geometries and sizes were simulated, to characterize the induced local strains. Additionally, a biofidelic computational model was developed to simulate the experimental results. This model would be valuable for prediction of localized overstraining and skin ruptures for common wound geometries in surgeries. Additionally, this knowledge may not only guide the development of better pre-surgical planning techniques, but can also be applied to robotic surgeries for selection of wound closure techniques which inhibit skin tissue overstraining and rupture. Section 2 will discuss the materials and methods of the experimental and DIC modeling framework, along with the development of the computational model. Section 3 will discuss the various results followed by conclusions in Section 4.

Test Specimen Fabrication and Digital Image Correlation (DIC) Framework
Two-part silicone material with a shore hardness of 30A (Mold Star 30 from Smooth-On, Inc., Macungie, PA, USA) was mixed in a 1:1 ratio by weight to fabricate biofidelic human skin phantoms [19,[24][25][26][27][28][29][30][31]. Ten test coupons were generated with 50 mm length, 9 mm width and 2.5 mm depth. Common elliptical wound geometries [10] of varying dimensions were created using a scalpel and regular scissors on 5 out of 10 coupons (Figure 2a). The direction along the length of the coupon was considered as the y-axis and the orthogonal direction along the width of the coupon was considered the x-axis. The wounds were classified based on their aspect ratios (AR, ratio of width along x-direction and length along y-direction, see Figure 2a) along the axes as Wound 1 (AR = 0.20, symmetric-short), Wound 2 (AR = 1.50, symmetric-long), Wound 3 (AR = 0.14, asymmetric-wide), Wound 4 (AR = 0.16, asymmetric-thin), and Wound 5 (AR = 0.07, transverse) respectively to study the effect of wound dimension on the strain fields. average skin response is similar to unwounded skin. Local changes around the wounds were not quantified and the influence of geometrical parameters of the wound on the surrounding skin straining was not characterized. To date, there is a lack of knowledge on localized straining of skin wounds. This inhibits the understanding of localized overstraining and possibilities of rupture in wounds during skin surgeries. To address this key gap in literature, a novel experimental framework was developed to study skin straining due to common wound geometries, using uniaxial mechanical testing and Digital Image Correlation (DIC). A biofidelic skin phantom [20][21][22][23] was employed, and common wound geometries and sizes were simulated, to characterize the induced local strains. Additionally, a biofidelic computational model was developed to simulate the experimental results. This model would be valuable for prediction of localized overstraining and skin ruptures for common wound geometries in surgeries. Additionally, this knowledge may not only guide the development of better pre-surgical planning techniques, but can also be applied to robotic surgeries for selection of wound closure techniques which inhibit skin tissue overstraining and rupture. Section 2 will discuss the materials and methods of the experimental and DIC modeling framework, along with the development of the computational model. Section 3 will discuss the various results followed by conclusions in Section 4.

Test Specimen Fabrication and Digital Image Correlation (DIC) Framework
Two-part silicone material with a shore hardness of 30A (Mold Star 30 from Smooth-On, Inc., Macungie, PA, USA) was mixed in a 1:1 ratio by weight to fabricate biofidelic human skin phantoms [19,[24][25][26][27][28][29][30][31]. Ten test coupons were generated with 50 mm length, 9 mm width and 2.5 mm depth. Common elliptical wound geometries [10] of varying dimensions were created using a scalpel and regular scissors on 5 out of 10 coupons (Figure 2a). The direction along the length of the coupon was considered as the y-axis and the orthogonal direction along the width of the coupon was considered the x-axis. The wounds were classified based on their aspect ratios (AR, ratio of width along x-direction and length along y-direction, see Figure 2a) along the axes as Wound 1 (AR = 0.20, symmetric-short), Wound 2 (AR = 1.50, symmetric-long), Wound 3 (AR = 0.14, asymmetric-wide), Wound 4 (AR = 0.16, asymmetric-thin), and Wound 5 (AR = 0.07, transverse) respectively to study the effect of wound dimension on the strain fields.
(a) In literature, DIC techniques have been used especially in micro and nano scale mechanical testing applications to understand the localization of strains in a material under varying loading conditions [32]. Especially, DIC have been used to understand crack initiation and propagation in hard composite materials [33]. Recently, DIC was used to understand anisotropic local skin tissue deformation [34]. However, to date, DIC has not been used to study soft materials under damage conditions.
In this work, a certain set of steps were followed to conduct the DIC experiments. Spray paint (KILZ White Interior/Exterior paint from Home Depot, Tuscaloosa, AL, USA) was dispersed on the surface of the test coupon in such a way that the speckle patterns are fine, dissimilar, and densely populates the surface to be studied ( Figure 2b). This process was repeated six times to obtain the desired speckle pattern [32]. A high-quality camera (20.1 Mega Pixel Sony DSC-H400 Single-Lens Reflex (SLR) camera, Tokyo, Japan) was used on a tripod to record the tensile tests. The videos were run through a video-to-image conversion software to obtain 12-14 images from each test. The average length of the videos recorded was 48 s. Each test was conducted and recorded three times to ensure repeatability.
Visual Image Correlation (VIC) 2D software (from Correlated Solutions Inc., Columbia, SC, USA) [32] was used for post processing of the image sets from the tensile tests. VIC-2D analyses subtle differences between two images taken before and after the deformation by a solid-state digital camera. Strains are calculated based on the spray particle displacements measured by VIC-2D throughout the different deformation phases during the tensile test. According to Bruck et al.'s numerical model of the deformations observed in a DIC study [32], the strain calculations along the length and width of a tensile test coupon is estimated using the equations presented in Figure 3. In literature, DIC techniques have been used especially in micro and nano scale mechanical testing applications to understand the localization of strains in a material under varying loading conditions [32]. Especially, DIC have been used to understand crack initiation and propagation in hard composite materials [33]. Recently, DIC was used to understand anisotropic local skin tissue deformation [34]. However, to date, DIC has not been used to study soft materials under damage conditions.
In this work, a certain set of steps were followed to conduct the DIC experiments. Spray paint (KILZ White Interior/Exterior paint from Home Depot, Tuscaloosa, AL, USA) was dispersed on the surface of the test coupon in such a way that the speckle patterns are fine, dissimilar, and densely populates the surface to be studied ( Figure 2b). This process was repeated six times to obtain the desired speckle pattern [32]. A high-quality camera (20.1 Mega Pixel Sony DSC-H400 Single-Lens Reflex (SLR) camera, Tokyo, Japan) was used on a tripod to record the tensile tests. The videos were run through a video-to-image conversion software to obtain 12-14 images from each test. The average length of the videos recorded was 48 s. Each test was conducted and recorded three times to ensure repeatability.
Visual Image Correlation (VIC) 2D software (from Correlated Solutions Inc., Columbia, SC, USA) [32] was used for post processing of the image sets from the tensile tests. VIC-2D analyses subtle differences between two images taken before and after the deformation by a solid-state digital camera. Strains are calculated based on the spray particle displacements measured by VIC-2D throughout the different deformation phases during the tensile test. According to Bruck et al.'s numerical model of the deformations observed in a DIC study [32], the strain calculations along the length and width of a tensile test coupon is estimated using the equations presented in Figure 3.

Mechanical Tensile Testing
Uniaxial mechanical testing was conducted on all biofidelic human skin phantoms using a universal testing machine (MTS Criterion Model 42, Eden Prairie, MN, USA) to measure the load versus deformation behavior under tension. A constant strain rate of 0.4 mm/s was used to conduct all experiments based on literature [34][35][36][37][38]. The standard size chosen for all test coupons was 34 mm × 9 mm × 2.5 mm (after clamping) [34][35][36][37][38]. Five coupons with no wounds were tested first to obtain baseline measurements and to calibrate the VIC 2D system. After baseline measurements, tensile tests ( Figure 4) were conducted on test specimens with wounds and strain changes with time were quantified using VIC 2D. Overstraining, defined as the occurrence of a high strain zone (i.e., strain ≥25% of average strain [4]), was identified during the tests. Possibility of rupture was indicated for incidences of high strain zones surrounded by low strain zones (i.e., strain ≤25% of average strain [4]).

Finite Element Modeling
The computational model of the coupon was developed based on the size of the actual coupons used in the experiments (50 mm length, 9 mm width and 2.5 mm depth). Figure 5 shows the transverse wound geometry observed under a microscope. The direction along the length of the coupon was considered to be the y-axis and along the width to be the x-axis. The geometries of the coupon and the wound were measured and modeled ( Figure 5) within the finite element software

Mechanical Tensile Testing
Uniaxial mechanical testing was conducted on all biofidelic human skin phantoms using a universal testing machine (MTS Criterion Model 42, Eden Prairie, MN, USA) to measure the load versus deformation behavior under tension. A constant strain rate of 0.4 mm/s was used to conduct all experiments based on literature [34][35][36][37][38]. The standard size chosen for all test coupons was 34 mm × 9 mm × 2.5 mm (after clamping) [34][35][36][37][38]. Five coupons with no wounds were tested first to obtain baseline measurements and to calibrate the VIC 2D system. After baseline measurements, tensile tests ( Figure 4) were conducted on test specimens with wounds and strain changes with time were quantified using VIC 2D. Overstraining, defined as the occurrence of a high strain zone (i.e., strain ≥25% of average strain [4]), was identified during the tests. Possibility of rupture was indicated for incidences of high strain zones surrounded by low strain zones (i.e., strain ≤25% of average strain [4]).

Mechanical Tensile Testing
Uniaxial mechanical testing was conducted on all biofidelic human skin phantoms using a universal testing machine (MTS Criterion Model 42, Eden Prairie, MN, USA) to measure the load versus deformation behavior under tension. A constant strain rate of 0.4 mm/s was used to conduct all experiments based on literature [34][35][36][37][38]. The standard size chosen for all test coupons was 34 mm × 9 mm × 2.5 mm (after clamping) [34][35][36][37][38]. Five coupons with no wounds were tested first to obtain baseline measurements and to calibrate the VIC 2D system. After baseline measurements, tensile tests ( Figure 4) were conducted on test specimens with wounds and strain changes with time were quantified using VIC 2D. Overstraining, defined as the occurrence of a high strain zone (i.e., strain ≥25% of average strain [4]), was identified during the tests. Possibility of rupture was indicated for incidences of high strain zones surrounded by low strain zones (i.e., strain ≤25% of average strain [4]).

Finite Element Modeling
The computational model of the coupon was developed based on the size of the actual coupons used in the experiments (50 mm length, 9 mm width and 2.5 mm depth). Figure 5 shows the transverse wound geometry observed under a microscope. The direction along the length of the coupon was considered to be the y-axis and along the width to be the x-axis. The geometries of the coupon and the wound were measured and modeled ( Figure 5) within the finite element software

Finite Element Modeling
The computational model of the coupon was developed based on the size of the actual coupons used in the experiments (50 mm length, 9 mm width and 2.5 mm depth). Figure 5 shows the transverse wound geometry observed under a microscope. The direction along the length of the coupon was considered to be the y-axis and along the width to be the x-axis. The geometries of the coupon and the wound were measured and modeled ( Figure 5) within the finite element software ANSYS (ANSYS 17.1, Canonsburg, PA, USA). Meshing was conducted using the higher order 3-D 20-node Solid 186 element, which perform well with large deflections, hyperelasticity, and large strain deformations [39,40]. The test specimen was constrained in all degrees of freedom at the bottom end, and pulled from the top end. Tensile test with a strain rate of 0.4 mm/s was simulated by transferring the Displacement-Time data from the experiments into ANSYS dynamic FE solver. The stress-strain plots obtained from experiments were compared with the results from the computational analysis. Also, the local strains were analyzed and the results were compared with the DIC based measurements.  [39,40]. The test specimen was constrained in all degrees of freedom at the bottom end, and pulled from the top end. Tensile test with a strain rate of 0.4 mm/s was simulated by transferring the Displacement-Time data from the experiments into ANSYS dynamic FE solver. The stress-strain plots obtained from experiments were compared with the results from the computational analysis. Also, the local strains were analyzed and the results were compared with the DIC based measurements.

Nonlinear Material Modeling
Soft tissue surrogates such as the human skin phantom exhibits nonlinear biomechanical behavior, which can be modeled using isotropic hyperelastic formulations such as Fung, Mooney-Rivlin, Ogden, Yeoh, and Veronda-Westmann [13,20,25,[41][42][43][44][45][46][47]. Hyperelastic models are based on the material dependent strain-energy function (denoted as  ) [48,49]. The strain energy function depends on either the principal stretches ( 1  , 2  and 3  ) along the x, y and z Cartesian coordinate axes, or the Cauchy-green tensor invariants namely 1 I , 2 I and 3 I (also functions of the principal stretches) [41] given by 1 2 3 ( , , ) In the current work, the Veronda-Westmann's hyperelastic model was used to characterize the nonlinear mechanical behavior of the human skin phantom. The strain energy functions of Veronda-Westmann's model is given by Equation (1).
The principal Cauchy stress is expressed in terms of the stretch and the strain energy function using Equation (2). Further, using Equations (1) and (2), Equation (3) is obtained. Using experimental stress-strain data from uniaxial tests in the current work, the hyperelastic curve fit parameters 1 c and 2 c were predicted using Equation (3) (see Table 1) with an R 2 correlation value of 0.982.
The principal Cauchy stress is expressed in terms of the stretch and the strain energy function using Equation (2). Further, using Equations (1) and (2), Equation (3) is obtained. Using experimental stress-strain data from uniaxial tests in the current work, the hyperelastic curve fit parameters c 1 and c 2 were predicted using Equation (3) (see Table 1) with an R 2 correlation value of 0.982.

Strain in Skin Phantom without Wound
Skin phantom specimens were tested in tension and DIC was used to model the strain distribution on the specimen surface ( Figure 6). To quantify the baseline strain versus time data in the skin phantom specimens, strains were measured using DIC in a skin phantom specimen without wound, and compared with the strains recorded using the experimental device. Figure 6 shows the average strain field post-processed using DIC. A total of 14 frames were extracted from the 40 s long video of the tensile test. Markers were placed at 6 different locations to obtain the average strain versus frames (Figure 7) for five coupons each tested twice. It should be mentioned here that the plot axes have been modified for all Strain (ε y ) versus Time Frame, from eyy [1]-Lagrange versus Index [1] for Figure 7 and following figures for clarity. A high R 2 correlation of 0.976 was estimated for the average strain versus time frame results and the data from the experimental device, and the overall error in DIC strain measurements was quantified to be 2.3 ± 0.3. The maximum difference in the local strains recorded at different locations of the specimen, and the average strain was 0.11 (14%), indicating no overstraining or rupture potential.

Strain in Skin Phantom without Wound
Skin phantom specimens were tested in tension and DIC was used to model the strain distribution on the specimen surface ( Figure 6). To quantify the baseline strain versus time data in the skin phantom specimens, strains were measured using DIC in a skin phantom specimen without wound, and compared with the strains recorded using the experimental device. Figure 6 shows the average strain field post-processed using DIC. A total of 14 frames were extracted from the 40 s long video of the tensile test. Markers were placed at 6 different locations to obtain the average strain versus frames (Figure 7) for five coupons each tested twice. It should be mentioned here that the plot axes have been modified for all Strain ( y  ) versus Time Frame, from eyy [1]-Lagrange versus Index [1] for Figure 7 and following figures for clarity. A high R 2 correlation of 0.976 was estimated for the average strain versus time frame results and the data from the experimental device, and the overall error in DIC strain measurements was quantified to be 2.3 ± 0.3. The maximum difference in the local strains recorded at different locations of the specimen, and the average strain was 0.11 (14%), indicating no overstraining or rupture potential.

Strains in Skin Phantom with Symmetric-Short Wound
Tensile tests were conducted on the specimen with symmetric and short wound geometry (Wound 1), and DIC was used to study the local strains around the wound perimeter. A significant strain build-up was observed around the wound periphery majorly concentrated on the left and right sides ( Figure 8). Markers were placed at the different locations around the wound periphery ( Figure 9) and the specific strain values were estimated. It was observed that the average strain value for the test specimen was lower than that at the points of strain build-up. The maximum strain recorded was approximately 31% greater than the overall average strain, indicating possible overstraining. The strain values recorded at the top and bottom points of the wound were lower than that of the average strain observed for the specimen. Also, the strain values recorded on the left and right side of the wound (at points P4 and P5 in the Figure 9) were similar, which could be due to the near symmetric elliptical geometry of the wound. No possible rupture sites, characterized by the incidence of high strain zones surrounded by low strain zones, were observed.

Strains in Skin Phantom with Symmetric-Short Wound
Tensile tests were conducted on the specimen with symmetric and short wound geometry (Wound 1), and DIC was used to study the local strains around the wound perimeter. A significant strain build-up was observed around the wound periphery majorly concentrated on the left and right sides ( Figure 8). Markers were placed at the different locations around the wound periphery ( Figure 9) and the specific strain values were estimated. It was observed that the average strain value for the test specimen was lower than that at the points of strain build-up. The maximum strain recorded was approximately 31% greater than the overall average strain, indicating possible overstraining. The strain values recorded at the top and bottom points of the wound were lower than that of the average strain observed for the specimen. Also, the strain values recorded on the left and right side of the wound (at points P4 and P5 in the Figure 9) were similar, which could be due to the near symmetric elliptical geometry of the wound. No possible rupture sites, characterized by the incidence of high strain zones surrounded by low strain zones, were observed.    . Average strain (y-axis) versus time frame (x-axis) measured using DIC from tensile test on the skin phantom specimen with symmetric-short wound 1.

Figure 9.
Average strain (y-axis) versus time frame (x-axis) measured using DIC from tensile test on the skin phantom specimen with symmetric-short wound 1.

Strains in Skin Phantom with Symmetric-Long Wound
The specimen with symmetric and long wound geometry (Wound 2) was subjected to tensile strain and 14 frames were extracted from 49 s visual data ( Figure 10). A high stress build-up was observed on the left and right side of the wound to start with, which shifted dominantly to the right side with increasing strains. Using 6 markers along the circumference of the wound as shown in the inset in Figure 11, the lowest strains were estimated at the top and bottom of the wound. All the markers along the left side of the wound (P1, P3 and P4) resulted in similar strain changes as seen in the Figure 11. The highest strains were recorded for the marker P2 on the right side of the wound, the maximum value of which was 34% higher than the average strain (indicating overstraining). Additionally, the strain at the bottom was observed to increase with time (with a maximum value of >30% over the average) with a zone of low strains (with a minimum value of <25% below the average strain) around it, which may be a possible rupture site, with high chances of skin tearing.

Strains in Skin Phantom with Symmetric-Long Wound
The specimen with symmetric and long wound geometry (Wound 2) was subjected to tensile strain and 14 frames were extracted from 49 s visual data (Figure 10). A high stress build-up was observed on the left and right side of the wound to start with, which shifted dominantly to the right side with increasing strains. Using 6 markers along the circumference of the wound as shown in the inset in Figure 11, the lowest strains were estimated at the top and bottom of the wound. All the markers along the left side of the wound (P1, P3 and P4) resulted in similar strain changes as seen in the Figure 11. The highest strains were recorded for the marker P2 on the right side of the wound, the maximum value of which was 34% higher than the average strain (indicating overstraining). Additionally, the strain at the bottom was observed to increase with time (with a maximum value of >30% over the average) with a zone of low strains (with a minimum value of <25% below the average strain) around it, which may be a possible rupture site, with high chances of skin tearing.

Strains in Skin Phantom with Asymmetric-Wide Wound
The test specimen with asymmetric and wide wound geometry (Wound 3) was subjected to tensile loads for 42 s and 14 frames were extracted for the DIC measurements. A high stress build-up was observed initially on the left-hand side of the wound which shifted eventually to the right-hand side (see Figure 12). On quantification of the strains along the circumference of the long wound with 7 markers (Figure 13), the strain fields were found to be close to each other all around the wound with the peak at point P2 initially which was on the left-hand side. After 8 s, the maximum strain was recorded for marker P5, which was within 25% of the average strain, indicating no overstraining. The significant differences in the strains at the left and right sides could be attributed to the asymmetry of the wound geometry. The lowest strains were observed at the top and bottom of the wound similar to other wounds. Additionally, no possible rupture sites were identified. Figure 11. Average strain (y-axis) versus time frame (x-axis) measured using DIC from tensile test on the skin phantom specimen 2 with wound.

Strains in Skin Phantom with Asymmetric-Wide Wound
The test specimen with asymmetric and wide wound geometry (Wound 3) was subjected to tensile loads for 42 s and 14 frames were extracted for the DIC measurements. A high stress build-up was observed initially on the left-hand side of the wound which shifted eventually to the right-hand side (see Figure 12). On quantification of the strains along the circumference of the long wound with 7 markers (Figure 13), the strain fields were found to be close to each other all around the wound with the peak at point P2 initially which was on the left-hand side. After 8 s, the maximum strain was recorded for marker P5, which was within 25% of the average strain, indicating no overstraining. The significant differences in the strains at the left and right sides could be attributed to the asymmetry of the wound geometry. The lowest strains were observed at the top and bottom of the wound similar to other wounds. Additionally, no possible rupture sites were identified.

Strains in Skin Phantom with Asymmetric-Thin Wound
Specimen with asymmetric and thin wound geometry (Wound 4) was subjected to the tensile test for 41 s and 14 frames were captured ( Figure 14). Eight markers were placed around the wound circumference and the strains were tracked (Figure 15) using DIC. High strain build-ups were observed at the top (P0 and P1), left side (P3) and bottom right (P5) locations along the wound geometry. The reason behind such arbitrary high strain zones could be the thin and asymmetric geometry of the wound. The maximum strain quantified for marker P5 was approximately 24% above the average strain, indicating possibility of overstraining at higher applied strains. The lowest strain was observed at the bottom most point of the wound. Strains estimated at the other locations were observed to be similar. Also, with time, the top locations P0 and P1 experienced high strains (with maximum values over 22% of average strain) surrounded by a low strain zone (with minimum values below 26% of average strain), which may be possible sites of rupture with further applied strains, cause progressive tissue tearing.

Strains in Skin Phantom with Asymmetric-Thin Wound
Specimen with asymmetric and thin wound geometry (Wound 4) was subjected to the tensile test for 41 s and 14 frames were captured (Figure 14). Eight markers were placed around the wound circumference and the strains were tracked (Figure 15) using DIC. High strain build-ups were observed at the top (P0 and P1), left side (P3) and bottom right (P5) locations along the wound geometry. The reason behind such arbitrary high strain zones could be the thin and asymmetric geometry of the wound. The maximum strain quantified for marker P5 was approximately 24% above the average strain, indicating possibility of overstraining at higher applied strains. The lowest strain was observed at the bottom most point of the wound. Strains estimated at the other locations were observed to be similar. Also, with time, the top locations P0 and P1 experienced high strains (with maximum values over 22% of average strain) surrounded by a low strain zone (with minimum values below 26% of average strain), which may be possible sites of rupture with further applied strains, cause progressive tissue tearing.

Strains in Skin Phantom with Transverse Wound
The specimen with the transverse wound geometry (Wound 5) was subjected to tensile test for 57 s and 12 frames were extracted for DIC measurements. It was observed that there was a symmetric strain build-up on the left and right side of the wound to start with, followed by a localized high strain build-up on the left size of the wound circumference, finally leading to tearing in that direction (see frame 12 of the Figure 16). Four markers were placed along the wound geometry to quantify the strain changes ( Figure 17). Point P1 at the left-hand corner had a steep rise in strain followed by the point P3 on the right side, P0 on the top, and P2 at the bottom. The maximum value of strain recorded for marker P1 was over 53% higher than the average strain, indicating massive overstraining. The average strain was found to be higher than the strains at the bottom point P2. Also, a low strain region (with minimum value <45% of average strain) was observed around the localized high strain zone at the location of marker P1, indicating high possibility of tissue rupture. The actual tearing observed in the experiment validated this hypothesis of high tearing possibility at the location of high strain surrounded by a low strain zone. Figure 15. Average strain (y-axis) versus time frame (x-axis) measured using DIC from tensile test on the skin phantom specimen 4 with wound.

Strains in Skin Phantom with Transverse Wound
The specimen with the transverse wound geometry (Wound 5) was subjected to tensile test for 57 s and 12 frames were extracted for DIC measurements. It was observed that there was a symmetric strain build-up on the left and right side of the wound to start with, followed by a localized high strain build-up on the left size of the wound circumference, finally leading to tearing in that direction (see frame 12 of the Figure 16). Four markers were placed along the wound geometry to quantify the strain changes ( Figure 17). Point P1 at the left-hand corner had a steep rise in strain followed by the point P3 on the right side, P0 on the top, and P2 at the bottom. The maximum value of strain recorded for marker P1 was over 53% higher than the average strain, indicating massive overstraining. The average strain was found to be higher than the strains at the bottom point P2. Also, a low strain region (with minimum value <45% of average strain) was observed around the localized high strain zone at the location of marker P1, indicating high possibility of tissue rupture. The actual tearing observed in the experiment validated this hypothesis of high tearing possibility at the location of high strain surrounded by a low strain zone.

Comparison of Strains in Skin Phantoms across Wound Geometries
The DIC results from the different wound geometries were compared. It was observed that the strain at the bottom most point of the wound was always the lowest and its value was lower than that of the average strains measured for the specimen throughout. For symmetric long wounds with the long edge along the y-axis (direction of tension), the maximum strains were observed at the left and right sides of the wound with similar values. The wider this long wound was, the greater was the difference between the strains on the left and right side. Asymmetry in the wound geometry also played a major role in the strain differences in the left and right side. For a short wound, the strain distribution was similar at the top, left and right due to the comparable aspect ratio (AR~1) of the wound length (along y-axis) and width (along x-axis). However, the geometry of the wound in case of short wounds was observed to play a significant role in deciding the superior strain build-up locations with tension, which may drive soft tissue rupture (as seen in case with wound 5). Also, for thin and long wounds (wound geometry 4), high strain build-up was observed on the top of the wound which may lead to rupture and increase in the length of the wound. All these findings provide valuable insights into the local wounded skin biomechanics for common wound geometries and sizes.
The average strain versus time plots of the specimens with no wound and the ones with the different wound geometries were compared as shown in Figure 18. From the study, the strains observed for wounded specimens were found to be lower than those without wounds. The lowest strains were observed for the specimen with the smallest size and the most symmetric wound (specimen 5). The strains were found to increase with decreasing AR. The strains for the wound with lowest aspect ratios (Wound 4, AR = 0.07) was found to be the closest to the strains observed in the specimens with no wounds.

Comparison of Strains in Skin Phantoms across Wound Geometries
The DIC results from the different wound geometries were compared. It was observed that the strain at the bottom most point of the wound was always the lowest and its value was lower than that of the average strains measured for the specimen throughout. For symmetric long wounds with the long edge along the y-axis (direction of tension), the maximum strains were observed at the left and right sides of the wound with similar values. The wider this long wound was, the greater was the difference between the strains on the left and right side. Asymmetry in the wound geometry also played a major role in the strain differences in the left and right side. For a short wound, the strain distribution was similar at the top, left and right due to the comparable aspect ratio (AR~1) of the wound length (along y-axis) and width (along x-axis). However, the geometry of the wound in case of short wounds was observed to play a significant role in deciding the superior strain build-up locations with tension, which may drive soft tissue rupture (as seen in case with wound 5). Also, for thin and long wounds (wound geometry 4), high strain build-up was observed on the top of the wound which may lead to rupture and increase in the length of the wound. All these findings provide valuable insights into the local wounded skin biomechanics for common wound geometries and sizes.
The average strain versus time plots of the specimens with no wound and the ones with the different wound geometries were compared as shown in Figure 18. From the study, the strains observed for wounded specimens were found to be lower than those without wounds. The lowest strains were observed for the specimen with the smallest size and the most symmetric wound (specimen 5). The strains were found to increase with decreasing AR. The strains for the wound with lowest aspect ratios (Wound 4, AR = 0.07) was found to be the closest to the strains observed in the specimens with no wounds.

Computational Modeling of Stress and Strain on Transverse Wound
A computational model of the biofidelic skin phantom with transverse wound geometry 5 was developed to model the local strains in tension and to compare it with the experimental results and

Computational Modeling of Stress and Strain on Transverse Wound
A computational model of the biofidelic skin phantom with transverse wound geometry 5 was developed to model the local strains in tension and to compare it with the experimental results and DIC measurements. Applying constant strain rates employed in the experiments to the computational model, the deformation of the specimen is shown in Figure 19 at the times 2 s, 24 s and 57 s respectively from the start of the test. and 57 s respectively from the start of the test.
The local strains around the wound studied using the computational analysis were compared with DIC measurements. Figure 20 captures the strain distribution around the wound for a constant strain rate of 0.4 mm/s at 2 s, 24 s and 57 s respectively. The maximum strain was observed at the left side of the wound in both the computational and DIC measurements. This localized strain build-up could be attributed to the geometry of the wound. Another common observation is the occurrence of minimum strains at the top and bottom of the wound. After 24 s, strain reduction was observed around the wound due to tearing in the DIC measurements, which was not captured by the computational model.   The local strains around the wound studied using the computational analysis were compared with DIC measurements. Figure 20 captures the strain distribution around the wound for a constant strain rate of 0.4 mm/s at 2 s, 24 s and 57 s respectively. The maximum strain was observed at the left side of the wound in both the computational and DIC measurements. This localized strain build-up could be attributed to the geometry of the wound. Another common observation is the occurrence of minimum strains at the top and bottom of the wound. After 24 s, strain reduction was observed around the wound due to tearing in the DIC measurements, which was not captured by the computational model. DIC measurements. Applying constant strain rates employed in the experiments to the computational model, the deformation of the specimen is shown in Figure 19 at the times 2 s, 24 s and 57 s respectively from the start of the test. The local strains around the wound studied using the computational analysis were compared with DIC measurements. Figure 20 captures the strain distribution around the wound for a constant strain rate of 0.4 mm/s at 2 s, 24 s and 57 s respectively. The maximum strain was observed at the left side of the wound in both the computational and DIC measurements. This localized strain build-up could be attributed to the geometry of the wound. Another common observation is the occurrence of minimum strains at the top and bottom of the wound. After 24 s, strain reduction was observed around the wound due to tearing in the DIC measurements, which was not captured by the computational model.   The maximum local strain variation with time around the wound, from the computational model, was compared with the experimental data obtained using DIC measurements ( Figure 21). It can be seen that there is good agreement between the computational data and the experimentally measured maximum strain. The maximum local strain variation with time around the wound, from the computational model, was compared with the experimental data obtained using DIC measurements ( Figure 21). It can be seen that there is good agreement between the computational data and the experimentally measured maximum strain. The average stress-strain responses of the computational and experimental models were compared for the transverse wound 5 geometry selected, and also with the normal unwounded skin. In Figure 22, the skin with the wound was observed to be more compliant than the normal skin. Comparing the experimental and computational plots of the wound under tension, an R 2 correlation index of 0.976 was estimated. The computational plot was observed to become stiffer compared to the experimental values after a strain of about 0.4. This observation may be attributed to skin phantom tearing which was observed during the experiments. However, the tearing was not simulated in the computational model, possibly resulting in higher stress values after a strain of 0.4.  The average stress-strain responses of the computational and experimental models were compared for the transverse wound 5 geometry selected, and also with the normal unwounded skin. In Figure 22, the skin with the wound was observed to be more compliant than the normal skin. Comparing the experimental and computational plots of the wound under tension, an R 2 correlation index of 0.976 was estimated. The computational plot was observed to become stiffer compared to the experimental values after a strain of about 0.4. This observation may be attributed to skin phantom tearing which was observed during the experiments. However, the tearing was not simulated in the computational model, possibly resulting in higher stress values after a strain of 0.4. The maximum local strain variation with time around the wound, from the computational model, was compared with the experimental data obtained using DIC measurements ( Figure 21). It can be seen that there is good agreement between the computational data and the experimentally measured maximum strain. The average stress-strain responses of the computational and experimental models were compared for the transverse wound 5 geometry selected, and also with the normal unwounded skin. In Figure 22, the skin with the wound was observed to be more compliant than the normal skin. Comparing the experimental and computational plots of the wound under tension, an R 2 correlation index of 0.976 was estimated. The computational plot was observed to become stiffer compared to the experimental values after a strain of about 0.4. This observation may be attributed to skin phantom tearing which was observed during the experiments. However, the tearing was not simulated in the computational model, possibly resulting in higher stress values after a strain of 0.4.

Conclusions
In this work, the mechanics of wounded skin was studied using experiments, DIC and computational modeling. Specifically, localized overstraining and possibilities of rupture in wounded tissue was investigated. A biofidelic human skin phantom was employed to develop skin tensile test specimens, and some common wounds geometries were created. The test specimens were subjected to uniaxial strains and DIC was used to investigate the strain fields developed due to different wound geometries. Local markers were placed around the wound circumference to closely monitor the strain changes around the wound during the tensile test. The symmetric wounds were found to result in maximum strains at the left and right sides and minimum strains at the top and bottom. Asymmetry and reducing wound thickness was observed to lead to significant differences in the strain values on the left and right sides of the wound, and also the development of localized high strains. Also, the wound aspect ratio was found to inversely affect the average strains developed around the wound periphery. Skin overstraining and possibilities of tearing were observed in all symmetric wounds, with the maximum effect in the transverse wound with the lowest aspect ratio. A novel computational model was developed based on the experiments with the common elliptical transverse wound geometry, exhibiting the most prominent overstraining effect and rupture. Hyperelastic material properties were adopted from the experiments and tensile strain was applied on the numerical model of the wounded skin. The results from the analyses were found to be in close agreement with experiments and DIC measurements. A few limitations of the computational work should be acknowledged. Only a few common elliptical wound geometries were considered in the study. Skin tissue tearing was not modeled in the current computational model which will be incorporated in future studies to improve the accuracy of the model. Also, the effect of biaxial pre-stress loading and skin tissue material anisotropy needs to be included to obtain more precise insights into the mechanics of wounded skin. To date, this study is the first attempt towards biomechanical modeling of wounded skin, with a focus on localized overstraining and rupture possibilities. The results of the study will help the selection of wound closure techniques which may cause minimal skin overstraining and rupture possibilities, informing pre-surgical planning and future applications in robotic surgeries.