Estudio por dinamica molecular de la resistencia a la falla por traccion de un modelo estadistico inspirado en el bambu Guadua angustifolia

Vista previa del material en texto

Estudio por dinámica molecular de la resistencia a la falla por tracción
de un modelo estad́ıstico inspirado en el bambú Guadua angustifolia
Universidad Nacional de Colombia
Facultad de Ciencias
Departamento de F́ısica
Trabajo Final de Maestŕıa presentado por:
Juan Alejandro Pérez Rangel
Asesor: Dr. Rer. Nat José Daniel Muñoz
2016
Índice general
7
Agradecimientos 8
Resumen 9
Abstract 10
Introducción 11
1. La guadua 16
1.1. Estructura de la guadua . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
1.2. Propiedades mecánicas de la guadua . . . . . . . . . . . . . . . . . . . . . . . . 18
2. Marco Teórico 22
2.1. Dinámica molecular . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2. Modelos estad́ısticos de fractura . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.2.1. Modelo de haces de fibra . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
3. Modelo del tejido parénquima en dirección longitudinal 30
3.1. Construcción del tejido de Parénquima . . . . . . . . . . . . . . . . . . . . . . . 30
3.2. Fuerzas y torques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.2.1. Fuerzas de contacto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.2.2. Fuerzas de los resortes . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
3.2.3. Fuerzas viscosas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
3.2.4. Torques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.3. Algoritmo de integración . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.4. El proceso de falla . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.5. Tracción axial y flexión . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
4. Resultados 38
4.1. Normalización de las constantes elásticas . . . . . . . . . . . . . . . . . . . . . . 39
4.2. Resultados de ensayos de tracción . . . . . . . . . . . . . . . . . . . . . . . . . . 40
4.2.1. Histogramas de avalanchas . . . . . . . . . . . . . . . . . . . . . . . . . . 40
2
Índice general 3
4.2.2. Escalamiento finito . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
4.2.3. Falla estructural en procesos de tracción . . . . . . . . . . . . . . . . . . 42
4.3. Resultados de ensayos de flexión . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
4.3.1. Histogramas de avalanchas . . . . . . . . . . . . . . . . . . . . . . . . . . 43
4.3.2. Escalamiento finito . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
4.3.3. Falla estructural en procesos de flexión . . . . . . . . . . . . . . . . . . . 45
Conclusiones 47
Bibliograf́ıa 50
Anexos 54
A. Unidades y parámetros de las simulaciones 55
B. Desarrollo del código para las simulaciones 61
Índice de figuras
1. Mapa de distribución del bambú en el mundo . . . . . . . . . . . . . . . . . . . 12
2. Estructura interna del tejido de parénquima . . . . . . . . . . . . . . . . . . . . 14
1.1. Culmmo de Guadua angustifolia según su edad . . . . . . . . . . . . . . . . . . 16
1.2. Partes de las plantas de guadua en la dirección longitudinal . . . . . . . . . . . . 17
1.3. Estructura anatómica transversal del culmo de Guadua angustifolia . . . . . . . 18
1.4. Ensayo de flexión a una probeta de guadua . . . . . . . . . . . . . . . . . . . . . 20
1.5. Módulo de ruptura y de elasticidad de la guadua ante ensayos de flexión según
la edad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.1. Interacción por contacto entre dos part́ıculas en dinámica molecular . . . . . . . 23
2.2. Poligonos 2D en DM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.3. Optimización de las fuerzas de contacto con esferorectánculos . . . . . . . . . . . 24
2.4. Esquema del modelo de haces de fibra . . . . . . . . . . . . . . . . . . . . . . . . 27
2.5. Modelo de contracción por secado de la parénquima . . . . . . . . . . . . . . . 27
2.6. Dinámica de fractura para un modelo de contracción de Guadua angustifolia . . 28
3.1. Modelo axial del tejido de parénquima . . . . . . . . . . . . . . . . . . . . . . . 31
3.2. Modelo de fuerzas de de contacto entre células . . . . . . . . . . . . . . . . . . . 32
3.3. Distribución de umbrales de ruptura en los enlaces celulares del tejido de
parénquima . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.4. Falla estructural del modelo bajo tracción axial y o flexión . . . . . . . . . . . . 35
3.5. Implementación de la deformación por flexión . . . . . . . . . . . . . . . . . . . 36
4.1. Configuraciones de tamaño en el modelo del tejido de parénquima . . . . . . . . 38
4.2. Influencia de una mayor rigidez entre tubos, como efecto de un mayor número
Nh de enlaces horizontales, sobre el promedio de la máxima avalancha 〈∆max〉 . 39
4.3. Influencia de una mayor rigidez entre tubos, como efecto de un mayor nmero Nh
de enlaces horizontales, sobre la máxima deformación bajo tracción 〈�max〉 . . . . 40
4.4. Histogramas de avalanchas ∆ para sistemas 32× 32 bajo pruebas de tracción axial 40
4.5. Escalamiento de 〈∆max〉, para pruebas de tracción . . . . . . . . . . . . . . . . . 41
4.6. Máxima deformación bajo tracción para igual rigidez . . . . . . . . . . . . . . . 42
4.7. Histograma de avalanchas: sistemas con células de igual tamaño . . . . . . . . . 43
4
Índice de figuras 5
4.8. Comparación de 〈∆max〉 según el tamaño de las células en el tejido de parénquima
para las pruebas flexión. El factor η corresponde con el promedio de enlaces
horizontales entre dos tubos adyacentes. (a) Sistemas 16×16. (b) Sistemas 32×32. 44
4.9. Escalamiento de 〈∆max〉, para esfuerzos de flexión . . . . . . . . . . . . . . . . . 45
4.10. Máxima deformación bajo flexión con sistemas de igual rigidez . . . . . . . . . . 45
4.11. Máxima deformación bajo flexión con sistemas de diferente rigidez . . . . . . . . 46
Indice de tablas
1.1. Comportamiento del MOE axial y MOR axial para diferentes regiones de la
pared del culmmo de la guadua Moso bambú . . . . . . . . . . . . . . . . . . . . 19
4.1. Parámetros del escalamiento de 〈∆max〉 para las pruebas de tracción axial y flexión. 44
A.1. Caracteŕısticas f́ısicas del tejido de parénquima . . . . . . . . . . . . . . . . . . . 55
A.2. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 8× 8 con células pequeñas . . . . . . . . . . . . . . . . 57
A.3. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 8× 8 con células grandes. . . . . . . . . . . . . . . . . 57
A.4. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 8× 8 con células de dos tamaños. . . . . . . . . . . . . 58
A.5. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 16× 16 con células pequeñas. . . . . . . . . . . . . . . 58
A.6. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 16× 16 con células grandes. . . . . . . . . . . . . . . . 58
A.7. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 16× 16 con células de dos tamaños. . . . . . . . . . . . 59
A.8. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 32× 32 con células pequeñas. . . . . . . . . . . . . . . 59
A.9. Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 32× 32 células grandes. . . . . . . . . . . . . . . . . . 59
A.10.Parámetros numéricos y unidades básicas de medida (UBM) para las
simulaciones de sistemas 32× 32 con células de dos tamaños. . . . . . . . . . . . 60
6
Dedicado a Lucy y Maŕıa, por todas las horas que invert́ı en este trabajo y que les
pertenećıan. Recuerden que lo único que nunca olvido, es que lasamo.
7
Agradecimientos
Al profesor José Daniel Muñoz, por su valioso tiempo, por sus aportes fundamentales
en cada reto que se asumió, por su paciencia y sobre todo por confiar en mis capacidades
como investigador. Sus enseñanzas y su ejemplo me quedan como insumos principales para la
construcción de mi carrera. A los compañeros del grupo de Simulación de Sistemas F́ısicos,
por aportar un punto de vista cŕıtico y respetuoso hacia mi trabajo. Agradezco también a la
Universidad Nacional y a Colciencias, por brindarme el escenario y los recursos para asumir
este proyecto.
A Lucy y Maŕıa, por apoyarme siempre y asumir mis dificultades como si fueran propias.
Este logro también es de ustedes. Recuerden que lo único que nunca olvido, es que las amo.
8
Resumen
El bambú Guadua angustifolia es un material de bajo costo y gran beneficio ambiental con
propiedades mecánicas excepcionales, lo que ha motivado su uso en construcción, diseño de
muebles y artesańıas. Su excelente resistencia ante esfuerzos de tracción y flexión se atribuye
principalmente su estructura fibrada. Sin embargo, el tejido blando de parénquima en el que se
insertan las fibras también presenta una estructura interesante, formada por tubos en dirección
axial en los que se células largas y cortas se alternan aleatoriamente. El objetivo del presente
trabajo es establecer si esta alternancia aleatoria de células largas y cortas aporta una ventaja
estructural a la resistencia mecánica de la G. angustifolia ante deformaciones de tracción axial
o de flexión. Para ello, se construyó un modelo bidimensional de elementos discretos (DEM)
en el que las células se representan por esferorectángulos organizados en tubos y unidos o a
través de resortes (que representan las lamelas entre células); cada resorte se rompe cuando
la fuerza que ejerce alcanza un valor umbral. Los umbrales de ruptura se asignan a partir
de una distribución amplia o estrecha de valores, que representan altos o bajos niveles de
desorden estructural respectivamente. Con este modelo se realizan simulaciones de pruebas de
tracción axial y de flexión para tres clases de sistemas: de células largas, de células cortas
o que alternan los dos tamaños al azar (como en el tejido real) para diferentes valores de
desorden en los enlaces. Los resultados muestran en todos los casos evidencia de la existencia
de una transición de fase entre un comportamiento tenaz del material a bajos desordenes en
los enlaces y un comportamiento dúctil a altos desordenes, con el tamaño de la avalancha
de fracturas más grande como parámetro de orden, coincidiendo con trabajos anteriores. Sin
embargo, el resultado más importante se obtiene al constatar que el modelo de dos tamaños de
células resiste deformaciones de tracción axial más grandes que los sistemas de células uniformes
cuando el desorden en los enlaces es muy bajo, alcanzando valores incluso mucho mayores que
los de sistemas con solo células largas. Los resultados sugieren, entonces, que la inserción al
azar de células cortas en un tejido originalmente compuesto solamente de células largas genera
un desorden estructural que compensa la uniformidad que se tiene cuando todos los umbrales
de ruptura son similares, y evita aśı la fragilidad t́ıpica de sistemas muy ordenados. En cuanto
a la flexión no se observa una ventaja evidente al incluir células de dos tamaños en el tejido.
Estos resultados refuerzan la suposición de que se requiere desorden estructural para mejorar
la resistencia del material a la falla e ilustra la potencia de los métodos de elementos discretos
para estudiar estos materiales.
9
Abstract
The bamboo Guadua angustifolia is a low-cost wood of great environmental significance and
exceptional mechanical properties that has been widely used in buildings, furniture and craft
works. Its excellent resistance in pulling and bending tests is usually ascribed to the fibres in the
wood, but the parenchimous soft tissue around the fibres also shows an interesting structure,
with alternating long and short cells along parallel tubes in axial direction. The main goal of
the present work is to establish if this random mixture of sizes gives any mechanical advantage
to avoid fail in pulling and bending test. A two-dimensional discrete element model (DEM)
was built to investigate for this hypothesis. The tissue is modelled as a set of shperorectangles
arranged in tubes and connected by springs (representing the lamellas among cells). Each spring
breaks when its force surpasses a threshold value, randomly chosen from either broad or narrow
intervals, representing high or low levels of structural disorder in the bonds, respectively. This
model was used to simulate pulling and bending tests for three kinds of systems: with long
cells, short cells or a random combination of both of them, for several values of structural
disorder. All results show evidence of a phase transition between a low-disordered fragile phase
and a high-disordered resistant phase, with the size of the largest avalanche as order parameter,
in agreement with previous studies. Nevertheless, the most valuable result is that the mixed
system is more resistant in pulling tests than the homogeneous ones when the disorder in
the breaking thresholds is low. This result suggests that inserting small cells at random in
a tissue of long cells creates a structural disorder that compensate for the typical fragility
of ordered systems when all bonds thresholds are similar (as in the real tissue). In bending
tests no similar difference is observed. These results support the general assumption that some
structural disorder is necessary to better resist against fail at high strains and, illustrates the
power and possibilities of DEM methods to investigate these materials.
10
Introducción
El bambú, junto al arroz y la caña de azúcar, es una de las especies mas importantes de
gramineas a nivel mundial. Al ser una planta de clima tropical, está presente en varios páıses
americanos. Entre ellos, Brasil presenta la mayor diversidad, con un total de 141 especies de
bambúes leñosos, y le siguen Colombia, con un total de 72 especies; Venezuela, con 60; Ecuador
con 44, y Costa Rica y México, con 39 especies leñosas, respectivamente [Lon05]. Uno de los
campos de mayor aplicación de la guadua es la construcción. Sus usos van desde estructuras en
varas, en placas de concreto como material de refuerzo, en laminados o como material principal
en construcciones de gran riqueza arquitectónica.
Todas estas aplicaciones son posibles gracias a las excepcionales propiedades mecánicas de la
guadua, como por ejemplo su excepcional resistencia ante esfuerzos de tracción y flexión. Esta
respuesta mecánica macroscópica debe ser consecuencia de su estructura a nivel microscópico,
pues los materiales que la forman no son suficientemente ŕıgidos como para explicar por śı
solos el comportamiento en bloque. En efecto, la guadua se destaca como otro de los materiales
naturales que se rompen con mucha dificultad, a pesar de que están hechos de substancias
que son consideradas frágiles, mas del 50 % de su composición es tejido blando denominado
parénquima. Otro ejemplo es el nácar, compuesto por carbonato de calcio (es decir, tiza, en
un 95 %), pero estructurado microscópicamente en forma de pequeñas columnas logra una
enorme resistencia [BB05]. La dentina, que es uno de los tejidos más duros del cuerpo humano,
está hecho también de arreglos especiales de carbonato de calcio. En una telaraña, los hilos
tienen resistencias mecánicas muy diversas, que junto con la respuesta viscoelástica del material,
permiten que la estructura se comporte de manera excepcional ante cargas distribuidas, como
las generadas por el viento [SWC12]. En todos estos materiales, el desorden en la estructura
o las inhomogeneidades en sus parámetros a nivel microscópico parecen lograr que el material
como un todo sea más resistente de lo que seespera. Explicar cómo se obtiene la respuesta
mecánica de este tipo de materiales a partir de su microestructura no es tarea fácil, pues su
comportamiento es el resultado de la interacción no lineal de sus muchos componentes, por lo
que estos materiales se suelen catalogar como complejos.
Al estudiar propiedades mecánicas de los materiales, uno de los comportamientos a nivel
microscópico de mayor interés es su resistencia a la falla, ya que una fractura que atraviese
el material y lo dañe es al final el efecto combinado de muchas fracturas microscópicas.
Con este objetivo se han construido modelos estad́ısticos que describen cómo se dan las
fracturas entre los elementos que forman el material y que interactúan entre śı de manera local.
Dichos modelos, denominados Modelos Estad́ısticos de Fractura (SMF por sus siglas en inglés),
representan una herramienta útil para establecer, de manera cualitativa, la relación entre la
11
12 Introducción
Figura 1: El bambú es nativo en 5 continentes: Africa, Asia, Sur America, Norte America y Australia,
imagen tomada de [Bam15].
micro y macromecánica del material, pero el camino no es fácil. Las dificultades quedan bien
expuestas en el trabajo de Kun et al. [FKH06]. Por un lado, es importante trabajar con modelos
realistas de los materiales que incluyan una representación detallada de la microestructura del
material, pues tales modelos hacen posible clarificar el efecto de los parámetros microscópicos
del material en su respuesta macroscópica; pero un modelo realista es demasiado pesado
computacionalmente, es decir, los tiempos de computo son considerables ya que el análisis
se basa en comportamientos estad́ısticos obtenidos a partir de miles de simulaciones. Por otro
lado, el estudio de procesos de daño y fractura de materiales desordenados aborda problemas
interesantes también para la f́ısica estad́ıstica como procesos cŕıticos, transiciones de fase, etc.
En este tipo de estudios, el interés se centra en simular sistemas suficientemente grandes, como
para contar con buena estad́ıstica. Esto exige simplificar el modelo a sus componentes esenciales,
pues de otra forma no será posible simularlos en un tiempo de cómputo razonable, con lo que los
resultados numéricos se limitan a encontrar la forma funcional de las dependencias, pero no los
valores exactos de los coeficientes para un material real. El arte consiste, entonces, en encontrar
sistemas estad́ısticos sencillos que den valores cuantitativos de orden general (como exponentes
en leyes de potencia, por ejemplo) que puedan ser comparados con datos experimentales.
Los modelos de fibras o fiber bundle models (FBMs) son ejemplos de dicha simplificación.
La estructura del modelo consiste en un arreglo paralelo de fibras, fijas en los extremos, y que
se someten a cargas paralelas a la fibra. Cada fibra obedece la ley de Hooke y por lo tanto
genera una respuesta elástica lineal ante el esfuerzo aplicado hasta que se rompe por superar
un esfuerzo ĺımite [Pie26]. La distribución de la carga de la fibra rota se redistribuye en las fibras
intactas, lo que puede ocasionar una nueva falla hasta que el sistema falla estructuralmente, es
decir, la avalancha de fallas produce que todas las fibras se rompan. Este tipo de modelos pueden
incluir caracteŕısticas de tipo estad́ıstico como la distribución probabilistica de los umbrales de
ruptura en las fibras [Dan45], lo que conlleva a comportamientos cŕıticos en los que la dinámica
de las avalanchas obedecen a una le de potencias [HH92] con exponentes caracteŕısticos. Este
13
tipo de comportamientos asemejan transiciones de fase continuas o abruptas en la dinámica de
fractura [FKH03, HP07].
Los SMFs, permiten modelar materiales diversos, desde estructuras sintéticas, como la fibra
de carbono y materiales compuestos, hasta tejidos naturales, en los que se tienen interacciones
celulares de diversos tipos. En general, si las fuerzas cohesivas que dan estabilidad estructural
a un material se modelan como enlaces con respuesta elástica, entonces se pueden aplicar
modelos estad́ısticos. La guadua, que es el material de interés en el presente trabajo, incluye
en su estructura interna tejido de fibra y enlaces celulares denominados lamelas, que dan
cohesión al tejido blando denominado parénquima. Entonces, sus comportamientos de falla ante
determinados esfuerzos se pueden estudiar a partir de FBMs. Estudiar la guadua con SMFs
es interesante no solo por los resultados de modelos anteriores, sino porque puede evidenciar
el papel que juega la estructura y el desorden de esa estructura celular sobre sus propiedades
mecánicas. Se cree que la causa de esa resistencia de la guadua reside en su estructura interna
pues, a diferencia de la madera, donde las fibras corren en todas direcciones, los tejidos de
la guadua se ordenan solamente en dirección axial, con haces de fibras fuertes dentro de un
tejido mucho más blando, la parénquima, que es el que se rompe cuando el material falla
[LEM06, LXR02].
Un referente fundamental para el presente trabajo es la investigación desarrollada por
Villalobos en su tesis doctoral [Vil12] para estudiar las fracturas producidas por la contracción
que se generan en procesos de secado. En su tesis, Villalobos construye varios modelos
bidimensionales inspirados en la estructura casi hexagonal que muestra el tejido parénquima
cuando se corta en un plano transversal al eje del tallo. Uno de los modelos trabajados por
Villalobos consiste en un arreglo casi hexagonal de masas puntuales (que representan las células)
unidas con una red de resortes elásticos (que representan las lamelas, o enlaces celulares). La
idea fundamental fue asignar a cada uno de los enlaces un umbral de ruptura, dejar fijas las
posiciones de todas las células en el peŕımetro y contraer la longitud natural de todos los
resortes, para ver cómo los enlaces se romṕıan en avalanchas. El análisis estad́ıstico de la
dinámica de falla encontró dos reǵımenes. Si los umbrales de ruptura son todos similares, se
obtiene un material tenaz, donde muchos enlaces intercelulares se rompen a la vez en una gran
falla, casi al inicio de la contracción, que prácticamente divide el sistema en dos desde el inicio.
Por el contrario, si los umbrales de ruptura se distribuyen sobre un intervalo amplio de valores,
se obtiene un material dúctil, donde los enlaces se destruyen aqúı y allá, generando pequeñas
fracturas que poco a poco se unen hasta atravesar el material y hacerlo fallar. El paso entre los
dos es una transición de fase, donde el parámetro de control es el nivel de desorden estructural
(entendido como el ancho relativo del intervalo de valores que toman los umbrales de falla) y
el parámetro de orden es el tamaño de la avalancha más grande.
Los resultados de Villalobos muestran la eficiencia de los FBMs para estudiar fractura
en guadua y destacan el acierto de modelar los enlaces celulares por medio de resortes con
respuesta elástica lineal y umbrales de ruptura, pero no agotan el estudio de la parénquima.
En efecto, cuando el corte se hace en un plano que contenga al eje del tallo, la parénquima
también muestra una estructura peculiar, pues alterna al azar células largas y cortas en tubos
paralelos que corren en dirección axial (Fig.2.c), todas unidas por medio de lamelas [Lie98].
La pregunta es si este desorden estructural contribuye también a la resistencia del material. El
14 Introducción
objetivo del presente trabajo es determinar si la alternancia aleatoria de células largas y cortas
en el tejido parénquima de la G. angustifolia aumenta la resistencia del material a la falla ante
deformaciones de flexión o tracción.
Figura 2: (a) Se muestra el culmo de gadua señalado por el recuadro rojo. (b) Se muestra una
sección transversal del culmo. Las regiones oscuras corresponden al tejido de fibra y la parte más
clara corresponde al tejido de parénquima, tomado de [Vil12].(c) Vista tridimensional de tejidodel
culmo con haces vasculares rodeados de parénquima, tomado de [Lie98]. (d) Fotograf́ıa del arreglo del
tejido de parénquima en la dirección axial en el que se dstinguen células largas y cortas alternadas
aleatoriamente, tomado de [Lie98].
Para dar respuesta a este interrogante, se construirá un modelo bi-dimensional de elementos
discretos que asemeje la estructura del tejido en el plano radial-axial, representando las células
con rectángulos que se ordenan en tubos paralelos que alternan aleatoriamente células largas y
cortas, como en el tejido real (Fig. 2.c y 2.d ). Al igual que en los modelos anteriores, las células
se unirán con resortes frágiles, que se rompen cuando la fuerza que ejercen supera cierto valor
umbral, escogido de manera aleatoria entre dos valores ĺımite. Luego, el sistema total se someterá
a deformaciones crecientes de tracción o flexión, integrando las ecuaciones de movimiento del
sistema, y se analizará la estad́ıstica de las avalanchas de rupturas y los valores de deformación
para los cuales ocurre una falla estructural. Finalmente, se compararán los datos obtenidos
con diferentes configuraciones de células: unas en la que el tejido se compone de células de dos
tamaños (largas y cortas) y, otras en las que todas las células tienen el mismo tamaño (todas
largas o todas cortas). De esa manera se espera poder determinar si esta alternancia en el
tamaño de las células aumenta o no la resistencia del material ante esfuerzos de tracción axial
y flexión.
El presente documento está organizado de la siguiente manera. En el caṕıtulo 1 se hace
una descripción de la estructura de la guadua, poniendo especial atención a la configuración
microscópica del tejido de parénquima en la dirección axial, y se muestran algunos resultados
experimentales que tienen que ver con las propiedades mecánicas de los tejidos de f́ıbra y
parénquima de la guadua al someterla a pruebas de flexión y tracción. El caṕıtulo 2 inicia
15
presentando las generalidades del método de dinámica molecular , para luego discutir diferentes
modelos de fractura, en especial el modelo de haz de fibras y el modelo de contracción inspirado
en G. angustifolia propuesto por G. Villalobos et. al [Vil12]. El caṕıtulo 3 presenta el modelo
propuesto para representar el tejido parénquima en un corte axial. El caṕıtulo 4 muestra
los resultados de los ensayos computacionales de flexión y tracción para tres combinaciones
diferentes de células: todas largas, todas pequeñas o alternadas al azar. Finalmente, se presentan
las conclusiones y recomendaciones más importantes.
Caṕıtulo 1
La guadua
El bambú Guadua angustifolia ha mostrado ser un material con propiedades mecánicas
excepcionales, como una enorme resistencia a esfuerzos de flexión, tracción y compresión. De
hecho, una viga de guadua laminada de sección cuadrada de 12cm y 120cm de largo se puede
flexionar en su centro hasta una cuarta parte de su longitud, y al soltarla la viga vuelve a su
estado original sin daños estructurales apreciables [Mon10]. Este tipo de ventajas ha motivado
su aplicación en campos como la construcción, el diseño de muebles y la creación de artesańıas.
A sus valiosas propiedades estructurales se suma el hecho de que en muchas regiones del mundo
el bambú se cultiva de manera sencilla, eficiente, a bajo costo y con beneficios ambientales
considerables (Fig. 1).
1.1. Estructura de la guadua
Figura 1.1: De izquierda a derecha un culmo de G. angustifolia al iniciar, uno joven (de 1 a 3 años),
uno adulto (de tres a seis años) y finalmente un culmo seco (mayor de 6 años) [XL02].
16
1.1. Estructura de la guadua 17
La planta de bambú tiene una estructura principal denominada culmo del que nacen
ramificaciones alternas [LXR02] (Fig. 1.2). En la dirección longitudinal, el culmo se desarrolla
en secciones ciĺındricas compuestas por nudos y entrenados, cuyos diámetros decrecen con la
altura. Sus propiedades mecánicas dependen de la sección de la vara que se considere. Por
ejemplo, la sección comprendida aproximadamente entre 0m y 4m de altura con diámetros de
alrededor de 18cm, conocida como cepa, es ideal para estructuras en las que los esfuerzos de
compresión son altos, la sección entre 4m y 10m de altura , conocida como basa, es la de mayor
uso en muebles, pues mantienen diámetros aproximadamente constantes a lo largo del culmo,
y la última sección, denominada sobrebasa, que constituye la sección más alta y más delgada,
es la de de mayor resistencia a la flexión. En términos de la edad, un culmo de guadua se
considera joven si tiene entre 1 y 3 años de edad y se considera maduro a los 6 años de edad. El
corte y uso del culmo es ideal entre los 2 y los 6 años [LXR02]. Después de los 6 años el culmo
se denomina sobremaduro y empieza a perder el contenido de humedad, lo que eventualmente
produce fracturas en su superficie exterior.
Figura 1.2: Dependiendo de la zona de la guadua, el material presenta caracteŕısticas espećıficas. Se
reconocen tres secciones principales: sobrepasa, basa y cepa. Tomado de [Bam].
Existen cerca de 1250 especies de Bambú, que se diferencian fundamentalmente en su
estructura interna. Sin embargo, las diferencias son consideradas mı́nimas en comparación con
la gran heterogeneidad estructural de las 20 000 especies de madera [Lie98]. Internamente, el
culmo de bambú presenta una estructura de haces vasculares rodeados de fibra, insertos dentro
del tejido blando denominado parénquima (Fig. 2). El culmo está formado por entre 47 % y 52 %
en peso de parénquima; entre 40 % y 42 % de fibra y, entre un 8 % y un 9 % de tejido conductivo
[Lie98, OVC07, LXR02] citados en [OSJ09]. Esta variedad de componentes permiten describir
la guadua como un material compuesto, reforzado con fibras, en las que el tejido de parénquima
tiene la función de distribuir las cargas axiales en el tejido de fibra [ZPST10].
18 Caṕıtulo 1. La guadua
Figura 1.3: Estructura anatómica del culmo de la Guadua angustifolia en un corte transversal: a) Se
indica la corteza y la parénquima en las regiones c y p . b) Detalle del haz vascular rodeado de fibra
en la región f, tomado de [XL02].
El arreglo tanto del tejido blando de parénquima como de la fibra y el tejido conductivo es
axial. Sin embargo, dependiendo de la región del culmo que se estudie se tienen caracteŕısticas
microestructurales diferentes. Por ejemplo, para las guaduas Moso bambú y G. angustifolia,
la distribución del tejido de fibra y de parénquima en la sección transversal del culmo no es
homogénea, sino que la concentración de tejido de fibra aumenta de la parte interior a la parte
exterior, aumentando también la densidad del material [DG14, OVC07, LAS10]. Para la Moso
bambú, la densidad promedio del tejido de fibra en las regiones más densas puede llegar a
800Kg/m3 y en regiones internas del culmo puede ser de 500Kg/m3; mientras que para el
tejido de parénquima la densidad promedio es de 350Kg/m3, [SA15]. En la G. angustifolia
se encuentra que el porcentaje de tejido de parénquima incrementa de las regiones de la cepa
externa a la interna con valores máximos y mı́nimos de 60 % y 35 %, respectivamente. Este
porcentaje también vaŕıa con la altura, pues la fracción de tejido parénquima es mayor en las
regiones de la cepa [OVC07].
1.2. Propiedades mecánicas de la guadua
Desde el punto de vista mecánico, la madera de la guadua presenta dos componentes: la
parénquima, que es un tejido blando, y estructuras aisladas de tejidos de fibra alrededor de
tejidos de conducción, que en conjunto podemos llamar las fibras, y que son más duras y
resistentes. Se ha encontrado que la parénquima tiene un el módulo de elasticidad aproximado
de 5,8GPa, lo que representa un 33 % de los valores obtenidos para la fibra, mientras que el
módulo de ruptura de sus paredes celulares toma valores aproximados de 230MPa, un 63 % de
los valores correspondientes a las fibras [YY07]. Cerca de lapared interna del culmo hay pocas
fibras insertas en la parénquima, y a medida que nos acercamos a la pared externa la densidad de
fibras va aumentando. Como consecuencia, la resistencia de la madera completa va aumentando
a medida que nos alejamos del eje del culmo. Pruebas de nano-indentación muestran que el
modulo de Young total crece desde 7,0GPa cerca de la superficie interna del culmo hasta
13GPa cerca de su borde externo [TRA+11]. Experimentos de micro-tensión también revelan
que el esfuerzo máximo de tensión alcanza valores aproximados de 200MPa para las regiones
internas del culmo y 700MPa para las regiones más externas [TRA+11]. Trabajos más recientes
1.2. Propiedades mecánicas de la guadua 19
indican que el módulo de Young y el esfuerzo de ruptura crecen linealmente en proporción con
la densidad de la madera como un todo [DG14] (Tabla. 1.1).
Distancia normalizada Densidad del tejido de fibra MOE axial MOR axial
0,10 580Kg/cm3 5GPa 40MPa
0,20 600Kg/cm3 7GPa 90MPa
0,80 900Kg/cm3 20Gpa 230MPa
Tabla 1.1: Comportamiento del MOE axial y MOR axial para diferentes regiones de la pared del
culmmo de la guadua Moso bambú. Los valores son tomados de [DG14]
La edad de la planta tiene incidencia en las propiedades mecánicas del culmo. El agua
que transporta el culmo desde la ráız viene cargada con sedimentos de śılice, que se van
acumulando en los espacios intersticiales entre las células de la parénquima. Luego de mucho
tiempo, se acumula tanto śılice que la parénquima ya no puede retener agua. El tejido, además,
se vuelve mucho más duro, razón por la cual se necesitan herramientas especiales (usualmente
recubiertas de diamante artificial) para cortarlas en la producción de papel. Una vara de más
de 6 años ya no puede retener un gran porcentaje de agua, y el culmo se torna seco (Fig. 1.1.c).
Investigaciones demuestran que este proceso de secado de la guadua en pie (antes de cosecharla)
induce fracturas debido al colapso de la estructura celular [OGT04, OGP05, SMO09] citados en
[Vil12]. La edad también afecta las propiedades mecánicas del material. En probetas de guadua
(compuestas de fibra y parénquima) sometidas a ensayos de flexión (Fig. 1.4), se encuentra que
las varas sobremaduras (mayores de 6 años) tienen mejores comportamientos para los módulos
de ruptura y elasticidad que varas más jóvenes (Fig. 1.5) [HAG07]. Es interesante notar que
estos cambios con la edad ocurren en el tejido parénquima, pero no en las fibras, que además
tampoco presentan diferencias significativas relacionadas con la sección del culmo a la que
pertenecen. En efecto, al caracterizar haces de fibra extráıdos de G. angustifolia, se encuentra
que no existe una diferencia significativa para los módulos de elasticidad en ensayos de tensión
axial respecto a la edad de la guadua (joven, madura o sobremadura) o respecto a la zona del
culmo del que se extrae la fibra (cercana al eje o cercana a la superficie externa del culmo)
[LEM06].
Existe evidencia experimental [ZPST10], de que cada uno de los dos tejidos (fibras y
parénquima) aportan al módulo de elasticidad axial (MOE) y al módulo de rotura (MOR) del
material en proporción a la fracción de volumen que ocupa cada tejido. De la extrapolación lineal
de los datos experimentales, se obtienen valores para ambos módulos si el material estuviese
formado de solo fibras o solo parénquima, como se ve en la Tabla. 1.2. Los valores de ensayos
realizados sobre haces de fibras resultan ser similares a los de la extrapolación, lo que valida
dichos resultados extrapolados. De acuerdo con ellos, el módulo de elasticidad axial del tejido
parénquima puro se situa alrededor de los 0.22GPa, y su módulo de rotura, alrededor de los
19MPa.
Aunque el buen comportamiento de la guadua ante esfuerzos de tracción axial y flexión
se atribuyen principalmente a su estructura fibrada, el tejido de parénquima presenta un
diseño interesante que también puede ayudar a sus propiedades mecánicas. Como ya se
20 Caṕıtulo 1. La guadua
mencionó, las células de parénquima se organizan en columnas orientadas en dirección axial,
alternando aleatoriamente células largas con células cortas (Fig. 2). Esta configuración sugiere
la existencia de un beneficio estructural, similar al que se obtiene en otros materiales en los que
la organización o diseño a nivel microscópico genera efectos en la macro mecánica del material
[BB05, SWC12]. El objetivo del presente trabajo es determinar si este diseño estructural del
tejido de parénquima efectivamente puede mejorar la resistencia a la falla del tejido ante
esfuerzos de tracción axial y flexión. Para ello, construiremos un modelo bidimensional de
elementos discretos que reproduzca ese diseño estructural y que permita estudiar la influencia
de dicho diseño en la resistencia del material a la falla.
(a)
(b)
Figura 1.4: (a) esquema de probeta para prueba de flexión y sitio de extracción de la guadua. (b)
ensayo de flexión de una probeta de guadua. Tomados de [HAG07].
(a) (b)
Figura 1.5: (a) Modulo de elasticidad MOE ( N
mm2
), promedio por edad de la guadua. (b) Modulo de
rotura MOR ( N
mm2
), promedio por edad de la guadua. El intervalo de confianza es del 95 %. Tomado
de [HAG07].
1.2. Propiedades mecánicas de la guadua 21
MOE axial MOR axial
Bloques de Moso Bambú Ef = 40,4GPa σf = 581,7MPa
Bloques de Moso Bambú Ep = 0,22GPa σp = 19,0MPa
Haces de fibra Ehf = 33,96GPa σhf = 482,2MPa
Tabla 1.2: Modulos de elasticidad y de ruptura para fibra y parénquima en la guadua Moso bambú.
Los valores son tomados de [ZPST10].
Caṕıtulo 2
Marco Teórico
En el caṕıtulo anterior se presentaron las principales caracteŕısticas del bambú. Por un lado,
se hizo una aproximación a la anatomı́a y caracteŕısticas de la varas de guadua en general, y
se logró identificar el tipo de tejidos que la componen: tejido de fibra y tejido conductivo, de
una parte, y parénquima, de otra; aśı como la forma en que éstos se distribuyen a través del
culmo. También se habló de su respuesta mecánica ante esfuerzos de tracción y flexión, y de
las amplias posibilidades de aplicación de la guadua; se listaron las propiedades mecánicas del
material , un campo de investigación activo en el que cada vez es más importante considerar
la incidencia de la estructura interna en su comportamiento mecánico. El objetivo de esta
trabajo es, precisamente, establecer mediante modelos numéricos del tejido la posible relación
entre la estructura microscópica de la parénquima en corte axial y su respuesta microscópica
ante esfuerzos de tracción y de flexión. En este caṕıtulo de presentarán los métodos numéricos
utilizados para simular el tejido.
2.1. Dinámica molecular
Cuando un sistema f́ısico está formado por multitud de cuerpos individuales que se ejercen
fuerzas entre si, como los granos en un suelo o las células en un tejido, los Métodos de elementos
discretos (DEM por sus siglas en inglés), son una excelente opción para su simulación. En estos
métodos los cuerpos se simulan como objetos geométricos (generalmente esferas o poliedros) que
pueden trasladarse y rodar, y las fuerzas y torques entre cada par de cuerpos se calculan a partir
de sus posiciones y velocidades relativas. Existen diferentes formas de calcular las interacciones,
lo que da lugar a diferentes variantes. Una de las más usadas es el método de dinámica molecular
(MD por sus siglas en inglés). Este método integra las ecuaciones de movimiento de Newton a
pasos de tiempo ∆t. En cada paso de tiempo se calculan las fuerzas sobre cada part́ıcula y se
integran para hallar sus nuevas posiciones y velocidades. Estas fuerzas pueden ser de dos tipos:
interacciones entre las part́ıculas, que pueden incluir potenciales de interacción a distancia o
fuerzas por contacto directo en choques, y fuerzas externas, como la gravedad o la interacción
viscosa con el medio. Las fuerzas de contacto se suelen calcular permitiendo quelas figuras
geométricas se intercepten parcialmente, midiendo la magnitud (profundidad, área o volumen)
de dicha intersección y calculando la fuerza a partir de estas medidas (Fig. 2.1). Entonces, la
22
2.1. Dinámica molecular 23
fuerza de contacto solo depende de las posiciones relativas entre las part́ıculas.
Figura 2.1: La fuerza de contacto se puede definir en términos de la longitud de penetración h. Las
posiciones y velocidades de cada part́ıcula en un instante posterior t+∆t, se obtiene con las posiciones,
velocidades y fuerzas en un instante inicial t
Si las part́ıculas son esferas o discos, la fuerza normal de contacto ~Fc entre dos de ellas
se puede modelar como proporcional una potencia de la distancia de interpenetración h,
Fc = Kch
3/2 o Fc = Kch
1, respectivamente. La constante de proporcionalidad Kc es función de
la rigidez del material del que están hechas las part́ıculas. Como esta fuerza es función de una
distancia, es fácil asignarle una enerǵıa potencial.
Sin embargo, en las situaciones reales, las part́ıculas (arenas, rocas, granos de arroz, etc.)
no son esféricas. Una alternativa que permite generar part́ıculas de formas aleatorias que se
acerquen a las formas reales, es usar los denominados poĺıgonos de Voronöı. El método consiste
en ubicar al azar un conjunto de puntos Pi, denominados puntos generadores, en una región
finita del espacio. Luego, para cada punto generador, se define el poĺıgono de Voronöı como el
lugar geométrico formado por todos los puntos del espacio que están mas cerca de Pi que de
cualquier otro punto generador (Fig. 2.2). Las formas aśı obtenidas son más cercanas a las de
los granos reales, pero al precio de complicar enormemente el cálculo numérico de las fuerzas
de contacto. Algunas de las expresiones más usadas consideran fuerzas que son proporcionales,
por ejemplo, al rea o (volumen) de intersección, ejercida en el punto medio de la linea que une
los puntos de intersección de los peŕımetros de los dos poĺıgonos, y perpendicular a ella. Estas
expresiones presentan inestabilidades numéricas y, en general, no permiten asignar de manera
sencilla una enerǵıa potencial elástica al contacto.
Una estrategia que concilia el uso de formas más reales con un tiempo de computo razonable,
es usar esferopoĺıgonos (Fig. 2.3), que son esencialmente poĺıgonos a los que se le han redondeado
sus esquinas mediante dos operaciones, llamadas dilatación y contracción [AM08]. Considere un
poĺıgono A y un disco B . La dilatación de A por B (indicada como A⊕B ) es el conjunto que
se obtiene al barrer con el disco alrededor del perfil del poĺıgono. Si el poĺıgono es un cuadrado,
el resultado es un cuadrado más grande, pero con las esquinas redondeadas (Fig. 2.3a). La
erosión también es una operación binaria entre conjuntos. Consiste básicamente en contraer
todos los lados del poĺıgono para convertirlo en uno más pequeño (Fig. 2.3b). La operación se
24 Caṕıtulo 2. Marco Teórico
Figura 2.2: (Izquierda) Poĺıgonos de Voronöı construidos a partir de puntos generadores, tomado de
[Gal10]. (Derecha) Sobrelamamiento de dos poĺıgonos en interacción por choque
(a) Dilatación (b) Erosión (c) Apertura
Figura 2.3: Para optimizar los cálculos de fuerzas de contacto se modelan las células con rectángulos
redondeados, figura tomada de [Gal10]
denota como A	B. Finalmente, la combinación de las dos operaciones, denominada apertura
y notada como A ◦ B = (A 	 B) ⊕ B, da como resultado una figura en la que se cambian
los vértices por esquinas redondeadas (Fig. 2.3c). Aśı es mucho mas sencillo caracterizar las
fuerzas de contacto para los choques entre dos part́ıculas. En efecto, si la constante elástica de
contacto Kc es suficientemente alta, los únicos choques entre dos esferopoĺıgonos serán vértice
contra vértice y vértice contra lado. Los dos corresponden a casos particulares de la fuerza de
contacto entre dos cilindros, que es conocida anaĺıticamente, y que puede aproximarse como
lineal. Esta fuerza corresponde ademas a una enerǵıa potencial cuadrática con la distancia de
interpenetración.
Los métodos de DM son eficientes para modelar sistemas de muchas part́ıculas. Aspectos
como interacciones f́ısicas, formas o la cantidad de part́ıculas se definen de forma sencilla en los
modelos y se pueden simular en tiempos de computo razonables. Estas ventajas han logrado
extender su aplicación al estudio de procesos de fractura en medios continuos, debido a la
2.2. Modelos estad́ısticos de fractura 25
importancia que se le confiere a la influencia de la estructura interna de los materiales en sus
comportamientos de fractura.
2.2. Modelos estad́ısticos de fractura
El estudio de procesos de fractura en materiales continuos es también uno de los campos
de aplicación de los DEM, en lo que se conoce como Modelos Estad́ısticos de Fractura o
(SMF por sus siglas en inglés). El material de estudio se representa formado por elementos
individuales (por ejemplo las células de un tejido) que se modelan como elementos discretos. Los
elementos se encuentran unidos por enlaces, que representan las fibras o moléculas que unen a
los elementos del material. Estos enlaces se suelen modelar como resortes o vigas, que se rompen
cuando alcanzan cierto umbral de esfuerzo o de deformación. Por lo tanto, en estos modelos
los elementos discretos no se rompen, sino que son los enlaces que estructuran el material los
que son propensos a fallar cuando sus esfuerzos o, equivalentemente, sus deformaciones superan
cierto umbral caracteŕıstico. Adicionalmente, cuando un enlace falla, los esfuerzos que soportaba
se ven redistribuidos en el resto de la estructura, unas veces de manera local, cuando la carga
del elemento roto se distribuye entre sus enlaces vecinos, y otras de manera global, cuando la
carga del elemento roto se redistribuye en todos los enlaces intactos. Este reordenamiento de las
cargas suele generar a su vez otras fracturas, casi de manera instantánea, y as sucesivamente,
en lo que se conoce como una avalancha. En este tipo de aproximaciones, una fractura que
atraviese el material y lo haga fallar es al final el efecto combinado de muchas fracturas
microscópicas locales. El estudio de estas avalanchas y de las fracturas que van generando
permiten caracterizar finalmente la resistencia del material como un todo, y establecer relaciones
entre los parmetros de resistencia microscópicos y la respuesta macroscópica del material.
Uno de los resultados más importantes, en lo que tiene que ver con los comportamientos de
fractura, es que el desorden estructural mejora la resistencia a la falla del material, entendida
como el esfuerzo máximo o la deformación máxima que puede alcanzar antes de romperse en
dos pedazos. El desorden estructural se puede implantar de diferentes maneras: por ejemplo,
usando elementos con formas y tamaños aleatorios [Vil12] (como los poĺıgonos de Voronöı) o
usando enlaces con umbrales [FKH06, Vil12] o constantes elásticas [SWC12, LdA85, PC03]
escogidos aleatoriamente a partir de una distribución probabiĺıstica. Este desorden estructural
aumenta la resistencia a la falla a través de diversos mecanismos. Por ejemplo, en el caso de un
material que se contrae por secado [Vil12], si todos los umbrales de ruptura de los enlaces son
iguales (poco desorden estructural), se obtiene que para un cierto valor de esfuerzo se produce
una fractura que crece por toda la estructura, y el material falla muy prontamente, como si
fuera frágil. Por el contrario, si los umbrales están distribuidos aleatoriamente sobre un intervalo
relativamente amplio de valores, el material comienza a fallar locamente aqúı y allá hasta que
en determinado momento falla completamente por la unión de muchas fallas locales. El hecho
más importante es que para este segundo caso la falla total se presenta para un valor de esfuerzo
mucho mayor que para el primer caso, lo que se puede asociar a un comportamiento dúctil.En efecto, las pequeñas fracturas aqúı y allá generan posibilidades de movimientos locales que
redistribuyen las cargas sin llegar a romper los enlaces, mientras que para el caso frágil todos
26 Caṕıtulo 2. Marco Teórico
los enlaces fallan al tiempo, sin darle oportunidad al sistema de reacomodarse. El cambio de un
comportamiento al otro se suele caracterizar como una transición de fase en donde el parámetro
de control es una medida del desorden estructural.
Este tipo de modelos estad́ısticos tienen dos ventajas importantes. Por un lado, es posible
estudiar el comportamiento mecánico del material de forma local, de manera que todas las
propiedades o caracteŕısticas microscópicas puedan ser implementadas en el modelo, al menos de
manera simplificada. Por otro lado, es sencillo incluir en los modelos caracteŕısticas estructurales
globales como la distribución de los umbrales de ruptura, lo que permite identificar la influencia
que puede tener la micro-mecánica del material en el comportamiento mecánico global. A
continuación describiremos dos ejemplos, a saber: el Fiber Bundle Model (FBM), que fue el
primer modelo SMF propuesto, y el Cell Model, propuesto por G. Villalobos, F. Kun y J.D.
Muñoz [Vil12], que es el antecedente directo del presente trabajo.
2.2.1. Modelo de haces de fibra
Un modelo de hace de fibra (FBM por sus siglas en inglés) consiste en un conjunto de fibras
paralelas entre si fijas a dos placas y sometidas a una carga paralela a la fibra (Fig. 2.4a, a).
Las primeras descripciones de este sistema fueron realizadas por Pierce en 1926 para estudiar
la resistencia de fibras de algodón [Pie26], pero el primer modelo estad́ıstico fue propuesto
por Daniels [Dan45]. En su modelo, todas las fibras obedecen la ley de Hooke con una misma
constante del resorte, pero cada una de ellas tiene un valor de esfuerzo ĺımite diferente, obtenido
a partir de una distribución de probabilidad conocida. Inicialmente cada una de las fibras
soporta una fracción de carga dada por Fi =
F
N
, donde F es la carga total y N es el número
de de fibras del sistema que permanecen intactas. Las placas se comienzan a alejar hasta que
la primera fibra se rompe. En ese momento, las fibras intactas se redistribuyen la carga, lo que
puede generar que otras fibras se rompan, y aśı sucesivamente, configurando una avalancha de
rupturas. El tamaño de la avalancha se suele definir como el número de fibras que se rompen.
Cuando la avalancha termina, la carga vuelve a aumentar y se produce una nueva avalancha,
y aśı sucesivamente, hasta que el material falla completamente. Este tipo de modelos en el
que la carga que soportaba la fibra rota se distribuye equitativamente en las fibras restantes
se denomina distribución global de carga (GLS por sus siglas en inglés). Se ha demostrado
anaĺıticamente que en este tipo de modelos el número P (∆) de avalanchas de tamaño ∆ se
distribuye como una ley de potencias, P (∆) ∼ ∆2,5, cuando el número de fibras N → ∞, lo
que es evidencia de un comportamiento cŕıtico [HH92]. En otro tipo de modelos, conocidos
como distribución local de carga (LLS por sus siglas en inglés) presentan otros exponentes
caracteŕısticos y, eventualmente, comportamientos que asemejan transiciones de fase continuas
o abruptas en los proceso de falla [FKH03, HP07, SZV99].
Trabajos posteriores implementan propiedades f́ısicas adicionales en las fibras con el fin de
aproximarse a modelos más realistas. Por ejemplo,aunque el modelo básico de fibras supone un
umbral de ruptura constante durante todo el proceso de falla [Dan45], se puede considerar una
dependencia temporal, es decir, se implementa un proceso de fatiga en el material disminuyendo
el umbral de ruptura de las fibras en el tiempo [Col56]. La fatiga o desgaste también se puede
implementar reduciendo el modulo de Young de las fibras intactas en cada proceso consecutivo
2.2. Modelos estad́ısticos de fractura 27
(a)
(b)
Figura 2.4: (a) Esquema de fibras paralelas en un FBM, tomado de [SPC10]. (b) Distribución de
avalanchas t́ıpicas bajo carga global en el haz. Tomado de [FKH06].
de aumento de carga, en lo que se denomina un FBM de daño continuo (CDFBM por sus siglas
en inglés) [FKH06].
FBM en modelo de Guadua angustifolia
Figura 2.5: (Izquierda) Foto de la sección transversal tejido de parénquima [VLV+03] . (Centro)
Ejemplo de una geometŕıa aleatoria inicial. (Derecha) Umbrales de ruptura uniformemente distribuidos
con centro en C. La cantidad de desorden es controlada variando el valor de W (semiancho de la
distribución). Tomado de [Vil12].
Uno de los trabajos de referencia en simulación de fracturas en bambú para la presente
investigación fue desarrollado por Villalobos en el 2012 [Vil12], quien en su tesis doctoral
28 Caṕıtulo 2. Marco Teórico
propuso varios modelos de fractura inspirados en la geometŕıa del tejido parénquima de la
guadua cuando se corta en un plano perpendicular al tallo. En este plano, la parénquima
aparece como un arreglo de celdas aproximadamente hexagonales (Fig. 2.5, izquierda), lo que
motivó que todos sus modelos tuviesen un diseño aproximadamente hexagonal. En uno de
ellos, el tejido se modela como una red triangular ligeramente desordenada de masas unidas
por resortes. Cada masa está unida a seis vecinas por medio de resortes frágiles y cada resorte
se rompe cuando su esfuerzo supera cierto umbral. Los umbrales se distribuyen uniformemente
entre los valores C−W y C+W (Fig. 2.5). Los nodos en la periferia se mantienen en posiciones
fijas, mientras que el sistema se contrae uniformemente reduciendo gradualmente la longitud
natural de todos los resortes por un factor α, hasta que uno de ellos se rompe. En ese momento,
los esfuerzos se redistribuyen de manera natural al integrar las ecuaciones de movimiento, lo
que puede provocar que otro resorte se rompa, y aśı sucesivamente, creando una avalancha de
fracturas que eventualmente cesa. Entonces, se reinicia la contracción hasta la siguiente falla,
y aśı sucesivamente, hasta que el material completo se rompe en dos.
Figura 2.6: (a) Distribución de tamaños de avalancha para diferentes valores de desorden estructural
ξ.(b) Parámetro de orden 〈∆max〉 de la transición, como función del desorden estructural ξ. (c)
Suceptibilidad
〈
m2
m1
〉
como función de ξ
Para los procesos de fractura, Villalobos encontró que el sistema se comporta de manera
diferente dependiendo de la razón ξ = W/C, que representa el grado de desorden estructural.
Para valores bajos de desorden estructural (es decir, cuando todos los umbrales de ruptura
son prácticamente iguales), se produce una gran avalancha, casi al comienzo del proceso, que
rompe el material en dos. En cambio, cuando el desorden estructural es alto, las fracturas
suceden aqúı y allá, produciendo pequeñas fisuras que eventualmente se unen en una gran
fisura que atraviesa el material. El ĺımite entre ambos comportamientos de fractura aparece
alrededor de ξc = 0,4 , punto en el que también la distribución de avalanchas cumple una ley
de potencias, P (∆) ∼ ∆2,6. El parámetro de orden de esta transición de fase, entre reǵımenes
de falla, resulta ser el tamaño de la avalancha máxima menos uno (que llamaremos ∆max), y
la susceptibilidad viene dada por el promedio de la razón m2/m1 entre el segundo y el primer
momento de la distribución de avalanchas (con mk = Σ∆
k − ∆kmax, y donde la suma recorre
todas las avalanchas). El colapso de las curvas de estos valores en función de ξ para diferentes
2.2. Modelos estad́ısticos de fractura 29
tamaños laterales L del sistema arrojan que estas dos cantidades escalan como
〈∆max〉 = Lβ/νf((ξ − ξc)L1/ν) , (2.1)
con valores de ξc = 0, 4 , β/ν = 1, 4 y 1/ν = 1, 0; y como〈
m2
m1
〉
= Lγ/νg((ξ − ξc)L1/ν) , (2.2)
con ξ = 0, 4, γ/ν = 1, 0 y 1/ν = 0, 9, muy cercanos a los anteriores. Estos resultados sugieren
que esta fractura por secado es diferente a las producidas por tracción o por colisión.
Para el presente trabajo,este modelo constituye un insumo fundamental en lo que tiene
que ver con la descripción del modelo de fractura para el tejido de parénquima de la guadua.
Espećıficamente, se implementará de la misma manera la distribución de umbrales de ruptura,
y el valor caracteŕısitico de desorden ξ = 0,4, para el cual el modelo presenta la transición
entre los comportamientos de fractura, representa un valor de referencia para los eventuales
resultados del nuevo modelo.
Caṕıtulo 3
Modelo del tejido parénquima en
dirección longitudinal
En los caṕıtulos anteriores hemos discutido con variados ejemplos cómo la estructura
microscópica del material influye en su resistencia global a la falla. En especial, se han descrito
diferentes estudios que muestran que el desorden estructural hace que el material como un todo
sea más resistente de lo esperado. Ese ha sido también el caso de los modelos inspirados en la
estructura del tejido parénquima de la G. angustifolia cuando se corta en un plano perpendicular
al eje del culmo [Vil12], pero el problema no termina alĺı. En efecto, el mismo tejido parénquima
también presenta una estructura particular cuando se corta con un plano axial, es decir que
contiene al eje del culmo. En efecto, en este caso se observa que las células se organizan en
tubos paralelos pegados entre śı, pero en los que se alternan células largas y cortas de manera
aleatoria (Fig. 3.1.a) [Lie98]. El objetivo del presente trabajo es determinar si esta aleatoriedad
mejora la resistencia del material, para lo cual construiremos un modelo de elementos discretos
que reproduzca la estructura del tejido. En lo que sigue, se describirá el modelo en detalle, aśı
como el proceso de simulación con el que se estudiará su resistencia.
3.1. Construcción del tejido de Parénquima
La configuración del tejido de parénquima en dirección axial consiste, pues, de un conjunto
de tubos en los que se alternan aleatoriamente células largas y cortas, unidas por lamelas. Para
modelarlo, se construyó un arreglo bidimensional de células rectangulares en una disposición
como la que se muestra en la Fig. 3.1. Las células se modelan como esferopoĺıgonos rectangulares
que, cuando chocan, se repelen por fuerzas de contacto elásticas, y que se encuentran también
unidas por resortes que representan los enlaces celulares establecidos por las lamelas. Cada
célula se enlaza con las demás a través de un conjunto de resortes, que salen aproximadamente
de sus vértices: cuatro se usan para conectarla con las células siguiente y anterior de su mismo
tubo, y los demás resortes se usan para conectarla con las células de los tubos vecinos. (Fig.
3.1.c). La condición para enlazar células horizontalmente es que el vértice de determinada célula
tenga dónde conectarse paralelamente a la célula del lado. Eventualmente puede suceder que
un vértice de una célula no se conecte a la derecha o a la izquierda, debido a las posiciones
30
3.2. Fuerzas y torques 31
relativas de los vértices de células adyacentes. Cada enlace tiene un umbral de ruptura (un
esfuerzo máximo) obtenido a partir de una distribución de probabilidad uniforme sobre un
intervalo que puede ser ancho o estrecho, como veremos más adelante. Finalmente, se agrega
también una fuerza viscosa proporcional a la velocidad que ayuda a disipar las oscilaciones que
se inducen al deformar el material.
(a) (b) (c)
Figura 3.1: (a) Imágen real del tejido de parénquima en la dirección axial, escala en µm, tomado de
[Lie98].(b) Modelo de arreglos tubulares del tejido de parénquima con células largas y cortas alternadas
aleatoriamente. (c) Union de las células mediante resortes que implementan los enlaces celulares.
3.2. Fuerzas y torques
Sobre las células actúan tres tipos de fuerzas: fuerzas de contacto, fuerzas de los enlaces y
una fuerza viscosa. Estas fuerzas, junto con los torques asociados a las dos primeras, determinan
el movimiento de cada célula.
3.2.1. Fuerzas de contacto
Las células en el modelo no son completamente ŕıgidas. Cuando dos células chocan se
produce un sobrelapamiento que permite definir la fuerza de contacto. Para nuestro modelo
solamente incluimos fuerzas de contacto normales, dadas por
~Fijc = Kchijn̂ , (3.1)
~Fjic = −~Fijc , (3.2)
donde n̂ es un vector normal a la linea de contacto (Fig. 3.2).
32 Caṕıtulo 3. Modelo del tejido parénquima en dirección longitudinal
(a) Contacto vértice-vértice (b) Contacto lado-vértice
Figura 3.2: Dado un par de células del tejido se dan dos tipos de choque: vértice-vértice y lado-vertice,
h es la distancia de sobrelapamiento y el vector director de la fuerza de contacto es normal a la linea
de contacto.
3.2.2. Fuerzas de los resortes
Los enlaces entre células se representan con resortes, que actúan según la Ley de Hooke,
~Fije = −Ke∆lê , (3.3)
~Fjie = −~Fije , (3.4)
donde ∆l = l− l0 es la longitud de elongación o compresión del resorte y ê es un vector unitario
paralelo al resorte.
En nuestro modelo, todos los resortes tienen la misma constante elástica, y sus longitudes
de reposo son tales que ninguno de ellos hace fuerza en la configuración inicial (es decir, cuando
la deformación del material es nula).
3.2.3. Fuerzas viscosas
Con el fin de atenuar las oscilaciones inducidas por los pasos sucesivos de deformación, se
añade sobre cada célula una fuerza viscosa proporcional a su velocidad,
~Fvis = −c~v , (3.5)
donde c es una constante. La constante c se escogió de tal manera que fuera suficientemente
grande para que las oscilaciones se amortiguaran rápidamente, pero no demasiado como para
alterar el proceso de fractura. Esto se verificó reduciendo la constante a la mitad y viendo que
los resultados permanećıan inalterados.
3.3. Algoritmo de integración 33
3.2.4. Torques
Para el movimiento de rotación de cada célula respecto a su centro, el torque resulta ser
I~α = ~rc × ~Fc + ~re × ~Fe , (3.6)
donde ~re y ~rc ubican los puntos de aplicación de la fuerza elástica de los resortes y la fuerza de
contacto sobre la célula en un choque.
3.3. Algoritmo de integración
Una vez definida la ecuación de movimiento para traslación y rotación de las células, se
resuelven para intervalos de tiempo constantes usando un algoritmo de integración, que a partir
de las posiciones y las velocidades en un tiempo t genere las nuevas posiciones y velocidades
en un tiempo (t + δt). Para seleccionar el algoritmo de integración mas adecuado se usan
criterios como el orden del error de truncamiento del algoritmo, el tiempo de computo y el
tamaño y caracteŕısticas del sistema. Un criterio sencillo consiste en simular un sistema aislado
y no disipativo, en el que la enerǵıa debe conservarse. El algoritmo de integración tiene un
error numérico que, acumulado en cada paso de integración, termina aumentando el valor de la
enerǵıa total del sistema. Este ensayo fija aproximadamente le tiempo máximo de integración.
En su trabajo de tesis te maestŕıa [Wil07] probó diferentes algoritmos con el fin de identificar
cual produćıa menos aumento de enerǵıa por errores numéricos, espećıficamente comparó los
algoritmos Verlet, Leapfrog, Velocidad Verlet-optimizado y algoritmo predictor corrector de
quinto orden. El algoritmo que presenta menor variación en la enerǵıa mecánica en todo el
tiempo de simulación es el Velocidad Verlet Optimzado [IPO02], y es el que se usará en este
trabajo, tanto para la traslución como para la rotación.
Sean ~R = y ~V = ~V (t) la posición y la velocidad, respectivamente, de una part́ıcula i en un
instante de tiempo t. El algoritmo de Velocidad Verlet optimizado consiste de cuatro pasos, a
saber:
1. Se calcula una primera estimación de la posición,
~R1 = ~R(t) + ~V (t)ξδt , (3.7)
donde ξ ' 0, 193183325037836 es el valor óptimo de un parámetro libre que permite
optimizar el algoritmo minimizando el error total de truncamiento.
2. Se calcula la fuerza en la nueva posición ~R1 (denotadapor ~F [~R1 ]), y con ella se obtiene
una primera estimación de la velocidad,
~V1 = ~V (t) +
1
m
~F [~R1]
δt
2
. (3.8)
34 Caṕıtulo 3. Modelo del tejido parénquima en dirección longitudinal
3. Se calcula una nueva posición auxiliar, ~R2, como
~R2 = ~R1 + ~V1(1− 2ξ)δt . (3.9)
4. La nueva posición, ~R2, genera una nueva fuerza, ~F [~R2], en esa posición. Finalmente, la
nueva posición, ~R(t+ δt), y la nueva velocidad, ~V (t+ δt), estarán dadas por
~R(t+ δt) = ~R2 + ~V (t+ δt)ξδt , (3.10)
~V (t+ δt) = ~V1 +
1
m
~F [~R2]
δt
2
. (3.11)
3.4. El proceso de falla
Cada enlace entre células, representado por un resorte, tiene asignado un umbral de ruptura,
que corresponde a la fuerza máxima que pueden soportar al estar elongado. La forma de
definir los umbrales en toda la extensión del tejido es uno de los procesos más importantes,
ya que permite implementar un desorden estructural en términos de una distribución de los
umbrales. El umbral σi para el enlace i-ésimo se escoge al azar en un intervalo (σmin, σmax),
que se puede caracterizar con dos valores: su valor central C = (σmax + σmin)/2 y su sem-
iancho W = (σmax − σmin)/2 (Fig.3.3). Si la distribución es amplia, los umbrales de ruptura
tienen valores muy diferentes, y el conjunto de enlaces es desordenado. Por el contrario, si la
distribución es angosta, los umbrales tienen valores similares a lo largo del tejido, y el conjunto
de enlaces es ordenado . Dicho desorden se puede caracterizar con el parámetro adimensional
ξ = W
C
, que ha mostrado caracterizar bien los procesos de falla de modelos inspirados en la
guadua [Vil12].
(a) (b)
Figura 3.3: (a) Distribución de umbrales de ruptura para los enlaces celulares, con W y C el ancho y el
centro de la distribución, respectivamente. Tomado de [Vil12]. (b) Disposición de los enlaces celulares
en un arreglo de células del tejido de parénquima.
3.5. Tracción axial y flexión 35
3.5. Tracción axial y flexión
(a)
(b) (c) (d)
Figura 3.4: Falla estructural del modelo bajo tracción axial o flexión. (a) sistema de 8 tubos con
una distribución de de dos tamaños. (b) sistema de 16 tubos con una distribución de dos tamaños.
(c) sistema de 16 tubos con todas las células del mismo tamaño. (d) sistema de 16 tubos con una
distribución de dos tamaños.
El modelo del tejido parénquima aśı construido se someterá a dos tipos de ensayos: de
tracción axial y de flexión.
La tracción axial se implementa mediante un proceso de deformación homogénea del tejido,
de la siguiente manera.
Cada célula incrementa su posición vertical y0i en una cantidad proporcional �ky0i hasta
alcanzar una nueva posición y1i = y0i(1 + �k), donde �k es el valor de deformación para el
cual se rompe un primer enlace.
Cuando se rompe el primer enlace se detiene la deformación, y se deja evolucionar el
sistema (es decir, se resuelve la segunda ley de Newton para cada célula, incluyendo una
fuerza viscosa proporcional a la velocidad) hasta que el sistema se relaja. Durante la
relajación puede ocurrir que más resortes se rompan, debido a la redistribución de cargas,
en lo que se conoce como una avalancha.
Finalizada la avalancha, se inicia una nueva deformación vertical, hasta lograr el inicio de
una nueva avalancha.
El proceso se repite hasta que las fallas provocadas por cada rutina deformación-relajación
generan una falla estructural, lo que significa que el sistema se rompe en dos (Fig.3.4).
En ese momento se alcanza el valor máximo de deformación �max =
∑n
k=1 �k.
Para cada valor de desorden ξ se analizan la estad́ıstica de las avalanchas y los valores de
deformación máxima para los cuales ocurre la falla estructural.
Para el caso de la flexión, se considera el tejido como la sección de una viga con eje neutro
inicialmente paralelo al eje axial del tejido y con radio de curvatura en un eje perpendicular (Fig.
36 Caṕıtulo 3. Modelo del tejido parénquima en dirección longitudinal
(a) (b)
Figura 3.5: Esquema de variables geométricas con las que se implementa la deformación por flexión. (a)
Variables geométricas para un sistema antes de la flexión. (b) Variables geométricas para un sistema
después de la flexión y cuando ha ocurrido la falla estructural
3.5). La deformación consiste en disminuir el radio de curvatura Rcurvatura1 en una cantidad
proporcional �kRcurvatura1 (Ec. 3.12), hasta que para alguna deformación �k algún enlace se
rompa. En este caso la deformación no solo provoca cambios en las posiciones yi de cada célula,
como en la tracción, sino que también aparecen un cambio en la posición xi y un cambio en la
orientación θi de cada célula.
Rcurvatura2 = Rcurvatura1(1− �k) (3.12)
Esta deformación se implementa como sigue:
Se define el eje neutro del tejido, que pase por la coordenada ycm del centro de masa y
sea paralelo al eje axial del sistema (Fig. 3.5).
Se define un centro de curvatura, y por lo tanto un radio de curvatura dado por la distancia
del centro de curvatura a la coordenada xcm del centro de masa del sistema. El radio de
curvatura debe ser suficientemente grande como para poder considerar el eje neutro como
un arco de circunferencia.
Cualquier punto en el eje neutro permite definir un arco de longitud S1 = Rcurvatura1θ1, y
la distancia del centro de masa de cualquier célula al eje neutro mide Ci (Fig. 3.5).
Se disminuye el radio de curvatura Rcurvatura1 en una cantidad proporcional �kRcurvatura1.
Esta deformación se hace garantizando que el eje neutro se curve, pero no se estire, tal
3.5. Tracción axial y flexión 37
cual sucede en el modelo de flexión de una viga ŕıgida.
Cada uno de los puntos del tejido que están por fuera del eje neutro se deforman respecto
a él de manera que a mayor distancia del eje neutro, mayor deformación. Por ejemplo,
los tubos que están a la derecha del eje se curvan y se estiran, mientras que los tubos que
están a la izquierda se curvan y se comprimen (Fig. 3.5). La deformación total del tubo
a lo largo de su eje está dada por Ciθ si se estira, o −Ciθ si se comprime.
Al igual que en la tracción, para cada valor de ξ se analiza la estad́ıstica de las avalanchas y
los valores de deformación para los cuales ocurre una falla estructural.
Caṕıtulo 4
Resultados
(a) (b) (c)
Figura 4.1: Las tres configuraciones a estudiar. (a) Un tejido construido con todas sus células de
igual tamaño, todas pequeñas. (b) Un tejido construido con todas sus células de igual tamaño, todas
grandes. (c) Un tejido construido con células de dos alturas posibles asignadas aleatoriamente. En
todos los casos se muestran sistemas de 16 tubos, L=16.
Como ya se comentó en caṕıtulos anteriores, el objetivo último de este trabajo es identificar
si la alternancia aleatoria de células largas y cortas en los tubos que forman el tejido parénquima
aumenta la resistencia a la falla del material ante esfuerzos de tracción o flexión. Para ello se
construyeron tres configuraciones diferentes de tejido: una en la que todas las células son del
mismo tamaño, pero pequeñas; otra en la que que todas las céulas son grandes y, finalmente,
una tercera configuración en la que las células tienen dos longitudes diferentes y se alternan de
manera aleatoria, ver (Fig.4.1). Los valores umbrales de los enlaces entre las células se obtienen
de distribuciones con diferentes niveles de desorden estructural, caracterizados por el parámetro
de desorden ξ definido en el caṕıtulo anterior (Figura 3.3). Cada realización de los tres modelos
se somete a esfuerzos de tracción y de flexión, hasta que el modelo falla (es decir, se rompe en
38
4.1. Normalización de las constantes elásticas 39
dos). El análisis se enfoca en tres aspectos principales: la distribución de avalanchas, el tamaño
de la máxima avalancha promedio que se presenta para diferentes valores de desorden y el valor
máximo de deformación que se tiene cuando se produce la falla estructural (Fig. 3.4). Cada
análisis se realizapara tres tamaños de sistema: 8× 8 células, 16× 16 células y 32× 32 células.
Cada realización del sistema (es decir, cada configuración de células con un valor de desorden
estructural ξ) se corrió 300 veces para sistemas de 8 × 8 células, 200 veces para sistemas de
16×16 células y 100 veces para sistemas de 32×32 células. El paquete completo de simulaciones
tomó 12 horas para los sistemas 8 × 8, 8 d́ıas para los sistemas 16 × 16 y 20 d́ıas de cómputo
para los 32× 32. Para ello se emplearon 10 equipos con procesadores a 2,0GHz, 4GB de Ram
y 8 núcleos cada uno.
4.1. Normalización de las constantes elásticas
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
ξ
10
1
10
2
10
3
〈∆
m
a
x
〉 
L=32 grandes
L=32 2 tamaños
(a)
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
ξ
10
1
10
2
10
3
〈∆
m
a
x
〉
L=32 pequeñas
L=32 grandes
L=32 2 tamaños
(b)
Figura 4.2: Influencia de una mayor rigidez entre tubos, como efecto de un mayor número Nh de
enlaces horizontales, sobre el promedio de la máxima avalancha 〈∆max〉 en función del desorden en
los umbrales de ruptura de los enlace ξ en ensayos de tracción axial para diferentes tamaños del
sistema. (a) Resultados cuando el efecto no se ha compensado, y todos los enlaces horizontales tienen
la misma constante elástica kh. (b) Resultados cuando el efecto se ha compensado, y el producto khNh
se mantiene constante.
Inicialmente se estudiaron dos configuraciones, una con todas las células grandes y otra con
células de dos tamaños, pero manteniendo una constante elástica k de los resortes igual a todos
los enlaces. Desde la construcción, el sistema con células de dos tamaños tenia mayor número
de enlaces entre un par de tubos adyacentes, lo que generaba un efecto de mayor rigidez entre
tubos. Para observar el efecto de esta mayor rigidez entre tubos se calculó el valor promedio de
la avalancha máxima 〈∆max〉 y el valor de máxima deformación al momento de la falla 〈�max〉
en función del desorden de enlaces ξ en ensayos de tracción axial para dos condiciones: una en
40 Caṕıtulo 4. Resultados
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
ξ
0.1
0.2
0.3
0.4
〈ε
m
ax
〉
L=16 grandes
L=16 2 tamaños 
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
ξ
0.1
0.2
0.3
0.4
〈ε
m
ax
〉
 L=32 grandes
