FR 3052891 A1 20171222 FR 1655705 A 20160620 1 Procédé d'estimation du facteur d'intensité des contraintes et procédé de calcul de durée de vie associé DOMAINE TECHNIQUE GENERAL L'invention concerne l'analyse de propagation de fissure dans des pièces mécaniques. Ces pièces sont essentiellement destinées à des aéronefs, mais peuvent être n'importe quel composant mécanique. La propagation est définie dans un cadre de propagation par fatigue, avec une série de cycles de chargement. En particulier, l'invention concerne les méthodes et systèmes de simulation, de détermination et d'interpolation de facteur d'intensité des contraintes (appelé FIC par la suite) d'une pièce modélisée numériquement que l'on souhaite analyser. L'invention concerne aussi les méthodes et systèmes de calcul de durée de vie de la pièce à analyser. Le FIC se décompose en trois grandeurs, notées KI, KII et KIII, correspondant aux modes respectivement d'ouverture de la fissure, de cisaillement plan et de cisaillement anti-plan. La présente description sera faite pour KI mais reste applicable aux autres grandeurs. Les méthodes numériques de propagation de fissure sont performantes. Par exemple, les méthodes des éléments finis étendus (XFEM) ou de fissuration conforme permettent la prévision des trajets de propagation de fissure de manière fiable ainsi que le calcul des FIC le long d'un front de fissure. Cependant, le calcul de la durée de vie en propagation de fissure n'est généralement pas intégré dans les codes du commerce ou incompatible avec les exigences industrielles. C'est pourquoi il est usuel que les industries de l'aéronautique développent leur propre code de calcul de la durée de vie dont le 3052891 2 fonctionnement se fonde sur des post-traitements particuliers des résultats issus des simulations numériques de propagation de fissure. Pour une dimension de fissure donnée et en un point donné du 5 front de fissure, le FIC est la grandeur caractérisant l'intensité du chargement au voisinage de la pointe de fissure. Selon l'amplitude de cette grandeur au cours d'un cycle de chargement de fatigue, la fissure va se propager plus ou moins vite. C'est le rôle des simulations numériques de propagation de fissure 10 de fournir les valeurs du FIC le long du front de fissure et pour différentes longueurs de fissure. A partir de cette liste de valeurs, les codes de calculs de la durée de vie doivent être en mesure de prédire la durée de vie en propagation de fissure. 15 Néanmoins, pour des raisons de capacité de calcul, les pas de propagation, c'est-à-dire d'itération du calcul des simulations numériques, ne sont pas pris trop petits, ce qui signifie que toutes les longueurs de fissure ne sont pas simulées. Afin de mener à bien le calcul de durée de vie, les FIC doivent donc être interpolés entre les pas de propagation simulés. Il est donc essentiel d'avoir une méthode d'interpolation efficace afin que les durées de vie prédites soient pertinentes. ETAT DE L'ART On rappelle que l'on se place dans le cadre d'une modélisation de fissure par fatigue, c'est-à-dire une modélisation de l'évolution d'une 30 fissure lorsque plusieurs cycles de charge sont appliqués. Les méthodes de fissuration 3D, généralement obtenues par éléments finis, permettent d'obtenir les FIC le long d'un front de fissure à chaque pas de propagation ainsi que les informations géométriques 20 25 3052891 3 qui la caractérisent (coordonnées des noeuds du front de fissure et faces et noeuds associés à l'une des faces de la fissure). En pratique, seul le FIC maximum le long du front de fissure pour chaque pas de propagation (voir en figure 1 l'unique point Pl, P2, P3, P4 sur chacun 5 des traits courbes qui représente chaque front de fissure Fl, F2, F3, F4 à un pas de propagation différent) intervient dans le calcul de la durée de vie. Pour pouvoir utiliser ces valeurs de FIC, il est utile de les mettre en correspondance avec une longueur de fissure correspondante. 10 Par la suite, les points entre chaque valeur tabulée sont obtenus à l'aide d'une interpolation linéaire entre les points de la tabulation. Il y a deux inconvénients majeurs concernant cette méthode d'interpolation des FIC. 15 Le premier concerne le fait que le FIC renseigné dans le formulaire est le FIC maximum le long du front de fissure pour chaque pas de propagation. En d'autres termes, ce n'est pas forcément le FIC correspondant au même point du front qui est renseigné (les points sur la figure 1 concernent différents endroits du front de fissure). On ne 20 tient ainsi pas compte, avec cette méthode, de la géométrie tridimensionnelle de la fissure et de son évolution au cours du temps qui impacte nécessairement l'évolution du FIC. Le second inconvénient concerne le fait que cette méthode ne permet pas de vérifier si le pas de propagation est suffisamment fin. Il a 25 été montré que si le pas de propagation est trop grossier, les durées de vie en propagation de fissure peuvent ne plus être conservatives. Il y a en outre une dernière faiblesse à cette approche relative à la géométrie de la fissure qu'il faut considérer, qui peut être complexe 30 (bifurcation par exemple). En effet, la fissure est une construction intrinsèquement tridimensionnelle et se pose le problème de la définition de la « longueur de la fissure. Cela limite fortement l'applicabilité de cette méthode. 3052891 4 Il existe des méthodes théoriques pour calculer la propagation de fissures. Par exemple, la méthode des fonctions de poids permet de 5 calculer des facteurs d'intensité des contraintes FIC à partir du champ de contrainte de la structure non-fissurée et de la géométrie de fissure (références [1], [2] et [3] par exemple - les références sont données en fin de description). Néanmoins, cette méthode est limitée à des cas simples ce qui la rend inutilisable en règle générale. 10 Les méthodes dites « perturbatives » sont plus complexes mathématiquement parlant, mais permettent de déterminer semianalytiquement le comportement en fissuration ainsi que la stabilité du front (voir référence [4]). Cependant, la majeure partie des développements ne s'applique que sous l'hypothèse de milieu au moins 15 semi-infini. PRESENTATION DE L'INVENTION L'invention propose ainsi un procédé d'estimation du facteur d'intensité 20 des contraintes dans une pièce modélisée numériquement, dans le cadre d'une modélisation de propagation de fissure par fatigue, comprenant les étapes suivantes : (E2) : obtention, à partir d'une modélisation numérique de la pièce à analyser, d'une pluralité de valeurs simulées en différents 25 points de la pièce à analyser, ladite pluralité de valeurs simulées comprenant, pour différents pas de propagation simulés d'une fissure tridimensionnelle qui apparait sur la modélisation numérique de la pièce, un ensemble des valeurs simulées du facteur d'intensité des contraintes en des 30 points associés, la position de ces points, et des données relatives à la surface fissurée (coordonnées des noeuds et faces des éléments tridimensionnels définissant l'une des faces de la 3052891 5 fissure, les faces des éléments pouvant être stockées sous la forme d'une table de connectivité), - (E3) : détermination, pour chaque pas de propagation et pour l'ensemble de valeurs simulées pour celui-ci, d'une part d'une 5 valeur convertie d'amplitude du facteur d'intensité des contraintes efficace global équivalent correspondant à celui d'une fissure plane à front droit et d'autre part d'une longueur convertie de la fissure équivalente considérée, lesdites valeurs converties étant déterminées par égalisation de l'énergie dissipée dans la fissure 10 tridimensionnelle de la modélisation numérique et de l'énergie dissipée dans la fissure d'un modèle standard de fissure plane à front droit, les énergies étant elles-mêmes déterminées en fonction des facteurs d'intensité des contraintes, - (E4) : Interpolation de valeurs converties d'amplitude de facteur 15 d'intensité des contraintes efficace global équivalent entre deux longueurs de fissure converties (a) successives, - (E5) : Stockage des valeurs converties d'amplitude de facteur d'intensité des contraintes efficace global équivalent ainsi interpolées avec les longueurs de fissure associées. 20 Ainsi, on se rapporte à un procédé global « thermodynamiquement » équivalent à celui d'une fissure plane à front droit. Une notion de longueur possédant un sens physique est générée, ce qui permet de s'affranchir des limitations citées précédemment. 25 L'invention concerne aussi les caractéristiques suivantes, prises seules ou en combinaison : - l'étape d'interpolation met en oeuvre une interpolation linéaire par morceaux, 30 -l'étape d'interpolation met en oeuvre une interpolation par minimisation de l'énergie de courbure, - l'étape d'interpolation comprend les sous-étapes suivantes : interpolation à l'aide d'une interpolation linéaire, 3052891 6 interpolation à l'aide d'une interpolation par minimisation de l'énergie de courbure, les deux interpolations étant interchangeables, calcul d'au moins une grandeur représentative d'une différence de valeurs d'amplitude du facteur d'intensité des contraintes entre les 5 deux interpolations, ladite différence étant calculée strictement entre deux valeurs de longueurs de la fissure correspondant à deux pas successifs, comparaison de ladite grandeur représentative de la différence avec un seuil prédéterminée, 10 si la grandeur représentative de la différence est supérieure au seuil, génération d'une consigne d'instruction de calcul à partir de la simulation numérique de la pièce à analyser des valeurs du facteur d'intensité des contraintes (FIC) le long du nouveau front simulé et des valeurs de la position de la fissure, entre les deux pas successifs, 15 - le procédé comprend préalablement à l'étape d'obtention (E2), une étape de calcul (El) mettant en oeuvre la simulation par éléments finis (ou éléments finis étendus) à des pas successifs de l'évolution de la fissure dans la pièce à analyser, ladite simulation étant effectuée par l'unité de traitement, 20 - les données numériques récupérées à d'obtention (E2) des données sont obtenues au moyen d'une simulation par éléments finis ou d'une simulation par éléments finis étendus, - pour l'étape de détermination (E3), le premier sous-ensemble de données permet de connaitre l'amplitude efficace du taux de restitution 25 d'énergie AG par les relations suivantes : 1 G = F (Ki+ Ki2i) + 1 - Ki2ii li E* = E en contraintes planes E* = E / (1 - u2) en déformations planes (AKeff (S))2 = A Ge f f (s). E* Où KI, Ku et Km sont, respectivement, les coefficient de FIC correspondant aux modes d'ouverture de la fissure, de cisaillement dans le plan et de cisaillement anti-plan, E est le module de Young, p est le 3052891 7 module de cisaillement, E* est le module de Young équivalent, v est le coefficient de Poisson, le deuxième sous-ensemble de données permet de connaitre l'incrément de surface fissuré par unité de longueur de front de fissure 5 dSurf, ces formulations permettent d'obtenir la valeur de l'énergie dissipée : Energie dissipée =f AG e f f (s) . dSurf (s). dl front de fissure dans lequel cette énergie dissipée est égalisée avec l'énergie dissipée de la fissure d'un modèle standard de fissure plane à front droit, 10 s'exprimant sous la forme : - glob glo b ,glob Energie dissipée = AGef f . dSurf ',front Où la notation glob est relative à ce modèle standard, dans lequel, grâce à une modélisation par une loi de Paris avec les paramètres de fermeture de fissure d'Elber, le lien entre dSurfeb et A Ggl° b est effectué ef f dSurfglob C. GLAGglob ' E*) n eff dN 15 Où C et n sont des coefficients de Paris, dans lequel, en égalisant les deux énergies dissipées, on obtient 1 AG e f f (s) . dSurf (s). dl 12+1 ( 1 Energie dissipée n+1 Meff glob 'front de fissure = C. (E*)nLfront C. (E*)nLfront dans lequel, au moyen des relations entre 3,Geff et AKeff, on alo détermine AK- et dans lequel on détermine Surfeb, c'est-à-dire la ef fb longueur a. 20 L'invention concerne également un procédé d'évaluation de la durée de vie d'une pièce modélisée numériquement, dans lequel on met en oeuvre un calcul de la durée de vie de la pièce faisant intervenir les valeurs converties d'amplitude de facteur d'intensité des contraintes efficace 3052891 8 global équivalent interpolées en fonction d'une longueur de fissure équivalente selon un procédé tel que décrit précédemment. L'invention concerne également un système comprenant une unité de 5 traitement comprenant des moyens de calcul et une mémoire, l'unité étant configurée pour mettre en oeuvre un procédé d'estimation ou un procédé d'évaluation de durée de vie en propagation de fissure par fatigue tel que décrit précédemment. 10 L'invention concerne également un produit programme d'ordinateur configuré pour être mis en oeuvre par un système tel que décrit précédemment, et comprenant des instructions pour provoquer la mise en oeuvre d'un procédé d'estimation ou un procédé d'évaluation de durée de vie tel que décrit précédemment. 15 PRESENTATION DES FIGURES D'autres caractéristiques, buts et avantages de l'invention ressortiront de la description qui suit, qui est purement illustrative et non limitative, 20 et qui doit être lue en regard des dessins annexés, sur lesquels : - La figure 1 illustre un front de fissure pour différents pas de propagation, - La figure 2 illustre un système permettant de mettre en oeuvre l'invention, 25 - La figure 3 illustre une pièce à analyser, - La figure 4 illustre deux procédés et des modes de réalisation associés selon l'invention, - La figure 5 illustre des points convertis en équivalence avec une fissure plane à front droit, à partir d'une fissure tridimensionnelle, 30 selon une théorie d'équivalence énergique par analyse thermodynamique, 3052891 9 - La figure 6 illustre ces mêmes points avec deux méthodes d'interpolation : linéaire par morceaux et par minimisation de l'énergie de courbure, - La figure 7 illustre l'évolution de l'interpolation polynomiale avec 5 le degré du polynôme, - La figure 8 illustre schématiquement une fissure plane à front droit, - Les figures 9a à 9c illustrent des fonctions de forme 1D issues d'une généralisation des éléments finis et servant à calculer 10 l'énergie de courbure, - Les figures 10 à 10c illustrent un changement de représentation pour lisser les métriques. DESCRIPTION DETAILLEE 15 En référence à la figure 2, on a représenté un système 10 d'interpolation de valeur de facteurs d'intensité des contrainte (FIC) K d'une pièce à analyser 20, qui est modélisée numériquement. Comme indiqué en introduction, le FIC K se décompose en trois données KI, KII, 20 KIII. La description sera donnée pour KI uniquement. La pièce à analyser 20 est une pièce destinée à l'aéronautique, pour laquelle la durée de vie doit pouvoir être estimée. La pièce 10 est typiquement une aube fan de compresseur ou de turbine, un disque moteur, une bride, une attache moteur, un carter telle qu'illustrée en 25 figure 3. Cette liste est donnée à titre illustratif, puisque, pour la mise en oeuvre du procédé, le type de pièce n'intervient cependant pas. Le système 10 comprend une unité de traitement de données 12, par exemple un ordinateur ou un serveur, disposant de moyens de 30 calcul 14 configurés pour mettre en oeuvre un procédé qui sera décrit plus en détail ci-après et, disposant avantageusement, d'une mémoire 16. Les moyens de calcul 14 peuvent par exemple être un calculateur de type processeur, microprocesseur, microcontrôleur, etc. La mémoire 16 3052891 10 peut être par exemple un disque dur, une mémoire « flash », ou un espace de stockage délocalisé, de type « cloud ». L'unité de traitement de données 12 peut également être adaptée 5 pour mettre en oeuvre des simulations numériques, telles qu'éléments finis, de la pièce à analyser 20. La simulation par éléments finis permet d'obtenir, pour chaque pas de simulation, des données relatives à ladite pièce et à la propagation de fissure. En particulier, il est possible d'obtenir pour chaque pas de 10 propagation, les valeurs du FIC en des points associés et la position de ces points, généralement le long du front de fissure. Les positions comprennent par exemple les coordonnées des noeuds du maillage définissant le front de fissure. Il est possible d'obtenir des données relatives à la surface fissurée : outre les coordonnées de noeuds, les 15 faces des éléments tridimensionnels définissant une parmi les faces de la fissure. Ces données sont stockées sous la forme d'une table de connectivité, classiquement connues en éléments finis. D'autres données, telles que le maillage de la surface fissurée ainsi que du front peuvent aussi être obtenues. 20 En référence à la figure 4, un procédé 100 d'interpolation des valeurs d'un FIC d'une fissure tridimensionnelle simulée par éléments finis va être décrit. Ce procédé 100 est aussi avantageusement utilisé pour améliorer 25 l'interpolation du FIC. Cela sera décrit par la suite. Dans une première étape El, une simulation numérique est effectuée sur la pièce à analyser 10. Cette simulation permet d'obtenir les données précitées. 30 Dans un mode de réalisation préférée, cette simulation est faite par éléments finis ou par éléments finis étendus. Cette étape El est effectuée soit par l'unité de traitement 12 du système, soit par un autre système. 3052891 11 Dans les deux cas, une étape E2 de réception desdites données (générées durant l'étape El) par l'unité de traitement 12 du système est mise en oeuvre. Plus précisément, on définit pour chaque front de fissure 5 simulé un premier sous-ensemble de données qui contiennent les valeurs du FIC en des points associés sur le front de fissure et leur position respective associée, et un deuxième sous-ensemble de données qui contient les données relatives à la surface fissuré (les coordonnées des noeuds définissant le front de fissure ainsi que le maillage surfacique 10 d'une des faces de la fissure, par exemple la table de connectivité des faces des éléments tridimensionnels concernés) telles qu'elles ont été définies précédemment. Ayant une pluralité de fronts simulés à des pas différents, on obtient ainsi une pluralité de premiers et de deuxièmes sous-ensembles de données. 15 Il s'agit de la première étape du procédé d'interpolation à proprement parler, dans la mesure où l'étape de simulation El n'est pas nécessairement mise en oeuvre expressément pour une interpolation ultérieure. 20 Dans l'étape E3, la problématique de la géométrie de la fissure en trois dimensions est ramenée à une problématique de fissure plane à front droit. La fissure plane à front droit est un modèle standard connu de l'homme de l'art. Pour cela, une équivalence thermodynamique fondée sur une égalité de l'énergie dissipée entre les deux modélisations 25 est effectuée. Le fondement physique est expliqué ci-dessous. Le premier principe de la thermodynamique, pour une transformation du front de fissure infinitésimale donne, exprimée en densité d'énergie par unité de longueur de front de fissure : SQ 6Wext densité d'énergie interne densité de chaleur densité de travail des ef forts extérieurs Le second principe de la thermodynamique permet d'exprimer 30 l'irréversibilité de la transformation : T. dS = SQ + G. dSurf densité d'énergie dissipée 3052891 12 La densité d'énergie dissipée est toujours positive. Le terme dSurf représente l'incrément de surface fissurée par unité de longueur de front de fissure au cours de la transformation infinitésimale. Cette quantité est homogène à une longueur : il s'agit en quelque sorte d'un incrément 5 de longueur de fissure au sens « thermodynamique » du terme. Il apparait ainsi une grandeur « physique » que l'on pourra considérer pour effectuer une interpolation qui présente un sens physique. Les éléments précédents permettent de rendre le problème de fissuration tridimensionnel initial équivalent à un problème de 10 propagation de fissure plane à front droit en élasticité plane. Pour plus de généralité, nous nous plaçons dans le cas de l'élasticité adaptée, c'est-à-dire que le matériau a déjà pu subir une déformation plastique initialement, mais la propagation de fissure se fait sous un chargement cyclique de fatigue en élasticité. 15 Le champ de température peut évoluer au cours du temps mais le champ de température ne varie pas spatialement. Par conséquent, à un cycle de fatigue est associée une température. Une hypothèse de modélisation consiste à considérer que la loi de propagation de fissure est une loi de Paris avec correction d'Elber 20 (références [5], [6]), donnée par l'équation suivante : dSurfn(T) n(T) = C (T). (AK. KR, n) = C(T). (AK, f f) avec KR , n = a(T). R + b(T) dN avec C(T) et n(T) sont les coefficients de la loi de Paris et a(T) et b(T) sont les paramètres de la loi de fermeture de fissure d'Elber, AK est l'amplitude du FIC, qui permet de s'affranchir du rapport de charge, (c'est-à-dire du rapport R=Kmin/Kmax) et AKeff est l'amplitude efficace 25 du FIC, qui prend en compte l'effet de fermeture de la fissure. Le lien entre AK et AKeff a été établi par Elber. L'incrément de surface dSurf fissuré par cycle est maintenant relié à l'amplitude du facteur d'intensité de contrainte efficace AKeff. Il existe une relation entre K (ou AKeff, par équivalence) et le 30 taux de restitution d'énergie G (ou l'amplitude dudit taux efficace AGeff par équivalence), donné par les relations suivantes : 3052891 13 1 ., ., 1 ., G = -E* (Kin + Kt') + - Kfii li E* = E en contraintes planes E* = E / (1 - u2) en déformations planes où E est le module de Young, E* est le module de Young équivalent, p est le module de cisaillement, y est le coefficient de Poisson. Le FIC est connu en chaque point du front de fissure et pour chaque pas de propagation à l'aide de la simulation numérique de 5 l'étape El et récupéré à l'étape E2. En remplaçant les valeurs par leur amplitude dans les relations précédentes, on obtient, par définition, une relation entre le taux de restitution d'énergie et l'amplitude du facteur d'intensité des contraintes efficaces (s étant l'abscisse curviligne le long du front de fissure) : (AKeff (s))2 = AG e f f(s).E* 10 L'énergie dissipée au cours du processus de fissuration s'obtient par la relation suivante (en restant dans l'hypothèse de transformation infinitésimale, donc en ne considérant qu'un cycle de chargement) : Energie dissipée =f AG e f f (s) . dSurf (s). dl front de fissure Le dSurf est connu à l'aide des données relatives à la surface 15 fissurée calculées par la simulation numérique à l'étape El et récupérées à l'étape E2. En faisant le postulat d'une équivalence du problème de fissuration 3D avec un problème de fissuration plane à front droit, on peut égaliser les énergies dissipées et les longueurs de front de fissure. 20 Par conséquent, il est légitime d'écrire : - glob glob ,glob Energie dissipée = AG e f f . dSur f . ',front La notation « glob » signifie que la donnée est propre à la modélisation en fissuration plane à front droit. 25 La fissuration plane à front droit respecte elle aussi la loi de propagation de Paris. On obtient par conséquent que 3052891 14 dSurfglob c. GLAGglob ' E*) n eff dN En égalisant les deux énergies dissipées et en réécrivant l'équation, on obtient la relation suivante 1 AGeff(SIdSlill(S).d/ n+1 lob = (Energie dissipée)n+1 ='front de fissure 1 "1j eff *)n-Lfront *)n-Lfront C.(E C.(E 5 Or dSurfeb/dN= dSurfeb (car dN = 1 ici, puisqu'on ne considère qu'un seul cycle) est homogène à une longueur de fissure. Par conséquent, son intégration au cours du temps entre le premier cycle et le dernier cycle permet d'obtenir une longueur de fissure équivalente Surfeb, notée longueur a, au sens thermodynamique. On obtient bien 10 une longueur qui a un sens physique. La relation précédente établie entre AKge7fb et AGglob permet ef f d'injecter l'amplitude du FIC efficace dans l'équation précédente On obtient ainsi une relation reliant l'amplitude du FIC efficace global équivalent dans la modélisation en fissure plane à front droit en 15 fonction de la longueur de la fissure Surfeb, c'est-à-dire la longueur a. Ces données-là pourront ensuite être interpolées. Cette construction physique est valide en toute généralité. Cependant, elle nécessite de connaître toutes les informations du problème pour chaque cycle de chargement ce qui est généralement 20 trop coûteux en termes de temps de calculs et d'espace mémoire pour la sauvegarde de l'intégralité des résultats (au moins un millier de cycle devraient être simulés). Les calculs par éléments finis se font souvent à pas de propagation maximal fixe le long du front de fissure, ou pré-intègrent sur un nombre de cycles donné (10, 100 ou 1000 par 25 exemple). Cela revient à dé-raffiner la discrétisation en pas de propagation pour réduire le temps de calcul. Par conséquent la loi de propagation est remplacée par : lob dSurq)neci dN = C. coe f f . GlAGegfIcib. E*) = Ceq. (lAGegr. E* 3052891 15 où coeff est le coefficient qui permet la pré-intégration de plusieurs cycles Cela revient à modifier les coefficients de la loi de Paris C et n par les coefficients Ceq = C. coeff et neq = n. 5 Le lien entre dSurf egqi°b et dSurf eb est donc le suivant : dSurfegqlob = coeff. dSurfglob De plus, le terme dSurf (s) est remplacé par la surface incrémentée entre deux pas de propagation simulés dSurfeq(s) (donc non infinitésimale). En conclusion, le raisonnement précédent s'applique approximativement en effectuant les remplacements précédents et la 10 longueur de fissure équivalente se déduit toujours de la même manière. Ainsi, dans l'étape E3, à partir de la pluralité de premiers et de deuxièmes sous-ensembles de données récupérées à l'étape E2, c'est-à-dire les données du FIC le long des différents fronts de fissures simulées 15 à des pas différents et des données relatives aux surfaces fissurées, et au moyen de ce qui a été expliqué précédemment, on détermine l'amplitude efficace globale du FIC AKgef.1°Ib et la longueur a de la fissure plane à front droit associée. Cette détermination est mise en oeuvre par l'unité de traitement 12. 20 Pour un front de fissure donné, les données du FIC simulées permettent de connaitre K (premier sous-ensemble de données) et au moyen des équations précédentes permettent de calculer le AKgef.1°Ib et les données de la surface fissurée dSurfeq (deuxième sous-ensemble de données) permettent de calculer la longueur a. Du fait de la pluralité de 25 données correspondant à chaque front, on obtient une pluralité de couple (a, Mgef7). L'équivalence obtenue peut aussi être appliquée à une fissure plane à front droit. La relation suivante montre que le système reste 30 invariant : 3052891 16 1 A &lob "front de fissure 'AGef f(S) . dSur feq(s). dl neq+1 Ceq- (Elnec I. ',front 1 AGef f. Ceq- (AGef f. rreci neq+1 Ceq. (E*)neci = AGeff De cela découle le fait que les surfaces fissurées sont les mêmes, donc le système final est identique au problème initial. Ceci prouve la cohérence de la méthode fondée sur une analyse thermodynamique du 5 problème. La figure 5 illustre dans le plan la représentation des couples (a, 3,Kgefiolb) . Les références F1, F2, F3, F4 sont données en référence à la figure 1, pour signifier à quel front de la fissure d'origine correspondent 10 les points du graphe. Il s'agit d'un abus d'écriture. On rappelle à présent que le but est de pouvoir déterminer des FIC pour des états de fissure qui n'ont pas été simulés numériquement, c'est-à-dire pour des états qui se produisent strictement entre deux pas successifs. 15 A cette fin, dans une étape E4, une interpolation des AKgef.1°Iben fonction de la longueur a est effectuée. L'interpolation est effectuée par les moyens de calcul 14 de l'unité de traitement 12. L'interpolation peut être effectuée selon plusieurs méthodes. On 20 privilégie les interpolations qui ne créent pas d'informations au problème ou qui en ajoute un minimum. Pour cela, une interpolation linéaire par morceaux IL ou une interpolation par minimisation de l'énergie de courbure IMC conviennent. Les deux interpolations peuvent aussi être effectuées, dans la mesure où elles peuvent avoir des applications 25 différentes. La figure 6 illustre ces deux interpolations. 3052891 17 L'interpolation par minimisation de l'énergie de courbure globale peut être vue comme une méthode fondée sur des considérations physiques propres au problème de fissuration. En effet, la littérature scientifique sur le sujet montre que si une fissure se propage sans 5 brancher, alors l'évolution du facteur d'intensité des contraintes avec l'extension de fissure peut être représentée par une fonction infiniment dérivable. Cela signifie que le fait de calculer la courbure d'une telle fonction présente un sens. De plus, minimiser une énergie de courbure intervient dans 10 quelques situations physiques comme, par exemple, dans la physique des bulles de savons ou plus généralement dans la physique des membranes. La figure 7 illustre une interpolation linéaire entre trois points. 15 Cet exemple nous montre qu'il est possible d'approximer de mieux en mieux une interpolation par des polynômes. L'approximation sera d'autant meilleure que le degré du polynôme sera élevé. Il est aussi possible de démontrer que l'écart maximum entre l'interpolation linéaire et l'approximation polynômiale tend vers zéro lorsque le degré du 20 polynôme tend vers l'infini. Ainsi, l'interpolation linéaire peut être vue comme le passage à la limite à l'infini que l'on appliquerait à une fonction infiniment dérivable. Cette approche confère un caractère aussi « physique » à l'interpolation linéaire par morceaux qu'à l'interpolation fondée sur la 25 minimisation de l'énergie de courbure. Réaliser une interpolation consiste à compléter l'information manquante à l'aide d'un principe qui doit au minimum respecter la « physique » du problème considéré, c'est-à-dire ne pas introduire d'éléments qui violeraient certaines équations fondamentales auxquelles 30 le problème obéit. C'est pour cette raison que l'interpolation par minimisation de l'énergie de courbure est une méthode « physique » car elle ne contredit pas la théorie de la mécanique de la rupture. La méthode consistant en 3052891 18 des interpolations linéaires par morceaux n'est pas physique au premier abord mais peut être vue comme le passage à la limite d'une méthode d'interpolation polynômiale globale qui, elle, serait physique et qui par extensions rendrait la méthode d'interpolation linéaire par morceaux 5 presque « physique » également. Un autre point de vue concernant la méthode d'interpolation linéaire serait de la voir comme une méthode de minimisation de la distance entre des points. Cette manière de voir les choses la place alors dans la catégorie des méthodes optimales. L'interpolation linéaire 10 minimise la distance tandis que l'autre minimise l'énergie de courbure sous les contraintes des positions des points à relier. La minimisation des distances peut être vue comme un moyen de relier des points sans introduire d'informations supplémentaires au problème initial et cela résulte en un champ continu. Dans le cas de la seconde méthode, on 15 introduit au moins une information supplémentaire qui est que le champ est régulier. Dans ce cas, la minimisation de l'énergie de courbure est une méthode optimale au sens où on minimise l'information physique à introduire alors que la méthode de l'interpolation linéaire est optimale au sens où on minimise l'information purement mathématique. 20 Enfin, dans une étape E5, après l'étape d'interpolation, les données interpolées sont stockées, afin de pouvoir y avoir accès pour d'autres applications. Le stockage se fait typiquement dans la mémoire 16. On parle généralement de « formulaire ». Le formulaire est un 25 tableau qui donne les facteurs d'intensité des contraintes en fonction de la longueur de la fissure. Par conséquent, l'étape d'interpolation E4 génère un fichier, sous forme d'un format texte ou d'un tableau contenant le formulaire, c'est à dire regroupant les données converties et les données converties 30 interpolées. L'étape E5 consiste à stocker ce fichier. En effet, comme indiqué en introduction, l'obtention du formulaire est généralement utilisée pour calculer la durée de vie d'une pièce. Par 3052891 19 conséquent, en référence à la figure 4, un procédé 200 de détermination de durée de vie de la pièce à analyser 20 est décrit. Dans une étape El, l'unité de traitement 12 reçoit le formulaire calculé à l'étape E4 et/ou stocké à l'étape E5 du procédé précédent. 5 Cette étape de récupération peut simplement consister à avoir accès à la mémoire 16. Dans une étape E'2, un procédé de calcul de la durée de vie en propagation de fissure par fatigue de la pièce est mis en oeuvre. Le document [7] décrit un tel procédé. 10 Il suffit d'utiliser un code existant et de remplacer le formulaire de FIC par celui proposé dans cette invention pour déterminer la durée de vie en propagation de fissure. Ce procédé 200 peut être mis en oeuvre par un autre système que celui décrit précédemment. 15 Raffinement de l'erreur Lorsque l'on dispose de la totalité des informations concernant le problème à résoudre, le résultat de l'interpolation sera toujours le même 20 quelle que soit la méthode d'interpolation utilisée et qui respecte la totalité des contraintes « physiques » du problème. Typiquement, sur les intervalles I1 de la figure 6, l'interpolation apparait intuitivement comme étant de qualité ; sur les intervalles 12, elle apparait intuitivement comme de moins bonne qualité. 25 Il est ainsi possible de voir les écarts entre deux interpolations optimales dites « physiques et ou mathématiques » comme un moyen d'identifier un manque d'information qui serait utile pour améliorer l'interpolation ou le degré de confiance en l'interpolation. Cette considération permet de mettre en place une étape de 30 vérification de la pertinence, illustré en figure 6. Dans deux sous-étapes E41 et E42 interchangeables de l'interpolation, on effectue une interpolation linéaire par morceaux et une interpolation par énergie de courbure. 3052891 20 On remarque que strictement entre deux valeurs de longueur correspondant à deux pas successifs simulés, c'est-à-dire strictement en zone interpolée, il existe des différences 15 de valeur AI(gi°b entre les e f f deux interpolations. 5 Dans une sous-étape E43, on calcule au moins une de ces différences 15, que l'on compare, dans une sous-étape E44, à une valeur seuil prédéterminée VSP. Le seuil prédéterminé VSP est choisi en fonction de la qualité de l'interpolation souhaitée. 10 Il est aussi possible de comparer une valeur de ces différences pour chaque zone entre deux pas successifs, ou de comparer une moyenne, ou une valeur maximale, etc. On parlera donc plus généralement d'une grandeur br représentative de la différence. Cette grandeur br indique qu'il existe un écart que l'on a pu quantifier par les 15 fonctions citées précédemment. Enfin, une sous-étape de comparaison E45 entre la grandeur br représentative et le seuil VSP est effectuée. Si la grandeur représentative br est inférieure au seuil VSP, on peut considérer que les deux interpolations sont de qualité et que les 20 données, si elles avaient été simulées, seraient proches des deux valeurs d'interpolation. Si la grandeur représentative br est supérieure au seuil VSP, alors on peut considérer qu'il y a une incertitude trop grande sur l'interpolation, et que l'obtention de nouvelles valeurs simulées est utile, 25 voire nécessaire. A cette fin, on génère une consigne d'instruction de calcul à partir de la simulation numérique, typiquement par éléments finis (ou éléments finis étendus) de la pièce à analyser 20 d'au moins une valeur du facteur d'intensité des contraintes FIC et la position du front de fissure entre les deux pas successifs. En d'autres termes, l'unité 30 de traitement 12 envoie une consigne pour effectuer une partie de l'étape El, avec un pas réduit pour au moins une simulation. On parle ainsi de raffinement de la discrétisation en longueur de fissure. 3052891 21 A titre d'exemple, une valeur seuil prédéterminée VSP de 2% peut convenir. Le critère dépend du cahier des charges. On peut par exemple décider que br =15 et que dès que br/max(mgef/ofb parmi interpolation linéaire et courbure) > 2% pour une 5 valeur de longueur a donnée, une consigne de re-calcul de la simulation est générée. Enfin, l'étape de stockage E5 peut comprendre le stockage des deux interpolations et des grandeurs représentatives des différences br et/ou des différences Ô. 10 ANNEXE 1 : exemple d'interpolation par minimisation de l'énergie de courbure L'interpolation par minimisation de l'énergie de courbure ne donne 15 des prédictions physiques cohérentes que dans le cadre d'une fissure plane à front droit en élasticité plane, telle qu'illustré en figure 8, où a représente la contrainte et a la longueur de la fissure. Or, le procédé décrit précédemment permet de considérer thermodynamiquement toute fissure tridimensionnelle comme une fissure vérifiant ces 20 propriétés-là. Pour effectuer cette interpolation, il faut déterminer l'énergie de courbure. 25 On définit la fonction f de la manière suivante (x est la longueur de fissure « virtuelle », appelé a auparavant) : di(' = f(x) dx La nécessité d'introduire des longueurs de fissure « virtuelles » est explicitée en ANNEXE 2 La fonction f s'écrit comme une combinaison linéaire des fonctions 30 polynômiales par morceau (de degré 5 pour assurer une régularité C2 sur toute la plage d'évolution de la longueur de fissure « virtuelle », il 3052891 22 s'agit d'une généralisation des éléments finis classiques) de référence données ci-dessous : cla(z) = 0,5 - 0.9375. z + 0,625. z3 - 0,1875. z5 C1,2(z) = 2. (0,3125 - 0.4375. z - 0.375. z2 + 0,625. z3 + 0,0625. z4 - 0,1875. z5) C1,3(z) = 4. (0,0625 - 0.0625. z - 0.125. z2 + 0,125. z3 + 0,0625. z4 - 0,0625. z5) C2,1(z) = 0,5 + 0.9375. z - 0,625. z3 + 0,1875. z5 C2,2(z) = 2. (-0,3125 - 0.4375. z + 0.375. z2 + 0,625. z3 - 0,0625. z4 - 0,1875. z5) C2,3(z) = 4. (0,0625 + 0.0625. z - 0.125. z2 - 0,125. z3 + 0,0625. z4 + 0,0625. z5) x - xi+1 x - xi z = . + xi+1 - xi+1 - xi Le premier indice indique le numéro du noeud sur lequel la valeur nodale est prise. 5 Le deuxième indice indique la nature de la valeur nodale : 1: Valeur nodale du facteur d'intensité des contraintes (FIC) 2: Valeur nodale de la dérivée par rapport à la variable z du FIC 3: Valeur nodale de la dérivée seconde par rapport à la variable z du FIC L'élément unidimensionnel associé à ces fonctions de forme est 10 définit de manière à ce que le premier noeud soit situé à z = -1 et le deuxième à z = 1. L'abscisse curviligne associée à la courbe donnant l'évolution du facteur d'intensité des contraintes en fonctions de la longueur de fissure « virtuelle » est donnée par : i-1 fx s(x) = ,\/(1+ (f(u))2) du = xi+1,\/(1 + (f(u))2) . du + ,\/(1 + (f(u))2) . du i=t xi X. 15 On notera : { x - xi+1 x - xi z = . . + . xi+1 - xi xi+1 - xi f(x) = fi (z) = fll.c1,1(z) + f12 . C1,2 (Z) + f13. C1,3 (Z) + f11.c2,1 (z) + f1+ 1- c2,2 (z) + f1+ 1- c2,3 (z) On notera que les degrés de liberté fll et flll = fLsont connus. L'optimisation portera sur les autres degrés de liberté. D'où l'on obtient (où Lx est la longueur de l'élément virtuel) : J-1f sx(x) = ,\/(1 + (fi(z))2). dz + -1 ,\/(1 + (42) . du) . 1=1 1 3052891 23 Le vecteur unitaire tangent s'obtient comme suit : t(Sx) = dsx I f Les figures 9a à 9c illustrent les valeurs des fonctions de formes C11, C12, C13, C21, C22 et C23 et leur dérivée et dérivée seconde. On remarque que ces courbes vérifient leurs contraintes spécifiques (0 ou 5 1) au niveau des noeuds. Nous pouvons désormais calculer le vecteur unitaire tangent à la courbe KI. Les formules données précédemment conduisent à : t = dYx = dYx dx 1 1 x dsx dx dsx (f(1 X)) dsx dz [ 21 2 dz ' dx ) ; X - Xj+1 X - Xi ) (f(lx)) . 1 + (fJ (xi+1 _ xi + xi+i _ xi Il ne reste plus que le calcul de la courbure. Pour cela, la dérivée 10 de tx par rapport à l'abscisse curviligne est calculée : dtx d idyx dx \ d2yx dx \ 2 dyx d2x ( 0 \ A 2 1 dzx dsx dsx dx dsx) dx2 dx d(sx)2 (f' (x)) (dsx) ( f(x)) d(sx)2 Pour déterminer les termes en dérivées par rapport à l'abscisse curviligne, les formules de passage différentielles de la coordonnée curviligne à la coordonnée virtuelle sont exprimées : 7 d \ dx 0 7 dsx dsx dx dz d2x dx )2 d2 \d(sx)2 / ci(sx)2 dex \ j / Diffsx Cette équation matricielle est inversible ce qui permet d'obtenir : 7 d \ dsx0 7 d dx dx dsx dz dzsx idsx \ 2 - d2 dX2 dx ) \d(sx)2/ Divffx 15 On remarque que Diffx-1 = Diffsx, ce qui permet de déterminer (dds:)2 et aXfacilement. Finalement, on trouve: dtx dsx = (,'(x)) . Diffx -1 (2,2) + (f(1-x)) . Diffx -1 (2,1) 3052891 Avec : 24 Ex x d Diffx = d2sx dx2 (dsx21 dx 0 La courbure s'en déduit automatiquement : x(x) = (Diffx-1(2,1))2 + (fr(x).Diffx-1(2,2) + f(x).Diffx-1(2,i))2 La matrice Diffx-1 doit désormais être explicitée. Pour cela, la formule classique suivante d'inversion de matrice 2x2 est utilisée: A = (a 1 ) detA ( d a ) = avec det(A) = a. d - b. c d() 5 On a donc : x - xj+1 x - Xi Z = Xi+1 - Xi + Xi+1 - Xi ds 3 2 3 det(Diffx) = dx I = [1 + (fi (z)) Diffx 1(x) det(Diffx) 4 fj(z)._dfi (z) 1 + (fi (z))2I Lx dz La fonctionnelle que l'on cherche à minimiser est la suivante (si le maillage est constitué de n-1 élément) : xn Wx = (f(U))2) (K(u))2. du On utilise préférablement alors un algorithme de Newton pour 10 déterminer les paramètres fk qui minimisent cette fonctionnelle. Il est ainsi nécessaire de connaître les vecteurs gradients et matrices Hessiennes (par rapport aux variables inconnues) des grandeurs intervenant dans le problème d'optimisation. 2 2 1 + fi (z) 0 3052891 25 2. (Diffx-1(2,1)) . Grad (Diffx-1(2,1)). +2. (f'(x). Diffx 1(2,2) + f(x). Diffx 1(2,1)) Grad(r (x))i. Diffx-1(2,2) + f'(x). Grad (Diffx-1(22)) ( +Grad(f(x))i. Diffx 1(2,1) + f(x). Grad (Diffx 1(2,1)) 2. Grad (Diffx-1(2,1)) . Grad (Diffx-1(2,1)) i i +2. (Diffx 1(2,1)) . H (Diffx 1(2,1)).. 1j 7Grad(fr (x))i. Diffx-1(2,2) + f'(x). Grad (Diffx-1(2,2)) \ i +2. Grad(f(x))i. Diffx-1(2,1) + f(x). Grad (Diffx-1(2,1)) / Grad(f ' (x))i. Diffx-1(2,2) + f'(x). Grad (Diffx-1(22))i +Grad(f(x))i. Diffx 1(2,1) + f(x). Grad (Diffx 1(2,1))i ( +2. (f'(x). Diffx 1(2,2) + f(x). Diffx 1(2,1)) 7 H(fr (x))ii. Diffx-1(2,2) + Grad(r (x))i. Grad (Diffx-1(22)) \ i +Grad(fr (x))i. Grad (Diffx 1(2,2)) + f'(x). H (Diffx 1(2,2)) i ii Grad(K(u)2)i = - Fi(x(u)2)i = H (f(x))ii. Diffx-1(2,1) + Grad(f(x))i. Grad (Diffx-1(2,1)) +Grad(f(x))i. Grad (Diffx 1(2,1))i + f(x). H (Diffx 1(2,1))ii 1 (1 + (f(u))2) . Grad(f(u))i. (K(u))2. du\ + xn 2. f ,\/ (1 + (f(u))2) Grad(K(u)2)i. du 3 1Xn (1 + (f(U))2) Z . Grad(f(u))i. Grad(f(u))i. (K(u))2. du\ H (Wx)ii = f xn 1 + JX1 (1 + (f(U))2) 2 . H(f(U)). (K(u))2. du xn + ,\/ (1 + (f(u))2) H(K(u)2)ii. du Cet algorithme peut être mis en oeuvre par l'unité de traitement 12. ANNEXE 2: Méthode des métriques variables Grad(Wx)i = f xn x1 5 3052891 26 La figure 10a illustre les points du front à différents pas de propagation. On remarque que la distance entre deux mêmes points du front n'est pas constante entre deux pas, ce qui signifie que l'échelle de longueur associée à chaque élément est différente. On parle alors de 5 métriques différentes. Le fait d'avoir des éléments avec des métriques différentes (discontinue donc au point d'intersection) tout en imposant un raccord de régularité C2 entre les éléments conduit à des instabilités 10 numériques. Un autre point de vue consiste à dire que les éléments n'ont pas les mêmes contributions à la fonctionnelle (du fait de leurs dimensions différentes) alors qu'ils véhiculent tous la même quantité d'information. 15 Quel que soit le point de vu, cela peut induire un biais numérique qu'il faut avantageusement traiter. Pour résoudre ce problème, il suffit de se placer dans un espace de représentation virtuel dans lequel les éléments d'interpolation ont 20 tous les mêmes dimensions. Par exemple, si les longueurs de fissure discrétisées sont x E 0,1; 0,2; 0,4; 0,6; 0,8}, on prendra xvirtue/ E {1; 2; 3; 4; 5; 6; 7; 8; 9; 10;11}. Pour passer de la représentation virtuelle à la représentation 25 réelle, on utilise une bijection entre les deux repères (voir figure 10b). On utilise les mêmes fonctions d'interpolation que précédemment. La correspondance entre les deux espaces de représentation se fait en imposant les degrés de liberté filet el- = fil à chaque noeud. Puis, pour avoir une métrique qui évolue de manière « douce » entre chaque 30 élément, le procédé présenté précédemment est mis en oeuvre. Cela permet d'obtenir une correspondance lisse entre les deux espaces de représentation comme représenté sur la figure 10c. 3052891 27 Références : [1] Bueckner HGZ (1970), "A novel principle for the computation of stress intensity factors", Angew, Math. Mech , Vol. 50, p. 529-546. 5 [2] Rice J.R. (1972), "Some remarks on elastic crack-tip stress fields", Int. J. Solids and Structures, 8, 751-758. [3] Rice J.R. (1989), "Weight function theory for three-dimensional 10 elastic crack analysis", Fracture Mechanics : Perspectives and Directions (Twentieth Symposium), ASTM STP 1020, R.P. Wei and R.P. Gangloff, Eds., American Society for Testing and Materials, Philadelphia, pp. 2957. 15 [4] Lazarus V. (1997), "Quelques problèmes tridimensionnels de la mécanique de la rupture fragile", Thèse Université Paris 6. [5] W. Elber. "The significance of fatigue crack closure". ASTM STP, 486 : 230-242, 1971. 20 [6] PC. Paris, F. Erdogan, 1963, "A critical analysis of crack propagation laws". J Basic Eng 85, pp 528-534. [7] « Nasgro - Fracture Mechanics and Fatigue Crack Growth Analysis 25 Software - Reference Manual», Version 6.2, 2011. REVENDICATIONS1. Procédé d'estimation (100) du facteur d'intensité des contraintes (FIC) dans une pièce modélisée numériquement, dans le cadre d'une modélisation de propagation de fissure par fatigue, comprenant les étapes suivantes : (E2) : obtention, à partir d'une modélisation numérique de la pièce à analyser (20), d'une pluralité de valeurs simulées en différents points de la pièce à analyser (20), ladite pluralité de valeurs simulées comprenant, pour différents pas de propagation simulés d'une fissure tridimensionnelle qui apparait sur la modélisation numérique de la pièce, un ensemble des valeurs simulées du facteur d'intensité des contraintes (FIC K) en des points associés, la position de ces points, et des données relatives à la surface fissurée de la fissure, - (E3) : détermination, pour chaque pas de propagation et pour l'ensemble de valeurs simulées pour celui-ci, d'une part d'une valeur convertie (AK glob) d'amplitude du facteur d'intensité des ef f contraintes efficace global équivalent correspondant à celui d'une fissure plane à front droit et d'autre part d'une longueur (a) convertie de la fissure équivalente considérée, lesdites valeurs converties étant déterminées par égalisation de l'énergie dissipée dans la fissure tridimensionnelle de la modélisation numérique et de l'énergie dissipée dans la fissure d'un modèle standard de fissure plane à front droit, les énergies étant elles-mêmes déterminées en fonction des facteurs d'intensité des contraintes, - (E4) : Interpolation de valeurs converties d'amplitude de facteur d'intensité des contraintes efficace global équivalent entre deux longueurs de fissure converties (a) successives, - (E5) : Stockage des valeurs converties d'amplitude de facteur d'intensité des contraintes efficace global équivalent ainsi interpolées avec les longueurs de fissure associées. 305 2 891 29 2. Procédé d'estimation selon la revendication 1, dans lequel l'étape d'interpolation (E4) met en oeuvre une interpolation linéaire par morceaux. 5 3. Procédé d'estimation selon la revendication 1, dans lequel l'étape d'interpolation (E4) met en oeuvre une interpolation par minimisation de l'énergie de courbure. 4. Procédé d'estimation selon la revendication 1, dans lequel l'étape 10 d'interpolation (E4) comprend les sous-étapes suivantes : - (E41) interpolation à l'aide d'une interpolation linéaire, - (E42) interpolation à l'aide d'une interpolation par minimisation de l'énergie de courbure, les deux interpolations étant interchangeables, - (E43) calcul d'au moins une grandeur représentative d'une différence 15 (6r) de valeurs d'amplitude du facteur d'intensité des contraintes efficace global équivalent (FIC) entre les deux interpolations, ladite différence étant calculée strictement entre deux valeurs de longueurs de la fissure (a) correspondant à deux pas successifs, - (E44) comparaison de ladite grandeur représentative de la différence 20 (6r) avec un seuil prédéterminée (VSP), - (E45) si la grandeur représentative de la différence (ôr) est supérieure au seuil (VSP), génération d'une consigne d'instruction de calcù1 à partir de la simulation numérique de la pièce à analyser (10) des valeurs du facteur d'intensité des contraintes (FIC) et des valeurs de la position de 25 la fissure, entre les deux pas successifs. 5. Procédé d'estimation selon l'une quelconque des revendications précédentes, comprenant, préalablement à l'étape d'obtention (E2), une étape de calcul (El) mettant en oeuvre la simulation par 30 éléments finis à des pas successifs de l'évolution de la fissure dans la pièce à analyser (20), ladite simulation étant effectuée par l'unité de traitement (12). 3052891 30 6. Procédé d'estimation selon l'une quelconque des revendications précédentes, dans lequel les données numériques récupérées à l'étape d'obtention (E2) des données sont obtenues par une simulation par éléments finis ou éléments finis étendus. 5 7. Procédé d'estimation selon l'une quelconque des revendications précédentes, dans lequel, pour l'étape de détermination (E3), le premier sous-ensemble de données permet de connaitre l'amplitude efficace du taux de restitution d'énergie AG par les relations 10 suivantes et leur équivalence en matière d'amplitude : 1 G = -* + KID 1 + Kii E* = E en contraintes planes E* = E/(1 - u2) en déformations planes (21Kef f(s))2 = ,dGeff (s). E* Où KI, KII et Km sont, respectivement, les coefficient de FIC correspondant aux modes d'ouverture de la fissure, de cisaillement dans le plan et de cisaillement anti-plan, E est le module de Young, E* est le module de Young équivalent, p est le module de cisaillement, v est le 15 coefficient de Poisson, le deuxième sous-ensemble de données permet de connaitre l'incrément de surface fissuré dSurf, ces formulations permettent d'obtenir la valeur de l'énergie dissipée : Energie dissipée =f il Ge f f(s).dSurf (s). dl front de fissure 20 dans lequel cette énergie dissipée est égalisée avec l'énergie dissipée de la fissure d'un modèle standard de fissure plane à front droit, s'exprimant sous la forme : Energie dissipéegl°b = AGee fi fb.dSurfglob .1 I,front Où la notation glob est relative à ce modèle standard, 3052891 31 dans lequel, grâce à une modélisation par une loi de Paris avec les paramètres de fermeture de fissure d'Elber, le lien entre dSurfeb et AGef fl°b est effectué dSurfeb C. GiAGglf o f b E*) dN e - où C et n sont des coefficients de Paris, 5 dans lequel, en égalisant les deux énergies dissipées, on obtient 1 . dSurf (s). dl n+1 Énergie dissipée n+1 glob f front de fissure ef f (s) AGeff = C. (E*YiLfront C. (E*PLfront dans lequel, au moyen des relations entre LGeff et zKeff, on détermine AKglob et dans lequel on détermine Surfe'', c'est-à-dire la eff longueur a. 10 8. Procédé d'évaluation (200) de la durée de vie d'une pièce (20) modélisée numériquement, dans lequel on met en oeuvre un calcul de la durée de vie de la pièce (20) faisant intervenir les valeurs converties d'amplitude de facteur d'intensité des contraintes efficace interpolées selon un procédé conforme à l'une des revendications 1 à 7. 15 9. Système comprenant une unité de traitement (12) comprenant des moyens de calcul (14) et une mémoire (16), l'unité étant configurée pour mettre en oeuvre un procédé d'estimation selon l'une quelconque des revendications 1 à 7, ou un procédé d'évaluation de durée de vie 20 selon la revendication 8. 10. Produit programme d'ordinateur configuré pour être mis en oeuvre par un système selon la revendication 9, et comprenant des instructions pour provoquer la mise en oeuvre d'un procédé d'estimation selon l'une 25 quelconque des revendications 1 à 7, ou un procédé d'évaluation de durée de vie selon la revendication 8. n