Reconstrucción de tensiones para el método de elementos finitos con mallas poligonales

 

Stress recovery for the polygonal finite element method

 

 

Octavio Andrés González-Estrada1, Sundararajan Natarajan2, Carlos Graciano3

 

 

1Escuela de Ingeniería Mecánica, Universidad Industrial de Santander, Ciudad Universitaria, Bucaramanga, Colombia.  Orcid: 0000-0002-2778-3389, Email: agonzale@uis.edu.co

2Department of Mechanical Engineering, Indian Institute of Technology – Madras, Chennai 600036, India.

Orcid: 0000-0002-0409-0096, Email: snatarajan@iitm.ac.in

3Departamento de Ingeniería Civil, Universidad Nacional de Colombia, Medellín, 050034, Colombia,

Orcid: 0000-0003-0659-7963, email: cagracianog@unal.edu.co

 

 

 

RESUMEN

 

El método de los elementos finitos es una de las herramientas numéricas más utilizadas para el diseño en ingeniería. En los últimos años se han desarrollado nuevas aproximaciones numéricas para extender el uso del método de elementos finitos a mallas poligonales. Dichas aproximaciones mejoran la precisión de la solución y aumentan la flexibilidad en el mallado. Sin embargo, como toda aproximación, es necesario cuantificar el valor del error inducido para poder validar los resultados obtenidos. En este trabajo se presenta el uso de una técnica de estimación del error de discretización para mallas de elementos finitos poligonales. La técnica está basada en la reconstrucción de la solución en tensiones mediante un procedimiento de mínimos cuadrados ponderados que considera la influencia de las ecuaciones de equilibrio. Se ha utilizado un problema con solución exacta para evaluar la efectividad del estimador, obteniendo buenos resultados a nivel local y global.

 

PALABRAS CLAVE: Método de elementos finitos, Elementos poligonales, Estimación del error, Adaptatividad, Reconstrucción de tensiones, Ajuste por mínimos cuadrados.

 

ABSTRACT

 

The finite element method is one of the most used numerical tools for engineering design. In recent years, novel numerical approximations have been proposed to extend the finite element method to meshes using arbitrary polygons. Such approximations are aimed to improve accuracy and increase the flexibility during the meshing process. However, as any approximation, they exhibit an error that requires to be quantified in order to validate the numerical results. In this paper, we present a technique to estimate the discretization error in energy norm for arbitrary polygonal finite element meshes. This recovery-based technique uses a moving least squares approach that considers constraints to represent the equilibrium equations.  We use two benchmark problems to evaluate the effectivity of the error estimator, with good results both locally and globally.

 

KEYWORDS: Finite element method, Polygonal elements, Error estimation, Adaptivity, Stress recovery, Moving Least Squares.

 

INTRODUCCIÓN

 

 

Este artículo puede compartirse bajo la licencia CC BY-ND 4.0 y se referencia usando el siguiente formato: O.A. González-Estrada, S. Natarajan,

C. Graciano, “Reconstrucción de tensiones para el método de elementos finitos con mallas poligonales”, UIS Ingenierías, vol. 16, no. 1, pp. 23-34, Enero-Junio 2017.

En el método de los elementos finitos (MEF) se han considerado tradicionalmente topologías con elementos triangulares o cuadriláteros en 2D, o sus extensiones a 3D. En los últimos años ha aumentado el interés por elementos finitos con formas pentagonales, hexagonales en [5] se presentó un proceso sistemático para la construcción de funciones de forma para elementos poligonales. Mucho más tarde, con los trabajos presentados en [6]–[8] se ha abierto una nueva era de elementos finitos llamada elementos finitos poligonales, en los cuales la topología no está limitada a usar cuatro


lados en 2D o seis caras en 3D. El uso de elementos finitos poligonales de más de cuatro lados ofrece una mayor flexibilidad en el mallado y una mejora en la precisión de los resultados [9]. Rashid y Gullet [10] propusieron el método de elementos finitos con topología de elemento variable (VETFEM), en el cual las funciones de forma para elementos convexos y no convexos son calculadas en el espacio físico (𝐱 Ω) usando un procedimiento de optimización. Ghosh et al. [11], [12] desarrollaron los modelos de elementos finitos con celdas Voronoï (VCFEM) para microestructuras heterogéneas de materiales compuestos o porosos. Las celdas en VCFEM son arbitrarias y no se ajustan a la definición estricta de Voronoï. Otros desarrollos incluyen elementos poligonales basados en métodos de nodos virtuales [13] y el método de elementos virtuales [14].

 

