Volkert Van Der Wijk
email: v.vanderwijk@utwente.nl

Sébastien Krut
email: sebastien.krut@lirmm.fr

François Pierrot
email: francois.pierrot@lirmm.fr

Just L Herder
email: j.l.herder@utwente.nl

Generic Method for Deriving the General Shaking Force Balance Conditions of Parallel Manipulators with Application to a Redundant Planar 4-RRR Parallel Manipulator

Keywords: Shaking Force Balancing, Linear Momentum, Parallel Manipulator, Machine Vibrations

This paper proposes a generic method for deriving the general shaking force balance conditions of parallel manipulators. Instead of considering the balancing of a parallel manipulator link-by-link or leg-by-leg, the architecture is considered altogether.

The first step is to write the linear momentum of each element. The second step is to substitute the derivatives of the loop equations, by which the general force balance conditions are obtained. Subsequently specific kinematic conditions are investigated in order to find advantageous, simple balance solutions.

As an example, the method is applied to a planar 4-RRR parallel manipulator, for which the force balance conditions and solutions are discussed and illustrated for each step respectively. By including the loop equations, linear relations of the motion among mechanism elements lead to an increase of balance possibilities. For specific kinematic conditions, additional linear relations among the motion of mechanism elements may be obtained, resulting in another increase of balance possibilities. For the latter, symmetric motion is an important feature for which a 4-RRR manipulator is advantageous.

I. Introduction

When mechanisms (i.e. manipulators or robots) are running at high speeds, shaking forces and shaking moments are a prominent cause of machine vibrations. Instead of applying damping to reduce the influence of these vibrations, with dynamic balancing the mechanism is designed to exert no shaking forces and shaking moments at all. Therefore dynamically balanced mechanisms exhibit, among others, reduced wear, noise, and fatigue [START_REF] Lowen | Survey of Investigations into the Balancing of Linkages[END_REF], increased accuracy [START_REF] Foucault | On the Development of a Planar 3-DOF Reactionless Parallel Mechanism Proceedings of IDETC[END_REF], [START_REF] Raaijmakers | zoekt snelheidslimiet pakken en plaatsen op (Besi attacks the speedlimit for pick and place motion) Mechatronica nieuws[END_REF], and increased payload capacity [START_REF] Lim | Payload Capacity of Balanced Robotic Manipulators Robotica[END_REF].

Disadvantages of dynamic balancing are the often considerable increase of mass, inertia, and complexity of the mechanism [START_REF] Van Der Wijk | Comparison of Various Dynamic Balancing Principles Regarding Additional Mass and Additional Inertia[END_REF]. These are the main reasons that there are few studies on dynamic balancing of multi-degree-of-freedom (multi-DoF) parallel mechanisms, which are already complex mechanisms by themselves because of their parallel structures.

To omit the complexity of using the loop equations, most studies on parallel mechanisms are involved with balancing of each mechanism link or each manipulator leg individually, including [START_REF] Foucault | On the Development of a Planar 3-DOF Reactionless Parallel Mechanism Proceedings of IDETC[END_REF], [START_REF] Gosselin | Synthesis and Design of Reactionless Three-Degree-of-Freedom Parallel Mechanisms[END_REF], [START_REF] Wu | Design of Reactionless 3-DOF and 6-DOF Parallel Manipulators Using Parallelepiped Mechanisms[END_REF], [START_REF] Wu | On The Dynamic Balancing of Multi-DOF Parallel Mechanisms With Multiple Legs[END_REF], [START_REF] Arakalian | Design of planar 3-DOF 3-RRR reactionless parallel manipulators[END_REF], [START_REF] Van Der Wijk | Dynamic Balancing of Clavel's Delta Robot In: Kecskemthy and Mller[END_REF]. A class of low-mass force balanced mechanisms of which the mechanism is considered altogether without the need of loop equations is presented in [START_REF] Van Der Wijk | On the development of low-mass force balanced manipulators[END_REF].

When the loop equations are not considered, the number of balance solutions that can be found is limited. With loop equations the linear relations of the motion among mechanism elements are included, resulting in an increased number of balance solutions. For example deriving the force balance conditions with the Linear Independent Vector Method (LIVmethod) [START_REF] Berkof | A New Method for Completely Force Balancing Simple Linkages J. of Engineering for Industry[END_REF], which considers the loop equations, results in a general description of the force balance from which the complete set of general force balance conditions can be derived. The LIV-method is well applicable to parallel mechanisms, for which [START_REF] Berkof | A New Method for Completely Force Balancing Simple Linkages J. of Engineering for Industry[END_REF] also shows a few examples. However the derivation of the balance conditions is cumbersome and specific for each mechanism.

Shaking force balancing by using linear momentum equations and shaking moment balancing by using angular momentum equations have shown to be systematic and intuitive approaches for finding the complete set of balance conditions [START_REF] Van Der Wijk | On the development of low-mass force balanced manipulators[END_REF], [START_REF] Van Der Wijk | Dynamic balancing of a single crankdouble slider mechanism with symmetrically moving couplers[END_REF]. A mechanism is force balanced when the linear momentum is conserved, and a mechanism is moment balanced when the angular momentum is conserved for any motion of the mechanism.

This paper proposes a generic method for deriving the general (i.e. complete set of) shaking force balance conditions of parallel manipulators by considering the loop equations. In addition, it will be shown that advantageous balance solutions can be found for specific kinematic conditions, an opportunity that is usually overlooked.

The method is applied to a planar 4-RRR parallel manipulator, which, if in addition to being kinematically redundant has actuation redundancy, can become a high performance robot [START_REF] Corbel | Towards 100G with PKM. Is actuation redundancy a good solution for pick-and[END_REF], [START_REF] Corbel | Actuation Redundancy as a Way for Improving the Acceleration Capabilities of 3T and 3T-1R Pickand-place PKM Journal of Mechanisms and Robotics[END_REF]. This is the case since actuation redundancy leads to superior kinetostatic performances. The kinematic redundancy of the 4-RRR manipulator will show to be advantageous for dynamic balancing as compared with a 3-RRR manipulator. q 42 q 41 q 5 q 31 q 32 q 22 q 21 q 12 q 42 q 32 q 22 q 31 q 21 q 5

x y First the linear momentum equations without loop equations are investigated (leg-by-leg approach). Subsequently the derivatives of the loop equations are substituted, gaining the general force balance conditions. Then specific kinematic conditions are investigated and a resulting prototype manipulator is presented.

II. Force balance solutions without loop equations

Figure 1 shows the general topology of the planar 4-RRR parallel manipulator with its parameters. It has four legs i with each two links with lengths l i1 and l i2 of which the angular rotations with the x-axis are θ i1 and θ i2 respectively. The platform has width 2c 5 , height 2d 5 , and orientation θ 5 . The fixed pivots at the base are A i , the joints between two links of a single leg are B i and the joints of the legs with the platform are C i . For having a maximal force transmission to the platform, joints C i are located at the platform as indicated for which links B i C i are crossing one another. The center-of-mass (CoM) of each link is defined by the parameters p ij and q ij along and perpendicular to the link respectively, with j being the number of the link of a leg. The CoM of the platform is defined by p 5 and q 5 with respect to the center of the platform.

In order to have a shaking force balanced manipulator (i.e. to have the base be without shaking forces), the linear momentum of all moving elements with respect to the base needs to be conserved (to be constant). The positions of the link CoMs can be written in [x, y] T notation as

r i1 = a i + p i1 cθ i1 -q i1 sθ i1 b i + p i1 sθ i1 + q i1 cθ i1 r i2 = a i + l i1 cθ i1 + p i2 cθ i2 -q i2 sθ i2 b i + l i1 sθ i1 + p i2 sθ i2 + q i2 cθ i2 with A i = [a i , b i ] T .
To simplify the equations by omitting θ 5 , the platform CoM is modeled by four masses µ 5i which distribute m 5 among the joints C i . The positions of C i write

r 5i = a i + l i1 cθ i1 + l i2 cθ i2 b i + l i1 sθ i1 + l i2 sθ i2
Deriving the linear momentum equations [START_REF] Van Der Wijk | Comparison of Various Dynamic Balancing Principles Regarding Additional Mass and Additional Inertia[END_REF], [START_REF] Van Der Wijk | Dynamic balancing of a single crankdouble slider mechanism with symmetrically moving couplers[END_REF] of the manipulator then results in

P O = 4 i=1 (m i1 ṙi1 + m i2 ṙi2 + µ 5i ṙ5i ) = -λ
with

λ i11 = m i1 p i1 + m i2 l i1 + µ 5i l i1 λ i12 = m i1 q i1 λ i21 = m i2 p i2 + µ 5i l i2 λ i22 = m i2 q i2
and µ 51 + µ 52 + µ 53 + µ 54 = m 5 (2)

µ 51 -µ 52 -µ 53 + µ 54 = m 5 p 5 c 5 (3)
-µ 51 -µ 52 + µ 53 + µ 54 = m 5 q 5 d 5

Equations (2-4) define the distributed masses at C i by which the CoM of the distributed masses and their summed value remain equal to the position of the CoM and the value of the mass of the platform. Due to the redundant leg, there are four distributed masses and with three equations then one of them is independent. When µ 54 is chosen as independent parameter, the other µ 5i are calculated with

µ 51 = m 5 1 2 + p 5 2c 5 -µ 54 (5)
µ 52 = m 5 - p 5 2c 5 - q 5 2d 5 + µ 54 (6)
µ 53 = m 5 1 2 + q 5 2d 5 -µ 54 (7)
Already various balance solutions can be obtained. For Eq. ( 1) to be constant (zero) for all motion of the manipulator, all

λ ijk 1x 8x 24x 28x a) b) c) d)
Fig. 2. Topologies of Force Balance Solutions obtained without using the Loop Equations (Leg-by-Leg Approach).

need to be zero. This results in the following force balance conditions.

q i1 = 0 (8) q i2 = 0 (9) m i1 p i1 + m i2 l i1 + µ 5i l i1 = 0 (10) m i2 p i2 + µ 5i l i2 = 0 (11)
Having q ij be zero means that the link CoMs are in line with their links. Four types of balance solutions are discussed individually.

A. Solutions for dependent link COM parameters

To satisfy Eqs. [START_REF] Van Der Wijk | Dynamic Balancing of Clavel's Delta Robot In: Kecskemthy and Mller[END_REF] and [START_REF] Van Der Wijk | On the development of low-mass force balanced manipulators[END_REF], p ij can be chosen to be the dependent balance parameters which then are calculated with

p i1 = -m i2 l i1 -µ 5i l i1 m i1 , p i2 = -µ 5i l i2 m i2
p ij will be negative when the other parameters are positive. This represents the balance topology shown in Fig. 2a which can be interpreted as adding a countermass to each link. According to Eqs. (2-4), the platform CoM then needs to be within the area bounded by its joints. When the platform CoM is outside this area, one or more p i2 will have to be positive.

B. Solutions for one independent link COM parameter

To satisfy Eqs. ( 10) and [START_REF] Van Der Wijk | On the development of low-mass force balanced manipulators[END_REF], one p ij can be chosen to be independent instead of one µ 5i which then is dependent. One of eight possibilities is when p 22 is independent, for which µ 52 is determined by

µ 52 = -p 22 m 22 l 22
The other parameters except p 22 then are calculated with

p i1 = -m i2 l i1 -µ 5i l i1 m i1 , p i2 = -µ 5i l i2 m i2
The resulting balance topology is shown in Fig. 2b. According to Eqs. [START_REF] Foucault | On the Development of a Planar 3-DOF Reactionless Parallel Mechanism Proceedings of IDETC[END_REF][START_REF] Raaijmakers | zoekt snelheidslimiet pakken en plaatsen op (Besi attacks the speedlimit for pick and place motion) Mechatronica nieuws[END_REF][START_REF] Lim | Payload Capacity of Balanced Robotic Manipulators Robotica[END_REF], the combined CoM of the platform and mass µ 52 at C 2 then needs to be within the area bounded by the three joints of the independent µ 5i .

C. Solutions for two independent link COM parameters

To satisfy Eqs. ( 10) and [START_REF] Van Der Wijk | On the development of low-mass force balanced manipulators[END_REF], two p ij can be chosen to be independent instead of two µ 5i which then are dependent. One of twenty-four possibilities is when p 22 and p 42 are independent, for which µ 52 and µ 54 are determined by

µ 52 = -p 22 m 22 l 22 , µ 54 = -p 42 m 42 l 42
The other parameters except p 22 and p 42 then are calculated with

p i1 = -m i2 l i1 -µ 5i l i1 m i1 , p i2 = -µ 5i l i2 m i2
The resulting balance topology is shown in Fig. 2c. According to Eqs. (2-4), the combined CoM of the platform and the masses µ 52 and µ 54 at their joints C 2 and C 4 respectively, then needs to be on the line between the two joints of the independent µ 5i .

