Ruike Shi

Stéphane Avril

Haitian Yang

Víctor A Acosta Santamaría

Yue Mei

Yiqian He
email: heyiqian@dlut.edu.cn

Adaptation of the virtual fields method for the identification of biphasic hyperelastic model parameters in soft biological tissues with osmotic swelling

Keywords: Biphasic hyperelasticity, parameter identification, chemoelasticity, virtual fields method, osmotic effects

Biphasic hyperelastic models have become popular for soft hydrated tissues and there is a pressing need for appropriate identification methods using full-field measurement techniques such as digital volume correlation. This paper proposes to address this need with the virtual fields method (VFM). The main asset of the proposed approach is that it avoids the repeated resolution of complex nonlinear finite-element models. By choosing special virtual fields, the VFM approach can extract hyperelastic parameters of the solid part of the biphasic medium without resorting to identifying the model parameters driving the osmotic effects in the interstitial fluid. The proposed approach is verified and validated through three different examples: first using simulated data and then using experimental data obtained from porcine descending thoracic aortas samples in osmotically active solution.

Introduction

The mechanical behavior of soft biological tissues, such as cartilage or arteries, can be commonly modeled as biphasic, accounting for the intertwined contributions of fiber networks and interstitial fluid [START_REF] Lanir | Osmotic swelling and residual stress in cardiovascular tissues[END_REF] . This biphasic behavior implies chemomechanical couplings regulating the osmotic pressure of the interstitial fluid [START_REF] Ehret | Inverse poroelasticity as a fundamental mechanism in biomechanics and mechanobiology[END_REF] . An example of such coupling is the swelling effect that occurs when there is a gradient of osmotic pressure, which can be attributed to the presence of proteoglycans [START_REF] Kovach | The importance of polysaccharide configurational entropy in determining the osmotic swelling pressure of concentrated proteoglycan solution and the bulk compressive modulus of articular cartilage[END_REF][START_REF] Lanir | Biorheology and fluid flux in swelling tissues. I. Bicomponent theory for small deformations, including concentration effects[END_REF] within collagenous fiber networks in the tissues [START_REF] Ehret | Inverse poroelasticity as a fundamental mechanism in biomechanics and mechanobiology[END_REF] .

The biphasic hyperelastic constitutive models of hydrated biological tissues were originally introduced by Eringen and Ingram in the 1960s [START_REF] Eringen | A continuum theory of chemically reacting media-I[END_REF][START_REF] Ingram | A continuum theory of chemically reacting media-II Constitutive equations of reacting fluid mixtures[END_REF] . These biphasic models have been widely used in the biomechanics research of soft biological tissues. For instance, -In the field of the biomechanics of cartilage, Mow et al. [START_REF] Mow | Biphasic creep and stress relaxation of articular cartilage in compression: theory and experiments[END_REF] and Lai et al. [START_REF] Lai | A triphasic theory for the swelling and deformation behaviors of articular cartilage[END_REF] established a biphasic model and a triphasic model. These models include mainly the collagenproteoglycan matrix and the interstitial fluid, respectively. Then Ateshian et al. [START_REF] Ateshian | The correspondence between equilibrium biphasic and triphasic material properties in mixture models of articular cartilage[END_REF][START_REF] Ateshian | Equivalence between short-time biphasic and incompressible elastic material responses[END_REF][START_REF] Ateshian | Modeling the matrix of articular cartilage using a continuous fiber angular distribution predicts many observed phenomena[END_REF] applied for the models to articular cartilage with osmotic loading. Besides, Wilson et al. [START_REF] Wilson | A comparison between mechanoelectrochemical and biphasic swelling theories for soft hydrated tissues[END_REF] studied the links between the biphasic swelling theories and mechano-electrochemical models. Wayne et al. [START_REF] Wayne | Application of the u-p finite element method to the study of articular cartilage[END_REF] developed a u-p finite element model for biphasic theory by coupling deformation and fluid pressure. Recently, Kandil et al. [START_REF] Kandil | Interlamellar-induced timedependent response of intervertebral disc annulus: A microstructure-based chemoviscoelastic model[END_REF] proposed a microstructure-based chemo-viscoelastic model considering the osmo-induced deformation and internal fluid variation.

-In the field of biomechanics of arteries, Lanir [START_REF] Lanir | Osmotic swelling and residual stress in cardiovascular tissues[END_REF][START_REF] Lanir | Biorheology and fluid flux in swelling tissues. I. Bicomponent theory for small deformations, including concentration effects[END_REF] developed the biocomponent theory in 1987. He applied this theory to residual stresses of cardiovascular tissues. Then Azeloglu et al. [START_REF] Azeloglu | Heterogeneous transmural proteoglycan distribution provides a mechanism for regulating residual stresses in the aorta[END_REF] proposed a numerical model combining the triphasic theory of Lai et al. [START_REF] Lai | A triphasic theory for the swelling and deformation behaviors of articular cartilage[END_REF] with the FEM to study the influence of proteoglycans on residual stresses in the aorta. Their model was eventually validated experimentally. Recently Santamaría et al. [START_REF] Acosta Santamaría | Three-dimensional full-field strain measurements across a whole porcine aorta subjected to tensile loading using optical coherence tomography-digital volume correlation[END_REF][START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] developed experimental methods for quantifying the chemoelastic effects in arterial tissues.

All those studies showed the importance of biphasic or triphasic theory to investigate the biomechanics of hydrated biological tissues, such as cartilage and arteries.

The identification of constitutive model parameters for soft biological tissues is still flourishing [START_REF]Material Parameter Identification and Inverse Problems in Soft Tissue Biomechanics[END_REF][START_REF] Avril | Overview of identification methods of mechanical parameters based on full-field measurements[END_REF] . Among the available identification approaches, the virtual fields method (VFM) has been widely used in solid mechanics [START_REF] Pierron | The Virtual Fields Method[END_REF][START_REF] Grédiac | Principe des travaux virtuels et identification. (Principle of virtual work and identification)[END_REF][START_REF] Grédiac | Special virtual fields for the direct determination of material parameters with the virtual fields method. 1--Principle and definition[END_REF][START_REF] Grédiac | Special virtual fields for the direct determination of material parameters with the virtual fields method. 2--Application to inplane properties[END_REF][START_REF] Grédiac | Special virtual fields for the direct determination of material parameters with the virtual fields method. 3. Application to the bending rigidities of anisotropic plates[END_REF][START_REF] Mei | General finite-element framework of the virtual fields method in nonlinear elasticity[END_REF][START_REF] Mei | Introducing regularization into the virtual fields method (VFM) to identify nonhomogeneous elastic property distributions[END_REF] , given its advantages, such as its insensitivity to the uncertainty of boundary conditions [START_REF] Grédiac | The virtual fields method for extracting constitutive parameters from full-field measurements: a review[END_REF] , robustness [START_REF] Avril | Sensitivity of the virtual fields method to noisy data[END_REF] , and fast convergence [START_REF] Avril | General framework for the identification of constitutive parameters from full-field measurements in linear elasticity[END_REF] .

Avril et al. [START_REF] Avril | Anisotropic and hyperelastic identification of in vitro human arteries from full-field optical measurements[END_REF] first applied VFM to arterial tissues to identify anisotropic hyperelastic material parameters. Kim et al. [START_REF] Kim | Experimental characterization of rupture in human aortic aneurysms using a full-field measurement technique[END_REF] used VFM to identify the material properties of human aneurysmal aortas. Bersi et al. [START_REF] Bersi | Novel methodology for characterizing regional variations in the material properties of murine aortas[END_REF] collected full-field biaxial data and used the VFM to estimate regional variations in material parameters for a microstructurally motivated constitutive model.

The identification of biphasic hyperelastic model parameters has been rather challenging.