En general, los métodos numéricos proveen una solución aproximada al problema físico en consideración [3], [15], [16]. La aproximación numérica tiene entonces un error que debe ser cuantificado para validar la calidad del modelo numérico. Las técnicas de estimación y control del error han sido ampliamente investigadas para FEM y otros métodos numéricos [17]–[19]. El estudio del error en las técnicas más recientes de aproximación numérica como, por ejemplo, los elementos finitos poligonales, resulta necesario para poder garantizar los resultados numéricos obtenidos para aplicaciones prácticas. En este trabajo se presenta el uso de una técnica de estimación del error de discretización, basada en reconstrucción de la solución, para aproximaciones de elementos finitos poligonales. La técnica de reconstrucción utiliza un procedimiento mejorado de mínimos cuadrados ponderados sobre la solución de tensiones de elementos finitos [20].

 

El artículo está organizado de la siguiente manera: En la sección de metodología se define en primer lugar el problema elástico lineal en consideración. Posteriormente, se presenta la técnica de elementos finitos poligonales. Se sigue con una introducción al concepto de estimación del error de discretización y se presenta la técnica de reconstrucción mediante mínimos cuadrados ponderados. En la sección de resultados numéricos se presentan los resultados para un do problema de referencia y, finalmente, en la última sección se dan las conclusiones.

 

 

2. METODOLOGÍA

 

2.1. Planteamiento del problema elástico lineal

 

En esta sección se presenta el problema 2D de elasticidad lineal para el cual se aplica la metodología propuesta. Considérese el vector de tensiones 𝛔 = {𝜎𝑥𝑥,𝜎𝑦𝑦, 𝜎𝑥𝑦}𝑻, los desplazamientos u y las deformaciones ε, definidos sobre el dominio 𝛺 𝟐. Sean b las cargas volumétricas, t las tracciones de Neumann y ̅𝐮 las condiciones de Dirichlet. El problema de elasticidad se define entonces como encontrar 𝐮 tal que satisfaga:

 

𝐋T𝛔 + 𝐛 = 0 en Ω (1)

𝐆𝛔 = 𝐭 en ΓN        (2)

𝐮 = ̅𝐮 en ΓD         (3)

𝛆(𝐮) = 𝐋𝐮    (4)

𝛔 = 𝐃𝛆(𝐮), (5)

 

donde L es el operador diferencial de elasticidad lineal, D es la matriz con la elasticidad del material y G es el operador que proyecta las tensiones para el equilibrio en el contorno.

 

El problema anterior puede ser expresado en su forma variacional como: Encontrar  𝐮 𝑽 tal que v 𝑽

 

∫      𝛆(𝐯)T𝐃𝛆(𝐮)dΩ = ∫𝐯T𝐛

Ω          Ω          (6) 

+ ∫ 𝐯T𝐭  dΓ, ΓN

 

donde 𝑽 = {𝐯|𝐯 [𝐇1(𝛀)]𝟐, 𝐯𝚪 = 𝟎}.

 

2.2. Resolución    mediante    elementos finitos

poligonales

 

Formulando la discretización del método de elementos finitos para el problema anterior, se busca la solución

𝐮 𝑽  tal que ∀𝐯 𝑽. Usando Galerkin se tiene:

 

∫      𝛆(𝐯)T𝐃𝛆(𝐮h)dΩ = ∫ 𝐯T𝐛

Ω          Ω          (7) 

+ ∫ 𝐯T𝐭  dΓ, ΓN

 

y se resuelve utilizando la teoría clásica de elementos finitos. El uso de elementos poligonales requiere el cálculo de funciones de forma que satisfagan las siguientes propiedades: no negativas, interpolación, partición de la unidad y completitud.

 

 

 

2.2.1. Aproximación sobre elementos poligonales arbitrarios.

Wachspress [5] introdujo el concepto de funciones de forma para geometrías arbitrarias, incluyendo funciones de interpolación para polígonos cóncavos o convexos. Las funciones de interpolación no convencionales Wachspress no corresponden a las funciones polinómicas usadas en el FEM. Funciones de interpolación basadas en subdivisiones Voronoï fueron introducidas por Sibson [21], y han sido usadas ampliamente en la mecánica. Las funciones de Laplace también son un esquema de interpolación basada en los elementos vecinos. Para avances recientes en elementos poligonales se recomienda la revisión presentada en [7]. En este trabajo se ha considerado el uso de funciones de forma de Laplace con un mapeo isoparamétrico. Las funciones de forma de tipo Laplace, también llamadas de funciones de vecinos naturales, proveen una función de ponderación natural para nodos distribuidos de forma arbitraria en el espacio. Para un punto P con n vecinos naturales, la función de forma de Laplace para el nodo PI  se escribe como:

 

