MECANICA DE MATERIALES 2 EXPERIMENTACÍON, MODELADO NUMERICO Y TEÓRICO
Octubre de 2014
Institución editora Facultad de Ingeniería – Universidad Tecnológica de Bolívar
Los conceptos y opiniones de los artículos contenido de esta publicación son responsabilidad de sus autores; En ningún momento comprometen las orientaciones y políticas de la facultad de Ingeniería de la Universidad Tecnológica de Bolívar.
Contacto: Prof. Jairo F. Useche, Ph.D Departamento de Ingeniería Mecánica y Mecatrónica Universidad tecnológica de Bolívar Parque Industrial Vélez – Pombo, km.1 Tel/fax: +575 6535337/6619240
[email protected]
MECANICA DE MATERIALES 2 EXPERIMENTACÍON, MODELADO NUMERICO Y TEÓRICO
Dirección Universitaria Jaime Bernal Villegas Rector William Arellano Cartagena Vicerrector Académico Maria del Rosario Gutierrez de Piñeres Vicerrectora Administrativa Jorge Del Rio Director de Investigaciones, Emprendimiento e Innovación
Facultad de Ingeniería Jairo Useche Vivero Decano Facultad de Ingeniería Raul Padrón Carvajal Secretario Académico Edgardo Arrieta Ortiz Director departamento de Ingeniería Mecánica y Mecatrónica
MECANICA DE MATERIALES 2 EXPERIMENTACÍON, MODELADO NUMERICO Y TEÓRICO Cuerpo Editorial Jairo Useche Vivero Editor General Comité Editorial Prof. Jairo Useche, Ph. D Universidad Tecnológica de Bolívar Prof. Juan Pablo Casas, Ph.D., M.Sc. Universidad de los Andes Prof. Alejandro Marañón, Ph.D., M.Sc. Universidad de los Andes
Prof. Jairo Useche, Ph. D Prof. Jairo Cabrera, Ph.D Prof. Edgardo Arrieta, M.Sc. Prof. E.L. Albuquerque, Ph.D. Prof. Renato Pavanello, Ph.D. Prof. Paulo Sollero, Ph.D. Prof. Alejandro Marañón, Ph.D. Prof. Jose Rafael Toro, M. Sc.
Comité Científico Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar Universidad Tecnológica de Bolívar
Diagramación y Diseño Sharicar Méndez Villamizar Erick Guerrero
CONTENIDO
A TRULY MESHLESS ELEMENT-FREE GUERIN METHOD BASED ON A NOVEL POINT IN POLYGON ALGORITHM. Autores: Kevin E. Patrón Hernández, Edgardo Arrieta. Universidad Tecnológica de Bolívar.…………………………………………….………………………………………………………………. 8 ANÁLISIS DE MODELOS DE DAÑO PARA LA CARACTERIZACIÓN DE DUCTOS DE GAS Y PETRÓLEO. Autores: Johnny Ortiz, Henry Ávila, Arístides Távara. Universidad Nacional de Trujillo.……………………………………………………………………………………………….…………… 18 ESTUDIO DE LAS PROPIEDADES BAJO IMPACTOS A BAJA VELOCIDAD DE UN HONEYCOMB CON ESFUERZO INTERNO. Autores: Pablo Andrés Montoya Macías, Juan Pablo Casas Rodriguez, Universidad de los Andes ………………………………………………………………………………………………………. 27 CONEXIÓN CON PERFORACIONES MÚLTIPLES EN MATERIALES ORTO TRÓPICOS. Autores: Cesar Echavarría, Beatriz Echavarría. Universidad Nacional Colombia……………..…………………………………………………………………………………..
de
35
ESTUDIOS DE LOS FENÓMENOS DE CONO DE AGUA Y GAS EN YACIMIENTOS DE PETRÓLEO. Autores: Gustavo Gontijo, A.V. Diaz Jr., E.L. Albuquerque, E.L.F. Fortaleza. Universidad de Brasilia.……………………………………………………………………………………………..………. 45 OBTENCIÓN DE CINTAS AMORFAS MEDIANTE LA TÉCNICA DEL MELT SPINNING VS. SIMULACIÓN POR VOLÚMENES FINITOS: ESTUDIO COMPARATIVO. Autores: Marcelo Barone, Soledad Gamarra, Marcelo Pagnola. Universidad de Buenos Aires………………………………………………………………………………………………………….. 53 SIMULACIÓN DEL COMPORTAMIENTO DE HONEYCOMB DE ALUMINIO SOMETIDA A CONDICIONES CUASI-ESTÁTICAS. Autores: Camilo Orrego, Juan Pablo Casas. Universidad de los Andes.…………………………………………………………………………………….…………….. 60 RESPUESTA ESTRUCTURAL DE PANELES TIPO SANDWICH CON NÚCLEO DE POLIURETANO SOMETIDO A IMPACTOS DE BAJA VELOCIDAD. Autores: Jesus Aguiar, Jairo Useche. U iversidad Te ológi a de Bolívar…………………………………………………………… 68 EVALUACIÓN DEL DESEMPEÑO DE ESPUMAS POLIMÉRICAS DE CELDAS ABIERTAS RELLENAS CON STF ANTE EVENTOS EXPLOSIVOS. Autores: Vicent Steve Robinson, Juan Pablo Casas. Universidad de los Andes………………………………………………………………………………………………………….. 75
MODELACIÓN DE UN MATERIAL ELASTÓMERO DE COMPORTAMIENTO NO-LINEAL APLICADO EN UN FENÓMENO DE INTERACCIÓN FLUIDO-ESTRUCTURA (F-E). Gustavo Suarez. Universidad Pontificia Bolivaria a……………………………………… 85 ANÁLISIS DE ESFUERZOS Y DEFORMACIONES DE GALERÍAS EN EXPLOTACIÓN SUBTERRÁNEA DE CARBÓN. Sebastián Jaramillo, Galin Quitian, Oswaldo Bustamante. Universidad Nacional de Colo ia………………………………………………………………………………………................ 92 ANÁLISIS DE FATIGA POR IMPACTO DE UNIONES ADHESIVAS. Juan Fernando Téllez, Juan Pablo Casas, Peter Alvarado. Universidad de los Andes…………………………………………………………………………………………………………… 103 REDES NEURONALES APLICADAS AL ENTRENAMIENTO DE RELACIONES CONSTITUTIVAS LINEALES Y NO-LINEALES. Fredy Mercado, Universidad de Buenos Aires. …………………………………………….. 114
FINITE ELEMENT ANALYSIS OF HIGH SPEED COMPOSITE CRAFT. Hugo Alvarez. Cotecmar……………………………………………………………………………….
127
ANÁLISIS NUMÉRICO DE IMPACTO HIDRODINÁMICO USANDO UN MÉTODO DE PARTÍCULAS. Autor: Jairo Cabrera. Universidad Tecnológica de Bolívar……………………………. 140 DYNAMIC CHARACTERIZATION OF HEXAGONAL CELL STRUCTURES. Autores: Diego Avendaño, Edgar Marañón, Juan Pablo Casas. Universidad de los Andes…………………………………………………………………………………………………………. 153 ANÁLISIS POR DEFORMACIÓN PLANA DE PELÍCULAS DELGADAS MULTICAPAS MEDIANTE EL MÉTODO DE LOS ELEMENTOS FINITOS. Emerson Escobar, Antonio Ramírez, Dubernoys Ramírez, Carlos Vidal. Universidad Autó o a del O ide te………………………………………………………………………… 167
ESTUDIO EXPERIMENTAL Y NUMÉRICO DE PANELES TIPO SANDWICH CON NÚCLEOS HONEYCOMB SOMETIDOS A IMPACTOS DE BAJA VELOCIDAD. Luis Gar ia, Jairo Use he, U iversidad Te ológi a de Bolívar…………………. 180
DEGRADACIÓN DE LAS PROPIEDADES MECÁNICAS DE COMPUESTOS LAMINADOS FABRICADOS CON FIBRA DE VIDRIO POR INMERSIÓN SALINA. Mary Arias Departamento de Ingeniería Mecánica Universidad Tecnológica de Bolívar, Colombia………………………………………………………………………………………………….. 185 MÉTODO DE INTEGRACIÓN POR CELDAS PARA ANÁLISIS DE LAMINADOS GRUESOS UTILIZANDO – BEM. Christian Harnish Departamento de Ingeniería Mecánica Universidad Tecnológica de Bolívar, Colombia……………………………………………………………………………………… 192 CREEP EN ACEROS HP MODIFICADOS. G. González Universidad Tecnológica Nacional; Bahía Arge ti a…………………………………………………………………………………………….......
Blanca, 201
UNA APROXIMACIÓN PRELIMINAR AL COMPORTAMIENTO DEL MIEMBRO INFERIOR FRENTE A CARGAS IMPULSIVAS MEDIANTE EL MÉTODO DE LOS ELEMENTOS FINITOS. J. calle Estudiante, Universidad de Los Andes, Bogotá Colombia………………….. 209 MÉTODO DE ENJAMBRE DE PARTÍCULAS Y EVOLUCIÓN DIFERENCIAL PARA LA OPTIMIZACIÓN DE PANELES ESTRUCTURALES. Guillermo E. Giraldo Grupo De Investigación En Materiales Y Estructuras GIMAT, Universidad Tecnológica De Bolívar…………………………………………………………………………………… 214
A truly meshless element-free Galerkin method based on a novel “point in polygon” algorithm Kevin E. Patrón , Edgardo W. Arrieta Department of Mechanical Engineering, Universidad Tecnológica de Bolívar, Colombia Abstract It is proposed a truly meshless element-free Galerkin method (EFGM) based on a novel algorithm that solves the “point in polygon” test efficiently for massive queries, and a fixed integration grid that contains a cloud of Gauss points and is completely independent of the nodal distribution. The two-dimensional physical domain is projected onto this fixed grid and numerical integration of the weak form is performed only over cells lying inside the domain and also over partially cut cells. Our novel point in polygon algorithm indirectly recognizes the shape of the domain boundaries and automatically selects the appropiate Gauss points whose contribution to the system matrices are significant. Thus, this approach turns the usual mesh-based EFGM into a “truly meshless” one. The integration accuracy may be easily modified by handling several parameters of this fixed grid and by recursively partitioning the partially cut cells. Numerical experiments exhibit a great performance by the implemented scheme in different applications, from solving geometric problems up to two-dimensional linear elasticity. The presented methodology may result very useful for movable boundary problems, such as those involving crack growth in fracture mechanics and structural shape optimization, since it avoids any remeshing process.
Keywords: search
1
meshless methods, element-free Galerkin method, fixed grid, point in polygon algorithm, k-nearest neighbor
Introduction
Over the last decade, a group of meshfree methods have been emerging to solve boundary value problems while trying to alleviate several disadvantages that were becoming increasingly evident in the well developed and mature Finite Element Method (FEM). The minimum requirement for these meshfree methods is that a predefined mesh is not necessary for the field variable approximation. The element-free Galerkin method [1] (EFGM) belongs to this category. Its shape functions create a “diffuse” (smooth) approximation which is not linked to integration cells or elements, although makes use of background cells (which can be independent of the field nodes) with the purpose of performing numerical integration of the weak form over the problem domain. Normally, the integrals over the domain are computed using Gaussian integration over a set of background cells that overlaps the domain exactly [2], so any typical FEM mesh is applicable under some requirements, and even the mesh vertices may be used as field nodes. It has been discovered that one main source of error associated with the cell-based integration scheme is precisely that the local support domains of the EFG shape functions may not align with the background cells1 . A natural path to more accurate integrations is the construction of integration cells which align with the shape function support domains. Thus, discontinuities of the integrand are avoided within quadrature domains and the integration of the weak form is considerably improved. However, when high-order Gauss quadrature is used in each cell, the integration error was reported to be insignificant [3]. Several approaches have been proposed with the aim of avoiding the background cells completely and turning the conventional EFGM into a very truly meshless method. For example, a nodal integration of the EFGM has been 1 In FEM, Gaussian quadrature is used in each element, and the supports of the shape functions corresponds to the cells used for integration
1
2 A novel “point in polygon” algorithm
2
proposed, in which the integrals of the weak form are evaluated only at the nodes [4]. This nodal integration suffered from spatial instability that resulted from underintegration of the weak form, whereby an additional stabilization term should be included in the problem functional to achieve a reasonable rate of convergence. In [5] a strain smoothing stabilization for nodal integration was proposed to eliminate spatial instability in nodal integration. On the other hand, in [6] it was presented a new body integration technique that relies on a partition of unity by a set of moving least squares shape functions, each defined on a small patch that belongs to a set of overlapping patches covering the domain (moving least squares quadrature). Anyway, numerical quadrature was used to perform the integrations over those patches and on their intersections with the domain (for those that were cut by the boundary), and an algorithm for the patches-intersection test had to be used. Unfortunately, completely removing the integration mesh makes the integration process, in general, much more difficult. The two-dimensional meshless EFGM that we propose here relies on a fixed integration grid totally independent of the nodal distribution, and a novel algorithm that solves the “point in polygon” test efficiently for massive queries, together working as follows. At first, the physical domain is projected onto this fixed grid, which is a cell-based structure that contains a cloud of Gauss points for integration purposes. In principle, this grid can be structured or unstructured, although for sake of simplicity a fixed structured grid would be prefered 2 . Some cells are expected to be outside the domain, while many other cells are completely inside the domain and a few cells partially cut by the boundary domain. Numerical integration of the weak form is performed only over the two latter. At last, a very simple integration scheme is carried out, which relies on taking into account only the contributions (to the system matrices) from those Gauss points that lie inside the physical domain. Our new point in polygon algorithm, based on an efficient k- nearest neighbor searching system, is able to speed up significantly this Gauss points selection3 . Our algorithm indirectly recognizes the shape of the domain boundaries by discarding integration points outside the physical domain. Thus, despite the presence of the fixed grid, this whole approach turns the conventional mesh-based EFGM into a truly meshless one. Movable boundary problems, such as crack growth in fracture mechanics and structural shape optimization, can be easily treated with this methodology thanks to its versatility in tracking the shape of the boundary without any remeshing process.
2
A novel “point in polygon” algorithm
The “point in polygon” test consists in determining whether a given point in the plane lies inside, outside, or on the boundary of a polygon. The simplest algorithm that solves this test in linear time O(n) is the “Ray casting algorithm” or “even-odd algorithm”, which solves the problem by casting a ray from the test point and counting how many edges of the polygon the ray intersects. Another typical linear algorithm is the “winding number algorithm” or “angle summation algorithm” [7]. Several variants of this two algorithms have been proposed, as seen in [8] and [9] . Many other more sophisticated algorithms have been presented, some of them make use of a background mesh and wall-sharing relationships, as in [10] and [11]. In the present section we propose a new algorithm that also solves this test efficiently and is very suitable for a great number of queries with the same polygon.
2.1
The main ideas
The underlying ideas behind this simple algorithm are mainly the so-called “orientation theorem”, and some criteria for classifying a point as inside or outside the polygon. Finding the nearest neighbor vertex is a crucial step in the latter idea, as will be discussed below. 2.1.1
The orientation theorem
This theorem is only applicable for strictly convex polygons, and states that, if a polygon is defined by a collection of vertices in anticlockwise order, then a point which is inside the polygon will always lie to the left of each side [12] . It’s easy to check wheter a point P (xp , yp ) lies to the left of a given side defined by the vertices Vi (xi , yi ) and Vi+1 (xi+1, yi+1 ) through the signed area of the triangle Vi Vi+1 P , defined as follows 2 3
the structured grid can be created very easily with the bounding box of the domain compared to brute-force algorithms
3
2 A novel “point in polygon” algorithm
1 (xi − xp ) AT = 2 (yi − yp )
(xi+1 − xp ) (yi+1 − yp )
If AT is greater than zero, P lies to the left of the side. If AT is equal to zero, P lies on the side. If AT is less than zero, P is to the right of the side.
P1 lies to the left of the side Vi Vi+1 , P2 is on Vi Vi+1 , and P3 is to the right of Vi Vi+1 . Fig. 1: A polygon side 2.1.2
Criteria for classifying a point
In the present subsection the criteria for classifying a point (as inside, outside or on the boundary of a polygon) will be presented. Let us consider a simple polygon formed by the union of line segments, each of these defined by two vertices. Let us call a “loop” to a pair of line segments of the polygon that are joined together by a central vertex. In the most general case in a simple polygon, there will be some “convex loops”, also there will be some “concave loops” and some others will be “straight loops”, depending on the shape that they describe.
V3 V4 V5 is a straight loop. V8 V9 V1 is a concave loop. V6 V7 V8 is a convex loop. Fig. 2: A simple polygon and several test points
4
2 A novel “point in polygon” algorithm
Generally, a test point will always be far from some loops, but at the same time close to another loops. It should be noted that, the fact that this point is inside or outside the polygon depends directly on certain geometric relationships between the point and its “ nearest loop”, specifically on the relative position between the point and its nearest loop, and on the shape of this loop (concave, convex or straight). Let us begin to solve the problem by identifying the nearest loop to this test point. One of the ways to achieve this is through the determination of the nearest neighbor vertex 4 , which would be the central node of the nearest loop. Once the nearest neighbor vertex has been found, its succesor and predecessor vertices can be easily identified with the vertex enumeration5 . Thus, these three vertices form the nearest loop. On the other hand, it is very simple to determine whether a loop is convex, concave or straight. Assuming the external boundary nodes are listed as mentioned above, the signed double area of the triangle formed by the vertices Vi−1 Vi Vi+1 indicates the shape of the loop , as follows xi−1 AP = xi xi+1
yi−1 yi yi+1
1 1 1
If Ap es greater than zero, the loop Vi−1 Vi Vi+1 is convex. If Ap is less than zero, the loop is concave. If Ap is equal to zero, the loop is straight.
Fig. 3: A convex loop (left) and a concave loop (right) Based on the orientation theorem, it is easy to check that, for a test point whose nearest loop is convex, the only condition for this point to be inside the polygon is that it lies to the left of both sides of the nearest loop. Otherwise, it may be out of the polygon or on its edge, depending on the values of the signed areas for each side of the loop. This statement is perfectly valid for a point whose nearest loop is straight. For example, in figure 2, the nearest loop for point P1 is the convex loop V1 V2 V3 . This point lies to the left of both sides , so it’s inside the polygon. On the other hand, the nearest loop for point P3 is the convex loop V9 V1 V2 . This point lies to the left of side V9 V1 but to the right of side V1 V2 , so it’s outside the polygon. The nearest loop for point P2 is the straight loop V3 V4 V5 . This points lies to the left of both sides, so it’s inside the polygon. The case where a point is on a side (or on the intersection of two sides, which means that coincides with a vertex) can be easily treated according to the orientation theorem. Now, for a point whose nearest loop is concave, if it lies to the left of at least one of the sides, then the point is inside the polygon. If the point lies on at least one side, then necessarily the point is on the edge. Otherwise, the point is out of the polygon. For example, in figure 2, the nearest loop for point P4 is the concave loop V5 V6 V7 . This point lies to the left of side V5 V6 , so it’s inside the polygon. On the other hand, the nearest loop for point P5 is the concave loop V8 V9 V1 . This point lies to the right of both sides, so it’s outside the polygon. By these criteria, any point can be easily classified as inside, outside, or on the boundary of a polygon. 4 5
finding the nearest neighbor vertex will be discussed in the next subsection vertices should be enumerated in increasing order in anticlockwise
3 A truly meshless element-free Galerkin method
2.1.3
5
Finding the nearest neighbor vertex
The simplest algorithm to find the nearest neighbor vertex consists in computing the Euclidean two-dimensional distances from the test point to each of the polygon vertex and then selecting the minimum distance. The vertex corresponding to this minimum distance will be the nearest neighbor vertex. This exhaustive search has the advantage of being trivial to implement and exhibits great performance when the polygon has few number of vertices. It uses O(n) storage and the nearest neighbor vertex search is done in linear time O(n). Unfortunately, this naive approach is not suitable for queries on polygons formed by a great number of vertices. An efficient data structure to perform the nearest neighbor vertex search is the k-d tree [13]. It is a space-partitioning data structure for organizing points in a k-dimensional space. The k-d tree is a binary tree in which every node is a k-dimensional point. With the aid of the distance bound and bounding box of each node, it is possible to quickly prune parts of the tree that could not include the nearest neighbor A kd-tree of a set of n nodes uses O(n) storage and can be constructed in O(n log n) time. With this data structure, the general k-nearest neighbor search can be done in O(k log n) time. In our meshless EFG implementation it has been included a Fortran 95 software called KDTREE2, developed by Matthew B. Kennel [14] , and licensed under the terms of the Academic Free Software License. It’s an open source package for neighbor searching in Euclidean space with k-d trees. This is our k-nearest neighbor nodes searching system, useful for the point in polygon algorithm.
2.2
The algorithm
The general algorithm has two stages, a preprocessing phase and the processing phase. First, a kd tree may be created (only once) with the polygon vertices, so the nearest neighbor vertex searching is speed up significantly. This process is even less time-consuming when many searchings are performed with the same polygon vertices. Secondly, it is necessary to determine (only once) the succesor and predecessor vertices for all of the polygon vertices, which can be done in linear time. This information can be stored easily in an array with the aid of the vertices enumeration. At this point, the preprocessing phase is finished. For a polygon formed by v vertices, the storage is O(v) and the preprocessing time is O(v log v). The processing phase consist in the test execution. First, given a test point, its closest loop is determined, and also its shape. Second, and lastly, a simple conditional statement (based on the shape of the nearest loop) can quickly solve the test. The expected query time with a polygon formed by v vertices and p test points is O(p log v). In contrast, with a linear “point in polygon” algorithm, the expected query time is O(pv). The inclusion of holes in the polygon requires a slight modification in the way the enumeration is carried out over those nodes belonging to internal borders. For these nodes, the enumeration should be in increasing order in clockwise. In this way, the previous criteria still holds and present no loss of generality. As can be seen, the present algorithm is very suitable for performing massive queries with the same polygon.
3
A truly meshless element-free Galerkin method
In this section it is presented a general overview of the element-free Galerkin method and our proposed truly meshless EFG scheme.
3.1
The element-free Galerkin method
The element-free Galerkin method (EFGM) is a meshfree method developed by T. Belytschko et. al. in 1994 [1], based on the diffuse elements (DEM) method proposed by Nayroles in 1992 [15]. The main features of the EFGM are as follows [2] : 1. A moving least square approximation [16] is employed for the construction of the approximation functions6 . 2. A Galerkin weak form is used to develop the discretized system equation. 3. Background cells for integration are required to carry out the intregrations involved in the weak form. 6
Also called shape functions
6
3 A truly meshless element-free Galerkin method
3.1.1
The EFG approximation scheme
The moving least square approximation uh (x) of the function u(x) is defined in the domain Ω by uh (x) =
m X
T
pj (x)aj (x) = p(x) a(x)
j
where m is the number of terms in the polynomial basis, p(x)is a basis functions (generally monomials) , and a(x) is a vector of coefficients. a(x) is obtained at any point x by minimizing a weighted, discrete L2 norm as follows J=
n X i
h i2 T w(x − xi ) p(xi ) a(x)−ui
where n is the number of points in the neighborhood of x for which the weight function w(x − xi ) 6= 0, and ui is the nodal value of u at x = xi . The stationariy of J with respect to a(x) leads to the following a(x) = A−1 (x)B(x)u where A and B are the matrices defined by A(x) =
n X i
T
w(x − xi )p(xi ) p(xi )
B(x) = [w(x − x1 )p(x1 ),w(x − x2 )p(x2 ), . . . , w(x − xn )p(xn )] Hence, the resulting EFG approximation to a function u(x) and its ith shape function are respectively as follows h
u (x) =
m n X X i
pj (x)(A
−1
φi (x)ui
i
j
φi (x) =
(x)B(x))ji ui =
n X
m X
pj (x)(A−1 (x)B(x))ji
j
3.1.2
Background cells for integration
In the conventional EFGM it is necessary to establish a cell structure, independent of the nodal configuration, in order to perform the integrations derived from the weak form system equation. In each cell, Gauss quadrature is used. These non-overlapping cells are generally arranged in a regular pattern, and provide a structure for the numerical integration. As stated by Beltyschko in [1], the number of quadrature points depends on the number of nodes √ in a cell: if m is the number of nodes in a cell, a nG × nG Gauss quadrature could be used, where nG = m + 2. Liu in [17] indicates that the sufficient requirement on the total number of quadrature points nQ is that nQ = (3 ∼ 9)n , where n is the total number of unfixed field nodes 7 , and plus, nQ should be at least 2/3 of n. It is well known that these cells should overlay the whole domain, so any mesh8 similar to an unstructured FEM mesh is applicable9 . In [2] the author recommends to employ the triangular background mesh generated on a triangulation tecnique, and the mesh vertices may be used as field nodes in the problem domain. Nevertheless, structured grids might be used, provided that the integrations are performed correctly over the real physical domain10 . 7 8 9 10
For 2D problems It’s possible to state that background cells constitute a mesh No connectivity is required, unlike FEM The following images were taken from [18] and [19]
3 A truly meshless element-free Galerkin method
7
Fig. 4: Classification of grids
3.2
Our meshless EFG approach
The proposed meshless EFG approach relies on a simple integration scheme that makes use of our “point in polygon” algoritm (see the previous section) and a fixed grid, which is created from the bounding box of the two-dimensional problem domain. The cell size initially may be set according to the requirements exposed above, although the cell density in the grid may be easily increased. The physical domain and its nodal representation is projected onto this grid, as shown in figure 5.
Fig. 5: An example of a fixed grid in the integration scheme 3.2.1
The 2D integration scheme
In order to solve a boundary value problem, its Galerkin weak form is used to develop the discretized system equation. The weak form involves performing integration over the physical domain and its boundaries. As the domain is overlapping a fix grid that contains a cloud of integration points, it is necessary to bound the integration limits to the physical domain. There will be cells completely inside the domain, other cells will be completely outside the domain, and a few cells will be cut by the boundary curve, as shown in figure 6. Having this in mind, the integration scheme relies on a simple idea : Gauss points outside the physical domain are left out from the computation. The contributions to the system matrices come only from Gauss points inside the domain. The boundary nodal representation is enough to perform the boundary integrals (by Gaussian quadrature) that are involved in the Galerkin weak forms. This integration scheme was discussed in detail in [20] . Evidently, an integration error was induced by cutting the integration cells and using only those Gauss points that fall inside the domain. In [20], convergence results
8
4 Partial numerical results
were presented by computing the relative error against the grid spacing. Superlinear to quadratic convergence is achieved in numerically integrating this cutting strategy. Our novel point in polygon algorithm is able to solve efficiently the Gauss points selection. A kd- tree is created with the boundary nodes in order to represent the domain, which may be treated as a polygon. The Gauss point density can be increased by reducing the cell size and by increasing the Gauss quadrature order. It is also possible to recursively partition the partially cut cells in order to obtain more accurate results in the weak form integration.
Fig. 6: The two-dimensional integration scheme This approach turns the usual mesh-based EFGM into a truly meshless EFGM because the novel algorithm indirectly recognizes the shape of the domain boundaries by discarding integration points outside the physical domain. Thus, despite the presence of the fixed grid, this whole approach turns the conventional mesh-based EFGM into a truly meshless one. The only requeriment is a typical nodal boundary representation, as shown in figure 6, the algorithm will do the rest of the process automatically.
4
Partial numerical results
A two-dimensional linear elasticity problem was solved with the proposed scheme. The weak form of this problem includes a penalization term (due to the lack of Kronecker delta property of the shape functions) in order to impose the boundary conditions. The typical cantilever beam subjected to a parabolic shear stress distribution is studied, because it is often used for benchmarking numerical methods since its analytic solution is known.
Fig. 7: Cantilever beam
9
5 Conclusions
The parameters for this cantilever beam are the following: loading (integration of the distributed traction) P = −1000, Young’s modulus : E = 3 × 107 , Poisson’s ratio : ν = 0.3, its height : D = 12, its length : L = 48, and a unit thickness. A total of 7 × 25 nodes were used for the domain and boundary representation. For the numerical integration a total of 4 × 10 cells were used. In each cell, a 4 × 4 Gaussian quadrature was used. A plane stress problem is considered. The cinematic constrains are considered as shown in the figure 7. The figure 8 shows the initial nodal representation employed for the problem domain. Besides, the analytical solution (final configuration of the beam) and the numerical solution obtained by the proposed scheme. Compared with the analytical result, it is observed that this meshless EFG approach produces very good results11 .
Fig. 8: Meshless EFG solution
5
Conclusions
A truly meshless element-free Galerkin method has been presented. A novel “point in polygon algorithm” has been developed, which indirectly recognizes the shape of the domain boundaries and automatically selects an appropiate number of Gauss points lying inside the physical domain, with the aim of performing numerical integration of the weak form over the two-dimensional domain, and line integrals over the bondaries. A flexible cloud of Gauss points is defined by a fixed grid completely independent of the nodal representation. In general, with the assistance of this algorithm it is possible to efficiently perform numerical integrations over any domain described by a meshfree nodal representation. The new meshless approach exhibits great performance in different applications, mainly on solving geometric problems and boundary value problems with meshless methods. The presented methodology may result very useful for movable boundary problems, such as those involving crack growth in fracture mechanics and structural shape optimization, since it avoids any remeshing process. A natural path to improve the proposed methodology would be to implement an algorithm that identifies the partially cut cells and recursively subdivides them12 to perform more detailed and accurate integrations.
References [1] Belytschko, T., Lu, Y.Y. and Gu, L. “Element-Free Galerkin Methods”. International Journal for Numerical Methods in Engineering, Vol. 37, 229-256 (1994) 11 12
at least in the displacements solution similar to the space partitioning carried out by quadtrees, which is easier than a refinement based on a triangulation algorithm
5 Conclusions
10
[2] Liu, G.R.(2010). “Meshfree Methods. Moving beyond the Finite Element Method”. (2nd. edition) . CRC Press, Taylor and Francis Group. [3] Dolbow, J., Belytschko, T. “Numerical Integration of the Galerkin Weak Form in Meshfree Methods”. Northwestern University. (1998) [4] Beissel, S. and Belytschko, T. “Nodal integration of the element-free Galerkin method”. Comput. Methods Appl. Mech. Engrg. 139 (1996) [5] Chen, J. , Wu, C., Yoon, S., You, Y. “A stabilized conforming nodal integration for Galerkin mesh-free methods”. International journal for numerical methods in engineering. (2001) [6] Duflot, M. , Nguyen-Dang, H. “A truly meshless Galerkin method based on a moving least squares quadratures”. Communications in Numerical Methods in Engineering (2000) [7] Alciatore, D. “A winding number and point-in-polygon algorithm”. Colorado State University. [8] Galetzka, M., Glauner, P. “A correct eved-odd algorithm for the point-in-polygon (PIP) problem for complex polygons”. (2012) [9] Gatilov, S. “Efficient angle-summation algorithm for point inclusion test and its robustness”. A. P. Ershov Institute of Informatic Systems, Novosibirsk, Russia. (2012) [10] Jin, Y. “An efficient inclusion test for a massive point distribution”. Lawrence Livermore National Laboratory, California, USA. (2013). [11] Zalik, B., Kolingerova, I. “A cell-based point in polygon algorithm suitable for large sets of points”. Computers & Geosciences. (2001) [12] Sloan, S. W. “A point in polygon program”. Advance Engineering Software, Vol. 7, No. 1. (1985) [13] Bentley, J.L. “Multidimensional Binary Search Trees Used for Associative Searching”. Stanford University. (1975) [14] Kennel, M. B. “KDTREE2: Fortran 95 and C++ software to efficiently search for near neighbors in a multidimensional Euclidean space”. Institute for Nonlinear Science, UC San Diego. (2004) [15] B. Nayroles, G. Touzot and P. Villon, “Generalizing the finite element method diffuse approximation and diffuse elements”, Comput. Mech., 10, 307-318 (1992). [16] Lancaster P. , Salkauskas K. “Surfaces Generated by Moving Least Squares Methods”. Mathematics of Computation, Vol. 37, No. 155. (1981) [17] Liu, G. R. , Gu. Y. (2005). “An Introduction to Meshfree Methods and Their Programming”. Springer [18] "Structured Grid" by Shyam2791 Own work. Licensed under Creative Commons Attribution-Share Alike 3.0 via Wikimedia Commons http://commons.wikimedia.org/wiki/File:Structured_Grid.PNG#mediaviewer/File:Structured_Grid.PNG [19] "Unstructured grid" by Shyam2791 Own work. Licensed under Creative Commons Attribution-Share Alike 3.0 via Wikimedia Commons http://commons.wikimedia.org/wiki/File:Unstructured_grid.PNG#mediaviewer/File:Unstructured_grid.PNG [20] Bobaru, Florin Ph.D. and Rachakonda, Srinivas, "E(FG)2: A NEW FIXED-GRID SHAPE OPTIMIZATION METHOD" (2006). Faculty Publications from the Department of Engineering Mechanics. Paper 75.
MODELACIÓN DE UN MATERIAL ELASTÓMERO DE COMPORTAMIENTO NO LINEAL APLICADO EN UN FENÓMENO DE INTERACCIÓN FLUIDO-ESTRUCTURA (F-E) Suárez G. Grupo de Matemáticas, Universidad Pontificia Bolivariana, Medellín, Colombia
[email protected]
Resumen: Se presenta el desarrollo matemático y computacional de un material elastómero de comportamiento hiperelástico polinomial el cual ha sido implementado para representar la dinámica estructural sometida a interacción fluido-estructura (F-E). Se calculó el comportamiento de la función (esfuerzo-deformación) para obtener los parámetros funcionales que componen el modelo hiperelástico a elegir. La aplicabilidad del material fue ensayada en un fenómeno de interacción F-E, donde la parte estructural se desarrolló con el material elastómero calculado. Se realizaron simulaciones computacionales que reflejaron funcionalidad mecánica de la estructura, con valores del desplazamiento del fluido de 1.5 m/s y una compresión estructural mayor al 40%. Palabras claves: Modelación matemática, material elastómero, interacción fluido-estructura (F-E), Ley fenomenológica, simulación computacional. INTRODUCCIÓN Entre estos materiales existentes, se destacan los materiales hiperelásticos, cuyo comportamiento no lineal requiere calcular el comportamiento de la función (esfuerzo-deformación) para obtener los parámetros funcionales que componen el modelo hiperelástico a elegir [1-4]. Las leyes constitutivas de los materiales hiperelásticos se derivan principalmente del potencial elástico de energía de deformación y presentan un comportamiento en grandes deformaciones con efecto reversible [1-4]. Esta recuperación es la memoria del material. MATERIALES Y MÉTODOS Se postula que existe cierta función W E de densidad de energía libre (o energía de deformación) de componentes del tensor de deformaciones, definida como [1-4]: ij
S ij
W E ij
2
W C ij
(1)
Dónde: S = Los componentes del segundo tensor de tensiones de Piola-Kirchhoff y que actúa como potencial de tensiones. W = La función de la energía de deformación por unidad de volumen no deformado. E = Los componentes del tensor de deformación Lagrangiana. C = Los componentes del tensor de deformación de Cauchy-Green. ij
ij
ij
En la presente investigación, las funciones esfuerzo-deformación para construir el modelo fueron desarrolladas a partir de datos experimentales de laboratorio, las cuales se graficaron obteniendo las funciones de regresión, de donde se obtienen parte de los parámetros, fig. 1. Se determinaron bases de datos experimentales de investigaciones referentes a la utilización y caracterización de elastómeros para representar fenómenos F-E [1-4]. La información requerida para el desarrollo de los cálculos se presenta en las tablas siguientes. Se obtuvieron los datos experimentales para establecer el comportamiento del material bajo estado uniaxial, ver tabla 1.
Tabla 1. Datos de Esfuerzo [Pa] vs. Deformación uniaxial. Deformación 0,1338 0,2675 0,3567 0,6242 0,8917 1,1592 1,4268 2,051 2,586 3,0318 3,7898 4,3694 4,8153 5,172 5,4395 5,707 5,9299 6,0637 6,1975 6,3312 6,465 6,5541 6,6433
Esfuerzo 10691 16800 21383 29019 36656 41238 47347 61093 73311 85530 1,11E+05 1,37E+05 1,62E+05 1,89E+05 2,14E+05 2,38E+05 2,64E+05 2,90E+05 3,15E+05 3,41E+05 3,67E+05 3,93E+05 4,43E+05
Se determinó el interpolador de la forma Polinomial N u 3 , que estableció el comportamiento del material bajo estado uniaxial, y se obtuvieron las constantes del modelo hiperelástico, parámetros que fueron requeridos dentro de la formulación de la ley constitutiva aplicada en las simulaciones, ver fig. 1: σ
Esfuerzo vs Deformación uni-axial
5,0E+05 4,5E+05 4,0E+05 3,5E+05
3,0E+05 2,5E+05 2,0E+05 1,5E+05 1,0E+05 5,0E+04
0,0E+00 0,00
1,00
2,00
3,00
Ecuación calculada: y =
ε
3495.6x3
4,00
-
24165x2
5,00
6,00
7,00
+ 68175x - 2707.9
-------y (Función tendencia) ------- vs σ (Función experimental)
Figura 1. Gráficos de Esfuerzo vs. Deformación uniaxial y descripción obtenida por métodos matemáticos de regresión Polinomial.
Se obtuvieron los datos para establecer el comportamiento del material bajo estado biaxial, ver tabla 2: Tabla 2. Datos de Esfuerzo [Pa] vs. Deformación biaxial. Deformación 2,00E-02 6,00E-02 0,11 0,14 0,2 0,31 0,42 0,68 0,94 1,49 2,03 2,43 2,75 3,07 3,26 3,45
Esfuerzo 6470 10963 16607 18078 22918 30529 35736 45522 53637 67470 87116 1,01E+05 1,20E+05 1,39E+05 1,55E+05 1,70E+05
Se determinó el interpolador de la forma Polinomial Nu 3 , que estableció el comportamiento del material bajo estado biaxial, y se obtuvieron las constantes del modelo hiperelástico, parámetros que fueron requeridos dentro de la formulación de la ley constitutiva aplicada en las simulaciones, ver fig. 2: σ
Esfuerzo vs Deformación biaxial
1,8E+05 1,6E+05 1,4E+05
1,2E+05 1,0E+05 8,0E+04 6,0E+04
4,0E+04 2,0E+04 0,0E+00 0
0,5
1
1,5
2
ε
2,5
3
3,5
4
-------y (Función tendencia) ------- vs σ (Función experimental)
Ecuación Calculada: y = 6390.9x3 - 28777x2 + 70539x + 8476.3
Figura 2. Gráficos de Esfuerzo vs. Deformación biaxial y descripción obtenida por métodos matemáticos de regresión Polinomial. Se obtuvieron los datos para establecer el comportamiento del material bajo esfuerzos cortantes, ver tabla 3:
Tabla 3. Datos de Esfuerzo [Pa] vs. Deformación cortante. Deformación 0,1034 0,1724 0,2828 0,4276 0,8483 1,3862 2,0000 2,4897 3,0345 3,4483 3,7793 4,0621
Esfuerzo 11032 16547 23166 28958 41369 53779 66190 76670 89356 1,03E+05 1,14E+05 1,25E+05
Se determinó el interpolador de la forma Polinomial Nu 3 , que estableció el comportamiento del material bajo esfuerzos cortantes, y se obtuvieron las constantes del modelo hiperelástico, parámetros que fueron requeridos dentro de la formulación de la ley constitutiva aplicada en las simulaciones, ver fig. 3: Esfuerzo vs Deformación efecto cortante 1,4E+05
1,2E+05
1,0E+05
8,0E+04
σ
6,0E+04
4,0E+04
2,0E+04
0,0E+00 0,00
0,50
1,00
1,50
2,00
2,50
ε
3,00
3,50
4,00
4,50
-------y (Función tendencia) ------ vs σ (Función experimental)
Ecuación Calculada: y = 2315.9x3 - 14202x2 + 48328x + 8806.5
Figura 3. Gráficos de Esfuerzo vs. Deformación cortante y descripción obtenida por métodos matemáticos de regresión Polinomial. Se obtuvieron los datos para establecer el comportamiento del material a efecto volumétrico, ver tabla 4:
Tabla 4. Datos de Esfuerzo [Pa] vs. Deformación del efecto volumétrico. Deformación Esfuerzo 0,8847 0,9127 0,9412 0,9703
1,59E+06 1,21E+06 8,15E+05 4,14E+05
Se determinó el interpolador de la forma Polinomial Nu 3 , que estableció el comportamiento del material a efecto volumétrico, y se obtuvieron las constantes del modelo hiperelástico, parámetros que fueron requeridos dentro de la formulación de la ley constitutiva aplicada en las simulaciones, ver fig. 4: σ
Esfuerzo vs Deformación efecto volumetrico
1,80E+06 1,60E+06 1,40E+06 1,20E+06 1,00E+06
8,00E+05 6,00E+05 4,00E+05 2,00E+05
0,00E+00 0,88
0,89
0,90
0,91
0,92
0,93
0,94
0,95
0,96
0,97
0,98
ε -------y (Función tendencia) ------- vs σ (Función experimental)
Ecuación Calculada: y = 8E+07x3 - 2E+08x2 + 2E+08x - 5E+07
Figura 4. Gráficos de Esfuerzo vs. Deformación a efecto volumétrico y descripción obtenida por métodos matemáticos de regresión Polinomial. La descripción del modelo de la forma polinomial se sustenta a partir del primer y segundo invariante de la deformación. El modelo fenomenológico se encuentra formulado mediante la densidad de energía de deformación. N
W
I
c ij
1
3
I i
2
3
i j 1
N
j
k 1
Donde: o 2 c1 0 c 0 1
J el
1
2k
(2) (3)
Y: ko
2 d
(4)
Con: d
1 2 c10
c 01
(5)
La descripción del modelo de la forma polinomial se sustenta a partir del primer y segundo invariante de la deformación. El modelo constitutivo fenomenológico utilizado se encuentra formulado mediante la densidad de energía de deformación. Los elementos físicos que la componen son: W: Energía de deformación potencial, donde µ 0: es el módulo de bulbo inicial, k: es el módulo cortante inicial o bulbo inicial; Ī1: Primer invariante de deformación desviador; Ī2 : Segundo invariante de deformación desviador; J: Determinante gradiente de la deformación elástica ( Fij ), c10, c01: Caracterización de las constantes de material de los invariantes de deformación; d: Parámetro de incompresibilidad del material; Nu: Número de constantes del material determinado por la suma de los Cij. Los demás parámetros que conforman el modelo constitutivo se obtienen de las regresiones, por medio de ajuste de los datos experimentales de las funciones del comportamiento del material a tensión y compresión uniaxial, a tensión y compresión biaxial, el efecto de tensión y compresión plana, además de la prueba de compresión volumétrica [4]. Para un modelo hiperelástico polinomial de 3o se debe realizar el ajuste de los datos experimentales con una ecuación de 3o y determinar los coeficientes de polinomio que harán parte de los parámetros que conforman la ecuación del modelo hiperelástico [1-3].
d c 0 1 c1 0
2
..... c ij
(6)
n
RESULTADOS Las bases de datos experimentales fueron determinadas de la caracterización de elastómeros para representar fenómenos estructurales. Los resultados evidencian velocidades apropiadas en el desplazamiento del fluido dada la compresión de la estructura, con valores de 1.5 m/s y una compresión estructural que entrega como resultado una eyección del fluido de 46% de su contenido (fracción de eyección), Fig. 5 y 6. Vol. *100% [ml3]
t [s] Figura 5. Gráfica del volumen desplazado por la compresión estructural obteniendo una eyección del fluido de 46% de su contenido (46 ml).
Figura 6. Reducción volumétrica de la estructura y apertura de dispositivos complementarios con velocidades de 1.5 m/s. Finalmente, se logró que el material representara la dinámica de una estructura con un comportamiento funcional tanto en la compresión como en la recuperación de dicho modelo. CONCLUSION Se realizaron simulaciones computacionales que reflejaron funcionalidad mecánica de la estructura, con valores del desplazamiento del fluido de 1.5 m/s y una compresión estructural mayor al 40%. Se concluye que el modelo de comportamiento hiperelástico desarrollado con la ley fenomenológica polinomial 3o garantizó apropiada compresión para construir diversidad de estructuras que confinen fluidos sometidos a eyección dinámica. Se logró que el material interpretará: estabilidad y equilibrio estructural en la interacción F-E, capacidad volumétrica y condiciones de trabajo, condición de las fuerzas externas, estado de deformación de la estructura en la interacción, resistencia geométrica de la estructura, cantidad de eyección del fluido, funcionalidad durante y después de la compresión. REFERENCIAS [1] G. Kluth, B. Després, “Discretization of hyperelasticity on unstructured mesh with a cell centered Lagrangian scheme”, Laboratoire Jacques-Louis Lions, Université Pierre et Marie Curie, 2009, France. [2] G. Xiao-Yan and Moe Riyad, “On Stress Analysis For A Hyperelastic Material”, Sulzer Carbomedics Inc. Austin, Texas, 2000, EEUU. [3] M. Kaliske and H. Rothert, “On the finite element implementation of rubber-like materials at finite strains”, Engineering Computations, Vol. 14, No. 2, p216-232. 1997. [4] Z. Svecová, V. Cúth, J. Slabeycius, M. Kopecký, “Evaluating of Hyperelastic Material Behavior”, Department of material engineering, Púchov, 2005. Slovak Republic.
Análisis de esfuerzos y deformaciones de galerías en explotación subterránea de carbón Oswaldo Bustamante.1 Sebastián Jaramillo. 2* Gali Quitian.3
1: Doctor en Ingeniería, Universidad Nacional de Colombia. Instituto de Minerales Cimex. Medellín, Colombia. 2: Ingeniero de Minas y Metalurgia, Universidad Nacional de Colombia. Instituto de Minerales Cimex. Medellín, Colombia. 3: Ingeniero Mecánico, Universidad Nacional de Colombia. Instituto de Minerales Cimex. Medellín Colombia. *
[email protected] Área temática: Modelado matemático y computacional Resumen El proyecto carbonífero El Bloque está ubicado en la vereda de Jonás del municipio de Fredonia, departamento de Antioquia y en la actualidad se adelanta una explotación subterránea con beneficio de los mantos de carbón por medio del método de cámaras y pilares variación: ensanche de tambores. Un macizo rocoso de carbón (o manto(s)), es un tipo de estructura sedimentaria que al ser beneficiada por la minería, altera sus condiciones y/o características geomecánicas naturales (iniciales), por lo que se hace necesario y pertinente conocer sus propiedades mecánicas y así dimensionar el comportamiento y respuesta de excavaciones múltiples bajo un campo de esfuerzos dinámicos. La mecánica de rocas tradicional tiene muy en cuenta las características físicas y geológicas particulares de cada tipo de roca para así poder dimensionar su comportamiento bajo la acción de esfuerzos y con base a la estadística edificar sus modelos. Sin embargo, existe la posibilidad de estudiar los fenómenos de lanzamientos de esfuerzos y de campo de deformación asumiendo el macizo rocoso como un medio continúo, pudiendo así, por medio del uso de herramientas y conceptos reológicos y de los elementos finitos (MEF) aproximar su comportamiento mecánico. En el acercamiento inicial se considera, estudia y analiza el comportamiento individual de la galería principal de manto dos (2), sin afectaciones laterales ni lanzamiento de esfuerzos provenientes de excavaciones vecinas; solo considerando la afectación de la carga efectiva que soporta el techo de la excavación. En segunda instancia, se analiza un galería vecina, bajo los mismos procedimientos y consideraciones anteriores; con el objetivo final de comparar como es su comportamiento conjunto bajo un escenario de lanzamiento de presiones. Se reportan los resultados obtenidos a través de la simulación de las excavaciones por el método de los elementos finitos trabajado con un programa comercial. La simulación fue alimentada con las propiedades y caracteres mecánicos obtenidos del macizo a partir de los estudios correspondientes y ensayos físico/mecánicos realizados tanto al manto de carbón como a la roca suprayacente (techo). La adopción inicial es la de un estado plano de esfuerzos y deformaciones con el fin de expandirlo a posteriori a modelos triaxiales. Palabras claves: macizo rocoso, geomecánica, explotación subterránea, medio continuo, elementos finitos.
Introducción “Los fenómenos físicos de la naturaleza no pueden ser descritos directamente y todo intento de hacerlo engendra un modelo. Existen dos de estos modelos de gran utilidad, el modelo corpuscular y el modelo fenomenológico” (Concha & Barrientos, 1996) La compleja naturaleza de los macizos rocosos y sus múltiples respuestas bajo un campo de esfuerzos dinámico lleva también a una utilización de métodos analíticos y numéricos muy distinta de la que se hace en otras ramas de la ciencia y las estructuras continuas. En la mecánica de rocas el macizo rocoso no se puede elegir y esta premisa es independiente del uso que se le valla a dar (macizo rocoso como elemento estructural, macizo rocoso con potencialidades de ser explotado etc.), por lo cual, un modelo físico / matemático para un macizo rocoso deberá estar adaptado a las condiciones de incertidumbre y aleatoriedad de sus propiedades, imperfección y limitación de las mediciones in situ, variabilidad en el tiempo de las cargas propias del macizo y el papel de la tridimensionalidad de los diseños en los mismos. Al introducir el concepto de medio continuo a la mecánica de rocas y aún más importante, considerar un macizo rocoso como tal, se hace necesario remitirse a la definición clásica, la cual establece que un medio continuo puede ser aproximado como el medio que no pierde sus propiedades al ser infinitamente dividido. Dicha definición puede ser interpretada bajo la concepción de la relación entre el volumen del medio y una propiedad cualquiera del mismo (siendo dicha propiedad una función del volumen). Al analizar el límite cuando el volumen tiende a cero, es muy probable encontrar lugares del espacio sin contener materia. Esto se debe a la naturaleza estadística del movimiento de las moléculas en esta región que, por ello, recibe el nombre de dominio de los efectos moleculares. Pero a partir de un punto del volumen, la propiedad puede ser considerada continua en el espacio. Si la escala de los fenómenos que son de nuestro interés es estrictamente mayor que dicho punto (denominado límite del contorno) la mecánica del medio continuo da resultados que pueden ser usados con seguridad. De esta manera podemos extrapolar el concepto y considerar el macizo rocoso como el ensamblaje de una serie de elementos estructurales interconectados mediante un número finito de nodos (Monge & Ramirez, 2004). Las rocas situadas a una cierta profundidad están sujetas a esfuerzos, resultado del peso de los estratos suprayacentes, tensiones tectónicas residuales, etc. Cuando se realiza una excavación subterránea, el campo de esfuerzos es alterado localmente y se produce una redistribución de las tensiones originales que existen en el medio. Las tensiones que originalmente actuaban en la roca extraída (porción de macizo), se redistribuyen y deben ser soportadas por la roca que se encuentra en las proximidades de la excavación. Si la roca, sin exceder su resistencia, puede soportar indefinidamente esta carga, no se hace necesaria la colocación de un sostenimiento; sin embargo, en la mayoría de los casos, la excavación tiende a cerrarse y sin la colocación de un sostenimiento adecuado llegaría a colapsar. Cuando se considera la interacción macizo-fortificación se suele asumir para el macizo un comportamiento elasto-plástico continuo controlado por el sostenimiento, representado en el contorno interior de la excavación por una presión radial uniforme (Rodriguez Garrido, 2003). Al analizar las condiciones actuales de respuesta del macizo rocoso en el cual el proyecto minero el Bloque adelanta labores (especialmente las evidencias en la galería principal de manto 2) surge el interrogante de como estimar y/o proyectar el campo de deformación en función del campo de esfuerzo activo (lanzamiento de esfuerzos) y de la interacción
macizo-fortificación debido a la respuesta que están brindando el techo y los hastiales de la excavación. Al considerar las condiciones iniciales en el problema que se está analizando, se debe tener presente que el comportamiento del macizo se caracteriza principalmente por presentar una respuesta (deformación y/o fractura) no lineal, anisotropía y heterogeneidad. Para encontrar una respuesta al campo de desplazamientos y esfuerzos generados por las condiciones anteriores se debe utilizar métodos numéricos que garanticen soluciones convergentes a la real y para dicho estudio se utiliza el método de los elementos finitos (MEF). Este tipo de método proporciona resultados con poco esfuerzo, sin embargo, su punto más débil se encuentra en las suposiciones iniciales que se deben realizar entorno al problema, por ejemplo, para el presente estudio se utiliza la suposición de medio continuo. No obstante, es de considerar que cuanto más precisa sea la conceptualización del problema, más precisa será la solución respecto a la respuesta observada en campo. (Bobet, 2010) El Método de elementos finitos es de por sí, uno de los más utilizados al momento de analizar sistemas con medios continuos o cuasi-continuos, este método consiste en la discretización de una región en una cantidad finita de pequeños elementos que se unen a través de la intersección de nodos, el principio de análisis de este método está basado en desplazamientos virtuales, los cuales establecen que, para un cuerpo en equilibrio, el total de trabajo interno asociado con el campo de desplazamientos virtuales aplicados sobre el cuerpo deben ser igual al trabajo externo virtual total. Cuerpo del trabajo Al analizar las condiciones geológico – estructurales del yacimiento, podemos concluir que la estructura está compuesta por un sistema de bloques angulares altamente perturbados, producto de la interacción de un considerable número de discontinuidades propias del sistema de fallas de Cauca – Romeral (Cárdenas, 2009); dentro del cual se destaca la falla inversa que cruza el área de concesión con actitud estructural N35W / 18NE. Al estar inmersos en una estructura sedimentaria; es más que prioritario considerar las diferentes unidades litológicas que la componen (columna estratigráfica), tal y como se muestra a continuación: Tabla 1. Columna estratigráfica proyecto carbonífero el Bloque en orden de profundidad
Fuente: Elaboración propia a partir de (Ramirez Álvarez, 2009) Conociendo las condiciones físicas y estructurales del macizo y las unidades litológicas que intervienen en la interacción macizo-distribución de cargas-fortificación (Arcillolita y Manto 2) de la zona a estudiar, podemos realizar una aproximación a los principales parámetros del macizo (como los son su resistencia a la compresión, resistencia a la tracción y módulo de elasticidad) con el uso del criterio de Hoek & Brown generalizado. Paralelamente a esto; podemos estimar los ángulos de fricción y las resistencias cohesivas para cada unidad litológica ajustando una relación lineal media a la curva generada por el criterio generalizado de Hoek & Brown (Hoek, Carranza-Torres, & Kurcum, 2002):
Para el dominio estructural de la galería de transporte, que corresponde al techo de la excavación, con una predominancia litológica de la Arcillolita (i.e litología que soporta las cargas) el criterio de Hoek & Brown queda como se muestra a continuación:
Ilustración 1. Criterios de ruptura de Hoek & Brown y de Mohr & Coulomb ajustados según el criterio de Hoek & Brown generalizado para la Arcillolita.
Fuente: Elaboración propia a partir del programa computacional ROCKLAB 1.032
Caracterizado por los siguientes parámetros:
Tabla 2. Parámetros mecánicos para la Arcillolita según el criterio de Hoek & Brown Clasificación según Hoek & Brown
Criterio de Hoek & Brown
Resistencia a la compresión simple (Mpa)
35
mb
GSI
40
mi
4
Criterio de MohrCoulomb Ajustado
0.112
Cohesión (Mpa)
0.201
s
0.0001
Ángulo de fricción (°)
23.02
a
0.511
Parámetros de macizo rocoso Resistencia a la tensión (Mpa) Resistencia a la compresión simple (Mpa) Módulo de elasticidad (Mpa)
-0.035
0.335
446.46
Factor de 0.8 perturbación (D) Módulo de elasticidad 8750 (Mpa) Fuente: elaboración propia a partir del programa computacional ROCKLAB 1.032
Para el dominio estructural de los pilares de protección de la cámara en desarrollo que corresponde al carbón, el criterio de Hoek & Brown queda como se muestra a continuación:
Ilustración 2. Criterios de ruptura de Hoek & Brown y de Mohr & Coulomb ajustados según el criterio de Hoek & Brown generalizado para el carbón. Fuente: Elaboración propia a partir del programa computacional ROCKLAB 1.032.
Tabla 3. Parámetros mecánicos para el carbón según el criterio de Hoek & Brown. Clasificación según Hoek & Brown
Criterio de Hoek & Brown
Criterio de MohrCoulomb Ajustado
Resistencia a la compresión simple (Mpa)
39
mb
1.387
Cohesión (Mpa)
0.316
GSI
40
s
0.0008
Ángulo de fricción (°)
49.95
mi
15
a
0.511
Factor de perturbación (D) Módulo de elasticidad (Mpa)
0.2 10725
Parámetros de macizo rocoso Resistencia a la tensión (Mpa) Resistencia a la compresión simple (Mpa) Módulo de elasticidad (Mpa)
-0.022
1.011
1276.12
Fuente: elaboración propia a partir del programa computacional ROCKLAB 1.032 A partir de la caracterización física / mecánica de las unidades litológicas que participan en la interacción macizo – distribución de cargas – fortificación se entra a analizar la geometría propuesta por el proyecto minero El Bloque en manto dos (2). Debido a las dificultades de lanzamientos de esfuerzo y campo de deformaciones que se han evidenciado en este sector de la mina; el proyecto minero El Bloque se vio obligado a realizar un modificación al método; dejando pilares de material de interés (carbón) como sostén permanente con el fin de evitar mayores dificultades a profundidad.
Ilustración 4. Geometría del método de explotación propuesto por el proyecto minero en manto dos (2). Vista en planta. Fuente: Elaboración propia Partiendo de la geometría propuesta por el proyecto minero en el sector de manto dos (2) y del conocimiento de la columna estratigráfica (responsable de la carga efectiva que soporta la excavación); se lleva a cabo la simulación de la galería de transporte en el programa computacional ANSYS R14.5 bajo los parámetros mecánicos hallados para la Arcillolita (dominio estructural que soporta las cargas en la galería. Ver tabla 2) y los siguientes parámetros geométricos: ancho del túnel: 2m; altura del túnel 2.5 m.
Ilustración 5. Modelo de desplazamientos totales en el eje “y” en el contorno de la excavación (metros) Fuente: elaboración propia a partir del programa computacional ANSYS R14.5
La distribución del campo de desplazamientos es congruente con el sentido de actuación de la carga suprayacente a la bóveda del túnel; puesto que a medida que se redistribuye el esfuerzo en el contorno de la excavación, el desplazamiento en el eje “y” va disminuyendo hasta el punto tal de ser casi nulo en la base del túnel. Considerando ahora el comportamiento del campo de esfuerzos en torno a la excavación tenemos que:
Ilustración 6. Simulación de la distribución del esfuerzo efectivo sobre el contorno de la excavación. (Newtons / m2). Fuente: elaboración propia a partir del programa computacional ANSYS R14.5 Al analizar el comportamiento de los esfuerzos en la bóveda de la excavación, se obtiene un aproximado de 1.95 Mpa en el punto donde comienza la redistribución debido a la forma circular de la misma. Cuando se realiza el comparativo entre el esfuerzo actuante y la resistencia a la compresión uniaxial de esta unidad litológica (0.335 Mpa) claramente se evidencia una superioridad del esfuerzo proveniente de la capa suprayacente que soporta la unidad; lo que intuitivamente llevaría a aproximar un escenario de falla. Sin embargo, el esfuerzo de confinamiento al cual se encuentra la unidad litológica debido a la profundidad de la excavación (con un valor máximo de 1.3 Mpa), le otorga energía al sistema generando así un desplazamiento de la zona limítrofe estabilidad / falla propia de la Arcillolita (Ver ilustración 1) configurando un escenario favorable para la galería de transporte. Considerando ahora la interacción de la galería de transporte con uno de los tambores que configuran el método de explotación tenemos que el campo de desplazamientos totales queda como se muestra a continuación:
Ilustración 7. Modelo de desplazamientos totales en “y” para la configuración galería – tambor (metros). Fuente: elaboración propia a partir del programa computacional ANSYS R14.5 Desde la perspectiva comparativa con la ilustración 5; tenemos una analogía en el aumento de izquierda a derecha del campo de desplazamiento en la galería de transporte, pero la interacción con el tambor, me genera una pérdida de simetría en este campo (en el contorno de la galería), ocasionando así un efecto viga (pandeamineto hacia la mitad del tambor) y una simetría respecto al eje central vertical del tambor, i.e una disminución progresiva del centro hacía en los desplazamientos. Considerando ahora el comportamiento del campo de esfuerzos en torno a la configuración planteada tenemos que:
Ilustración 8. Simulación de la distribución de los esfuerzos en el eje “y”. (Newtons / m2). Fuente: elaboración propia a partir del programa computacional ANSYS R14.5 Al realizar el análisis conjunto entre la ilustración 8 y la ilustración 2 (Criterio de ruptura de Hoek & Brown y de Mohr & Coulomb ajustados según el criterio de Hoek & Brown generalizado para el carbón) podemos apreciar que el esfuerzo que me genera el efecto viga a lo largo del tambor (aproximadamente 6.8 Mpa) supera la resistencia a la compresión simple del domino litológico del tambor (i.e manto de carbón) lo que conlleva a configurar un escenario potencial de falla. Si se tiene en cuenta el esfuerzo de confinamiento al cual se encuentra el tambor (con un valor máximo de 0.74 Mpa) y el esfuerzo de confinamiento necesario para estar exactamente en la zona limítrofe estabilidad / falla propia del carbón (0.72 Mpa) vemos que el desplazamiento hacia la
zona de estabilidad es mínimo; presentando así una alta probabilidad de falla en la zona central del tambor. Conclusiones
A través de esta investigación, se encuentra que existe la posibilidad de analizar medios discontinuos y heterogéneos con el método de los elementos finitos; (ampliamente utilizado en la simulación de medios continuos) esto a través de la integración de un software como ANSYS 14.5R y un criterio empírico de caracterización de macizos rocosos realizado a través del software ROCKLAB 1.032 A partir de los resultados arrojados por la ilustración 6 y las evidencias de campo; podemos concluir que aunque el desplazamiento hacia la zona de estabilidad en la gráfica esfuerzo principal mayor VS esfuerzo principal menor tiene una magnitud considerable las evidencias de campo demuestra una tendencia marcada de un desplazamiento hacia la zona de inestabilidad Se hace necesario realizar mediciones in situ de los desplazamientos o deformaciones de la clave y los hastiales de las galerías con el fin de poder considerar modelos reológicos y de fluencia (sobre todo en la galería de transporte); esto en función del dominio litológico y sus propiedades físicomecánicas (Arcillolita) Se concluye que la interacción de galerías múltiples bajo escenarios de esfuerzos dinámicos, genera un cambio local en el campo de deformaciones y una redistribución del mismo hacia sectores que puedan sufrir fenómenos como el efecto viga
Agradecimientos Agradecimiento especial al instituto de minerales CIMEX de la Universidad Nacional de Colombia sede Medellín por brindar el acompañamiento, conocimientos teóricos y aportar el software licenciado para el desarrollo de este trabajo. Al proyecto carbonífero El Bloque, por tener siempre las puertas abiertas y la buena disposición de brindar la perspectiva desde las experiencias y dificultades que día a día tienen en la mina. Al laboratorio de pavimentos por brindar el espacio y el personal necesario para llevar a cabo la caracterización física de las muestras Bibliografía Bobet, A. (2010). NUMERICAL METHODS IN GEOMECHANICS, 35(1), 27–48. Cárdenas, H. (2009). Perfil del suroeste antioqueño. Retrieved July 26, 2014, from http://antioquia.gov.co/antioquiav1/organismos/planeacion/descargas/perfiles/Perfil_subregional_Suroeste.pdf Concha, F., & Barrientos, A. (1996). Mecánica Racional Moderna.pdf. Concepción. Hoek, E., Carranza-Torres, C., & Kurcum, B. (2002). Hoek & Brown failure criterion – 2002 Edition (pp. 267–273). Monge, L. A., & Ramirez, P. (2004). Mecánica de Rocas : Fundamentos de la Ingeniería de Taludes.
Ramirez Álvarez, Y. D. (2009). Diseño y evaluación de los ademes de madera, en la empresa C.I CARMINALES en el municipio de Fredonia. Universidad Nacional de Colombia - Sede Medellín. Rivera Pinzón, L. K. (2004). Mejoramiento del rendimiento del sistema de transporte en la mina el bloque de Fredonia Antioquia. Universidad Nacional de Colombia - Sede Medellín. Rodriguez Garrido, M. (2003). Evaluación del coeficiente de seguridad del sostenimiento de galerías y túneles en función de su rigidez y distancia al frente en diferentes macizos rocosos y caracterización mediante el método impacto-eco. Escuela técnica superior de ingenieros de minas.
ANÁLISIS DE FATIGA POR IMPACTO EN UNIONES ADHESIVAS Téllez Martínez, Juan Fernando*; Casas Rodríguez Juan Pablo+, Alvarado Prieto Peter• *Departamento de Ingeniería Mecánica, Universidad de los Andes, Bogotá, Colombia
[email protected] +Departamento de Ingeniería Mecánica, Universidad de los Andes, Bogotá, Colombia
[email protected] •Departamento de Ingeniería, Fuerza Aérea Colombiana, Bogotá, Colombia
[email protected]
1. ResumenEste trabajo busca analizar la respuesta mecánica de uniones adhesivas de traslape simple ante impactos repetitivos de baja velocidad y energía, generados por un dispositivo de proyección vertical de masas, diseñado en la Universidad de los Andes. Dicho análisis se realiza empleando materiales compuestos como sustrato, los cuales utilizan resina epóxica como matriz y fibra de vidrio como refuerzo, y adhesivos estructurales epóxicos para fabricar las uniones. Adicionalmente se analiza la influencia del espesor de la capa de adhesivo en la resistencia mecánica de dichas uniones. Los resultados serán analizados utilizando modelos de vida a fatiga y métodos de análisis óptico para estudiar la falla. 2. Palabras claveFatiga por impacto; Impacto de baja velocidad; Uniones adhesivas; Traslape simple; Energía de impacto; Espesor. 3. IntroducciónLas cargas dinámicas en una estructura generan esfuerzos repetidos que resultan en el deterioro de las propiedades mecánicas de los materiales y eventualmente en la falla de los mismos, dichos esfuerzos pueden tomar la forma de un patrón sinusoidal cuya amplitud es constante y relativamente baja. Parámetros como el esfuerzo máximo ( ), esfuerzo mínimo ( ), y el esfuerzo medio ( ) se suman a su descripción, este tipo de fatiga puede también llamarse fatiga estándar. Sin embargo, recientemente ha crecido el interés por analizar los efectos de impactos de baja velocidad en componentes y estructuras, los cuales pueden ser causados por cargas vibrantes o por otro tipo de eventos. Dicho fenómeno es denominado fatiga por impacto [1]. En la industria aeronáutica los vehículos se encuentran expuestos al fenómeno de fatiga por impacto durante toda su vida útil en diversas situaciones, algunos ejemplos de dichas situaciones son: La caída de herramientas durante el mantenimiento de las aeronaves, golpes de aves, partículas, granizo u otro tipo de escombros durante el vuelo, decolaje y aterrizaje, impactos repetitivos recibidos durante los aterrizajes, entre otras [2]. En la industria automotriz puede presentarse en situaciones donde pequeños escombros del camino son proyectados hacia diferentes componentes como al tanque de gasolina (o al tanque de hidrógeno en automóviles híbridos, particularmente) [3], Sin embargo estos son
solo algunos ejemplos, de muchas otras situaciones que pueden ser perjudiciales en la vida de operación de diferentes componentes. Uno de los principales objetivos en la industria del transporte siempre ha sido la optimización del peso y es por esto que las uniones adhesivas han tomado mayor importancia en los últimos tiempos. Gracias a su menor peso, alta resistencia, reducción de concentradores de esfuerzos y capacidad de unión de diferentes adherentes, en comparación con otro tipo de uniones no permanentes, como remaches y tornillos, su utilización en diferentes aplicaciones ha crecido rápidamente [4], sin embargo, se debe tener especial cuidado con su sensibilidad y deterioro ante diversos factores ambientales [1]. Debido a este crecimiento se requiere un mejor entendimiento de su respuesta mecánica ante esfuerzos de pelado, de corte y a la combinación de ellos (entre otros) cuando se presentan condiciones de fatiga como los mencionados anteriormente. Adicionalmente se tiene la necesidad de analizar la influencia (desempeño y costo) de la capa de espesor del adhesivo en la resistencia de las uniones bajo diferentes regímenes de carga. Algunos autores han analizado la influencia del espesor de la capa de adhesivo en la resistencia de diferentes tipos de uniones, por ejemplo Kinloch y Moore [5], hacen una recopilación de diferentes estudios bajo condiciones cuasiestáticas donde se estudia dicha influencia. En el primer estudio Kinloch y Shaw aplican principios de mecánica de la fractura para estudiar la tenacidad a la fractura de las uniones adhesivas utilizando configuraciones de doble viga en cantiléver, compuestas por acero de bajo carbono y adhesivos epóxicos endurecidos. La Figura 1 muestra algunos de sus resultados, donde se observa como las gráficas de tenacidad a la fractura vs espesor de la capa de adhesivo llegan a un máximo y luego decaen hasta un valor de tenacidad a la fractura que es independiente del espesor de la capa de adhesivo.
Figura 1. Tenacidad a la fractura vs espesor de capa de adhesivo [5].
La forma de la Figura 1 se puede explicar en términos de la zona plástica en la punta de la grieta. Las tres secciones de la gráfica de tenacidad a la fractura vs espesor de capa de adhesivo se pueden explicar de la siguiente forma:
i.
Para espesores de capa de adhesivo pequeños, la zona plástica es relativamente grande. Es por esto que no puede desarrollarse completamente dentro del adhesivo y por lo tanto la tenacidad a la fractura del mismo no alcanza su máximo.
Figura 2. Espesor de capa de adhesivo menor a la zona plástica [5].
ii.
Cuando la zona plástica en la punta de la grieta se puede desarrollar completamente de forma ajustada dentro de la capa de adhesivo, la tenacidad a la fractura del alcanzará su máximo. Por la presencia de los sustratos puede presentarse una distorsión del campo de esfuerzos en la punta de la grieta y achatar un poco la forma de la zona plástica convirtiéndola en una elipse. El espesor en el cual sucede todo esto es igual a 2 veces el radio de la zona plástica. Por esta razón el máximo de la tenacidad a la fractura obtenido excederá el valor te tenacidad a la fractura del adhesivo bajo condiciones de deformación plana.
Figura 3. Zona plástica igual al espesor de la capa de adhesivo [5].
iii.
Una vez el espesor de la capa de adhesivo resulta relativamente mayor que el tamaño de la zona plástica, la tenacidad a la fractura del adhesivo se vuelve independiente del mismo y aproximadamente igual a la tenacidad a la fractura del adhesivo bajo condiciones de deformación plana.
Figura 4. Espesor de capa de adhesivo mayor a la zona plástica [5].
Kinolch y Moore comparan el comportamiento de la tenacidad a la fractura vs el espesor de la capa de adhesivo en diferentes configuraciones de uniones adhesivas, concluyendo que dicho comportamiento es independiente de la geometría, debido a que el resultado es recurrente en todas las pruebas realizadas. Por otro lado existen algunos investigadores que han realizado estudios en el tema de fatiga, un ejemplo de ello es el trabajo realizado por Casas-Rodríguez et al [1] [4], donde se analiza el fenómeno de fatiga por impacto de baja velocidad en uniones de traslape simple, utilizando materiales como aluminio, fibra de carbono y adhesivos epóxicos.
Dichos estudios presentan gráficas F-N y E-N comparando fatiga de amplitud constante con fatiga por impacto, mostrando la gran influencia de este fenómeno en el deterioro de las propiedades mecánicas de los diferentes materiales al reportar una disminución del número de ciclos para falla de dos órdenes de magnitud en fatiga por impacto. Adicionalmente se presentan gráficas comparativas de la velocidad de crecimiento de grieta (siendo mucho mayores las velocidades de fatiga por impacto), demostrando lo crítico que puede llegar a ser este fenómeno. En cuanto a los modos de falla, la fatiga por impacto presenta fallas más evidentes y una superficie más irregular en comparación con la fatiga de amplitud constante. Finalmente otros autores [6] han integrado el estudio de la influencia de espesor de la capa de adhesivo en condiciones dinámicas de carga, utilizando uniones de traslape simple compuestas por sustratos de fibra de vidrio y adhesivos epóxicos. Dicha investigación reporta una disminución de dos órdenes de magnitud en los ciclos para falla al tener un aumento en la capa de espesor de adhesivo de 2,5 a 5 mm. Esto sigue demostrando la importancia de continuar con el estudio de la influencia del espesor de la capa de adhesivo bajo diferentes regímenes de carga. 4. MetodologíaLa metodología de la presente investigación está basada en la experimentación y se complementa con modelos analíticos que permiten analizar la vida a fatiga de las uniones adhesivas. Las uniones de traslape simple son fabricadas a partir de materiales compuestos (resina epóxica y fibra de vidrio) y adhesivos epóxicos, empleados en la fabricación de aeronaves de entrenamiento. Se analizarán uniones con 1, 5 y 10 mm de espesor en la capa de adhesivo con el fin de determinar la influencia de dicho espesor en la resistencia mecánica de la unión. Para cada uno de los espesores mencionados se realizan pruebas cuasiestáticas (Esfuerzo vs Deformación) y de fatiga por impacto (Energía de impacto vs Número de ciclos) con el fin de comparar la tendencia encontrada en cada régimen. 4.1. Materiales Fibra de vidrio pre-impregnada 7781-550 E-761LT. Fibras tejidas (0° - 90°).
Resina epóxica RhinoTM 1307 - LV
Endurecedor epóxico RhinoTM 3102 (30 min)
“Flocked cotton fiber” (Cotton flox). Relleno utilizado para el adhesivo (21%).
4.2.
Manufactura
Para las uniones de traslape simple se manufacturan paneles de 101.6 mm (Dirección 0° de las fibras) * 177.8 mm. A partir de los cuales se obtienen probetas de 177.8 mm * 25.4 mm con 25.4 mm de traslape.
El proceso de manufactura de las probetas está compuesto por los siguientes pasos: -
Medición y corte de la fibra de vidrio según las dimensiones de los paneles.
-
-
Formación de los laminados (16 láminas por panel). Proceso de curado: Alistamiento (Tela superficial, Bolsa de vacío, Cinta de sellado) y Curado (Presión de vació y curva de calentamiento en el horno). Unión de los paneles según los tres espesores establecidos (1, 5 y 10 mm). Se utiliza dispositivo para control del espesor. El filete utilizado no incluye el espesor del sustrato. Las uniones se curan a temperatura ambiente (20°C) durante 7 días. Corte de los paneles curados y traslapados, utilizando chorro de agua de alta presión para obtener las probetas finales.
Cada lámina de material tiene 0.009 pulgadas (0.2286 mm) de espesor. Se utilizan 16 láminas por panel que otorgan un espesor final que oscila entre 4 y 5 mm. La manufactura de las probetas se lleva a cabo en el laboratorio de materiales compuestos ubicado en el comando aéreo de mantenimiento (CAMAN) de la fuerza aérea colombiana a una temperatura de 22°C y 55% de humedad relativa. 4.3.
Ensayos mecánicos
4.3.1. Pruebas cuasiestáticas Las uniones adhesivas en configuración de traslape simple se prueban a tensión en la máquina INSTRON 3367 de la Universidad de los Andes con el fin de determinar el esfuerzo último de las mismas y así comparar el comportamiento del material sometido a cargas cuasiestáticas y dinámicas. Esta prueba se realiza a una velocidad de 12.7 mm/min según la norma ASTM D3807. 4.3.2. Pruebas de fatiga por impacto Se realizan pruebas de fatiga por impacto en las uniones adhesivas en configuración de traslape simple, utilizando el dispositivo de impacto por proyección vertical de masas (DIPVM) o “drop weight impact tester” de la Universidad de los Andes (Figura 5). Dicho dispositivo cuenta con sensores láser para medir el desplazamiento, sensores piezoeléctricos para medir la fuerza de impacto y una fotocelda que permite su acople con una cámara de alta velocidad. La máquina (DIPVM) cuenta con un martillo de 13.9 kg y una altura máxima de disparo es de 1.25 m, lo cual genera una energía de impacto de 170 J y una velocidad de 4.95 m/s aproximadamente [7]. Dichos parámetros pueden ser variados según la adición de diferentes masas y la masa del impactor. Las pruebas se realizan utilizando energías entre 0,6 y 3 J que corresponden a alturas de impacto entre 4,5 y 21 mm, y velocidades de impacto entre 0,3 y 0.6 m/s. La adquisición de datos para su posterior análisis se lleva a cabo con la ayuda del software LabView® sincronizado con máquina.
Figura 5. Dispositivo de impacto por proyección vertical de masa.
Resultados4.4. Pruebas cuasiestáticas La Figura 6 y la Figura 7 muestran los resultados obtenidos en las pruebas cuasiestáticas de las probetas de traslape simple. Es evidente que un aumento en el espesor de la capa de adhesivo genera una disminución considerable en la fuerza de ruptura de la unión. Según lo explicado por [5], el espesor de la capa de adhesivo de 1 mm supera dos veces el radio de la zona plástica en la punta de la grieta del adhesivo y por esto al aumentar dicho espesor solo se observa una disminución en la resistencia mecánica de la unión. Al comparar los resultados para los espesores de 1 y 5 mm, el esfuerzo cortante máximo (τ máximo) disminuye en un 47,16%, para este último. De forma similar, existe una disminución de 67,05% al comparar los espesores de 1 y 10 mm.
Figura 6. Comparación pruebas cuasiestáticas 1, 5 y 10 mm.
La Figura 7 muestra la disminución del esfuerzo cortante máximo (τ máximo) al aumentar el espesor de la capa de adhesivo teniendo en cuenta el error porcentual presente en las pruebas realizadas. La pendiente de dicha disminución es más pronunciada entre los espesores de 1 y 5 mm, mientras que para los espesores de 5 a 10 mm se atenúa un poco, incluso tendiendo a un posible valor asintótico que resulta ser el valor crítico de la tenacidad a la fractura del adhesivo bajo condiciones de deformación plana. Este comportamiento concuerda con lo reportado en [5], cuya explicación se basa en cirterios de mecánica de la fractura.
Figura 7. Comparación resultados pruebas cuasiestáticas 1, 5 y 10 mm.
4.5. Pruebas de fatiga por impacto La Figura 8 muestra la onda característica de un impacto característico para una energía de 0,62 J aproximadamente. El análisis de las pruebas de impacto se realiza graficando la energía de impacto, que es la variable controlable (eje vertical primario), y la fuerza máxima de impacto (representada por el pico de la onda), normalizada con respecto a la carga cuasiestática máxima (eje vertical secundario), versus el número de impactos para falla (fractura total de la unión adhesiva).
Figura 8. Onda de impacto.
La Figura 9 presenta los resultados obtenidos al someter a las uniones adhesivas en configuración de traslape simple a impactos repetitivos. Para un solo impacto se puede observar que la fuerza es mayor a la cuasiestática por el efecto visco elástico de los polímeros. Sin embargo, se puede observar falla en las uniones adhesivas sometidas a 10 impactos aproximadamente con el 60% de la carga cuasiestática, para los espesores de 5 y 10 mm, lo cual representa claramente el fenómeno de fatiga por impacto y una característica de daño acumulado en la unión adhesiva. Una vez detectado el fenómeno se realiza la comparación de los diagramas E-N bajo condiciones de fatiga por impacto para los tres espesores de capa de adhesivo. Como en los resultados observados en las pruebas cuasiestáticas, la variación entre los espesores de 1 y 5 mm es considerablemente mayor que la evidenciada entre los espesores de 5 y 10 mm, lo cual demuestra congruencia en el comportamiento de las uniones adhesivas sometidas a diferentes regímenes de carga. La constante de deterioro de la unión, que corresponde a la pendiente del diagrama E-N, aumenta con el espesor de la capa de adhesivo. Por el contrario la energía necesaria para la ruptura de la unión en un solo impacto disminuye (intercepto con el eje Y,
diagrama E-N). Ambas características describen una disminución en la resistencia de la unión bajo condiciones de fatiga. En general, para la mayoría de los niveles de energía de impacto las diferencias en el número de impactos para los diferentes espesores de capa de adhesivo están entre 1 y 2 órdenes de magnitud, siendo el espesor de 1 mm el de mejor desempeño. Este comportamiento concuerda con los resultados de las pruebas cuasiestáticas analizadas previamente y demuestra lo crítico que resulta el aumento del espesor de la capa de adhesivo después del valor óptimo, igual a dos veces el radio de la zona plástica en la punta de la grieta. 1,6
3,5 1 mm
1,4
5 mm
2,5
1,2
10 mm
2
1,0 0,8
1,5
0,6
1
F impacto / F qstat.
Energìa de impacto (J)
3
0,4
0,5
0,2
0 1
10
100
1.000
0,0 10.000
Número de ciclos (Nf) Figura 9. Comparación diagrama E-N 1, 5 y 10 mm.
4.6. Análisis óptico de falla La Figura 10 muestra las micrografías de la falla de las probetas de traslape simple sometidas a fatiga por impacto, obtenidas utilizando un microscopio electrónico de barrido JEOL-JSM 6490 LV. A partir de estas imágenes se puede observar una falla interlaminar por el sustrato, debido a que la sección mostrada en la figura A muestra la ausencia de las fibras en la matriz de resina epóxica y algunos rastros de fractura de las mismas, mientras que la figura B muestra las fibras pegadas al adhesivo, las cuales fueron arrancadas de la matriz, evidenciando la falla cohesiva por el sustrato mencionada anteriormente.
Figura 10. MEB falla uniones traslape simple. A) Sustrato-falla interlaminar. B) Adhesivo-fibras sustrato.
Se encuentra el fenómeno de fatiga por impacto como identificable y repetible según lo reportado en [1] y [4]. Además de esto el comportamiento de las uniones adhesivas es congruente en condiciones cuasiestáticas y dinámicas, teniendo una explicación basada en los conceptos de mecánica de la fractura para las primeras y generando las bases para poder generar hipótesis acerca de las últimas. 5. Conclusiones Los resultados de fatiga por impacto se aproximan con una función logarítmico-lineal, donde el eje X es el número de impactos para fractura y el eje Y corresponde a la energía de impacto, según el modelo de vida a fatiga de Johnson ( ). El diagrama E-Nf muestra un aumento considerable en la constante de deterioro (A) al aumentar el espesor de la capa de adhesivo de 1 a 5 mm. Sin embargo, el aumento de la constante es significativamente menor al variar dicho espesor de 5 a 10 mm. Aumentar el espesor de la capa de adhesivo en las uniones de traslape simple, después de cierto valor, resulta sumamente nocivo para su resistencia tanto en condiciones cuasiestáticas como de fatiga por impacto. En la presente investigación tanto en el régimen cuasi-estático como en fatiga por impacto se encuentra el mismo tipo de falla, la cual se observa en el sustrato de manera interfacial y demuestra alta resistencia en la unión adhesiva. Además de esto la tendencia de la reducción de la resistencia al aumentar el espesor de la capa de adhesivo es congruente en ambos regímenes. 6. AgradecimientosLos autores agradecen enormemente a la Universidad de los Andes quien brindó su apoyo en la realización de las pruebas de laboratorio, a la Fuerza Aérea Colombiana quien suministró los materiales y laboratorios para la fabricación de los especímenes utilizados en la presente investigación.
7. Bibliografía-
[1] J. P. Casas Rodriguez, I. A. Ashcroft y V. V. Silberschmidt, «Delamination in adhesively bonded CFRP joints: Standard fatigue, impact-fatigue and intermittent impact,» Composites Science and Technology, nº 68, pp. 2401-2409, 2008. [2] M. Meo, R. Vignjevic y G. Marengo, «The response of Honeycomb sandwich panels under low velocity impact loading,» International Journal of Mechanical Sciences, pp. 1301-1325, 2005. [3] K. Azouaoui, Z. Azari y G. Pluvinage, «Evaluation of impact fatigue damage in glass/epoxy composite laminate,» International Journal of Fatigue, nº 32, pp. 443-452, 2010. [4] J. P. Casas Rodríguez, I. A. Ashcroft y V. V. Silberschmidt, «Damage in adhesively bonded CFRP joints: Sinusoidal and impact-fatigue,» Composites Science and Technology, nº 68, pp. 2663-2670, 2008. [5] A. J. Kinloch y D. R. Moore, «The influence of adhesive bond line thickness on the toughness of adhesive joints,» The Application of Fracture Mechanics to Polymers, Adhesives and Composites, 2004. [6] I. Sridhar, J. H. Tang y N. Srikanth, «Static and fatigue failure analysis of adhesively bonded thick composite single lap joints,» Composites Science and Technology, 2013. [7] D. F. Avendaño, «Caracterización dinámica de estructura celular hexagonal,» Universidad de los Andes, Bogotá D.C., 2013. [8] I. Maekawa, «The influence of stress wave on the impact fracture strength of cracked member,» International Journal of Impact Engineering, nº 32, pp. 351-357, 2005. [9] W. M. Banks, D. H. Nash y O. S. David West, «An experimental study of damage accumulation in balanced CFRP laminates due to repeated impact,» Composite Structures, pp. 247-258, 2008. [10] B. Khan, R. M. Rao y N. Venkataraman, «Low velocity impact fatigue studies on glass epoxy composite laminates with varied material and test parameters-effect of incident energy and fibre volume fraction,» Journal of Reinforced Plastics and Composites, nº 14, pp. 1150-1159, 1995.
Redes Neuronales aplicadas al entrenamiento de Relaciones Constitutivas lineales y no lineales Mercado N., Fredy. Facultad de Ingenier´ıa, Universidad de Buenos Aires, Av. Paseo Col´on 850, CABA, Buenos Aires, Argentina. 30 de julio de 2014
Resumen
de unidades o neuronas. Cada neurona est´a vinculada a las neuronas de las capas siguientes por medio En este trabajo se emplean Redes Neuronales Ar- de una conexi´on (pesos sin´apticos), la cual se caractetificiales (RNAs) con entrenamiento supervisado pa- riza por modificarse progresivamente conforme se enra aprender las relaciones constitutivas de materiales se˜ nan patrones a la red. A este proceso le llamaremos el´ asticos lineales y no lineales is´ otropos con miras a aprendizaje. El objetivo propuesto en este trabajo es emplearlas como sustitutos de la relaciones consti- lograr, mediante un m´etodo iterativo (Backpropagatutivas matriciales cl´ asicas en c´ odigos de elementos tion), la modificaci´on de los pesos sin´apticos de toda finitos. Los datos de entrenamiento son generados de la red tras la ense˜ nanza repetida de patrones o estaforma sint´etica debido a la ausencia de datos de en- dos (vectores) de deformaciones y sus correspondiensayos experimentales. Las RNAs elegidas son el Per- tes tensiones. Lo anterior con el fin de reemplazar la ceptr´ on Simple y el Perceptr´ on Multicapa, mejorados RC cl´asica con una RC definida por la topolog´ıa de mediante el algoritmo de Propagaci´ on del error hacia la RN y el valor que tengan sus pesos sin´apticos, de atr´ as (Backpropagation). Los resultados muestran un tal forma que: aprendizaje satisfactorio del comportamiento deforσ = CRN ε (2) maci´ on-tensi´ on de los diferentes materiales estudiaDonde CRN es la RC dada por la Red Neuronal. Las dos. RN podr´ıan ser u ´tiles para aprender RCs de materiales con comportamientos altamente no lineales, don1. Introducci´ on de las tensiones (salidas de la red) est´ an en funci´ on de variables independientes (entradas de la red) que Las Relaciones Constitutivas (RC) relacionan tenintroducen no linealidades en el comportamiento del siones, energ´ıa libre y flujo de calor con las deformamaterial, como son la tasa de deformaci´on y el enciones del cont´ınuo y la temperatura [3]. En este tradurecimiento con la deformaci´ on para el caso de los bajo limitaremos el alcance a la relaci´ on entre tensiometales. Se propone entonces entrenar redes a partir nes y deformaciones dentro de un material is´otropo. de datos sint´eticos que consideraremos producto de Las RC pueden entenderse tambi´en desde el punto de ensayos experimentales, y una vez entrenada procevista funcional: deremos a comprobar la calidad de la relaci´ on consσ = Cε (1) titutiva encontrada y a tratar de validarla para un Donde σ es un vector que almacena las componen- conjunto de datos de entrada diferentes a los datos tes que definen el estado de tensiones, C es la RC de entrenamiento, es decir, desconocidos para la red. matricial cl´ asica de 6 x 6 componentes y ε el vector En la primera parte del art´ıculo se presentar´ a un que almacena las componentes que definen el estado resumen de la teor´ıa del Perceptr´on con Backpropade deformaciones. Los dos estados descritos y la RC gation. Luego, se tratar´ a la tem´atica de las RC y se est´ an definidos para un punto de un material. introducir´an el modelo el´astico lineal 3D y el modelo La Red Neuronal Artificial (RNA) que ser´a em- el´astico lineal 2D para deformaci´ on plana, para luego pleada es el Perceptr´ on, el cual est´ a formado por ca- mencionar detalles sobre el pre y post-procesamiento pas, y donde cada capa posee determinado n´ umero de los datos de entrenamiento, el criterio de conver1
2
3 RELACIONES CONSTITUTIVAS
gencia y la topolog´ıas de la red. Por u ´ltimo se mostrar´ an tres aplicaciones del Perceptr´ on con Backpropagation a dos problemas lineales y uno no lineal (RC elasto-pl´ astica).
g(h) = tanh βh
(7)
para rango − 1 ≤ g(h) ≤ 1 Su derivada es
2.
Perceptr´ on Multicapa con Backpropagation
El perceptr´ on es una Red Neuronal de aprendizaje supervisado que transforma una entrada en una salida. Una topolog´ıa ilustrativa de un perceptr´on multicapa puede ser observada en la Figura 1. El funciona-
g ′ (h) = β(1 − g 2 ) = β(1 − (tanh βh)2 ) El par´ametro β a menudo es 1 [8].
2.2.
Backpropagation
Consiste en corregir los pesos wij para cada capa a partir del error entre la salida deseada y la salida actual de la red para un patr´on determinado. Los pasos para realizar las correcciones son bastante conocidos y se pueden hallar en la referencia [8]. El algoritmo de Backpropagation se aplica para cada patr´on presentado a la red hasta presentar todos los patrones. A cada presentaci´on de todos los patrones se le llama un “epoch”. El aprendizaje se extender´a por tantos “epoch” como sea necesario hasta que se cumpla el criterio de convergencia o bien hasta un n´ umero de “epochs” determinado por el usuario de la red. Para todos los pesos de la red el algoritmo es:
Figura 1: Red Neuronal Feed-Forward de 2 capas. Topara EPOCH=1 hasta MAXEPOCH mada de [8]. para PATRON=1 hasta NPATRONES ∆Wij = ∆Wij + ηδi Vj miento del Perceptr´ on, capa por capa, es como sigue: fin dado un patr´ on µ, la unidad escondida j recibe una Wij = Wij + ∆Wij entrada neta ξkµ fin X hµj = wjk ξkµ (3) k
3.
Relaciones constitutivas
la cual produce la salida Vjµ = g(hµj ) La unidad de salida i recibe la entrada Vjµ X Wij Vjµ hµi =
(4)
(5)
3.1.
Modelo general - 3D
El modelo cl´asico de la relaci´ on constitutiva lineal entre tensiones y deformaciones para el caso m´ as general (3D) est´a dado por la ecuaci´on 8:
j
que produce la salida final Oiµ = g(hµi )
2.1.
(6)
Funci´ on de activaci´ on
Para la funci´ on de activaci´ on g(h) normalmente se utiliza una funci´ on sigmoidea. Esta funci´ on debe ser diferenciable y lo normal es desear que se sature en ambos extremos [8]. La funci´ on empleada para el caso de estudio no lineal es:
1 σ ν 11 1 − ν σ22 ν 1 − ν σ33 = A 0 τ12 τ23 0 τ13 0
ν 1−ν 1 ν 1−ν 0
ν 1−ν ν 1−ν 1 0
0
0
0
0
0
0
0
0
0 1 − 2ν 2(1 − ν)
0
0
0
0 1 − 2ν 2(1 − ν) 0
0
ε 11 0 ε22 0 ε 33 γ 0 12 γ 0 23 1 − 2ν γ13 2(1 − ν)
(8)
A=
E(1 − ν) (1 + ν)(1 − 2ν)
´ 4 ENTRENAMIENTO DEL PERCEPTRON
3
donde σ y τ son tensiones normales y de corte, y ε y γ son deformaciones normales y de corte, respectivamente.
3.2.
Modelo para deformaci´ on plana 2D
Este modelo puede obtenerse a partir de la RC tridimensional de la ecuaci´ on 8. En este caso las deformaciones est´ an sobre un plano 2D, m´ as no las tensiones. Esto indica que: ε33 = ε13 = ε23 = 0
(9)
σ13 = σ23 = 0
(10)
1. Se selecciona la topolog´ıa de la RNA que pueda representar el problema. El art´ıculo de la Referencia [9] limita esta opci´on s´olo a redes de m´ as de dos capas, sin embargo el art´ıculo de la Referencia [5] logra emplear una red sin capas escondidas para modelar una relaci´ on constitutiva el´astica lineal. 2. Adquisici´on de datos que contengan la informaci´on que caracteriza el problema muestreando el comportamiento de un material mediante modelaci´on num´erica o experimentaci´on real. 3. Inclusi´on de la RNA como una subrutina en un programa de computadora.
Simplificando a partir de la Ecuaci´ on 8 y teniendo Est´a comprobado que las RNA son aproximadores en cuenta que para el caso de deformaci´ on plana la universales para cualquier funci´on con muchas variacomponente ε33 = 0 obtenemos (con γ12 = 2ε12 ): bles independientes, donde la informaci´ on sobre la dependencia funcional es transferida a la RNA por medio de un n´ umero de datos discretos que describan 1−ν ν 0 σ11 ǫ el valor de la funci´ on para un n´ umero suficientemente 11 ν σ22 1−ν 0 largo de puntos de muestra en el espacio de las va ǫ (11) =M riables de inter´es. El n´ umero de capas escondidas, el (1−2ν) 22 0 τ12 0 2 n´ umero de neuronas y la repartici´on de pesos entre 2ε12 capas escondidas deben ser definidas por el usuario. σ33 ν ν 0 Estas caracter´ısticas son dependientes del problema E y no es posible dar una prescripci´on general, por lo M= (12) cual la pr´actica com´ un es determinar la mejor confi(1 + v)(1 − 2ν) guraci´on topol´ogica tratando de minimizar el n´ umero de capas ocultas y el n´ umero de neuronas por cada 4. Entrenamiento del Per- capa [1].
ceptr´ on
Para entrenar la red se modifican los pesos de las interconexiones mediante un m´etodo iterativo para obtener una se˜ nal de salida lo m´ as pr´ oxima posible a la salida deseada como respuesta a un patr´on de entrada. A este proceso se le llama entrenamiento. Para ello deben conocerse de antemano datos de entrada y salida de la red. Una parte de estos datos puede ser utilizada para entrenarla, mientras que la otra parte puede usarse para realizar pruebas que confirmen un correcto aprendizaje [9].
4.1.
Procedimiento
4.2.
Normalizaci´ on de datos
Antes de entrenar la red tanto las variables de entrada como las de salida deber´an ser normalizadas en cierto rango. Dependiendo de la funci´on de activaci´on que utilicemos podemos necesitar normalizar en un rango de 0 a 1, tal como recomiendan las referencias [14] y [10]. Esto nos llevar´ıa a usar un m´etodo de normalizaci´on ampliamente utilizado, que es X′ =
X − 0,95Xmin 1,05Xmax − 0,95Xmin
(13)
donde X son los datos originales, y Xmin y Xmax son los valores m´ınimos y m´aximos de X, respectivaEl grupo de datos de entrada y salida puede obmente. X ′ son los datos unificados de los correspontenerse a partir de evaluaciones num´ericas o ensayos dientes X. Este procedimiento de pre-procesamiento experimentales. Una vez finalizado el entrenamiento, puede hacer m´as eficiente el entrenamiento de la red la capacidad de generalizaci´ on de la red nos permineuronal [13]. te evaluar tensiones a partir de deformaciones desconocidas para la red. El procedimiento para modelar Otra forma de normalizar los datos de entrada es la relaciones constitutivas se divide en tres pasos impor- propuesta por Hashash, Jung y Ghaboussi en la Refetantes: rencia [7]. Se emplea cuando los datos corresponden a
´ 5 EJEMPLOS NUMERICOS
4
deformaciones que ser´ an empleadas en la capa de entrada de la red neuronal. El vector de deformaciones ser´ıa escalado entonces por un factor apropiado: N εN i
εi = ε Si
tal que
−1
1 mm hasta la rueda giratoria. En el PFC, en cambio, el Gap es inferior a 1 mm. En este último proceso (utilizado a escala industrial) la zona de metal líquido está limitada entre la boquilla y la rueda, esto permite la producción de cintas con anchos del orden de 150 a 200 mm [3]. Se ha demostrado que muchas de las propiedades físicas de las cintas dependen de los valores de los parámetros del proceso [4; 5; 6] tales como la velocidad de la rueda (VR) [7], presión de eyección de gas (PE) [6] , temperatura de fusión (Tf) [8; 9], Gap [10]. A su vez, también se ha establecido que las variaciones de estos parámetros cambian las dimensiones geométricas y sobre todo el espesor de las cintas. Aumentando la velocidad tangencial de la rueda, se incrementa la tasa de enfriamiento del material, lo que influye también en el espesor de la cinta obtenida [6]. Cabe destacar que, para alcanzar las mejores propiedades magnéticas, es necesario controlar el espesor de las cintas obtenidas por melt spinning [11].
2. Metodología: En este trabajo se estudia la técnica de CBMS, explicada en párrafos anteriores. En la Figura 1 se presenta un esquema del equipamiento utilizado para la obtención de cintas, generalmente, amorfas debido a la tasa de enfriamiento involucrada (~106 K/s) [12]. La bobina de inducción del horno calienta la aleación sobre el punto de fusión, y mediante la aplicación de argón a sobrepresión se expulsa el material fundido a través de la boquilla, con una velocidad de eyección que permite la formación de la cinta sobre la rueda giratoria. Se utilizó una aleación madre de Fe78B13Si9 (%at), obtenida a partir de ferroaleaciones comerciales. Dicha aleación se logra luego de varias fundiciones en atmósfera de vacío. La aleación madre se tritura en piezas de aproximadamente 5 mm, y posteriormente es refundida a la temperatura de fusión del metal 1300 °C en un crisol de cuarzo que posee la boquilla de Nitruro de Boro en su extremo. Procedimientos explicados en: Obtención del lingote de aleación madre, Obtención de las cintas [13]. Para nuestro ensayo los set ups del proceso se muestran en la Tabla I. Como producto se obtuvo una cinta que tiene un espesor medio de 32 µm y un ancho promedio (A) de 1,18 mm, Figura 2. Como se mencionó en párrafos precedentes, la producción de materiales por esta técnica tiene gran cantidad de variables de suma importancia y, la elevada velocidad del proceso de eyección es la principal limitación para poder ajustar dichos parámetros y proyectarlos a escala industrial. Por lo expuesto, el proceso productivo se filmó con una cámara de alta velocidad para tener registros reales de la experiencia, utilizando un equipo de filmación Visionresearch, Phantom-HD con 4700 pps y una resolución 320 x 640. Luego de la captura del procedimiento, el material fílmico se editó con el software Phantom Camera Control [14], que permitió procesar todo el ensayo cuadro por cuadro, logrando identificar cada paso de forma sistemática. El período considerado para este trabajo tiene que ver con el tiempo transcurrido desde el inicio de la eyección t0 hasta la formación inicial de la cinta amorfa t8. La subdivisión en ocho intervalos se debe a la captación de los cuadros y responde únicamente al funcionamiento del equipo de filmación. En la Tabla II se muestra cómo se dividieron los intervalos de la eyección con el correspondiente incremento relativo y absoluto de tiempo referido a t0. Velocidad tangencial de la rueda VR [m/s] Gap [mm] Presión de eyección [bar] Tabla I: Set Ups del proceso
Tiempo t relativo [ms] t absoluto [ms] t0
0,000
0,000
t1
2,127
2,127
t2
0,213
2,340
t3
0,213
2,553
t4
0,425
2,978
t5
0,426
3,404
t6
1,064
4,468
t7
1,702
6,170
t8
4,467
10,637
Tabla II: Intervalos de tiempo
40 3 0,3
Para el modelo matemático se utilizó el programa OpenFOAM®, y se replicó la geometría comprendida por el crisol, metal fundido, boquilla y superficie de la rueda. Se realizó un mallado de 259600 celdas y se establecieron: la presión de eyección del gas argón y la velocidad de giro de la rueda como condiciones de contorno. Se resolvió con el solver interFOAM que, basado en el método de volumen de fluido (VOF), permitió trabajar con dos fluidos inmiscibles (aleación metálica y gas argón) definiendo la fracción de volumen de cada fase involucrada en las celdas [15]. Luego de la confección del “meshing” y el procesamiento con interFOAM, la simulación se visualizó con el software ParaView [16].
Figura 1: Esquema de equipo utilizado
Figura 2: Cinta obtenida a partir de aleación madre de Fe78 Si9 B13
3. Resultados: Para poder contrastar las muestras fotográficas con los perfiles obtenidos por OpenFOAM®, se tomaron las escenas con el software ParaView para los mismos intervalos propuestos para las muestras anteriores (t0…t8). En las siguientes figuras se observan, a modo comparativo, las muestras fotográficas (a) obtenidas durante la producción de la cinta amorfa y (b) obtenidas con ParaView.
(a)
(b)
Figura 3: t0=0. Comienzo de eyección, este es el tiempo inicial considerado para contrastar con la simulación. En (a) se puede observar el destello característico que se produce cuando el material alcanza el punto óptimo para comenzar la eyección. En este instante todo el material se encuentra en fase líquida, tal cual se observa en (b).
(a)
(b)
Figura 4: t1= 2,127 ms. En este momento se puede observar la salida desde la boquilla. Todo el material permanece en estado líquido.
(b)
(a)
Figura 5: t2= 2,34 ms. Conforme el proceso sigue avanzando, se ve el ensanchamiento de la columna líquida en la parte inferior a la misma.
(a)
(b)
Figura 6: t3= 2,553 ms. Persiste el ensanchamiento inferior, este es el instante previo a que el jet de metal líquido haga contacto con la rueda giratoria.
(a)
(b)
Figura 7: t4= 2,978 ms. En este instante se pone en contacto el metal líquido con la rueda giratoria. El detalle de (b) muestra la formación de un plano en la parte inferior del jet líquido. Es el comienzo de la transferencia de calor a la taza de 106 Ks-1.
(DSM) (USM)
(b)
(a)
Figura 8: t5= 3,404 ms. Comienza a formarse el Downstream Meniscus (DSM) y el Upstream Meniscus (USM).
(DSM) (USM)
(DSM)
(USM)
Puddle
Tri‐junction
(a)
Puddle
Tri‐junction
(b)
Figura 9: t6= 4,468 ms. En este instante ya se observa el perfil característico de esta técnica. Se consolidaron los DSM y USM, el Puddle y el Tri-junction (frontera límite donde comienza a formarse la cinta) [17].
(a)
Solidificación
Solidificación
(b)
Figura 10: t7= 6,170 ms. En este instante se observa el inicio de la solidificación, es el comienzo de la formación de la cinta.
Cinta
(a)
Cinta
(b)
Figura 11: t8= 10,637 ms. En este instante el proceso alcanza continuidad.
4. Conclusiones: Tanto en las muestras fotográficas como en las simulaciones se pueden observar los componentes más relevantes en la conformación de las cintas, a saber: Downstream Meniscus (DSM), Upstream Meniscus (USM), Puddle y el Tri-junction (frontera límite donde comienza a formarse la cinta). Estas zonas dan sostén a la columna ferrostática de metal fundido, además de ser las fronteras de intercambio térmico. La similitud de los perfiles de conformación en los intervalos comparados nos permite presentar los primeros resultados alentadores de esta metodología y suponer que las condiciones de contorno planteadas en el modelo matemático son confiables, con la posibilidad de escalar este análisis a procesos semi – industriales. 5. Referencias: [1] http://www.openfoam.com/ [2] Paul H. Steen, Christian Karcher - FLUID MECHANICS OF SPIN CASTING OF METALS School of Chemical Engineering, Cornell University, Ithaca, New York 14853-5204 USA, Institute for Fluid Mechanics, Dresden University of Technology, 01062 Dresden, Germany. [3] Kurokawa et al. - US Patent No: 5,908,068 (1999) [4] M. Srinivas, B. Majumdar, D. Akhtar, A. P. Srivastava, D. Srivastava. Influence of wheel speed during planar flow melt spinning on the microstructure and soft magnetic properties of Fe 68.5 Si18.5 B9Nb3Cu1 ribbons. [5] Tetsuji Saito. Electrical resistivity and magnetic properties of Nd-Fe-B alloys produced by meltspinning technique. Jun 19, 2010. [6] G. Pozo Lopez, L.M. Fabietti, A.M. Condo, S.E. Urreta. Microstructure and soft magnetic properties of Finemet-type ribbons obtained by twin-roller melt-spinning. Jun 1, 2010.
[7] Victor I. Tkatch *, Alexander I. Limanovskii, Sergey N. Denisenko, Sergey G. Rassolov. The effect of the melt-spinning processing parameters on the rate of cooling. Physics and Engineering Institute of National Academy of Sciences of Ukraine, Rose Luxemburg street, 72, Donetsk 83114, Ukraine Received 5 June 2000. [8] K. H. Cho and M. C. Kim, Numerical Modeling of Planar Flow Casting Process (2006) [9] Keiyu Nakagawa, Teruto Kanadani, Yasuyuki Mori*2 and Yuto Ishii*3. The Effect of Jetting Temperature on the Fabrication of Rapidly Solidified Fe-Si-B Systems Alloys Using SingleRoller Melt Spinning*1. Faculty of Engineering, Okayama University of Science, Okayama 7000005, Japan. [10] Chunbai Wang, Numerical modeling of free surface and rapid solidification for simulation and analysis of melt spinning. [11] Zhaoyang Wu, Xi’an Fan, Guangqiang Li, Zhanghua Gan, Jian Wang, Zhan Zhang, "Evolution from amorphous to nanocrystalline and corresponding magnetic properties of Fe-Si-B-Cu-Nb alloys by melt spinning and spark plasma sintering", Materials Science and Engineering B 187 (2014), p. 61–66.[12] Heping Liu, Numerical Simulation of Initial Development of Fluid Flow and , Wenzhi Chen, Shengtao Qiu, Guodong LiuHeat Transfer in Planar Flow Casting Process [13] M., Pagnola, M., Barone M.,”Estudio Magnético de Cintas de FeSiB Obtenidas mediante Melt Spinning”, SAM-CONAMET, 13er Congreso internacional en Ciencia y Tecnología de Metalurgia y Materiales, Puerto Iguazú, 2013, p. 143. [14] http://www.visionresearch.com/Products/Accessories--Options/PCC-Software/ [15] Márquez Damián S., "Description and utilization of interFoam multiphase solver", FinalWork, Computational Fluid Dynamics. [16] http://www.paraview.org/ [17] E.A. Theisen( a), M.J.Davis (b), S.J.Weinstein (c), P.H.Steen (d). Transient behavior of the planar-flow melt spinning process. (a) Metglas Inc. (b) Engineering Sciencesand Applied Mathematics,NorthwesternUniversity. (c) Department of Chemical and Biomedical Engineering,Rochester Institute of Technology. (d) School of Chemical and Biomolecular Engineering.
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS
SIMULACION POR MEDIO DE MEDIO DE ELEMENTOS FINITOS DEL COMPORTAMIENTO DE HONEYCOMB DE ALUMINIO SOMETIDOS A COMPRESION CUASI-ESTATICA Orrego Caicedo, Camilo Jose. Casas Rodriguez, Juan Pablo
[email protected],
[email protected] Departamento de ingeniería mecánica Universidad de los Andes
Resumen—La investigación de cómo afectan diversos parámetros en la respuesta de los honeycombs a bajas tasas de deformación unitaria compresivas se hacen mediante el uso de programas de elementos finitos. Estos se validan por medio de experimentaciones realizadas anteriormente y usando las mismas condiciones de material y geometría. A partir de resultados comparados, se varían parámetros geométricos como espesor de la lámina y la altura de la celda unitaria para ver la importancia de estos. Además, también se estudiara como afectara el uso de adhesivos en la respuesta tanto en la región elástica como la plástica compresiva. Los resultados obtenidos para las simulaciones de comprobación validan las pruebas realizadas al ser similares a los datos experimentales. También se obtiene una tendencia positiva entre el modulo de Young y el esfuerzo compresivo plástico con los parámetros antes mencionado, lo cual tiene coherencia con modelos desarrollados anteriormente. Palabras claves—Honeycomb, Simulación, Esfuerzo Plástico Compresivo, Modulo Elástico, Esfuerzo Elástico Máximo, Espesor, Altura, Tamaño de Malla, Adhesivo. Introducción—Los honeycombs, desde el punto de vista de manufactura, son un tipo de configuración en la cual se crean láminas de un mismo material y se unen de tal forma que
parece un panal de abeja, de ahí su nombre. También se puede ver como un material celular bidimensional que es relativamente fuerte y rígido en el plano normal a la microestructura pero débil en plano [1]. Estos se pueden utilizar en una estructura tipo sándwich, lo cual se refiere a que se encuentran en medio de 2 láminas, las cuales mantienes unido al honeycomb y ayuda al proceso de dispersión de energía [2]. Estos sistemas tienen alta resistencia al pandeo, buena rigidez relativa al peso y estabilidad, lo cual lo hace importante en muchas industrias, donde la durabilidad estructural y la tolerancia al daño son vitales en el momento de diseño [3]. Uno de los principales problema al estudiar daños en honeycombs en estructura tipo sándwich es que presenta tipos de fallas dinámicos altamente complejos y la imposibilidad de analizar el comportamiento del material en tiempo real [3]. Por eso, se estudian a partir de resultados experimentales y se crean modelos matemáticos que se aproximen en lo posible a estos. Con el desarrollo de la tecnología, se visualizó el uso de los nuevos softwares de elementos finitos como herramientas útiles en el estudio de los honeycombs. La principal ventaja de este método es la relación costosresultados, ya que existen una gran variedad
El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS de simuladores comerciales que pueden obtener resultados a distintas condiciones cercanos a los obtenidos en un laboratorio de pruebas.
Figura 1. Honeycomb. (Zhang, Zhang, Wen. 2014) Resultados de simulaciones realizadas anteriormente han dado resultados importantes, como por ejemplo Chen y Ozaki [4] mostraron la importancia de la altura de los honeycombs en el modulo elástico del elemento, y Yamashita y Gotoh analizaron el efecto de la forma de la celda y su espesor en las propiedades del honeycomb [5]. Aun con estos avances, se podrían hacer más pruebas para seguir estudiando el efecto de diversos parámetros en las propiedades de los honeycombs, a partir de los resultados de las simulaciones. El propósito de este proyecto es usar un programa de modelación con el cual se pueda analizar el efecto que tiene el cambio de propiedades geométricas del honeycomb, como lo puede ser su altura o el ancho de la lámina con la cual está hecho, en su comportamiento en el diagrama esfuerzodeformación. Para eso, se realizaran simulaciones variando estos parámetros y se analizaran los resultados reportados por el simulador. También, se verá los efectos del
uso de un adhesivo epoxi, que es la forma usada en la manufactura, en la simulación. Para tener un mayor criterio, otro parámetro a analizar será el tamaño de la malla con la cual se hará el análisis de elementos finitos. Metodología—El programa de elementos finitos que fue usado para realizar las simulaciones fue Marc Mentat versión estudiantil 2013 y profesional 2012. El elemento a analizar fue la celda unitaria del honeycomb, que se modeló como 3 láminas que entre ellas forman un ángulo de 120°, de las cuales 2 tienen espesor normal y la otra tiene espesor doble. Para el análisis del comportamiento del honeycomb con el adhesivo, se agrega otro elemento, el cual tendrá la misma malla que las laminas, en el centro de lo que sería la lamina de espesor doble, volviéndolo 2 laminas separadas. La longitud de cada lámina es de 1,74 mm. El espesor y la altura son variables para el análisis.
Figura 2. Imagen de la celda unitaria con espesor 0.025 mm.
El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS El enmallado se hace manualmente, haciéndolo más fino en la zona donde la lamina se dobla. Esto se hace para evitar distorsiones que afecten los resultados. El tamaño de la malla en los ejes X y Y se mantiene constante en todas los modelos y es igual a 1/3 del espesor de los elementos y del largo de las laminas. La subdivisión realizada para el eje Z es variable, para el estudio de su efecto. El límite de este está condicionado a la convergencia de los resultados. El material usado para la simulación es el aluminio 5052, para así poder ser comparado con experimentos anteriores. Para las simulaciones con adhesivos, se incluye un material adicional con las propiedades de un epoxi. Además, se define que los elementos del enmallado se modelan como sólidos. Propiedad del Valor Aluminio Módulo de Young 70.3 (GPa) Coeficiente de 0.33 Poisson Esfuerzo de fluencia 292 (MPa) Esfuerzo Ultimo 313 (MPa) Deformación de fallo 12 (%) Tabla 1. Propiedades del aluminio usado en la simulación [6] Propiedad del epoxi Valor Módulo de Young 2.01 (GPa) Esfuerzo de fluencia 29 (MPa) Esfuerzo Ultimo 57.7 (MPa) Deformación de fallo 10.4 (%) Tabla 2. Propiedades del epoxi usado en la simulación [7].
Adicionalmente, se agregan las condiciones de fronteras necesarias. Las restricciones de desplazamiento usadas son de 2 tipos: es de los ejes X y Y y se aplica en un único nodo, que se encuentra en la parte inferior del honeycomb, en el lugar donde se intercepta las laminas; el otro es del eje Z y se ubica en todos los nodos inferiores del honeycomb. Se simula la compresión cuasi-estática mediante el uso de un nodo exterior a la geometría analizada. A este, se le asigna que se traslade únicamente en el eje Z por medio de restricciones en los demás ejes. Posteriormente, se une a los nodos superiores del honeycomb mediante enlaces, lo cual hace que todos presenten el mismo desplazamiento en el sistema fuera del plano. En el postprocesamiento, a este único nodo se examina la fuerza y el desplazamiento que reporta, lo que facilita el análisis. Se convierten los datos en esfuerzo y deformación por medio de ecuaciones. Véase ecuaciones 1 y 2.
Donde en la formula (1), el área se define a partir de un rombo que recubra todo el elemento y su área a partir de teoría desarrollada [8]. Véase ecuación 3.
Las parámetros utilizados para las pruebas fueron espesores de 0.025, 0.05, 0.075, 0.1, 0.25 y 0.5 mm para una misma altura de 12 mm, y alturas de 5, 10, 12, 15, 20 y 25 mm para un mismo espesor de 0.025 mm. La combinación de prueba es de un espesor de 0.025 mm y una altura de 12 mm.
El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS Se usó el método arc-lenght del software Marc, el cual indica que las posibles inestabilidades en el sistema son producto de las cargas mecánicas y no las reporta como error, para realizar las simulaciones. El incremento de la fuerza se define por el método Crisfield, que lo hace por medio de una función cuadrática, la cual es resuelta usando la iteración de Newton, y a partir de estos, buscar un punto de equilibrio en el sistema.
Figura 3. Muestra del sistema usado para generar la deformación en el modelo. Resultados y Análisis 3,5 Esfuerzo ingenieril (MPa)
El método arc-lenght usado para esta prueba fue definida con un número máximo de incrementos de 200, un numero deseado de reciclajes sobre incrementos de 5 y un máximo radio de arc-lenght sobre el inicial de 10. Esto se definió a partir de las pruebas realizadas en los modelos con un espesor de 0.025 mm, una altura de 12 mm y adhesivo, y comparándolas con los resultados experimentales en un honeycomb con las mismas características. Pero, debido a que este método es de incremento automático, se dificulta la obtención de muchos resultados [9]. El módulo de Young se calculó dividiendo el esfuerzo elástico máximo entre la deformación a la que ocurre y el esfuerzo plástico compresivo se calculó promediando los valores en la meseta presentada en las gráficas.
3 Con Adhesivo (1/25)
2,5 2
Con Adhesivo (1/50)
1,5
Con Adhesivo (1/75)
1 0,5
Experimental
0 0
0,02
0,04
0,06
Deformacion unitaria (mm/mm)
Grafica 1. Comparación de modelo con adhesivo con resultados experimentales En la grafica 1, se puede ver como el resultado de la simulación de prueba se aproxima a los resultados experimentales, especialmente para la malla de 1/50 y, hasta cierto punto, la malla de 1/25, lo cual indica El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS que el modelo computacional es bueno y se puede proceder con las otras pruebas.
700 Modulo Elastico (MPa)
9000 8000 Modulo Elastico (MPa)
800
7000 6000 5000 4000
600 500 400 300 200 100
3000
0
2000
0
10
1000
20
30
Altura (mm)
0 0
0,2
0,4
0,6
Grafica 3. Modulo elástico del honeycomb vs altura de honeycomb con distintas mallas
Grosor (mm)
Grafica 1. Modulo elástico del honeycomb vs espesor de las paredes con distintas mallas
70 60 50 Esfuerzo plastico de compresion
40 30
Esfuerzo elastico maximo
20 10
Esfuerzo ingenieril (MPa)
Esfuerzo ingenieril (MPa)
80
4,5 4 3,5 3 2,5
Esfuerzo plastico de compresion
2 1,5
Esfuerzo elastico maximo
1 0,5 0 0
10
20
30
Altura del honeycomb (mm)
0 0
0,2
0,4
0,6
Grosor de lamina (mm)
Grafica 2. Esfuerzo elástico máximo del honeycomb y esfuerzo plástico de compresión vs espesor de las paredes con distintas mallas
Grafica 4. Esfuerzo elástico máximo del honeycomb y esfuerzo plástico compresivo vs altura del honeycomb con distintas mallas Estas graficas muestran que el tamaño de la malla en los rangos evaluados no afectan significativamente los resultados, por lo tanto se usara el promedio de estos datos a cada parámetro evaluado para los siguientes análisis. Se pueden presentar errores en la integración numérica debido a la geometría
El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS
Modulo Elastico (MPa)
12000 10000 8000 Simulaciones
6000
800 700
Modulo Elastico (MPa)
compleja del modelo, ya que produce extrapolación de los datos y nodos saliendo de la estructura [10].
600 500 400
Simulaciones
300
Chen
200
Merghani
100 0
Ashby
4000
0
Merghani
20
30
Altura (mm)
2000 0 0
0,2
0,4
Grafica 7. Comparación de los resultados promedios de modulo elástico de las simulaciones con los modelos teóricos según varia la altura del honeycomb.
0,6
Grafica 5. Comparación de los resultados promedios de modulo elástico de las simulaciones con los modelos teóricos según varia el espesor de las laminas. 250
200 150 Simulaciones
100
Wierzbicki 50
Esfuerzo Plastico de Compresion (MPa)
Grosor (mm)
Esfuerzo Plastico de Compresion (MPa)
10
3,5 3 2,5 2 Simulaciones
1,5
Wierzbicki
1 0,5 0 0
10
20
30
Altura (mm)
0 0
0,2
0,4
0,6
Grosor (mm)
Grafica 6. Comparación de los resultados promedios de esfuerzo plástico de compresión de las simulaciones con los modelos teóricos según varía el espesor de las láminas.
Grafica 8. Comparación de los resultados promedios de esfuerzo de compresión plástica de las simulaciones con los modelos teóricos según varía la altura del honeycomb.
Los resultados muestran una tendencia positiva en la relación entre los parámetros estudiados versus el espesor, algo El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS consistente con los modelos analíticos. También se ve como la altura afecta a los datos, aunque en menor medida y no siempre de la misma manera: negativa para el esfuerzo máximo elástico, positiva para el modulo elástico y no afecta al esfuerzo plástico compresivo. Estos 2 últimos resultados son concordantes con ciertos modelos teóricos hallados. Se puede ver también como los modelos analíticos desarrollados para el modulo elástico no se aproximaron a los resultados de la simulación, mientras que los del esfuerzo plástico compresivo se aproximan a los resultados para grosores pequeños, por lo cual es un indicio que las asunciones realizadas para desarrollar estos modelos no son válidas para espesores mayores después de cierto valor. Aun con esto, las tendencias mostradas son las mismas entre los modelos analíticos y los resultados computacionales, aunque varía que tan marcados se presenta. Errores en las simulaciones variando la altura generaron valores no tan cercanos a los valores analíticos para ciertos casos, pero los demás son acordes a lo esperado. Conclusiones Los resultados de las simulaciones de comparación, que son las de adhesivos con 0.025 mm de espesor y 12 mm de altura, se acercan bastantes a los experimentales, lo cual los valida. El ancho de las láminas usadas fue el factor más importante en el comportamiento de los honeycombs simulados, presentando una tendencia positiva significativa en todos los parámetros analizados.
presentando distintas tendencias en cada parámetro analizado. El uso de adhesivo aumenta los valores de los parámetros analizados en la relación a los resultados sin este, principalmente el esfuerzo compresivo plástico.
Los modelos teóricos utilizados para el cálculo del esfuerzo de compresión se presentaron validos para espesores de laminas pequeños Los modelos teóricos utilizados para el cálculo del modulo elástico de honeycombs no son los apropiados para compararse con los valores experimentales y de las simulaciones La variación del tamaño de malla puede afectar los resultados de las simulaciones, presentando divergencia en estos, pero los convergentes son relativamente cercanos entre ellos. Agradecimientos—A mi padre Nicolás y a mi madre Dilia, a mi hermano Pablo, a Gladys, a Brayan, a mi tía Mercedes y sus hijos José Antonio y Matías, a mis amigos, y al profesor Juan Pablo Casas. Referencias [1] Wilbert. A, Jang. W. Y, Kyriakides. S, Floccari, J. F. Buckling and progressive crushing of laterally loaded honeycomb. International Journal of Solids and Structures, Volumen 48, Ejemplar 5. Marzo 2011. [2] Bitzer. T. N. Honeycomb technology: Materials, Design, Manufacturing, Applications and Testing. Chapman & Hall. 1997
La altura, en menor medida, también afecta al comportamiento del elemento, pero El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Universidad de los Andes. Orrego, Casas. SIMULACION POR MEDIO DE ELEMENTOS FINITOS [3] Meo. M, Vignjevic. R, Marengo. G. The response of honeycomb sandwich panels under low-velocity impact loading. International Journal of Mechanics Sciences, Volumen 47 Ejemplar 9. Septiembre 2005 [4] Chen. D.H, Ozaki. S. Analysis of in-plane elastic modulus for a hexagonal honeycomb core: Effect of core height and proposed analytic method. Composite Structures. Volumen 88, Ejemplar 1. Marzo 2009. [5] Petrone. G, Rao. S, De Rosa. A, Mace. B.R, Franco. F, Bhattacharya. D. Behavior of fibrereinforced honeycomb core under low velocity impact loading. Composites Structures, Volumen 100. Junio 2013. [6] Avendaño. D.F, Casas-Rodriguez. J. P, Maranon. A. Dynamic characterization of
hexagonal cell structures. 14 th Pan-American Congress of Applied Mechanics. Marzo 2014. [7] Casas-Rodriguez. J.P, Ashcroft. I.A, Silberschmidt. V.V. Delamination in adhesively bonded CFRP joints: Standard fatigue, impact-fatigue and intermittent impact. Composites Science and Technology. Volumen 68, Ejemplar 12. Septiembre 2008. [8] Aaron-Jeyasingh, V.M. Analytic modeling of metallic honeycomb for energy absorption and validation with FEA. Wichita State University. Mayo 2005 [9] Volume A: Theory and User Information. Marc 2012. [10] Errors arising In FEA. Recuperado de http://www.tech.plym.ac.uk/sme/mech335/fea errors.htm el 10 de julio de 2014.
El CUARTO SIMPOSIO NACIONAL SOBRE MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTINUAS – SMEC 2014
Respuesta estructural de paneles tipo sándwich con núcleo de poliuretano sometido a impactos de baja velocidad Aguiar J., Useche, J. Universidad Tecnológica de Bolívar Cartagena, Colombia Resumen Este trabajo presenta los resultados de un estudio experimental y numérico del comportamiento de paneles tipo sándwich con núcleo de poliuretano y láminas de fibra de vidrio con tejido Woven Roving sometidas a impactos de baja velocidad utilizando indentadores planos, semiesféricos y cónicos, con distintos niveles de energía. La variación en el tiempo de la carga de impacto es monitoreada utilizando una celda de carga diseñada para tal propósito. El desplazamiento y velocidad del indentador fueron medidas utilizando una técnica de correlación no destructivas basadas igualmente en CDI para caracterizar la zona de daño después del impacto. Igualmente se utilizó seccionamiento de las probetas para reconstruir la zona afectada por el daño. Los efectos de la forma del indentador y energía inicial sobre la respuesta del panel y tipo de daño generado fueron estudiados. Buscando complementar los resultados experimentales, fueron modelados numéricamente utilizando el método de los elementos finitos a través del software LS-DYNA®. Fue utilizado el método de daño MAT169 implementado en este paquete. Los resultados numéricos obtenidos concuerdan en alto grado con los resultados experimentales hallados.
1. Introducción Los ingenieros se encargan de diseñar la mayoría de los productos manufacturados y los sistemas de elaboración necesarios y óptimos para su producción. Hay que tener en cuenta que los materiales son necesarios para fabricar productos, por lo tanto, los ingenieros deben conocer la estructura interna y las propiedades de los materiales, de tal manera que puedan elegir los más adecuados para cada aplicación y crear los mejores métodos para procesarlos. Al fin de cuentas, la producción y elaboración de los materiales hasta convertirlos en productos terminados, constituyen una parte importante de la economía actual. Expertos en investigación y desarrollo crean nuevos materiales o modifican las propiedades de los existentes. Los ingenieros de diseño usan materiales actuales, modificados o nuevos para diseñar y crear nuevos productos y sistemas. En ocasiones los ingenieros requieren de un nuevo material para su diseño, y la tarea de crearlo será encomendada a científicos e ingenieros especialistas en investigación [1]. El uso de los materiales compuestos ha revolucionado la industria, evolucionando a través del tiempo y continúan haciéndolo así como el hombre y la ingeniería. Remplazando otros tipos de
materiales en distintas aplicaciones, más que todo en las industrias aeronáutica, electrónica de la aviación, marítimas y navales, automotriz, de estructuras civiles y de uso deportivo. El énfasis de este trabajo se enfoca en aplicaciones para embarcaciones pequeñas ya que son más económicos, livianos, y con la configuración adecuada pueden llegar a tener resistencias bastante altas, incluso se prevé un aumento del 5% en el uso futuro de este tipo de materiales. Los ensayos de alta velocidad, el objetivo es simular el impacto que produce la colisión de una masa pequeña sobre una estructura. En este tipo de impactos, el efecto está muy localizado en la zona de alrededor del impacto, debido a la duración del impacto, por lo que las condiciones de contorno de la estructura no son determinantes en el comportamiento de la misma. Los impactos a baja velocidad, tienen como objetivo simular el impacto que produce la colisión de un objeto de masa considerable con una baja velocidad sobre una estructura. Debido a las duraciones de estos impactos, las ondas de tensión se propagan hasta el contorno del elemento y se reflejan varias veces durante el proceso de impacto por lo tanto la respuesta de la estructura es global, influyendo tanto su geometría como sus condiciones de contorno [2]. 2. Materiales y métodos La estructura sándwich es fabricada con paredes o pieles de fibra de vidrio y núcleo de poliuretano, de la cual se ha fabricado una lámina para dividirla y extraer varias probetas de 20x20 cms para posteriormente ser instaladas en la máquina. Para la realización de está lámina se necesitaron los siguientes componentes: 2Kg fibra MAT 600 1Kg fibra Woven Roving 4Kg de resina 1060 Catalizador 1Kg de Poliuretano (A y B) ¼ de lámina de Madecor Construcción de la lámina Para la construcción de la lámina lo primero que se necesitó construir fue un molde en madera tipo Madecor, de unos 55x55 cms y con altura de 2.5 cms. Las paredes de la lámina sándwich fueron hechas de fibra de vidrio con una configuración MatRoving-Mat. Luego que ya se tuviera el poliuretano (Solución A y B) y las dos láminas de fibra de vidrio se procede a realizar la lámina.
El molde fue cubierto con bolsa plástica y cinta adhesiva para evitar filtraciones y el daño del molde por el control de la expansión de la espuma.
Figura 1. Molde para la lámina hecho en madera tipo madecor.
Se colocan las dos paredes de fibra en la base del molde y en la tapa para que una vez se deposite el poliuretano, este quedará adherido a las paredes de fibra luego de la reacción.
Figura 2 - 3. Muestra del molde de madera con su tapa
Gracias a la colaboración de la empresa COTECMAR, se pudo obtener una muestra del material tipo sándwich que se utiliza en la empresa como material de refuerzo, la diferencia radica en que es usada una espuma especial para núcleos de materiales compuestos estructurales llamada Divinycell serie H100.
Figura 4. Muestra de estructura sándwich.
Figura 5. Muestra de estructura sándwich.
Figura 6. Panel estructural brindado por la empresa COTECMAR.
Construcción de la máquina de impacto La máquina de impacto ha sido construida para un proyecto anterior, y en estos momentos está siendo rediseñada para una mejor instrumentación.
Figura 7. Máquina de impacto a modificar Entre las modificaciones a realizar están: -
Diferentes tipos de intentador (plano, semiesférico y cónico) Cabezal con nuevo diseño con la posibilidad de añadir masa adicional Columnas más gruesas
Modelo numérico y experimental Para estas pruebas se utilizará un software de análisis numérico llamado LS-DYNA®, teniendo claro que son ensayos de impacto a baja velocidad que se realizarán con distintos tipos de intentador. Se realizarán tablas tipo: Carga- Desplazamiento, Energía de Impacto- Energía Absorbida, Índice de absorción- Energía de impacto, entre otras; todo esto con el fin de validar estos resultados con los obtenidos numérica y analíticamente. Como aclaración principal cabe resaltar que lo que se realizará en el en análisis numérico es con respecto a una lámina ideal de estructura sándwich con núcleo de espuma de poliuretano, de modo que se buscará la mejor adaptación del modelo a los resultados experimentales. Es importante conocer las condiciones reales de la lámina al realizar la comparación con el modelo numérico: -
La cantidad de soportes que se puedan colocar el momento de utilizar el software debe ser la misma que componen la máquina de impacto. No se usa adhesivo en este método, en reemplazo se usa la reacción misma de la espuma de poliuretano. En ventaja con otros núcleos para este tipo de estructuras, la homogeneidad de la espuma de la probeta es muy buena. La fibra de vidrio es hecha a mano, esto afecta directamente en la homogeneidad de la misma, por la posible creación de burbujas.
Figura 8. Ilustración de las dimensiones de las probetas para los ensayos de impacto a baja velocidad.
Figura 9. Análisis por FEM de los ensayos de impacto a baja velocidad en estructuras sándwich con núcleo de espuma de poliuretano [3].
REFERENCIAS [1] William F. Smith. Fundamentos de la ciencia e ingeniería de materiales. McGraw-Hill, México, D.F. (2006). Capítulo 1 Pg. 3. [2] Daniel rueda Lastres. Modelización analítica de paneles sándwich sometidos a impactos de baja velocidad. Universidad Carlos III de Madrid. Departamento de mecánica de medios continuos y teoría de estructuras. Leganés-Madrid, España. (2012). Cap 2.3 Impacto sobre estructuras sándwich. Pagina 9. [3] Jin Zhou, Mohamad Zaki Hassan, Zhongwei Guan, Wesley J. Cantwell. The low velocity impact response of foam-based sandwich panels. University Of Liverpool. Liverpool, UK. (2012). Pg 1784.
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
Evaluación del desempeño de espumas poliméricas de celda abierta rellenas con STF ante eventos explosivos Robinson-Luque, V.S.1, Casas Rodriguez, J.P.2 Departamento de Ingeniería Mecánica Universidad de los Andes, Bogotá, Colombia 1
[email protected],
[email protected]
Resumen En la presente investigación se evalúa la influencia del STF en el impulso generado por una carga explosiva (tasas de deformación >106 s-1 [1]) de 0.8 gr de PENT (≈1 gr equivalente de TNT) por medio del uso del péndulo balístico, el cual permite conocer (a través de su oscilación) el impulso generado por la energía residual de la onda de choque a tres distancias de separación (DS) diferentes (100, 150 y 200mm). Para esto, se elaboraron paneles tipo sándwich con placas (frontales y posteriores) de aluminio de 1mm de espesor con diferentes núcleos, i.e. espuma polimérica, espuma polimérica rellena de STF, espuma rígida SAN (estireno acrilonitrilo) y espumas SAN en configuración de densidad ascendente y descendente. Por otro lado, para la calibración se utilizó una placa acero CR de ½” de espesor con el propósito de encontrar el impulso máximo transmitido al péndulo. Los impulsos encontrados en la experimentación para los diferentes paneles fueron comparados encontrando que al rellenar la espuma polimérica con STF se puede lograr atenuar hasta un 35% más de energía para un DS de 100 mm con respecto a los otros núcleos estudiados. Sin embargo, al aumentar el DS dicha diferencia se pierde, encontrando que a un DS de 200 mm la espuma sin impregnar absorbe un 16% más energía que al ser impregnada con STF. Palabras Clave: STF (Shear Thickening Fluid), Espumas poliméricas, Espumas de celda abierta rellenas con fluido (FFOCF), Espumas rígidas (SAN) altas deformaciones, explosiones, péndulo balístico.
Abstract This paper evaluates the influence of STF over the impulse generated by an explosive charge (strain rate >106 s-1 [1]) of 0.8 gr of PENT (≈1 gr equivalent of TNT) employing a ballistic pendulum, which allows to identify (by its oscillation) the impulse generated by the residual energy of the shock wave at three different standoff distances (DS = 100, 150 and 200mm). In order to develop the experiments, different sandwich panels were developed by using aluminum plates of 1mm of thickness with different cores, i.e. polymeric foam, STF filled polymeric foam, rigid SAN foam (Styrene-acrylonitrile) and ascendant and descendant density SAN foam. On the other hand, to achieve the calibration of the pendulum a ½” CR steel plate was used in order to find the maximum transmitted impulse. The obtained impulses for the different tested panels were compared to each other; finding that STF filled foam can attenuate up to 35% more energy for a DS of 100mm in relation to the other studied cores. However, by increasing the DS that difference decrease obtaining that at DS of 200 mm the non-impregnated foam absorbs 16% more energy than the STF filled foam. Key Words: STF (Shear Thickening Fluid), Polymeric Foams, Fluid Filled Open Foam (FFOCF), Rigid Foams (SAN), High strains, Blast, ballistic pendulum
Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
1
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
1. INTRODUCCIÓN La necesidad de desarrollar diferentes materiales y mecanismos con la capacidad de absorber energía ha cobrado gran interés en las últimas décadas debido a la cantidad de aplicaciones que estos presentan, especialmente en la protección tanto personal como de infraestructuras contra fenómenos dinámicos tales i.e. impactos y explosiones. Entre los materiales más utilizados para este fin se encuentran las espumas, las cuales gracias al comportamiento de su microestructura bajo cargas compresivas, pueden mantener un esfuerzo constante en un rango amplio de deformaciones; este fenómeno resulta ser de gran utilidad para absorber energía, además de funcionar como amortiguador atenuando las aceleraciones a un cuerpo. Sin embargo, aunque las espumas presentan un gran potencial para absorber energía (entre otras aplicaciones), existen numerosos escenarios prácticos en los que el espesor, el peso y/o la rigidez de la espuma limita la implementación de espumas rígidas en estructuras tipo sándwich [2-4]. Generalmente, una espuma de celda abierta se encuentra rellena de aire, el cual no genera una contribución considerable al comportamiento mecánico de la misma frente a diferentes cargas [2]. Por esta razón diferentes autores han propuesto rellenar estas espumas con diferentes fluidos (Newtonianos y No-Newtonianos) [5, 6] con el fin de mejorar la respuesta dinámica de este tipo de materiales; además, este relleno líquido también permite lograr un mayor control sobre la disipación de energía por medio del fluido viscoso ante una deformación (compresiva) mucho menor. Uno de los fluidos que ha generado mayor interés para este tipo de aplicaciones es el STF (Shear Thickening Fluid) debido principalmente a su habilidad para cambiar de un fluido con una baja viscosidad a uno de alta viscosidad llegando a un estado casi sólido. Es precisamente este comportamiento el que hace que este tipo de materiales sea deseado en las aplicaciones de absorción de energía [6-8], pues mientras éstos se encuentran en condiciones normales (baja viscosidad) le permiten al usuario tener una gran movilidad y flexibilidad, pero cuando el STF se ve sometido a impactos (altos esfuerzos) empieza la transición de su viscosidad, por medio de la cual se absorbe parte de la energía generada por el impacto y ayuda a desviar la energía restante. Además, mediante la selección de los diferentes componentes del STF, i.e. el tipo y el tamaño de las partículas o el disolvente utilizado es posible seleccionar y sintonizar el inicio de la transición de la viscosidad para una aplicación específica [7]. Sin embargo, la efectividad en la absorción de energía por parte de este tipo de fluidos depende de la interacción que presenta con la geometría que lo contiene (o en la que se encuentra impregnado). Las espumas elásticas de celda abierta son uno de los mejores contenedores para el STF debido a que el daño en su microestructura es mínimo incluso cuando éstas se encuentran sometidas a grandes deformaciones, motivo por el cual su capacidad para contener un fluido no se ve afectada de manera significativa. Además, tanto la espuma elástica como el STF poseen un proceso altamente reversible motivo por el cual pueden ser implementados en aplicaciones de cargas cíclicas [4]. Con el fin de estudiar el comportamiento que presentan los materiales ante este fenómeno se han desarrollado diferentes dispositivos, entre los cuales se encuentran el péndulo balístico, el cual, es una técnica de medición ampliamente utilizada por distintos autores [9-13] para obtener el impulso generado por cualquier carga, pues este responde al momento total entregado al área de la superficie de una masa Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
2
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
movible. Dicho momento incluye la explosión en el aire, cualquier fragmento y cualquier partícula en el fluido [14]. Entre los estudios elaborados se destaca el de Humphreys [9] quien desarrolla una serie de pruebas a placas de acero de diferentes tamaños y propiedades las cuales son fijadas a un péndulo balístico. Para esta prueba el autor instaló una cámara de alta velocidad sobre el péndulo con el objetivo de registrar las deformaciones sufridas por las placas en tiempo real (Figura 1). En este estudio el autor reportó que la explosión tiene una duración de 150 μs empezando en el segundo cuadro (desde la derecha) y terminando en el cuarto sin presentar deformación alguna en la placa. Para el sexto cuadro ya se ha alcanzado una deformación plástica final en la placa en un lapso de tiempo de alrededor de 0.0005s, el cual es extremadamente menor comparado a los 4.43s del periodo natural del péndulo. De esta forma se evidencia que la deformación plástica en la placa ha ocurrido mucho antes de que el péndulo balístico alcance su máxima oscilación (balanceo).
Figura 1. Movimiento actual de la placa (85 μs por cuadro) [9]
De la misma manera cabe destacar el trabajo desarrollado por en la unidad de investigación de impactos, explosiones y supervivencia (BISRU por sus siglas en ingles) de la Universidad de Ciudad del Cabo (UCT) encabezada por el profesor Gerald Nurick quienes han sido los líderes en la implementación de este dispositivo en el análisis del comportamiento de materiales ante cargas explosivas. Entre los trabajos desarrollados se encuentra el presentado por Jacob et al [10] quienes evalúan la influencia de la distancia de separación entre el explosivo y una placa circular de 53 mm de diámetro, encontrando que para distancias de separación inferiores al radio de placa (13-40 mm), la carga explosiva es considerada como localizada, mientras que al ser la distancia de separación superior al diámetro (100-300 mm) la carga se considera como uniformemente distribuida sobre toda la placa. Adicionalmente los autores concluyen que las deflexiones en el punto medio decrecen significativamente a medida que la distancia de separación aumenta de 13 a 50 mm para una misma carga, pero cuando las distancias se encuentran entre 75 y 300 mm las deflexiones en el medio son similares.
Figura 2. Fotografía de las capas en ordenadas de manera secuencial con respecto a la distancia de separación (DS) [10] Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
3
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
Hassan et al [12] investiga la influencia de la variación de la densidad de núcleo de espumas de PVC con placas de aluminio con un espesor de 1.7 mm adheridas en la parte frontal y posterior de la espuma utilizando un adhesivo de curado rápido (Timbond). En su estudio Hassan et al presentan que el daño al interior de los paneles tipo sándwich es mucho más severa cuando la densidad del núcleo aumenta, a diferencia de la espuma de menor densidad (60 kg/m3) en la que no se evidencia desprendimiento entre las capas o de fractura en el núcleo en los rangos de impulso trabajados.
(a)
(b)
Figura 3. Sección transversal de las estructuras tipo sándwich con una densidad de núcleo de (a) 60 kg/m3 y (b) 200 kg/m3 para 4 diferentes impulsos [12].
El presente artículo tiene como objetivo evaluar el desempeño de diferentes núcleos, en especial la espuma de celda abierta rellena de STF en términos del impulso transmitido y peso del compuesto con el fin de observar sus beneficios y limitaciones al ser utilizado como equipo de protección para atenuar el impulso generado por la detonación de una carga explosiva. 2. MATERIALES Los paneles tipo sándwich constan de dos placas de aluminio los cuales cubre el núcleo. El aluminio trabajado es 6061 con un espesor de 1 mm, mientras que los núcleos utilizados se presentan en la Tabla 1. Tabla 1. Núcleos utilizados
Nomenclatura
Matriz
Observaciones
Espuma
Espuma elastomérica de poliuretano
Espuma rosada
STF
Espuma elastomérica de poliuretano
SAN
Espuma rígida SAN
Rellena de STF (fracción de volumen del 58.13%) Referencia A800
SANas
Espuma rígida SAN (Ascendente)
Referencias A600-A800-A1200
SANdes
Espuma rígida SAN (Descendente)
Referencias A1200-A800-A600
Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
4
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
En este punto cabe mencionar que la nomenclatura presentada en la Tabla 1 será utilizada a lo largo del presente documento. 2.1. Preparación de las muestras Primero, se cortó la lámina de aluminio en cuadrados 120mm de lado, los cuales fueron mecanizados para abrir los agujeros con el fin de sujetar las muestras al péndulo. Posteriormente se cortaron cuadrados 72mm de lado para los diferentes núcleos utilizados. Los espesores de cada sándwich y su configuración son presentados en la Figura 4.
(a)
(b)
(d)
(c)
(e)
Figura 4. Sección transversal de las estructuras tipo sándwich utilizadas con un núcleo de (a) espuma, (b) STF, (c) SAN (A800), (e) SANas y (d) SANdes.
Por otro lado, para adherir las capas del compuesto se utilizó Sikaflex 221 con un espesor de 0,5 mm (aproximadamente) con el fin de asegurar que la deflexión encontrada en las láminas de aluminio sea la misma que la sufrida por el núcleo. 3. PROCEDIMIENTO EXPERIMENTAL Como se mencionó anteriormente, para el presente estudio se decidió trabajar con una carga explosiva constante de 0,8 gr de PENT para todas las pruebas desarrolladas, la cual es la carga de un detonador comercial N°8. El impulso transmitido será determinado por medio del péndulo balístico. 3.1. Péndulo balístico El péndulo balístico consiste generalmente de una masa suspendida por cuatro cables de acero (ver Figura 5). En este caso particular, el péndulo es impulsado al detonar una carga explosiva al interior del tubo de choque, ubicado en la parte frontal del módulo de pruebas. El principio de funcionamiento parte de la conservación de momentum, es decir, durante la explosión el impulso generado por dicho fenómeno causará que el péndulo se balancee; ese balanceo permite conocer el impulso conociendo la geometría del péndulo. Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
5
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
Figura 5. Péndulo balístico convencional desarrollado por GIE Uniandes
Con el fin de garantizar el correcto funcionamiento del péndulo balístico, es necesario garantizar que los 4 cables de acero que lo soportan se encuentren sometidos a la misma tensión, con esto se garantizará que el recorrido del péndulo solo presente translación y se pueda calcular el impulso. Para esto se utiliza un nivel ubicado en el centro geométrico del rectángulo generado por los 4 puntos de anclaje tanto en dirección horizontal como vertical. Por otro lado, de ser necesario, la longitud de los cables de acero puede ser modificada para nivelar el péndulo por medio de los tensores de tipo comercial que se conectan directamente con el péndulo. De esta manera el impulso es transmitido a través del centroide del péndulo. 3.2. Módulo de pruebas En la Figura 5 se presenta el módulo de pruebas del péndulo balístico, el cual consta de unas placas de soporte en medio de las cuales se sujetarán las muestras/especímenes a ser probados. Dichas placas de soporte se conectan al péndulo por medio de 4 barras espaciadoras las cuales deben ser lo suficientemente largas para permitir la deformación de las muestras. 3.3. Carga explosiva Con el propósito de garantizar que el detonador se ubique en el centro del tubo de choque se utiliza un disco de poliestireno expandido de ½” de espesor y un diámetro igual al diámetro interno del tubo; en el centro de dicho disco se realiza un agujero del mismo diámetro que el detonador utilizado. Se decidió utilizar poliestireno expandido debido a que este material no afecta las propiedades de la onda de choque por su baja densidad y a que este material se desintegra totalmente debido a las altas temperaturas generadas producto de la detonación.
(a)
(b)
Figura 6. Esquema del montaje para (a) la calibración y (b) los compuestos tipo sándwich. Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
6
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
En la Figura 6 se presentan los esquemas de los montajes utilizados para las pruebas desarrollados. La distancia de separación (DS) se encuentra representada en la misma figura como “x”. En la Figura 6a se presenta el esquemático del montaje para la calibración del péndulo. Para la calibración se utilizará una placa de acero estructura de ½” de espesor la cual será la encargada de recibir toda la energía producto de la detonación del explosivo; dado que dicha placa es rígida se asegura que toda la energía producida por el explosivo sea transmitida al péndulo sin presentar perdidas asociadas a la deformación de la muestra. Dicha placa se asegura al módulo de pruebas mediante 6 pernos de 3/8”. La masa total del péndulo es de 73.04 kg incluyendo la placa de calibración. En la Figura 6b se presenta el esquemático del montaje de pruebas para el análisis del comportamiento de los diferentes núcleos. Como se mencionó en la sección 2.1 los compuestos tipo sándwich fueron desarrollados con placas de aluminio tanto frontal y posterior con un espesor de 1mm. Adicionalmente, se puede observar la presencia de placas de separación con las cuales se fija la distancia entre las placas de sujeción cuando las muestras tienen muy poca rigidez o cuando se dispone de poco material para ser atornillado a dichas placas, de esta manera se garantiza que la muestra soporte toda la carga de la detonación sin verse afectada por presiones producidas al apretar demasiado el espécimen. 3.4. Matriz de pruebas En la Tabla 2 se presenta la matriz de pruebas utilizada. En esta se presenta la distancia de separación (DS) y el número de repeticiones realizadas para cada caso, según la nomenclatura presentada en la Tabla 1. Tabla 2. Matriz de pruebas
Núcleo
Repeticiones
DS (mm)
2
200
2
150
2
Calibración
Espuma
STF
Núcleo
Repeticiones
DS (mm)
2
200
2
150
100
2
100
2
200
2
200
2
150
2
150
2
100
1
100
2
200
2
200
1
150
1
150
2
100
2
100
SAN
SANas
SANdes
4. ANÁLISIS Y RESULTADOS. En la presente sección se presentan los resultados obtenidos de las pruebas experimentales realizadas con el péndulo balístico. 4.1. Análisis de impulsos En la Figura 7 se presentan los impulsos obtenidos para las pruebas efectuadas. Los valores encontrados fueron normalizados con respecto al impulso promedio más alto con el fin de tener una comparación directa entre todas las pruebas realizadas. Las Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
7
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
pruebas con repetición fueron ponderadas y la media encontrada es el valor reportado en las gráficas suministradas. Al comparar los valores de calibración se observa una caída en el impulso a medida que el DS aumenta, lo cual es consecuente con la teoría de propagación de ondas acústicas en el aire; esta respuesta demuestra una debida elaboración de las pruebas. En la Figura 7a se puede observar que la espuma rellena de STF presenta una reducción del impulso de casi el 35% a comparación de los otros núcleos estudiados. Sin embargo a medida que el DS aumenta esta porcentaje se reduce hasta el punto de encontrar que a un DS de 200mm (Figura 7c) la espuma polimérica sin impregnar absorbe un 16% más del impulso generado que su contraparte rellena de STF. Este fenómeno puede deberse a que el STF al ser un fluido dilatante responde en función de su tasa de deformación, es decir, al estar el detonador más cerca del compuesto dicha tasa crece y los mecanismos de disipación de energía del STF se vuelven mucho más efectivos.
(a)
(b)
(c) Figura 7. Impulsos normalizados a diferentes DS para la placa de calibración los paneles tipo Sándwich estudiados.
4.2. Análisis de impulsos en función del peso del compuesto. Por otro lado, resulta de gran interés el evaluar la respuesta de dichos paneles en función de su masa, puesto que en muchas de las aplicaciones de protección el peso resulta ser de gran interés, i.e. cuando se desea realizar protección personal, el peso del equipo puede limitar el movimiento del usuario y es incómodo. De esta manera, en la Figura 8 se presentan los impulsos normalizados con respecto a la masa del compuesto. En la Figura 8a se puede apreciar que la energía que se transfiere al péndulo por cada kg de compuesto es 58.2% menor para el caso de la espuma con STF con respecto a su contraparte (espuma sin relleno), lo cual resulta ser muy importante si se tiene en cuenta que al rellenar la espuma con STF su masa aumenta alrededor de 1.6 veces.
Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
8
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
Por otro lado, al observar la Figura 8b y la Figura 8c se observa que la diferencia presentada anteriormente disminuye de manera drástica por efectos de la distancia de separación, siendo siempre la espuma con relleno de STF la que mejor relación impulso-peso presenta.
(a)
(b)
(c) Figura 8. Impulsos normalizados en función del peso del panel estudiado a diferentes DS.
5. CONCLUSIONES En la presente investigación se evaluaron diferentes núcleos de un compuesto tipo sándwich en función del impulso transferido encontrado a través de la implementación de un péndulo balístico. Los resultados fueron presentados en la Figura 7 y en la Figura 8, en las que se encontró que al impregnar una espuma polimérica de poliuretano con STF la capacidad para absorber energía incrementa hasta un 35%, cuando la carga explosiva se encuentra a una distancia de separación de 100mm. Sin embargo, como la activación del STF depende de la tasa de deformación a la que se ve sometida, a medida que se aleja la carga explosiva el porcentaje de energía absorbida disminuye hasta presentar un aumento en la energía transferida del 16% con respecto a su contraparte sin impregnar. Finalmente, a pesar del aumento en la energía transferida por parte de la espuma impregnada con STF este núcleo presenta la menor relación de energía transferida por masa lo cual es un criterio importante en el diseño de equipos de protección. Por otro lado la configuración ascendente de la espuma tipo SAN también presenta una gran capacidad para absorber energía, además de tener la capacidad para brindarle a las estructuras en las que se utilice mayor rigidez. 6. AGRADECIMIENTOS Los autores desean agradecer a la industria militar Colombiana (INDUMIL) por el apoyo brindado a esta investigación.
Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
9
CUARTO SIMPOSIO EN MECÁNICA DE MATERIALES Y ESTRUCTURAS CONTÍNUAS SMEC 2014
7. REFERENCIAS [1] M. A. Meyers, Dynamic behavior of materials. New York: John Wiley & Sons, Inc., 1994. [2] L. J. Gibson and M. F. Ashby, Cellular Solids: Structure and Properties: Cambridge University Press, 1999. [3] K. Kitagawa, K. Takayama, and M. Yasuhara, "Attenuation of shock waves propagating in polyurethane foams," Shock waves, vol. 15, pp. 437-445, July 2006. [4] M. A. Dawson, "Modeling the Dynamic Response of Low-Density, Reticulated, Elastomeric Foam Impregnated with Newtonian and Non-Newtonian Fluids," Doctor of Philosophy in Mechanical Engineering PhD, Department of Mechanical Engineering, Massachusetts Institute of Technology., 2008. [5] M. A. Dawson, G. H. McKinley, and L. J. Gibson, "The Dynamic Compressive Response of Open-Cell Foam Impregnated With a Newtonian Fluid," Journal of Applied Mechanics, vol. 75, p. 11, 2008. [6] M. A. Dawson, G. H. McKinley, and L. J. Gibson, "The Dynamic Compressive Response of an Open-Cell Foam Impregnated With a Non-Newtonian Fluid," Journal of Applied Mechanics, vol. 76, p. 8, 2009. [7] G. Bettin, "Energy Absorption of Reticulated Foams Filled with Shear-Thickening Silica Suspensions," Master of Science in Mechanical Engineering MSc, Department of Mechanical Engineering, Massachusetts Institute of Technology, 2005. [8] A. S. Lim, S. L. Lopatnikov, N. J. Wagner, and J. W. Gillespie, "An experimental investigation into the kinematics of a concentrated hard-sphere colloidal suspension during Hopkinson bar evaluation at high stresses," Journal of Non-Newtonian Fluid Mechanics, vol. 165, p. 9, 2010. [9] J. S. Humphreys, "Plastic Deformation of Impulsively Loaded Straight Clamped Beams," Journal of Applied Mechanics, vol. 32, pp. 7-10, 1965. [10] N. Jacob, G. N. Nurick, and G. S. Langdon, "The effect of stand-off distance on the failure of fully clamped circular mild steel plates subjected to blast loads," Engineering Structures, vol. 29, pp. 2723-2736, 2007. [11] I. M. Snyman, "Impulsive loading events and similarity scaling," Engineering Structures, vol. 32, pp. 886-896, 2010. [12] M. Z. Hassan, Z. W. Guan, W. J. Cantwell, G. S. Langdon, and G. N. Nurick, "The influence of core density on the blast resistance of foam-based sandwich structures," International Journal of Impact Engineering, vol. 50, pp. 9-16, 2012. [13] F. de Borbón and D. Ambrosini, "Dynamic response of composites sandwich plates with carbon nanotubes subjected to blast loading," Composites Part B: Engineering, vol. 45, pp. 466-473, 2013. [14] C. E. Needham, Blast Waves: Springer Verlag, 2010.
Cuarto Simposio en Mecánica de Materiales y Estructuras Continuas - SMEC 2014
10