# Contrôle de la puissance de chauffe du procédé PFR à l'aide d'abaques dynamiques

Nicolas Bur, Pierre Joyot, Pierre Villon, Francisco Chinesta

► To cite this version:

Nicolas Bur, Pierre Joyot, Pierre Villon, Francisco Chinesta. Contrôle de la puissance de chauffe du procédé PFR à l'aide d'abaques dynamiques. 12e Colloque national en calcul des structures, CSMA, May 2015, Giens, France. hal-01515061

HAL Id: hal-01515061

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

Submitted on 27 Apr 2017

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 CC0 1.0 - Universal - International License
# Contrôle de la puissance de chauffe du procédé PFR à l'aide d'abaques dynamiques

N. Bur<fnref n="1" />, P. Joyot<fnref n="1" />, P. Villon<fnref n="2" />, F. Chinesta<fnref n="3" />

**Résumé** — Afin d'améliorer le procédé de Placement de Fibre Robotisé (PFR), la puissance de la source de chauffe utilisée doit être contrôlée précisément. Nous proposons de construire des abaques liant, entre autres, la puissance et l'accélération du robot. De la sorte nous établissons une loi de chauffe à partir d'un profil de vitesse. Nous validons notre approche en résolvant l'équation de la chaleur à partir de la vitesse donnée et de la loi de chauffe précédemment calculée.

**Mots clés** — Composites, Placement de Fibre Robotisé, Simulation numérique, Réduction de modèle, PGD, Contrôle dynamique.

## 1 Introduction

Parmi les procédés de fabrication de structures en composites thermoplastiques, le Placement de Fibre Robotisé (PFR) connaît un essor considérable depuis les années 80 [5]. Le principe de fonctionnement schématisé à la Figure 1 est le suivant : Les fibres, qui se présentent sous forme de bandelettes, sont acheminées jusqu'à la tête de dépose fixée en bout d'un bras robot. Elles sont guidées par la tête et déposées sur le moule. Un flux laser, porté par fibre optique et focalisé à l'endroit idoine, permet de fondre la matrice plastique qui enrobe les fibres ; simultanément, un rouleau de compactage, comme son nom l'indique, assure la cohésion entre le substrat et la fibre nouvellement déposée et fondue. Lorsque le robot parvient à l'extrémité de la pièce, les fibres sont sectionnées et le robot se repositionne pour déposer une nouvelle bande de matériau.

<figure><img src="image_2.png" /><figcaption>FIGURE 1 – Schématisation du PFR</figcaption></figure>

La modélisation numérique de ce procédé est sujette à d'importants travaux de recherche. Dans [3], un système thermo-mécanique a été résolu, fournissant non seulement le champ de température au sein de la pièce étudiée, mais également les contraintes résiduelles induites par l'histoire thermique.

Cependant, le modèle mis en place se base sur l'hypothèse du régime quasi-statique, la vitesse utilisée étant considérée constante dans le temps. Nous proposons ici d'étudier le régime transitoire afin de connaître l'évolution de la température lorsque la vitesse varie. Cette approche vise à répondre à la nécessité de pouvoir piloter la source de chaleur en fonction du déplacement du robot, afin d'atteindre en permanence la température optimale du procédé.

Nous allons donc tirer parti de la méthode PGD afin d'établir des abaques numériques solution de

<footnote n="1">¹ Estia-Recherche, {n.bur, p.joyot}@estia.fr</footnote>

<footnote n="2">² GSM, Université de Technologie de Compiègne, pierre.villon@utc.fr</footnote>

<footnote n="3">³ GeM, École Centrale Nantes, francisco.chinesta@ec-nantes.fr</footnote>
modèles multi-paramétriques de l'équation de la chaleur.

## 2 PGD et minimisation du résidu

Dans cette section, nous donnons l'expression tensorielle de la PGD, dans sa variante basée sur la minimisation du résidu (voir [1]). Pour plus de détails sur cette méthode numérique, introduite par P. Ladevèze dans [4], le lecteur intéressé pourra se référer à [2].