D. Solutions for three independent link COM parameters

To satisfy Eqs. ( 10) and ( 11), three p ij can be chosen to be independent instead of three µ 5i which then are dependent. One of twenty-eight possibilities is when p 12 , p 22 , and p 42 are independent, for which µ 51 , µ 52 , and µ 54 are determined by

µ 51 = -p 12 m 12 l 12 , µ 52 = -p 22 m 22 l 22 , µ 54 = -p 42 m 42 l 42
The other parameters then are calculated with

p i1 = -m i2 l i1 -µ 5i l i1 m i1 , p 32 = -µ 53 l 32 m 32
The resulting balance topology is shown in Fig. 2d. According to Eqs. (2-4), with only one independent µ 5i the combined CoM of the platform and the three dependent masses µ 51 , µ 52 , and µ 54 at their joints C 1 , C 2 , and C 4 respectively, needs to be exactly at the joint of the independent µ 5i , C 3 .

III. General force balance solutions by including loop equations

For the 4-RRR manipulator there are three independent loop equations which can be written with

r 51 = r 52 + 2c 5 cθ 5 sθ 5 r 51 = r 53 + 2c 5 cθ 5 sθ 5 + 2d 5 sθ 5 -cθ 5 r 51 = r 54 + 2d 5 sθ 5 -cθ 5
The time derivatives of these equations are

l
These loop equations can be substituted in the linear momentum equation Eq. ( 1) in various ways. A practical choice is to first eliminate θ 5 and θ5 by which two equations remain.

Choosing the CoMs of links l 12 and l 32 to be independent parameters, the loop equations are rewritten as

sθ
U l12 q 12 + c 5 l 42 d5U l12 p 12 + m 32 l42 l32 p 32 -l42 U l32 q 32 -c5l42 d5U l32 p 32 + c5l42 d5U (µ 51 -µ 53 ) + (µ 54 + µ 53 )l 42 = 0 m 21 q 21 + m 12 l 21 U l 12 p 12 + d 5 l 21 c 5 U l 12 q 12 + m 32 l 21 l32 q 32 -l 21 U l32 p 32 -d 5 l 21 c5U l32 q 32 + l21 U (µ 51 -µ 53 ) = 0 m 41 q 41 + m 12 -l41 U l12 p 12 + c5l41 d5U l12 q 12 + m 32 l41 l 32 q 32 + l41 U l 32 p 32 -c5l41 d 5 U l 32 q 32 + l41 U (µ 53 -µ 51 ) = 0 m 22 q 22 + m 12 l22 U l 12 p 12 + d5l22 c 5 U l 12 q 12 + m 32 l 22 l 32 q 32 -l 22 U l 32 p 32 -d 5 l 22 c 5 U l 32 q 32 + l 22 U (µ 51 -µ 53 ) = 0 m 42 q 42 + m 12 -l 42 U l 12 p 12 + c 5 l 42 d 5 U l 12 q 12 + m 32 l 42 l32 q 32 + l 42 U l32 p 32 -c 5 l 42 d5U l32 q 32 + l42 U (µ 53 -µ 51 ) = 0 (17)
These are the general force balance conditions of the planar 4-RRR parallel manipulator. In Fig. 3 a variety of force balance topologies is shown that can be obtained from these balance conditions. One important difference of these results with the results of Fig. 2 is that the link CoMs are not restricted to be along the lines through the joints (q ij can be nonzero). In addition, more mass parameters are independent. For instance when the mass of each element is known, for the solutions of Fig. 2 there are three independent mass position parameters while for the solutions of Fig. 3 there are six independent mass position parameters. For the latter generally holds that there are two independent mass position parameters per loop equation, while for the former two of them are because of the mass distribution of the platform (by which in fact a linear relation between the motion of the platform and of the links was already introduced) and one is because of the redundant leg. If for the six mass position parameters the position of the CoM of the platform is defined together with q i2 = 0 or q i1 = 0, then the solutions of Fig. 2 are obtained.

Figure 4a shows a balance topology when the CoMs of l 12 and l 32 are on the line through their joints. Figure 4b shows a fully symmetric topology. The advantage of this topology with respect to the topology of Fig. 2a is that there is more design freedom.

The general force balance conditions are valid for any condition of the kinematics. Each link and the platform can have any dimension and each pivot can have any position. In addition, the force balance conditions also hold for planar 1-RRR, 2-RRR, and 3-RRR (parallel) manipulators. Legs of the 4-RRR can be 'taken out' by simply filling out zero mass values. As long as two legs remain, any mass arrangement of the parallel manipulator is possible. For a 1-RRR mechanism, which is serial, the position of the CoM of each element is determined.

IV. Force balance solutions for specific kinematic conditions

The previous section showed that linear relations among the motion of mechanism elements enhance the balance possibilities. Since the 4-RRR manipulator has three-DoF only, the lin- q 42 q 41 q 5 q 31 q 32 q 22 q 21 q 12 q 42 q 11 q 32 q 22 q 31 q 21 A 4 q 41 l 11 q 42 q 41 q 5 q 31 q 32 q 22 q 21 q 12 q 42 q 11 q 22 q 32 q 21 q 31 a) b)

A 4 q 41 l 11 A 1 A 2 A 3 B 1 B 2 B 3 B 4 C 1 C 2 C 4 C
A 1 A 2 A 3 B 1 B 2 B 3 B 4 C 1 C 2 C 4 C
Fig. 5. Specific Kinematic Conditions for which for a Nonrotating Platform Velocity terms a) θ11 and θ42 , θ12 and θ41 , θ21 and θ32 , and θ22 and θ31 ; and b) θ11 and θ41 , θ12 and θ42 , θ21 and θ31 , and θ22 and θ32 , become linearly dependent.

ear momentum equations Eq. ( 22) could be rewritten to solely depend on three angular velocities. However, the six velocity terms are generally nonlinearly related and therefore this will not influence the general force balance conditions of Eq. ( 17). But for some specific kinematic conditions, relations among the velocity terms do become linear.

Two specific kinematic conditions are shown in Fig. 5 for which for each case four velocity terms become linearly dependent. For the configuration of Fig. 5a the velocity terms of θ11 and θ42 , θ12 and θ41 , θ21 and θ32 , and θ22 and θ31 are linearly dependent for a nonrotating platform (θ 5 = 0 and θ5 = 0 for all motion). This is the case when for the link lengths hold: l 11 = l 42 , l 12 = l 41 , l 21 = l 32 , and l 22 = l 31 and the pivots A i are located such that the following links remain parallel for all motion of the mechanism: l 11 l 42 , l 12 l 41 , l 21 l 32 , and l 22 l 31 . This implies that the line through pivots A 1 and A 4 and the line through A 2 and A 3 are parallel to the line through C 1 and C 4 and the line through C 2 and C 3 , respectively.

The force balance conditions for these specific kinematic conditions can be derived by substituting the relations

θ 11 = θ 42 θ11 = θ42 θ 12 = θ 41 θ12 = θ41 θ 21 = θ 32 θ21 = θ32 θ 22 = θ 31 θ22 = θ31 θ 5 = 0 θ5 = 0 (18)
in the linear momentum equations Eq. ( 22). The loop equations already included some of these relations. The resulting linear momentum equations Eq. ( 23) now have become dependent on four velocity terms. The conditions for which these equations

l 32 q 32 + l 41 U l 32 p 32 -c 5 l 41 d 5 U l 32 q 32 + l 41 U (µ 53 -µ 51 ) = 0 (19)
The twelve general force balance conditions of Eq. ( 17) have reduced to eight force balance conditions, which in fact are a combination of the twelve force balance conditions. Therefore the balance solutions of Fig. 3 and Fig. 4 are still valid. With respect to Eq. ( 17), there are four additional independent mass position parameters (p ij and q ij ), yielding new balance possibilities. Figures 6a and6b show possible resulting balance topologies. For these topologies the positions of the CoM of the platform and the CoMs of all links l i2 can be chosen arbitrary, while the CoMs of links l i1 are determined for balance. This means that if countermasses are used, they only need to be applied on the links that are directly pivoted to the base, which is an advantageous feature for low mass addition [START_REF] Van Der Wijk | Guidelines for Low Mass and Low Inertia Dynamic Balancing of Mechanisms and Robotics In: Kröger and Wahl[END_REF]. It is also possible to have some links l i2 contain countermasses instead of links l i1 . However, floating masses increase the complexity of the structural design of the manipulator considerably.

For the configuration of Fig. 5b the velocity terms of θ11 and θ41 , θ12 and θ42 , θ21 and θ31 , and θ22 and θ32 are linearly dependent for a nonrotating platform (θ 5 = 0 and θ5 = 0 for all motion). This is when for the link lengths hold: l 11 = l 41 , l 12 = l 42 , l 21 = l 31 , and l 22 = l 32 and the pivots A i are located such that the following links remain parallel for all motion of the mechanism: l 11 l 41 , l 12 l 42 , l 21 l 31 , and l 22 l 32 . This implies that the line through pivots A 1 and A 4 and the line through A 2 and A 3 are parallel to the line through C 1 and C 4 and the line through C 2 and C 3 , respectively.

The force balance conditions for this configuration can be obtained from the linear momentum equations Eq. ( 22) in an equivalent way as for the configuration of Fig. 5a l 12 q 12 + m 41 q 41 + m 12 -l41 U l12 p 12 + c5l41 d5U l12 q 12 + m 32 l 41 l 32 q 32 + l 41 U l 32 p 32 -c 5 l 41 d 5 U l 32 q 32 + l 41 U (µ 53 -µ 51 ) = 0 m 31 q 31 -m 32 l 31 l 32 q 32 + m 21 q 21 + m 12 l 21 U l12 p 12 + d 5 l 21 c5U l12 q 12 + m 32 l21 l32 q 32 -l21 U l32 p 32 -d5l21 c5U l32 q 32 + l21 U (µ 51 -µ 53 ) = 0

m 22 q 22 + m 12 l22 U l 12 p 12 + d5l22 c 5 U l 12 q 12 + m 32 l 22 l 32 q 32 -l 22 U l 32 p 32 -d 5 l 22 c 5 U l 32 q 32 + l 22 U (µ 51 -µ 53 ) = 0 m 42 q 42 + m 12 -l 42 U l 12 p 12 + c 5 l 42 d 5 U l 12 q 12 + m 32 l 42 l32 q 32 + l 42 U l32 p 32 -c 5 l 42 d5U l32 q 32 + l 42 U (µ 53 -µ 51 ) = 0 (21)
Also for this configuration eight force balance conditions are obtained which are a combination of the twelve conditions of Eq. ( 17). Figure 6c shows a possible balance topology in which the positions of the CoMs of the platform, of two links l i1 , and of two links l i2 can be chosen arbitrary. When applying countermasses, a countermass is needed for each set of links l 11 and l 41 , links l 12 and l 42 , links l 21 and l 31 , and links l 22 and l 32 . A disadvantage of this topology is that countermasses at links l i2 cannot be omitted.

For many high performance pick and place manipulators the platform does not need to rotate. The rotational degree of freedom, however, may be a necessity to be able to move the platform with accurate translations, compensating for tolerances and production inaccuracies.

Figure 7 shows a prototype manipulator, currently being developed and tested, which is derived from the balance topology of Fig. 6b. The configuration of Fig. 6b is in a singular position, but this is solved by having platform joint C 1 be coincident with C 2 and C 3 be coincident with C 4 . An advantage of this configuration with respect to the other configurations is that no arms need to cross one another. Four countermasses are applied near the base pivots, which are the rotational axes of direct drive motors.

The specific kinematic conditions found for the 4-RRR manipulator could be partly valid for the 3-RRR manipulator when two of the three legs are moving symmetrically. However, the third leg cannot be included and the kinematic conditions to have two legs move symmetrically will cause serious problems for the force transmission to the platform. Applying the specific kinematic conditions to a planar 2-RRR parallel manipulator results in a well balanced configuration when the two pivots at the base are coincident and the two joints at the platform are coincident, as illustrated in Fig. 8.

V. Discussion

The positions of the CoMs of the mechanism elements were described along the shortest way towards the base by which the least time dependent parameters are involved for each of them, which is common practice. It is also possible to write the position of for example the CoM of link l 12 along leg three. This, however, increases the complexity of the calculations considerably. For the linear momentum equations with the loop equations included, the choice of description will not affect the obtained results. For the linear momentum equations without the loop equations this choice will lead to balance solutions which 13th World Congress in Mechanism and Machine Science, Guanajuato, México, 19-25 June, 2011 A12-523 are less compact, less symmetric, and with increased countermass addition.