Nadeen et al. [START_REF] Chahine | Direct measurement of osmotic pressure of glycosaminoglycan solutions by membrane osmometry at room temperature[END_REF] designed an apparatus for measuring osmotic pressure and estimated how it contributes to cartilage stiffness. More recently, Santamaría et al. [START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] measured 3D strain fields induced by chemoelastic effects in arteries. These studies showed the significance of the biphasic constitutive model in arteries. Nevertheless, to the best of the authors' knowledge, there has been no existing work dedicated to the identification of biphasic constitutive models using VFM. Besides, although soft tissues (e.g. arterial wall) could be multi-layered structures, there is no work on the identification of layer-specific parameters for soft tissues in the framework of biphasic hyperelastic model. The objective of this work is to address this lack. Considering that biphasic constitutive models are computationally challenging due to coupling between elasticity and osmotic effects, we designed a novel specific VFM framework for chemoelasticity.

The paper is organized as follows: Section 2 introduces the materials and methods, including constitutive relationship for biphasic hyperelasticity, general principle for VFM and the application of VFM for biphasic hyperelasticity. Section 3 reports the results of the verification of the proposed model and its application into practical problems, first conducted with simulated data and then with experimental data. Finally, the conclusion and discussion are given in Section 4.

Materials and methods

The introduction for biphasic hyperelasticity in arterial tissues

The framework of the biphasic constitutive relationship used in this paper is proposed by Azeloglu et al. [START_REF] Azeloglu | Heterogeneous transmural proteoglycan distribution provides a mechanism for regulating residual stresses in the aorta[END_REF] , in which the residual stresses are hypothesized to arise from the fixedcharge density of the proteoglycans present in the arterial tissues, inducing the Donnan osmotic pressure relative to the external environment as only equilibrium conditions are considered for simplicity, when the solid and fluid velocities have reduced to zero. Therefore, the Cauchy stress T is assumed to have two components:

T = T i + T s (1)
where T i is the interstitial fluid component and T s is the solid matrix component.

As illustrated in Fig. 1, the total deformation gradient F is realized by the multiplying a prestretched part F pre and a mechanical part F m [START_REF] Kandil | Interlamellar-induced timedependent response of intervertebral disc annulus: A microstructure-based chemoviscoelastic model[END_REF] . The F pre is a stress-free swelling due to osmotic effects without any mechanical loading. Therefore, the total deformation gradient

F = F m F pre . ( 2
)
The interstitial fluid component of Cauchy stress T i takes the form based on the Donnan equilibrium theory [START_REF] Azeloglu | Heterogeneous transmural proteoglycan distribution provides a mechanism for regulating residual stresses in the aorta[END_REF][START_REF] Donnan | The theory of membrane equilibria[END_REF][START_REF] Overbeek | The donnan equilibrium[END_REF] , that is,

T i = -(p -p * )I = -Rθ (c F ) 2 + (c * ) 2 -c * I ( 3
)
where I is the identity tensor, p * is the ambient pressure in the external bath (henceforth considered to be zero, p * = 0), c * is the external bath salt osmolarity, R is the universal gas constant, and θ is absolute temperature, c F is the proteoglycan fixed-charge density relating to the solid matrix relative volume change (J = det F ) via

c F = φ w 0 c F 0 J -1 + φ w 0 ( 4
)
where c F 0 and φ w 0 are the fixed-charge density and water content, respectively. Note that in this paper the above chemical parameters are assumed to be homogeneous.

For the solid matrix component, consider an isotropic neo-Hookean constitutive relationship [START_REF] Bower | Applied Mechanics of Solids[END_REF] :

T s = µJ -5/3 B - 1 3 I 1 I + κ(J -1)I (5)
where B = F F T is the left Cauchy-Green deformation tensor, I 1 is the first invariant of the right Cauchy-Green deformation tensor C, and µ and κ are the shear modulus and bulk modulus of the solid matrix, respectively.

In this paper, the reference configuration is set as the zero osmotic loading stage, and a conceptual sequence of configurations is assumed as shown in Fig. 1 according to Kandil et al. [START_REF] Kandil | Interlamellar-induced timedependent response of intervertebral disc annulus: A microstructure-based chemoviscoelastic model[END_REF] , in which the stress-free chemical configuration is an intermediate virtual state due to osmotic swelling. In this way, the total deformation is divided into chemical-induced and mechanical-induced parts. Note that a non-zero prestretched deformation F pre is produced to maintain the Cauchy stress T = 0 in the stress-free chemical configuration. By combining the Eq.(3), Eq.( 4) and Eq.( 5), and note that the prestretch F pre has the form

F pre = λ pre I = (J pre ) 1/3 I ( 6
)
where J pre = det F pre , we have the equation for solving F pre as

-p + κ(J pre -1) = 0 . ( 7
)
Finally, the prestretch can be calculated as

F pre = 1 + p κ 1/3 I . ( 8
)
4

General principle for VFM

The virtual fields method (VFM) is based on the principle of virtual work, for quasi-static conditions with neglecting the acceleration and the body force [START_REF] Pierron | The Virtual Fields Method[END_REF] , which may be written such as,

- Ω T : ε * dΩ + ∂Ωt t • u * d∂Ω t = 0 ( 9
)
where T is the actual Cauchy stress tensor across the domain of interest Ω and is related to the strains through constitutive equations, strains are derived from gradients of the measured displacement field u, t are the tractions applied on a part of the boundary ∂Ω t , u * is a kinematically admissible virtual displacement field, and ε * is the associated virtual strain field.

First, a constitutive model is chosen, one can write in the general case

T = G (ε), ( 10
)
where G is a given function of the actual strain components. The constitutive parameters are also involved in G . Therefore Eq.( 9) can be written as

- Ω G (ε) : ε * dΩ + ∂Ωt t • u * d∂Ω t = 0 . ( 11
)
Then, virtual fields are chosen. Each virtual field introduced into equation Eq.( 26) yields one scalar equation. The constitutive parameters are then sought as solutions to a set of such equations. The construction of the virtual fields is a key issue of the method. An important feature is the fact that the above equation is verified for any kinematically admissible virtual field u * [START_REF] Pierron | The Virtual Fields Method[END_REF] . Kinematically admissible means that u * must be continuous across the whole volume and it must be equal to the prescribed displacement on the boundary. The virtual field u * can be constructed analytically or solving equations created automatically [START_REF] Mei | On improving the accuracy of nonhomogeneous shear modulus identification in incompressible elasticity using the virtual fields method[END_REF] . The detailed procedure of constructing virtual field for the proposed model in this paper will be introduced later.

The adaptation of the VFM to biphasic hyperelasticity 2.3.1. Governing equations

Assuming quasi-static conditions and no body forces [START_REF] Avril | Anisotropic and hyperelastic identification of in vitro human arteries from full-field optical measurements[END_REF] , the VFM equations for the biphasic hyperelasticity can be written as

- Ω (T i + T s ) : ε * dΩ + ∂Ωt t • u * d∂Ω t = 0 (12)
where T s and T i denotes the solid matrix component and interstitial fluid component for actual Cauchy stress across the domain of interest, respectively, t are the tractions applied on the boundary ∂Ω t , Ω is the domain volume in the current configuration, and ε * is the virtual strain field.

Definition of virtual strain fields

In this paper, the virtual displacement fields are defined by solving a serial of equations.

Firstly it is assumed that the osmotic pressure terms vanish under the principle of virtual work with the use of some special virtual fields, i.e. a special virtual field with tr(ε * ) = 0,

yielding Ω -p(J)I : ε * dΩ = 0 . ( 13
)
Further if we consider a three-layers arterial wall tissue for example, including intima, media and adventitia. Then, Eq.( 12) can be written as

- 3 i=1 Ω i T (µ i ) : ε * dΩ + ∂Ωt t • u * d∂Ω t = 0 ( 14
)
where Ω i is the domain of the i-th layer, and µ i is the shear modulus in Ω i .

