Application of the Stiffness-Force method in the analysis of the fiber-matrix adherence in composite laminates

1 Introducción

El método de la rigidez (SM) en elementos finitos puede interpretarse como una forma aproximada de resolver las condiciones de equilibrio, derivadas éstas del principio de los trabajos virtuales. En el método de la rigidez la ecuación de gobierno es

${\textstyle \left[K\right]\left\{x\right\}=}$$\left\{P\right\}$ (1)

El método de las fuerzas integrado dual (IFMD) tiene la misma ecuación de gobierno de la ecuación (1), sin embargo, en este método además de las ecuaciones de equilibrio, se tienen en cuenta las ecuaciones de compatibilidad, derivadas éstas del principio del trabajo virtual complementario [1].

2 El método Rigidez-Fuerza (SFM)

Las ecuaciones del IFMD para un elemento continuo discretizado en elementos finitos se obtienen a partir del principio del trabajo virtual y del trabajo virtual complementario.

${\textstyle \int _{V}^{}{\left\{\delta \epsilon \right\}}^{T}\left\{\sigma \right\}dv=}$${\left\{\delta {a}^{i}\right\}}^{T}\left\{{P}^{i}\right\}$ (2)

${\textstyle \int _{V}^{}{\left\{\delta \sigma \right\}}^{T}\left\{\epsilon \right\}dv=}$${\left\{\delta {P}^{i}\right\}}^{T}\left\{{P}^{i}\right\}$ (3)

{aⁱ} y {Pⁱ} son los vectores de desplazamientos y de fuerzas nodales respectivamente. Con los desplazamientos nodales, se conocen los desplazamientos en cualquier punto del elemento [2]:

${\textstyle \left\{u\right\}=\left[N\right]\left\{{a}^{i}\right\}}$

(4)

donde [N] es la matriz de funciones de interpolación de desplazamientos A partir de los desplazamientos nodales se conocen también las deformaciones en el elemento,

${\textstyle \left\{\epsilon \right\}=}$$\left[L\right]\left[N\right]\left\{{a}^{i}\right\}=\left[B\right]\left\{{a}^{i}\right\}$ (5)

L] es la matriz de operadores diferenciales y [B] es la matriz de forma. En este método las tensiones y los desplazamientos se interpolan de forma independiente [3]:

${\textstyle \left\{\sigma \right\}=}$$\left[Y\right]\left\{{F}^{j}\right\}$ (6)

En el que [Y] es la matriz de funciones de interpolación de fuerzas que satisfacen las ecuaciones de equilibrio. El vector de tensiones se expresa en función del vector deformaciones mediante la matriz de coeficientes de flexibilidad del material [S]:

${\textstyle \left\{\sigma \right\}=}$$\left[S\right]\left\{\epsilon \right\}$ (7)

Sustituyendo las expresiones (5) y (6) en el principio del trabajo virtual (2) se obtiene la ecuación de equilibrio del elemento:

${\textstyle \left(\int _{V}^{}{\left[B\right]}^{T}\left[Y\right]dv\right)\left\{{F}^{j}\right\}=}$$\left\{{P}^{i}\right\}$ (8)

Donde la integral de volumen es la matriz de equilibrio del elemento ${\textstyle \left[E\right]=}$$\int _{v}^{}{\left[B\right]}^{T}\left[Y\right]dV$ . Así la ecuación de equilibrio (8)se escribe como:

${\textstyle \left[E\right]\left\{{F}^{j}\right\}=}$$\left\{{P}^{i}\right\}$ (9)

Ahora, sustituyendo las expresiones (6), (7) y (9) en el principio del trabajo virtual complementario (3), éste se escribe como:

${\textstyle \left(\int _{V}^{}{\left[Y\right]}^{T}\left[S\right]\left[Y\right]dv\right)\left\{{F}^{j}\right\}=}$$\left[{E}^{e}\right]\left\{{a}^{i}\right\}$

(10)

Donde la integral de volumen es la matriz de flexibilidad ${\textstyle \left[{G}^{e}\right]=}$$\int _{v}^{}{\left[Y\right]}^{T}\left[S\right]\left[Y\right]dV$ del elemento [4].

${\textstyle \left[{G}^{e}\right]\left\{{F}^{j}\right\}=}$$\left[{E}^{e}\right]\left\{{a}^{i}\right\}$ (11)

De la ecuación (11) se tiene el vector de fuerzas independientes del elemento {F^j}:

${\textstyle \left\{{F}^{j}\right\}={\left[{G}^{e}\right]}^{-1}\left[{E}^{e}\right]\left\{{a}^{i}\right\}}$

