# Control of a cantilever array by Periodic Networks of Resistances

H. Hui, Y. Yakoubi, Michel Lenczner, N. Ratier

► To cite this version:

H. Hui, Y. Yakoubi, Michel Lenczner, N. Ratier. Control of a cantilever array by Periodic Networks of Resistances. 11th International Conference on Thermal, Mechanical and Multiphysics Simulation and Experiments in Micro/Nanoelectronics and Systems (EuroSimE 2010), Apr 2010, Bordeaux, France. pp.1-6, 10.1109/ESIME.2010.5464558 . hal-00677432

HAL Id: hal-00677432

https://hal.science/hal-00677432v1

Submitted on 14 Apr 2021

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_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Control of a Cantilever Array by Periodic Networks of Resistances

H. Hui¹,³ Y. Yakoubi² M. Lenczner³ and N. Ratier³

¹School of Mechatronic Northwestern Polytechnical University, 127 Youyi Xilu, 710072 Xi’an Shaanxi, China

²Laboratoire Jacques-Louis Lions, Université Pierre et Marie Curie, Boîte courrier 187, 75252 Paris Cedex 05, France

³FEMTO-ST, Département Temps-Fréquences Université de Franche-Comté,
26, Chemin de l'Epitaphe, 25030 Besançon Cedex, France

Email: hui.hui@femto-st.fr, yyakoubi@ann.jussieu.fr, michel.lenczner@utbm.fr and nicolas.ratier@femto-st.fr

## Abstract

In this paper, we present a two-scale model including an optimal active control for a one-dimensional cantilever array with regularly spaced actuators and sensors. With the purpose of implementing the control in real time, we propose an approximation that may be realized by an analog distributed electronic circuit. More precisely, our analog processor is made by Periodic Network of Resistances (PNR). The control approximation method is based on two general concepts, namely on functions of operators and on the Dunford-Schwartz representation formula. We conducted validations of the control approximation method as well as of its effect in the complete control loop.

## 1. Introduction

In the past decade, a number of papers have been focused on semi-decentralized distributed optimal control for systems with distributed actuators and sensors. Most of them are dealing with infinite length systems, see [1] and [10] for systems governed by partial differential equations, and [3] for discrete systems. In the papers [4] and [5] the authors have introduced an approximation of an optimal control to a finite length beam endowed with a periodic distribution of piezoelectric sensors and actuators. Even if it was giving satisfactory results, it was suffering from some limitations. In [9] it has been extended so that to cover a larger range of systems and to increase its precision and robustness. Indeed, the new method does not require that each operator of the state equation and of the cost functional be functions of a same operator but they must be only functions of a same operator up to some change of variable operators. Regarding precision, the Taylor series approximating a function of an operator has been replaced by the use of the Dunford-Schwartz representation formula followed by a quadrature rule for the contour integral.

Here we apply our new method to a recently developed and validated two-scale model of cantilever arrays, submitted in the paper [8]. It is rigorously justified thanks to an adaptation of the two-scale approximation method introduced in [6] and detailed in [7]. Its main advantage is that in the same time it requires little computing effort and it is reasonably precise.

This paper presents results from an implementation of the new semi-decentralized optimal control strategy on the two-scale model of cantilever arrays. We provide results regarding precision and cost. However our calculations have been carried out using the simplest optimal control strategy, namely a Linear Quadratic Regulator. As in [5], we

also provide a realization of the semi-decentralized control scheme through a Periodic Network of Resistances (PNR), implementing a finite difference scheme for the partial differential operator in the Dunford-Schwartz formula. Finally, we quote that the entire approach can be extended to other linear optimal control problems, i.e. LQG or $H_\infty$ controls as well as to more physical actuating and sensing principles.

## 2. A Two-Scale Model of Cantilever Arrays

We consider a one-dimensional cantilever array comprised of an elastic base, and a number of clamped elastic cantilevers with free end, see Figure 1. Assuming that the number of cantilevers is sufficiently large, a homogenized model was derived using a two-scale approximation method. This is reported in the detailed paper [7] devoted to static regime. The corresponding model extended to dynamic regime is introduced in the letter [6]. The modelling papers were written in view of Atomic Force Microscopy application.

