# Two-time-scale inverse problems: formulation and solution

Guillaume Puel, Denis Aubry

<figure><img src="image_1.png" /></figure>

## To cite this version:

Guillaume Puel, Denis Aubry. Two-time-scale inverse problems: formulation and solution. Journal of Physics: Conference Series, 2015, 657 (012005), 10.1088/1742-6596/657/1/012005. hal-01338270

HAL Id: hal-01338270

https://hal.science/hal-01338270

Submitted on 13 Jul 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
PAPER • OPEN ACCESS

# Two-time-scale inverse problems: formulation and solution

To cite this article: G Puel and D Aubry 2015 *J. Phys.: Conf. Ser.* **657** 012005

View the [article online](#) for updates and enhancements.

## Related content

- *Effect of surface tension on the mode selection of vertically excited surface waves in a circular cylindrical vessel*
Jian Yong-Jun, E Xue-Quan, Zhang Jie et al.

- *Generalized Zakharov equations for a streaming or drifting plasma*
S Heuraux and G Leclert

- *Derivation of the low-frequency generalized Zakharov equations*
S Heuraux and G Leclert

<figure><img src="image_3.png" /></figure>

IOP ebooks™

Bringing together innovative digital publishing with leading authors from the global scientific community.

Start exploring the collection-download the first chapter of every title for free.
# Two-time-scale inverse problems: formulation and solution

G Puel and D Aubry

Laboratoire MSSMat (CentraleSupélec – CNRS UMR8579 – Université Paris-Saclay), Grande Voie des Vignes, 92290 Châtenay-Malabry, France

E-mail: {guillaume.puel,denis.aubry}@centralesupelec.fr

**Abstract.** The aim here is to study two-time-scale models and their associated parameter identification. When it is possible to consider two well-separated time scales, and when the fast component of the applied loading is periodic, a periodic time homogenization scheme, similar to what exists in space homogenization, can be used to derive a homogenized model. A parameter identification process for this latter is then proposed, consisting in homogenizing with respect to time a classical identification strategy based on the use of adjoint state formulations; it is then applied to an academic example showing the benefits of such a strategy.

## 1. Introduction

As far as the numerical calculation of a time-dependent model is concerned, the question of the computational cost can be of utmost relevance, especially when the considered model deals with fast phenomena, which require the use of very small time steps, when compared with the length of the time interval of study. In order to drastically reduce the computational cost, a periodic time homogenization method, such as [1], can be used when well-separated time scales can be defined, and when the fast component of the applied loading is periodic. The resulting homogenized model is then cheaper computationally, for it can be solved using time steps related to the slow time scale only, whereas the fast time scale is taken into account in an average way in the homogenization scheme.

In order them to give accurate predictions, such time-homogenized models have to be compared with experimental data. The key point is to define an identification strategy able to deal with such models in such a way that the process remains cheap and efficient. The aim of this paper is to analyze on a specific academic example the different choices made throughout the whole identification process, and what is their impact on the identification results.

## 2. Periodic time homogenization method

Periodic time homogenization, as it was initially proposed in [1], can be seen as a transposition to time of the classical periodic space homogenization methods, such as the techniques described in [2, 3]. It consists in separating two time scales, a slow one $t$ and a fast one $\tau = t/\xi$, by introducing, for every time-dependent quantity in the equations of the reference problem, an asymptotic expansion, in terms of the ratio $\xi$:

$$ \alpha(t, \tau) = \alpha_0(t, \tau) + \xi\alpha_1(t, \tau) + \xi^2\alpha_2(t, \tau) + O(\xi^3) \quad (1) $$

<figure><img src="image_4.png" /></figure>

Content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
When this ratio $\xi$ is very small, it is possible to consider that the two time scales are independent, and that any derivative with respect to time has to use the partial derivatives with respect to the two time scales:

$$d_t \alpha = \partial_t \alpha + \frac{1}{\xi} \partial_{\tau} \alpha \qquad (2)$$

where $d_t$, $\partial_t$, et $\partial_{\tau}$ stand for the total time derivative, the partial derivative with respect to the slow time and the partial derivative with respect to the fast time respectively. Moreover, if the applied loading has a component, which is periodic with respect to the fast time $\tau$, it is reasonable to assume that every variable is quasiperiodic, meaning that it is periodic with respect to $\tau$ (with associated period $T_f$).

Using the asymptotic expansion (1) for every quantity in the reference equations, and balancing the different orders of $\xi$, the time-homogenized equations are determined by averaging over a fast period $T_f$ the different quantities:

$$<\alpha> = \frac{\xi}{T_f} \int_0^{T_f/\xi} \alpha(t, \tau) \, d\tau \qquad (3)$$

allowing to separate slow-evolving phenomena from fast-time periodic components using the quasi-periodicity assumption. The residual quantities associated with this average are then denoted as $\alpha^* = \alpha - <\alpha>$, and depend on both time scales $t$ and $\tau$ a priori: they usually have to verify a very simple problem (e.g. linear elastic). Eventually, these homogenized equations are solved relatively to the slow time scale only, by introducing at each slow time step the averaged influence of the fast cycles corresponding to the solution of the fast problem, allowing to solve for all the time-homogenized variables with a drastically reduced computational cost.

References about the periodic time homogenization still tend to be scarce, as showed in [4], even if the scope of applications can be quite large: see for example [5, 6, 7, 8, 9]. Recently, we focused on validating the method for different cases of simulations of structures withstanding fatigue loads with two periodic components:

* material fatigue with a viscoplastic law defining two hardening variables in [10];

* material fatigue with an isotropic damage law in [11];

* extension of the method to three different time scales in [4, 12];

* extension of the method to the case of resonant excitation in [4].

### 3. Gradient-based parameter identification with adjoint state formulations

In this section, a general identification strategy is proposed, which is relevant for both reference and time-homogenized models: it will thus be illustrated in this latter case in Section 4. The forward problem is considered as an implicit formulation with a vector function $\mathcal{F}$ considering $\mathbf{u}(t)$ over a time interval $[0, T]$:

$$\begin{gathered} \mathcal{F}(\mathbf{u}(t), d_t \mathbf{u}(t), d_t^2 \mathbf{u}(t), \mathbf{p}, t) = \mathbf{0} \qquad (4) \\ \mathbf{u}(0) = \mathbf{U}_0 \qquad d_t \mathbf{u}(0) = \mathbf{V}_0 \end{gathered}$$

where $d_t$ and $d_t^2$ are the first- and second-order time derivatives respectively. $\mathbf{U}_0$ and $\mathbf{V}_0$ stand for the initial conditions of the dynamic problem. Whereas $\mathbf{u}$ is the state vector of size $N$, composed of all the time-dependent degrees of freedom (DOFs) describing the studied problem, $\mathbf{p}$ stands for the vector containing the $P$ scalar parameters associated with the differential equation (4).

The identification problem consists in finding the parameter vector $\mathbf{p}_{opt}$ such that the solution $\mathbf{u}(t; \mathbf{p}_{opt})$ of (4) obtained with the parameters $\mathbf{p}_{opt}$ is as close to the available experimental data as possible. These latter are indeed compared with the corresponding quantities $\mathbf{A}\mathbf{u}(t; \mathbf{p})$,
where $\mathbf{A}$ is a projection operator allowing to select, for each quantity, the closest DOF to the experimental measurement point. In order to use consistent notations, the corresponding experimental quantity is denoted $\mathbf{A}u_{exp}(t)$; however, it does not mean that such a vector $u_{exp}(t)$ actually exists.

The following misfit function is then introduced: it consists of a norm measuring the discrepancy between the quantities predicted with the forward model (4) and experimental data:

$$ \mathcal{J}(\mathbf{p}) = \frac{1}{2} \int_{0}^{T} |\mathbf{A}(\mathbf{u}(t; \mathbf{p}) - \mathbf{u}_{exp}(t))|^2 dt + \frac{1}{2} |\mathbf{R}(\mathbf{p} - \mathbf{p}_0)|^2 \quad (5) $$

where $\mathbf{u}(t; \mathbf{p})$ satisfies Equation (4). The $L^2$-norm proposed here is completed with a Tikhonov regularization term, allowing to deal with the ill-posedness of the identification problem, by bounding the magnitude of the parameter vector $\mathbf{p}$ to be identified: this regularization term uses a vector $\mathbf{p}_0$ containing nominal values corresponding to a priori experience, and a diagonal weighing matrix $\mathbf{R}$. Eventually, the solution of the identification problem can be sought as the parameter vector minimizing the misfit function $\mathcal{J}(\mathbf{p})$:

$$ \mathbf{p}_{opt} = \arg \min_{\mathbf{p}} \mathcal{J}(\mathbf{p}) \quad (6) $$

The determination of this minimum is achieved using gradient-based minimization methods, therefore the question of avoiding local minima by means of an appropriate regularization process should be carefully addressed. In some cases, rather than using the classical Tikhonov regularization term, the a priori experience can be introduced in some specific ways, as in [13]. Similarly, the fact of using a homogenized model in the parameter identification process can introduce a regularizing effect, just as explained in [14]. However, we will not address here this specific question, but rather focus on the identification process itself.

To estimate the gradient of the misfit function, we solve an adjoint state problem. A typical example in mechanical engineering is given in [15], where the parameters of an elastoplastic material law are identified with indentation tests. In the strategy proposed here, the generic form of the adjoint state problem is as follows:

$$ \begin{gathered} \nabla_{\mathbf{u}} \mathcal{F}^{\mathrm{T}} \mathbf{z} - d_t (\nabla_{d_t \mathbf{u}} \mathcal{F}^{\mathrm{T}} \mathbf{z}) + d_t^2 (\nabla_{d_t^2 \mathbf{u}} \mathcal{F}^{\mathrm{T}} \mathbf{z}) = \mathbf{A}^{\mathrm{T}} \mathbf{A} (\mathbf{u} - \mathbf{u}_{exp}) \quad (7) \\ (\nabla_{d_t^2 \mathbf{u}} \mathcal{F}^{\mathrm{T}} \mathbf{z})|_{t=T} = \mathbf{0} \qquad (\nabla_{d_t \mathbf{u}} \mathcal{F}^{\mathrm{T}} \mathbf{z})|_{t=T} = \mathbf{0} \end{gathered} $$

where $\nabla_{\mathbf{u}} \mathcal{F}$, $\nabla_{d_t \mathbf{u}} \mathcal{F}$ and $\nabla_{d_t^2 \mathbf{u}} \mathcal{F}$ stand for the directional derivatives of $\mathcal{F}$ with respect to $\mathbf{u}$, $d_t \mathbf{u}$ and $d_t^2 \mathbf{u}$ respectively. The adjoint state problem is then a time-backward differential equation with final conditions, and where the first-order sensitivities of the forward problem are concerned.

Once the adjoint state problem (7) is solved, it can be shown that the misfit function gradient with respect to the parameter vector $\mathbf{p}$ can be expressed as:

$$ \nabla_{\mathbf{p}} \mathcal{J}(\mathbf{p}) = \mathbf{R}^{\mathrm{T}} \mathbf{R}(\mathbf{p} - \mathbf{p}_0) - \int_{0}^{T} \nabla_{\mathbf{p}} \mathcal{F}^{\mathrm{T}} \mathbf{z}(t) dt \quad (8) $$

This specific way of estimating the misfit function gradient can be compared with a classical finite difference formula, such as the central finite difference scheme: in this latter case, when the parameter vector is of size $P$, the gradient calculation is obtained by evaluating the misfit function in $2P$ additional ‘points’, each couple of points corresponding to two symmetrical perturbations of the misfit function associated with each parameter in the vector $\mathbf{p}$. The resulting computational cost for each gradient evaluation consists of $2P$ solutions of the forward
problem (4) and 2P time integral evaluations. By contrast, when the adjoint state solution is used, only two differential equation solutions are required: one for the forward problem (4) and one for the adjoint problem (7). The associated computational cost for each gradient evaluation is then 2 differential equation solutions, and P time integral evaluations. The resulting gain is as high as the number of parameters to be identified. Moreover, it is easier to control the accuracy of the gradient estimate with the adjoint state method than with finite difference formulas, for which the choice of the discretization steps has a strong influence on the final estimate.

## 4. Parameter identification of a model with two time scales

### 4.1. Forward reference and time-homogenized problems

Here an academic example is proposed to discuss the different steps of the parameter identification problem associated with a time-homogenized model. This latter consists of a straight bar of length L withstanding at one end a normal force $f_s$ with two periodic components.

$$ \partial_x \sigma + c_K \mathrm{d}_t \partial_x \sigma = \rho \mathrm{d}_t^2 u \qquad (9) $$

$$ u|_{x=0} = 0 \qquad (\sigma + c_K \mathrm{d}_t \sigma)|_{x=L} = f_s \qquad (10) $$

$$ \sigma = E (\partial_x u - \varepsilon^p) \qquad (11) $$

$$ \mathrm{d}_t \varepsilon^p = \left( \frac{|\sigma|}{K} \right)^n \mathrm{sign} \sigma \qquad (12) $$

where $\partial_x$ is the partial space derivative, and zero initial conditions are assumed for the displacement $u$, $\mathrm{d}_t u$ and the plastic strain $\varepsilon^p$. $K$ and $n$ are the two parameters to be identified, using the measured displacement $u_{exp}$ at $x = L$, whereas $E$, $\rho$ and $c_K$ are assumed known.

The time-homogenized equations come from the zeroth-order expression of the forward problem (9)-(10)-(11)-(12), where the fast-time average (3) has been previously applied:

$$ \partial_x \langle \sigma_0 \rangle = 0 \qquad (13) $$

$$ \langle u_0 \rangle_{x=0} = 0 \qquad \langle \sigma_0 \rangle_{x=L} = \langle f_s \rangle \qquad (14) $$

$$ \langle \sigma_0 \rangle = E (\partial_x \langle u_0 \rangle - \varepsilon_0^p) \qquad (15) $$

$$ \partial_t \varepsilon_0^p = \left\langle \left( \frac{|\sigma_0|}{K} \right)^n \mathrm{sign} \sigma_0 \right\rangle \qquad (16) $$

where it has been assumed that $\rho L^2/(ET_f^2) \le O(1)$ and $c_K/T_f \le O(1)$, and where it can be shown that $\varepsilon_0^p(x, t, \tau) = \varepsilon_0^p(x, t)$. The fast-time average term associated with the evolution equation (16) is evaluated by means of a numerical integration formula, such as the trapezoidal rule, and uses the solution of the residual system defined as the difference between the reference problem and the zeroth-order time-homogenized system (13)-(14)-(15)-(16).

### 4.2. Identification strategy

The first step consists in describing, through the misfit function, the discrepancy between the time-homogenized model's predictions and experimental data: indeed, on the one hand, the model has been solved on slow time steps only, whereas, on the other hand, the experimental data can be available on a much finer scale.

The most efficient choice in terms of computation cost is to use time-homogenized quantities in the misfit function, since it allows to address the time integral on slow time steps only:

$$ \mathcal{J}_0(K, n) = \frac{1}{2} \int_0^T |\langle u_0 \rangle (L, t; K, n) - \langle u_{exp} \rangle (t)|^2 \mathrm{d}t \qquad (17) $$

Whereas $\langle u_0 \rangle (L, t; K, n)$ is the solution of the zeroth-order time-homogenized forward problem (13)-(14)-(15)-(16), $\langle u_{exp} \rangle (t)$ stands for the corresponding experimental quantity,
which is obtained by fast-averaging the experimental data for each slow time step $t_k$ of the time-homogenized displacement:

$$ \langle u_{exp} \rangle (t_k) = \frac{1}{T_f} \int_{t_k}^{t_k+T_f} u_{exp}(t) dt \quad (18) $$

The gradient of this misfit function is evaluated as in Section 3, using the solution of an adjoint state problem. This latter is as follows, with $\lambda_0(x, t)$ only:

$$ \partial_t \lambda_0 = \langle u_0 \rangle_{|x=L} - \langle u_{exp} \rangle \quad (19) \\ \lambda_0|_{t=T} = 0 $$

and it can be shown that it corresponds to the zeroth-order time-homogenized version of the time-backward adjoint state problem associated with the identification process for the reference problem. Indeed, Figures 1 and 2 show the comparison between the two corresponding adjoint state solutions for the identification process detailed previously: the adjoint state solution corresponding to the time-homogenized problem is homogeneous along the bar, and is very close to the (homogeneous) fast-time average of the (heterogeneous) adjoint state solution associated with the reference problem. This equation can be solved using the slow time steps $t_k$ only, which allows to derive the solution in a way as efficient as for the time-homogenized forward solution (13)-(14)-(15)-(16).

<figure><img src="image_5.png" /><figcaption>**Figure 1.** Adjoint state solutions for reference (in red) and time-homogenized (in blue) problems, at $x = 0$ and $x = L$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 2.** Zoom of Figure 1 (same color conventions).</figcaption></figure>

The misfit function's gradient then consists of the two following partial derivatives:

$$ \frac{\partial \mathcal{J}_0}{\partial K} = \int_0^T \int_0^L \left\langle \frac{n}{K} \left( \frac{E |\partial_x u_0 - \varepsilon_0^P|}{K} \right)^n \mathrm{sign} (\partial_x u_0 - \varepsilon_0^P) \right\rangle \lambda_0 \, dx \, dt \quad (20) $$

$$ \frac{\partial \mathcal{J}_0}{\partial n} = - \int_0^T \int_0^L \left\langle \left( \frac{E |\partial_x u_0 - \varepsilon_0^P|}{K} \right)^n \mathrm{sign} (\partial_x u_0 - \varepsilon_0^P) \log \left( \frac{E |\partial_x u_0 - \varepsilon_0^P|}{K} \right) \right\rangle \lambda_0 \, dx \, dt \quad (21) $$

Once again, it can be shown that these relations correspond to the zeroth-order time-homogenized estimates of the two misfit function's gradient components obtained for the identification problem associated with the reference problem. Actually, this is a result that already occurs in periodic space homogenization, as shown for example in [16].

In order to evaluate how the identification process performs, synthetic data $u_{exp}(t)$ are created by solving the reference forward problem (13)-(14)-(15)-(16), using $K_{exp} = 100 \cdot 10^6$ SI units and $n_{exp} = 10$ as parameter values. The bar withstands a two-periodic loading with $F = 0.05$ Hz and $F/\xi = 500$ Hz. $K_0 = 50 \cdot 10^6$ SI units and $n_0 = 5$ are chosen as initial parameter values for
the identification process, which is based on an interior-reflective Newton method [17] in order to minimize the misfit function $J^0$. Since we do not want to address specifically the question of regularization, no Tikhonov term is added. Figure 3 shows the comparison between the identified model ($K_{id} = 100.89 \cdot 10^6$ SI units and $n_{id} = 9.94$) and the synthetic reference, more precisely the variations of the longitudinal plastic strain, which is not directly observable. When compared with what is obtained when the inverse problem related to the reference problem (9)-(10)-(11)-(12) is considered, the computational cost associated with the identification process is significantly reduced ($10^5$ time steps instead of $10^7$ here).

<figure><img src="image_7.png" /><figcaption>**Figure 3.** Parameter identification of a time-homogenized model: experimental (in red) and predicted (in blue) longitudinal plastic strain.</figcaption></figure>

**5. Conclusion**

Here we have proposed a first preliminary study of a two-time-scale parameter identification process, using time-homogenized models: the adaptation of a classical identification strategy based on an adjoint state formulation to estimate the misfit function's gradient can be used in this specific framework: this leads, on the proposed example, to the determination of the time-homogenized counterpart of the adjoint solution associated with the reference identification problem. Despite its simplicity, the academical example studied here showed the relevance of the strategy and its reduced computational cost: these results can be viewed as a first step before dealing with more complex cases of study.

**References**

[1] Guennouni T and Aubry D 1986 *Comptes rendus de l'Académie des sciences. Série II* **303** 1765-1768

[2] Bensoussan A, Lions J L and Papanicolaou G 1978 *Asymptotic Analysis for Periodic Structures* (Burlington, MA: Elsevier)

[3] Sanchez-Palencia E 1980 *Non-homogeneous media and vibration theory* (Berlin Heidelberg New York: Springer-Verlag)

[4] Puel G and Aubry D 2014 *International Journal for Multiscale Computational Engineering* **12** 291-318

[5] Yu Q and Fish J 2002 *Computational Mechanics* **29** 199-211

[6] Shishkina E V, Blekhman I I, Cartmell M P and Gavrilov S N 2008 *Nonlinear Dynamics* **54** 313-331

[7] Lévy A, Le Corre S, Poitou A and Soccard E 2011 *International Journal for Multiscale Computational Engineering* **9** 53-72

[8] Dontsov E V and Guzina B B 2012 *Journal of the Mechanics and Physics of Solids* **60** 1791-1813

[9] Papon A, Yin Z Y, Riou Y and Hicher P Y 2013 *International Journal for Numerical and Analytical Methods in Geomechanics* **37** 1470-1491

[10] Aubry D and Puel G 2010 *Procedia Engineering* **2** 787-796

[11] Devulder A, Aubry D and Puel G 2010 *Computational Mechanics* **45** 637-646

[12] Puel G and Aubry D 2012 *European Journal of Computational Mechanics* **21** 312-324

[13] Puel G, Bourgeteau B and Aubry D 2013 *Mechanical Systems and Signal Processing* **36** 354-369

[14] Chung C B and Kravaris C 1988 *Inverse Problems* **4** 973-994

[15] Constantinescu A and Tardieu N 2001 *Inverse Problems in Engineering* **9** 19-44

[16] Kesavan S and Saint Jean Paulin J 1997 *SIAM Journal on Control and Optimization* **35** 1557-1573

[17] Coleman T F and Li Y 1994 *Mathematical Programming* **67** 189-224