Cálculo de los factores de intensificación de tensiones generalizados en uniones adhesivas metal-composite

MATERIALES COMPUESTOS 05 627 Cálculo de los factores de intensificación de tensiones generalizados en uniones adhesivas metal-composite A. Barroso, P. Toro, V. Mantiþ, F. París y J. Ortiz. Grupo de Elasticidad y Resistencia de Materiales, Universidad de Sevilla, España. RESUMEN Las uniones adhesivas entre materiales metálicos y compuestos generan puntos críticos donde la geometría y/o las propiedades de los materiales cambian bruscamente. El estado tensional que predice la Teoría de la Elasticidad lineal en estos puntos es singular (el valor de las tensiones no está acotado). La caracterización completa de estos estados tensionales no acotados se consigue con el cálculo de los exponentes característicos (órdenes de singularidad de tensiones), las funciones características y los Factores de Intensificación de Tensiones Generalizados (FITG). En este trabajo se presenta un procedimiento para el cálculo de los FITG basado en un ajuste por mínimos cuadrados y el empleo de modelos numéricos de las uniones utilizando el Método de los Elementos de Contorno (MEC). Se presentan comparaciones para problemas con solución analítica y ejemplos de aplicación. 1. INTRODUCCIÓN Las uniones adhesivas entre materiales de distinta naturaleza (p.e. metal-composite) generan una serie de puntos críticos, en el entorno de los cuales la teoría de la elasticidad lineal predice tensiones no acotadas. Este estado singular de tensiones se debe tanto a las geometrías de esquinas re-entrantes que aparecen como a la discontinuidad en las propiedades elásticas de los materiales. En el caso particular de las uniones adhesivas entre láminas de metal y laminados de material compuesto, estos puntos críticos se identifican con claridad en la Fig.1, que representa una unión adhesiva a solape simple entre una lámina de aluminio y un laminado de carbono. Metal adherent 90o 45o 0o Aluminio Laminado de carbono 90o 45o 0o Adhesivo Metal adherent 90o 45o 0o Aluminio Laminado de carbono 90o 45o 0o Adhesivo 1 i N T 0 T N T i T i-1 x 1 x 2 1 i N T 0 T N T i T i-1 x 1 x 2 Fig. 1 - Esquema de un encolado metal-composite y el sistema de referencia local. Tomando un sistema de referencia polar (r, T ) centrado en el vértice de la esquina (ver Fig.1), el estado singular de tensiones en el entorno de estos puntos se puede aproximar mediante un desarrollo en serie de la forma: 628 MATERIALES COMPUESTOS 05 )(),( 1 1 TTV O ijk K k k ij f r K rk ¦  |, (i,j=r, T  (1) donde O k son los exponentes característicos (0< O k<1), fijk( T ) son funciones suaves dentro de cada material que contienen la dependencia angular de la componente tensional V ij asociadas al modo singular O k y Kk son los Factores de Intensificación de Tensiones Generalizados (FITG en adelante) que miden el peso de cada término dentro del desarrolllo en serie asintótico. La caracterización del estado tensional singular es de vital importancia si se considera que dicho estado tensional puede ser el responsable de la aparición de una grieta o fallo cuya propagación arruine globalmente la estructura. Al igual que la Mecánica de la Fractura Elástica Lineal (MFEL) puede predecir con cierta fiabilidad la propagación de una fisura cuando los FIT (Kk) alcanzan un cierto valor admisible (la tenacidad a fractura del material), parece razonable considerar la posibilidad de que el estado singular de tensiones, en estas configuraciones de esquinas multimateriales (Fig. 1) controlado por los FITG (Kk), juegue un papel importante en la iniciación del fallo. Para analizar la viabilidad de este enfoque a través de la Mecánica de la Fractura (sin grietas) es necesario disponer de las herramientas que permitan el cálculo de los parámetros involucrados en la ecuación (1): O k, fijk( T ) y Kk. En el apartado 2 se revisará la determinación de O k y fijk( T ), en el apartado 3 se abordará el cálculo de Kk., por último en el apartado 4 se verificará la eficiencia del procedimiento propuesto con trabajos previos y se abordarán nuevos problemas de encolados metal-composite. 2. CARACTERIZACIÓN DE LOS MODOS SINGULARES DE TENSIÓN Los exponentes característicos O k y las funciones características fijk( T ) dependen exclusivamente de la geometría y condiciones de contorno locales de la esquina así como de las propiedades elásticas de los materiales. En la actualidad y bajo las siguientes hipótesis: materiales homogéneos elásticos lineales, adhesión perfecta entre los materiales, condiciones de contorno homogéneas y estados planos (tensión plana, deformación plana o deformación plana generalizada) se dispone de herramientas versátiles y eficientes para el cálculo de valores precisos de los O k, ver Dempsey y Sinclair (1979, 1981), Pageau et al (1994), Ting (1997), Mantiþ et al (1997) y Poonsawat et al (1998, 2001) entre otros. En este trabajo se emplea una herramienta semianalítica desarrollada por los mismos autores, Barroso et al (2003), que proporciona los órdenes de singularidad de tensiones en esquinas multimateriales. 3. FACTORES DE INTENSIFICACÓN DE TENSIONES GENERALIZADOS En el cálculo de los FITG sí influye la geometría y condiciones de contorno globales del problema, por lo tanto es necesario el uso de modelos numéricos y/o técnicas experimentales para el cálculo de los mismos. Adicionalmente, en la mayoría de los casos, es necesario un postprocesado de los datos obtenidos, quedando la precisión final de los resultados de los FITG estrechamente ligada a la "calidad" de ambos procesos. Por ejemplo, dado el carácter asintótico de la representación del estado tensional, MATERIALES COMPUESTOS 05 629 ecuación (1), FITG extraídos de la solución lejos del vértice de la esquina (y por lo tanto fuera del alcance de representación de dicha aproximación asintótica) pueden no ser representativos del estado tensional real, mientras que FITG extraídos de la solución excesivamente cerca de la esquina pueden estar afectados por errores numéricos de la discretización del modelo. Existen técnicas que tratan de evitar estas dificultades, mediante la evaluación de integrales realizadas en el entorno de la esquina que son teóricamente independientes del camino, ver por ejemplo los trabajos de Qian y Akisanya (1999), Wu (2001), Banks-Sills y Sherer (2002) y Cisilino y Ortiz (2005) entre otros. Una revisión de los procedimientos de cálculo de los FITG se puede encontrar en Helsing y Jonsson (2002a). En el caso particular de los FITG se echa en falta un mayor número de problemas con resultados fiables que sirvan de referencia para el ajuste de nuevos métodos de cálculo. Comentarios muy interesantes acerca de la validez de los resultados numéricos publicados en la literatura se pueden leer en Helsing y Jonsson (2002b). En el presente trabajo se ha puesto a punto un sencillo procedimiento de cálculo de los FITG que está basado en un ajuste por mínimos cuadrados. Los modelos numéricos utilizados se han realizado utilizando un código de desarrollo propio que implementa el Método de los Elementos de Contorno (MEC en adelante) dada su mayor precisión en problemas de esta índole (tensiones singulares y resultados en los contornos e intercaras). 3.1 Cálculo de los FITG mediante un ajuste por mínimos cuadrados. De forma similar a la utilizada por Yang y Munz (1995) quienes utilizaron Elementos Finitos, en el presente trabajo se utiliza un procedimiento de cálculo de FITG basado en un ajuste por mínimos cuadrados. Se utilizan los valores numéricos de los desplazamientos ur(r, T ) y u T (r, T ) extraídos del modelo numérico realizado con el MEC (unum) y los desplazamientos de la representación asintótica de la fórmula (1) (useries). Se define una función de error J, ver (2), en la que las únicas incógnitas son los FITG, Kk.  >@ 2 101 11 ),...,,,(),(,..., ¦¦¦  AN j M n kjn series jn num kKKruruKKJ D DD TT (2) En la función J se incorporan tres sumatorios: el sumatorio en D considera sólo una componente de los desplazamientos si $ =1 y las dos componentes si $ =2, el sumatorio en j marca el número de aristas que intervienen en el ajuste (para N materiales existen N+1 caras/intercaras donde hay resultados del modelo MEC) y por último el sumatorio en n establece el conjunto de nodos en cada cara/intercara que se incorpora al ajuste. El conjunto de valores de Kk que minimiza la función J se obtiene de la resolución del siguiente sistema de ecuaciones lineales: 0 ),...,( 1 w w i k K KKJ , (i=1,...,k) (3) La solución de este sistema de ecuaciones depende del valor de los índices $ , N y M. Algunos comentarios respecto de esta dependencia se verán en el apartado siguiente. 630 MATERIALES COMPUESTOS 05 4. EJEMPLOS Para la evaluación de la bondad del método implementado se han tomado dos problemas sencillos resueltos en trabajos previos: una placa a tracción de un material isótropo con una entalla a 90º y una esquina bimaterial isótropa libre-libre. Si bien el procedimiento se puede aplicar a un espectro de configuraciones más amplio, dos son la razones por las cuales no se han tomado ejemplos de mayor complejidad (presencia de materiales de comportamiento anisótropo, mayor número de materiales, etc.): a) una vez que se dispone de la información de O k y fijk( T ) para un problema concreto, independientemente del número de materiales, naturaleza de los mismos y las condiciones de contorno, el método de ajuste se aplica de idéntica manera, b) las expresiones analíticas (si están disponibles) de problemas más complejos son muy poco manejables y no aportan dificultad añadida al método presentado. 4.1 Placa a tracción con una entalla a 90º. El primer problema que se ha resuelto es un problema sencillo de una placa de material isótropo que presenta una entalla a 90º en uno de sus laterales sometida a una tracción uniforme en las caras superior e inferior. Se han empleado valores de a/w=0.5 y h/w=1, siendo a=5cm. El modelo MEC tiene 709 nodos y 709 elementos lineales, con un tamaño de 0.1 cm en las caras horizontales y verticales. En el vértice de las caras inclinadas el primer elemento tiene un tamaño de 10-7 cm, los cuales crecen progresivamente, con un factor de 1.5, hasta alcanzar un tamaño de 0.1cm que permanece constante en el resto de la cara, dando un total de 100 elementos en dichas caras. a 2h 270º w x y V V a 2h 270º w x y V V Fig. 2Placa a tracción. Esta configuración de cargas simétrica, que genera un modo I de apertura en el entorno del vértice, ha sido estudiado con detalle por Helsing y Jonsson (2002a) de donde se ha tomado el valor de referencia para evaluar los errores: KI=4.295886967699. Utilizando ur y u T ( $ =2) y las dos caras inclinadas (N=2) para la evaluación de J en (2), se representa en la Fig. 3 el valor absoluto del error relativo entre la solución obtenida y el valor de referencia para dos aproximaciones diferentes: en la Fig. 3a) se han tomado 3 términos del desarrollo asintótico, k=3 en (1), para evaluar KI, mientras que en la Fig. 3b) se han tomado los primeros 7 términos. Se han representado todas las combinaciones que resultan de tomar un conjunto de nodos consecutivos, siendo el eje derecho el nodo inicial desde el que se realiza el ajuste y el eje izquierdo el nodo final (la mitad inferior de la gráfica carece de sentido). En la Fig. 3 se aprecian los efectos de borde, con errores significativos a medida que el ajuste se realiza utilizando un pequeño conjunto de nodos muy alejado del vértice (zona izquierda de las gráficas) y los errores asociados a la discretización del modelo cuando el ajuste se realiza tomando pocos nodos muy cerca del vértice (zona derecha de la gráfica). A medida que se toman más términos del desarrollo en serie para representar la solución asintótica del problema en el vértice, estos efectos de borde son menos acusados, por lo que de manera evidente se aprecia una significativa reducción de los MATERIALES COMPUESTOS 05 631 errores en la zona izquierda de la gráfica entre a) 3 términos y b) 7 términos. En cualquier caso en ambos ejemplos se observa para un amplio rango de conjuntos "razonables" de nodos (evitando grupos de nodos pequeños muy cerca y muy lejos del vértice) que los resultados son, desde un punto de vista ingenieril, bastante razonables (<4%). Los errores relativos obtenidos con el método desarrollado en Ortiz et al (2005) utilizando un método de integral conservativa se encuentran, para este problema, siempre por debajo del 1%. 20 40 60 80 NODO FINAL 20 40 60 80 100 NODO INICIAL 0 1 2 3 4 5 ERROR 20 40 60 80 NODO FINAL 0 1 2 3 4 5 ERROR 20 40 60 80 NODO FINAL 20 40 60 80 100 NODO INICIAL 0 1 2 3 4 5 ERROR 20 40 60 80 NODO FINAL 0 1 2 3 4 5 ERROR a) b) Fig. 3Errores relativos (en %) con: a) 3 términos y b) 7 términos. 4.2 Esquina bimaterial isótropa. La herramienta desarrollada se ha verificado también en una esquina bimaterial isótropa estudiada por Qian y Akisanya (1999) y que se muestra en la Fig. 4. La solución de referencia dada por Qian y Akisanya se obtuvo para unos parámetros de Dundurs D =0.8 y E =0.2 lo que hace que el material de la parte superior de la Fig. 4 sea unas 10 veces más rígido que el situado en la parte inferior. Las expresiones utilizadas por Qian y Akisanya para tensiones y desplazamientos son: V 90º 120º h h h V 90º 120º h h h Fig. 4Esquina bimaterial. )(),( 1 1 TTV O m ijk N k k m ij frHr k ¦  )(),( 1 TT O m ik N k k m igrHru k ¦ (4) ),,,,(· 21 1 TTOEDV O kkk ahH k  (5) donde m es el material, Hk es el FITG asociado a O k, ak es el FITG adimensional (igual a Hk si h=1 y V =1) y f( T ) y g( T ) son funciones conocidas, ver Qian y Akisanya (1999). Utilizando integrales de contorno independientes del camino, finalmente transformadas en integrales de dominio, Qian y Akisanya obtienen para los tres primeros modos (sólo el primero de ellos singular) los siguientes valores de O k y ak: O 1=0.6747 (a1=0.6301), O 2=1.1637 (a2=-0.3671) y O 3=1.5938 (a2=0.5443). Se ha realizado un modelo MEC del problema con un tamaño de elementos 0.025h lejos 632 MATERIALES COMPUESTOS 05 del vértice y 10-8h en el vértice usando elementos lineales. Los órdenes de singularidad de tensiones se calcularon utilizando la herramienta desarrollada en Barroso et al (2003) obteniéndose: O 1=0.673473, O 2=1.167477, O 3=1.589147 muy similares a los de Qian y Akisanya. El cálculo de los FITG se ha llevado a cabo utilizando (2) y (3) con ur y u T (A=2), las tres aristas que confluyen al vértice (N=3) y un grupo de nodos "razonable", es decir, evitando nodos muy cercanos y muy lejanos a la esquina. Considerando tres términos del desarrollo asintótico y un grupo de nodos entre el 15 y el 55, que se corresponde con la zona más "plana" de la representación, se obtuvo a1=0.673688 que difiere un 6.92% del resultado de Qian y Akisanya. Ortiz et al (2005) con una técnica de integral conservativa, han obtenido un valor medio de a1=0.67829, que difiere un 7.65% del resultado de Qian y Akisanya y un 0.68% del valor obtenido en el presente trabajo. 4.3 Unión adhesiva metal-composite. Tras la verificación de la bondad del método desarrollado, se ha analizado una unión encolada metal-composite como la que se muestra en la Fig. 5. Se trata de una unión de doble solape entre una chapa de aluminio (de 3.2 mm de espesor) y un laminado unidireccional de fibra de carbono de 8 capas a 0º (de 1.6 mm de espesor) utilizando un adhesivo epoxy que tras el curado en prensa a 180ºC presenta un espesor de 0.1 mm. La unión se va a someter a un estado de cortadura a tracción. El modelo MEC de la unión utiliza el plano de simetría y se presenta parcialmente (sólo aparece la zona de solape) en la Fig. 5 incluyendo ya la geometría deformada (amplificada x50). En las esquinas multimateriales del modelo se ha utilizado un tamaño de 10-8 mm y un factor de crecimiento de 1.5, de forma que el modelo completo tiene un total de 1484 nodos. 50 mm 50 mm 12.5 mm Aluminio espesor=3.2 mm [0º]8espesor = 1.6 mm [0º]8espesor = 1.6 mm V Espesor de la capa de adhesivo = 0.1 mm 0.00 0.50 1.00 1.50 2.00 2.50 3.00 3.50 42.00 47.00 52.00 57.00 62.00 geometría deformada (X 50) geometría indeformada CFRP 8 láminas 0º Aluminio Fig. 5Geometría y modelo MEC de la unión adhesiva metal-composite. Con la herramienta desarrollada en Barroso et al (2003) se han obtenido los órdenes de singularidad de tensiones O k y las funciones f( T ) y g( T ) de todas las esquinas multimateriales que aparecen en este modelo. Los FITG se han calculado utilizando las expresiones (2) y (3) una vez resuelto el modelo numérico. En la Fig. 6 se muestran los valores del primer modo singular para dos de las esquinas del modelo, mostrándose en trazo continuo la geometría indeformada y con puntos la deformada (x10). MATERIALES COMPUESTOS 05 633 1.50 1.55 1.60 1.65 1.70 1.75 1.80 47.00 47.50 48.00 48.50 49.00 49.50 50.00 50.50 51.00 Esquina 1 Esquina 2 ALUMINIO ADHESIVO LAMINADO DE CARBONO Deformada (x10) Indeformada Esquina 1: 1O 1 = 0.236764 K1 = 0.01576 Esquina 2: 1O 1 = 0.0210332 K1 = -0.003086 Fig. 6Detalle de dos esquinas del modelo y los resultados obtenidos. La esquina 1 la forma el laminado de carbono en contacto con el adhesivo, mientras que la esquina 2 es la que forma el rebose de adhesivo con la chapa de aluminio. Los valores de K que se presentan en la Fig. 6 están normalizados, de forma que VTT _ T =0º=K/(2 S r)1-O. En el caso de la esquina 1 y dado el pequeño espesor de adhesivo (0.1 mm) se ha realizado el ajuste entre r=10-6 y r=0.03 mm (medidos desde la esquina) de forma que se evitan aquellos nodos muy cerca de la esquina (entre 10-8 y 10-6 mm) y aquellos nodos muy lejos del campo asintótico (>1/3 del espesor de adhesivo). Tomando otros grupos de nodos en el entorno (r=10-6, r=0.03 mm) los valores de K obtenidos son razonablemente similares (<4%). En el caso de la esquina 2, cuatro grupos de nodos diferentes entre 10-7 y 0.33 mm han dado variaciones de K inferiores al 2%. En ambos casos se han tomado los tres primeros términos del desarrollo asintótico además del término de giro como sólido rígido, que es un término independiente de r, con O k=1 en la fórmula (1). 5. CONCLUSIONES Y DESARROLLOS FUTUROS En el presente trabajo se ha desarrollado un procedimiento de cálculo de los Factores de Intensificación de Tensiones Generalizado (FITG) en esquinas multimateriales anisótropas basado en un ajuste de mínimos cuadrados. Se han utilizado en dicho ajuste los desplazamientos de la solución asintótica y la solución numérica de un modelo de Elementos de Contorno. La herramienta se ha mostrado eficiente y robusta en problemas ya resueltos por otros autores y sencilla de uso en problemas nuevos. El procedimiento desarrollado, junto con otro trabajo previo de los mismos autores (que permite la obtención de los órdenes de singularidad de tensiones en esquinas multimateriales anisótropas) completa una potente herramienta para la completa caracterización de esquinas multimateriales anisótropas. Con el estado tensional singular caracterizado en detalle se puede definir a medida un programa de ensayos experimentales para evaluar si el valor de los FITG juega algún papel en una posible propuesta de fallo de este tipo de uniones. AGRADECIMIENTOS Estudio financiado por el Ministerio de Ciencia y Tecnología (MAT2003-03315). 634 MATERIALES COMPUESTOS 05 REFERENCIAS BANKS-SILLS, L. AND SHERER, A. (2002) A conservative integral for determining stress intensity factors of a bimaterial notch, Int. J. Fracture 115, 1-26. BARROSO, A., MANTIý, V. AND PARÍS, F. (2003) Singularity analysis of anisotropic multimaterial corners, Int. J. Fracture, 119/1, 1-23. CISILINO, A.P. Y ORTIZ, J. (2005) Boundary element analysis of three-dimensional mixed-mode cracks via the interaction integral, Comput. Method Appl. M. 194, 935956. DEMPSEY, J.P. AND SINCLAIR, G.B. (1981)On the singular behaviour of a bimaterial wedge, J. Elasticity, 11/3, 317-327. DEMPSEY, J.P. AND SINCLAIR, G.B. (1979) On the stress singularities in the plane elasticity of the composite wedge, J. Elasticity, 9/4, 373-391. HELSING, J. AND JONSSON, A. (2002b) On the accuracy of benchmark tables and graphical results in the applied mechanics literature, ASME J. Appl. Mech, .69/1, 88-90. HELSING, J. AND JONSSON, A. (2002a) On the computation of stress fields on polygonal domains with V-notches, Int. J. Numer. Meth. Eng., 53, 433-453. MANTIý, V. PARÍS, F. AND CAÑAS, J. (1997) Stress singularities in 2D orthotropic corners, Int. J. Fracture 83, 67-90. ORTIZ, J.E., MANTIý, V. AND PARÍS, F. (2005) A domain-independent integral for computation of Stress Intensity Factors along three-dimensional crack front and edges by BEM, (Submitted for publication). PAGEAU, S.S., JOSEPH, P.F. AND BIGGERS, JR., S.B. (1994) The order of stress singularities for bonded and disbonded three-material junctions, Int. J. Solids Struct., 31, 2979-2997. POONSAWAT, P., WIJEYEWICKREMA, A.C. AND KARASUDHI, P. (1998) Singular stress fields of an anisotropic composite wedge with a frictional interface, ASCE 12th Engineering Mechanics Conference, La Jolla, San Diego, CA, 578-581. POONSAWAT, P., WIJEYEWICKREMA, A.C. AND KARASUDHI, P. (2001) Singular stress fields of angle-ply and monoclinic bimaterial wedges, Int. J. Solids Struct., 38, 91-113. QIAN, Z.Q. AND AKISANYA, A.R. (1999) Wedge corner stress behaviour of bonded dissimilar materials, Theor. Appl. Fract. Mec., 32, 209-222. SEWERYN, A. AND MOLSKI, K. (1996) Elastic stress singularities and corresponding generalized stress intensity factors for angular corners under various boundary conditions, Eng. Frac. Mech., 55/4, 529-556. TING, T.C.T. (1997) Stress singularities at the tip of interfaces in polycristals, Damage and Failure of Interfaces (Ed. by Rossmanith). Balkema, Rotterdam, 75-82. WU, K.C. (2001) Near-tip field and the associated path-independent integrals for anisotropic composite wedge, The Chinese Journal of Mechanics, 17/1, 21-28. YANG, Y.Y. AND MUNZ, D. (1995) Stress intensity factor and stress distribution in a joint with an interface corner under thermal and mechanical loading, Computers & Structures, Vol 57. No.3, pp:467-476.