Soit $\Omega$ un domaine multidimensionnel. Considérons la formulation faible d'un problème linéaire donnée par

$$a(\Psi(x_1, \dots, x_N), \Psi^*(x_1, \dots, x_N)) = b(\Psi^*(x_1, \dots, x_N)), \quad (1)$$

pour laquelle nous recherchons une solution approchée écrite sous la forme continue

$$\begin{align}
\Psi(x_1, \dots, x_N) &\approx \sum_{i=1}^{M} \alpha^i \psi_1^i(x_1) \cdot \psi_2^i(x_2) \cdot \dots \cdot \psi_N^i(x_N); \nonumber \\
&= \sum_{i=1}^{M} \left( \alpha^i \prod_{j=1}^{N} \psi_j^i(x_j) \right) = \Psi^m(x_1, \dots, x_N). \tag{2}
\end{align}$$

La représentation séparée de cette approximation est construite à l'aide d'un schéma itératif d'enrichissement-projection. En utilisant les valeurs nodales de chaque fonction, l'équation précédente permet d'établir la forme discrète

$$\begin{align}
\Psi(x_1, \dots, x_N) &= \sum_{i=1}^{M} \alpha^i \psi_1^i(x_1) \otimes \psi_2^i(x_2) \otimes \dots \otimes \psi_N^i(x_N); \nonumber \\
&= \sum_{i=1}^{M} \left( \alpha^i \prod_{j=1}^{N} \psi_j^i(x_j) \right) = \Psi^m(x_1, \dots, x_N). \tag{3}
\end{align}$$

où nous avons conservé la même notation pour ne pas alourdir le texte. De même, pour simplifier, le symbole $\prod$ désigne ici le produit tensoriel quand il indiquait le produit classique dans l'équation (2).

Afin de traiter le cas d'opérateur symétrique, nous nous intéressons à la résolution de

$$\mathcal{A}^T \mathcal{A} \Psi = \mathcal{A}^T \mathcal{B}, \quad (4)$$

avec une représentation séparée de $\Psi$ comme indiqué dans (2), où $\mathcal{A}$ et $\mathcal{B}$ désignent les opérateurs gauche et droit de l'équation (1). Ils s'expriment comme des tenseurs :

$$\mathcal{A} = (A_k^i)_{\substack{1 \le k \le N \\ 1 \le i \le n_A}}, \quad \mathcal{B} = (B_k^j)_{\substack{1 \le k \le N \\ 1 \le j \le n_B}} \quad (5)$$

Partant du postulat que nous connaissons les $m$ premiers modes donnant l'approximation $\Psi^m$ de $\Psi$, nous cherchons un nouveau mode $\psi^{m+1} = \prod_{j=1}^{N} \psi_j^{m+1}$ pour enrichir cette approximation en résolvant

$$\mathcal{A}^T \mathcal{A} \Psi^m + \mathcal{A}^T \mathcal{A} \prod_{j=1}^{N} \psi_j^{m+1}(x_j) = \mathcal{A}^T \mathcal{B}. \quad (6)$$

De façon plus explicite, afin de calculer la $k^e$ fonction de ce mode $\psi^{m+1}$, on a

$$\sum_{l=1}^{n_A} \sum_{j=1}^{n_A} \sum_{i=1}^{m} \left( \alpha^i \underline{A}_k^l \underline{A}_k^j \psi_k^i \prod_{\substack{h=1 \\ h \neq k}}^{N} \psi_h^T \underline{A}_h^l \underline{A}_h^j \psi_h^i \right) + \sum_{l=1}^{n_A} \sum_{j=1}^{n_A} \left( \underline{A}_k^l \underline{A}_k^j \psi_k^{m+1} \prod_{\substack{h=1 \\ h \neq k}}^{N} \psi_h^{m+1} \underline{A}_h^l \underline{A}_h^j \psi_h^i \right) = \sum_{l=1}^{n_A} \sum_{j=1}^{n_B} \left( \underline{A}_k^l \underline{B}_k^j \prod_{\substack{h=1 \\ h \neq k}}^{N} \psi_h^{m+1} \underline{A}_h^l \underline{B}_h^j \right). \quad (7)$$