L=32 2 tamaños
Figura 4.3: Influencia de una mayor rigidez entre tubos, como efecto de un mayor nmero Nh de enlaces
horizontales, sobre la máxima deformación bajo tracción 〈�max〉. (a) Resultados cuando el efecto no se
ha compensado, y todos los enlaces horizontales tienen la misma constante elástica kh. (b) Resultados
cuando el efecto se ha compensado, y el producto khNh se mantiene constante.
la que todos los enlaces (horizontales y verticales) tienen la misma constante de resorte (Fig.
4.2a), y otra donde el valor de la constante elástica de los resortes horizontales kh se ajustó
de manera que el producto khNh (con Nh el número de enlaces horizontales) se mantuviese
constante (Fig. 4.2b), para diferentes tamaños del sistema. Se observa en efecto que la mayor
rigidez disminuye el tamaño de la avalancha máxima. Por lo tanto, en todo lo que sigue se
utiliza kh ajustada para que el producto khNh se mantenga constante.
4.2. Resultados de ensayos de tracción
4.2.1. Histogramas de avalanchas
10
0
10
1
10
2
10
3
∆
10
0
10
1
10
2
10
3
10
4
N
(∆
)
ξ=0.01
ξ=0.1
ξ=0.2
ξ=0.3
ξ=0.4
ξ=0.5
ξ=0.6
ξ=0.8
10
0
10
1
10
2
10
3
∆
10
0
10
1
10
2
10
3
10
4
N
(∆
)
ξ = 0.2
ξ = 0.8
ξ = 0.5
(a)
10
0
10
1
10
2
10
3
∆
10
0
10
1
10
2
10
3
10
4
N
(∆
)
ξ=0.1
ξ=0.2
ξ=0.3
ξ=0.4
ξ=0.5
ξ=0.6
ξ=0.8
10
0
10
1
10
2
10
3
∆
10
0
10
1
10
2
10
3
10
4
N
(∆
)
ξ = 0.2
ξ = 0.8
ξ = 0.5
(b)
10
0
10
1
10
2
10
3
∆
10
0
10
1
10
2
10
3
10
4
N
(∆
)
ξ=0.01
ξ=0.1
ξ=0.2
ξ=0.3
ξ=0.4
ξ=0.5
ξ=0.6
ξ=0.8
10
0
10
1
10
2
10
3
∆
10
0
10
1
10
2
10
3
10
4
N
(∆
)
ξ = 0.2
ξ = 0.8
ξ = 0.5
(c)
Figura 4.4: Histogramas de avalanchas ∆ para sistemas 32×32 bajo pruebas de traccin axial (a) células
de un solo tamaño pequeñas. (b) células de un solo tamaño grandes. (c) células de dos tamaños.
4.2. Resultados de ensayos de tracción 41
Las primeras pruebas desarrolladas fueron los ensayos de tracción para las diferentes
configuraciones de tamaños . La figura 4.4 muestra la distribución de los tamaños de avalanchas
para diferentes valores de desorden estructural ξ en las tres configuraciones. En todos los casos se
observa que para valores pequeños de ξ se presentan avalanchas muy grandes, evidenciadas por
la aparición de un segundo pico al final de los histogramas (Fig. 4.4). En cambio, para sistemas
más desordenados son más frecuentes las avalanchas pequeñas. Este resultado cualitativo
coincide con los encontrados en anteriores modelos de fractura para guadua [Vil12].
4.2.2. Escalamiento finito
-80 -60 -40 -20 0 20 40
(ξ − ξ
c
)L
1/ν
10
0
10
1
10
2
〈∆
m
a
x
〉 
 /
 L