In order to construct the virtual strain fields ε * in Eq.( 13)a simple approach is to construct a linear virtual displacement field to produce a constant virtual strain field,

u * = (η 1 x)e x + (η 2 y)e y + (η 3 z)e z . ( 15
)
where it is required to make η 1 + η 2 + η 3 = 0 for meeting Eq.( 13).

For the case of identification for multi-layers in Eq.( 14), the virtual fields are required to be linearly independent when there are many unknown parameters to be identified simultaneously, and in order to realize tr(ε * ) = 0, these more complicate virtual fields can be constructed by defining the curl of a potential vector field H [START_REF] Mei | On improving the accuracy of nonhomogeneous shear modulus identification in incompressible elasticity using the virtual fields method[END_REF] u * = ∇ × H , (16) since the divergence of the curl above always yields

div(∇ × H) = 0 , ( 17
)
where H = H i e i (i = x, y, z) and e i is the unit vector in the global cartesian coordinate system.

Thus, we can reconstruct virtual displacement fields such as

u * = ∇ × H = ∂H z ∂y - ∂H y ∂z e x + ∂H x ∂z - ∂H z ∂x e y + ∂H y ∂x - ∂H x ∂y e z . ( 18
)
In this work, under the kinematically admissible condition, we calculated the virtual displacement field by assuming a special case in an uniaxial tension for a plate as shown in Fig. 2, where displacement boundary condition is

u x | x=0 = u y | y=0 = u z | z=0 = 0 which is
usually used, and define the virtual fields as follow:

u * x = ∂H z ∂y - ∂H y ∂z = ηx k f (z) (19)
u * y = ∂H x ∂z - ∂H z ∂x = 0 (20) u * z = ∂H y ∂x - ∂H x ∂y = ηz k+1 g(x) , ( 21
)
where k is the exponent parameter of the power functions, η is a scaling factor parameter which is usually set to a small value to prevent excessive virtual deformation, f (z) and g(x)

are undetermined functions. Note that power functions are used for the virtual displacements in Eq.( 19) and ( 21). The purpose of power functions is to generate virtual strain fields that avoid linear dependent VFM equations.

In order to simply the problem, we set H y = 0, in this way H x and H z can be directly calculated by integration leading to

H x = -ηy • z k+1 g(x) ( 22
)
H z = ηy • x k f (z) . ( 23
)
Then, by substituting Eqs.( 22) and ( 23) into Eq.( 20), we can deduce that

g(x) = -kx (k-1) (24)
f (z) = (k + 1)z k . ( 25
)
Finally, combining Eqs.( 19), ( 20), ( 21), ( 24) and ( 25), we obtain the general equations for these special virtual fields such as

     u * x = η(k + 1)x k z k u * y = 0 u * z = -ηkx k-1 z k+1 . ( 26
)

Summary of the complete procedure

Based on the theory and methods introduced above, the procedure for the implementation of the proposed model is briefly shown in Table 1, mainly including three parts in the following:

(1)Full-field measurements: The displacement fields used in this study have three different possible origins, namely analytical solutions, FE analyses and real experiments. For simple homogeneous (monolayer) uniaxial tension, the analytical solution can be used. For a multilayered structure, the FEM is used to provide simulated measurements and derive the deformation gradient F and the Jacobian J. For the experimental case, the displacements are measured using optical coherence tomography (OCT) and digital volume correlation (DVC) before deducing F , E and J, The traction t for the VFM in numerical verifications is pointwise on the boundary, but as for experiment case, the measured stress is an average force by assuming cross-section perfectly homogeneous as it is hard to measure the heterogeneous traction distribution on the boundary in the experiment.

(2)Building the VFM equations: From the obtained F fields, we determine T according to Eq.( 1),( 3),( 4), (5). Special virtual fields based on Eq.( 26) are applied to Eq.( 14).

Since T depends on the unknown parameters µ, the virtual work equation Eq.( 14) is not satisfied exactly, and a residual can be derived, as defined in Eq.( 27). This residual is then minimized. The cost function can be written such as

f (µ 1 , µ 2 , µ 3 ) = 3 v=1 - 3 i=1 Ω i T : ε * v dΩ + ∂Ω t • u * v d∂Ω 2 , ( 27
)
which represents the quadratic deviation between the internal virtual work (IVW) and the external virtual work (EVW), with v labelling the virtual field and i labelling the layer.

Using Gauss quadrature in isoparametric elements, integrals in the cost function can be changed into discrete sums, such as

f (µ 1 , µ 2 , µ 3 ) = 3 v=1 - 3 i=1 ele T (µ i , X ele ) : ε * v (X ele ) + eleb t eleb A eleb u * v (X eleb ) 2 (28)
where ele and eleb denote the element and boundary element respectively, A eleb is the area of one element side, ε * v (X ele ) denotes the value of the v-th virtual strain field at the Gauss integration point in the ele-th element, and ε * v is derived from u * v .

(3) Solving the VFM equations: We use the MATLAB f minsearch function to search for µ, in which the derivative-free Nelder-Mead simplex optimization algorithm is used for finding the minimum value of an unconstrained multivariable scalar function [START_REF] Lagarias | Convergence Properties of the Nelder-Mead Simplex Method in Low Dimensions[END_REF] .

The problem can be expressed as

     Find µ = {µ 1 , µ 2 , µ 3 } f (µ) = 3 v=1 -W In v + W Ex v 2 min f (µ) (29)
where W In v represents the IVW under the v-th virtual field, and W Ex v represents the EVW under the v-th virtual field. Eventually, the cost function is minimized through an iterative scheme and yields the unknown material parameters. 5. Choice of virtual field that satisfies tr(ε * ) = 0 to avoid the effect of the chemo term. 6. Computation of the value of cost function f . 7. Use of the optimization method to find a vector of µ by minimizing f . 8. Checking whether the identified value equals the parameters used in forward simulation.

Noise investigation method

In order to study the influence of measurement errors on identifications, we add artificial white noise to the displacement data obtained by the FEM simulations as [START_REF] Thakolkaran | NN-EUCLID: Deep-learning hyperelasticity without stress data[END_REF] .

u a,n i = u a,fem i + δ a,∆n i with δ a,∆n i ∼ N (0, σ 2 u ) ∀(a, i) ∈ Ω (30)
where u a,fem i and δ a,∆n i denote the true displacement from FEM and random noise, respectively, the δ a,∆n i represents the Gaussian white noise with zero mean and standard deviation σ u ranging between 10 -4 (low noise) to 10 -3 (high noise) [START_REF] Flaschel | Unsupervised discovery of interpretable hyperelastic constitutive laws[END_REF] . For each fixed value of σ u , 40 realizations of Gaussian white noise δ a,∆n i were generated. The average identified results for all the samples were calculated [START_REF] Yiqian | Solving inverse couple-stress problems via an element-free Galerkin (EFG) method and Gauss-Newton algorithm[END_REF] .

The relative error was defined as

error = µ i -μi µ i × 100% (31)
where µ i and μi are the exact shear modulus and average identified value under 40 samples respectively.

Verification of the new VFM based on analytical or FEM simulations

In this paper, the proposed model as established in previous section will be verified in theoretical cases before applying into the practical problems. In the verification, the fullfiled measurements are assumed to be 'virtually' obtained from analytical solutions or FE analyses. The derivations for analytical solutions and FE schemes are introduced in the following, respectively.

Analytical solution for a plate under tension

For some cases with simple geometry and boundary conditions, the analytical solution of the biphasic hyperelastic problem can be derived. For example, consider a plate under uniaxial tension, as shown in Fig. 2, where λ is the prescribed axial stretch, and α and β are