Ainsi, pour obtenir le mode $m+1$, il faut résoudre l'équation précédente dans une boucle de type point-fixe par exemple :

1. les différentes fonctions $\psi_h^{m+1}$ sont connues (par initialisation ou d'après les itérations), $\forall h \in \{[1, \dots, N], h \neq k\}$;

2. $\psi_k^{m+1}$ est calculé avec l'équation (7);

3. la liste des fonctions $\psi_h^{m+1}$ est mise à jour avec le $\psi_k^{m+1}$ tout juste calculé;

4. on retourne au premier point pour calculer la fonction $\psi_k^{m+1}$ suivante;

5. le calcul s'arrête lorsqu'une certaine précision est atteinte.
# 3 Loi de chauffe

## 3.1 Mise en place du modèle

Nous nous focalisons ici sur une phase d'accélération du robot, comprise entre deux paliers de vitesse (voir Figure 2).

<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Données du modèle</figcaption></figure>

Sur les intervalles $I_0$ et $I_2$, nous faisons l'hypothèse du régime permanent, ce qui permet de ne résoudre qu'une équation de convection-diffusion plutôt que l'équation de la chaleur. Cette dernière n'est étudiée que sur l'intervalle de temps $I_1$.

Les modèles développés tiennent compte des différents phénomènes thermiques : convection avec l'air ambiant sur les surfaces exposées, conduction avec le moule, mais aussi gestion du contact non parfait aux inter-plis à l'aide de Résistances Thermiques de Contact (RTC).

## 3.2 Abaque stationnaire

Le problème à résoudre ici est donné par la formulation faible

$$ egin{cases} 
\text{trouver } \tilde{u}_S \in H_0^1(\Omega) \text{ tel que } \forall u^* \in H_0^1(\Omega), \\ 
\int_{\Omega} (k_{\parallel} \tilde{u}_{S,x} u_x^* + k_{\perp} \tilde{u}_{S,z} u_z^* - \rho C_p v \tilde{u}_{S,x} u^*) \, dx + \int_{\Gamma_a} h_a \tilde{u}_S u^* \, ds + \int_{\Gamma_M} h_m \tilde{u}_S u^* \, ds \\ 
+ \sum_k \int_{I_1} h_i (\tilde{u}_S^- - \tilde{u}_S^+) u^* \, ds = \int_{\Gamma_L} \Phi u^* \, ds + \int_{\Gamma_a} h_a (u_a - T_0) u^* \, ds + \int_{\Gamma_M} h_m (u_m - T_0) u^* \, ds, 
\end{cases} 
\quad (8) $$

Puisque l'on cherche à déterminer la meilleure puissance pour obtenir une température optimale $T_{opt}$ dans une zone de contrôle, et ce pour les deux paliers $I_0$ et $I_2$, la vitesse et la puissance sont ajoutées comme coordonnées du système. Ainsi, en post-traitant la solution obtenue (normalisation sur la zone de contrôle), nous obtenons un abaque retournant la température en fonction de la vitesse et de la puissance.

Ceci permet de déterminer pour les vitesses $v_e$ et $v_s$ les puissances $q_0$ et $q_2$ à fournir sur les deux paliers $I_0$ et $I_2$ afin d'atteindre une température avoisinant $T_{opt}$ dans la zone de contrôle.

## 3.3 Abaque transitoire

Ici, la formulation faible du système à résoudre est

$$ egin{cases} 
\text{trouver } \tilde{u}_{Dp} \in H_0^1 (\Omega \times I) \text{ tel que } \forall u^* \in H_0^1 (\Omega \times I) \\ 
\int_{\Omega \times I} (k_{\parallel} \tilde{u}_{Dc_x} u_x^* + k_{\perp} \tilde{u}_{Dc_z} u_z^* + \rho C_p (\tilde{u}_{Dc_x} - v(t) \tilde{u}_{Dc_x}) u^*) \, dx \, dt + \int_{\Gamma_a \times I} h_a \tilde{u}_{Dp} u^* \, ds \, dt \\ 
+ \int_{\Gamma_M \times I} h_m \tilde{u}_{Dp} u^* \, ds \, dt + \sum_i \int_{I_1 \times I} h_i (\tilde{u}_{Dp}^- - \tilde{u}_{Dp}^+) u^* \, ds \, dt = \\ 
\int_{\Gamma_L \times I} (\Phi(t) - \Phi(t=0)) u^* \, ds \, dt + \int_{\Omega \times I} \rho C_p (v(t) - v(t=0)) \tilde{u}_{S_x} u^* \, dx \, dt 
\end{cases} 
\quad (9) $$

où $\tilde{u}_S$ est donnée par l'abaque statique et correspond à la solution initiale.

Dans un souci de simplicité, puisque le profil de vitesse est donné ($v_e$ et $v_s$ sont connus), nous avons choisi pour cet abaque de ne prendre comme coordonnées que le temps ($t$), l'espace ($x$ et $z$) et la puissance ($q$). Pour être plus précis, ce $q$ correspond en fait à la puissance à l'instant $t_2$ : nous avons choisi de paramétrer la puissance sous forme d'une droite, à partir de $q_0$ déterminé ci-dessus pour l'instant $t_1$.
La figure donne l'évolution de la température maximale au cours du temps pour différentes valeurs de la puissance q.

<figure><img src="image_4.png" /><figcaption>FIGURE 3 – Évolution de la température pour différentes valeurs de q</figcaption></figure>

À partir de la solution calculée, nous définissons alors la meilleure puissance comme solution de

$$ \min_{q} \left( \sqrt{(T - T_{opt})^2} \right) \qquad (10) $$

en prenant auparavant la norme sur la zone de contrôle (en espace) mais aussi la norme sur une durée caractéristique $τ$. Par exemple, avec $τ = 0.25 * (t_2 - t_1)$, nous obtenons la courbe permettant d'obtenir la meilleure puissance.

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

## 3.4 Conclusion

L'utilisation de cet abaque dynamique permet de déterminer la meilleure puissance afin d'atteindre la température optimale dans la zone de contrôle en chaque instant. Une simulation prenant en compte le profil de vitesse et la loi de chauffe associée est en cours afin de vérifier la pertinence de notre modèle.

## Remerciements

Ces travaux s'inscrivent dans le cadre du projet FUI Impala, financé par OSEO, le Conseil régional d'Aquitaine et le Conseil général des Pyrénées Atlantiques.

## Références

[1] Ammar, A. : The proper generalized decomposition : a powerful tool for model reduction. International Journal of Material Forming 3(2), 89–102 (2009)

[2] Chinesta, F., Keunings, R., Leygue, A. : The Proper Generalized Decomposition for Advanced Numerical Simulations. Springer (2014)

[3] Chinesta, F., Leygue, A., Bognet, B., Ghnatios, C., Poulahon, F., Bordeu, F., Barasinski, A., Poitou, A., Chatel, S., Maison-Le-Poec, S. : First steps towards an advanced simulation of composites manufacturing by automated tape placement. International Journal of Material Forming (2012)

[4] Ladevèze, P. : New algorithms : mechanical framework and development (in french). Technical report, LMT-Cachan 57 (1985)

[5] Lukaszewicz, D.H.J.A., Ward, C., Potter, K.D. : The engineering aspects of automated prepreg layup : History, present and future. Composites Part B : Engineering 43(3), 997–1009 (2012)