,

𝐽=1 𝛼𝐽       

(8) 

𝑠 (𝐱)

 

 

donde 𝛼𝐼(𝐱) es la función de ponderación de Laplace, 𝑠𝐼(𝐱) es la longitud del borde de Voronoï asociado con P y el nodo PI y 𝑰(𝐱) es la distancia Euclidiana entre los puntos P y PI, ver Figura 1.

 

2.2.2. Integración numérica.

Las funciones de forma sobre elementos poligonales son polinomios racionales, por lo tanto, la definición de una cuadratura de integración eficiente constituye una dificultad. Una posible solución es la subdivisión del elemento físico en triángulos sobre los cuales se puede aplicar posteriormente una de las cuadraturas usadas para dicha topología. La subdivisión se utiliza únicamente para facilitar la integración y no introduce grados de libertad adicionales. Este procedimiento implica un mapeo isoparamétrico de dos niveles y se sustenta en la positividad del Jacobiano de la transformación. 

 

2.3. Estimación del error en norma energética

 

En la ausencia de otros tipos de errores, el error de discretización de elementos finitos se define como 𝐞 = 𝐮𝐮. Para cuantificar la calidad del campo 𝐮 se suele evaluar la norma energética del error e, la cual se define como:

 

.        (9) 

 

Figura 1. Diagrama Voronoï de un punto P.

 

Introduciendo el error en el campo de tensiones como 𝐞𝛔 = 𝛔𝛔𝒉 y utilizando la ecuación constitutiva se puede reescribir:

 

.   (10) 

 

Siguiendo el trabajo de Zienkiewicz y Zhu [22], es posible obtener una estimación 𝜺 del error exacto ǁ𝐞ǁ  mediante la aproximación

 

,    (11) 

 

siendo 𝐞𝛔 = 𝛔 𝛔 el error aproximado en tensiones y 𝛔 el campo reconstruido de tensiones. Si 𝛔 converge a la solución exacta a una velocidad más alta que la solución de elementos finitos, el estimador es asintóticamente exacto. Es decir, el error aproximado tiende al error exacto con el refinamiento de la malla. La calidad del estimador depende entonces de la técnica usada para la reconstrucción de tensiones.

 

 

2.4. Reconstrucción de tensiones mediante mínimos cuadrados

 

La técnica moving least squares (MLS) está basada en un ajuste de mínimos cuadrados ponderados hacia el punto de prueba donde se reconstruye el valor de la tensión. La técnica considera una expansión polinómica para cada una de las componentes del campo de tensiones

reconstruidas, en la forma

 

𝜎𝑖(𝐱) = 𝐩(𝐱)𝐚i(𝐱),      𝑖 = 𝑥𝑥, 𝑦𝑦, 𝑥𝑦, (12)

 

donde p representa una base polinómica y a son los

coeficientes a evaluar

 

𝐩(𝐱) = {1 𝑥 𝑦 𝑥2 𝑥𝑦 𝑦2 … }      (13)

    𝐚𝑖(𝐱)

𝑇              (14)

= {𝑎0𝑖(𝑥) 𝑎1𝑖(𝑥) 𝑎2𝑖(𝑥) 𝑎3𝑖(𝑥)…  }

 

Para 2D, la expresión para evaluar los campos de

tensiones reconstruidas se escribe

 

𝜎𝑥𝑥∗ (𝐱)

𝛔(𝐱) = {𝜎𝑦𝑦 (𝐱)} = 𝐏(𝐱)𝐀(𝐱)

𝜎𝑥𝑦∗ (𝐱)

(15)

𝐩(𝐱)   0      0      𝐚𝑥𝑥(𝐱)

 

 

El formato de la ecuación (15), que considera las tres componentes del campo de tensiones en una sola ecuación, resulta útil para imponer las ecuaciones de restricción utilizadas para imponer las condiciones de equilibrio en la reconstrucción.

 

Supongamos que 𝛘 es un punto en 𝛀𝐱, siendo 𝛀𝐱 el soporte de un punto 𝐱 definido por el radio 𝑹𝜴𝒙. La aproximación MLS para cada componente de tensión en

𝛘 está dada por

 

σ𝑖 (𝐱, 𝛘) = 𝐩(𝛘)𝐚𝒊(𝐱),     ∀𝛘 Ωx.     (16)

 