F =   λ 0 0 0 α 0 0 0 β   . ( 32
)
T 22 = -p + µ(λαβ) -5/3 ( 2 3 α 2 - 1 3 λ 2 - 1 3 β 2 ) + κ(λαβ -1) = 0 (33)
T 33 = -p + µ(λαβ) -5/3 ( 2 3 β 2 - 1 3 λ 2 - 1 3 α 2 ) + κ(λαβ -1) = 0 . ( 34
)
The relationship α = β can be derived from the equations above, therefore the Eq.( 34) can be rewritten as

T 33 = -p + µ(λβ 2 ) -5/3 1 3 β 2 - 1 3 λ 2 + κ(λβ 2 -1) = 0 (35)
where p based on Eq.( 3) and ( 4), can be elaborated as

p =   Rθ   φ w 0 c F 0 λβ 2 -1 + φ w 0 2 + (c * ) 2 -c *   + p *   . ( 36
)
To solve this nonlinear problem, the Newton method is applied, yielding

β n+1 = β n - T 33 (β n ) T ′ 33 (β n ) (37)
where β n is the n-th iteration of β. For any λ, β and α, hence F can be deduced with this iterative resolution.

FEM solutions

The FEM is an alternative approach when analytical resolution is not possible. Let x = φ(X, t) describe the motion of a body from the initial reference configuration X ∈ k(0)

to its current configuration x ∈ k(t). Define the deformation gradient as

F = ∇ X φ . ( 38
)
The general total potential energy for the biphasic hyperelasticity model is

Π = Ω Ψ dΩ - Ω b • φdΩ - ∂Ω t • φd∂Ω . ( 39
)
where Ψ is the strain energy, b is the body force vector per unit current volume of Ω and t are the tractions on the boundary ∂Ω.

Minimization of the potential energy with respect to φ results in a nonlinear system of equation, and the variational forms in the spatial description are

Ω T : ∇ x δφdΩ - Ω b • δφdΩ - ∂Ω t • δφd∂Ω = 0 . ( 40
)
In finite element discretization, isoparametric interpolations of geometry variables X and field variables φ are written as

X = n I=1 N I (ξ)X I , φ = n I=1 N I (ξ)φ I ( 41
)
where ξ denotes the coordinates in the reference element, and n is the number of nodes in each element.

Substituting Eq.( 41) into Eq.( 40) and using the Newton-Raphson scheme to solve the nonlinear system of equations [START_REF] Waffenschmidt | A gradient-enhanced largedeformation continuum damage model for fibre-reinforced materials[END_REF] , we obtain

R φ i + K φ i • ∆φ i+1 = 0 ( 42
)
where R φ is the global residual vector, and the element stiffness matrix is given by

K φ = Ω ∇ T x N • C • ∇ x N dΩ + Ω [∇ T x N • T • ∇ x N ]IdΩ (43)
and the spatial elasticity tensor C is given by

C = C osmo + C e ( 44
)
where C e is the spatial elasticity tensor for the solid matrix, which can be written as

C e = 4J -1 (F ⊗F ) : ∂ 2 W ∂C 2 : F T ⊗F T (45)
in which C = F T F is the right Cauchy-Green deformation tensor, and W is a strain energy density function for the solid. Here, a hyperelastic constitutive model [START_REF] Bower | Applied Mechanics of Solids[END_REF] is used and is expressed as

W = 1 2 µ J -2/3 I 1 -3 + 1 2 κ(J -1) 2 (46)
and C osmo in Eq.( 44) is the tensor of the osmotic modulus,

C osmo = J -1 (F ⊗F ) : 2 ∂ (-JpC -1 ) ∂C : F T ⊗F T = -J ∂p ∂J I ⊗ I + p(2I⊗I -I ⊗ I) ( 47
)
which results in part from the change in osmotic pressure with a change in tissue volume [START_REF] Ateshian | The correspondence between equilibrium biphasic and triphasic material properties in mixture models of articular cartilage[END_REF][START_REF] Chahine | Direct measurement of osmotic pressure of glycosaminoglycan solutions by membrane osmometry at room temperature[END_REF] .

The detailed expression of C osmo can be found in Azeloglu et al. [START_REF] Azeloglu | Heterogeneous transmural proteoglycan distribution provides a mechanism for regulating residual stresses in the aorta[END_REF] .

Application to experimental data

Uniaxial tensile tests were previously performed to stimulate and investigate chemoelastic effects in arteries [START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] . In these experiments, porcine descending thoracic aortas were excised and used as test samples. Propylene glycol (PG) and phosphate-buffered saline (PBS) were mixed to prepare an osmotically active solution.

Three rectangular samples (10 × 58mm), named Aorta1, Aorta2 and Aorta3, were immersed in an 80% PG solution and tested in uniaxial tension (stepwise stress-relaxation tests). The loading direction was aligned with (Y-axis). Strain fields were measured during these tests using optical coherence tomography (OCT) and digital volume correlation (DVC) [START_REF] Acosta Santamaría | Three-dimensional full-field strain measurements across a whole porcine aorta subjected to tensile loading using optical coherence tomography-digital volume correlation[END_REF] .

Before carrying out the stepwise uniaxial tensile tests, a 1.15 stretch preconditioning (displacement increment of 5.2mm) was applied. The initial distance between clamps was 35.1mm, as shown in Fig. 3.

Figure 3: Setup used for obtaining experimental datasets and acquisitions of experimental data using the OCT-DVC technique [START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] OCT volumetric images were acquired under a field of view (FOV) of 2 × 4 × 2.51mm

(X-, Y-and Z-axes, respectively, see Fig. 3(b)). Four specific regions of interest were defined on the OCT volume images: the entire aortic wall (global), the intima region, the media region and the adventitia region. 3D displacement fields were measured using DVC with a local correlation algorithm (LA-DVC). Ideal correlation conditions were reached directly in the 3D OCT image sequence acquisition and using a rigid body translation test. The implemented correlation parameters can be found in our previous work [START_REF] Acosta Santamaría | Three-dimensional full-field strain measurements across a whole porcine aorta subjected to tensile loading using optical coherence tomography-digital volume correlation[END_REF][START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] . The voxel size was a key parameter to reduce the uncertainties of the full-field measurements. The subset size parameter was used as (128, 96, 64, 32, 16, and 8 voxels). After measuring and fitting the displacement fields with tricubic functions, the components of the Green-Lagrange strain tensor were calculated with MatLab ® . Here the error of measurement is not considered and only the mean value of displacement measurement is used for simplification.

Results

First case study: homogeneous material with uniaxial tension

This first case study is a plate under uniaxial tension, as introduced in Section 2.4.1,

where an analytical solution is derived. A stretch λ is applied in the x direction at the right-side of the plate, and the constraints are applied as shown in Fig. 2. The parameters of the chemo-mechanically coupled constitutive model are reported in Table 2. 3),( 5), while the chemo part p will be increased with a lower c * as shown in Fig. 4(c), therefore, a lower total Cauchy stress is produced.

Parameter identification

Since the material considered in the first case study is homogeneous, there is only one unknown parameter, µ. Thus, one single linear virtual displacement field based on Eq.( 15)

is considered, that is, The shear modulus in solid part is well identified as shown in Table 3 regardless the existence of chemo term and the convergence process of µ is shown in Fig. 5. Due to many uncertainties in the chemo terms, we verify that the proposed VFM is immune to chemoelastic coupling effects by simulating three different external bath salt osmolarity environments and one different initial fixed-charge density in the tissue. The identified results of shear modulus µ are always equal to the exact shear modulus µ, as reported in Table 3. Moreover, we also test the identification results under different initial guess values and target(exact) values, the results are provided in Table 4, in which the initial values of shear modulus µ changes from 1000 Pa to 10000 Pa and the exact value of shear modulus µ changes from 130000 Pa to 150000 Pa, it is shown that the identified results maintain a good accuracy with relative error less than 1%.