With the substitution of the loop equations care must be taken that the parameters of the platform c 5 and d 5 remain fully coupled. For instance when the loop equations Eqs. [START_REF] Berkof | A New Method for Completely Force Balancing Simple Linkages J. of Engineering for Industry[END_REF][START_REF] Van Der Wijk | Dynamic balancing of a single crankdouble slider mechanism with symmetrically moving couplers[END_REF][START_REF] Corbel | Towards 100G with PKM. Is actuation redundancy a good solution for pick-and[END_REF] are substituted in Eq. ( 1) for θ 12 , θ12 , θ 42 , and θ42 , c 5 and d 5 become partly decoupled. This reduces the complexity of the linear momentum equations, but it also restricts the balance solutions that can be derived.

Although this paper does not investigate the shaking moment balancing, it is noted that the 4-RRR manipulator has also advantageous moment balance features. If the manipulator is symmetrically arranged, for example equivalent to the prototype manipulator of Fig. 7, and links l i1 have equal inertia, links l i2 have equal inertia, and the platform is symmetric, then for the specific kinematic conditions the manipulator is completely moment balanced for motion along the lines of symmetry. Motion aside these lines will result in relatively low shaking moments.

VI. Conclusion

A generic method for deriving the general shaking force balance conditions of parallel manipulators was proposed which considers the robot architecture as a whole, rather than link by link or leg by leg. The focus of the method is to find and include linear relations among velocity terms of the linear momentum equations. Therefore the loop equations and specific kinematic conditions were considered.

The method was applied to a planar 4-RRR parallel manipulator, for which the force balance conditions and solutions were discussed and illustrated. A prototype was presented, featuring a simple balance solution for perfect force balance and advantageous moment balancing, taking advantage of kinematic redundancy, symmetric arrangements, and the intended manipulator task.

13th World Congress in Mechanism and Machine Science, Guanajuato, México, 19-25 June, 2011 A12-523

Appendix

The linear momentum equation Eq. ( 1) after substitution of the loop equations Eq. ( 15) and ( 16

After substitution of Eqs. (18), the linear momentum for the specific configuration of Fig. 4a becomes

P O =    

Fig. 1 .

 1 Fig. 1. Redundant Planar 4-RRR Parallel Manipulator with its Parameters. Joints C i are located at the Platform for having a Maximal Force Transmission to the Platform.

13thFig. 3 .Fig. 4 .

 34 Fig. 3. A Selection of Force Balance Topologies from the General Force Balance Conditions.

Fig. 6 .

 6 Fig. 6. Force Balance Topologies for a-b) the Kinematic Conditions of Fig. 5a; c) the Kinematic Conditions of Fig. 5b.

  by substituting the relations θ 11 = θ 41 θ11 = θ41 θ 12 = θ 42 θ12 = θ42 θ 21 = θ 31 θ21 = θ31 θ 22 = θ 32 θ22 = θ32 θ 5 = 0 θ5 = 0 (20)

Fig. 7 .

 7 Fig. 7. Prototype Manipulator derived from Fig. 6b with Joint C 1 Coincident with C 2 and C 3 Coincident with C 4 (Patent Pending).

2 Fig. 8 .

 28 Fig. 8. Planar 2-RRR Parallel Manipulator derived from Fig. 6b with only Two Legs and having the Base Pivots be Coincident and the Platform Joints be Coincident. When the CoM of the Platform is at the Joint, the Platform can also rotate without influencing the Force Balance.

  13th World Congress in Mechanism and Machine Science, Guanajuato, México, 19-25 June, 2011

															A12-523
														B 1
	p 21	m 21	l 21	B 2	m 32	l 22 p 22	m 22	C 3	d 5 m 5 p 5 c 5 d 5	C 4 c 5	C 1	l 12 q 12	m 12 p 12	m 11 q 11	l 11 p 11
	p 31 l 31	m 31	p 32	l 32		C 2			m 42	l 42 p 42	l 41	p 41 m 41
													B 4
			B 3

  111 sin θ 11 -λ 112 cos θ 11 λ 111 cos θ 11 -λ 112 sin θ 11 θ11 + -λ 211 sin θ 21 -λ 212 cos θ 21 λ 211 cos θ 21 -λ 212 sin θ 21 θ21 + -λ 311 sin θ 31 -λ 312 cos θ 31 λ 311 cos θ 31 -λ 312 sin θ 31 θ31 + -λ 411 sin θ 41 -λ 412 cos θ 41 λ 411 cos θ 41 -λ 412 sin θ 41 θ41 + -λ 121 sin θ 12 -λ 122 cos θ 12 λ 121 cos θ 12 -λ 122 sin θ 12 θ12 + -λ 221 sin θ 22 -λ 222 cos θ 22 λ 221 cos θ 22 -λ 222 sin θ 22 θ22 + -λ 321 sin θ 32 -λ 322 cos θ 32 λ 321 cos θ 32 -λ 322 sin θ 32 θ32 + -λ 421 sin θ 42 -λ 422 cos θ 42 λ 421 cos θ 42 -λ 422 sin θ 42 θ42

  The resulting linear momentum equations Eq. (22) after substituting these equations in Eq. (1) are presented in the Appendix. When the λ ijk in Eq. (22) are replaced with their original values, the twelve conditions for which the linear momentum equations are zero for all motion become m 11 p 11 + m 12 l 11 (1 -p 12 32 p 32 -l 41 U l 32 q 32 -c 5 l 41 d 5 U l 32 p 32 +

	12 cθ 12 sθ 32 cθ 32 with U = c 5 θ12 = -l 22 U l 12 l 11 l 12 l 42 U l 12 d 5 l 22 c 5 U l 12 θ32 = -l 31 l 32 l 22 l 32 l 42 l 32 sθ 42 sθ 11 cθ 11 cθ 22 -sθ 22 -cθ 42 sθ 42 sθ 22 cθ 22 sθ 31 cθ 31 sθ 22 cθ 22 cθ 42 θ42 + θ11 + θ22 + θ42 + θ22 + θ31 + θ22 + l 21 c 5 U l 12 l 21 U l 12 l 41 U l 12 d 5 l 21 c 5 l 41 d 5 U l 12 c 5 l 42 d 5 U l 12 -sθ 21 cθ 21 -cθ 41 sθ 41 sθ 21 cθ 21 sθ 41 cθ 41 sθ 42 cθ 42 l 21 l 32 sθ 21 cθ 21 θ21 + θ21 + θ41 + θ21 + θ41 + θ42 (15) l 41 l 32 sθ 41 cθ 41 θ41 + U l 32 -cθ 21 cθ 21 θ21 + l 22 U l 32 -cθ 22 cθ 22 θ22 + l 41 U l 32 cθ 41 -sθ 41 θ41 + l 42 U l 32 cθ 42 -sθ 42 θ42 -d 5 l 21 c 5 U l 32 sθ 21 cθ 21 θ21 -d 5 l 22 c 5 U l 32 sθ 22 cθ 22 θ22 -c 5 l 41 d 5 U l 32 sθ 41 cθ 41 θ41 -c 5 l 42 d 5 U l 32 sθ 42 cθ 42 θ42 (16) d5 + d 5 c5 . l 12 ) m 31 p 31 + m 32 l 31 (1 -p 32 l32 ) m 11 q 11 -m 12 l11 l 12 q 12 m 31 q 31 -m 32 l 31 l 32 q 32 m 21 p 21 + m 22 l 21 + m 12 -l21 U l 12 q 12 + d5l21 c 5 U l 12 p 12 + m 32 l 21 l 32 p 32 + l 21 U l 32 q 32 -d 5 l 21 c 5 U l 32 p 32 + d 5 l 21 c 5 U (µ 51 -µ 53 ) + (µ 52 + µ 53 )l 21 m 41 p 41 + m 42 l 41 + m 12 l41 U l12 q 12 + c5l41 d5U l12 p 12 + m 32 l 41 c 5 l 41 d 5 U (µ 51 -µ 53 ) + (µ 54 + µ 53 )l 41 m 22 p 22 + m 12 -l22 U l12 q 12 + d5l22 c5U l12 p 12 + m 32 l22 l 32 p 32 + l22 U l 32 q 32 -d5l22 c 5 U l 32 p 32 +	= 0 = 0 = 0 = 0 = 0 = 0

l d5l22 c 5 U (µ 51 -µ 53 ) + (µ 52 + µ 53 )l 22 = 0 m 42 p 42 + m 12 l 42

  In this case the force balance conditions result in m 11 p 11 + m 12 l 11 (1 -p12 l 12 ) + m 41 p 41 + m 42 l 41 + m 12 U (µ 51 -µ 53 ) + (µ 54 + µ 53 )l 41 = 0 m 31 p 31 + m 32 l 31 (1 -p 32 l 32 ) + m 21 p 21 + m 22 l 21 + m 12

	l 41 U l12 q 12 + c 5 l 41 d5U l12 p 12 +
	m 32	l41 l32 p 32 -l41 U l32 q 32 -c5l41 d5U l32 p 32 +
	c5l41
	d 5 -l21 U l12 q 12 + d5l21 c5U l12 p 12 +
	m 32 d5l21 c 5 U (µ 51 -µ 53 ) + (µ 52 + µ 53 )l 21 l21 l 32 p 32 + l21 U l 32 q 32 -d5l21 c 5 U l 32 p 32 +	= 0
	m 22 p 22 + m 12	-l 22 U l12 q 12 + d 5 l 22 c5U l12 p 12 +
	m 32	l22 l32 p 32 + l22 U l32 q 32 -d5l22 c5U l32 p 32 +
	d5l22 c5U (µ 51 -µ 53 ) + (µ 52 + µ 53 )l 22	= 0
	m 42 p 42 + m 12	l 42 U l 12 q 12 + c 5 l 42 d 5 U l 12 p 12 +
	m 32	l 42 l32 p 32 -l 42 U l32 q 32 -c 5 l 42 d5U l32 p 32 +