(12)

Multiplicando en ambos lados de dicha expresión por la matriz [E^e] de equilibrio, el primer término de la ecuación corresponde al vector de fuerzas externas aplicadas en los nodos, {Pⁱ},

${\textstyle \left\{{P}^{i}\right\}=\left[E\right]{\left[G\right]}^{-1}{\left[E\right]}^{T}\left\{{a}^{i}\right\}}$

(13)

Así tenemos la ecuación de gobierno IFMD [5] y la matriz de rigidez del elemento:

${\textstyle \left\{{P}^{i}\right\}=\left[K\right]\left\{{a}^{i}\right\}}$

(14)

${\textstyle \left[K\right]=\left[E\right]{\left[G\right]}^{-1}{\left[E\right]}^{T}}$

(15)

3 Elemento cuadrilátero axisimétrico

3.1 Formulación SFM

A continuación se va a aplicar el método SFM para formular el elemento cuadrilátero axisimétrico de 4 nodos de la figura 1.

Figura 1. Elemento cuadrilátero axisimétrico.

Como se ha visto en el apartado anterior, para obtener la matriz de rigidez del elemento es necesario calcular previamente las matrices de equilibrio y de flexibilidad.

En la formulación de este elemento se van a utilizar las mismas funciones de interpolación lineales que en el cuadrilátero general [N_i] [6]:

${\textstyle \left[{N}_{i\,}\right]={\frac {1}{4}}\left(1+\xi {\xi }_{i}\right)\left(1+\right.}$$\left.\eta {\eta }_{i}\right)$

(16)

En el caso de un elemento axisimétrico la matriz de operadores diferenciales [L] es la siguiente:

$\left[L\right]=\left[{\begin{array}{cccc}{\frac {\partial }{\partial r}}&{\frac {1}{r}}&0&{\frac {\partial }{\partial z}}\\0&0&{\frac {\partial }{\partial z}}&{\frac {\partial }{\partial r}}\end{array}}\right]$

(17)

La matriz de deformación del elemento [B] se obtiene derivando la matriz de funciones de interpolación es la siguiente:

(18)

3.3 Matriz de Equilibrio [E^e]

Para el caso de un cuadrilátero general axisimétrico la matriz de equilibrio se calcula mediante la siguiente expresión:

${\textstyle \left[E\right]=2\pi \int _{V}^{}{\left[B\right]}^{T}\,\left[Y\right]{r}_{C}\left|J\right|dV}$

(19)

Donde r_c es la coordenada de cada punto de integración y ${\textstyle \left|J\right|}$ el determinante del Jacobiano. En el caso de un cuadrilátero general la matriz de forma [B] varía con las coordenadas de cada punto, y como matriz de funciones de interpolación de fuerzas [Y] se ha tomado la siguiente expresión [7,8]:

${\textstyle \left[Y\right]=\left[{\begin{array}{ccccccc}1&0&0&J_{11}^{2}\eta &J_{21}^{2}\xi &0&0\\0&1&0&J_{12}^{2}\eta &J_{22}^{2}\xi &0&0\\0&0&0&0&0&1&J_{12}^{}\xi +J_{22}^{}\eta \\0&0&1&J_{11}J_{12}\eta &J_{21}J_{22}\xi &0&0\end{array}}\right]}$

(20)

3.4 Matriz de Flexibilidad [Ge]

La matriz de flexibilidad del elemento se calcula, en un cuadrilátero axisimétrico mediante la siguiente expresión:

${\textstyle \left[G\right]=2\pi \int _{V}^{}{\left[Y\right]}^{T}\,\left[S\right]\left[Y\right]{r}_{C}\left|J\right|dV}$

(21)

Donde [S] es la matriz de flexibilidad del material:

$\left[S\right]=\left[{\begin{array}{cccc}{\frac {1}{E}}&-{\frac {\nu }{E}}&-{\frac {\nu }{E}}&0\\-{\frac {\nu }{E}}&{\frac {1}{E}}&-{\frac {\nu }{E}}&0\\-{\frac {\nu }{E}}&-{\frac {\nu }{E}}&{\frac {1}{E}}&0\\0&0&0&{\frac {2\left(1+\nu \right)}{E}}\end{array}}\right]$

(22)

3.5 Matriz de Rigidez del elemento [K^e]

Conocidas las matrices de equilibrio y de flexibilidad del elemento, lse obtiene la matriz de rigidez del elemento:

${\textstyle \left[{K}^{e}\right]=\left[{E}^{e}\right]{\left[{G}^{e}\right]}^{-1}{\left[{E}^{e}\right]}^{T}}$

(23)

Una vez calculadas las matrices de rigidez de los elementos, la matriz de rigidez del sistema se obtiene ensamblando las matrices de rigidez de los elementos.

4 Análisis mediante SFM del ensayo pull-out de fibra única

El ensayo de tracción de fibra única es una técnica experimental muy utilizada para la caracterización de las propiedades de la interfase entre la fibra y la matriz en materiales compuestos [9]. Existen en la bibliografía aproximaciones analíticas [10,11], basadas en funciones de tensión que cumplen las condiciones de equilibrio, y aproximaciones numéricas [12,13] que, con diferentes tipos de elementos y altos grados de discretización, proporcionan la distribución de tensiones en la fibra y en la matriz.

En este trabajo, se ha analizado por elementos finitos mediante el método SFM el modelo axisimétrico representado en la figura 2 del ensayo pull-out de fibra única.

Figura 2. Ensayo pull-out de fibra única.

Para el análisis mediante el método SFM se ha realizado una subrutina de usuario UEL de un elemento cuadrilátero de cuatro nodos axisimétrico, y se ha implementado en el programa comercial Abaqus.

El estado de cargas y las condiciones de contorno del modelo axisimétrico analizado se indican en la figura 3, que en un primer paso se correponden con el estado intacto o previo al despeque.

Figura 3. Modelo axiximétrico y condiciones de contorno.

En las figuras 4 y 5 se muestra la distribución de tensiones axiales σ^f_z y de cortadura σ^f_rz de la fibra en la unión fibra-matriz.

Figura 4. Distribución de tensiones σ_zen la unión fibra matriz.

Figura 5. Distribución de tensiones τ_rz en la unión fibra matriz.

5 Conclusiones

Se ha formulado mediante el método SFM un cuadrilátero general axisimétrico y se ha implementado en una subrutina de usuario (UEL) de Abaqus.

6 Agradecimientos

Los autores agradecen a la Universidad del País Vasco (UPV/EHU) la financiación del Grupo de Investigación Mecánica de Materiales GIU 16/51 en la convocatoria de 2016.

Referencias

[1] S.N. Patnaik et al. Integrated Force Method Solution to Indeterminate Structural Mechanics Problems. National Aeronautics and Space Administration (2004).

[2] E. Oñate. Cálculo de Estructuras por el Método de Elementos Finitos” Centro Internacional de Métodos Numéricos en Ingeniería (1995).

[3] Kaljevic, S.N. Patnaik, D.A. Hopkins. Three dimensional structural analysis by the integrated force method. Computers and Structures 59, 4,pag 691-706 (1996).

[4] J. Robinson, Integrated theory of finite element methods. Wiley and sons (1973).

[5] S.N.Patnaik, M.S. Nagraj, Analysis of continuum by the integrated force method. Computers and Structures, 26: 899-905. (1987).

[6] E. Oñate. Cálculo de Estructuras por el Método de Elementos Finitos. Centro Internacional de Métodos Numéricos en Ingeniería, Barcelona (1995).

[7] THH Pian, K Sumihara. Rational approach for assumed stress finite elements, Int J Numer Methods Eng. 20: 1685-1695 (1984).

[8] C. Zhang, D. Wang, J. Zhang, W. Feng, Q. Huang. On the equivalence of various hybrid finite elements and a new orthogonalization method for explicit element stiffness formulation, Finite Elements Anal.Des. 43: 321-332 (2007).

[9] E. Graciani, V. Mantič, F. París, J. Varnab. Numerical analysis of debond propagation in the single fibre fragmentation test. Composites Science and Technology. 69,15–16: 2514–2520 (2009).

[10] M.Y.Quek, C.Y. Yue. Axisymmetric stress distribution in the single filament pull-out test. Materials Science and Engineering, A198: 105-116 (1994).

[11] Y. Wang, L. Zhou, Z. Wang, H. Huang, L. Ye. Stress distributions in single shape memory alloy fiber composites. Material and Design:32: 3783-3789 (2011).

[12] W. Sun, F. Lin. Computer modeling and FEA simulation for composite single fiber pull-out. Journal of thermoplastic composite materials. 14: 327-343 (2001).

[13] V. Bheemreddy, K. Chandrashekhara, L.R. Dharani, G.E. Hilmas. Modeling of fiber pull-out in continuous fiber reinforced ceramic composites using finite element method and artificial neural networks. Computacional Materials Science. 79: 663-673 (2013)