u * =   0.2x -0.1y -0.1z   ε * =   0.2 -0.1 -0.1   . ( 48
)
Based on the noise analysis method in Section 2.3.4, we tested the influence of noise on the identification, as reported in Table 5. It is shown that the error remains less than 3%, with a relatively low noise level (σ u ≤ 5 × 10 -4 ) but gradually increases with the increase of noise level.

Second case study: multilayer material under biaxial tension

The second case study relies on the chemomechanically coupled constitutive model that we implemented within the Abaqus/Standard commercial FE software by means of the user subroutine UMAT, based on Eqs.( 3),( 4),( 5),( 45) and (47).

The geometry of this second case study is shown in Fig. 6 to simulate the layered structure of arteries. The size of the three-layers plate is 0.08×0.04×0.03m (i.e. l=0.08, m=0.04 and n=0.03) and the geometry is discretized with a fine mesh of 768 cubic elements. The boundary conditions are the same as the first case study in Fig. 2, with

u x | x=0 = u y | y=0 = u z | z=0 = 0.
Biaxial displacement loading is applied to the X and Y directions, including 1.30 stretch in X direction and 1.05 stretch in Y direction. Accordingly, the hyperelasitic constitutive model is used for arterial tissues in this example. The values for the constitutive parameters are those reported in Tables 6. Note here we use a large bulk modulus to make the material nearly incompressible.

FEM solutions

Firstly, we verify the accuracy of the FE solutions for direct analysis. As it is not possible to find the analytical solution for the three-layers examples, so we use the FE model introduced in Section 2.4.2 to solve the previous first case example with different bath osmolarity values. As it is shown in Fig. 7, it is found that there is a very good agreement between FE and analytical solutions.

Parameter identification

For the second case study, we can still consider the constant virtual strain field of Eq.(48), which is represented in Fig. 8(a). However, two other heterogeneous virtual fields are required to be introduced to ensure independence among the three VFM equations. Based on Eq.( 26), two other virtual fields can be determined as illustrated in Fig. 8(b) and 8(c) where parameters setting details can also be found. Displacement loading / pre The identification for each layer is verified by considering with different chemoelastic terms. As reported in Table 7, the same values of µ are identified even if we assume different values for the external salt osmolarity. This demonstrates that the proposed model is effective in shear modulus identification without knowing the chemoelatic terms. The process of minimizing the cost function f (µ 1 , µ 2 , µ 3 ) is provided in Fig. 9.

The influence of noise on the identification is reported in Table 8. The material parameters are well recovered with a low noise level (σ u < 5 × 10 -4 ) and it is also shown that the error of identified results is increased with the increase of noise level. However, when the noise level reaches a higher level (σ u = 1 × 10 -3 ), the results deviate significantly from the reference values.

Third case study: uniaxial tensile test of descending thoracic aortas

The third example consider the parameter identification based on the measurements in a uniaxial tensile test of descending thoracic aortas described in Section 2.5. Three samples with the same size were tested. The thickness ratio of each layer was set as 1:8:1 (adventitia:media:intima), according to the size of experimental sample [START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] . The geometry for the sample is shown in Fig. 10. The measurement for the average Green-Lagrange strains with the standard deviation(std) and the corresponding applied force in the different deformed configurations for each sample are reported in Tables 9, 10, and 11. on in these samples, the identified shear moduli show only marginal variations with the applied force, indicating that the identification method is able to eliminate the non-purely elastic effects. Indeed, it is assumed that these effects only affect the hydrostatic pressure to which the VFM is immune. This assumption permitted the separate identification of the shear modulus of each layer, as shown in Fig. 12. The largest shear moduli are obtained in the media layer, with values of approximately 1 MPa. The intima and adventitia appeared to be at least half as stiff. Although the obtained absolute values may be specific to the conditions used for the purpose of this study (immersion in PG to induce osmotic effects), the relative values between each layer indicate a significant gradient of elastic properties across the thickness.

Conclusion and discussion

The mechanical behavior of soft biological tissues can be significantly affected by the chemical potentials of interstitial fluids, resulting in obvious chemomechanical coupling phenomena. And chemoelastic constitutive models are widely used for soft biological tissues.

However, to the best of our knowledge, there have been very few studies related to parameter identification of these models using VFM combined with optical full-field measurement techniques.

In this paper, VFM was applied for the parameter identification of biphasic chemoelastic models using full-field measurements. The proposed approach avoids iterative resolutions Therefore, the proposed approach will be useful for further parameter identifications in soft biological tissues.

The proposed approach was assessed to identify the layer-specific stiffness properties of arteries based on strain fields. The identification results based on numerical simulations showed that the proposed model is immune to chemoelastic effects. Moreover, the shear modulus of the media of porcine descending thoracic aortas immersed in an 80% PG solution was identified with values of approximately 1 MPa, however, the identified shear modulus is smaller than the instantaneous shear modulus, which was identified to be approximately 3 MPa in a previous work [START_REF] Acosta Santamaría | Characterization of chemoelastic effects in arteries using digital volume correlation and optical coherence tomography[END_REF] . Indeed, the instantaneous shear modulus was found by assuming incompressibility of the solid part in the biphasic hyperelastic model, whereas identification with the VFM did not require such an assumption. Another interesting finding is that the layer-specific shear moduli can be identified simultaneously with the VFM. The shear modulus of the media layer is nearly twice as stiff as that of the intima and adventitia. The ratios between interlayer shear moduli of arteries are similar to those obtained in the work of Peña et al. [START_REF] Peña | Layer-specific residual deformations and uniaxial and biaxial mechanical properties of thoracic porcine aorta[END_REF] , in which they were approximately 1:2:1. Therefore, the proposed approach was proved to provide an effective identification of interlayer gradients of shear moduli despite the chemoelastic effects. A potential application of the proposed VFM approach is the multiscale mechanical characterization of biological tissues since both the effective material parameters at the macroscopic scale (as in Example 1) and the multilayer or heterogeneous material parameters at the mesoscopic scale (as in Examples 2 and 3) can be well identified, showing effective multiscale identification for soft biological tissues even in the presence of chemoelastic effects.

A main limitation in this study is that the VFM was only tested with the neo-Hookean model. Due to the obvious anisotropic behavior of aortic wall, future work will extend the approach to anisotropic hyperelastic models, such as the Holzapfel [START_REF] Holzapfel | A new constitutive framework for arterial wall mechanics and a comparative study of material models[END_REF] and MA-HGO [START_REF] Nolan | A robust anisotropic hyperelastic formulation for the modelling of soft tissue[END_REF] models which are used extensively to model collagen fibre reinforced biological materials, and the use of neo-Hookean model could produce a systematic variation of the apparent shear modulus with increasing strain. Another limitation of the work is that we ignore the stiffening effect of soft tissues for the sake of simplicity, a more complex VFM for nonlinear problems is required to consider this case.

In summary, we have presented a new VFM approach for parameter identification that is well suited for chemoelastic constitutive relationships in soft biological tissues. The identification of layer-specific elastic parameters can be applied even in the presence of complex chemoelastic effects. Important developments will extend the approach to multiscale constitutive models, and anisotropic effects of soft tissues.

FFigure 1 :F

 1 Figure 1: Multiplicative decomposition of the deformation gradient F

Figure 2 :

 2 Figure 2: The geometry and boundary conditions for a plate under tension

  Fig. 4(b), in which a smaller external bath salt osmolarity produces a relatively smaller total T 11 , it can be explained that as the T 11 has hyperelastic part and chemo part according to Eqs.(3),(5), while the chemo part p will be increased with a lower c * as shown in Fig.4(c),

Figure 4 :Figure 5 :

 45 Figure 4: Analytical results under different external bath salt osmolarities and initial fixed-charge density c F 0 = 200 meq/L in uniaxial tension. (a) Stretch β variation in the Z-direction.(b) Cauchy stress T 11 . (c) Osmotic pressure p.