Para obtener los coeficientes A utilizamos la aproximación continua de MLS descrita en [23], donde se minimiza el siguiente funcional:

 

 

 

Evaluando 𝜕𝐽/𝜕𝐀 = 0 se obtiene el sistema lineal

𝐌(𝐱)𝐀(𝐱) = 𝐆(𝐱) usado para calcular A, donde

 

       (18)

 

 

Asumiendo que existen 𝑛 puntos de muestreo de coordenadas 𝛘𝒍(𝑙 = 1 … 𝑛) dentro del soporte de 𝐱, con peso 𝑯𝒍 y siendo |𝐉𝑙)| el determinante del jacobiano, las expresiones de las ecuaciones (17,18) pueden ser evaluadas numéricamente como

 

𝑛

𝐽(𝐱) = ∑ 𝑊(𝐱𝛘𝑙)[𝛔(𝐱, 𝛘𝑙)

𝑙=1

𝛔(𝛘𝑙)]2|𝐉(𝛘𝑙)|𝐻𝑙

𝑛

𝑇(𝛘𝑙)𝐏(𝛘𝑙)|𝐉(𝛘𝑙)|𝐻𝑙 (19)

𝐌(x) = ∑ 𝑊(𝐱𝛘𝑙)𝐏

𝑙=1

𝑛

𝐆(𝐱) = ∑ 𝑊(𝐱𝛘𝑙)𝐏𝑇(𝛘𝑙)𝛔(𝛘𝑙)|𝐉(𝛘𝑙)|𝐻𝑙

𝑙=1

 

Los puntos de integración para la evaluación numérica de las integrales mostradas corresponden a los puntos de integración en el dominio 𝛀𝐱, para los cuales el campo de tensiones ya está calculado. En (18), 𝑊 es la función de ponderación de MLS, la cual ha sido tomada como una

función de cuarto orden tal que 𝑠

 

𝑊(𝐱𝛘)