<figure><img src="image_2.png" /><figcaption>Fig. 1: Array of Cantilevers</figcaption></figure>

After a number of simplifications, the approximate homogenized model expressed in the two-scale referential, which is a rectangle $\Omega = (0, L_B) \times (0, L_C^*)$. The parameters $L_B$ and $L_C^*$ represent respectively the base length in the macroscale $x$-direction and the scaled cantilever length in the microscale $y$-direction. The base is modelled by the line $\Gamma = \{(x, y) | x \in (0, L_B) \text{ and } y = 0\}$, and the rectangle $\Omega$ is filled by an infinite number of cantilevers. We describe the system motion by its bending displacement only. So, the base is governed by an Euler-Bernoulli beam equation with two kinds of distributed forces, one exerted by the attached cantilevers and the other, denoted by $u(t, x, 0)$, originating from an actuator distribution. The bending displacement, the mass per unit length, the bending coefficient and the scaled cantilever width being denoted by $w(t, x, 0)$, $\rho^B$, $R^B$ and $\ell_C^*$, the base governing equation states

$$ \rho^B \partial_{tt}^2 w + R^B \partial_{x \cdots x}^4 w + \ell_C^* R^C \partial_{yyy}^3 w = u \text{ in } \Gamma. \quad (1) $$

The base is assumed to be clamped, so the boundary con-
ditions are

$$w = \partial_x w = 0, \qquad (2)$$

at its ends. Each cantilever is oriented in the y-direction, and its motion is governed by the Euler-Bernoulli equation distributed along the y-direction. It is subjected to a control force $u(t, x, y)$ taken as distributed along each whole cantilever. It can be replaced by any other realistic force distribution. Denoting by $w(t, x, y)$, $\rho^C$ and $R^C$ cantilever bending displacements, the mass per unit length, and the bending coefficient, the governing equation in $(x, y) \in \Omega$ is

$$\rho^C \partial_{tt}^2 w + R^C \partial_{y\dots y}^4 w = u, \qquad (3)$$

endowed with the boundary conditions