Figure 6 :

 6 Figure 6: The geometry and displacement loading for the three-layers plate

Figure 7 :Figure 8 :

 78 Figure 7: Verification of the FE model against the analytical model on the results of Cauchy stress

Figure 9 :

 9 Figure 9: The identification process in the second case study. (a) Convergence plot of the novel VFM approach for each layer parameter. (b) Cost function f descent.

X

Figure 10 :

 10 Figure 10: Geometry and boundary condition for an aortic sample immersed in a PG solution and subjected to a uniaxial tensile test

  that the identified results for different layers of arteries are the same for ten different group of initial guesses. Then, the identification results for three samples using the proposed VFM are shown in Figs.11(a), 11(b), and 11(c), respectively. Despite the chemoelastic effects going

Figure 11 :

 11 Figure 11: Identified shear moduli of tested samples with average strains across each layer and the entire aortic wall. (a) Aorta1, (b)Aorta2,(c)Aorta3.

Figure 12 :

 12 Figure 12: Identified shear moduli (mean value and deviation bar) for each layer of the tested samples

Table 1 :

 1 Steps of the identification procedure for biphasic hyperelastic parameters 1. Forward simulation of chemomechanically coupled model and experimental tests in osmotically active solution. 2. Acquisition of full-field displacement from numerical calculation or experimental measurements. 3. Computation of deformation gradient F through the given full-field displacement. 4. Computation Cauchy stress T from chemomechanically constitutive equations.

Table 2 :

 2 Parameters of the chemo-mechanically coupled constitutive model

	Type	Description	Symbol Values	Units
		Universal gas constant	R	8.314	kg • m 2 s -2 K -1 mol -1
		Absolute temperature	θ	298	K
	Chemoelastic	Ambient pressure	p *	0	kPa
	term	Initial fixed-charge density	c F 0

Table 3 :

 3 Identified shear modulus in the first case study

	Type	Symbol			Value
	External bath salt osmolarity	c *	50 mosM 100 mosM 500 mosM 500 mosM
	Initial fixed-charge density	c F 0	40 meq/L 40 meq/L 40 meq/L 60 meq/L
	Tensile stretch	λ	1.2	1.2	1.2	1.2
	Initial guess shear modulus	µ	1000Pa	1000Pa	1000Pa	1000Pa
	Identified modulus	µ	140000Pa 140000Pa 140000Pa 140000Pa
	Exact modulus	µ	140000Pa 140000Pa 140000Pa 140000Pa

Table 4 :

 4 The identification results of different exact values and initial guesses

	Identified µ(Pa)		Exact values µ(Pa)
	and relative error(%)	130000 135000 140000 145000
		1000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
		2000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
		3000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
	Initial guess µ(Pa)	4000 5000	129533 (-0.36) 129533 (-0.36)	135192 (0.14) 135192 (0.14)	140000 (0) 140000 (0)	144500 (-0.34) 144500 (-0.34)	(0.08) (0.08)
		6000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
		7000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
		8000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
		9000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)
		10000	129533 (-0.36)	135192 (0.14)	140000 (0)	144500 (-0.34)	(0.08)

Table 5 :

 5 The relative error of identified results under different levels of measurement noises Standard deviation σ u 1 × 10 -4 3 × 10 -4 5 × 10 -4 7 × 10 -4 1 ×-3

	Exact µ(Pa)	140000	140000	140000	140000	140000
	Mean identified μ(Pa) 140349	139443	136382	131183	119993
	error	-0.25%	0.4%	2.58%	6.30%	14.29%

Table 6 :

 6 Parameters of the chemo-mechanically coupled constitutive model for second case study

	Type	Description		Symbol Values	Units
		Universal gas constant	R	8.314	kg • m 2 s -2 K -1 mol -1
		Absolute temperature	θ	298	K
	Chemoelastic	Ambient pressure	p *	0	kPa
	term	Initial fixed-charge density	c F 0	40	meq/L
		Water content	φ w 0	70	%
		External bath salt osmolarity	c *	[20,100,500,1000] mosM
		Material 1	Shear modulus Bulk modulus	µ 1 κ	100 6000	kPa kPa
	Hyperelastic term	Material 2	Shear modulus Bulk modulus	µ 2 κ	200 6000	kPa kPa
		Material 3	Shear modulus Bulk modulus	µ 3 κ	300 6000	kPa kPa

Table 7 :

 7 The identified results for different layers under different external bath osmolarities (%) 1.388 × 10 -3 5.788 × 10 -4 1.247 × 10 -3

	c *	Type	Material 1	Material 2	Material 3
		Initial value of µ(Pa)	150000	150000	150000
	20 mosM Identified µ(Pa)	100000.95	199999.81	299998.34
	100 mosM Identified µ(Pa)	99999.02	200001.77	300003.04
	500 mosM Identified µ(Pa)	100001.24	199998.37	299996.26
	1000 mosM Identified µ(Pa)	99997.62	200001.04	300006.52
		Exact µ(Pa)	100000	200000	300000
		Average relative error

Table 8 :

 8 The relative error of identified results in three layers under different levels of measurement noises

	σ u	μ1	Mean identified values μ(Pa) error(%) μ2 error(%)	μ3	error(%)
	1 × 10 -4	99622.35 0.38	200185.09 -0.09	298801.27 0.40
	3 × 10 -4	97340.56 2.66	197091.80 1.45	292038.95 2.65
	5 × 10 -4	92950.14 7.05	189106.28 5.45	279275.88 6.91
	7 × 10 -4	86117.79 13.88	175407.71 12.30	259909.74 13.36
	1 × 10 -3	63550.56 36.45	142574.59 28.71	209265.30 30.24
	Exact values 100000.00 -	200000.00 -	300000.00 -

Table 9 :

 9 The average Green-Lagrange strains from the experiment for Aorta1

	Type	Loading step	Measured stress(MPa)	E xx (%) (std) E yy (%) (std) E zz (%) (std)
		3rd	0.07	-1.510 (0.42) 4.352 (0.26) 0.145 (1.18)
	Global	5th 7th	0.1 0.13	-2.124 (0.52) 7.090 (0.37) 0.108 (2.32) -2.596 (0.71) 9.554 (0.51) 0.170 (3.56)
		9th	0.17	-3.125 (0.76) 11.847 (0.60) 0.353 (4.76)
		3rd	0.07	-1.558 (0.42) 4.381 (0.35) -2.271 (4.57)
	Intima	5th 7th	0.1 0.13	-2.136 (0.44) 7.157 (0.49) -3.528 (6.99) -2.600 (0.65) 9.603 (0.67) -4.768 (8.89)
		9th	0.17	-3.129 (0.74) 11.949 (0.79) -5.566 (10.18)
		3rd	0.07	-1.489 (0.39) 4.376 (0.23) 1.003 (0.34)
	Media	5th 7th	0.1 0.13	-2.127 (0.54) 7.112 (0.32) 1.603 (0.63) -2.582 (0.68) 9.570 (0.45) 2.263 (1.60)
		9th	0.17	-3.114 (0.72) 11.857 (0.52) 3.007 (2.67)
		3rd	0.07	-1.434 (0.67) 4.273 (0.35) -2.120 (3.15)
	Adventitia	5th 7th	0.1 0.13	-2.043 (0.85) 6.958 (0.46) -3.525 (4.96) -2.509 (1.05) 9.407 (0.65) -4.795 (6.73)
		9th	0.17	-3.020 (1.24) 11.630 (0.79) -6.137 (7.91)