c5l42 d5U (µ 51 -µ 53 ) + (µ 54 + µ 53 )l 42 = 0 13th World Congress in Mechanism and Machine Science, Guanajuato, México, 19-25 June, 2011 A12-523 m 11 q 11 -m 12 l 11

  ) writes P O = (-λ 111 + l 11 l 12 λ 121 )sθ 11 + (-λ 112 + l 11 l 12 λ 122 )cθ 11 (λ 111 -l11 l12 λ 121 )cθ 11 + (-λ 112 + l11 l12 λ 122 )sθ 11 211 + l 21 U l 12 λ 122 -d 5 l 21 c 5 U l 12 λ 121 -l 21 l 32 λ 321 -l 21 U l 32 λ 322 + d 5 l 21 c 5 U l 32 λ 321 )sθ 211 -l 21 U l 12 λ 122 + d 5 l 21 c 5 U l 12 λ 121 + l 21 l 32 λ 321 + l 21 U l 32 λ 322 -d 5 l 21 c 5 U l 32 λ 321 )cθ 21 + (-λ 212 -l 21 U l12 λ 121 -d 5 l 21 c5U l12 λ 122 -l 21 l32 λ 322 + l 21 U l32 λ 321 + d 5 l 21 c5U l32 λ 322 )sθ 21 311 + l 31 l 32 λ 321 )sθ 31 + (-λ 312 + l 31 l 32 λ 322 )cθ (λ 311 -l 31 l32 λ 321 )cθ 31 + (-λ 312 + l 31 l32 λ 322 )sθ 31 411 -l 41 U l 12 λ 122 -c 5 l 41 d 5 U l 12 λ 121 -l 41 l 32 λ 321 + l 41 U l 32 λ 322 + c 5 l 41 d 5 U l 32 λ 321 )sθ 41 + (-λ 412 + l 41 U l12 λ 121 -c 5 l 41 d5U l12 λ 122 -l 41 l32 λ 322 -l 41 U l32 λ 321 + c 5 l 41 d5U l32 λ 322 )cθ 41 (λ 411 + l41 U l 12 λ 122 + c5l41 d 5 U l 12 λ 121 + l41 l 32 λ 321 -l41 U l 32 λ 322 -c5l41 d 5 U l 32 λ 321 )cθ 41 + (-λ 412 + l 41 U l 12 λ 121 -c 5 l 41 d 5 U l 12 λ 122 -l 41 l 32 λ 322 -l 41 U l 32 λ 321 + c 5 l 41 d 5 U l 32 λ 322 )sθ 41 221 + l 22 U l 12 λ 122 -d 5 l 22 c 5 U l 12 λ 121 -l 22 l 32 λ 321 -l 22 U l 32 λ 322 + d 5 l 22 c 5 U l 32 λ 321 )sθ 22 + (-λ 222 -l22 U l12 λ 121 -d5l22 c5U l12 λ 122 -l22 l32 λ 322 + l22 U l32 λ 321 + d5l22 c5U l32 λ 322 )cθ 22 (λ 221 -l 22 U l 12 λ 122 + d 5 l 22 c 5 U l 12 λ 121 + l 22 l 32 λ 321 + l 22 U l 32 λ 322 -d 5 l 22 c 5 U l 32 λ 321 )cθ 22 + (-λ 222 -l 22 U l12 λ 121 -d 5 l 22 c5U l12 λ 122 -l 22 l32 λ 322 + l 22 U l32 λ 321 + d 5 l 22 c5U l32 λ 322 )sθ 22 421 -l 42 U l12 λ 122 -c 5 l 42 d5U l12 λ 121 -l 42 l32 λ 321 + l 42 U l32 λ 322 + c 5 l 42 d5U l32 λ 321 )sθ 42 + (-λ 422 + l42 U l 12 λ 121 -c5l42 d 5 U l 12 λ 122 -l42 l 32 λ 322 -l42 U l 32 λ 321 + c5l42 d 5 U l 32 λ 322 )cθ 42 (λ 421 + l 42 U l 12 λ 122 + c 5 l 42 d 5 U l 12 λ 121 + l 42 l 32 λ 321 -l 42 U l 32 λ 322 -c 5 l 42 d 5 U l 32 λ 321 )cθ 42 + (-λ 422 + l42 U l12 λ 121 -c5l42 d5U l12 λ 122 -l42 l32 λ 322 -l42 U l32 λ 321 + c5l42 d5U l32 λ 322 )sθ 42

		θ11 +
	    	(-λ 21 + (-λ 212 -l21 U l12 λ 121 -d5l21 c5U l12 λ 122 -l21 l32 λ 322 + l21 U l32 λ 321 + d5l21 c5U l32 λ 322 )cθ 21 (λ     	θ21 +
		(-λ θ31 +
	    	(-λ     	θ41 +
	    	(-λ     	θ22 +
	    	(-λ     	θ42



  (-λ 111 + l 11 l 12 λ 121 -λ 421 -l 42 U l 12 λ 122 -c 5 l 42 d 5 U l 12 λ 121 -l 42 l 32 λ 321 + l 42 U l 32 λ 322 + c 5 l 42 d 5 U l 32 λ 321 )sθ 11 + (-λ 112 + l 11 l12 λ 122 -λ 422 + l 42 U l12 λ 121 -c 5 l 42 d5U l12 λ 122 -l 42 l32 λ 322 -l 42 U l32 λ 321 + c 5 l 42 d5U l32 λ 322 )cθ 11 (λ 111 -l11 l 12 λ 121 + λ 421 + l42 U l 12 λ 122 + c5l42 d 5 U l 12 λ 121 + l42 l 32 λ 321 -l42 U l 32 λ 322 -c5l42 d 5 U l 32 λ 321 )cθ 11 + (-λ 112 + l 11 l 12 λ 122 -λ 422 + l 42 U l 12 λ 121 -c 5 l 42 d 5 U l 12 λ 122 -l 42 l 32 λ 322 -l 42 U l 32 λ 321 + c 5 l 42 d 5 U l 32 λ 322 )sθ 11 211 + l 21 U l 12 λ 122 -d 5 l 21 c 5 U l 12 λ 121 -l 21 l 32 λ 321 -l 21 U l 32 λ 322 + d 5 l 21 c 5 U l 32 λ 321 )sθ 21 + (-λ 212 -l21 U l12 λ 121 -d5l21 c5U l12 λ 122 -l21 l32 λ 322 + l21 U l32 λ 321 + d5l21 c5U l32 λ 322 )cθ 21 (λ 211 -l 21 U l 12 λ 122 + d 5 l 21 c 5 U l 12 λ 121 + l 21 l 32 λ 321 + l 21 U l 32 λ 322 -d 5 l 21 c 5 U l 32 λ 321 )cθ 21 + (-λ 212 -l 21 U l12 λ 121 -d 5 l 21 c5U l12 λ 122 -l 21 l32 λ 322 + l 21 U l32 λ 321 + d 5 l 21 c5U l32 λ 322 )sθ 21 (-λ 311 + l 31 l32 λ 321 -λ 221 + l 22 U l12 λ 122 -d 5 l 22 c5U l12 λ 121 -l 22 l32 λ 321 -l 22 U l32 λ 322 + d 5 l 22 c5U l32 λ 321 )sθ 31 + (-λ 312 + l31 l 32 λ 322 -λ 222 -l22 U l 12 λ 121 -d5l22 c 5 U l 12 λ 122 -l22 l 32 λ 322 + l22 U l 32 λ 321 + d5l22 c 5 U l 32 λ 322 )cθ 31 (λ 311 -l 31 l 32 λ 321 + λ 221 -l 22 U l 12 λ 122 + d 5 l 22 c 5 U l 12 λ 121 + l 22 l 32 λ 321 + l 22 U l 32 λ 322 -d 5 l 22 c 5 U l 32 λ 321 )cθ 31 + (-λ 312 + l31 l32 λ 322 -λ 222 -l22 U l12 λ 121 -d5l22 c5U l12 λ 122 -l22 l32 λ 322 + l22 U l32 λ 321 + d5l22 c5U l32 λ 322 )sθ 31 321 )sθ 41 + (-λ 412 + l 41 U l 12 λ 121 -c 5 l 41 d 5 U l 12 λ 122 -l 41 l 32 λ 322 -l 41 U l 32 λ 321 + c 5 l 41 d 5 U l 32 λ 322 )cθ 41 (λ 411 + l 41 U l12 λ 122 + c 5 l 41 d5U l12 λ 121 + l 41 l32 λ 321 -l 41 U l32 λ 322 -c 5 l 41 d5U l32 λ 321 )cθ 41 + (-λ 412 + l41 U l 12 λ 121 -c5l41 d 5 U l 12 λ 122 -l41 l 32 λ 322 -l41 U l 32 λ 321 + c5l41 d 5 U l 32 λ 322 )sθ 41

			
			   	θ11 +
	    	(-λ     	θ21 +
			
	   		   	θ31 +
	    	(-λ 411 -l41 U l12 λ 122 -c5l41 d5U l12 λ 121 -l41 l32 λ 321 + l41 U l32 λ 322 + c5l41 d5U l32 λ     	θ41	(23)

13th World Congress in Mechanism and Machine Science, Guanajuato, México, 19-25 June, 2011

13th World Congress in Mechanism and Machine Science, Guanajuato, México, 19-25 June, 2011 A12-523


G Leugering
email: leugering@am.uni-erlangen.de

A A Novotny
email: novotny@lncc.br

G Perla Menzala

J Sokołowski
email: jan.sokolowski@iecn.u-nancy.fr

ON SHAPE OPTIMIZATION FOR AN EVOLUTION COUPLED SYSTEM

Keywords: piezoelectricity, electromechanical interaction, shape sensitivity analysis

A shape optimization problem in three spatial dimensions for an elasto-dynamic piezoelectric body coupled to an acoustic chamber is introduced. Well-posedness of the problem is established and first order necessary optimality conditions are derived in the framework of the boundary variation technique. In particular, the existence of the shape gradient for an integral shape functional is obtained, as well as its regularity, sufficient for applications e.g. in modern loudspeaker technologies. The shape gradients are given by functions supported on the moving boundaries. The paper extends results obtained by the authors in [9] where a similar problem was treated without acoustic coupling.

INTRODUCTION

Shape optimization for coupled models is an emerging field of research required for applications in modern key-technologies. In the present paper a model for interactions between elastic, piezo-electric and acoustic fields with non stationary partial differential equations is proposed and analyzed. The geometrical domain is decomposed into regions with different physical properties, and the sub-domains are coupled by means of appropriate transmission conditions for the equations under considerations. The problem is chosen in such a way, that the results can be applied for a broad class of models, with the appropriate modifications, if necessary. The configuration is viewed e.g. as a loudspeaker in an acoustic chamber. The question asked in applications concerning loudspeakers, beepers or energy harvesters is about the shape and the topology of the material components involved. See [START_REF] Wein | On the effect of self-penalization of piezoelectric composites in topology optimization[END_REF][START_REF] Wein | Topology optimization of piezoelectric actuators using the simp method[END_REF][START_REF] Wein | Topology optimization of a piezoelectric loudspeaker coupled with the acoustic domain[END_REF][START_REF] Wein | Topology optimization of a piezoelectric-mechanical actuator withsingle-and multiple-frequency excitation[END_REF] for the original engineering problem formulation along with topology optimization results based on the classical SIMP method. Indeed, a main objective is e.g. to maximize the acoustic pressure in the chamber by choosing appropriately shaped elasto-piezo-systems. However, in these articles the problem was concerned with optimizing the topology of the piezo-patches only. Moreover, only a time-harmonic solution was considered. Time dependent piezo-electric coupled systems have been investigated in the literature before, e.g. in [START_REF] Melnik | Existence and uniqueness theorems of the generalized solution for a class of non-stationary problem of coupled electroelasticity[END_REF][START_REF] Mercier | Existence, uniqueness and regularity results for piezoelectric systems[END_REF]. Multilayered piezo-actuator devices have been studied e.g. in [START_REF] Geis | Piezoelectricity in multi-layer actuators modelling and analysis in two and three dimensions[END_REF]. In [START_REF] Leugering | Shape sensitivity analysis of a quasi-electrostatic piezoelectric system in multilayered media[END_REF], the dynamic problem without acoustic coupling was first studied with respect to well-posedness and shapesensitivity analysis. See the references in [START_REF] Leugering | Shape sensitivity analysis of a quasi-electrostatic piezoelectric system in multilayered media[END_REF] for further information about the literature in this context. In this paper the same authors consider the fully coupled dynamic system involving also the acoustic chamber. For the mathematical theory concerning the evolution problems the reader may refer to e.g., [START_REF] Lasiecka | Mathematical control theory of coupled PDEs[END_REF].

In order to avoid additional difficulties with respect to geometrical singularities, and in order to have a simpler presentation of the results, we decide to use a layered system as in fig. 1.

Topological sensitivity analysis is not performed in the paper, however we can refer to the related papers which include the topological derivatives for the stationary models. The shape and topological sensitivity analysis of partial differential equations is an efficient tool in numerical solution of optimum design problems for distributed parameter systems. We are interested in shape sensitivity analysis in three spatial dimensions of the complete model of the interaction between elastic, piezoelectric and acoustic fields. There is a difference between stationary problems and evolution problems in this respect. To be more precise, the difference concerns the singular domain perturbations, e.g. the analysis of the influence of nucleation of small voids on the solutions of the mathematical model. Such an analysis can be performed in the framework compound and matched asymptotic expansions for stationary models, and it is unknown in the case of full evolution model. The asymptotic analysis is not, however, our principal activity in the present paper, we refer the reader to [START_REF] Cardone | Asymptotic analysis, polarization matrices, and topological derivatives for piezoelectric materials with small voids[END_REF] for some results in this direction for piezoelectric elastic bodies.

The techniques of boundary variations, which we employ in the paper, is the speed method. By this method, material and shape derivatives are determined for the partial differential equations involved, and the Hadamard structure theorem for shape gradients is used in order to identify the boundary density function of the shape gradient which, in turn, can be used in numerical methods for shape optimization.

This means that the first part of our analysis in the framework of shape sensitivity analysis is devoted to the so-called material derivatives of solutions to the boundary value problems in the stationary case or to the evolution initial-boundary value problems in the evolution case. The analysis which leads to the material derivatives is usually performed in the fixed domain setting by an application of the implicit function theorem. To this end transport mapping for the family of admissible domains is constructed, and by construction the mapping is a diffeomorphism in three spatial dimensions between admissible domains. We need some regularity of domains and of the mappings to assure all necessary properties of the diffeomorphism.

In optimum design of elastic structures the topological derivatives can be determined by asymptotic analysis with respect to the small parameter which governs perturbations of coefficients in a regular case or singular domain perturbations in limit cases of small voids and/or rigid inclusions. We point out that for evolution problems that case of singular domain perturbations is still out of the reach, however regular perturbations in coefficients make no additional difficulties compared with the shape sensitivity analysis.

In the paper the shape gradient (5.48) is obtained for shape functional (2.8) defined for the model introduced in Section 2.1. We need the expression of the shape gradient to be given by a function, for the purposes of numerical methods of shape optimization. Therefore, the regularity issue we adress in the paper can be described as follows. Under minimal regularity assumptions for the model and for the shape functional, determine the expressions for Eshelby tensors (5.36) and (5.37) in such a way, that the traces of tensors on moving boundaries are given by functions. Therefore, the shape gradient is given by a function and the levelset methods of shape optimization can be applied in order to solve numerically the associated optimization problems.

THE PROBLEM FORMULATION

Let us consider an open bounded domain Ω of R 3 with smooth boundary ∂Ω. We assume that Ω has the form Ω = D D 0 , where D and D 0 are open bounded domains with D 0 ⊂ D and Ω denotes the closure of Ω. In addition, let B i , with i = 0, 1, 2, 3, be open subsets with smooth boundary Γ i , such that, for j = 0, 1, 2, B j ⊂ B j+1 , with B 0 = D 0 and B 3 = D. We set

Ω P = B 1 B 0 , Ω M = B 2 B 1 and Ω A = B 3 B 2 .
In summary, as shown in figure 1, the mutually disjoints open domains Ω P , Ω M , Ω A have boundaries

∂Ω P = Γ 0 ∪ Γ 1 , ∂Ω M = Γ 1 ∪ Γ 2 and ∂Ω A = Γ 2 ∪ Γ 3 ,
respectively. We remark that the order of domains can be chosen in reverse order such that the acoustic part is inside and represents an acoustic chamber. According to our above motivation, Ω M and Ω P represent the regions where mechanical and piezoelectric devices are located, respectively, and Ω A represents the acoustic chamber.

2.1. The model. We are interested in the following system

       1 c 2 ϕ tt -∆ϕ = f in Ω A × (0, T ) w tt -divS = g in Ω M × (0, T ) u tt -divσ -divψ = = h 0 in Ω P × (0, T ) (2.1)
where the first equation describes the acoustic wave propagation, the second one is the linear elasticity system and the last coupled system represents the electromechanical interaction phenomenon. The equations are coupled at layers Γ j (j = 1, 2). In particular, ϕ is the acoustic potential scalar field, S is the mechanical stress tensor, σ is the electromechanical stress tensor and ψ the electrical displacement field. The constitutive laws describing the elastic behavior and piezoelectric effects, both in the linearised case of small mechanical deformations and electric fields, are    S(w) = Aε(w) , σ(u, q) = Cε(u) -P e(q) , ψ(u, q) = P ⊤ ε(u) + De(q) , (2.2) where w = w(x, t) and u = u(x, t) are the mechanical and electromechanical displacements, respectively, and q = q(x, t) is the electric potential. In addition, A and C are the elasticity fourth-order tensors respectively associated to the elastic and electromechanical parts, P the piezoelectric coupling third-order tensor and D the dielectric second-order tensor. As usual A, C and D satisfy the symmetry conditions A ijkl = A jikl = A klij , C ijkl = C jikl = C klij and D ij = D ji , whereas P satisfies P ijk = P jik . Furthermore, there exist nonnegative constants a 0 , c 0 and d 0 such that

A ijkl X ij X kl ≥ a 0 X 2 ij , C ijkl Y ij Y kl ≥ c 0 Y 2 ij , D ij z i z j ≥ d 0 z 2 i ,
where Einstein's summation convention is used. It is assumed for simplicity that all constitutive tensors are piecewise constant, i.e., constant in each layer. The mechanical strain tensors ε(u), ε(w) and the electric vector field e(q) are given by

ε(u) = ∇ s u := 1 2 (∇u + ∇u ⊤ ), ε(w) = ∇ s w := 1 2 (∇w + ∇w ⊤ ) and e(q) = -∇q , (2.3)
We complement the system (2.1) with the following initial conditions

   ϕ(x, 0) = ϕ 0 (x), ϕ t (x, 0) = ϕ 1 (x), w(x, 0) = w 0 (x), w t (x, 0) = w 1 (x), u(x, 0) = u 0 (x), u t (x, 0) = u 1 (x), (2.4)
and boundary conditions of the form

ψ • n = 0 u = 0 on Γ 0 × (0, T ), ∂ϕ ∂n = - 1 c ϕ t on Γ 3 × (0, T ), (2.5)
where n is the outward unit normal vector pointing toward the exterior of Ω. Finally, we consider the following transmission conditions    u = w σn = Sn q = q P on Γ 1 × (0, T ) and

w t • n = - ∂ϕ ∂n Sn = -ϕ t n on Γ 2 × (0, T ), (2.6)
where n = n (i) = -n (i-1) is the unit normal vector pointing toward the exterior of B i . We also assume the compatibility condition q P (x, 0) = q P t (x, 0) = 0.

2.2. Shape functional. We consider a shape functional of the form

J Ω (ϕ t , w) = T 0 J Ω (ϕ t , w) , (2.7)
with J Ω (ϕ t , w) defined as

J Ω (ϕ t , w) := α 1 2 Ω A (ϕ t -p ⋆ ) 2 -β Ω M (div(w)η + w • ∇η) , (2.8)
where p ⋆ is a target acoustic pressure, η is an arbitrary scalar function, α = 1β and β ∈ [0, 1]. By taking η| Γ 1 = 0 and η| Γ 2 = 1, we have

J Ω (ϕ t , w) = α 1 2 Ω A (ϕ t -p ⋆ ) 2 -β Γ 2 w • n , (2.9)
where w • n is the normal component of the mechanical displacement on the interface Γ 2 between the acoustic chamber and the mechanical device, respectively represented by Ω A and Ω M . It means that we want to maximize the mechanical displacement and the acoustic pressure by taking p ⋆ large enough.

STATE EQUATIONS

In this section the existence and the regularity of weak solutions to the model of coupled equations in multilayered domain is established by Theorem 1. The same results are valid for other coupled systems introduced in the paper, including the material and shape derivatives as well as the adjoint state equations.

Weak solutions.

In order to derive a weak formulation of the piezoelectric problem (2.1)-(2.6) we introduce the following bilinear forms

a A (ϕ, ϕ) := ∇ϕ, ∇ϕ Ω A , a M (w, w) := A∇ s w, ∇ s w Ω M , a M M (u, u) := C∇ s u, ∇ s u Ω P , a EE (q, q) := D∇q, ∇q Ω P ,
a M E (u, q) := P ⊤ ∇ s u, ∇q Ω P , a EM (q, u) := P ∇q, ∇ s u Ω P , and spaces

W A = H 1 (Ω A ), W M = [H 1 (Ω M )] 3 , W P = [H 1 (Ω P )] 3 , W E = H 1 (Ω P ), (3.1)
as well as

W = {(ϕ, w, u, q)(t) ∈ W A × W M × W P × W E : u = 0 on Γ 0 , w = u on Γ 1 and q = q P (t) on Γ 1 , for each t ∈ (0, T )} , (3.2)
W = {( ϕ, w, u, q) ∈ W A × W M × W P × W E : u = 0 on Γ 0 , w = u on Γ 1 and q = 0 on Γ 1 } . (3.3)
Then the weak formulation of (2.1)-(2.6) is obtained by multiplying the equations with test functions ( ϕ, w, u, q)∈ W(Ω), respectively, followed by integration by parts. It reads: for each t ∈ (0, T ) and any ( ϕ, w, u, q)∈ W(Ω), find the acoustic potential ϕ, the mechanical displacement w, the electromechanical displacement u and the electric potential q, with (ϕ, w, u, q) ∈ W, such that

       1 c 2 ϕ tt (t), ϕ Ω A + a A (ϕ(t), ϕ) -w t (t) • n, ϕ Γ 2 + 1 c ϕ t (t), ϕ Γ 3 + w tt (t), w Ω M + a M (w(t), w) + ϕ t (t), w • n Γ 2 + u tt (t), u Ω P + a M M (u(t), u) + a EM (q(t), u) +a EE (q(t), q) -a M E (u(t), q) = 0 . (3.4)
In order to put this into a more convenient format, we introduce the variable W := (ϕ, w, u, q) and the bilinear forms

A(W, W ) := a A (ϕ, ϕ) + a M (w, w) + a M M (u, u) + a EE (q, q) + a EM (q, u) -a M E (u, q), (3.5)
B(W, W ) := -w • n, ϕ Γ 2 + ϕ, w • n Γ 2 + 1 c ϕ, ϕ Γ 3 , (3.6)
where the symbol •, • K denotes the usual inner product for elements of functional spaces defined in a domain K. Notice that

A(W, W ) = a A (ϕ, ϕ) + a M (w, w) + a M M (u, u) + a EE (q, q), B(W, W ) = 1 c ϕ, ϕ Γ 3 .
The space W can be seen as the form-domain of A(•, •). The weak system (3.4) can be rewritten as

MW tt , W Ω + B(W t , W ) + A(W, W ) = 0, ∀W ∈ W , (3.7)
where M = diag ( 1 c 2 I, I, I, 0). Still, (3.7) is not a standard vectorial dissipative wave equation in weak form, the mass matrix-operator M is singular. Therefore, a proof of well-posedness seems to be at order. 3 × {0} and compatibility condition q P (x, 0) = q P t (x, 0) = 0, then, there exists a unique weak solution to (3.4) belonging to the class

Theorem 1. Given f ∈ L 2 (0, T ; L 2 (Ω A )), g ∈ L 2 (0, T ; [L 2 (Ω M )] 3 ), h ∈ L 2 (0, T ; [L 2 (Ω P )] 3 ), (ϕ 0 , w 0 , u 0 , 0) ∈ W, (ϕ 1 , w 1 , u 1 , 0) ∈ L 2 (Ω A ) × [L 2 (Ω M )] 3 × [L 2 (Ω P )]
ϕ ∈ L ∞ (0, T ; H 1 (Ω A )) , ϕ t ∈ L ∞ (0, T ; L 2 (Ω A )) , w ∈ L ∞ (0, T ; [H 1 (Ω M )] 3 ) , w t ∈ L ∞ (0, T ; [L 2 (Ω M )] 3 ) , u ∈ L ∞ (0, T ; [H 1 (Ω P )] 3 ) , u t ∈ L ∞ (0, T ; [L 2 (Ω P )] 3 ) , q ∈ L ∞ (0, T ; H 1 (Ω P )) . (3.8)
In addition, if we assume that q P (t) ∈ C 2 (Γ 1 ) and the initial data satisfy the compatibility conditions (ϕ 0 , w 0 , u 0 , ϕ 1 , w 1 , u 1 ) ∈ W, with

W = {(ϕ 0 , w 0 , u 0 ) ∈ H 2 (Ω A ) × [H 2 (Ω M )] 3 × [H 2 (Ω P )] 3 , (ϕ 1 , w 1 , u 1 , 0) ∈ W (Ω) : u = 0, ψ 0 • n = 0 on Γ 0 , w = u, σ 0 n = S 0 n on Γ 1 and q = 0 on Γ 1 w 1 • n = - ∂ ∂ n ϕ 0 , S 0 n = -ϕ 1 n on Γ 2 , ∂ ∂n ϕ 0 = - 1 c ϕ 1 on Γ 3 } , (3.9)
then the solution belongs to the (more regular) class

ϕ ∈ L ∞ (0, T ; H 2 (Ω A )) , ϕ t ∈ L ∞ (0, T ; H 1 (Ω A )) , ϕ tt ∈ L ∞ (0, T ; L 2 (Ω A )) , w ∈ L ∞ (0, T ; [H 2 (Ω M )] 3 ) , w t ∈ L ∞ (0, T ; [H 1 (Ω M )] 3 ) , w tt ∈ L ∞ (0, T ; [L 2 (Ω M )] 3 ) , u ∈ L ∞ (0, T ; [H 2 (Ω P )] 3 ) , u t ∈ L ∞ (0, T ; [H 1 (Ω P )] 3 ) , u tt ∈ L ∞ (0, T ; [L 2 (Ω P )] 3 ) , q ∈ L ∞ (0, T ; H 2 (Ω P )) . (3.10)
Proof. The proof of Theorem 1 is relegated to the Appendix. Remark 2. We can replace ϕ with c 2 ϕ in (3.4) which amounts to multiplying the first equation in (2.1) by c 2 . This is the form used in the sequel.

3.2.

Outlines of the shape sensitivity analysis. Theorem 1 implies the existence and the regularity of solutions to the model as well as to the systems which are obtained for material and shape derivatives as well as for the adjoint state.

• If the solution of state equation belongs to the class (3.10), then all boundary conditions for the shape derivatives are well defined. • For the shape functional under consideration the shape differentiability is achieved for the material derivatives belonging to the class (3.8). • Once the existence of the material derivatives is established for the model, the existence of shape derivatives follows from the relation (4.10). • If material derivatives belong to the class (3.10), then the shape derivatives belong to the class (3.8).

• If material derivatives belong to the class (3.8) then, in view of (4.10), the shape derivatives are given by very weak solution of the system.

Therefore, what we really need for the proof of shape differentiability of the functional, is the existence of the regular solution to the model, and the existence of sufficiently smooth material derivatives which can be used in order to obtain the shape differentiability of the functional. The adjoint state allows us to simplify the form of the shape gradient, but there is no implication of the adjoint state on the shape differentiability of the functional. The existence of material derivatives implies the existence of the shape derivatives as well as the differentiability of the shape functional by means of the Hadamard structure Theorem [START_REF] Sokołowski | Introduction to shape optimization -shape sensitivity analysis[END_REF] for the shape gradient.

SHAPE SENSITIVITY ANALYSIS

Formal derivation of the coupled equations for shape derivatives of solutions to the model under considerations leads to the shape gradient of the cost functional. In this derivation the transmission conditions on the interfaces should be taken into account, it means that the derivatives with respect to the shape parameter τ → 0 are evaluated from both sides of the interface. In our model one exterior boundary Γ 0 and one interface Γ 1 move according to the boundary perturbations rule defined by the speed velocity method. In formal derivation no attention is payed to the regularity of solutions, however we are interested in the resulting shape gradient regularity since the regularity has the important implications on the numerical methods. If the shape gradient is given by a distribution which lives on the moving boundaries or interfaces, this property should be taken into account when computing numerically the descent direction for gradient type numerical methods of shape optimization. On the other hand the levelset methods for shape optimization require the shape gradient of the cost given by a function, the shape gradient becomes the coefficient of the associated Hamilton-Jacoby equations for the levelset function.

On the other hand, the proof of the shape differentiability of the cost functional relies on the material derivatives of solutions to the model. The stability analysis of the model which results in the material derivatives is performed in the fixed domain setting. In this way the shape gradient of the continuous shape functional is precisely determined and it can be used for numerical computations.

For sake of simplicity, in this section we consider that the f , g and h in (2.1) are identically zero. We also consider that the initial conditions (2.4) are homogeneous. We observe that the only source in the system is given by q = q P (x, t) on Γ 1 × (0, T ), which satisfies the compatibility condition, namely, q P (x, 0) = q P t (x, 0) = 0.

The perturbed domain, parameterized by τ ∈ R + small enough, is denoted as

Ω τ = {x τ ∈ R 3 : x τ = x + τ V, x ∈ Ω, τ ≥ 0} , (4.1)
where V is a smooth vector field defined in Ω that represents the shape change velocity. Thus, the original domain is retrieved by setting τ = 0, that is Ω 0 ≡ Ω. In particular, we are interested in the perturbations of the boundary Γ 0 of the electromechanical device and of the interface Γ 1 between the mechanical and electromechanical devices. It means that the shape change velocity field can be defined as

V = 0 on Γ 2 ∪ Γ 3 = ∂Ω A . (4.2)
The shape functional defined in the perturbed domain reads

J Ωτ (ϕ τ,t , w τ ) = T 0 J Ωτ (ϕ τ,t , w τ ) , (4.3)
where ϕ τ = ϕ τ (x τ , t) and w τ = w τ (x τ , t), together with u τ = u τ (x τ , t) and q τ = q τ (x τ , t), are solutions of the following variational problem defined in the perturbed domain Ω τ : for each t ∈ (0, T ) and any

( ϕ, w, u, q) ∈ W(Ω τ ), find (ϕ τ , w τ , u τ , q τ ) ∈ W(Ω τ ), such that            ϕ τ,tt , ϕ Ω A + c 2 ∇ϕ τ , ∇ ϕ Ω A -c 2 w τ,t • n, ϕ Γ 2 + c ϕ τ,t , ϕ Γ 3 + w τ,tt , w Ω M + A∇ s τ w τ , ∇ s τ w Ω M + ϕ τ,t , w • n Γ 2 + u τ,tt , u Ω P + C∇ s τ u τ , ∇ s τ u Ω P + P ∇ τ q τ , ∇ s τ u Ω P = 0 , D∇ τ q τ , ∇ τ q Ω P -P ⊤ ∇ s τ u τ , ∇ τ q Ω P = 0 . (4.4)
with homogeneous initial conditions. In addition, the sets W(Ω τ ) and W(Ω τ ) are defined analogously as before.

4.1. Material derivatives of solutions. We are going to evaluate material and shape derivatives for the state system, and two formulae for the shape gradient including the distributed representation and the boundary representation. Before start, let us introduce the following notation for material derivative of a function ξ(x)

ξ(x) = d dτ ξ τ (x τ ) τ =0 . (4.5)
We assume for the sake of simplicity that the only source in the system is given by q = q P (x, t) on Γ 1 × (0, T ), which satisfies the compatibility condition q P (x, 0) = q P t (x, 0) = 0. In addition, we have the nonhomogeneous initial conditions for all functions.

For each t ∈ (0, T ) and any ( ϕ, w, u, q) ∈ W, find the acoustic potential ϕ, the mechanical displacement w, the electromechanical displacement u and the electric potential q, with (ϕ, w, u, q) ∈ W, such that

           ϕ tt , ϕ Ω A + c 2 ∇ϕ, ∇ ϕ Ω A -c 2 w t • n, ϕ Γ 2 + c ϕ t , ϕ Γ 3 + w tt , w Ω M + A∇ s w, ∇ s w Ω M + ϕ t , w • n Γ 2 + u tt , u Ω P + C∇ s u, ∇ s u Ω P + P ∇q, ∇ s u Ω P = 0 , D∇q, ∇ q Ω P -P ⊤ ∇ s u, ∇ q Ω P = 0 . (4.6)
Beside the above system, for the state equation the initial and boundary conditions are imposed, and the potential q P is prescribed on Γ 1 × (0, T ).

The material derivatives are given by the solutions to the following system

                                   φtt , ϕ Ω A + c 2 ∇ φ, ∇ ϕ Ω A -c 2 ẇt • n, ϕ Γ 2 + c φt , ϕ Γ 3 + ẇtt , w Ω M + A∇ s ẇ, ∇ s w Ω M + φt , w • n Γ 2 + utt , u Ω P + C∇ s u, ∇ s u Ω P + P ∇ q, ∇ s u Ω P = ∇w ⊤ (A∇ s w) + ∇ w ⊤ (A∇ s w), ∇V Ω M -w tt • w + A∇ s w • ∇ s w, divV Ω M + ∇u ⊤ (C∇ s u) + ∇ u ⊤ (A∇ s u) + ∇q ⊗ P ⊤ ∇ s u + ∇ u ⊤ P ∇q, ∇V Ω P -u tt • u + C∇ s u • ∇ s u + P ∇q • ∇ s u, divV Ω P ,
D∇ q, ∇ q Ω P -P ⊤ ∇ s u, ∇ q Ω P = ∇q ⊗ D∇ q + ∇ q ⊗ D∇q -∇u ⊤ P ∇ q -∇ q ⊗ P ⊤ ∇ s u, ∇V Ω P -D∇q • ∇ q -P ⊤ ∇ s u • ∇ q, divV Ω P .

(4.7)

The system becomes closed provided we complement the system (4.7) with the following initial conditions

   φ(x, 0) = ∇ϕ 0 (x) • V (x, 0), φt (x, 0) = ∇ϕ 1 (x) • V (x, 0), ẇ(x, 0) = ∇w 0 (x)V (x, 0), ẇt (x, 0) = ∇w 1 (x)V (x, 0), u(x, 0) = ∇u 0 (x)V (x, 0), ut (x, 0) = ∇u 1 (x)V (x, 0), (4.8)
and boundary conditions

ψ • n = -ψ • ṅ u = 0 on Γ 0 × (0, T ), ∂ φ ∂n = - 1 c φt on Γ 3 × (0, T ) . (4.9)
In addition, the potential ∇q P (x, t) • V (x, 0) is prescribed on Γ 1 × (0, T ) for the material derivative of the electric potential q.

Theorem 3. The material derivatives of solutions for the system (4.6) are given by (4.7) along with the initial conditions (4.8) and boundary conditions (4.9).

Proof. The proof of Theorem 3 is given in Section 5.3.

Shape derivatives of solutions.

The system of equations with the initial and boundary conditions is derived for the shape derivatives of solutions to the model. The shape derivatives lead to the shape gradient of the cost functional. By the Hadamard representation Theorem of the shape gradient, it follows that it is a distribution which lives on the moving boundary. From the point of view of numerical methods of shape optimization, it is preferable to have the shape gradient given by a function. The shape derivatives are given by solutions to the linearized equations with respect to the shape by using the speed method. The initial and boundary value problem for the linearized equations, in view of the shape functional under considerations, lead to appropriate adjoint state equations. All together the obtained system defines the regularity of the shape gradient which is expressed in terms of the shape derivatives, the adjoint state and the integrand of the shape functional. By the regularity assumptions on the data, the sufficient regularity of the shape gradient can be achieved. In fact, the regularity of the data is also required for derivation of the shape gradient using the material derivatives. Roughly speaking, the proof of shape differentiability is performed in the material derivatives framework in the fixed domain setting. However, in general, the shape gradient identification is possible with the shape derivatives.

Condition 4. In this section the normal component

v n := V • n
of the velocity vector field is nonnull only on the boundary Γ 0 and on the interface Γ 1 . It means that only Γ 0 and Γ 1 are perturbed by an action of the shape velocity field V .

We have the following relation between material and shape derivatives, since in general case the material derivative of a function ξ can be written as

ξ = ξ ′ + ∇ξ, V .
(4.10)

From relation (4.10) it follows that the shape derivatives looses the spatial regularity compared to the material derivatives. For hyperbolic problems this property should be taken into account in order to assure the regularity of shape derivatives in terms of the regularity of the data to the state equation.

The shape derivatives satisfy the homogeneous system (2.1) of the form

       ϕ ′ tt -c 2 ∆ϕ ′ = 0 in Ω A × (0, T ) w ′ tt -divS ′ = 0 in Ω M × (0, T ) u ′ tt -divσ ′ -divψ ′ = = 0 0 in Ω P × (0, T ) (4.11)
along with the homogeneous initial conditions,

   ϕ ′ (x, 0) = 0, ϕ ′ t (x, 0) = 0, w ′ (x, 0) = 0, w ′ t (x, 0) = 0, u ′ (x, 0) = 0, u ′ t (x, 0) = 0 . (4.12)
and nonhomogeneous boundary and interface conditions obtained below from (2.5) on Γ 0 and from (2.6) on Γ 1 , respectively. The constitutive relations (2.2) are in the same form for the shape derivatives, therefore are not repeated here. Boundary conditions for shape derivatives on Γ 0 . Now, we derive the boundary conditions on Γ 0 .

• The homogeneous Dirichlet boundary condition for the displacement field u = 0 leads to the homogeneous boundary condition for the material derivative, and in view of (4.10) becomes the nonhomogeneous boundary condition for the shape derivative

u ′ = - ∂u ∂n V • n = -v n ∂u ∂n on Γ 0 × (0, T ) , (4.13)
• The homogeneous Dirichlet boundary condition for the normal component of the vector field ψ written in the form ψ τ (x τ ) • n τ (x τ ) = 0 becomes the nonhomogeneous boundary condition for the normal component of the shape derivative vector field after differentiation with respect to τ ,

ψ ′ • n + v n n • Dψ • n -ψ Γ • ∇ Γ v n = 0 on Γ 0 × (0, T ) , (4.14)
where we denote by ψ Γ := ψ-(ψ•n)n the tangential component of the field ψ on the moving boundary Γ 0 × (0, T ) • The third condition in (2.5) is just repeated for ϕ ′ since the boundary Γ 3 × (0, T ) is independent of the shape parameter τ .

Boundary conditions for shape derivatives on Γ 1 . Now, we derive the transmission conditions on the interface Γ 1 .

• The transmission condition for displacement fields u = w leads to nonhomogeneous transmission condition for the shape derivatives obtained in the same way as for homogeneous Dirichlet boundary condition, actually

u ′ + v n ∂u ∂n = w ′ + v n ∂w ∂n on Γ 1 × (0, T ) , (4.15)
• In the similar way the boundary value for the shape derivative q ′ of the potential q is obtained

q ′ + v n ∂q P ∂n = 0 on Γ 1 × (0, T ) , (4.16)
• The equality of normal stresses σn = Sn on the interface Γ 1 × (0, T ) leads to the nonhomogeneous transmission conditions for normal stresses of shape derivatives σ ′ n, S ′ n,

σ ′ n -v n (h + 2κSn) + div Γ (v n σ Γ ) = S ′ n -v n (g + 2κσn) + div Γ (v n S Γ ) on Γ 1 × (0, T ) , (4.17)
where κ is the mean curvatrure of Γ 1 , σ Γ = σn -(σn • n)n is the tangential stress on Γ 1 , div Γ is the tangential divergence on Γ 1 , and S Γ = Sn -(Sn • n)n is the tangential stress on Γ 1 . Therefore, we complement the system (4.11) with the following boundary and transmission conditions

ψ ′ • n = -v n n • Dψ • n + ψ Γ • ∇ Γ v n u ′ = -v n ∂u ∂n on Γ 0 × (0, T ) , (4.18)
         u ′ + v n ∂u ∂n = w ′ + v n ∂w ∂n σ ′ n -v n (h + 2κSn) + div Γ (v n σ Γ ) = S ′ n -v n (g + 2κσn) + div Γ (v n S Γ ) q ′ = -v n ∂q P ∂n on Γ 1 × (0, T ) (4.19) w ′ t • n = - ∂ϕ ′ ∂n S ′ n = -ϕ ′ t n on Γ 2 × (0, T ) , (4.20)
∂ϕ ′ ∂n = - 1 c ϕ ′ t on Γ 3 × (0, T ) . (4.21)
where n is the outward unit normal vector pointing toward the exterior of Ω.

Theorem 5. For the shape derivatives of the solutions to the coupled model described in Section 2.1, we have:

• The shape derivatives ϕ ′ , w ′ , u ′ , q ′ , of the solutions ϕ, w, u, q, for the system (4.6) are given by (4.11), (4.12), (4.18)-(4.21) in the strong formulation. • For the regularity of the weak solutions to this system it is required that the following assumption q D ∈ L ∞ (0, T ; H 2 (Ω P )) is satisfied, which implies the regularity of the nonhomogeneous Dirichlet boundary condtion for the shape derivative q ′ ,

∂q D ∂n V • n ∈ L ∞ (0, T ; H 1/2 (Γ 1 )) . (4.22)
• According to (3.10), (4.10) and (4.22), there exist shape derivatives of the solutions to the system (3.4) with the following regularity given by weak solutions to the following system:

ϕ ′ ∈ L ∞ (0, T ; H 1 (Ω A )) , ϕ ′ t ∈ L ∞ (0, T ; L 2 (Ω A )) , w ′ ∈ L ∞ (0, T ; [H 1 (Ω M )] 3 ) , w ′ t ∈ L ∞ (0, T ; [L 2 (Ω M )] 3 ) , u ′ ∈ L ∞ (0, T ; [H 1 (Ω P )] 3 ) , u ′ t ∈ L ∞ (0, T ; [L 2 (Ω P )] 3 ) , q ′ ∈ L ∞ (0, T ; H 1 (Ω P )) ,
-equations are given by (4.11); -initial conditions are homogeneous (4.12); -boundary and and transmission conditions are given by (4.18)-(4.21).

Proof. The proof is standard, taking into account the specificity of the hyperbolic systems, the simplest case of the wave equation is covered in details e.g., by Cagnol and Zolésio [START_REF] Cagnol | Shape derivative in the wave equation with dirichlet boundary conditions[END_REF], see also Sokolowski and Zolésio [START_REF] Sokołowski | Introduction to shape optimization -shape sensitivity analysis[END_REF] as well as Delfour and Zolésio [START_REF] Delfour | Shapes and Geometries[END_REF]. Formally, the equations for the shape derivatives are derived by an application of the Reynolds' Transport Theorem to the variational formulation of the model in variable domain setting.

Then, the boundary conditions on moving boundary and moving interface are found from the results given in [START_REF] Sokołowski | Introduction to shape optimization -shape sensitivity analysis[END_REF] for the shape derivatives of the elasticity boundary value problems. The initial conditions are derived from the assumption that the initial conditions for the model are shape independent i.e., the shape derivatives of initial conditions are null.

SHAPE DIFFERENTIABILITY OF A FUNCTIONAL

5.1. Adjoint system. In order to simplify further calculations, let us introduce the adjoint states ϕ a , w a , v and p, which are solutions of the following variational system: For each t ∈ (0, T ) and any ( ϕ, w, v, p) ∈ W(Ω), find the adjoint acoustic potential ϕ a , the adjoint mechanical displacement w a , the adjoint electromechanical displacement v and the adjoint electric potential p, with (ϕ a , w a , v, p) ∈ W, such that

               ϕ a tt , ϕ Ω A + c 2 ∇ϕ a , ∇ ϕ Ω A -w a t • n, ϕ Γ 2 -c ϕ a t , ϕ Γ 3 + w a tt , w Ω M + A∇ s w a , ∇ s w Ω M + c 2 ϕ a t , w • n Γ 2 + v tt , v Ω P + C∇ s v, ∇ s v Ω P -P ∇p, ∇ s v Ω P = α ϕ tt -p ⋆ t , ϕ Ω A + β( η, div( w) Ω M + ∇η, w Ω M , D∇p, ∇ p Ω P + P ⊤ ∇ s v, ∇ p Ω P = 0, (5.1)
with the following final conditions

ϕ a (x, T ) = 0 and ϕ a t (x, T ) = -α(ϕ t (x, T ) -p ⋆ (x, T )), w a (x, T ) = w a t (x, T ) = 0, v(x, T ) = v t (x, T ) = 0. (5.2)
From the above system, we can define the adjoint mechanical stress tensor S a , the electromechanical stress tensor σ a and the adjoint electrical displacement ψ a as following    S a (w a ) = Aε(w a ), σ a (v, p) = Cε(v) + P e(p), ψ a (v, p) = -P ⊤ ε(v) + De(p).

(5.

3)

The strong system associated to the adjoint problem reads as follows

       ϕ a tt -c 2 ∆ϕ a = α(ϕ tt -p ⋆ t ) in Ω A × (0, T ) w a tt -divS a = 0 in Ω M × (0, T ) v tt -divσ a -divψ a = = 0 0 in Ω P × (0, T ) (5.4)
with final conditions given by (5.2), boundary conditions

ψ a • n = 0 v = 0 on Γ 0 × (0, T ), ∂ϕ a ∂n = 1 c ϕ a t on Γ 3 × (0, T ), (5.5)
and transmission conditions of the form v = w a (σ a -S a )n = -βηn on Γ 1 × (0, T ) and

w a t • n = -c 2 ∂ϕ a ∂n S a n = (βη -c 2 ϕ a t )n on Γ 2 × (0, T ). (5.6)
In addition, we have p(x, t) = 0 on Γ 1 × (0, T ), which naturally satisfies the compatibility condition.

Remark 6.

It is important to observe that the adjoint system is a time reversal problem, which should be solved by taking t ← Tt. In this case the boundary condition on Γ 3 becomes dissipative, namely

∂ϕ a ∂n = - 1 c ϕ a t on Γ 3 × (0, T ) , (5.7)
and for the adjoint system holds Theorem 1.

Proposition 7.

There is a unique weak solution (ϕ a , w a , v, p) satisfying the regularity (3.8) of Theorem 1 for the adjoint system (5.4), (5.2), (5.5), (5.6) 5.2. Shape derivative calculation. We are going to denote by ϕ τ,t := ∂ϕ τ ∂t the time derivative of the function ϕ τ which is defined in Ω τ .

Let us perform the shape sensitivity analysis of the functional J Ωτ (ϕ τ,t , w τ ). Thus, we need to calculate its derivative with respect to the parameter τ at τ = 0, that is

T 0 JΩ (ϕ t , w) = JΩ (ϕ t , w) := d dτ J Ωτ (ϕ τ,t , w τ ) τ =0
.

(5.8)

In order to proceed, it is convenient to introduce an analogy to classical continuum mechanics [START_REF] Gurtin | An introduction to continuum mechanics[END_REF] whereby the shape change velocity field V is identified with the classical velocity field of a deforming continuum and τ is identified as an artificial time parameter (we refer to [START_REF] Sokołowski | Introduction to shape optimization -shape sensitivity analysis[END_REF] for analogies of this type in the context of shape sensitivity analysis). In this case, by making use of the concept of material derivative of a spatial field [START_REF] Gurtin | An introduction to continuum mechanics[END_REF][START_REF] Gurtin | Configurational forces as basic concept of continuum physics[END_REF] and considering the Reynolds' Transport Theorem, the shape derivative of the functional J Ω (ϕ t , w) is given by

JΩ (ϕ t , w) = D Ω (J Ω (ϕ t , w)), V + D ϕt (J Ω (ϕ t , w)), φt + D w (J Ω (ϕ t , w)), ẇ , (5.9)
where

T 0 D Ω (J Ω (ϕ t , w)), V = β T 0 Ω M (∇w ⊤ η + ∇η ⊗ w) • ∇V -β T 0 Ω M (div(w)η + w • ∇η)divV = β T 0 ∇w ⊤ η + ∇η ⊗ w, ∇V Ω M -β T 0 div(w)η + w • ∇η, divV Ω M , (5.10)
and

T 0 D ϕt (J Ω (ϕ t , w)), φt = α T 0 Ω A (ϕ t -p ⋆ ) φt = α Ω A (ϕ t -p ⋆ ) φ T 0 -α T 0 Ω A (ϕ tt -p ⋆ t ) φ = α Ω A (ϕ t (T ) -p ⋆ (T )) φ(T ) -α T 0 Ω A (ϕ tt -p ⋆ t ) φ = α (ϕ t (T ) -p ⋆ (T )), φ(T ) Ω A -α T 0 (ϕ tt -p ⋆ t ), φ Ω A , (5.11)
T 0 D w (J Ω (ϕ t , w)), ẇ = -β T 0 Ω M (div( ẇ)η + ẇ • ∇η) = -β T 0 ( η, div( ẇ) Ω M + ∇η, ẇ Ω M ) .
(5.12)

Thus, since the acoustic chamber remains fixed, we have

JΩ (ϕ t , w) = β T 0 ∇w ⊤ η + ∇η ⊗ w, ∇V Ω M -β T 0 div(w)η + w • ∇η, divV Ω M - T 0 α (ϕ tt -p ⋆ t ), φ Ω A - T 0 β( η, div( ẇ) Ω M + ∇η, ẇ Ω M + α (ϕ t (T ) -p ⋆ (T )), φ(T ) Ω A .
(5.13) 5.3. Proof of Theorem 3. Let us now calculate the derivative of the state system (4.4) with respect to the parameter τ at τ = 0. Thus, by making use again of the concept of material derivative of a spatial field [START_REF] Gurtin | An introduction to continuum mechanics[END_REF][START_REF] Gurtin | Configurational forces as basic concept of continuum physics[END_REF] and considering the Reynolds' Transport Theorem, we obtain:

• For the acoustic chamber

ϕ tt , ϕ • Ω A = φtt , ϕ Ω A , (5.14) ∇ϕ, ∇ ϕ • Ω A = ∇ φ, ∇ ϕ Ω A , (5.15) w t • n, ϕ • Γ 2 = ẇt • n, ϕ Γ 2 , (5.16) ϕ t , ϕ • Γ 3 = φt , ϕ Γ 3 .
(5.17)

• For the mechanical device

w tt , w • Ω M = ẇtt , w Ω M + Ω M (w tt • w)divV , (5.18) ϕ t , w • n • Γ 2 = φt , w • n Γ 2 , (5.19) A∇ s w, ∇ s w • Ω M = A∇ s ẇ, ∇ s w Ω M + Ω M (A∇ s w • ∇ s w)divV - Ω M (∇w ⊤ (A∇ s w) + ∇ w ⊤ (A∇ s w)) • ∇V .
(5.20)

• For the piezoelectric device D∇q, ∇ q • Ω P = D∇ q, ∇ q Ω P + Ω P (D∇q • ∇ q)divV -Ω P (∇q ⊗ D∇ q + ∇ q ⊗ D∇q) • ∇V , (5.24) P ⊤ ∇ s u, ∇ q • Ω P = P ⊤ ∇ s u, ∇ q Ω P + Ω P (P ⊤ ∇ s u • ∇ q)divV -Ω P (∇u ⊤ P ∇ q + ∇ q ⊗ P ⊤ ∇ s u) • ∇V .

u tt , u • Ω P = utt , u Ω P + Ω P (u tt • u)divV , (5.21
(5.25)

where we have used the fact that the admissible variations ϕ, w, u and q do not depend on the parameter τ . Thus, the derivative with respect to the shape parameter τ of the state system, after some rearrangements, becomes (4.7).

Distributed Shape Gradient.

Theorem 8. The form of distributed gradient of shape functional (4.3) defined in variable domain setting, is given by (5.35), (5.36), (5.37). In addition, for the strong solutions we have the divergence free Eshelby tensors (5.46).

Proof. By setting in (4.7) ϕ = ϕ a , w = w a , u = v and q = p in the above, we obtain

                                               ϕ a tt , φ Ω A + c 2 ∇ φ, ∇ϕ a Ω A -c 2 ẇt
• n, ϕ a Γ 2 + c φt , ϕ a Γ 3 + w a tt , ẇ Ω M + A∇ s ẇ, ∇ s w a Ω M + φt , w a • n Γ 2 + v tt , u Ω P + C∇ s u, ∇ s v Ω P + P ∇ q, ∇ s v Ω P = ∇w ⊤ (A∇ s w a ) + (∇w a ) ⊤ (A∇ s w), ∇V Ω M w tt • w a + A∇ s w • ∇ s w a , divV Ω M + ∇u ⊤ (C∇ s v) + ∇v ⊤ (A∇ s u) + ∇q ⊗ P ⊤ ∇ s v + ∇v ⊤ P ∇q, ∇V Ω P u tt • v + C∇ s u • ∇ s v + P ∇q • ∇ s v, divV Ω P + ϕ a tt , φ Ω A -φtt , ϕ a

Ω A + w a tt , ẇ Ω M -ẇtt , w a Ω M
+ v tt , u Ω P -utt , v Ω P , D∇ q, ∇p Ω P -P ⊤ ∇ s u, ∇p Ω P = ∇q ⊗ D∇p + ∇p ⊗ D∇q -∇u ⊤ P ∇p -∇p ⊗ P ⊤ ∇ s u, ∇V Ω P -D∇q • ∇p -P ⊤ ∇ s u • ∇p, divV Ω P .

(5.26)

where we have introduced the terms ± ϕ a tt , φ Ω A , ± w a tt , ẇ Ω M , ± v tt , u Ω P in the left hand side of the first equality. Using integration by parts, we have (5.32) Thus, identity (5.26) can be re-written as

T 0 ϕ a tt , φ Ω A - T 0 φtt , ϕ a Ω A = ϕ a t , φ Ω A | T 0 -φt , ϕ a Ω A | T 0 = ϕ a t (T ), φ(T ) Ω A = -α(ϕ t (T ) -p ⋆ (T )), φ(T ) Ω A , (5.27)
                                           ϕ a tt , φ Ω A + c 2 ∇ϕ a , ∇ φ Ω A -w a t • n, φ Γ 2 -c ϕ a t , φ Γ 3 + w a
tt , ẇ Ω M + A∇ s w a , ∇ s ẇ Ω M + c 2 ϕ a t , ẇ • n Γ 2 + v tt , u Ω P + C∇ s v, ∇ s u Ω P -P ∇p, ∇ s u Ω P = ∇w ⊤ (A∇ s w a ) + (∇w a ) ⊤ (A∇ s w), ∇V Ω M w tt • w a + A∇ s w • ∇ s w a , divV Ω M + ∇u ⊤ (C∇ s v) + ∇v ⊤ (A∇ s u) + ∇q ⊗ P ⊤ ∇ s v + ∇v ⊤ P ∇q, ∇V Ω P u tt • v + C∇ s u • ∇ s v + P ∇q • ∇ s v, divV Ω P -∇ s v, P ∇ q Ω P -P ∇p, ∇ s u Ω P , D∇p, ∇ q Ω P + P ⊤ ∇ s v, ∇ q Ω P = ∇q ⊗ D∇p + ∇p ⊗ D∇q -∇u ⊤ P ∇p -∇p ⊗ P ⊤ ∇ s u, ∇V Ω P -D∇q • ∇p -P ⊤ ∇ s u • ∇p, divV Ω P + ∇p, P ⊤ ∇ s u Ω P + P ⊤ ∇ s v, ∇ q Ω P ,

(5.33)

FIGURE 1 .

 1 FIGURE 1. Layered domain represented by Ω.

)

  C∇ s u, ∇ s u • Ω P = C∇ s u, ∇ s u Ω P + Ω P (C∇ s u • ∇ s u)divV -Ω P (∇u ⊤ (C∇ s u) + ∇ u ⊤ (A∇ s u)) • ∇V ,(5.22)P ∇q, ∇ s u • Ω P = P ∇ q, ∇ s u Ω P +Ω P (P ∇q • ∇ s u)divV -Ω P (∇q ⊗ P ⊤ ∇ s u + ∇ u ⊤ P ∇q) • ∇V ,(5.23)

  a • n Γ 2 = -T 0 w a t • n, φ Γ 2 .

We then subtract weak-(⋆) convergent subsequences and pass to the limit in the equations. The fulfillment of the initial data is proved by a standard argument. Note that also non-homogenous boundary conditions for q (and u ) can be easily handled.

Author G.L. acknowledges support by the DFG-Cluster of Excellence: Engineering of Advanced Materials (grant). Author G.P.M acknowledges support by the Brazilian Agency CNPq under grant 301134/2009-0. Authors A.A.N. and J.S. acknowledge support by the Brazilian-French research programme CAPES/COFECUB under grant 604/08 beetwen LNCC (Petrópolis) and IECN (Nancy). The research of J.S. was partially supported by the Projet ANR GAOS Geometrical analysis of optimal shapes.

where we have introduced the terms ± P ∇p, ∇ s u Ω P and ± P ⊤ ∇ s v, ∇ q Ω P . In the same way, let us set ϕ = φ, w = ẇ, v = u and p = q in the adjoint system (5.1), to obtain

D∇p, ∇ q Ω P + P ⊤ ∇ s v, ∇ q Ω P = 0 , (5.34)

By comparing (5.33) with (5.34), and using (5.13) we observe that

where the last term of (5.13) is absorbed by (5.27) and we have used the fact that ∇ s v, P ∇ q Ω P = P ⊤ ∇ s v, ∇ q Ω P and P ∇p, ∇ s u Ω P = ∇p, P ⊤ ∇ s u Ω P . In addition, the Eshelby tensors [START_REF] Eshelby | The elastic energy-momentum tensor[END_REF] Σ M and Σ P are respectively given by

(5.36)

with σ, ψ and σ a , ψ a given, respectively, by (2.2) and (5.3).

Boundary Shape Gradient.

Theorem 9. By the structure theorem for a shape differentiable shape functionals [START_REF] Sokołowski | Introduction to shape optimization -shape sensitivity analysis[END_REF], from (5.35) the boundary formulae of the shape gradient is obtained. In general, the shape gradient on the boundary is given by a distribution. However, for the strong solutions, in view of (5.46), the boundary formula for the shape gradient takes the form (5.47), and in such a case the shape gradient on the moving boundary is given by a function.

Proof. After applying the divergence theorem in (5.35), we observe that

(5.38) 1) is the unit normal vector pointing toward the exterior of B i . Let us calculate the divergence of the tensors Σ M and Σ P given by (5.36) and (5.37), respectively

since ∇(div(w)η + w • ∇η) = div(∇w ⊤ η + ∇η ⊗ w). Integration by parts yields

and after some arrangements, we obtain divΣ M = (∇w) ⊤ (w a tt -divS a ) + (∇w a ) ⊤ (w tt -divS) , (5.44) divΣ P = (∇u) ⊤ (v tt -divσ a ) + (∇v) ⊤ (u tt -divσ) + ∇qdivψ a + ∇pdivψ .

(5.45)

Finally, by taking into account the strong systems (2.1) and (5.4), we have the following important results

divΣ M = divΣ P = 0 .

(5.46)

In addition, since V = 0 on Γ 2 , and from these last results together with (5.38,5.39), we obtain the final expression for the shape derivative of the functional JΩ (ϕ t , w), namely

(5.47) with Σ M and Σ P given respectively by (5.36) and (5.37). The above form of shape derivative of the distributed functional can serve us to identify the shape gradient.

Since the shape functional in question is differentiable in the sense of the shape sensitivity analysis in [START_REF] Sokołowski | Introduction to shape optimization -shape sensitivity analysis[END_REF], we can apply the structure theorem to this end. In particular, from the boundary and transmission conditions, namely, (2.5), (5.5) and (2.6), (5.6), respectively, it is straightforward to verify that the above equation holds the structure theorem. Therefore, it is sufficient to take into consideration the speed vector fields normal to the boundaries and the interfaces. This observation influences only two boundary integrals with the Eshelby tensor, and the result is the following.

Corollary 10. The density g of the boundary shape gradient of the distributed shape functional is given by the following expression

As it is indicated before, in order to apply the level-set strategy of shape optimization, it is required that the density g of the boundary shape gradient is given by functions supported on the boundaries and on the interfaces.

CONCLUSIONS

In this paper the shape optimization problem for coupled non stationary partial differential equations is analysed. Beside the existence of an optimal shape under realistic conditions, the form of the shape gradient is established in usual expressions necessary for applications of numerical methods, say, for boundary integrals. This means, that the shape optimization problem can be solved by the discretization of the continuous shape gradient and the appropriate finite elements in spatial variables and the finite differences in time variable, for example. The numerical realization, however, will be subject to a forthcoming publication.

APPENDIX A. PROOF OF THEOREM 1

Without lost of generality, in the proof we can assume that the boundary condition for q (i.e. q P ) as well as the nonhomogeneous terms f , g and h are identically equal to zero. We can use the Galerkin procedure. Thus, we introduce sequences {( ϕ m , w m , u m , q m ) m∈N } in W and the finite dimensional spaces W m = span{( ϕ 1 , w 1 , u 1 , q 1 ), ..., ( ϕ m , w m , u m , q m )} so that the union over all such spaces is dense in W. Clearly, if we take the test functions (ϕ m , w m , u m , q m ), with ϕ m = ϕ = ϕ, w m = w = w, u m = u = u and q m = q = q, and initial conditions ϕ m (0

{0} respectively, then, it follows that problem (3.4) has a local solution in an interval [0, t m ). In order to extend the solution to [0, +∞), the finite dimensional system of ordinary differential equations associated to (3.4) now reads as follows. For any ( ϕ, w, u, q

Taking as test functions (ϕ m t (t), w m t (t), 0, 0) in (A.1) and integrating over [0, t] we obtain

Now as test functions in (A.1) we take (0, 0, u m t (t), 0), (0, 0, 0, q m t (t)) and (0, 0, 0, q m (t)) and to obtain the identities 1 2

and a EE (q m , q m t ) = a M E (u m , q m t ) and a EE (q m , q m ) = a M E (u m , q m ) , (A.4) respectively. Then, it follows that

by symmetry. Using (A.4) in (A.3) and integrating over [0, t] we obtain

In order to obtain an initial condition for q m we need to solve

for any ξ ∈ {span{ξ 1 , ξ 2 , ..., ξ m }, ξ i ∈ H 1 (Ω) : ξ i = 0 on Γ 1 }. Since we know the regularity of u m 0 we can apply the Lax-Milgran lemma and obtain a unique solution q m (0) belonging to H 1 (Ω) with q m (x, 0) = 0 on Γ 1 . Furthermore q m (0) ≤ c u m 0 [H 1 (Ω P )] 3 . (A.9) Using the coercivity of the bilinear forms a M M and a EE in W P and W E respectively to obtain from (A.3) and (A.7). In case we consider f , g, h and q D different from zero, we use Gronwall's inequality at this point.

for some positive constant C. A standard argument shows that also

Using the a priori energy estimates (A.10) and (A.11) we can then extract subsequences {ϕ m }, {ϕ m t }, {ϕ m tt }; {w m }, {w m t }, {w m tt }; {u m }, {u m t }, {u m tt }, which we relabel by original indices converging for K := A, M, P weak-(⋆) in L ∞ (0, T ; W K (Ω)), and weak in L 2 (0, T ; W K (Ω) * ), respectively, to elements ϕ * , ϕ * t , ϕ * tt ; w * , w * t , w * tt ; u * , u * t , u * tt . Standard arguments reveal that these elements solve the weak system (3.4) and that the initial data are matched in the corresponding spaces as well.

As for the second part of the theorem, we first differentiate the approximate weak system and take (ϕ m tt (t), 0, 0, 0) and then (0, w m tt (t), 0, 0). One adds the results and integrates with respect to time to obtain

As for the piezoelectric part, we take (0, 0, u m tt (t), 0) as test functions and then, after another differentiation of the weak form with respect to t, (0, 0, 0, q m t (t). We obtain 1 2

Integration with respect to time leads to:

We need estimates on ϕ m tt (0) , w m tt (0) [L 2 (Ω M )] 3 , u m tt (0) [L 2 (Ω P )] 3 and q m t (0) W in terms of our data. As now u m t (0) ∈ W P we can uniquely solve the second equation of (3.4) to obtain q m t (0) W E ≤ C{ u m t (0) W P }. Moreover, for given u m (0) ∈ H 2 (Ω P ) let q m (0) ∈ H 2 (Ω P ) be such that    div(D∇q) = div(P ⊤ ∇ s u m (0)) in Ω P D∇q • n = P ⊤ ∇ s u • n on Γ 1 q = 0 on Γ 0 (A.15)

Then q m (0) H 2 (Ω) P ≤ { u m (0) H 2 (Ω P ) 3 } and after evaluating the strong solution at t = 0 we obtain u m tt (0) ≤ C u m (0) H 2 (Ω P ) 3 We can now proceed as before, in order to obtain the a priori estimate u m tt (t)