$$\left\{ \begin{array}{ll} \partial_y w = 0 & \text{at } y=0, \\ \partial_{yy}^2 w = \partial_{yyy}^3 w = 0 & \text{at } y = L_C^* \end{array} \right. \qquad (4)$$

representing an end clamped in the base, and a free end. The weak formulation associated to (1-4) states as

$$\int_0^{L_B} (\rho^B \partial_{tt}^2 w + R^B \partial_{xx}^2 w \partial_{xx}^2 v)|_{\Gamma} dx \\ + l_C^* \int_{\Omega} \rho^C \partial_{tt}^2 w + R^C \partial_{yy}^2 w \partial_{yy}^2 v dydx \qquad (5) \\ = \int_0^{L_B} (u v)|_{\Gamma} dx + l_C^* \int_{\Omega} u v dydx,$$

for any regular function $v$, satisfying in particular the conditions: $v = \partial_x v = 0$ at both end of the base and $\partial_y v = 0$ at $y = 0$ at the junction.

### 3. Model Reformulation

To simplify the model, but keeping its distributed feature, we discretize in the y-direction projecting on a basis $K_n(y) = \int_0^y y T'_n(y) dy$, where $T_n(y)$ is the basis of Chebyshev polynomial. We define the approximations of the displacement and of the control

$$\left\{ \begin{array}{l} w(t, x, y) \approx \sum_{n=1}^{N} w_n(t, x) K_n(y), \\ u(t, x, y) \approx \sum_{n=1}^{N} u_n(t, x) K_n(y), \end{array} \right. $$

where $w_n(t, x)$ and $u_n(t, x)$ are the polynomial coefficients in the approximation of $w$ and $u$ respectively. We also choose $v \approx \sum_{m=1}^N v_m(t, x) K_m(y)$, so we find that $(w_n(t, x))_{n=1,2,\dots,N}$ are the solutions to a set of equations posed on $\Gamma$,

$$\sum_{n,m=1}^{N} M_{m,n} \partial_{tt}^2 w_n + K_{m,n}^B \partial_{x\dots x}^4 w_n \\ + K_{m,n}^C w_n = \sum_{n,m=1}^{N} \tilde{B}_{m,n} u_n \text{ in } [0, \infty) \times \Gamma. \qquad (6)$$

The boundary conditions are $\sum_{n=1}^N w_n(t, 0) K_n(0) = \sum_{n=1}^N \partial_x w_n(t, 0) K_n(0) = 0$ and $\sum_{n=1}^N w_n(t, L_B) K_n(0) =$

$\sum_{n=1}^N \partial_x w_n(t, L_B) K_n(0) = 0$. In (6), we use the notations

$$\begin{aligned} M_{m,n} &= \rho^B (K_m K_n)|_{\Gamma} + l_C^* \rho^C \int_0^{L_C} K_m K_n \, dy, \\ K_{m,n}^B &= R^B (K_m K_n)|_{\Gamma}, \\ K_{m,n}^C &= l_C^* R^C \int_0^{L_C} \partial_{yy}^2 K_m \partial_{yy}^2 K_n \, dy, \\ \tilde{B}_{m,n} &= (K_m K_n)|_{\Gamma} + l_C^* \int_0^{L_C} K_m K_n \, dy. \end{aligned}$$

The LQR problem is set for control variables $(u_n)_{n=1,2,\dots,N} \in L^2(\Gamma)^N$ and for the cost functional

$$\mathcal{J} = \int_{0}^{+\infty} \sum_{n=1}^{N} \| \partial_{xx}^2 w_n(t, x) \|_{L^2(\Gamma)}^2 \qquad (7) \\ + \| u_n(t, x) \|_{L^2(\Gamma)}^2 \, dt.$$

The choice of the functional is related to vibration stabilization of the microcantilever array.

### 4. Classical Formulation of the LQR Problem

Now, we write the above LQR problem in a classical abstract setting, see [2], even if we do not detail the functional framework. We set $z^T = (w_n \ \partial_t w_n)_{n=1,2,\dots,N}$ the state variable, $u^T = (u_n)_{n=1,2,\dots,N}$ the control variable, $A = \begin{pmatrix} 0_{N \times N} & I_{N \times N} \\ -(M^{-1}(K^B \partial_x^4 + K^C))_{N \times N} & 0_{N \times N} \end{pmatrix}$ the state operator, $B = \begin{pmatrix} 0_{N \times N} \\ (M^{-1}\tilde{B})_{N \times N} \end{pmatrix}$ the control operator, $C = \begin{pmatrix} \partial_{xx}^2 I_{N \times N} & 0_{N \times N} \\ 0_{N \times N} & 0_{N \times N} \end{pmatrix}$ the observation operator, and $S = I_{N \times N}$ the weight operator. Consequently, the LQR problem, consisting in minimizing the functional under the constraint (6), may be written under its usual form as

$$\begin{array}{l} \partial_t z(t, x) = Az(t) + Bu(t) \\ \text{for } t > 0 \text{ and } z(0) = z_0, \end{array} \qquad (8)$$

with the minimized cost functional (7). We know that $(A, B)$ is stabilizable and that $(A, C)$ is detectable, in the sense that the system is controllable and observable. It follows that for each $z_0$, the LQR problem (8) admits a unique solution

$$u^* = -Kz, \qquad (9)$$

where $K = S^{-1}B^*P$, and $P$ is the unique self-adjoint non-negative solution to the operational Riccati equation

$$A^*P + PA - PBS^{-1}B^*P + C^*C = 0. \qquad (10)$$

### 5. Semi-Decentralized Approximation

This Section is devoted to formulate the approximation method. The mathematical derivation has been introduced in a paper [9]. We denote by $\Lambda$, the mapping: $\Lambda: f \to w$, where $w$ is the unique solution of $\partial_x^4 w = f$ in $\Gamma$ with the boundary conditions $w = \partial_x w = 0$ for $x = \{0, L_B\}$. The spectrum $\sigma(\Lambda)$ is discrete and made up of real eigenvalues $\lambda_k$. They are solutions to the eigenvalue problem $\Lambda\phi_k = \lambda_k\phi_k$ with $\|\phi_k\|_{L^2(\Gamma)} = 1$. In the sequel, $I_\sigma = (\sigma_{\min}, \sigma_{\max})$ refers to an open interval that includes the complete spectrum.
## 5.1 Factorization of K by a Matrix of Functions of Λ

In this part, we introduce the factorization of the controller K under the form of a product of a matrix of functions of Λ. To do so, we introduce the change of variable operators $\Phi_Z = \begin{pmatrix} \Lambda^{\frac{1}{2}} & 0 \\ 0 & I \end{pmatrix}$, $\Phi_U = I$ and $\Phi_Y = \begin{pmatrix} \partial_x^2 \Lambda^{\frac{1}{2}} & 0 \\ 0 & I \end{pmatrix}$, from which we introduce the matrices of functions of Λ, $a(\Lambda) = \Phi_Z^{-1} A \Phi_Z$, $b(\Lambda) = \Phi_Z^{-1} B \Phi_U$, $c(\Lambda) = \Phi_Y^{-1} C \Phi_Z$ and $s(\Lambda) = \Phi_U^{-1} S \Phi_U$, simple to implement on a semi-decentralized architecture. A straightforward calculation yield

$$ 
\begin{aligned}
a(\lambda) &= \begin{pmatrix} 0 & I \\ \widetilde{M} & 0 \end{pmatrix}, b(\lambda) = \begin{pmatrix} 0 \\ M^{-1}\widetilde{B} \end{pmatrix}, \\
c(\lambda) &= \begin{pmatrix} I & 0 \\ 0 & 0 \end{pmatrix}, \text{ and } s(\lambda) = I,
\end{aligned}
$$

where $\widetilde{M} = -M^{-1}(K^B\lambda^{-1/2} + K^C\lambda^{1/2})$. From (9), the optimal controller K admits the factorization

$$ K = k(\Lambda) = \Phi_U q(\Lambda) \Phi_Z^{-1}, \quad (11) $$

where $q(\lambda) = s^{-1}(\lambda)b^T(\lambda)p(\lambda)$, and where for all $\lambda \in \sigma$, $p(\lambda)$ is the unique self-adjoint nonnegative matrix solving the algebraic Riccati equation

$$ a^T(\lambda)p + pa(\lambda) - pb(\lambda)s^{-1}(\lambda)b^T(\lambda)p + c^T(\lambda)c(\lambda) = 0. $$

## 5.2 Approximation of the Functions of Λ

We build the approximation in two steps. Firstly, we use a rational approximation $k_R(\Lambda)$ of $k(\Lambda)$, then it is approximated by another function $k_{R,M}$ which is simple to discretize, and yields an accurate approximation. To do so, we use the Dunford-Schwartz formula, see [12], representing a function of an operator, because it involves only the operator $(\zeta I - \Lambda)^{-1}$ which may be simply and accurately approximated. Since the function $k(\Lambda)$ is not known, the spectrum $\sigma(\Lambda)$ cannot be easily determined, so we approximate $k(\lambda)$ by a highly accurate rational approximation $k_R(\Lambda)$, then the Dunford-Schwartz formula is applied to $k_R(\Lambda)$ with a path tracing out ellipses including $I_\sigma$ but no poles. Since the interval $I_\sigma$ is bounded, for each function $k_{ij}(\lambda)$ have a rational approximation over $I_\sigma$, we write under a global formulation, (which may be understood component wise)

$$ k_R(\lambda) = \frac{\sum_{m=0}^{R^N} d_m \lambda^m}{\sum_{m'=0}^{R^D} d'_{m'} \lambda^{m'}}, \quad (12) $$

where $d_m, d'_{m'}$ are matrices of coefficients and $R = (R^N, R^D)$ is the couple comprised of the matrices $R^N$ of numerator polynomial degrees and the matrices $R^D$ of denominator polynomial degrees. The path $C$, in the Dunford-Schwartz formula,

$$ k_R(\Lambda) = \frac{1}{2i\pi} \int_C k_R(\zeta) (\zeta I - \Lambda)^{-1} d\zeta, $$

is chosen to be an ellipse parameterized by $\zeta(\theta) = \zeta_1(\theta) + i\zeta_2(\theta)$, with $\theta \in [0, 2\pi]$. The parametrization is used as a change of variable, so the integral can be approximated by a quadrature formula involving $M$ nodes $(\theta_l)_{l=1,...,M} \in [0, 2\pi]$, and $M$ weights $(\omega_l)_{l=1,...,M}$, $I_M(g) = \sum_{l=1}^{M} g(\theta_l) \omega_l$.

In the following equations, we state that the matrices $k_R(\zeta)$ associated to the rational approximation with the numerator polynomial degrees $R^N$ and the denominator polynomial degrees $R^D$. So, for each $z \in L^2(\Gamma)^{2N}$ and $\zeta \in C$, we introduce the 2N-dimensional vector field

$$ v^{\zeta} = -i\zeta'k_R(\zeta)(\zeta I - \Lambda)^{-1}z. $$

Decomposing $v^\zeta$ into its real part $v_1^\zeta$ and its imaginary part $v_2^\zeta$, the couple $(v_1^\zeta, v_2^\zeta)$ is solution of the system

$$ \left\{ \begin{array}{l} \zeta_1 v_1^\zeta - \zeta_2 v_2^\zeta - \Lambda v_1^\zeta = \operatorname{Re}(-i\zeta' k_R(\zeta)) z, \\ \zeta_2 v_1^\zeta + \zeta_1 v_2^\zeta - \Lambda v_2^\zeta = \operatorname{Im}(-i\zeta' k_R(\zeta)) z. \end{array} \right. \quad (13) $$

Thus, combining the rational approximation $k_R$ and the quadrature formula yields an approximate realization $k_{R,M}(\Lambda)$ of $k(\Lambda)$,

$$ k_{R,M}(\Lambda) z = \frac{1}{2\pi} \sum_{l=1}^{M} v_1^{\zeta(\theta_l)} \omega_l. \quad (14) $$

This formula is central in the method, so it is the center of our attention in the simulations. A fundamental remark is that, a "real-time" realization, $k_{R,M}(\Lambda)$ $z$, requires solving $M$ systems like (13) corresponding to the $M$ quadrature nodes $\zeta(\theta_l)$. The matrices $k_R(\zeta(\theta_l))$ could be computed "off-line" once and for all, and stored in memory, so their determination would not penalize a rapid real-time computation. In total, the ultimate parameter responsible of accuracy in a real-time computation, apart from spatial discretization discussed in next Section, is $M$ the number of quadrature points.

## 6. Circuit Implementation

To realize an optimal control by a set of distributed circuits, we introduce a spatial discretization and synthesis of Equation (13). The interval $\Gamma$ is meshed with regularly spaced nodes separated by a distance $h$, we introduce $\Lambda_h^{-1}$ the finite difference discretization of $\Lambda^{-1}$, associated with the clamping boundary condition. In practice, the discretization length $h$ is chosen small compared to the distance between cantilevers. Then, $z_h$ denoting the vector of nodal values of $z$, for each $\zeta$ we introduce $(v_{1,h}^\zeta, v_{2,h}^\zeta)$, a discrete approximation of $(v_1^\zeta, v_2^\zeta)$, solution of the discrete set of equations,

$$ \zeta_1 v_{1,h}^\zeta - \zeta_2 v_{2,h}^\zeta - \Lambda_h v_{1,h}^\zeta = \operatorname{Re}(-i\zeta' k_R(\zeta)) z_h, \quad (15) $$

$$ \zeta_2 v_{1,h}^\zeta + \zeta_1 v_{2,h}^\zeta - \Lambda_h v_{2,h}^\zeta = \operatorname{Im}(-i\zeta' k_R(\zeta)) z_h. \quad (16) $$

Finally, an approximate optimal control, intended to be implemented in a set of spatially distributed actuators, could be estimated from the nodal values,

$$ k_{R,M,h} z_h = \frac{1}{2\pi} \sum_{l=1}^{M} v_{1,h}^{\zeta_l} \omega_l, $$
<figure><img src="image_3.png" /><figcaption>Fig. 2: Analog computation of $\Lambda_h v_1$.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 3: Five adjacent interior cells.</figcaption></figure>

estimated at mesh nodes in the following. We propose a synthesis of (15–16) by a distributed electronic circuit. The system is rewritten under the manageable form (17–18) and for the sake of simplicity, we use the notations $\alpha = Re(-i\zeta'k_R(\zeta))z_h$, $\beta = Im(-i\zeta'k_R(\zeta))z_h$, $v_1 = v_{1,h}^{\zeta}$, and $v_2 = v_{2,h}^{\zeta}$.

$$v_1 = \frac{\zeta_1}{\zeta_1^2 + \zeta_2^2} (\alpha + \Lambda_h v_1) + \frac{\zeta_2}{\zeta_1^2 + \zeta_2^2} (\beta + \Lambda_h v_2), \quad (17)$$

$$v_2 = \frac{\zeta_1}{\zeta_1^2 + \zeta_2^2} (\beta + \Lambda_h v_2) - \frac{\zeta_2}{\zeta_1^2 + \zeta_2^2} (\alpha + \Lambda_h v_1). \quad (18)$$

## 6.1 Analog computation of $\Lambda_h v_1$ and $\Lambda_h v_2$

The analog computation of $\Lambda_h v_1$ and $\Lambda_h v_2$ are made by Periodic Network of Resistances(PNR) circuits [11]. These electronic circuits have been developed to solve a large class of PDEs by analog computation. More exactly, PNR circuits compute the finite difference solution of a PDE.

PNR circuits are gathering of cells (Figure 2), the interior cells are indexed by $k = 1, \dots, N-1$, while the boundary cells correspond to $k = -1, 0, N$ and $N+1$. We will show that the circuits solve the equations $\Lambda_h^{-1}u_1 = i_1$. If the current sources $i_1$ are replaced by a voltage controlled current sources defined by $i_1 = gv_1$ (with $g$ is a real number), the voltage outputs of the circuits $u_1$ solve $g(\Lambda_h v_1)$ and so $\Lambda_h v_1$. The computation of $\Lambda_h v_2$ is done in the same way.

The interior cell $k$ which compute $(\Lambda_h v_1)_k$ is represented on Figure 3 with its two adjacent cells on each side.

We call $\rho_1$ the resistances between the potentials $u_1^{(k)}$ and $u_1^{(k\pm2)}$, and $\rho_2$ the resistances between the potentials $u_1^{(k)}$ and $u_1^{(k\pm1)}$. By applying the Kirchhoff Current Law (KCL) at node $u_1^{(k)}$, rearranging some terms and dividing by $h^4$, the equation of the cell $k$ can be written under the form:

$$\frac{1}{h^4} \left[ -\frac{1}{\rho_1} u_1^{(k-2)} - \frac{1}{\rho_2} u_1^{(k-1)} + 2u_1^{(k)} \left( \frac{1}{\rho_1} + \frac{1}{\rho_2} \right) - \frac{1}{\rho_2} u_1^{(k+1)} - \frac{1}{\rho_1} u_1^{(k+2)} \right] = \frac{1}{h^4} i_1^{(k)}.$$

If one choose the negative potential $\rho_1 = -h^4\rho_0$ and positive potential $\rho_2 = h^4\rho_0/4$, then the potential at node $u_1^{(k)}$ is expressed as a function of its neighbor voltages as

$$\frac{1}{h^4} [u_1^{(k-2)} - 4u_1^{(k-1)} + 6u_1^{(k)} - 4u_1^{(k+1)} + u_1^{(k+2)}] = \rho_0 i_1^{(k)},$$

which is the stencil of the differential operation $\Lambda^{-1}$. Consequently, the whole electronic circuit composed of $N-1$ cells computes the finite difference approximation $u_1 = \Lambda_h i_1 = g(\Lambda_h v_1)$. The numerical value of $\rho_0$ only changes the magnitude of the voltages $u_1^{(k)}$. The values of the resistances inside a cell depend only on the circuit topology and are easily expressed as a function of $\rho_1$ or $\rho_2$, here the resistances of the cells can be taken as $r_1 = r_3 = r_4 = r_6 = \rho_1/4$ and $r_2 = r_5 = \rho_2/2$.

The VCCS (Voltage Controlled Current Source) $i_1^{(k)}$ of Figure 3 is controlled by the voltage $v_1^{(k)}$ through the equation $i_1^{(k)} = gv_1^{(k)}$. The four boundary cells are represented in Figure 4. The imposed values of the voltages correspond
<figure><img src="image_5.png" /><figcaption>Fig. 4: Four boundary cells.</figcaption></figure>

to the clamping boundary condition. Remark that the terminals denoted by a cross are not connected, so the resistances are linked by one side at them can be removed without changing the behavior of the circuits. They are saved to show clearly the real difference between interior cells and boundary cells.

## 6.2 Analog computation of equation (17)

The analog computation of Equation (17) can be made by an array of classical non inverting summing amplifiers of Figure 5. Notice that there is no current exchange between these circuits and PNR inputs and outputs, see buffers in Figure 3. Analysis of the circuit of Figure 5 leads to (19).

<figure><img src="image_6.png" /><figcaption>Fig. 5: Analog computation of the k-th equation (17).</figcaption></figure>

With a proper choice of resistances, Figure 5 solve (17),

$$v_1^{(k)} = \frac{R_1 + R_2}{R_1} \left[ \frac{R_u}{R_b} \alpha + \frac{R_u}{R_b} g(\Lambda_h v_1)_k + \frac{R_u}{R_c} \beta + \frac{R_u}{R_d} g(\Lambda_h v_2)_k \right], \quad (19)$$

where $\frac{1}{R_u} = \frac{1}{R_a} + \frac{1}{R_b} + \frac{1}{R_c} + \frac{1}{R_d}$.

## 6.3 Analog computation of equation (18)

In a very similar way, the analog computation of Equation 18 can be made by an array of classical difference summing amplifiers of Figure 6. Analysis of the circuit of Figure 6 leads to (20). With a proper choice of resistances, Figure 6 solve (18),

$$v_2^{(k)} = \frac{R_u}{R_w} \frac{R'_1}{R'_a} \beta + \frac{R_u}{R_w} \frac{R'_2}{R'_b} g(\Lambda_h v_2)_k - \frac{R'_2}{R'_c} \alpha - \frac{R'_2}{R'_d} g(\Lambda_h v_1)_k, \quad (20)$$

where $\frac{1}{R_v} = \frac{1}{R'_a} + \frac{1}{R'_b} + \frac{1}{R'_1}$ and $\frac{1}{R_w} = \frac{1}{R'_c} + \frac{1}{R'_d} + \frac{1}{R'_2}$.

<figure><img src="image_7.png" /><figcaption>Fig. 6: Analog computation of the k-th equation (18).</figcaption></figure>

## 7. Numerical Simulation

In this Section, we validate the approximation method, established in Section 5, by a numerical simulation. We consider a silicon array comprised of an elastic base clamped of 10 elastic cantilevers, with base dimensions $L_B \times l_B \times h_B = 500 \mu m \times 16.7 \mu m \times 10 \mu m$, and one cantilever dimensions $L_C \times l_C \times h_C = 41.7 \mu m \times 12.5 \mu m \times 1.25 \mu m$. The model parameters of base and cantilever: the bending coefficient $R^B = 1.09 \times 10^{-5} N/m$, $R^C = 2.13 \times 10^{-4} N/m$ the mass per unit length $\rho^B = 0.0233 \text{kg/m}$, $\rho^C = 0.00291 \text{kg/m}$. In the rational approximation, the numerator polynomial degrees $R^N$ and the denominator polynomial degrees $R^D$ can be chosen sufficiently large (namely $R^N = R^D = 20$) so that the relative errors between the exact solution k and its rational approximation $k_R$, $e = \frac{\|k_R - k\|_{L^2(I_\sigma)}}{\|k\|_{L^2(I_\sigma)}}$, can be in the order of $10^{-8}$. This choice has no effect on real-time computation time.

<figure><img src="image_8.png" /><figcaption>Fig. 7: The relative error between the exact solution and the final approximation</figcaption></figure>
Numerical integrations have been performed with a
standard trapezoidal quadrature rule. The relative error,
$E = \frac{\|k_{R,M} - k\|_{L^2(I_\sigma)}}{\|k\|_{L^2(I_\sigma)}}$, between the exact control function
and final approximation are shown in Figure 7, for the num-
ber of nodes *M* varying from 5 to 20. It may be easily tuned
without changing spatial complexity associated with the fi-
nite difference discretization of $\Lambda^{-1}$.

We also present the ratio of the computation time of
solving the whole system for varying number of nodes *M*
to the reference computation time of solving the whole sys-
tem for *M* = 20, see Figure 8.

<figure><img src="image_9.png" /><figcaption>Fig. 8: The ratio of computation time</figcaption></figure>

**8. Conclusion**

In this paper, we have presented a semi-decentralized
approximation of an optimal control operator applied to a
two-scale model of microcantilever arrays. This model is
discretized in y-direction projecting on a transformed ba-
sis of Chebyshev polynomials. It has been shown that the
semi-decentralized optimal controller can be implemented
by a set of distributed electronic circuits. Numerical sim-
ulations have been carried out to validate the method and
study its performances. This method can be extended to
other optimal control theories, such as LQG or H∞.

**References**

1. B. Bamieh, F. Paganini, and M.Dahleh. Distributed control of spatially invariant systems. *IEEE Transactions on Automatic Control*, 47(7):1091–1107, 2002.

2. R. F. Curtain and H. Zwart. *An introduction to infinite-dimensional linear systems theory*. Texts in Applied Mathematics. Springer-Verlag, 1995.

3. R. D'Andrea and G. E. Dullerud. Distributed control design for spatially interconnected systems. *IEEE Trans. Automat. Control*, 48(9):1478–1495, 2003.

4. M. Kader, M. Lenczner, and Z. Mrcarica. Approximation of an optimal control law using a distributed electronic circuit: application to vibration control. 328(7):547–53, 2000.

5. M. Kader, M. Lenczner, and Z. Mrcarica. Distributed optimal control of vibrations: a high frequency approximation approach. 12(3):437–446, 2003.

6. M. Lenczner. A multiscale model for atomic force microscope array mechanical behavior. *Applied Physics Letters*, 90:091908, 2007.

7. M. Lenczner and R. C. Smith. A two-scale model for atomic force microscopes arrays in static operating regime. *Maths. and. Compt. Modelling*, 46:776–805, 2007.

8. Michel Lenczner, Emmanuel Pillet, Scott Cogan, and Hui Hui. A multiscale model of cantilever arrays and its updating. 7:125–135, 2009.

9. Michel Lenczner and Youssef Yakoubi. Semi-decentralized approximation of optimal control for partial differential equations in bounded domains. *Comptes Rendus Mécanique*, 337:245–250, 2009.

10. F. Paganini and B. Bamieh. Decentralization properties of optimal distributed controllers. 2(9):1877 – 1882, 1998.

11. Nicolas Ratier. Towards 2d electronic circuits in the spatial domain. *Proceedings of the 13th WSEAS international conference on Circuits*, pages 212–218, 2009.

12. K. Yosida. Functional analysis. Classics in Mathematics. Springer-Verlag, reprint of the sixth edition edition, 1980.