Table 10 :

 10 The average Green-Lagrange strains from the experiment for Aorta2

	Type	Loading step	Measured stress(MPa)	E xx (%) (std) E yy (%) (std) E zz (%) (std)
		3rd	0.06	-0.628 (2.26) 3.572 (1.80) 0.493 (0.38)
	Global	5th 7th	0.1 0.13	-1.135 (3.60) 6.412 (2.92) 0.627 (0.87) -1.475 (4.93) 9.210 (4.10) 0.523 (1.63)
		9th	0.17	-1.982 (6.02) 11.708 (5.09) 0.385 (2.53)
		3rd	0.06	-0.370 (2.93) 4.134 (2.01) -2.470 (10.00)
	Intima	5th 7th	0.1 0.13	-0.717 (4.68) 7.579 (3.22) -3.383 (14.95) -1.105 (6.36) 10.846 (4.38) -4.156 (18.03)
		9th	0.17	-1.688 (8.12) 13.758 (5.22) -5.080 (19.37)
		3rd	0.06	-0.733 (2.55) 3.658 (1.88) 0.645 (1.60)
	Media	5th 7th	0.1 0.13	-1.309 (3.94) 6.555 (2.96) 0.919 (1.50) -1.726 (5.18) 9.555 (3.94) 1.152 (1.22)
		9th	0.17	-2.225 (6.19) 12.021 (5.15) 1.627 (1.19)
		3rd	0.06	-1.233 (2.91) 1.912 (1.63) -2.372 (5.29)
	Adventitia	5th 7th	0.1 0.13	-1.511 (4.12) 3.393 (2.70) -3.950 (7.58) -1.775 (4.62) 4.830 (4.07) -4.421 (8.65)
		9th	0.17	-2.076 (5.25) 5.926 (5.35) -4.922 (11.71)

Table 11 :

 11 The average Green-Lagrange strains from the experiment for Aorta3

	Type E Global Loading step Measured stress(MPa) 3rd 0.08 -1.315 (1.28) 4.940 (0.62) -0.269 (0.71) 5th 0.12 -2.303 (1.92) 8.377 (1.03) -0.708 (1.66) 7th 0.17 -2.923 (2.54) 11.174 (1.24) -1.020 (2.76)
		9th	0.21	-3.368 (3.08) 14.004 (1.60) -1.124 (3.80)
		3rd	0.08	-1.506 (1.13) 4.929 (0.58) -1.655 (3.22)
	Intima	5th 7th	0.12 0.17	-2.642 (1.82) 8.405 (0.89) -3.531 (5.50) -3.432 (2.45) 11.206 (1.02) -5.844 (8.40)
		9th	0.21	-3.978 (2.94) 14.030 (1.34) -7.480 (10.18)
		3rd	0.08	-1.235 (1.55) 4.909 (0.69) 0.156 (0.71)
	Media	5th 7th	0.12 0.17	-2.210 (2.17) 8.342 (1.06) 0.279 (0.70) -2.836 (2.88) 11.166 (1.26) 0.687 (1.08)
		9th	0.21	-3.287 (3.36) 14.069 (1.57) 1.472 (1.40)
		3rd	0.08	-1.367 (1.47) 4.860 (1.01) -2.878 (4.84)
	Adventitia	5th 7th	0.12 0.17	-2.213 (2.13) 8.097 (1.66) -5.542 (7.43) -2.590 (2.60) 10.910 (2.24) -8.711 (11.56)
		9th	0.21	-3.072 (3.25) 13.703 (2.71) -11.764 (13.49)

xx (%) (std) E yy (%) (std) E zz (%) (std)

Table 12 :

 12 The identification results for experimental measurements from different initial guesses

	Initial guesses µ(MPa)	Identified values µ(MPa)
	Intima Media Adventitia -Intima Media Adventitia
	0.1	0.1	0.1	-0.4746 1.0040	0.4831
	0.2	0.2	0.2	-0.4746 1.0040	0.4831
	0.3	0.3	0.3	-0.4746 1.0040	0.4831
	0.4	0.4	0.4	-0.4746 1.0040	0.4831
	1.0	1.0	1.0	-0.4746 1.0040	0.4831
	2.0	2.0	2.0	-0.4746 1.0040	0.4831
	3.0	3.0	3.0	-0.4746 1.0040	0.4831
	10.0	10.0	10.0	-0.4746 1.0040	0.4831
	1.0	2.0	3.0	-0.4746 1.0040	0.4831
	3.0	2.0	1.0	-0.4746 1.0040	0.4831

Firstly, we investigate the influence of initial guess on the identification results as shown in Table

12

, in which one group of experimental measurements in Table

9

is used, it is shown

Acknowledgements

The research leading to this paper was funded by the NSFC Grants [12072063, 11972109], the Fundamental Research Funds for the Central Universities [DUT21YG129], the Grants of State Key Laboratory of Structural Analysis for Industrial Equipment [S22403, GZ21104], and the Alexander von Humboldt Foundation [1217594].

Data availability statement

The data discussed in this work are available on request.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.


Pierre-Jean Benghozi

Recent information of a very different nature has shown the challenge posed by search engines in the internet economy.

-60% of people seeking health-related information consult a search engine. -On average, Americans spend 17 minutes a day on Google. -40% of online advertising revenue is directed towards search engines.

-Google captures 95% of European requests to search engines (75% in the USA).

Introduction

In a digital economy which is increasingly becoming an information and attention-based economy, search engines are occupying a central role shaped by the symbolic status of Google both in terms of its renown and its dominance of internet traffic.

Despite their importance, economic models and competition structures for search engines are not widely understood. Thanks to their top-heavy IT infrastructures, they come under the network economy. And because they supply service applications, they can also be classified as part of the content economy. Google's disruptive success is down to the very fact that it has been able to innovate in the technical, service, use and business model spheres, all at the same time.

This diversity explains why it is so difficult to grasp the search engine phenomenon on a global scale. The aim of this workshop is to describe the market by systematic consideration of all the problems posed by search engine development, in terms of regulation, economic models and relations with access suppliers.

Content access tools

Search engines are at the crossroads of information technology, telecommunications, marketing and the cultural sector. On the one hand, they orchestrate information exchanges on reading uses in return for information on readers, whilst on the other hand they orchestrate information exchanges on readers in return for advertising. Their technical dimension has considerable consequences on the range of content offered. In the first place, search engines have physical constraints linked to the time required to scan the entire web (3 or 4 days minimum). This means that they are not particularly suitable for "red hot" information, unlike portals based on publishing and selection activities (like Yahoo in its early days). Then again, despite the size of the web which they process, in order to be effective good general search engines mustn't surpass the expectations of internet users. They all tend to supply the same results, which isn't the case with specialized search engines, as their users expect them to come up with unexpected results.

Search engines are thus hard to pin down to an intermediary role between supply of content and demand for information. They are tools generating reading contents and part of an "access industry" based on the computerization of reading media. In this light, we could say that search engines represent a new stage of library industrialization, based on the same tool structures and organized according to library techniques, the index system, responses to specific requests and reading interfaces.

In economic terms, search engines differ from libraries in their relationship with the additional costs of little-used content matter: the "long tail". The negligible nature of these costs in the digital economy makes it possible to create a commercial model which is able to break free from the subsidies on which libraries are dependent. Search engines bring together supply and demand for content in a fractal sense; they can insert advertising -both charged and specific -at each interaction node because there is no pre-defined referencing system as you would find in a traditional library, with an existing encyclopedic knowledge organization system. Against this background, relationships between libraries and search engines form an interesting topic for study and in certain cases, it has even been suggested that readers' rights should be established to respond to this industrialization of the reading environment. These relationships exist in largescale digitization projects (Google Library, Quaero) or more global schemes within digital archives (accessibility, diversity, exclusion, property rights, fair use, investment costs, etc.).