β
/ν
L=8
L=16
L=32
0 0.2 0.4 0.6 0.8
ξ 
10
0
10
1
10
2
10
3
〈∆
m
a
x
〉 
(a)
-200 -100 0 100 200
(ξ − ξ
c
)L
1/ν
10
0
10
1
〈∆
m
a
x
〉 
 /
 L
β
/ν
L=8
L=16
L=32
0 0.2 0.4 0.6 0.8
ξ 
10
1
10
2
10
3
〈∆
m
a
x
〉 
(b)
-40 -30 -20 -10 0 10 20
(ξ − ξ
c
)L
1/ν
10
0
10
1
〈∆
m
a
x
〉 
 /
 L
β
/ν
L=8
L=16
L=32
0 0.2 0.4 0.6 0.8
ξ 
10
1
10
2
10
3
〈∆
m
a
x
〉 
(c)
Figura 4.5: Comportamiento de 〈∆max〉 , según ξ. (a) sistemas con células de un solo tamaño pequeñas.
Los parámetros del escalamiento son, ξc = 0,5,
1
ν = 1,47,
β
ν = 0,8. (b) sistemas con células de un solo
tamaño grandes. Los parámetros del escalamiento son, ξc = 0,5,
1
ν = 1, 7,
β
ν = 1,1. (c) sistemas con
células de dos tamaños. Los parámetros del escalamiento son, ξc = 0,5,
1
ν = 1, 2,
β
ν = 1, 0. Las gráficas
insertas muestra el comportamiento antes del escalamiento.
Ahora, se analiza cuál seŕıa el valor cŕıtico de desorden ξc para el cual el sistema pasa
de tener avalanchas at́ıpicas muy grandes a tener avalanchas pequeñas mas frecuentes. Una
manera aproximada de encontrar esta valor cŕıtico consiste en buscar cuál valor genera un
histograma de avalanchas que distribuya aproximadamente como una ley de potencias (una
tendencia lineal en la gráfica Log-Log, del histograma). Para las tres configuraciones de células
(todas pequeñas, todas grandes, y de dos tamaños), se encuentra que valores ξ < 0,5 siempre
generan picos at́ıpicos del histograma, y para desordenes ξ > 0,5, las curvas siempre caen a
avalanchas más pequeñas al final del histograma, por lo que el valor de ξ ∼ 0,5 es nuestra
primera estimación para el desorden cŕıtico.
El siguiente paso consiste en intentar un escalamiento finito del parámetro de orden
utilizando simulaciones realizadas para diferentes tamaños del sistema. Como parámetro de
orden hemos escogido el tamaño de la máxima avalancha (menos 1), que es el mismo utilizado
por Villalobos et. al. [Vil12]. Para el escalamiento, se corrieron tres tamaños diferentes:
8 × 8, 16 × 16 y 32 × 32 células. Estas simulaciones tomaron tiempos crecientes de cómputo:
sistemas de 8 × 8, con un total de 300 pruebas para cada configuración de tamaños, tomaron
aproximadamente 12 horas; sistemas de 16 × 16, con un total de 200 pruebas para cada
42 Caṕıtulo 4. Resultados
configuración, tomaron aproximadamente 8 d́ıas, y sistemas de 32×32 células, con 100 pruebas
por configuración, tomaron aproximadamente 20 d́ıas.
En la figura 4.5 se muestran los resultados del escalamiento tomado el valor cŕıtico de
desorden como ξc = 0,5. Los exponentes con los que escala 〈∆max〉, en función del tamaño del
sistema son ν = 0,68 y β = 0,54 para los sistemas con células pequeñas, ν = 0,59 y β = 0,65 en
sistemas con células grandes y ν = 0,83 y β = 0,83 para sistema con células de dos tamaños.
Los escalamientos no son de la calidad suficiente como para concluir absolutamente que estos
son los exponentes cŕıticos para el sistema, pero al menos dan una estimación del rango de
valores posibles.
4.2.3. Falla