= {1 − 6𝑠2 + 8𝑠3 − 3𝑠4  si  |𝑠| ≤ 1   (20)

0      si  |𝑠| > 1

 

Donde 𝑠 representa la función de distancia normalizada dada por

 

     (21)

 

En el enfoque continuo de MLS cada uno de los puntos de muestreo 𝛘𝑙 está ponderado usando su área de influencia dada por |𝐉(𝛘𝑙)|𝐻𝑙. La literatura indica que la aproximación continua es más precisa que la aproximación discreta [23], especialmente cuando la distribución de puntos no es uniforme dentro del soporte. La continuidad del campo 𝛔 está garantizada por el procedimiento de MLS continuo descrito, ya que la función de ponderación asegura que los puntos de muestreo entren o salgan del soporte de una manera gradual y suave cuando 𝐱 se mueve.

 

2.5. Satisfacción de la ecuación de equilibrio en la frontera

 

La ecuación de equilibrio en el contorno debe satisfacerse en cada punto a lo largo del contorno de Neumann. En [24]–[26], considerando una técnica tipo SPR (superconvergent patch recovery), se imponen condiciones de equilibrio en el contorno utilizando multiplicadores de Lagrange. Sin embargo, este enfoque produce discontinuidades en una formulación MLS cuando pasamos de un soporte en el interior del dominio a un soporte que se interseca con el contorno.

 

Con el fin de evitar discontinuidades en el campo reconstruido, utilizamos un enfoque basado en el punto más cercano para introducir la imposición de la ecuación de equilibrio de manera suave en el funcional. Para un punto 𝐱 Ω, cuyo soporte Ωx se interseca con Γ, las ecuaciones de restricción para el equilibrio se consideran solamente en los puntos más cercanos  𝛘𝑗 Γ, como se muestra en la Figura 2.

 

 

 

Figura 2. Soporte MLS con las condiciones de contorno aplicadas en los puntos más cercanos.

 

Nótese que para un soporte dado se puede tener más de un punto más cercano, como ocurre para un punto 𝐱 que se acerca a una esquina, donde se toma un punto por cada lado de la esquina con el fin de evitar problemas de continuidad.

 

Expresamos el vector de tensión 𝛔(𝐱, 𝛘) en un sistema coordenado 𝑥̃𝑦̃ alineado con el contorno en 𝛘𝑗, tal que 𝑥̃ es el vector normal al contorno, rotado un ángulo 𝛼 con respecto a 𝑥:

 

𝛔̃(𝐱, 𝛘) = 𝐑(𝛼)𝛔(𝐱, 𝛘), (22)

 

donde 𝐑 es la matriz de rotación de tensiones

 

𝐫𝑥̃𝑥̃

𝐑 = [𝐫𝑦̃𝑦̃ ]

𝐫𝑥̃𝑦̃

cos2α sin2α  sin(2α)       (23)

= [    sin2α  cos2α −sin(2α)]

−sin(2α)/2  sin(2α)/2    cos(2α)

 

El funcional de MLS expresado en su versión continua e incorporando las restricciones en los contornos se escribe:

 

𝑛

𝐽(𝐱) = ∑ 𝑊(𝐱𝛘𝑙)[𝛔(𝐱𝑙,𝛘𝑙)

𝑙=1

𝛔(𝛘𝑙)]2|𝐉(𝛘𝑙)|𝐻𝑙

𝑛𝑏𝑐

+ ∑ 𝑊̃ (𝐱

𝑗=1

𝛘𝑗)[𝜎𝑖̃ (𝐱, 𝛘𝑗)

𝑛              𝜎𝑖𝑒𝑥̃   (𝛘𝑗)]2 (24)

𝐽(𝐱) = ∑ 𝑊(𝐱𝛘𝑙)[𝐏(𝐱𝑙,)𝐀(x)

𝑙=1

𝝈(𝐱𝑙)]2|𝐉(𝐱𝑙)|𝐻𝑙

𝑛𝑏𝑐

+ ∑ 𝑊̃ (𝐱

𝑗=1

𝛘𝑗)[𝑟𝑖̃(𝛼)𝐏(𝛘𝑗)𝑨(𝐱)

𝜎𝑖𝑒𝑥̃  (𝛘𝑗)]2𝑖̃ = 𝑥̃𝑥̃, 𝑥̃𝑦̃

 

Donde 𝑛𝑏𝑐 es el número de puntos 𝛘𝑗 sobre el contorno donde se consideran las condiciones de Neumann conocidas . Evaluando 𝜕𝐽/𝜕𝑨 = 0 se obtiene el sistema lineal 𝐌(𝐱)𝐀(𝐱) = 𝐆(𝐱) que se resuelve para 𝐀,

donde

 

 

 

En las expresiones anteriores 𝑊̃ es una función de ponderación definida como

 

(26)

 

 

 

 

 

 

 

 

La función  ̃𝑊 tiene dos características principales:

 

1.                    Incluye a 𝑊 de manera que la restricción en la frontera del dominio se introduzca de manera suave en el funcional 𝐽(𝐱). Así se asegura que el campo reconstruido mantenga la continuidad.

2.                    Incluye el término 𝑠−1 de forma que el peso de la restricción en 𝐽(𝐱) aumente a medida que nos acercamos al contorno (cuando  𝐱 𝛘𝑗, 𝑠 0), por lo tanto, 𝛔 tiende a cumplir la condición de equilibrio en el contorno. En el caso de que  𝑠 = 0, usamos multiplicadores de Lagrange.

 

2.6. Satisfacción de la ecuación de equilibrio interno

 

Con el fin de mejorar la calidad del campo reconstruido se trata de imponer adicionalmente la satisfacción del equilibrio interno dado por

 

𝛔 + 𝐛 = 0      (27)

 

Las derivadas espaciales de 𝛔, considerando (15), se pueden expresar como 

 

𝛔 = ( 𝐏)𝐀 + 𝐏( 𝐀)     (28)

 

El primer término en (28) puede evaluarse derivando la base polinómica, quedando definido un pseudoequilibrio [27], [28]. Para evaluar el segundo término derivamos el sistema lineal 𝐌𝐀 = 𝐆:

 

( 𝐌)𝐀 + 𝐌( 𝐀) = 𝐆    (29)

 

Evaluando 𝐀 de (29) y reemplazando en (28) se tiene 

(30)

,

(31)

 

donde las derivadas parciales de 𝐌 y 𝐆 con respecto a  𝑥

son

 

    (32)

   (33)

 

Derivando (20) se tiene

 

(34)

        (35)

 

Las ecuaciones (34,35) son funciones de 𝐀, por lo tanto, se pueden escribir los términos de la ecuación de equilibrio interno como función del vector de incógnitas

𝐀:

 

 

 

 

donde 𝑖𝑗 (𝑖 = 𝑥𝑥, 𝑦𝑦, 𝑥𝑦 y 𝑗 = 𝑥, 𝑦) representa la fila j correspondiente a la i-ésima componente de tensiones. Estas expresiones definen las restricciones para el equilibrio interno en 𝐱 entre los coeficientes 𝐀. Se utilizan multiplicadores de Lagrange para imponer dichas restricciones, resultando en el sistema lineal

 

,       (38)

 

donde 𝐂 y 𝐃 son los términos utilizados para imponer la restricción y 𝛌 es el vector de multiplicadores de Lagrange. Cabe señalar que, en el procedimiento propuesto en este trabajo no se ha considerado el término 𝐂𝑇𝛌 al evaluar las derivadas parciales de 𝐀  La técnica MLS utilizada produce así un campo de tensiones reconstruidas que es continuo y cuasi-estáticamente admisible.

 

3. RESULTADOS NUMÉRICOS

 

En esta sección se presenta un problema elástico lineal en 2D con solución exacta, utilizado para evaluar el desempeño del estimador de error propuesto para mallas con elementos poligonales. Para analizar los resultados se utilizan el error global en norma energética y los índices de efectividad global y local que se definen a continuación.

 

La efectividad del estimador de error se calcula con el índice global

 

.       (39)

 

Szabó y Babuška [15] describen como un buen estimador del error, capaz de ser utilizado para guiar procesos adaptativos, a aquel con efectividades dentro del rango de .

 

Para la efectividad local definimos un índice a nivel de

elemento D como

 

        (40)

,

 

donde es la efectividad en el elemento calculada en

 

3.1. Placa cuadrada

 

Se modela una placa cuadrada de , como se muestra en la Figura 3, con propiedades del material dadas por  para el módulo de Young y  para el coeficiente de Poisson. El problema se define de forma que la solución exacta del problema elástico en

desplazamientos está dada por

 

      (41)

 

 

Los valores exactos de las componentes de tensión se aplican sobre el contorno de Neumann, indicado por la línea discontinua en la Figura 3. Las tensiones a imponer como condición de contorno del problema de elementos finitos, se obtienen a partir de la solución en desplazamientos bajo condiciones de deformación plana, y se escriben:

 

       (42)

 

 

Figura 3. Placa cuadrada de .

 

Se imponen las siguientes fuerzas volumétricas en el problema para garantizar el equilibrio interno:

 

 

   (43)

.

 

Los resultados se obtienen para secuencias de cuatro mallas de elementos poligonales, polyFEM, y se comparan con mallas de elementos finitos convencionales con elementos cuadriláteros lineales, FEM-Q4, y triángulos lineales, FEM-T3. En la Figura 4 se muestran las primeras mallas de la secuencia de refinamiento uniforme utilizada para cada topología. Para un mismo número de grados de libertad (GDL), las mallas de polígonos de -lados tienen más elementos a menores valores de , siendo máximo para triángulos, con . Las mallas se construyen de forma que se mantiene un número similar de gdl para los tres tipos de elemento a medida que se refina. Para evaluar el error se utiliza la técnica de reconstrucción MLS propuesta.

 

En la Tabla 1 se muestra la evolución del error global estimado en norma energética, , la efectividad global, , y la velocidad de convergencia, s, a medida que aumenta el número de grados de libertad. La velocidad de convergencia total del error en norma energética es cercana a 0.5 para los tres tipos de elementos utilizados, indicando que el estimador de error MLS utilizado para polyFEM tiene un comportamiento asintóticamente exacto. Este mismo comportamiento se había verificado anteriormente para FEM convencional, con elementos T3 y Q4 [20]. El índice de efectividad global muestra valores muy cercanos al valor ideal de , indicando un buen desempeño de la técnica de reconstrucción de tensiones para los tres tipos de elemento. Se observa que el estimador tiende a subestimar ligeramente el error exacto

                         

(a) FEM-T3 con 102 GDL. (b) FEM-Q4 con 102 GDL.             (c) polyFEM con 148 GDL.

 

Figura 4. Malla 1 para cada secuencia de mallas uniformemente refinadas. 

 

para elementos triangulares lineales (FEM-T3), y tiende a sobreestimar el error para elementos cuadriláteros lineales (FEM-Q4) y poligonales (polyFEM).

 

En la Figura 5 se muestra el índice de efectividad local, D, para las mallas utilizadas de elementos poligonales. Se observa que los valores locales se acercan al valor ideal de 𝐷 = 0. Se nota para este problema un balance de zonas donde el error exacto se sobreestima (azul) y zonas donde el error se subestima (rojo).  poligonales. Los valores de efectividad global son cercanos a los valores teóricos de 𝜃 = 1, y de un orden de magnitud similar al obtenido para mallas de elementos finitos convencionales. Asimismo, se ubican dentro del rango óptimo de 𝜃 = [0.8, 1.2].

 

El procedimiento de reconstrucción considera ecuaciones de restricción para introducir la satisfacción de las ecuaciones de equilibrio, con el fin de mejorar la calidad del campo reconstruido.

 

Tabla 1. Error estimado en norma energética 𝛆MLS, efectividad 𝜃 y velocidad de convergencia s. Resultados para mallas poligonales (polyFEM), mallas convencionales de elementos cuadriláteros (FEM-Q4) y triángulos (FEM-T3) lineales.

 

 

polyFEM

 

 

FEM-T3

 

 

FEM-Q4

 

GDL

𝛆MLS

𝜃

s

GDL

𝛆MLS

𝜃

s

GDL

𝛆MLS

𝜃

s

148

34.9524

1.0028

 

102

65.3052

0.9923

 

102

53.0947

1.0012

 

260

29.4995

1.0027

0.30

366

34.2888

0.9932

0.50

366

26.6543

1.0011

0.54

580

17.7977

1.0052

0.63

1374

17.5159

0.9944

0.51

1374

13.3733

1.0016

0.52

7060

5.2154

1.0004

0.49

5310

8.8539

0.9974

0.50

5310

6.6906

1.0011

0.51

 

 

4. CONCLUSIONES

 

En este trabajo se ha presentado una técnica de estimación del error de discretización para aproximaciones de elementos finitos poligonales, basada en reconstrucción de la solución. La técnica utiliza un procedimiento mejorado de MLS para la obtención del campo de tensiones reconstruido.   elementos Se utilizó un ejemplo numérico con solución exacta para poder evaluar la efectividad del estimador, y compararlo con el observado para otros tipos de aproximaciones numéricas.

 

 

 

Figura 5. Índice de efectividad local (D) en las mallas analizadas.


Los resultados obtenidos muestran un buen desempeño del estimador del error propuesto para mallas de  Los valores de los índices de efectividad local son cercanos al valor óptimo de 𝐷 = 0. La distribución de los valores de 𝐷 indica que no existen zonas con altas subestimaciones o sobrestimaciones del error a nivel local.

 

El estimador propuesto dio buenos resultados y se presenta como una buena opción para la extensión de este trabajo a problemas con singularidades, en 3D y para adaptatividad de mallas poligonales.

 

5. AGRADECIMIENTOS

 

OA González-Estrada agradece el apoyo a través del proyecto Capital Semilla 1742 de la Universidad Industrial de Santander.

 

 6. REFERENCIAS

 

M. Kraus, A. Rajagopal, and P. Steinmann, “Investigations on the polygonal finite element method: Constrained adaptive Delaunay tessellation and conformal interpolants,” Comput. Struct., vol. 120, pp. 33–46, 2013.

K. Jayabal, A. Menzel, A. Arockiarajan, and S.

M. Srinivasan, “Micromechanical modelling of switching phenomena in polycrystalline piezoceramics: application of a polygonal finite element approach,” Comput. Mech., vol. 48, pp. 421–435, 2011.

P.      Pavankumar,      K.     Jayabal,     and   A.

Arockiarajan, “A comparative study between finite element and polygonal finite element approaches for electromechanical coupled linear problems,” Intergrated Ferroelectr., vol. 120, pp. 90–101, 2010.

K. Y. Sze and N. Sheng, “Polygonal finite element method for nonlinear constitutive modeling of polycrystalline ferroelectrics,” Finite Elem. Anal. Des., vol. 42, no. 2, pp. 107–129, 2005.

E. L. Wachspress, “A rational basis for function approximation,” Lect. notes Math., 1971.

G. Dasgupta, “Interpolants within Convex

Polygons: Wachspress’ Shape Functions,” ASCE - J. Aerosp. Eng., vol. 16, no. 1, pp. 1–8, 2003.

N. Sukumar and E. A. Malsch, “Recent Advances in the Construction of Polygonal Finite Element Interpolants,” Arch. Comput. Methods Eng., vol. 13, no. 1, pp. 129–163, 2006.

M. Arroyo and M. Ortiz, “Local MaximumEntropy Approximation Schemes,” Lect. notes Comput. Sci. Eng., vol. 57, p. 1, 2006.

N. Sukumar and A. Tabarraei, “Conforming polygonal finite elements,” Int. J. Numer. Methods Eng., vol. 61, no. 2045–2066, pp. 2045–2066, 2004. [10] M. M. Rashid and P. M. Gullet, “On a finite element method with variable element topology,” Comput. Methods Appl. Mech. Eng., vol. 190, no. 11–12, pp. 1509–1527, 2000.

S. Ghosh and Y. Liu, “Voronoi cell finite element model based on micropolar theory of thermoelasticity for heterogeneous materials,” Int. J. Numer. Methods Eng., vol. 38, no. 8, pp. 1361–1398, Apr. 1995.

S. Ghosh and S. Moorthy, “Elastic-plastic analysis of arbitrary heterogeneous materials with the Voronoi Cell finite element method,” Comput. Methods Appl. Mech. Eng., vol. 121, no. 1–4, pp. 373–409, 1995.

X. H. Tang, S. C. Wu, C. Zheng, and J. H. Zhang, “A novel virtual node method for polygonal elements,” Appl. Math. Mech., vol. 30, no. 10, pp. 1233– 1246, 2009.

L. da Veiga, F. Brezzi, A. Cangiani, G. Manzini, L. D. Marini, and A. Russo, “Basic principles of virtual element methods,” Math. Model. Methods Appl. Sci., vol. 23, p. 199, Nov. 2013.

B. A. Szabó and I. Babuška, Finite Element Analysis. New York: John Wiley & Sons, 1991.

O. A. González-Estrada, J. Leal Enciso, and J. D. Reyes Herrera, “Análisis de integridad estructural de tuberías de material compuesto para el transporte de hidrocarburos por elementos finitos.,” UIS Ing., vol. 15, no. 2, pp. 105–116, 2016.

J. J. Ródenas, O. A. González-Estrada, J. E. Tarancón, and F. J. Fuenmayor, “A recovery-type error estimator for the extended finite element method based on singular+smooth stress field splitting,” Int. J. Numer. Methods Eng., vol. 76, no. 4, pp. 545–571, 2008.

O. A. González-Estrada, S. Natarajan, J. J. Ródenas, H. Nguyen-Xuan, and S. P. A. Bordas,

“Efficient recovery-based error estimation for the smoothed finite element method for smooth and singular linear elasticity,” Comput. Mech., vol. 52, no. 1, pp. 37– 52, Sep. 2013.

M. Ainsworth and J. T. Oden, A posteriori Error Estimation in Finite Element Analysis. Chichester: John Wiley & Sons, 2000.

J. J. Ródenas, O. A. González-Estrada, F. J. Fuenmayor, and F. Chinesta, “Enhanced error estimator based on a nearly equilibrated moving least squares recovery technique for FEM and XFEM,” Comput. Mech., vol. 52, no. 2, pp. 321–344, Aug. 2013.

R. Sibson, “A vector identity for the Dirichlet tessellation,” Math. Proc. Cambridge Philos. Soc., vol. 87, no. 01, p. 151, 1980.

O. C. Zienkiewicz and J. Z. Zhu, “A simple error estimator and adaptive procedure for practical engineering analysis,” Int. J. Numer. Methods Eng., vol. 24, no. 2, pp. 337–357, 1987.

G.     R.     Liu,   “MFree       Shape        Function

Construction,” in Mesh Free Methods. Moving beyond the Finite Element Method, Boca Ratón, Florida: CRC Press, 2003, p. 693.

P. Díez, J. J. Ródenas, and O. C. Zienkiewicz, “Equilibrated patch recovery error estimates: simple and accurate upper bounds of the error,” Int. J. Numer. Methods Eng., vol. 69, no. 10, pp. 2075–2098, 2007.

J. J. Ródenas, M. Tur, F. J. Fuenmayor, and A. Vercher, “Improvement of the superconvergent patch recovery technique by the use of constraint equations: the SPR-C technique,” Int. J. Numer. Methods Eng., vol. 70, no. 6, pp. 705–727, 2007.

J. J. Ródenas, O. A. González-Estrada, P. Díez, and F. J. Fuenmayor, “Accurate recovery-based upper error bounds for the extended finite element framework,” Comput. Methods Appl. Mech. Eng., vol. 199, no. 37–40, pp. 2607–2621, 2010.

A. Huerta, Y. Vidal, and P. Villon, “Pseudodivergence-free element free Galerkin method for incompressible fluid flow,” Comput. Methods Appl. Mech. Eng., vol. 193, no. 12–14, pp. 1119–1136, Mar. 2004.

Q. Z. Xiao and B. L. Karihaloo, “Statically


Reconstrucción de tensiones para el método de elementos finitos con mallas poligonales