Caracterización numérica de esquinas multimateriales 3-D

MATERIALES COMPUESTOS 05 869 CARACTERIZACIÓN NUMÉRICA DE ESQUINAS MULTIMATERIALES 3-D J.E. Ortiz, A. Barroso, V. Mantiþ y F. París Grupo de Elasticidad y Resistencia de Materiales. Escuela Superior de Ingenieros. Universidad de Sevilla. RESUMEN Se presenta un método eficiente para calcular los Factores de Intensidad de Tensiones Generalizados (FITG) asociados a términos de tensión singulares y no-singulares en aristas de esquinas 3-D formadas por la unión de varios materiales. Este método está basado en el concepto de la integral-H. Se modelan 2 esquinas bimateriales en las cuales se calculan los FITGs asociados a los 2 primeros términos del desarrollo asintótico de la solución elástica lineal en el entorno de la arista. 1. INTRODUCCIÓN Las esquinas formadas por la unión de materiales, generalmente 2 o 3 materiales, son frecuentes en uniones adhesivas de materiales homogéneos o compuestos. El empleo de este tipo de materiales se ha incrementado en forma notable durante las últimas décadas. Sin embargo, el uso de este tipo de uniones en estructuras primarias es limitado por la falta de conocimiento sobre los mecanismos de fallo asociados a la unión adhesiva y criterios de fallo macroscópicos que permitan predecir en forma fiable su resistencia. Los criterios de fallo de esquinas reentrantes que se investigan usualmente están basados en el concepto del FITG, véase por ejemplo Seweryn (1994), Qian y Akisanya (1998) y Labossiere y Dunn (2001) para esquinas en medios homogéneos 2-D, bimateriales 2-D y 3-D, respectivamente. Por todo ello el desarrollo de las metodologías para calcular el FITG de forma eficiente sobre todo en una esquina 3-D es de interés actual. No todos los métodos para el cálculo del FIT en grietas pueden extenderse para el cálculo del FITG en esquinas, y aunque algunos puedan extenderse, no necesariamente permiten calcular los FITG asociados a términos de tensión de orden superior (Chen et al. 2004). La integral-H que permite calcular el FITG en esquinas 2-D, se deriva a partir del segundo teorema de reciprocidad de Betti y de la condición de biortogonalidad entre las familias de autovalores positivos y negativos. Este método ha sido aplicado para el cálculo del FITG en esquinas en medios homogéneos 2-D (Szabó y Babuška, 1991) o bimateriales 2-D (Banks-Sills et al. 2002). Para problemas de grietas 3-D una expresión en forma de una integral de superficie fue deducida por Meda y colaboradores (Meda et al. 1998) a partir de la integral-H. Todas las metodologías basadas en las integrales conservativas necesariamente requieren de uso de un campo auxiliar. Por su parte, la integral-H requiere solamente evaluar las variables naturales como desplazamientos y tensiones, lo cual representa una inherente ventaja de este método y a su vez puede implicar una mejor precisión en comparación con otras integrales de conservación que requieren gradientes de desplazamientos y/o de tensiones. Por otro lado, el Método de los Elementos de Contorno (MEC) es adecuado para la evaluación de las integrales de 870 MATERIALES COMPUESTOS 05 conservación. Los desplazamientos y tensiones en puntos internos pueden ser obtenidos por el MEC con alta precisión a partir de las representaciones integrales de contorno en comparación con otras técnicas como elementos finitos. El presente trabajo ha sido motivado por la necesidad de desarrollar una herramienta numérica para calcular de forma directa y precisa los FITGs asociados a los términos de tensiones singulares y algunos de orden superior a lo largo de la arista en una esquina multimaterial 3-D. Para este propósito se presenta una nueva formulación para la integral-H. La expresión independiente del dominio para esta integral es presentada y ha sido implementada usando las ventajas de la solución del MEC. Se han modelado 2 placas bimateriales conteniendo una esquina reentrante cada una de ellas. Se presentan los resultados del FITG asociados a los 2 primeros términos del desarrollo de la solución elástica en la zona circundante a la arista. 2. CAMPO ELÁSTICO EN UNA ESQUINA 3-D Considérese un sólido 3-D compuesto de dos materiales isótropos elásticos lineales distintos adheridos perfectamente (cumpliéndose las condiciones de equilibrio de tensiones y compatibilidad de desplazamientos en la interfaz bimaterial). El sólido contiene una esquina recta con un ángulo sólido diedro de valor I  I , ver Figura 1. Se define el sistema de coordenadas cilíndricas (r, T ,z) asociado al sistema de coordenadas cartesianas (x,y,z), ubicado sobre un punto fijo de la arista tal como se muestra en la Figura 1. La expresión general del desarrollo en serie del campo de desplazamiento cercano a la arista de la esquina puede ser escrita en la siguiente forma (Andersson et al., 1995 y Costabel et al., 2004): ),0()()(),,( 10 ¦¦ tt dd pj qp m pj j p pj mqprzKzr OOTT O fu Syyyjjjjy (1) d d d. D D. D D D D D d d d d d d d donde um es el vector desplazamiento, el índice m (m=1,2) es el número del material y O p+j denota el autovalor de la arista que determina el nivel de singularidad y su correspondiente autofunción )( T m pj f. El pj K representa el FITG de la arista, asociado al autovalor O p+j, que depende de la geometría y de la carga externa aplicada. Los autovalores O p y las autofunciones )( 0 T m p f son extraídas de los problemas 2-D de la esquina bimaterial. Los FITGs y las autofunciones para 0 j son denotadas usualmente como p K y ),( T m p frespectivamente. x Vértice Arista y z O I 1 I 2 Figura 1: Esquina bimaterial 3-D. MATERIALES COMPUESTOS 05 871 3. INTEGRAL INDEPENDIENTE DEL DOMINIO 3-D La integral-H independiente del camino para una esquina 2-D bimaterial establece una relación entre dos soluciones lineales elásticas del problema. Usando esta integral y la condición de biortogonalidad entre las familias de las autofunciones asociadas a los autovalores Op positivos y negativos, respectivamente, se demuestra que el FITG asociado a Op puede ser determinado. Para el caso 3-D la condición de ortogonalidad no siempre prevalece. Sin embargo, una integral independiente de dominio para el cálculo del FITG a lo largo de la arista de la esquina 3-D ha sido deducida por Ortiz et al. (2005). La expresión final de esta integral adopta la siguiente forma:  >@ ³ ³  cc  L Lz Vjiijiij n p zd)z(q dVquu )(zK , '' VV , (2) donde ui y Ti son los vectores desplazamiento y tensión del problema bajo estudio, y i u' y i T' los correspondientes vectores desplazamiento y tensión de una solución auxiliar. La evaluación de las integrales se lleva a cabo en un sistema de coordenadas locales (x’, y’, z’) con origen en la posición global en un punto fijo zn sobre la arista donde el FITG  np zK es evaluado, y con los ejes x’, y’ y z’ paralelos a los ejes x, y y z de la Figura 1, respectivamente. El volumen V es un dominio cilíndrico de longitud 2L alrededor de la arista conteniendo el punto zn. La función auxiliar q está definida en el dominio V como el producto de qr(r) y qz(z’). La primera de estas funciones tiene valor unitario sobre la arista, cero sobre la superficie cilíndrica del volumen V y es suave en el resto del volumen, mientras que la segunda es una función suave entre -L y +L y cero en estos puntos. La aproximación vía (2) para obtener el FITG asociado a Op resulta eficiente, si bien se requiere que los llamados términos de “sombra” de carácter tridimensional de (1) no afecten la aplicación de la condición de biortogonalidad entre las familias de autofunciones 2-D por lo que se debe verificar que 0<Op<Osing.+1 donde Osing. representa el menor exponente singular de Op que aparece en la expansión de la serie asintótica de la solución elástica en la vecindad de la arista. 4. EJEMPLOS NUMÉRICOS: ESQUINAS BIMATERIALES El procedimiento para obtener  np zK vía la expresión (2) ha sido implementado en la etapa de posprocesamiento de un código del MEC 3-D, los detalles pueden seguirse en Ortiz et al. (2005). Se demuestra la eficiencia de la metodología presentada evaluando los dos primeros FITGs asociados a Op en 2 esquinas bimateriales que se detallan a continuación. 4.1 Descripción de los problemas Se propone calcular los dos primeros FITGs asociados a Op en una esquina bimaterial y estudiar la incidencia del cambio del ángulo sólido de la esquina y las propiedades mecánicas sobre el FITG. Para ello se han considerado 2 probetas cuyos ángulos sólidos entre la interfaz bimaterial y la superficie libre son (90º, 120º) y (120º, 120º), que de aquí en adelante se denominaran probetas A y B, respectivamente. Las dimensiones y condiciones de contorno de las probetas se muestran en las Figuras 2 (a) 872 MATERIALES COMPUESTOS 05 y (c). Los módulos elásticos de los materiales son: E1=208GPa y E2=22GPa, y los coeficientes de Poisson Q1=0.3 y Q2=0.33. Los parámetros de Dundurs D y E del bimaterial son 0.8 y 0.2, respectivamente. Para esta combinación de materiales y ángulos sólidos, los dos primeros exponentes Op en la arista de la probeta A son 0.67347 y 1.16748 mientras que estos valores correspondientes a la arista de la probeta B son 0.60619975 y 1.13610937, respectivamente. Adicionalmente, se ensayó una segunda combinación de materiales usando los modelos anteriores. Las propiedades mecánicas son: E1=70GPa y E2=3GPa, y Q1=0.33 y Q2=0.38 (D=0.915089 y E=0.17455). Con esta nueva combinación de materiales los Op son 0.61262 y 1.24756 en la arista de la probeta A, y 0.58078097 y 1.23198660 en la probeta B. Los modelos del MEC consisten en 227 elementos cuadráticos y 1003 nodos y se muestran en las Figuras 2 (b) y (d). Figura 2: (a) Representación esquemática de la probeta A. (b) Malla en elementos de contorno de la probeta A. (c) Representación esquemática de la probeta B. (d) Malla en elementos de contorno de la probeta B. L 120º 120º Material 1 Material 2 L L 0.7L (c) (d) L L 120º 90º Material 1 Material 2 L 0.7L (a) (b) V V MATERIALES COMPUESTOS 05 873 4.2 Resultados y discusión Con el propósito de estimar la precisión de los resultados obtenidos, inicialmente se calculó el FITG asociado al término de tensión singular en la probeta A bajo la condición de deformación plana con D=0.8 y E=0.2. Estos valores se presentan en la Tabla 1 para 3 dominios de integración conjuntamente con la solución publicada por Barroso et al. (2004) usando el método de ajuste por mínimos cuadrados. La diferencia relativa porcentual máxima que se alcanza entre ambas soluciones es 0.8%. Se define el factor adimensional  tzap/ igual a )/()/( 1p LtzK p O V , donde V es la carga externa aplicada y L es una longitud característica de la probeta. Los valores de  tza / 1 a lo largo de la arista para las probetas A y B se presentan en las Figuras 3 y 5, respectivamente, mientras que los valores de  tza / 2 en las Figuras 4 y 6, respectivamente. Obsérvese que )/( 1tza es casi constante a lo largo de la arista para ambas combinaciones de materiales en ambas probetas, y de valor cercano al correspondiente en deformación plana (ver Tabla 1 y Fig. 3, para la probeta A con D=0.8 y E=0.2). En cambio )/( 2tza presenta un comportamiento casi constante para pequeños valores de z/t y tiende a f cuando z/to1, ver Figuras 4 y 6. Según Leguillon y Sanchez-Palencia (1999) el FITG )/( 1tzK cerca del vértice de una esquina puede verse afectado por el exponente de singularidad /1 que ocurre en dicho vértice, en general, el FITG )/( tzK p será afectado por todos los /q con q d p (Hallström et al., 1997). Cuando algún /1 es menor que Op, el )/( tzK p tiende a valores infinitos para posiciones cercanas al vértice o en caso contrario a cero. Debido a que  tza / 1 es casi constante a lo largo de la arista se puede suponer que 11 O #/ , ver Figura 3. En cambio el comportamiento de la curva  tza / 2 está controlado en mayor medida por la diferencia (/1-O2). Este valor es -0.494 y -0.635 en la probeta A para el primer y segundo bimaterial, respectivamente, suponiéndose la aproximación antes mencionada para /1. Así, )( 21 O / es menor para el bimaterial con D=0.8 y E=0.2 por lo que su )/( 2tzK asociado tiene un valor cada vez mayor cerca del vértice. El mismo análisis puede realizarse en la arista de la probeta B. K1/( V .L1-O) z/t r/L = 0.2 0.3 0.4 Barroso et. al. (2004) 0.000 0.67700 0.67877 0.67938 0.673688 0.500 0.67697 0.67876 0.67938 1.000 0.67669 0.67853 0.67918 Tabla 1: Factor de intensidad de tensiones generalizado K1 en la esquina de la probeta A en deformación plana. 874 MATERIALES COMPUESTOS 05 Es importante hacer notar que cuando )/( 1tzK tiende a infinito cerca del vértice, la severidad del orden de la singularidad asociada a dicho vértice es mayor que el correspondiente valor en la arista, por lo que en este caso el vértice de la esquina puede jugar un papel preponderante en la iniciación de grietas. En los problemas de esquinas bimaterial 3-D aquí estudiados, )/( 1tzK es casi constante a lo largo de la arista. Sin embargo, los )/( 2tzK tienen un valor cada vez más pronunciado para posiciones más cercanas al vértice haciendo que estos términos puedan incidir en mayor grado en el despegue o rotura de este tipo de materiales. Figura 3: Valores de )/( 1tza en la esquina bimaterial 3-D, Probeta A. 0.0 0.2 0.4 0.6 0.8 0.0 0.2 0.4 0.6 0.8 1.0 a 1 (z/t) z / t Serie1 Serie2 D E O   D E O   0 10 20 30 40 50 60 70 0.0 0.2 0.4 0.6 0.8 1.0 z / t a 2 (z/t) mat 1 mat D    E   O   D    E   O   Figura 4: Valores de )/( 2tza en la esquina bimaterial 3-D, Probeta A. MATERIALES COMPUESTOS 05 875 5. CONCLUSIONES Se ha presentado la formulación de una integral independiente del dominio para el cálculo de los Factores de Intensidad de Tensiones (FITG) a lo largo de la arista de una esquina 3-D. La formulación propuesta ha sido implementada en la etapa de posprocesamiento en un código del MEC. En comparación con otros métodos basados en integrales conservativas, solamente se requieren valores de desplazamientos y tensiones en puntos internos y del contorno del sólido. El presente procedimiento es válido para calcular los FITGs )(zK p asociados al término singular de tensiones así como a los términos de orden superior que verifican la condición 0<Op<Osing.+1, donde Osing representa el exponente singular más bajo en el desarrollo en serie de la solución elástica lineal en el entorno de la arista. El procedimiento desarrollado ha sido empleado para caracterizar esquinas bimateriales. Se han calculado los FITGs )(zK p asociados a los primeros exponentes Op en dos esquinas bimateriales, estableciéndose en 0.0 0.2 0.4 0.6 0.8 1.0 0.00.20.40.60.81.0 z/t a 1 (z/t ) Serie1 Serie2 D E O   D E O   Figura 5: Valores de )/( 1tza en la esquina bimaterial 3-D, Probeta B. Figura 6: Valores de )/( 2tza en la esquina bimaterial 3-D, Probeta B. -4 0 4 8 12 16 0.0 0.2 0.4 0.6 0.8 1.0 z/t a 2 (z/t) mat1 mat2 D E O   D E O   876 MATERIALES COMPUESTOS 05 general que )(zK p asociados a exponentes superiores de Op tienden a valores infinitos para posiciones cerca de un vértice en la esquina 3-D. Debe enfatizarse que los FITG fueron obtenidos sin usar elementos especiales o mallas muy finas. Así, en primer lugar, ningún requerimiento especial es necesario en el código MEC, la aplicación del MEC estándar es suficiente y, en segundo lugar, se ahorran recursos computacionales así como de tiempo de cálculo. AGRADECIMIENTOS El presente trabajo ha sido desarrollado durante la estancia de J.E.O. en la Universidad de Sevilla financiado por el Ministerio Español de Educación, Cultura y Deporte. A.B., V.M. y F.P. agradecen igualmente la financiación del proyecto No. MAT2003-03315 al Ministerio Español de Ciencia y Tecnología. REFERENCIAS 1. ANDERSSON, B., FALK U., BABUŠKA I. AND PETERSDORFF T.B. (1995). Reliable stress and fracture mechanics analysis of complex components using a h-p version of FEM. International Journal for Numerical Methods in Engineering 38, 2135-2163. 2. BANKS-SILLS L. (2002). A conservative integral for determining stress intensity factors of a bimaterial notch. International Journal of Fracture 115, 1-26. 3. BARROSO A., TORO P., MANTIý V. AND PARÍS F. (2004). Evaluation of generalized stress intensity factors in anisotropic elastic multimaterial corners. 11th European Conference on Composite Materials, Rhodas, Grecia. Pag. 1-10. 4. CHEN Y.H. AND LU T.J. (2004). On the path dependence of the J-integral in notch problems. International Journal of Solids and Structures 41, 607–618. 5. COSTABEL M., DAUGE M. AND YOSIBASH Z. (2004). A quasi-dual function method for extracting edge stress intensity functions. SIAM Journal Math. Anal. 35, 1177-1202. 6. HALLSTRÖM S. AND GRENESTEDT J.L. (1997). Mixed mode fracture of cracks and wedge shaped notches in expanded PVC foam. International Journal of Fracture 88, 343-358, 7. LABOSSIERE P.E.W. AND DUNN M.L. (2001). Fracture initiation at threedimensional bimaterial interface corners. Journal of the Mechanics and Physics of Solids 49, 609-634. 8. LEGUILLON D. AND SANCHEZ-PALENCIA E. (1999). On 3D cracks intersecting a free surface in laminated composites. International Journal of Fracture 99, 25-40. 9. MEDA G., MESSMER T. W. (1998). Sinclair G.B. and Solecki J.S. Pathindependent H-integral for three-dimensional fracture mechanics. International Journal of Fracture 94; 217-234. 10. ORTIZ J.E., MANTIý V. AND PARÍS F. (2005). A domain-independent integral for computation of stress intensity factors along three-dimensional crack fronts and edges by BEM. Enviado a publicación. 11. QIAN Z.Q. AND AKISANYA A. R. (1998). An experimental investigation of failure initiation in bonded joints. Acta Material 46, 14, 4895-4904. 12. SEWERYN A. (1994). Brittle fracture criterion for structures with sharp notches. Engineering Fracture Mechanics 47, 673-681. 13. SZABÓ, B.A. AND BABUŠKA, I. (1991). Finite Element Analysis, John Wiley & Sons, New York.