The example set by a library such as the one in Lyons is an interesting case in point. The benefits that Google can offer the largest library in France's second city show that even now the largest libraries don't have the resources to digitize their catalogues. This process seems however crucial in order to reach a broader environment, as well as enhancing the availability, access tools and quality of the reading material supplied locally. This outsourcing of specific IT resources plays a part in adding value to document collections, but also creates new problems in its own right. Every web application turns now to be conceived as a 'library' these days, and traditional libraries have to compete with other operators. With digital technology, it is no longer enough for a library to hold external documents for consultation by internal readers; the trend is increasingly towards finding internal documents for distribution/promotion in the outside world.

New market forms

Search engines have many different business models and present a wide variety of economic configurations. Most are based on a "free" business model which makes the search engine market very similar to the media market. It strives to reach the widest possible audience and places great importance on the creation of commercial brands to guarantee a good reputation and appeal.

Search engines do not restrict themselves to redefining the supply and economy of contents. They also design new market forms for information and communication where values have shifted from the applications and the information to the audience. In other words, search engines are "third party" prescribers who influence both the content producers and the consumers. They form part of a fragmented system comprising a primary market in which consumers choose the content, a prescribers' market in which consumers select a portal or a search engine and a referral market in which search engines "sell" their promotion capacity to advertisers.

Search engines are thus causing large-scale changes in the intermediary economy and multi-sided platform by extending the limits of their economic model. Some disruptive innovations concerning organization or business models thus have the ability to challenge established competitive positions: e.g. when existing software publishers are challenged by online applications. This then creates more competition between alternative business models than between rivals within the same economic model. So whilst traditional internet participants target primary product markets by encouraging the integration of different functionalities within a specific technology to enable them to control consumers better, a large part of Google's success is down to the fact that it bases its economic model and remuneration system on the web referencing market. Supplying sponsorised contextual information alongside results is the cornerstone of Google's income and this is thus very different from an everyday advertising support function and more akin to an advertising consultancy role and a "placement" service. Customization and contextual adverts provide the "perfect market" in some ways for advertising by adapting the advertising level and cost precisely to each market niche.

This search engine business model, primarily advertising-based (98% of Google's income comes from Adwords), can also be regarded as one of its weaknesses. Others should not go unmentioned. The relevant search engine market is hard to define because there is considerable interaction with other markets (the advertisers' market) and it combines a number of advertising forms, online advertising, banners, sponsored links and search engines. It also has to contend with language barriers (note importance of Chinese, Russian or Arab markets). Finally, it has to deal with the emergence of information and communication forms which define new ways of distributing value, and especially with the development of RSS flows, widgets or closed social networks. From an economic standpoint, however, market analysis reveals that entry barriers are starting to be formed, based especially on IT storage capacity and information processing within servers.

So we need to ask ourselves whether the concentration trend we are noticing nowadays will still benefit Google and, in general, other search engines in the future, rather than other portal forms. In the history of search engines, the dominant players have always been in a virtual monopoly situation, but these leaders changed very regularly. Operators such as Exalead are an interesting case in point, as they illustrate a number of alternative dynamics. The company (website creator, consultant to technology suppliers) sells its search engine to companies and the Exalead search engine available on the web is simply a "showcase" which is not based on an economic model in its own right.

Regulations -what regulations?

The considerable growth in online advertising and the dominant role played by search engines in this sector have led to financial issues regarding the bottlenecks represented by the advertisers' market and audience control. Competition thus tends to arise not only between other digital players, but also in associated markets centred on advertising (e.g. the press). These economic and legal pressures lead to new questions on the regulation methods which apply in this sector: analysis of the search engine market, effects on associated sectors both upstream (content suppliers) and downstream (advertisers), new forms of competition emerging in the software and infrastructure market, monitoring competition and current financial trends (Yahoo/Microsoft merger talks). A number of major issues also seem to be affected: vertical relations between content publishers and access suppliers, copyright and intellectual property, trademark law, distribution and protection of personal data, the efficiency of the process of creating cultural programs and innovation.

The main question we need to ask is whether or not the search engine market should come under a generic set of regulations, and how market criteria can be applied to non-market objectives (freedom of expression, protection of privacy, etc.). Opinions often differ dramatically on this subject. First, we should remember that specific sector regulations (especially those of the telecoms sector) do not cover the search engine market, even though certain concepts might be useful such as the notion of a Minimum Quality of Service. Consequently, some believe there are no real problems with competition in this sector and that the regulatory authorities therefore don't need to get involved; even if problems do arise, a specific regulation would not be needed if traditional competition law was able to be applied.

On the other hand, others feel that the internet economy accentuates and reformulates some traditional questions on the subject of competition such as selectivity and discrimination between content suppliers (more or less well ranked and referenced by the search engines) or the risk of the main players blocking innovation. The Open Source strategy is sometimes quoted by way of example, as those players who get involved in publishing standards may seek to benefit from the community of developers to control radical innovations and prevent them from emerging.

From a competition viewpoint, relations between downstream and upstream of the value chain, between electronic communications operators and content suppliers (audiovisual in particular) are still the major source of potential conflict. Just as sites such as Google, eBay, Yahoo and Skype defend the creation of "open" platforms, including mobiles, to offset the control of subscribers by telecoms operators, the latter are simultaneously advocating the notion of passing on some of the investment needed in infrastructures to these sites by invoicing them for use of the network according to the applications used. This question has come directly to the fore in the USA, especially thanks to the debate on 'netneutrality' or on exclusivity clauses or tariff discrimination.

Intellectual property

Compared to traditional information aggregators and databases (directories, television guides, etc.), search engines raise specific legal issues regarding the protection of intellectual property rights. However, these questions arise in a radically different way because of the technology of search engines. Questions are no longer just about referencing and aggregation, but also the new status of tags and hyperlinks, the existence of "cache" files or the responsibilities of hosts. The justice system has already turned its attention to these legal issues surrounding search engines, but the debate extends way beyond a number of significant recent decisions such as the one involving Belgian press publishers.

Collecting and indexing information from the web represents the first source of conflict because search engines develop their market activities by exploiting content which they haven't paid for. They are then in direct competition with the information suppliers as regards advertisers, readers and audiences.

A second major source of conflict, which is both little known and often underestimated, concerns trademark law (affecting 20% of Google search requests). This trademark and advertising law is little used nowadays, but search engines could be attacked more than they are in this respect, since the referencing strategies used by many sites lead to illegal situations involving deception or overly similar advertising on the part of competitors. In this case, the responsibility and the role -active or passive -of the search engine have yet to be determined, especially if they are directly involved in suggesting key words.

Protection of personal data

Search engines need to know the background of internet users' requests in order to provide pertinent and personalized responses and, at the same time, target consumers by advertising using sponsored links. This is meant to improve the traceability of search requests and consultation processes by diversifying within applications so as to be able to enhance the process of personalizing the supply of information. Such investments are only justified if the personal data collected can be exploited and used for added value which will inevitably be at odds with the rights of individuals: the more or less anonymous nature of the collected information, transparency and length of time for which the information is kept, marketing of the acquired data.

The individual protection policies in force conflict, in material terms, with the diverse approaches and creation of public data in various countries, even though the networks are international and transcend borders, and despite the fact that the European Data Protection Directive does not apply to search engines.

In terms of keeping search requests on file, systematic scanning of e-mails by search engines, and indeed the exploitation of identity profiles as put forward in social networks, consumer protection requirements have led to heated controversies concerning the best solution to adopt. Some defend strict regulations imposing very short limits and data storage times with a view to defending individuals and their right to keep certain pieces of information confidential. Others feel that freedom of expression and personal choices can be more successfully defended by improving service performance within a transparent context, whilst all public documents should be available by free access, without any censorship as a general principle.

Conclusion

To date, the pressures described above have only come to light occasionally, since search engines cover very different areas of the legislative system such as competition, personal data and intellectual property in particular. Against the background of a developing and highly innovative internet, the originality of search engines calls for new approaches which the regulatory bodies seem poorly equipped to handle. This is the topic of current debates, especially in Europe.