Métodos Generalizados Para El Cálculo Estático De Estructuras

Transcript

Proyecto Final de Carrera Métodos generalizados para el cálculo estático de estructuras de cables y simulación de la interacción dinámica catenaria pantógrafo según la norma europea EN50318 D. Miguel Such Taboada Director Dr. D. Alberto Carnicero López Madrid, 25 de mayo de 2008 Índice general 1. Introducción 1 2. Objetivos 3 3. Historia de la ecuación de la catenaria 5 4. Clasicación de las estructuras de cables 9 4.1. 4.2. 4.3. Estructuras de cables lineales . . . . . . . . . . . . . . . . . . . . . 4.1.1. Líneas de transmisión de energía eléctrica 4.1.2. Catenarias de trenes de alta velocidad 9 . . . . . . . . . . 10 . . . . . . . . . . . . 12 4.1.3. Puentes colgantes . . . . . . . . . . . . . . . . . . . . . . . . 15 4.1.4. Arcos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 4.1.5. Sistemas de transporte por cables . . . . . . . . . . . . . . . 19 Estructuras de cables planas . . . . . . . . . . . . . . . . . . . . . . 20 4.2.1. 21 Cubiertas de edicios . . . . . . . . . . . . . . . . . . . . . . Estructuras de cables tridimensionales . . . . . . . . . . . . . . . . 22 I Equilibrio estático de estructuras de cables 24 5. Métodos de cálculo. Estado del arte 26 5.1. Método de desplazamientos no lineales 5.1.1. 5.2. . . . . . . . . . . . . . . . . 29 Redes de cables . . . . . . . . . . . . . . . . . . . . . . . . . 29 El método de la rejilla . . . . . . . . . . . . . . . . . . . . . . . . . I 31 II ÍNDICE GENERAL 5.3. Método de la densidad de fuerza . . . . . . . . . . . . . . . . . . . . 34 5.4. Método de determinación de tensiones por mínimos cuadrados 37 . . . 6. Desarrollo teórico del método propuesto 40 6.1. Formulación en coordenadas locales . . . . . . . . . . . . . . . . . . 41 6.2. Formulación en coordenadas globales 44 6.2.1. . . . . . . . . . . . . . . . . . Consideraciones sobre el cable elástico . . . . . . . . . . . . 47 6.3. Generalización a 3D . . . . . . . . . . . . . . . . . . . . . . . . . . 48 6.4. Ensamblado y resolución . . . . . . . . . . . . . . . . . . . . . . . . 49 6.4.1. Referencias teóricas del problema . . . . . . . . . . . . . . . 52 6.4.2. Familia de métodos Gauss-Newton . . . . . . . . . . . . . . 52 6.4.3. Familia de métodos de región de conanza . . . . . . . . . . 7. Vericación de la implementación del modelo 54 57 7.1. Contrastación con el método de elementos nitos (MEF) . . . . . . 57 7.2. Simulación de sistema de transporte triangular . . . . . . . . . . . . 59 7.3. Comparación de un sistema de cables en 3D 62 7.4. Comparativa de cálculo de rigidez de una catenaria ferroviaria 7.5. Sistemas de transporte por cables conectados por poleas 7.6. . . . . . . . . . . . . . . . . 64 . . . . . . 66 Cálculo del pendolado de una catenaria de tren de velocidad alta . . 68 8. Ejemplo de aplicación 8.1. Creación de una malla de elementos nitos . . . . . . . . . . . . . . 73 73 9. Conclusiones 77 II Interacción Dinámica Catenaria-Pantógrafo 79 10.Estado del Arte 81 11.Formulación del problema dinámico en cables 84 11.1. Formulación del elemento co-rotacional . . . . . . . . . . . . . . . . 89 ÍNDICE GENERAL III 12.Formulación del contacto catenaria-pantógrafo 98 13.Integración temporal 107 13.1. La familia β -Newmark 13.2. El método α-Generalizado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110 112 14.Validación con la norma EN50318 114 15.Conclusiones 119 III Reducción Dinámica mediante Física Multicuerpo 120 16.Estado del arte 122 17.Frecuencias naturales y modos de vibración 124 17.1. Frecuencias propias en catenarias ferroviarias . . . . . . . . . . . . . 18.El método de la superposición modal 18.1. Condiciones iniciales 126 130 . . . . . . . . . . . . . . . . . . . . . . . . . . 19.Mecánica multicuerpo 132 134 19.1. Acoplamiento de modelos físicos . . . . . . . . . . . . . . . . . . . . 134 19.2. Aplicación a catenarias con modelos FEM 137 . . . . . . . . . . . . . . 20.Modelo multicuerpo jerárquico para la reducción del sistema 140 20.1. Formulación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 20.2. Resultados y vericación . . . . . . . . . . . . . . . . . . . . . . . . 144 20.3. Análisis de sensibilidad . . . . . . . . . . . . . . . . . . . . . . . . . 146 21.Conclusiones 156 ÍNDICE GENERAL IV IV Conclusiones y Aportaciones Originales 158 V Bibliografía 163 Índice de guras 3.1. Ejemplo de tienda romana . . . . . . . . . . . . . . . . . . . . . . . 6 3.2. Comparación entre una parábola y una catenaria . . . . . . . . . . 7 4.1. Línea de transporte de Energía eléctrica . . . . . . . . . . . . . . . 11 4.2. Partes de una catenaria . . . . . . . . . . . . . . . . . . . . . . . . . 13 4.3. Detalle de la sustentación de una catenaria . . . . . . . . . . . . . 14 4.4. Puente sobre el río Min, China . . . . . . . . . . . . . . . . . . . . . 16 4.5. Puente sobre el río Ródano . . . . . . . . . . . . . . . . . . . . . . . 16 4.6. Simulación numérica de una arcada [AGR06] . . . . . . . . . . . . . 18 4.7. Primera página del libro de De Ulloa . . . . . . . . . . . . . . . . . 19 4.8. Foto aérea de un teleférico . . . . . . . . . . . . . . . . . . . . . . . 20 4.9. Simple estructura de tensegridad . . . . . . . . . . . . . . . . . . . 23 5.1. Nodo de una red de cables . . . . . . . . . . . . . . . . . . . . . . . 28 5.2. Método de los desplazamientos no lineales . . . . . . . . . . . . . . 31 5.3. Red de cables con proyección ortogonal . . . . . . . . . . . . . . . . 32 6.1. Sistema de Coordenadas Locales . . . . . . . . . . . . . . . . . . . . 42 6.2. Relación entre coordenadas locales y globales . . . . . . . . . . . . . 44 6.3. Diagrama de cuerpo libre . . . . . . . . . . . . . . . . . . . . . . . . 46 6.4. Simple estructura de cables . . . . . . . . . . . . . . . . . . . . . . 50 7.1. Validación con MEF . . . . . . . . . . . . . . . . . . . . . . . . . . 58 7.2. Esquema de la disposición del sistema de transporte triangular . . . 60 V VI ÍNDICE DE FIGURAS 7.3. Situación inicial y nal de la estructura . . . . . . . . . . . . . . . . 63 7.4. Cálculo del equilibrio de un cable 64 7.5. Catenaria utilizada por Wu y Brennan 7.6. Comparación en la distribución de rigidez . . . . . . . . . . . . . . 66 7.7. Contraste gráco de los resultados obtenidos . . . . . . . . . . . . . 68 8.1. Posición de equilibrio de la catenaria 74 8.2. Vano central de la catenaria con malla MEF . . . . . . . . . . . . . 75 8.3. Desplazamientos desde el equilibrio de los nodos . . . . . . . . . . . 76 11.1. Prisma diferencial sometido a esfuerzo axil . . . . . . . . . . . . . . 85 11.2. Deformación de green . . . . . . . . . . . . . . . . . . . . . . . . . . 86 11.3. Deformación plana . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 11.4. Deformación del elemento corrotacional . . . . . . . . . . . . . . . . 90 12.1. Problema de contacto generalizado . . . . . . . . . . . . . . . . . . 98 . . . . . . . . . . . . . . . . . . . . . 100 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103 12.4. Penetración en arista . . . . . . . . . . . . . . . . . . . . . . . . . . 104 14.1. Catenaria de referencia EN50318 (10 vanos) . . . . . . . . . . . . . 115 14.2. Geometría y desplazamiento en los vanos centrales a 250 km/h . . . 116 14.3. Geometría y fuerza de contacto en los vanos centrales a 250 km/h . 117 14.4. Geometría y desplazamiento en los vanos centrales a 300 km/h . . . 117 14.5. Geometría y fuerza de contacto en los vanos centrales a 300 km/h . 118 17.1. Catenaria denida por la norma EN50318 . . . . . . . . . . . . . . 127 17.2. Modo de vibración 1 (1.0182 Hz) . . . . . . . . . . . . . . . . . . . 127 17.3. Modo de vibración 3 (3.0555 Hz) . . . . . . . . . . . . . . . . . . . 128 17.4. Modo de vibración 5 (5.0938 Hz) . . . . . . . . . . . . . . . . . . . 129 17.5. Modo de vibración 7 (7.1341 Hz) . . . . . . . . . . . . . . . . . . . 129 12.2. Sistemas de referencia locales 12.3. Penetración . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 VII ÍNDICE DE FIGURAS 19.1. Sistema multicuerpo . . . . . . . . . . . . . . . . . . . . . . . . . . 136 19.2. Catenaria ferroviaria EN50318 de 3 vanos . . . . . . . . . . . . . . . 138 19.3. Descomposición de la catenaria por vanos . . . . . . . . . . . . . . . 138 y B 139 20.1. Descomposición de la catenaria por vanos . . . . . . . . . . . . . . . 141 20.2. Paso de vano modal a vano FEM . . . . . . . . . . . . . . . . . . . 141 20.3. Paso de vano FEM a vano Modal . . . . . . . . . . . . . . . . . . . 141 19.4. Ampliación de la ligadura en el hilo de contacto entre los vanos A 20.4. Fuerza de contacto: FEM vs. multicuerpo FEM-Modal . . . . . . . 145 20.5. Desplazamientos: FEM vs. multicuerpo (FEM+Modal) . . . . . . . 146 20.6. Fuerza de contacto con 15 metros de análisis FEM . . . . . . . . . . 148 20.7. Fuerza de contacto con 20 metros de análisis FEM . . . . . . . . . . 149 20.8. Fuerza de contacto con 30 metros de análisis FEM . . . . . . . . . . 149 20.9. Fuerza de contacto con 50 metros de análisis FEM . . . . . . . . . . 151 20.10.Fuerza de contacto con análisis modal de 30 modos de vibración . . 152 20.11.Desplazamiento con 15 metros de análisis FEM . . . . . . . . . . . 152 20.12.Desplazamiento con 20 metros de análisis FEM . . . . . . . . . . . 153 20.13.Desplazamiento con 30 metros de análisis FEM . . . . . . . . . . . 153 20.14.Desplazamiento con 50 metros de análisis FEM . . . . . . . . . . . 154 20.15.Desplazamiento con análisis modal de 30 modos de vibración . . . . 154 20.16.Análisis de tiempos . . . . . . . . . . . . . . . . . . . . . . . . . . . 155 Índice de tablas 6.1. Ensamblado del sistema de ecuaciones . . . . . . . . . . . . . . . . 50 6.2. Algoritmo de Gauss-Newton . . . . . . . . . . . . . . . . . . . . . . 53 6.3. Algoritmo de la región de conanza . . . . . . . . . . . . . . . . . . 55 7.1. Comparación de resultados (Caso I) . . . . . . . . . . . . . . . . . . 59 7.2. Comparación de resultados (caso II.a) . . . . . . . . . . . . . . . . . 61 7.3. Comparación de resultados (Caso II.b) . . . . . . . . . . . . . . . . 62 7.4. Desviación del punto de equilibrio . . . . . . . . . . . . . . . . . . . 64 7.5. Comparación con los resultados de Peyrot . . . . . . . . . . . . . . 65 7.6. Comparativa de cálculo de rigidez . . . . . . . . . . . . . . . . . . . 67 7.7. Contraste numérico de los resultados . . . . . . . . . . . . . . . . . 69 7.8. Datos de la catenaria CRU 220 . . . . . . . . . . . . . . . . . . . . 70 7.9. Validación de la catenaria CRU 220 con 1 vano . . . . . . . . . . . 71 7.10. Validación O.Lopez-Garcia - Catenaria CRU220 - 4 vanos . . . . . . 72 14.1. Validación con el modelo de referencia . . . . . . . . . . . . . . . . 116 17.1. Sensibilidad del mallado de las frecuencias naturales . . . . . . . . . 128 20.1. Comparativa de resultados en fuerzas . . . . . . . . . . . . . . . . . 150 20.2. Comparativa de resultados en desplazamientos . . . . . . . . . . . . 151 VIII Capítulo 1 Introducción Los antecedentes de este trabajo se encuentran en la formulación e implementación de métodos numéricos ecaces para el análisis de las catenarias ferroviarias más habituales en Europa que venía realizando un grupo de investigación en la ETSIICAI de la Universidad Ponticia Comillas dirigido por el director de este proyecto. El autor del mismo contactó con dicho grupo y se planteó la posibilidad de realizar una formulación e implementación de un método general para el análisis estático, no sólo de catenarias ferroviarias, sino de cualquier tipo de estructura de cables. Además, se planteó la posibilidad utilizar dicho método para continuar con el desarrollo de un modelo de la interacción dinámica entre el pantógrafo y la catenaria de trenes de alta velocidad. Para ello, se estudió el trabajo realizado anteriormente y se ha tratado de superar para cumplir con la normativa europea EN50318 relativa a la validación modelos de simulación dinámica de la interacción catenaria-pantógrafo. Un código certicado permite validar catenarias para que puedan ser instaladas en las lineas ferroviarias europeas cumpliendo con lo establecido en las normas de interoperabilidad. La mayoría de los métodos de análisis estáticos de estructuras de cables están formulados de forma especíca para cada problemática. En el caso de cambiar alguno de los parámetros que denen el problema o tratar de extenderlos a otros tipos dejan de ser válidos. Por otro lado, actualmente existen pocos códigos de simulación de la interacción catenaria-pantógrafo aptos para certicar catenarias 1 CAPÍTULO 1. INTRODUCCIÓN 2 ferroviarias. El trabajo realizado en este campo es el que se resume en el documento presentado. El documento se dividirá en tres partes claramente diferenciadas. En la primera parte se procederá al desarrollo, formulación y validación del método para calcular el equilibrio estático de cualquier estructura de cables. En la segunda parte, se desarrollará, formulará y validará el modelo de interacción dinámica pantógrafo-catenaria y se validará según la normativa europea EN50318. Para ello se implementará una herramienta nueva de elementos nitos en MATLAB para poder implementar los nuevos modelos numéricos. La tercera parte, aprovechando la exibilidad del código de elementos nitos y la potencia del método presentado en la segunda parte, trata de obtener un modelo simplicado que permita optimizar el diseño de catenarias ferroviarias en un tiempo razonable. Cada parte consta de una breve revisión del estado del arte sobre cada uno de los temas a tratar, la formulación teórica de cada uno de los métodos, la validación de los respectivos métodos y diferentes casos de estudio y una breve conclusión. Capítulo 2 Objetivos Los objetivos planteados en la realización de este proyecto son: 1. Realizar una profunda revisión sobre los métodos empleados hasta el momento para la resolución de problemas de estructuras con cables. 2. Recopilar todos los trabajos posibles, experimentales o teóricos, que presenten resultados con los que validar los modelos que se desarrollen. 3. Desarrollar un modelo general para el cálculo de la posición de equilibrio estático de estructuras tridimensionales de cables basado en la ecuación exacta de la catenaria. 4. Implementar dichos modelos en un código ampliamente utilizado de propósito general como es Matlab. La implementación debe ser lo sucientemente exible para permitir la reproducción de cualquier problema de estructuras con cables. 5. Vericar, empleando la información recogida en la literatura cientíca, la validez del método desarrollado comprobando su exactitud, robustez y exibilidad. 6. Desarrollar un modelo de la interacción dinámica pantógrafo-catenaria basado en el método de los elementos nitos que tenga la precisión requerida por la norma europea EN50318 [CEN99]. 3 CAPÍTULO 2. OBJETIVOS 4 7. Implementar en un código de propósito general una herramienta que permita resolver problemas mediante el método de los elementos nitos. Debe hacerse de una manera lo sucientemente exible como para introducir el nuevo modelo dinámico de interacción catenaria-pantógrafo. 8. Vericar la validez del método mediante los requisitos especicados en la norma europea EN50318 y comprobando su exactitud, robustez y exibilidad. 9. Desarrollar un modelo reducido de la interacción catenaria-pantógrafo mediante la aplicación de técnicas multicuerpo jerárquicas con asignación dinámica de modelos. 10. Introducir dicho modelo en la herramienta de elementos nitos desarrollada en este proyecto. 11. Extraer las conclusiones oportunas en cuanto a la validez de los modelos. Capítulo 3 Historia de la ecuación de la catenaria Desde que el hombre aprendió a anudar y tejer bras naturales, formando así las primeras cuerdas, las ha utilizado para construir diferentes estructuras. En un principio, éstas tan solo servían como herramientas de caza y pesca. Posteriormente comenzaron a utilizarse con nes constructivos; los barcos de antiguas civilizaciones como la vikinga o la egipcia ofrecen una de las primeras referencias de estos usos, pues estaban provistos de redes para soportar y fortalecer sus velas [CCH84]. Sin embargo, el ámbito náutico no fue el único beneciado: a nivel más cotidiano, las primeras civilizaciones también se ayudaban del uso de cuerdas en tensión para levantar tiendas, así como para dotar de más estabilidad a las carpas una vez levantadas (sirvan como ejemplo las tiendas que solían transportar las legiones romanas durante las largas campañas de guerra). Por otro lado, también se hizo necesario salvar desniveles para poder desplazarse con más comodidad y velocidad. Ya en las civilizaciones del mundo antiguo, chinos e incas necesitaron, al aumentar las relaciones sociales y económicas de la época, cruzar ríos y montañas con mayor velocidad. Con este n se construyeron los primeros puentes colgantes. Estos puentes tenían la virtud de ser fáciles de fabricar y requerían un material muy ligero. Los primeros eran muy rudimentarios. 5 CAPÍTULO 3. HISTORIA DE LA ECUACIÓN DE LA CATENARIA 6 Figura 3.1: Ejemplo de tienda romana No pasaban de cuerdas o cadenas anudadas, pero la técnica de fabricación se fue perfeccionando con el tiempo, obteniéndose los precursores de los cables de acero tan usados hoy en día. Con el avance de la ciencia y la tecnología empezaron a surgir nuevas estructuras de cables. La electricación de las ciudades hizo necesario el transporte de electricidad a través de grandes líneas aéreas. Asimismo, los ferrocarriles abandonaron progresivamente el motor de vapor, y empezaron a estudiarse nuevos métodos para transmitir energía a los trenes. Los edicios, entregados al arte, empezaron a diseñarse con cubiertas curvas utilizando entramados de cables en tensión. La complejidad creciente de este tipo de estructuras hizo necesario entender mejor el comportamiento mecánico de los materiales. Ya en el siglo XV, Leonardo da Vinci había empezado a preguntarse cómo se comportaría un cable en tensión. En alguno de sus bocetos, Da Vinci fue el primero en dibujar una catenaria. En 1615 Beeckman diseñó un puente colgante suponiendo que la curva que éste adoptaba era una parábola. No obstante, esta solución no fue ampliamente conocida hasta que, dos siglos después, volviera a ser redescubierta por el ingeniero ruso Fuss, ahijado de Euler, a quien se encargó que diseñara un puente sobre el río Neva en San Petersburgo. Galileo, en matematiche, intorno à due nuove scienze, Discorsi e dimostrazioni publicado en 1638, armó que la forma CAPÍTULO 3. HISTORIA DE LA ECUACIÓN DE LA CATENARIA 7 que debe adoptar una cadena al ser colgada entre dos puntos debe ser parabólica, conclusión a la que llegó tomando como modelo el vuelo de un proyectil [Irv81]. A mediados del siglo XVII el astrónomo, físico y matemático holandés Christiaan Huygens ya sabía que Galileo estaba equivocado. No obstante, como dijo Huygens, la diferencia entre las dos curvas no es muy grande tal y como se ve en la gura 3.2. Parabola Catenaria 2.4 2.2 2 1.8 1.6 1.4 1.2 1 0.8 −1.5 −1 −0.5 0 0.5 1 1.5 Figura 3.2: Comparación entre una parábola y una catenaria En 1690 Jacob Bernoulli publicó Acta Erudiatorum, documento en el que se ex- plica por primera vez el concepto de integral. Para mostrar la potencia de la nueva herramienta de cálculo, Jacob propuso utilizarla para resolver denitivamente el problema al que Galileo no supo dar la solución correcta. Este reto fue resuelto facto de por tres personas: John Bernoulli (hermano de Jacob), Leibnitz y Huygens. Bernouilli y Leibniz aplicaron el cálculo diferencial, por aquel entonces recién descubierto. Huygens, por su parte, utilizó un método gráco. Es difícil saber quién lo hizo primero, ya que las respuestas se publicaron en un corto espacio de tiempo y la mala relación entre los autores no facilitó la tarea. Los hermanos Bernoulli además formularon la ecuación diferencial de equilibrio de una cadena sometida a diferentes estados de carga. Dentro de sus análisis llegaron a incluir la deformación elástica de los cables aplicando la ley de Hooke a sus ecuaciones. CAPÍTULO 3. HISTORIA DE LA ECUACIÓN DE LA CATENARIA Huygens fue quien le dio el nombre de catenaria 8 a la curva. Este nombre proviene de la palabra latina  catenarius , que signica cadena. También se le llamó funicular, basado en la denominación latina para cuerda. Hoy en día se reserva esta denominación para los vehículos o artefactos cuya tracción se realiza por medio de una cuerda, cable o cadena. El incremento de la complejidad de los problemas estructurales continuó planteando nuevos retos similares al de la forma de la catenaria. Un profundo estudio en el estudio de la historia de la resistencia de materiales se encuentra en el excelente libro de Timoshenko Historia de la resistencia de materiales [Tim83]. Por otro lado, no se debe olvidar la estrecha relación de este tema con el núcleo central del presente trabajo: la resistencia de materiales ha desempeñado un papel fundamental en el diseño y construcción de un sinnúmero de obras de la ingeniería, cuya belleza aún hoy nos sigue sobrecogiendo. Capítulo 4 Clasicación de las estructuras de cables Hoy en día el uso de los cables para la formación de estructuras se halla ampliamente extendido. Este fenómeno se explica al comparar el coste que suponen las estructuras rígidas con el desembolso, signicativamente menor, que demandan las estructuras de cables. Atendiendo a su conguración espacial, éstas se pueden dividir en tres grandes grupos: las estructuras de cables lineales, usadas generalmente en transporte, ya sea de energía o de objetos; las estructuras de cables planas, que gozan de una creciente popularidad debido a su belleza artística, y que se usan principalmente en edicaciones a modo de cubiertas (deben incluirse en este grupo las estructuras en forma de membrana); y, por último, las estructuras tridimensionales, las menos usuales y quizás las de menor interés práctico en la actualidad, a pesar de que en la naturaleza se encuentran muy a menudo mallas tridimensionales, compuestas por bras exibles con una innidad de utilidades. 4.1. Estructuras de cables lineales Las estructuras de cables lineales tienen la característica de avanzar en una dimensión. En general, cada cable de la estructura sólo conecta con otro cable en 9 CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 10 un punto llamado habitualmente nudo, si bien en algunos casos, como en el de las catenarias ferroviarias, se forman mallas verticales para aumentar la rigidez, con lo que conectaría con más de un cable. Para avanzar sin contactar con el suelo, la estructura está soportada con unos apoyos cuya distancia depende de la tensión del cable, su peso y la caída permitida. A los tramos de cable conectados entre dos apoyos se les llama vanos. La tensión de la línea suele transmitirse a través de poleas situadas en los apoyos. Sin embargo, debido al rozamiento que aparece en estas poleas no es posible tener un cable continuo con una sola tensión, sino que se deben formar diferentes tramos independientes mecánicamente. A continuación se presentan algunas de las tipologías más habituales. 4.1.1. Líneas de transmisión de energía eléctrica La red de transporte de energía eléctrica es la parte del sistema de suministro eléctrico constituida por los elementos necesarios para llevar hasta los puntos de consumo, y a través de grandes distancias, la energía generada en las centrales hidroeléctricas, eólicas, térmicas, de ciclo combinado o nucleares. Para ello, la energía eléctrica producida debe ser transformada previamente a un nivel superior de tensión. Esto es necesario, ya que, para un determinado nivel de potencia a transmitir, al elevar el voltaje se reduce la corriente y, por lo tanto, se reducen las pérdidas por efecto Joule. Parte fundamental de la red de transporte de energía eléctrica son las líneas de transporte. Se llama línea de transporte de energía eléctrica o línea de alta tensión al medio físico mediante el cual se realiza la transmisión de la energía eléctrica a grandes distancias. Está constituida tanto por el elemento conductor, usualmente cables de aleaciones de cobre o aluminio, como por sus elementos de soporte, las torres de alta tensión. Al estar éstas formadas por estructuras hechas de perles de acero, como medio de sustentación del conductor se emplean aisladores de disco, y herrajes para soportarlos. El proceso de tendido de una línea para transporte de energía eléctrica es una CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 11 técnica bien conocida. Se colocan unas poleas ancladas a las cadenas de aisladores que cuelgan de las crucetas de las torres y se pasa el cable. Posteriormente se procede al tensado y al engrapado del cable a las cadenas de aisladores; las compañías suelen exigir que las cadenas de aisladores queden en posición vertical. El proceso para obtener esta disposición se denomina engrapado. Figura 4.1: Línea de transporte de Energía eléctrica Existen diversos métodos de cálculo para determinar la posición de las grapas; sin embargo, la mayoría de ellos son muy simplicados, como lo demuestra la gran variedad de resultados que se obtiene para cálculos realizados sobre un mismo conjunto de vanos. A un conjunto de vanos unidos por poleas se le llama cantón. Con el método presentado en este trabajo sería posible calcular de manera exacta, teniendo en cuenta tanto la deformación elástica debido a la tensión como a la provocada por una distribución de temperaturas en los cables, la longitud y la tensión de cada tramo, así como la distancia entre el suelo y el cable conductor. Además es posible CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 12 calcular el punto en el que se deben anclar las grapas. Esta información es crítica para el diseño de la línea ya que la capacidad de la misma dependerá de la altura a la que estén los conductores en su punto mínimo. Algunas herramientas actuales ofrecen cálculos aproximados, pero generalmente no incorporan el efecto producido por la temperatura del cable. Gracias al modelo desarrollado no sólo es posible realizar dicho cálculo de una forma rápida y able, sino que además es posible aplicar métodos de optimización de estructuras para mejorar el diseño de este tipo de líneas, minimizando así coste, consumo y riesgo de fallo. 4.1.2. Catenarias de trenes de alta velocidad En el sector ferroviario, con la palabra catenaria se denomina a todo el conjunto de elementos que constituye la línea aérea de transporte y suministro de energía eléctrica a los trenes. Está situada sobre los raíles y avanza mayoritariamente en su misma dirección, aportando la energía eléctrica necesaria mediante un elemento de frotación denominado pantógrafo. El elemento fundamental de la catenaria es el cable de frotación con el pantógrafo de la locomotora; a este cable se le denomina hilo de contacto(ver gura 4.2). Para que el rozamiento entre el pantógrafo de la locomotora y el hilo de contacto sea lo más homogéneo posible, es necesario que el hilo de contacto mantenga constante su altura respecto a los carriles. Cuando las velocidades a las que se desplazan los trenes son relativamente bajas, de hasta 50 km/h aproximadamente, es suciente en el montaje de los hilos de contacto que la diferencia de altura entre los apoyos y el centro del vano sea del 1 por 1000 de la longitud del vano, y con un máximo de 20 cm, valores que se pueden conseguir mediante el propio tense mecánico del hilo de contacto. Sin embargo, cuando la velocidad aumenta, esta diferencia de alturas entre el apoyo y el centro del vano se vuelve más crítica, siendo necesaria una mayor uniformidad en las alturas. Como el tense mecánico del hilo de contacto no puede aumentar indenidamente, es necesario tender otro cable, denominado sustentador, y sujetar el hilo de contacto al nuevo cable tendido mediante unas retenciones, CAPÍTULO 4. 13 CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES Hilo Sustentador Hilo de contacto Pendola Pendola en Y 1.2 1 Altura [m] 0.8 0.6 0.4 0.2 0 −0.2 0 20 40 60 Distancia [m] 80 100 120 Figura 4.2: Partes de una catenaria denominadas péndolas, situadas longitudinalmente cada cierta distancia. De esta forma, y mediante la mayor o menor longitud de las péndolas, se consigue mantener constante la altura del hilo de contacto sobre los raíles. A la hora de montar estas estructuras es necesario conocer la longitud de las péndolas antes de ensamblar la catenaria. Aunque se realizan aproximaciones para minimizar el gasto de material, debido a la inexactitud de los métodos se tiene que realizar un calibrado manual midiendo cada péndola. Estos tendidos cubren grandes distancias y el proceso de calibrado supone un gran gasto de tiempo y de dinero. Como se muestra más adelante, la metodología tratada resuelve el problema con suciente precisión como para acelerar dicho proceso. Otro problema de especial interés desde un punto de vista cientíco y tecnológico es la interacción dinámica entre el pantógrafo y la catenaria. Adquiere gran importancia en las líneas de alta velocidad ya que, para que el tren funcione con normalidad, es necesario que el CAPÍTULO 4. 14 CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES Figura 4.3: Detalle de la sustentación de una catenaria pantógrafo no se despegue del cable y que no se produzca una vibración excesiva. Los métodos actuales para resolver problemas de dinámica del sólido deformable, como por ejemplo el método de los elementos nitos, pueden tardar entre 8 y 10 horas en calcular una respuesta de los que entre un 10 y 15 % de este tiempo se consume en el cálculo de la conguración de referencia. O. Lopez-Garcia et al. utilizaron una metodología que permitía reducir el tiempo de cálculo de la conguración de equilibrio incial, tal y como explican en [LGCT06]. Pese a que su modelo da tan buenos resultados como el mostrado en este trabajo, resulta demasiado rígido para su aplicación en problemas más generales y de mucho interés como, por ejemplo, el estudio de la zona de contacto entre un cantón y el siguiente, donde se produce una leve discontinuidad en la interacción entre la catenaria y el pantógrafo. El modelo que se presenta en este trabajo, aunque utiliza una idea similar, resulta mucho más exible y permite resolver los diferentes problemas que presentan las estructuras de cables, siendo por tanto una ecaz herramienta para el diseño de dichas estructuras. CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 15 4.1.3. Puentes colgantes Una de las construcciones que más ha impulsado el avance de este tipo de estructuras han sido los puentes colgantes. El primer puente colgante del que se tiene constancia es el construido en Yunnan, China, alrededor del 65 a.C. (si bien la identidad de su constructor constituye una incógnita) [Ron78]. En el Imperio del centro, tanto este puente como los que lo sucedieron se caracterizaban por colgar suspendidos de cadenas de hierro, algo que aún tardaría siglos en llegar a Europa. Por su parte, los Incas ya habían comunicado los Andes, antes de la llegada de Cristóbal Colón, por medio de puentes colgantes. Éstos estaban pensados para el tránsito a pie en cualquier época del año, y se construían con cuerdas tejidas a base de una hierba muy común en Sudamérica, el ichu (algunos siguen en pie hoy día, gracias a sucesivas restauraciones efectuadas con las mismas técnicas tradicionales que emplearon los primitivos artíces)[Wri00]. El diseño de estos puentes, junto a los que se construirían siglos más tarde en Europa, mejoró con el paso de los años. En estos primeros puentes colgantes el tablero estaba soportado directamente sobre los cables, por lo que tenía la forma de una catenaria cuya caída aumentaba conforme la cadena o cable se destensaba. Se añadían, además, otros cables o cadenas a mayor altura para usarlos a modo de barandilla (un ejemplo de este tipo de puentes lo se puede encontrar en la gura 4.4). Más adelante el diseño incorporaría cables secundarios, unidos al principal, que lograrían mantener la plataforma en posición horizontal. Este esquema mejorado, que cuenta con una ejemplar representación en el famoso puente sobre la Gate Golden (puerta dorada) de San Francisco, ha perdurado hasta nuestros días. Se conoce que, ya en el siglo XVII, había puentes hechos con cuerdas Europa. Muchos de ellos se construyeron con nes bélicos y hay constancia de ello en diversas crónicas de la época [Dre32]. Se cree que el primer puente europeo hecho con cadenas se construyó en Inglaterra en 1741. Contaba con 60 m de luz y, al estar destinado al uso diario de los trabajadores de las minas inglesas, su tosquedad lo situaba a años CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 16 Figura 4.4: Puente sobre el río Min, China luz del renamiento alcanzado por las estructuras chinas. Por lo demás, en Inglaterra no se vuelve a tener constancia de la existencia de ningún otro puente construido a base de cadenas hasta el año 1814. En lo que respecta al continente, la introducción en él de puentes colgantes de cables contó entre sus pioneros con los señores Sequin d'Annonay [Dre32], quienes, en 1823 propusieron al gobierno francés un diseño para la construcción de un puente de grandes dimensiones en Tournon, atravesando el río Ródano, cuyo boceto es el de la gura 4.5. Empezaron construyendo un modelo de 19 m de largo y 60 cm de ancho sobre el río Galore en Saint Vallier para obtener datos experimentales. El puente se abrió en agosto de 1825. Figura 4.5: Puente sobre el río Ródano Tras este breve recorrido histórico por la evolución de los puentes colgantes, cabe adelantar algunos comentarios relativos a las dicultades y problemas que presentan. En los puentes colgantes se dan dos tipos de problemas diferentes: el análisis de la CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 17 respuesta no lineal de los cables y el análisis de los pilares. La metodología que se desarrolla en este trabajo hace posible la resolución estática de ambos problemas, si bien es preferible el trabajo conjunto con elementos nitos para obtener un análisis más detallado de las tensiones en los pilares. Gran parte de la rigidez de los puentes colgantes proviene de la tensión de los cables. Debido a la naturaleza geométrica de esta rigidez, el sistema modica de forma no lineal sus propiedades frente a cargas externas. Cuanto mayor es la tensión a la que se está sometida la estructura, más se puede aproximar a un modelo lineal. Sin embargo, en estructuras menos rígidas estos modelos responden peor. Este problema es importante estudiarlo, ya que resulta crítico frente a la respuesta dinámica ocasionada por agentes externos. 4.1.4. Arcos A través de los siglos, los arcos se han revelado no sólo como un indispensable elemento estructural en todo tipo de construcciones, sino también como reejo de la evolución de las técnicas arquitectónicas, a menudo revelando con su forma la pertenencia de un edicio a uno u otro periodo histórico: el progresivo perfeccionamiento de su diseño ha permitido evolucionar hacia la construcción de edicios cada vez más esbeltos. Un recorrido por la evolución de los arcos debe comenzar con la inevitable mención al estilo románico, caracterizado por la omnipresencia de los arcos de medio punto. El origen de estas estructuras data de los tiempos de esplendor de la antigua Mesopotamia, pasando con posterioridad a Roma (de donde procede la particular denominación del estilo románico). Los creadores de las antiguas catedrales románicas infundieron en éstas la capacidad de transmitir quietud y recogimiento dotándolas de paredes gruesas, compactas y sin apenas ventanas para poder levantar naves que, a pesar de todo, eran, en comparación, bastante estrechas. Esto era debido en gran parte a la inecacia de los arcos de medio punto que estaban situados sobre puertas y columnas. El uso de los arcos apuntados u ojivales se introdujo en la arquitectura de la mano de los árabes y más adelante surgiría otro tipo de arcos apuntados que sería CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 18 característico del estilo gótico. El uso de este tipo de arcos, además de modicar estéticamente los edicios, aumentó la ecacia de la estructura, pues, gracias a su verticalidad, las presiones laterales se reducían considerablemente respecto a las producidas con la utilización del arco de medio punto, permitiendo así salvar mayores espacios. La evolución, a lo largo de los años, de esa idea que generó la transición a los arcos apuntados, llevó a considerar la introducción de los arcos con forma de catenaria en la construcción. Al verse sometida a una fuerza distribuida vertical, la catenaria, por razones geométricas, tan sólo soporta tensión axial. Aplicando esta idea a los arcos se obtiene una estructura que sólo se verá sometida a este tipo de esfuerzos, aumentando considerablemente la altura a la que se pueden elevar las columnas, así como la resistencia de las mismas. Utilizando esta idea, y ayudándose por modelos experimentales de cuerdas, Antoni Gaudí diseñó la Sagrada familia en 1883, iglesia que, como es sabido, sigue en construcción hoy en día. Figura 4.6: Simulación numérica de una arcada [AGR06] Con la metodología presentada en este trabajo se podrían reproducir los análisis realizados por Gaudí e incluso obtener curvas nuevas conociendo unos pocos datos como, por ejemplo, los puntos máximos deseados o la longitud de los arcos. Una vez obtenida la geometría podrían introducirse en programas de cálculo de estructuras para conrmar que el diseño tiene las propiedades deseadas. CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 19 4.1.5. Sistemas de transporte por cables Al igual que puentes colgantes, ya existían teleféricos hechos con cuerdas en Sudamérica antes de que fuera descubierta por los europeos. De Ulloa, un escritor español, describe en su libro Viaje histórico por la américa meridional, tal y como se cuenta en [Dre32], un tipo de puente llamado tarabita usado para cruzar los valles de la cordillera de los Andes. Figura 4.7: Primera página del libro de De Ulloa Un cable hecho de bambú se enviaba de un lado del valle, donde quedaba atado a un poste, a la otra ladera del valle, donde una polea servía para tensar el sistema. Elevando un extremo por encima del otro, y utilizando una canasta sucientemente grande para que un hombre se pudiera sentar en ella, era posible cruzar sin dicultad. Para facilitar el regreso se colocaba otro artilugio similar inclinado en dirección opuesta. De esta forma a los habitantes de la zona les era posible cruzar en relativamente poco tiempo la Cordillera de los Andes. Por otro lado, los teleféricos son sistemas muy utilizados en la actualidad pa- CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 20 ra el transporte de pasajeros. Generalmente se construyen con nes turísticos, ya que permiten transitar por terreno muy abrupto y en condiciones desfavorables sin necesidad de instalar una gran cantidad de postes. Esto los hace especialmente interesantes para remontes de alta montaña; de hecho, es uno de los sistemas más populares en las estaciones de esquí abiertas al público. Figura 4.8: Foto aérea de un teleférico Utilizando una conguración de tres cables es posible transportar objetos entre dos puntos cualesquiera en un área, en lugar de hacerlo entre dos puntos jos, como ocurre en los funiculares tradicionales. Una de las dicultades que entraña este método es el control de dicho sistema, ya que la posición del objeto a transportar depende de las tensiones aplicadas en los dos extremos libres del sistema. Con la metodología presentada se pueden conocer tanto la posición del objeto conocidas las tensiones como las tensiones necesarias para transportar el objeto al punto deseado. 4.2. Estructuras de cables planas Una estructura de cables es plana cuando tiene forma de malla o membrana. Este tipo de estructuras, debido a su ligereza, cuenta con una rigidez especíca bastante elevada; como en los casos unidimensionales, gran parte de la rigidez inherente del sistema viene dada por la tensión de los materiales. CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 21 4.2.1. Cubiertas de edicios Las primeras estructuras de este tipo, tal y como las se conocen hoy en día, fueron las cuatro cubiertas de los pabellones construidos por el ingeniero ruso V.G.Shookhov para una exposición en Nizjny-Novgorod en 1896. Durante los años treinta algunas cubiertas de tamaño medio se construyeron en Estados Unidos y en Europa, si bien ninguna gozó de relevancia signicativa. Cuando, en 1950, Matthew Nowicki diseño la State Fair Arena se dió un gran paso en el desarrollo de este tipo de cubiertas. Por desgracia, ese mismo año Nowicki murió en un accidente aéreo, pero su trabajo fue continuado por el arquitecto William Henry Deitrick y el ingeniero Fred Severud, quienes en 1953 completaron el edicio. Durante una visita a Estados Unidos, un estudiante alemán de arquitectura, llamado Frei Otto, vio los dibujos del Raleigh Arena en la ocina de Nueva York de Fred Severud. Otto se dio cuenta de que el proyecto aunaba muchas de sus mismas ideas para conseguir una construcción con la mínima cantidad de material. Tras su graduación en 1952, Otto comenzó a investigar sobre cubiertas colgantes. Su investigación, que fue presentada en su tesis doctoral Das Hängende Dach (La cubierta colgante), se convirtió en el primer documento dedicado exclusivamente a este tipo de estructuras. Tras interesarse por el trabajo realizado por Otto, Peter Stromeyer, dueño de una de las mayores compañías de fabricación de tiendas de campaña del mundo, contactó con el arquitecto, con lo que comenzó una fructífera relación. En 1957 Otto abrió un centro de investigación sobre construcción de estructuras ligeras en Berlín para optimizar el proceso. En 1964 añadió dicho centro de investigación al homónimo de la Universidad de Stuttgart, cuyo trabajo, desarrollado entre los años 1957 y 1965, fue publicado en los dos volúmenes de Tensile Structures [OTS67]. Frei Otto fue el responsable de la construcción y el desarrollo de gran cantidad de las estructuras tensadas construidas durante los 60 y los 70. Entre ellas, la primera gran estructura fue la del pabellón de Alemania de la exposición universal de 1967 en Montreal. CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 22 La creciente complicación de este tipo de estructuras ha sido la responsable de gran parte del desarrollo de las técnicas de cálculo para la obtención de la geometría de equilibrio durante la etapa de diseño. Estos métodos se han desarrollado desde mediados del siglo XX cuando la potencia de los ordenadores no era aún comparable a la que existe hoy en día. Por ello, estos métodos suelen ser muy rígidos y están sujetos a diferentes restricciones con el n de simplicar los cálculos. Con el procedimiento que se describe en este trabajo es posible estudiar el comportamiento estático de este tipo de estructuras, tanto las formadas por cables como aquéllas formadas por membranas con rapidez y precisión lo que permitiría diseñar estructuras aún más complejas. 4.3. Estructuras de cables tridimensionales Utilizando esta metodología también es posible resolver situaciones en las que intervengan cables interconectados formando redes tridimensionales de cables. Las aplicaciones más relevantes de este tipo de disposiciones son las denominadas Estructuras de Tensegridad. Surge así el concepto de tensegridad como principio estructural basado en el uso aislado de componentes en compresión dentro de una red de componentes en tensión, de forma que los elementos de compresión no se toquen y los elementos en tensión denan el sistema espacialmente. Estas estructuras son CAPÍTULO 4. CLASIFICACIÓN DE LAS ESTRUCTURAS DE CABLES 23 Figura 4.9: Simple estructura de tensegridad muy utilizadas en arte, ya que la forma que describen una vez completadas son muy estilizadas. En un ámbito más práctico, este tipo de estructuras establece el comportamiento mecánico de células y moléculas, así como el del ADN [Ing93]. Además, a nivel molecular, diversos estudios han analizado su utilidad para desentrañar el movimiento de los organismos unicelulares [CI99]. Parte I Equilibrio estático de estructuras de cables 24 25 En esta primera parte se propone un nuevo método para el cálculo de la posición de equilibrio estático de estructuras tridimensionales de cables. Este método se basa en las ecuaciones analíticas de la catenaria y supone una generalización de la aplicación previa para el cálculo de equilibrio inicial de catenarias realizado por el equipo de investigación en mecánica computacional del ICAI coordinado por el director de este proyecto. En los siguientes capítulos se profundizará en el método y estarán estructurados de la siguiente manera: En primer lugar se expone una revisión de los principales métodos de cálculo utilizados para su resolución, capítulo 5. A continuación, el capítulo 6 presenta el método teórico propuesto para la resolución de estructuras de cables a partir del desarrollo de las ecuaciones analíticas de la catenaria. El capítulo 7 presenta diferentes casos que permiten comprobar la robustez, precisión y exibilidad del modelo teórico y su implementación práctica, contrastando los resultados con otros publicados en diversas revisas cientícas. El capítulo 8 muestra una de las aplicaciones prácticas para las que se está empleando el modelo en la actualidad. Por último, el Capítulo 15 presenta brevemente las conclusiones del trabajo. Las referencias empleadas en el desarrollo del trabajo serán presentadas en orden alfabético al nal del documento. Capítulo 5 Métodos de cálculo. Estado del arte La mecánica de los medios continuos trata de predecir el comportamiento de los cuerpos cuando sobre ellos actúan fuerzas externas, comportamiento éste que depende de una serie de parámetros divisos en dos grandes grupos: por un lado, parámetros intrínsecos, basados en las propiedades del cuerpo o sistema que se estudia (geometría, masa o elasticidad), y, por otro lado, parámetros circunstanciales, que dependen del estado en que se encuentre el sistema (fuerzas externas, velocidad o posición). El comportamiento, pues, viene regido por un conjunto de ecuaciones en derivadas parciales acopladas, que tiene solución analítica en los casos más sencillos. Sin embargo, cuando se trata de aproximar una realidad más compleja habitualmente se emplean métodos numéricos de integración, como el método de los elementos nitos, el de las diferencias nitos, métodos espectrales, elementos de contorno, etc. Utilizando un método numérico es posible encontrar solución al problema de equilibrio inicial de sistemas de cables. En la mayoría de las estructuras, la conguración de referencia es conocida ya que esta no depende de la distribución de las tensiones internas. En las estructuras tensadas, como son las formadas por cables, la conguración inicial depende de las tensiones internas, que son a priori desconocidas, y que deben ser determinadas. La resolución de este problema constituye lo que se denomina problema de equilibrio inicial y es el paso previo a la obtención de la respuesta (ya sea estática o dinámica) de una estructura tensada frente a una 26 CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 27 acción exterior. Una manera de clasicar los diferentes métodos de resolución consistiría en la diferenciación entre los parámetros especicados por el diseñador y los que son tratados como incógnitas [HA82]. Los parámetros involucrados en un problema de equilibrio inicial son los siguientes: La topología de la estructura, que dene las conectividades de los miembros que la forman. Las cargas externas. Incluir éstas suele complicar el problema de equilibrio inicial, ya que la magnitud y la dirección de las cargas pueden depender de la conguración inicial de referencia La geometría de la estructura, uno de los dos parámetros clave del problema de equilibrio inicial, y especialmente importante para calcular las tensiones que actuarán en la estructura en cada momento: para una estructura en tensión, la curvatura es el parámetro que más afecta al comportamiento estructural; La distribución de las fuerzas internas, que se revela como el segundo parámetro clave, pues para conseguir un diseño seguro y económico es fundamental encontrar una distribución de fuerzas apropiada. El problema de equilibrio inicial es un problema estático puro, por lo que no es necesario introducir ecuaciones dinámicas. Sin embargo, algunos métodos, como, por ejemplo, el método de desplazamiento no lineal, utilizan ecuaciones cinemáticas para resolver el problema tal y como se comentará posteriormente. Este método en concreto requiere la especicación de ciertas propiedades del material, si bien dicha especicación no tiene por qué referirse necesariamente a las propiedades reales: pueden usarse propiedades cticias para controlar la solución de la conguración de referencia [HA82]. Como se ha mencionado anteriormente, las cargas externas pueden complicar el problema de equilibrio inicial, por lo que se suele asumir que los miembros de la CAPÍTULO 5. 28 MÉTODOS DE CÁLCULO. ESTADO DEL ARTE estructura no tienen peso y que ninguna carga actúa en los nodos. Sin embargo, para obtener una solución completa, las fuerzas externas estarán presentes en muchas de las ecuaciones expuestas en este capítulo, aunque normalmente sean despreciadas. Inicialmente, el único requisito sobre la conguración de referencia es que debe estar en equilibrio. Considérese un nodo i en un red de cuatro cables, como se puede observar en la gura 5.1. Las ecuaciones de equilibrio en las direcciones x ,y y z en el nodo se pueden escribir como: xk − xi xl − xi xm − xi xj − xi + Tik + Til + Tim + Fxi = 0, Lij Lik Lil Lim yj − yi yk − yi yl − yi ym − yi Tij + Tik + Til + Tim + Fyi = 0, Lij Lik Lil Lim zj − zi zk − zi zl − zi zm − zi Tij + Tik + Til + Tim + Fzi = 0, Lij Lik Lil Lim Tij (5.1) (5.2) (5.3) Como el equilibrio inicial es un problema estático, cualquier conguración con la que se satisfagan las ecuaciones anteriores en cada nodo será una solución del problema. Dependiendo de cual de los métodos de resolución que se exponen a continuación se utilice, las incógnitas de estas ecuaciones pueden ser las tensiones, las longitudes o las posiciones obteniéndose diferentes soluciones. No obstante, algunas soluciones son mejores ya que no todas responden a la realidad física que se busca. Figura 5.1: Nodo de una red de cables CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 29 A continuación se describirán los diferentes métodos que han sido utilizados por diferentes autores para obtener estas soluciones anteriormente. Se realizará un breve resumen de cada uno de ellos, resaltando a la vez tanto sus ventajas como sus inconvenientes. 5.1. Método de desplazamientos no lineales Entre los primeros métodos aplicados en la resolución de problemas de equilibrio inicial cobra especial relevancia el método del desplazamiento no lineal. Éste se basa en la técnica de los elementos nitos para el análisis del comportamiento estructural con grandes desplazamientos. Con frecuencia, el mismo algoritmo se aplica en la resolución tanto de problemas de equilibrio inicial como de problemas en los que aparezcan cargas externas. Sin embargo, este método se ve aquejado de grandes desventajas ya que es preciso tensar previamente la estructura para aproximarse al equilibrio lo que ralentiza mucho el proceso ya que requiere varios cálculos consecutivos. El método de los desplazamientos no lineales se puede resumir de la siguiente forma: primero, se establece una malla de elementos en equilibrio con una distribución distribución de fuerzas jada por el diseñador. Se crea una forma tridimensional de la malla desplazando los puntos de soporte de forma casi vertical a partir de sus posiciones iniciales hasta los puntos en los que estará anclada la estructura, y, por último, se usa un algoritmo iterativo para obtener la conguración de equilibrio de la estructura deformada. 5.1.1. Redes de cables Argyris fue uno de los primeros investigadores en utilizar el método de los desplazamientos no lineales para resolver problemas de equilibrio inicial en redes de cables, tal como describe en [AAB74]. Su método fue desarrollado para encontrar la forma de las cubiertas usadas en el estado olímpico de Munich, construido pa- CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 30 ra las olimpiadas de 1972 usando barras para representar los cables en su modelo numérico. Barnes detalla en [Bar88] un método similar en el que a una estructura inicialmente desequilibrada se le permite experimentar una vibración amortiguada hasta estabilizarse en una posición de equilibrio. Cualquiera de estos métodos permite conocer la geometría en equilibrio del problema. Sin embargo, los desplazamientos en los nodos jos pueden aumentar hasta que aparezca una distribución de fuerzas desfavorable. Por eso, una vez que los nodos jos han llegado a sus posiciones nales se realiza un ajuste de fuerzas modicando las longitudes iniciales de los elementos mediante el siguiente procedimiento. Un elemento tipo barra con un comportamiento elástico que cumpla la ley de Hooke conserva la longitud total y, por lo tanto, considerando L0 + ∆L0 = L0 + ∆L0 , donde L0 es la longitud antes del reajuste y L0 es la obtenida tras el mismo, se obtiene la relación de la ecuación 5.4. L0 = L0 + ∆L0 L0 + ∆L0 = . T 1+ 1 + AE (5.4) Después de este paso de ajuste la estructura ya no está en equilibrio, por lo que se necesitan algunas iteraciones para establecer el equilibrio nal tanto de longitudes como de fuerzas. El mencionado paso de ajuste modica el valor de las fuerzas respecto a las establecidas inicialmente. No obstante, dada la levedad de esta variación, se obtendrá una solución en la que la distribución de fuerzas será cercana a la jada inicialmente. Otra forma de mantener el control sobre las fuerzas es usar un módulo de elasticidad muy pequeño para los cables a costa de perder el control sobre la longitud de los mismos. El método de los desplazamientos no lineales se puede resumir de la siguiente manera. Las variables especicadas por el ingeniero son: CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 31 Figura 5.2: Método de los desplazamientos no lineales • la topología de la estructura • las condiciones de contorno, y • las propiedades de los materiales. Las incógnitas del problema son: • la geometría de la estructura, y • la distribución de fuerzas internas. La solución está restringida por la siguiente condición: • se debe especicar una distribución de fuerzas inicial. 5.2. El método de la rejilla Los métodos de resolución de problemas de equilibrio inicial han sido desarrollados para evitar los problemas asociados al método de los desplazamientos no lineales. Con el n de obtener un problema lineal equivalente muchos de estos métodos imponen ciertas restricciones sobre la solución. En particular, el desarrollo original de Siev y Eidelmann de 1962, pionero entre estos métodos, permite resolver la posición de equilibrio inicial de redes de cables asumiendo una condición de ortogonalidad sobre las mismas. CAPÍTULO 5. 32 MÉTODOS DE CÁLCULO. ESTADO DEL ARTE Figura 5.3: Red de cables con proyección ortogonal Su método usa las ecuaciones 5.1-5.3 a las que aplica las restricciones sobre la geometría, las condiciones de contorno y la distribución de esfuerzos internos de la red obteniendo como resultado un problema lineal cuya única incógnita es la altura de cada nodo[SE64]. Siev y Edelmann proponen asumir que la proyección horizontal del cable es ortogonal; es decir, una malla de tamaño xi = xk = xm y yi = yj = yl (véase gura 5.1), con ∆l. Aplicando esta modicación sobre las ecuaciones 5.1- 5.3: ∆l ∆l + Til =0 Lij Lil ∆l ∆l Tik + Tim =0 Lik Lim Tij Puesto que Tij L∆lij y Til L∆lil los cables en la dirección la dirección y x, (5.5) (5.6) son las componentes horizontales de las fuerzas de y Tik L∆lik y Tim L∆l im las componentes horizontales en y que no se introducen cargas externas en el plano horizontal, se demuestra que las fuerzas en dicho plano son constantes. Llamando fuerzas horizontales en el nudo reescribir las ecuaciones i en las direcciones x e y 5.3 como sigue: Hix y Hiy a las respectivamente, se pueden CAPÍTULO 5. 33 MÉTODOS DE CÁLCULO. ESTADO DEL ARTE Hix (zj − 2zi + zl ) + Hiy (zk − 2zi + zm ) + Fiz = 0 (5.7) Si se especican las componentes horizontales, la ecuación 5.7 es lineal y sus únicas incógnitas son las coordenadas z de los nodos libres. La ecuación 5.7 es la forma discreta de la ecuación de equilibrio vertical de una membrana, como se demuestra en [TWK59]: δ2z δ2z H˜x 2 + H˜y 2 + F˜z = 0, δx δy donde H˜x y H˜y (5.8) son las componentes horizontales de la distribución de fuerzas de tensado (N/m) en las direcciones x e y, respectivamente, y F˜z es la intensidad de 2 carga vertical (N/m ). El método de la rejilla se puede resumir de la siguiente manera: Las variables especicadas por el diseñador son: • topología estructural, y • condiciones de contorno. Las incógnitas del problema son: • geometría de la estructura, y • distribución de fuerzas internas. Las restricciones para la solución del problema son las siguientes: • limitado a cables rectos, • fuerzas horizontales constantes a lo largo de los cables, y • limitado a redes de cables con proyecciones planas rectas. CAPÍTULO 5. 34 MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 5.3. Método de la densidad de fuerza En la sección 5.2 se ha obtenido una solución del problema de equilibrio inicial mediante la resolución de un sistema de ecuaciones lineales equivalentes. Sin embargo, debido a las restricciones impuestas sobre los problemas, con el método anterior sólo se pueden resolver algunos de ellos. El método de la densidad de la fuerza nos permite abordar aquellos problemas sobre los que no se pueden aplicar todas las restricciones. Para obtener el sistema lineal equivalente, este método utiliza el articio matemático, desarrollado en [GB88], que se detalla a continuación. Inicialmente se parte de las ecuaciones de equilibrio de fuerzas 5.1 5.3 que son no lineales ya que la longitud de cada elemento es una función de las coordenadas de los nodos. Especicando las fuerzas y las longitudes, a partir de ahora denominadas q, en lugar de especicar las fuerzas de cada elemento las ecuaciones anteriores se ven modicadas de la siguiente manera: qij (xj − xi ) + qik (xk − xi ) + qil (xl − xi ) + qim (xm − xi ) = 0 (5.9) qij (yj − yi ) + qik (yk − yi ) + qil (yl − yi ) + qim (ym − yi ) = 0 (5.10) qij (zj − zi ) + qik (zk − zi ) + qil (zl − zi ) + qim (zm − zi ) = 0 (5.11) Con este cambio de variables se ha conseguido obtener un sistema de ecuaciones lineales cuyo estado de equilibrio tiene la densidad de fuerza indicada en cada elemento sin necesidad de imponer ninguna otra restricción. Este método es apropiado para obtener una primera aproximación; pero, si se desea estudiar más a fondo la estructura, es necesario aplicar un análisis posterior como los detallados en las secciones 5.1 y 5.2. La diferencia entre los métodos anteriores es que el método de los desplazamientos no lineales utiliza un número de ecuaciones igual al número de grados de libertad, mientras que el número de ecuaciones usado por el método de la densidad de fuerza es igual al número de restricciones adicionales impuestas, que en la mayoría de los casos suele ser menor que el número de grados de libertad, CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 35 como se demuestra en [Sch74]. Esta nueva metodología introducida por Schek ha suscitado, sin embargo, el estudio y desarrollo posterior de este método para su uso en aplicaciones diversas solventado las dicultades propias del método. Mollaert aplicó el método de la densidad de fuerza a estructuras compuestas tanto por cables como por elementos rígidos trabajando a compresión como detalla en [Mol84]. Para obtener la solución fuera del plano de los nodos jos, separó los miembros en tensión de aquellos en compresión cambiando, a continuación, los elementos substraídos de cada subestructura por fuerzas externas equivalentes diseñando, de esta forma cada parte por separado. Asimismo, el método de la densidad de fuerza se usó de forma conjunta con el de optimización por mínimos cuadrados para generar el patrón de corte de estructuras compuestas por membranas, tal y como se sugiere en [MT90]. Gracias a la simplicidad de la formulación de este método como la del de optimización por mínimos cuadrados, se pueden resolver problemas muy complejos en poco tiempo aunque se usen mallas muy nas. Estas propiedades hacen que el método de la densidad de fuerza sea preferible ante otros métodos, como por ejemplo el de la relajación dinámica, a la hora de obtener estos patrones. En la formulación de Schek se asume que la directriz de los cables es recta, lo cual va dejando de ser cierto a medida que la densidad de fuerza de los cables disminuye. No obstante, en la referencia [HA82] se extiende el método a uno más general donde esta directriz es curva, además de añadir elementos que reejan la física de una membrana. Esta ampliación se basa en asumir la matriz de rigidez geométrica como KG xg = 0 donde (5.12) KG es la matriz de rigidez geométrica de la estructura y xg el vector de coor- denadas nodales (x, y y z ). La ecuación 5.12 se puede aplicar a cualquier modelo de elementos nitos estructural y, aunque parece una ecuación de rigidez normal, las incógnitas son las coordenadas nodales en lugar de sus desplazamientos. Para es- CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 36 tructuras compuestas únicamente por elementos barra, el conjunto de ecuaciones de 5.12 es idéntico al que se obtiene utilizando el método de la densidad de fuerza. Para elementos simples estas matrices se pueden calcular analíticamente; sin embargo, al implementar muchos otros elementos éstas deben ser calculadas por integración numérica. Incluso tras conociendo la geometría la determinación de las tensiones en elementos complejos puede ser problemática. Christou implementó un elemento catenaria elástica en el método de la densidad de fuerza, considerando así la carga distribuida por los cables como reeja [Chr96]. Con la matriz de rigidez obtenida se puede resolver la geometría de equilibrio del problema, tras lo que se requiere un proceso iterativo para hallar la tensión en los cables la cual esta regida por una ecuación no lineal. No obstante, en estructuras muy tensas es común despreciar las cargas distribuidas. Más recientemente, Lai et al. han empleado el método de la densidad de fuerza para diseñar la forma de un reector desplegable con aplicaciones espaciales como describen en [LYP98]. Para ello, transformaron la membrana original en una red de cables equivalentes para, de esta forma, utilizar el conjunto de ecuaciones 5.9 - 5.11. Estos trabajos muestran como el método de la densidad de la fuerza no ha perdido su vigencia al pasar las décadas, pues, aunque fue introducido hace más de treinta años, en [LS71], aún hoy surgen nuevas áreas de aplicación. El método de la densidad de fuerza se puede resumir de la siguiente manera: Las variables especicadas por el ingeniero son: • Topología estructural, y • condiciones de contorno. Las incógnitas del problema son: • Geometría de la estructura, y • distribución de fuerzas internas. Las siguientes son las restricciones adicionales del método: CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE • se encuentra limitado a elementos cable rectos, y • la densidad de fuerza ha de estar jada para cada elemento. 37 5.4. Método de determinación de tensiones por mínimos cuadrados Todos los métodos descritos anteriormente consideran la geometría de la estructura como una incógnita del problema; ahora bien, puede darse el caso de que la geometría de la estructura sea conocida de antemano. En estas situaciones se debe determinar la distribución de fuerzas que satisface el sistema de ecuaciones, para lo cual se presentan dos métodos que derivan de la ecuación: At=f, (5.13) representando 5.13 la forma matricial de las ecuaciones 5.1 - 5.3 se obtiene Aγ,αβ =   α = x, y, z    donde β = j, k, l, m     γ = 1, 2, ..., N αβ −αγ Lγβ tγ,α = Tγ,β f = Fγ,α siendo α los grados de libertad y minado por β los nodos adyacentes al nodo de estudio deno- γ. El problema se resuelve diferente manera según la discusión del sistema de ecuaciones 5.13. Si el sistema es incompatible, es decir, que tiene más ecuaciones que incógnitas, entonces es necesario utilizar un método para buscar la solución que se aproxime más al equilibrio. Esto se puede hacer empleando el método de los mínimos cuadrados que se resume en la siguiente expresión. AT At = AT f, (5.14) No obstante es preciso remarcar que con este método tan solo obtiene el equilibrio de la estructura en el sentido de los mínimos cuadrados. Una de sus mayores CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 38 desventajas es que no se tiene ningún control sobre la distribución de fuerzas, no se pueden restringir las fuerzas de compresión y, además, la distribución de las mismas puede resultar muy irregular aunque jando algunas de ellas se puede controlar mejor. Frente a estas desventajas cabe señalar un rasgo positivo: la solución se obtiene con mucha rapidez puesto que viene dada por la resolución de un sistema simétrico de ecuaciones lineales. Si, por el contrario, el sistema es compatible indeterminado, con un número innito de soluciones para obtener el equilibrio, entonces se dene y resuelve una distribución exacta de tensiones t*. En general, este sistema de fuerzas no llevará al equilibrio por lo que las tensiones se pueden expresar como la suma de un conjunto de fuerzas ideales y de sus respectivas desviaciones del equilibrio, t = t* + ∆f. Como las fuerzas ideales se especican directamente, (5.15) ∆f se convierte en la incóg- nita del problema. La ecuación 5.13 se puede reescribir como en la ecuación 5.16, y su solución óptima se dene como el conjunto de desviaciones que tengan la menor norma euclídea. A∆t = f − ∆t*. (5.16) Recurriendo a una formulación clásica de multiplicadores de Lagrange, como reeja la expresión 5.17, es posible resolver el problema de optimización cuya solución óptima es la que se escribe de forma explícita en 5.18 ∆tT ∆t − 2kT [A∆t − (f − At*)] → m´ın .
−1 ∆t = AT AAT (f − At*) . (5.17) (5.18) Si la geometría de la estructura y la distribución de fuerzas jadas son compatibles la distribución de fuerzas de 5.15 varía muy poco de la especicada. Por ello, dado que cumple el equilibrio de manera exacta debe obtenerse una solución bastante suave. Sin embargo, en caso de que sean incompatibles, pueden aparecer grandes CAPÍTULO 5. MÉTODOS DE CÁLCULO. ESTADO DEL ARTE 39 desviaciones en la fuerzas. Una ventaja del sistema indeterminado con más incógnitas que ecuaciones es que el diseñador tiene algún control sobre la distribución de fuerzas aunque la geometría se especique de manera exacta. Por el contrario, si el sistema es indeterminado o hiperestático el procedimiento sólo está compuesto por matrices simétricas. El método se puede resumir como sigue: Las variables especicadas por el ingeniero son: • topología de la estructura, • condiciones de contorno, y • geometría de la estructura. La incógnita del problema es • la distribución de fuerzas internas. Se concluye aquí la presentación de las técnicas de resolución de estructuras de cables más empleadas. Junto a estas técnicas existen otros métodos particulares que permiten resolver multitud de problemas especícos. Algunos de estos métodos se utilizarán para contrastar los resultados. En el capítulo siguiente se presenta el desarrollo teórico del modelo propuesto que no puede ser encuadrado en ninguno de los grupos anteriores ya que todos estos tratan de resolver el equilibrio mediante la proyección geométrica de las tensiones tal y como muestra la ecuación 5.3. Sin embargo, el método propuesto obtiene las tensión horizontal y vertical de forma analítica lo que le permite mayor exibilidad al diseñador de la estructura que puede jar valores geométricos, topológicos y físicos indistintamente. Capítulo 6 Desarrollo teórico del método propuesto El estudio de las aplicaciones de la sección 3 revela que algunas de éstas, como pueden ser las líneas de transporte de energía eléctrica, las catenarias de los ferrocarriles o los funiculares, requieren un tratamiento más exacto de su comportamiento. Por ejemplo, pequeñas desviaciones entre el modelo y la realidad pueden ser muy relevantes tanto para cálculos estáticos como dinámicos, críticos en el dimensionamiento de las estructuras. Ciertamente todas estas aplicaciones requieren un modelado que responda a su realidad física pero, en algunos casos, las hipótesis que usualmente realizan los modelos desarrollados en la sección anterior no permiten reproducir el comportamiento real de cada cable. La verosimilitud de estas hipótesis se apoya en la analogía entre estructuras de cables y estructuras de barras, es decir, en el tratamiento de estructuras de cables muy tensos como elementos de directriz recta. Esto es cierto cuando despreciar las cargas distribuidas es plausible. Aunque algunos de los métodos analizados en la sección 5 utilizan elementos curvos para modelar los cables, tarde o temprano han de aplicarse ciertas simplicaciones con las que se pierde la precisión necesaria. Otra forma en la que este tipo de métodos aborda casos en los que la tensión en los cables sea pequeña es discretizar el cable continuo en elementos rectos más pequeños asimilables a los 40 CAPÍTULO 6. DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO 41 eslabones de una cadena. Sin embargo, esto conlleva un considerable aumento del número de incógnitas que resulta, especialmente crítico para aquellos casos en los que concurran en grandes desplazamientos. Este capitulo presenta la formulación analítica de un método de cálculo de estructuras de cables empleando la ecuación exacta de la catenaria. La implementación de éste método en el lenguaje de programación MATLAB R para la resolución de los distintos casos expuestos más adelante se ha consolidado en una herramienta llamada CALESCA. En primer lugar se desarrollará en coordenadas locales la ecuación de la catenaria para más adelante extenderla a un sistema generalizado de coordenadas globales que, nalmente, serán transformadas a una base tridimensional. Mediante la correcta manipulación de las mismas se podrá estudiar con precisión la física del problema de equilibrio inicial. El cálculo analítico de las tensiones de los cables permite conocer su valor exacto. Este hecho hace que el método presentado diera de los métodos usados en la actualidad, los cuales emplean proyecciones geométricas de las tensiones de los cables, consideradas constantes en muchos casos erróneamente. La resolución a partir de un punto de vista teórico ofrece la posibilidad de obtener de forma precisa la posición y la tensión de cualquier punto de la estructura, así como la longitud exacta entre dos puntos de la misma. 6.1. Formulación en coordenadas locales El comportamiento estático de un cable se rige, como se anticipó en la sección 3, por un conjunto de ecuaciones diferenciales y algebraicas. Para su deducción, considérese un cable simplemente colgado entre los puntos A y B, tal y como se muestra en la gura 6.1. El origen de coordenadas de vertical de M, ˜ Oηξ está situado sobre la siendo éste el punto mínimo de la curva. Tratando el cable como inextensible, el peso de la porción de cable desde el punto M hasta un punto cualquiera P será ws, donde s es la longitud del cable CAPÍTULO 6. 42 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO 2 B 1 A 2 M c 1 ˜ O Figura 6.1: Sistema de Coordenadas Locales entre dichos puntos. Esta fuerza actuará sobre el centro de gravedad de la porción de cable MP, la cual será compensada por la tensión del cable en denominadas, respectivamente, T y H. Al ser M P y en el punto mínimo de la curva, M H tan solo tiene componente horizontal. Utilizando el equivalente reducido al centro de gravedad, se aplican las ecuaciones de equilibrio de fuerzas en los ejes horizontal y vertical según el diagrama de cuerpo libre de la 6.3, con lo que se obtiene con lo que, escribiendo H = wc T cos ϕ = H (6.1) T sin ϕ = ws (6.2) y dividiendo 6.1 por 6.2, resulta s = c tan ϕ (6.3) Ésta es la forma intrínseca de la ecuación de la catenaria, en la cual la constante c es conocida como la constante de la catenaria. La forma cartesiana de la ecuación 6.3 es c· dξ =s dη CAPÍTULO 6. 43 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO y derivando de nuevo se obtiene 2 c· dξ ds = = 2 dη dη  1+ dξ dη 2 ! 21 Integrando esta ecuación resulta c sinh−1 donde K1 dξ = η + K1 dη es la constante de integración. Dado que el origen de coordenadas está situado bajo el punto mínimo de la curva, M, en K1 = 0, ξ = 0 dξ dη = 0, y por lo tanto resultando dξ η = sinh dη c (6.4) Integrando de nuevo se obtiene ξ = c cosh donde que K2 η + K2 c es otra constante de integración. Si en el punto K2 = 0, η = 0, ξ = c, se deduce y por lo tanto ξ = c cosh η c (6.5) que es la ecuación cartesiana de la catenaria respecto del sistema ˜ . Oηξ A partir de las ecuaciones 6.3 y 6.4 se obtiene que la longitud de arco de cable es s = c sinh η c (6.6) La tensión del cable se obtiene por la suma de cuadrados de las ecuaciones 6.1 y 6.2, con lo que se obtiene T 2 = w 2 s2 + c 2
de lo que, utilizando 6.5 y 6.6, se deduce T = wξ y, por lo tanto, T = wc cosh η c (6.7) CAPÍTULO 6. 44 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO Éstas son las ecuaciones básicas de la catenaria, que se pueden encontrar en diferentes libros de mecánica clásica, así como en parte de la bibliografía referente al estudio de estructuras de cables [Pug68, Irv81]. 6.2. Formulación en coordenadas globales Las ecuaciones anteriores son útiles en el caso de estudiar un único cable. Pero si se trata de estudiar varios cables, es necesario obtener las ecuaciones en un sistema más general que el empleado anteriormente. Para ello se de ha desarrollado una metodología similar a la empleada en [Cel06], que será generalizada más adelante para estructuras más complejas. 2 x2 B xm i n x1 1 A 2 y2 M y1 c 1 ˜ O O Figura 6.2: Relación entre coordenadas locales y globales Para ello, considerando el sistema Oxy se aplica la ecuación 6.5 a los extremos, los cuales aparecen en la gura 6.2, de lo que resulta: ξ2 = c cosh x2 − xmin c (6.8) ξ1 = c cosh xmin − x1 c (6.9) y CAPÍTULO 6. 45 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO Aplicando la ecuación 6.6 a los mismos puntos se obtiene: s2 = c sinh x2 − xmin c (6.10) s1 = c sinh xmin − x1 c (6.11) y Restando 6.8 de 6.9: ξ2 − ξ1 = y2 − y1 = c cosh xmin − x1 x2 − xmin − c cosh c c (6.12) Sumando 6.10 de 6.11: l = c sinh x2 − xmin xmin − x1 + c sinh c c (6.13) Aplicando la siguiente identidad hiperbólica a 6.12 cosh a − cosh b = 2 sinh 1 1 (a + b) · sinh (a − b) 2 2 (6.14) se obtiene:  y2 − y1 = 2c · sinh Despejando xmin x2 − x1 2c   sinh x1 + x2 xmin − 2c c  (6.15) de 6.15 se obtiene: xmin x1 + x2 = − c · asinh 2 y2 − y1
−x1 2c sinh x22c ! (6.16) Se puede apreciar que el mínimo de la curva depende exclusivamente de la posición de los extremos y de la constante. Suponiendo que éstos son conocidos, puede deducirse fácilmente la tensión vertical y horizontal en los extremos utilizando un sistema de coordenadas locales como el presentado en la sección 6.1. El diagrama de cuerpo libre de la gura 6.3 representa un cable desde el punto mínimo hasta uno de sus extremos en un sistema de coordenadas local ˜ . Oηξ Aplicando un balance de fuerzas en la dirección horizontal se obtiene: CAPÍTULO 6. DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO 46 Figura 6.3: Diagrama de cuerpo libre Th = T cos θ Sustituyendo en 6.17 T por su valor, hallado en Th = wc cosh En el punto mínimo (6.17) η = 0, cosh η c
=1 6.7, se obtiene η cos θ c y (6.18) cos θ = 1; por lo tanto Th = wc (6.19) Esta tensión horizontal permanece constante a lo largo de todo el cable, ya que no hay ninguna otra fuerza externa horizontal. La componente vertical de la tensión en los extremos puede hallarse compensando las fuerzas verticales: Tv = T sinθ = ws Sustituyendo en 6.20 s por su valor, hallado en Tv = wc sinh η  c (6.20) 6.6, resulta (6.21) Dada la posición de los extremos del cable, y, por tanto, siendo conocidas la distancia horizontal, d, y vertical, v, entre sus extremos el estado de dicho cable CAPÍTULO 6. 47 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO queda denido por las variables c, l y T. Asimismo se pueden utilizar las ecuacio- nes 6.7 y 6.13 para que, dada una de ellas, puedan calcularse las otras dos. La forma explícita y adimensionalizada de éstas referida al primer extremo es: El (c, l) = x2 − xmin xmin − x1 l − sinh − c sinh =0 c c c Ec (c, T ) = T xmin − x1 − c cosh =0 w c (6.22) (6.23) De añadirse otro cable quedaría perdida la información de la posición de uno de los extremos, con lo que aparecerían tres nuevas incógnitas. Éstas serían resueltas con las tres ecuaciones de balance de fuerzas en dichos puntos. Por lo tanto, para todo sistema bien denido, existirá un número nito de soluciones. 6.2.1. Consideraciones sobre el cable elástico Si se considera la elasticidad de los materiales, la longitud nal del cable, lf , dependerá de la tensión del cable y de las propiedades del mismo. La ley de Hook en su forma integral es Z ∆l = L T (x) ds EA y, sustituyendo en ella la tensión obtenida en la ecuación 6.7, se obtiene wc ∆L = 2EA        c 2 (xmin − x1 ) 2 (x2 − xmin ) senh + senh + (x2 − x1 ) 2 c c Por lo tanto, la longitud real del cable para cada estado de carga sería Lf inal = Lreposo + ∆L Al ser muy pequeño este incremento de longitud, la variación del peso por unidad de longitud se puede considerar constante, wLreposo ≈ wLf inal . Ésta sería la formulación completa para un solo cable perfectamente exible o para una cadena sin rozamiento colgada entre dos puntos. Un sistema más complejo estaría compuesto de más elementos interconectados en nodos. Para resolver este CAPÍTULO 6. 48 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO tipo de estructuras es necesario la generalización de este método a un sistema de referencia tridimensional. 6.3. Generalización a 3D Tal y como formuló Newton, las fuerzas en cualquier punto de un sistema mecánico en reposo deben estar equilibradas. Por ello, en cada nodo de unión de dos o más cables las fuerzas deben equilibrarse simultáneamente. Para aplicar las ecuaciones de equilibrio, estas tensiones deben descomponerse en un sistema de ejes ortogonales. Ya que la única fuerza distribuida aplicada sobre el cable es la gravedad, y solo está contenida en el plano vertical, las tensiones de los mismos también pertenecen al plano vertical que pasa por los extremos. En el sistema local para un cable contenido en un plano, la tensión en los extremos se puede dividir en y Tv . Para aplicar una matriz de transformación ortogonal se denominará Tp Th a la tensión que aparecería en el plano perpendicular al que une los extremos del cable que siempre será dirección Tp = 0. Utilizando la matriz de giro alrededor de un eje en la z,  √ ∆x ∆y 2 −∆x2   ∆y T =  √∆y2 −∆x2  0 −√ ∆y ∆y 2 −∆x2 √ ∆x ∆y 2 −∆x2 0 0    0    1 es posible obtener las tensiones en el sistema global de coordenadas a partir de las expresiones deducidas en la sección anterior. Estas tensiones tan sólo dependen de las posiciones relativas de los extremos del cable, la densidad lineal, de la catenaria, c. Tv correspondería a la componente en el eje z w, y la constante de la tensión y Th estará en el plano horizontal en la dirección que une a los extremos, o lo que es lo mismo: CAPÍTULO 6. 49 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO    √ ∆x ∆y 2 −∆x2 T  x      √ ∆y  Ty  =      ∆y2 −∆x2 Tz 0 −√ √ ∆y ∆y 2 −∆x2 ∆x ∆y 2 −∆x2 0 0     Th    0  Tp      Tv 1 (6.24) Gracias a ello es sencillo descomponer las fuerzas de cada nudos en sus componentes verticales y horizontales para obteniendo las siguientes ecuaciones de equilibrio estático: X Donde T = (Tx , Ty , Tz )T y T + Fp = 0 (6.25) Fp = (Fpx , Fpy , Fpz )T . Estas ecuaciones se pueden aplicar a cada nodo, siempre que todas las fuerzas que conuyen en el sean incógnitas. Si algún nodo es de contorno las reacciones (fuerzas exteriores) serán incógnitas. 6.4. Ensamblado y resolución Una vez conocidas las ecuaciones que denen la física del problema, se dene un sistema de ecuaciones no lineales de la siguiente forma. G (x) = 0 donde   G = (G1 , G2 , ..., Gn )  x = (x1 , x2 , ..., xn ) Dependiendo del tipo de problema las funciones objetivo, (6.26) Gi , pueden ser de dos tipos: de equilibrio o de compatibilidad. Las primeras responden al equilibrio estático en los tres grados de libertad (ecuación 6.25). Las otras, más especícas de este tipo de estructuras, ajustan la compatibilidad entre la forma, la longitud y la tensión. La forma de la catenaria se rige por el valor de la constante de la catenaria, c, fundamental para conocer la longitud y la tensión en el cable, relacionadas por las ecuaciones 6.22 y 6.23. Por otro lado, en el vector incógnitas: posiciones nodales, catenaria, c, x (cuadro 6.1) pueden diferenciarse dos tipos de (xi , yi , zi ), y propiedades del cable, la constante de la la longitud del cable, l y la tensión del cable, T. Para que el problema esté correctamente denido, deben denirse tantas incógnitas como ecuaciones. En CAPÍTULO 6. 50 DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO cada nudo interno existen tres ecuaciones de equilibrio, una en cada dirección; en cada nodo de contorno en los que se desee conocer o utilizar las reacciones también aparecerán tres ecuaciones de equilibrio; y por cada cable se deberán cumplir dos ecuaciones de compatibilidad más. Por otro lado, el número de incógnitas será el mismo ya que: en cada nudo interno se deseará conocer su localización geométrica; en cada nudo de contorno se calcularán las reacciones respecto a los ejes coordenados; y en cada cable se podrán averiguar dos de las tres propiedades características: l y T. c, Esta sería la disposición más simple ya que, sin modicar la cantidad de ecuaciones y de incógnitas, es posible aportar información adicional como podría ser alguna coordenada nodal pudiendo conocer así más información sobre los cables. A modo de resumen se muestra la siguiente tabla: Nudos internos Nudos de contorno Cables Ecuaciones P¯ Fi = 0 P¯ Fi = 0 Ec, El Incógnitas xi , yi , zi ext , Fyiext , Fziext Fxi Dos entre c, l, T Tabla 6.1: Ensamblado del sistema de ecuaciones Sirva como ejemplo la estructura que se muestra en la gura 6.4. Ésta esta compuesta por los cables 1, 2 y 3 de longitudes conocidas conectados en el nodo en el que se aplica una fuerza externa, 1 Fp = (Fpx , Fpy , Fpz )T . Figura 6.4: Simple estructura de cables Además del nodo 1 hay otros tres nodos que pertenecen al contorno y por lo tanto no formarán parte de las incógnitas ya que no se desea conocer las reacciones CAPÍTULO 6. DESARROLLO TEÓRICO DEL MÉTODO PROPUESTO 51 en los soportes. Las incógnitas que presenta este problema son las constantes de los cables c1 , c2 x1 , y1 y z1 y c3 , las coordenadas de la posición de equilibrio del nodo de conexión, y la tensión, T1 , T2 y T3 , en uno de los extremos de cada cable. Por lo tanto, resultan un total de 9 incógnitas en el problema. En cada cable hay denidas dos ecuaciones que equilibran la constante de la catenaria con la tensión y la longitud. Estas ecuaciones son en la expresión 6.23, y El1 , El2 y El3 , Ec1 , Ec2 y Ec3 , denidas según la ecuación 6.23. En el nodo interno se aplicarán las tres ecuaciones de equilibrio de fuerzas tal y como se muestra en la ecuación 6.25 cuyos valores se obtienen a partir de ecuación 6.24. La necesidad de sistematizar este tipo de cálculos de cara a su implementación precisa recurrir a la expresión matricial de las ecuaciones que rigen el comportamiento de cada cable. Así se obtiene una matriz por cada elemento que posteriormente se ensamblarán en un único sistema global de la forma 6.26. Para ello se debe relacionar las ecuaciones nodales de cada cable con las ecuaciones correspondientes en el sistema global y añadir las fuerzas internas calculadas para dicho nudo. Por último se sumarían las fuerzas externas que actúen sobre el mismo obteniendo el primer tipo de funciones objetivo. Las de segundo tipo se obtienen directamente de cada cable ya que sólo dependen de la posición de sus extremos. El sistema de ecuaciones que resulta está compuesto por un conjunto de ecuaciones no lineales, para cuya resolución existe una gran cantidad de métodos. Dada la gran variedad de problemas que es posible resolver, en algunos de ellos la no linealidad puede resultar crítica. Para abordar los diferentes casos que se presentan en este trabajo se han empleado las dos familias de métodos que se detallan a continuación. Ambas familias tratan de resolver un problema no lineal de mínimos cuadrados, el cual consiste en la búsqueda de un mínimo global a la suma de los cuadrados de las n funciones objetivo; es decir, n 1X 2 G (x) minimizar G (x) = 2 i=1 i donde x ∈ G (xk )). Para más referencias se pueden consultar los siguiente artículos [CV01] o [Sor97] que desarrollan el tema con más profundidad. Capítulo 7 Vericación de la implementación del modelo Habiendo concluido la formulación y generalización del método propuesto para la resolución de la estática de estructuras de cables, este capítulo pretende llevar a cabo una validación exhaustiva según dos vertientes claramente diferenciadas. La primera de ellas abordará la validación numérica del método propuesto respecto a otros recogidos en la literatura cientíca relativa al tema que ocupa este trabajo. Por otra parte, la exibilidad y robustez del método ha de ser igualmente vericada mediante la aplicación a distintos problemas de muy diversas áreas empleando una metodología común a todos, como es la expuesta en este trabajo. 7.1. Contrastación con el método de elementos nitos (MEF) Como primer caso susceptible de validación el método propuesto en este trabajo se ha contrastado con una formulación clásica de elementos nitos mediante la resolución de un ejemplo sencillo. Hasta el momento el MEF es reconocido como uno de los métodos más potentes para resolver este tipo de problemas debido a que, gracias a su exibilidad, resulta sencillo la sistematización y resolución de muchos 57 CAPÍTULO 7. 58 VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO problemas. Sin embargo, debido a que para obtener una precisión aceptable se debe aumentar la densidad del mallado, los tiempos de cálculo requeridos para resolver problemas relativamente sencillos crecen considerablemente. El caso que se va a analizar es un cable soportado en un extremo y al que se le aplica una fuerza de tensado horizontal en el otro extremo como se ve en la gura 7.1. Además de la gravedad, en el punto medio del cable se aplica una fuerza transversal,que sacará el cable de su plano de 100N . 100 N 0 -2 Eje z -4 -6 0 -8 -10 -2 -12 -14 0 -4 100 N 5 -6 10 15 Eje x 20 25 Eje y -8 Figura 7.1: Validación con MEF Para un problema tan simple, el tiempo de cálculo para la herramienta que se ha desarrollado es de alrededor de un segundo. El tiempo requerido para resolver este problema utilizando el modelo de elementos nitos depende críticamente del número de elementos que se haya utilizado. Hay que tener en cuenta que, mientras que la metodología que introduce este trabajo tan solo utiliza un elemento por cable, el MEF utiliza un conjunto de elementos tipo barra enlazados en los extremos. CAPÍTULO 7. 59 VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Debido a la geometría, este caso es no lineal y puede producir graves problemas de convergencia si las cargas no se aplican en el orden adecuado. Para obtener la solución por elementos nitos fue necesario aplicar la tensión horizontal en primer lugar y, una vez calculada la posición de equilibrio, introducir la gravedad y la fuerza transversal. Para realizar las simulaciones se ha utilizado una masa por unidad de longitud w = 1,34kg/m, de un área transversal de E = 2,1 · 1011 N/m2 A = 2 · 10−4 m2 y un módulo de elasticidad Para comprobar la importancia del número de elementos se ha realizado el mismo estudio con 20, 50 y 150 elementos por cable para, de esta manera, comprobar si los resultados convergen a la misma solución de equilibrio. Nodo 1 Nodo 2 x y z x y z Obtenidos 12.1003 -13.2282 6.05015 24.2006 0 0 MEF (20e) 12.1008 -13.233 6.0505 24.202 0 0 0.0042 0.0360 0.0058 0.0058 0 0 12.1004 -13.229 6.0502 24.201 0 0 0.0009 0.0058 0.0009 0.0017 0 0 MEF (150e) 12.1003 -13.228 6.0501 24.201 0 0 ∆% 4.3538 10 0.0018 0.0008 0.0017 0 0 ∆% MEF (50e) ∆% −5 Tabla 7.1: Comparación de resultados (Caso I) Se puede comprobar que conforme aumenta el número de elementos en la malla de elementos nitos, la solución tiende a la solución obtenida con el método propuesto obteniéndose diferencias muy pequeñas ya con 50 elementos. 7.2. Simulación de sistema de transporte triangular El TRS (Triangular Skyline System ) es un sistema mediante el cual se pueden transportar elementos en una zona de terreno sin dañar otros objetos que se en- CAPÍTULO 7. 60 VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO cuentren en el camino. El sistema está compuesto por tres cables suspendidos desde tres torres y conectados en un punto móvil del que podrá colgar una carga. Su principal utilidad es el traslado de objetos en áreas forestales o de alto valor ecológico. Utilizando este sistema, se podría realizar una tala controlada de un área de bosque sin dañar, tan solo, aquellos árboles que se deseen. El desplazamiento de la carga se ceñirá al área triangular formada por las bases de las torres de la que cuelgan los cables. Figura 7.2: Esquema de la disposición del sistema de transporte triangular En 1972, Kanzaki y Sakai presentaron un primer análisis estático en [KS72]. Su desarrollo no tenía en cuenta la deformación elástica de los cables y además no se pudieron obtener resultados para algunas de las posiciones de los cables. No obstante, presentaron datos contrastados que se utilizarán para vericar el correcto funcionamiento del método propuesto en este trabajo. El problema se modela como 3 cables de longitud desconocida y tensión conocida anclados en los puntos (15, 680, 771). P 1 = (260, 210, 786), P 2 = (320, 685, 790) y P3 = Para contrastar correctamente los resultados con el modelo de refe- rencia, no se tiene en cuenta la deformación elástica de los cables, si bien, ya que las tensiones son relativamente pequeñas, se pueden despreciar para cables con un CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO 61 modulo de elasticidad medio. El peso por unidad de longitud de los cables es de w = 0,1N/m y el peso de la carga es de Tension [N] Posición[m] Longitud[m] W = 100N . Kanzaki Obtenidos ∆% T1 3200 3200 - T2 2994 2994 - T3 4069 4069 - x 145,5 145,5 0 y 610,3 610,29 0,013 z 751,6 751,58 0,003 s1 418 418,06 0,06 s2 193,7 193,7 0,0048 s3 149,3 149,23 0,069 Tabla 7.2: Comparación de resultados (caso II.a) Como se puede observar las diferencia, ∆, respecto a los resultados publicados anteriormente son mínimos. Con el método que se ha desarrollado es posible resolver otro tipo de problemas asociados. Por ejemplo se puede jar una tensión y el punto del plano al que debe desplazarse el carro. En el cuadro 7.3 se muestran los resultados al resolver de nuevo el problema utilizando T3 , x e y como datos. Se puede comprobar como los resultados son prácticamente idénticos a los calculados en el caso anterior, así como a los resultados publicados anteriormente en [KS72]. Además se ha comprobado la exibilidad que ofrece esta metodología ya que permite jar aquellos parámetros que sean más convenientes para el diseño sin necesidad de desarrollar una nueva metodología. Por otro lado la robustez del método queda patente al resolver de forma adecuada diferentes problema en distintos tipos de variables. CAPÍTULO 7. 62 VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Tension[N] Posición[m] Longitud[m] Kanzaki Obtenidos ∆% T1 3200 3199,9 0,003 T2 2994 2994,1 0,002 T3 4069 4069 - x 145,5 145,5 - y 610,3 610,29 - z 751,6 751,58 0,003 s1 418 418,06 0,060 s2 193,7 193,7 0,005 s3 149,3 149,23 0,069 Tabla 7.3: Comparación de resultados (Caso II.b) 7.3. Comparación de un sistema de cables en 3D Este ejemplo contrasta los resultados obtenidos al analizar una estructura de cables tridimensional con los encontrados en [HL06]. Se trata de una estructura sencilla que consta de tres cables unidos en un punto en el que están soportados por un muelle vertical. En el punto de unión se aplica adicionalmente una carga de 1000N en el plano horizontal. Los resultados publicados anteriormente presentan el desplazamiento del punto de unión a partir de un punto arbitrario inicial. Con el método presentado es posible conocer la posición de equilibrio todos los puntos de la estructura, además de la del punto de unión. En la gura 7.3 se puede ver un esquema de la estructura en la que se aprecian los cables 1, 2 y 3 unidos en el punto A' de equilibrio. También está señalado el punto inicial A. En el cuadro 7.4 se presentan los resultados para los desplazamientos desde el punto arbitrario A al de equilibrio A0 con las diferencias, ∆, u ,v y w respecto a los resultados en [HL06]. Se puede ver como las diferencias son del mismo orden o menores que las ob- CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO 63 Figura 7.3: Situación inicial y nal de la estructura tenidas en las comparativas anteriores. El método utilizado por Huang y Lan fue contrastado con los resultados obtenidos por Peyrot y Goulois y que ha sido referenciado más de 20 veces. Éstos utilizaron un método similar para resolver la posición de equilibrio de la estructura y, por lo tanto, obtienen resultados del mismo orden. Peyrot y Goulois publican resultados de las tensiones horizontales y verticales de un cable colgado dependiendo de su tensión. Para comprobar las diferencias entre los métodos se ha contrastado con estos resultados que aparecen en [PG79] y que se muestran en la gura 7.4. A partir de los resultados presentados en el cuadro se observa como las diferencias se mantienen pequeñas sea cual sea la tensión de los cables. No obstante, crece lentamente conforme aumenta la tensión en los mismos. Esto puede ser debido al método con el que se calcula la deformación elástica. En la mayoría de los casos, se CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO ∆ [HL06] Propuesto u 26.47 26.527 0.2146 v 41.14 41.105 0.0860 w -2.87 -2.8833 0.4632 64 ( %) Tabla 7.4: Desviación del punto de equilibrio 60 m 19:945 N 19:225 N 15:744 N 9:178 N 504:43 N 22:16 N {329:08 N 3:0625 N 20 m 40 m 60 m 80 m Figura 7.4: Cálculo del equilibrio de un cable supone una tensión constante a lo largo del cable para calcular la deformación, sin embargo en el método que se presenta la deformación elástica se calcula utilizando la tensión en cada punto. 7.4. Comparativa de cálculo de rigidez de una catenaria ferroviaria El cálculo de la rigidez estática de una catenaria es una aplicación de gran interés en ingeniería ferroviaria. A día de hoy la Especicación Técnica de Interoperabilidad para catenarias hace referencia al cálculo de esta rigidez ya que es bastante fácil obtener de forma aproximada en comparación con un cálculo dinámico. CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Propuesto Peyrot y Goulois Tx Ty Tx Ty ∆x ( %) ∆y ( %) 3.0625 19.945 3.061 19.93 -0.0144387 0.00992054 9.178 19.225 9.172 19.24 0.000950228 0.0105126 22.16 15.744 22.15 15.73 -0.0180386 0.0269166 504.43 -329.08 504 -328.8 0.085245 -0.085086 65 Tabla 7.5: Comparación con los resultados de Peyrot Este cálculo tiene un fundamento completamente estático. La rigidez se dene como la relación entre la fuerza aplicada en un punto y el desplazamiento en dicho punto debido a ésta. Utilizando esta denición es posible obtener un modelo simplicado masa-muelle para realizas cálculos dinámicos aproximados como se explica en [LGCM07]. 1.6 1.4 1.2 Altura [m] 1 0.8 0.6 0.4 0.2 0 −0.2 0 50 100 150 Distancia [m] Figura 7.5: Catenaria utilizada por Wu y Brennan Wu y Brennan publicaron en [WB98a] un análisis de rigidez sobre la catenaria que se muestra en la gura 7.5. Se han contrastado los resultados publicados con los obtenidos utilizando el método propuesto en este trabajo y los resultados se muestran en la gura 7.6. Esta gura muestra la rigidez en el vano central de la CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO 66 catenaria calculada utilizando ambos métodos y la diferencia entre ellos. Figura 7.6: Comparación en la distribución de rigidez Se puede comprobar como las diferencias son pequeñas, especialmente en aquellos puntos donde hay péndolas. Esto es debido a que los datos publicados por Wu y Brennan no ofrecían tantos valores como los que se han usado para calcular la distribución de rigidez con el método propuesto y han sido interpolados. No obstante los resultados son satisfactorios. 7.5. Sistemas de transporte por cables conectados por poleas Bruno y Leonardi desarrollaron en [BL99] un método para el cálculo no lineal de sistemas de transporte a través de cables, como por ejemplo funiculares o remontes de esquí. En él desarrollan varios ejemplos en los que buscan la posición de equilibrio de dos cables conectados a través de una polea móvil. Contrastan su modelo utilizando una aproximación numérica de elementos nitos. Gracias a la versatilidad de la implementación presentada en este trabajo es posible introducir elementos móviles como las poleas e incluso obtener posiciones CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Distancia Rigidez (N/m) Diferencia (m) CALESCA Wu Brennan ∆% 0 3960.5 4094.26 3.37748 1.13636 4160.64 4299.18 3.32975 2.43506 5379.85 5672.13 5.43284 3.24675 4028.27 4299.18 6.72522 3.8961 3342.47 3459.02 3.4869 5.19481 2726.63 2905.74 6.56891 6.33117 2608.29 2659.84 1.97621 7.46753 2615.2 2659.84 1.70676 9.09091 2074.07 2106.56 1.56659 10.7143 1831.49 1881.15 2.71125 12.3377 1881.88 1881.15 0.0387361 15.4221 1522.65 1532.79 0.66552 17.5325 1585.57 1614.75 1.8405 20.4545 1374.33 1409.84 2.58364 22.5649 1470.02 1491.8 1.48188 25 1335.49 1368.85 2.49843 67 Tabla 7.6: Comparativa de cálculo de rigidez de equilibrio en las que alguna de las reacciones sea una condición impuesta, como es el caso de una polea móvil ya que se puede modelar como un nodo de contorno situado de tal forma que no aparezca reacción en el plano horizontal. Asimismo sería posible añadir el efecto del rozamiento estático de dicha polea o del sistema que utilizase para desplazarse. En este caso no se ha tenido en cuenta ningún tipo de rozamiento para poder contrastar los resultados con la referencia. El caso de vericación está compuesto por un cable de longitud módulo de Young longitud E = 2 · 106 kg/cm2 , w = 6,327kg/m. P 2 = (300, 0, 50) área A = 8,05cm2 y peso por unidad de Dicho cable está jado en los puntos unidos por una polea situada a z3 = 100m. L0 = 500m, P 1 = (0, 0, 0) y CAPÍTULO 7. 68 VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO x3=282.8 m x3=47.6 m T L2 T CALESCA Brunno 100 m L2 P1 T T P2 L1 L1 300 m Figura 7.7: Contraste gráco de los resultados obtenidos En la gura 7.7 se pueden ver las curvas obtenidas utilizando el método propuesto(en colores) superpuestas a las obtenidas por Bruno y Leonardi (lineas negras continuas y discontinuas). Asimismo, el cuadro 7.7 reeja los resultados numéricos de ambos cables. En ella, x3 es la posición de equilibrio de la polea, del cable desde el soporte izquierdo hasta la polea, la polea hasta el soporte derecho, T L2 L1 la longitud la longitud del mismo desde la tensión del cable en la polea y ∆, en tan- to por ciento. Todas las longitudes están expresadas en metros y las tensiones en kilogramos al igual que se encuentra en la referencia. 7.6. Cálculo del pendolado de una catenaria de tren de velocidad alta Como siguiente caso de validación se abordará el cálculo de la longitud de las péndolas de una catenaria de velocidad alta. Este es un claro ejemplo de problema CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO CALESCA Referencia ∆( %) x3 283.124 282.81 0.111175 L1 447.225 446.37 0.191514 L2 52.7751 52.86 0.160537 T 1833.53 1830 0.193029 x3 47.3323 47.67 0.708307 L1 110.96 111.07 0.0992642 L2 389.04 388.18 0.221612 T 1481.23 1478 0.218694 69 Tabla 7.7: Contraste numérico de los resultados de equilibrio inicial ya que la longitud de dichas péndolas es desconocida y, por lo tanto, también la posición de los hilos de contacto y sustentado como se ha explicado en la sección 4.1.2. Anteriormente, O.Lopez-Garcia et al. desarrollaron un método bidimensional especíco para este cálculo en [LGCT06] partiendo de la ecuación exacta de la catenaria, aproximando las ecuaciones de compatibilidad descritas en la sección 6.2 y añadiendo el efecto de las péndolas independientemente utilizando un doble algoritmo de resolución de sistemas de ecuaciones no lineales Newton-Raphson para corregir la longitud de las mismas. Estos autores comparan sus resultados con los obtenidos por métodos aproximados pero ampliamente utilizados en el sector ferroviario descritos por Montesinos y Carmona en [MC02]. El cálculo del pendolado se realizará sobre la catenaria de velocidad alta CRU 220 diseñada por RENFE cuyas características se denen en la tabla 7.8. A partir de estos datos, se han sido resuelto y contrastado dos problemas diferentes. El primero de ellos es el cálculo de la longitud de cada una de las péndolas, Np, para un cantón de un solo vano de la catenaria CRU 220. En este caso los resultados obtenidos se han contrastado con dos modelos aproximados que se pueden encontrar en [MC02] (Ref.1 y Ref.2) y el modelo publicado en [LGCT06]. Notar que ∆ representa la diferencia relativa en tanto por ciento entre el modelo propuesto CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Longitud del vano 52 m Sección del sustentador 2 153 mm Masa lineal del sustentador 1,414 kg/m Sección del hilo de contacto 2 120 mm Masa lineal del hilo de contacto 1,07 kg/m Número de hilos de contacto 2 Tensión mecánica en el sustentador 15892 N Tensión mecánica en hilo de contacto 15009 N Flecha inicial 20 cm Grifa sustentador 120 g Grifa hilo de contacto: 350 g Tipo de péndola Cable de cobre exible 70 2 de 16 mm -0.1 kg/m Tabla 7.8: Datos de la catenaria CRU 220 y el obtenido utilizando los resultados encontrados en [LGCT06], ya que es el más preciso de los anteriores. Como reeja el cuadro 7.9, la diferencia entre los resultados que aparecen en la literatura y los obtenidos con esta metodología es pequeña, aumentando la misma en la comparación con los modelos de Montesinos y Carmona (Ref.1 y Ref.2). No obstante, se justica ya que el método empleado por los autores aproxima el modelo y acepta algunos errores. Con el n de evidenciar la versatilidad de CALESCA sin que la precisión del método se vea resentida, se ha analizado un segundo caso en el que se calcula la longitud de las Np péndolas en un cantón de 4 vanos, contrastando los resultados posteriormente con los que aparecen en [LGCT06] (Ref. 3). Los resultados del cuadro 7.10 evidencian la precisión del modelo propuesto. La máxima diferencia encontrada es de 0,217 %, el cual es justicable debido a que en el método de O. Lopez-Garcia la constante de la catenaria se aproxima a partir de CAPÍTULO 7. VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Np Ref. 1 Ref. 2 Ref. 3 Propuesto ∆( %) 1 1.081 1.084 1.0859 1.08543 0.0437 2 0.935 0.939 0.942048 0.941346 0.0745 3 0.823 0.827 0.831808 0.830931 0.105 4 0.743 0.746 0.752262 0.75126 0.133 5 0.695 0.698 0.703388 0.702309 0.153 6 0.681 0.681 0.688061 0.686957 0.160 7 0.695 0.698 0.703388 0.702309 0.153 8 0.743 0.746 0.752262 0.75126 0.133 9 0.823 0.827 0.831808 0.830931 0.105 10 0.935 0.939 0.942048 0.941346 0.074 11 1.081 1.084 1.0859 1.08543 0.044 71 Tabla 7.9: Validación de la catenaria CRU 220 con 1 vano la tensión del cable en lugar de la tensión horizontal como dene la expresión exacta c= Tx . w En virtud de las dos comparativas llevadas a cabo en este apartado, se verica que los resultados obtenidos son correctos y corresponden a la posición de equilibrio. CAPÍTULO 7. 72 VERIFICACIÓN DE LA IMPLEMENTACIÓN DEL MODELO Np Ref. 3 Propuesto ∆( %) Np Ref. 3 Propuesto ∆( %) 1 1,08317 1,08269 0,0449 23 1,06223 1,06158 0,061 2 0,937504 0,936783 0,0768 24 0,918425 0,917494 0,101 3 0,825445 0,824544 0,109 25 0,808226 0,807073 0,142 4 0,744079 0,743047 0,139 26 0,728709 0,727397 0,180 5 0,693383 0,69227 0,1609 27 0,679853 0,678444 0,207 6 0,676231 0,675091 0,169 28 0,664532 0,663091 0,217 7 0,689731 0,688613 0,162 29 0,679853 0,678444 0,207 8 0,736774 0,735731 0,142 30 0,728709 0,727397 0,180 9 0,814485 0,813567 0,113 31 0,808226 0,807073 0,143 10 0,922885 0,922143 0,080 32 0,918425 0,917494 0,101 11 1,06489 1,06438 0,0483 33 1,06223 1,06158 0,060 12 1,06223 1,06158 0,0607 34 1,06495 1,06432 0,059 13 0,918425 0,917494 0,101 35 0,922969 0,922059 0,099 14 0,808226 0,807073 0,143 36 0,814588 0,813463 0,138 15 0,728709 0,727397 0,180 37 0,736891 0,735613 0,173 16 0,679853 0,678444 0,207 38 0,689856 0,688486 0,197 17 0,664532 0,663091 0,217 39 0,676359 0,674962 0,207 18 0,679853 0,678444 0,207 40 0,693508 0,692144 0,197 19 0,728709 0,727397 0,180 41 0,744195 0,74293 0,170 20 0,808226 0,807073 0,143 42 0,825546 0,824442 0,134 21 0,918425 0,917494 0,101 43 0,937584 0,936702 0,0941 22 1,06223 1,06158 0,0607 44 1,08323 1,08263 0,0550 Tabla 7.10: Validación O.Lopez-Garcia - Catenaria CRU220 - 4 vanos Capítulo 8 Ejemplo de aplicación 8.1. Creación de una malla de elementos nitos Una aplicación muy interesante del método desarrollado, que se está utilizando en la actualidad, es la generación de mallas de elementos nitos que cuenten con toda la información necesaria para estar en equilibrio inicial. Para ello se le debe indicar tanto la posición de los nodos como la deformación elástica de cada elemento para, de esta manera, compensar las fuerzas externas y las fuerzas internas. Como se ha comentado anteriormente, el cálculo de la respuesta dinámica de una estructura de cables, como puede ser una catenaria ferroviaria, puede tardar entre 8 y 10 horas de los cuales entre un 10 % y un 15 % se consume en el cálculo del equilibrio inicial por lo labroioso que resulta utilizando los métodos actuales tal y como se comentó al hablar del método del desplazamiento no lineal. Para aplicar esta metodología al equilibrio inicial de una estructura de cables en general, y de una catenaria ferroviaria en particular, se procedería de la siguiente manera. En primer lugar se resuelve la posición de equilibrio utilizando el método presentado, en este caso implementado en CALESCA, con lo que el resultado se obtiene en apenas unos segundos. La representación de dicha posición de equilibrio se puede ver en la gura 8.1. Esta disposición se obtiene a partir de los datos facilitados al diseñador. En este 73 CAPÍTULO 8. 74 EJEMPLO DE APLICACIÓN 1.4 1.2 1 Altura [m] 0.8 0.6 0.4 0.2 0 −0.2 0 20 40 60 80 100 Distancia [m] 120 140 160 180 Figura 8.1: Posición de equilibrio de la catenaria caso, se jarían las propiedades de los cables que forman la estructura, la disposición de las péndolas, la tensión de los hilos sustentador y de contacto y la posición del hilo de contacto. Si en lugar de jar esta posición se preere jar la longitud de las péndolas también es posible. Sin embargo, en la etapa de diseño estas longitudes aún son desconocidas y es preciso calcularlas. Fijado el tamaño de la malla de elementos nitos con que se desea discretizar la estructura, se procede a calcular posiciones intermedias de los cables para unirlos por elementos nitos tipo viga o barra según corresponda. Esto se hace utilizando la ecuación de la catenaria. Cuanto mayor sea el número de elementos mayor será la precisión del modelo MEF que se utilice. No obstante, cuanto mayor sea dicho número mayor será el tiempo necesario para su cálculo. La malla resultante se muestra en la gura 8.2. Como se puede ver los nuevos nodos coinciden con la gura anterior. Con esta malla, ya en equilibrio, el tiempo de cálculo se reduce a apenas 10 segundos utilizando una malla con 1120 grados de libertad y una precisión de 10−8 . CAPÍTULO 8. 75 EJEMPLO DE APLICACIÓN 1.2 1 Altura [m] 0.8 0.6 0.4 0.2 0 −0.2 60 70 80 90 Distancia [m] 100 110 120 Figura 8.2: Vano central de la catenaria con malla MEF En la gura 8.3 se muestra la diferencia entre el equilibrio encontrado utilizando el método propuesto y el mismo tras el equilibrado por elementos nitos. Como se puede ver, la máxima diferencia es mínima. Esta se debe fundamentalmente al método de los elementos nitos aunque es posible reducir el error tanto como se desee aumentando el número de elementos de la malla. A partir de esta disposición se pueden realizar diferentes cálculos dinámicos como por ejemplo la interacción dinámica con el pantógrafo o la respuesta frente al viento. CAPÍTULO 8. 76 EJEMPLO DE APLICACIÓN −7 3 x 10 Diferencia 2.5 Diferencia [m] 2 1.5 1 0.5 0 0 20 40 60 80 100 Distancia [m] 120 140 160 180 Figura 8.3: Desplazamientos desde el equilibrio de los nodos Capítulo 9 Conclusiones En esta primera parte se propone un nuevo método para el cálculo estático de estructuras tridimensionales de cables. Este método se basa en las ecuaciones analíticas de la catenaria y supone una generalización de una aplicación previa para el cálculo de equilibrio inicial de catenarias presentado en [LGCT06]. La validación del método se ha llevado a cabo mediante la comparación con problemas de estructuras de cables publicados en la literatura cientíca. La similitud de estos resultados con los publicados es muy elevada, no superando las diferencias el casos y el 1% 7% en el peor de los en valor medio. Las aportaciones más relevantes de este método respecto al anteriormente citado son las siguientes: El método anterior consideraba sólo cables estructuras de cables bidimensionales en el plano vertical. En el nuevo método, aprovechando que las catenarias están contenidas en un plano vertical, se incorpora una formulación semianalítica completamente tridimensional. Además se tiene en cuenta la elasticidad de los cables. A pesar de que esta suele ser despreciable algunas aplicaciones requieren la consideración de deformación elástica en estructuras de cables. Las restricciones propias de las estructuras de cables son incorporadas en la 77 CAPÍTULO 9. CONCLUSIONES 78 formulación de una forma natural. Además pueden emplearse para incluir modelos de dispositivos mecánicos como poleas, pretensores de cables, soportes, muelles, etc. Respecto a otros métodos de simulación de estructuras de cables las principales ventajas que éste presenta son las siguientes: Debido a que este método está basado en las ecuaciones exactas de la catenaria, las ventajas propias del método propuesto en [LGCT06] le son inherentes. Desde un punto de vista numérico, el método presenta una alta eciencia, no dependiendo de la discretización espacial, determinandose el tamaño del problema de la topología de la estructura. Por otra parte la mayoría de algoritmos de resolución de sistemas de ecuaciones no lineales algebraicos requieren el cálculo de la matriz jacobiana. A pesar de que su cálculo es laborioso, las expresiones analíticas de la matriz jacobiana son intrínsecas al método lo cual le reporta una mayor rapidez y precisión. El problema de equilibrio inicial y el cálculo del equilibrio estático bajo cargas se resuelven siguiendo el mismo algoritmo. Gracias al tratamiento de los parámetros conocidos y desconocidos ambos tipos de problemas se pueden resolver fácilmente. El método tratado en esta parte es perfectamente válido para la mayoría de las aplicaciones ingenieriles tales como catenarias ferroviarias, sistemas de transporte de energía eléctrica, redes de metro, funiculares, etc. Parte II Interacción Dinámica Catenaria-Pantógrafo 79 80 El método desarrollado en la primera parte permite obtener una malla de elementos nitos de una catenaria en equilibrio inicial. Aprovechando dichos resultados, en la segunda se propone un nuevo método para simulación de la interacción dinámica pantógrafo-catenaria basado en el método de los elementos nitos. La resolución de este problema utilizando dicha técnica tiene varios problemas desde el punto de vista numérico y matemático. En la presente parte se detalla el método propuesto para resolverlo y se estructurará de la siguiente manera: En primer lugar se expone una revisión de los principales métodos de cálculo utilizados en la actualidad para resolver problemas dinámicos estructurales, capítulo 10. A continuación, los capítulos 11, 12 y 13 presentan la formulación utilizada para la resolución de la interacción dinámica catenaria-pantógrafo. El capítulo 14 presenta la validación del método con la norma EN50318. Por último, el Capítulo 15 presenta brevemente las conclusiones del trabajo. Las referencias empleadas en el desarrollo del trabajo serán presentadas en orden alfabético al nal del documento. Capítulo 10 Estado del Arte El conocimiento de la dinámica del sistema catenaria-pantógrafo es crucial para el correcto diseño de líneas de alta velocidad. En los últimos años se ha hecho un gran esfuerzo en comprender el comportamiento del sistema con el tren en marcha: aproximación analítica de la fuerza de contacto (Ockendon, [OT71]), análisis por elementos nitos de la variación de la fuerza de contacto (Vinayagaligam, [Vin83]), determinación de frecuencias y modos propios (Wormley, [EOSW88]), etcétera. Uno de los trabajos más interesantes que se ha publicado en este contexto es el realizado por T.X. Wu y M. J. Brennan, tem Dynamics Basic Analytical Study of Pantograph-catenary Sys- [WB98b], en el que se propone un modelo simplicado del sistema catenaria-pantógrafo de un grado de libertad. Parte de un modelo de pantógrafo de dos grados de libertad. Con una masa superior representa las pletinas del pantógrafo y con otra inferior introduce el efecto de la inercia del armazón articulado. Por su parte, la catenaria se sustituye por una serie de resortes de distinta rigidez contra los que las pletinas del pantógrafo contactan. Si la rigidez del pantógrafo es muy superior a la que presenta la catenaria pueden condensarse las masas del pantógrafo en una única masa. De este modo, la catenaria quedaría representada mediante un resorte de rigidez variable. Es importante observar que en este modelo simplicado se desprecia la masa de la catenaria. La distribución de rigidez a lo largo de la catenaria no varía de forma suave (aparecen picos de rigidez al paso 81 CAPÍTULO 10. 82 ESTADO DEL ARTE por las péndolas), lo que diculta el análisis matemático del comportamiento del sistema -respuesta temporal, estabilidad del pantógrafo, etc.- al no disponerse de una expresión analítica que proporcione la rigidez. Representar la variación de rigidez en la catenaria por medio de la envolvente puede ser una buena aproximación si los incrementos de rigidez al paso por las péndolas son reducidos. Sin embargo, cuando dichos incrementos alcanzan un valor apreciable, el empleo de la envolvente como distribución de rigidez deja de ser apta debido a los errores introducidos. Este modelo, aunque simple, es capaz de describir gran parte de las características dinámicas del sistema, permitiendo efectuar análisis interesantes, como son la estabilidad del conjunto catenaria-pantógrafo o la obtención analítica tanto del desplazamiento vertical de la catenaria como de la fuerza de contacto. Estos mismos autores publicaron en 1999 una continuación del estudio anterior (Dynamic stiness of a railway overhead wire system and its eect on pantograph- catenary system dynamics [WB99]) donde se analiza el efecto que la propagación de ondas tiene sobre la rigidez de la catenaria. Para ello desarrollaron un nuevo modelo de catenaria que contemplara aspectos no recogidos en el modelo anterior, contando con la masa de la catenaria y permitiendo introducir el efecto de las ondas de perturbación sobre el comportamiento del sistema. O.Lopez-García y J.L.Maroño profundizaron en este modelo en Inuence of stiness and contact modelling on catenarypantograph system dynamics [LGCM07] donde realizan un análisis de la sensibilidad de la fuerza de contacto a la rigidez de la catenaria. Otro trabajo relevante es el publicado por G. Poetsch et al., ry dynamics and control Pantograph/Catena- + [PEM 97]. Este interesante documento aborda gran parte de la problemática que plantea el contacto entre catenaria y pantógrafo, especialmente de cara al análisis matemático de la dinámica del sistema. Resulta especialmente interesante el epígrafe 2, Modelling and simulation of the system dynamics, en el que se hace un repaso de los principales métodos de modelado tanto de la catenaria como del pantógrafo. Asimismo, se analizan una serie de métodos numéricos CAPÍTULO 10. 83 ESTADO DEL ARTE apropiados para la resolución de las ecuaciones dinámicas resultantes. Para nalizar, aparece un completo análisis de las ventajas que podría ofrecer la implantación de pantógrafos activos con ánimo de conseguir un mejor comportamiento del sistema. Por último, el trabajo publicado por Tong-Jin Park et al., analysis for the pantograph of a high speed rail vehicle Dynamic sensitivity [PHJ03] analiza las caracte- rísticas dinámicas del sistema mediante elementos nitos, de forma que sea posible efectuar un análisis paramétrico de sensibilidad y mejorar las prestaciones del conjunto. Concretamente se presenta la catenaria como una malla de elementos nitos tipo viga y el pantógrafo como un sistema discreto masa-muelle-amortiguador de tres grados de libertad. El modelo resultante del acoplamiento catenaria-pantógrafo se ensayó de forma que se apreciaran las diferencias en la fuerza de contacto según se hacían variar los parámetros del pantógrafo. De esta forma pudieron obtenerse conclusiones interesantes acerca de los valores que han de tener los distintos elementos del pantógrafo de cara a optimizar el comportamiento del sistema. Capítulo 11 Formulación del problema dinámico en cables En el cálculo estructural, las ecuaciones de campo que rigen el comportamiento del medio continuo se obtienen mediante la aplicación de las ecuaciones de equilibrio, compatibilidad y comportamiento. Éstas generan un sistema de ecuaciones diferenciales que relaciona las variables de campo con funciones conocidas que recogen el efecto de los distintos parámetros involucrados en el problema. La descripción general del movimiento en dirección longitudinal de una barra puede describirse, por tanto, partiendo de estas ecuaciones. Considerando el sólido de la gura 11.1 sometido a esfuerzos en dirección axil, es posible obtener la primera ecuación diferencial del equilibrio. ∂N (x, t) ∂2u = ρA 2 − n(x, t) ∂x ∂t siendo N (x, t) y n(x, t) (11.1) cargas en dirección axil, concentrada y distribuida res- pectivamente. La ecuación de comportamiento (11.2) es la que relaciona, en segundo lugar, la carga impuesta dulo elástico E y el área N (x, t) A con la deformación ε(x) del sólido mediante el mó- de la barra. N (x, t) = AEε(x) 84 (11.2) CAPÍTULO 11. 85 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES n(x) N ( x) N (x)+@N(x) d(x) Figura 11.1: Prisma diferencial sometido a esfuerzo axil Por último, la ecuación de compatibilidad plasmará la relación de la deformación ε(x) y los desplazamientos. Para el caso unidimensional y de acuerdo a una formulación lagrangiana, dado que el objeto de estudio es el propio sólido, se partirá de la medida de la deformación de Green: ε(x) = donde lo y l l2 − lo2 2lo2 (11.3) las longitudes inicial y deformada de la barra respectivamente. Esta medida de la deformación converge a la denición general para pequeñas deformaciones como se demuestra haciendo el desarrollo en serie de Taylor para l ≈ lo de (11.4). (l + ∆l)2 − l2 2l2 2 1 l + ∆l2 + 2l∆l − l2 ∆l = ≈ 2 l2 l ε(l ≈ lo ) = (11.4) Con la perspectiva de la formulación dinámica de cables, es preciso contemplar la no linealidad geométrica y su inuencia en la deformación de los mismos. Esto exige extender la denición de Green (11.3) al estado de deformación bidimensional que reeja la gura 11.2: CAPÍTULO 11. 86 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES z |)dx ( 1+du dx | | dx (dw dx ) dx x Figura 11.2: Deformación de green Considerando para este propósito el segmento elemental al eje x, dx inicialmente paralelo si su estado deformado es el de la gura 11.2, su longitud nal ds puede ser evaluada mediante el teorema de Pitágoras sobre los desplazamientos, ecuación (11.5): 2 ds =  ∂u dx + dx ∂x 2  + ∂w dx ∂x 2 (11.5) sustituyendo (11.5) en (11.3), la componente horizontal (la componente de interés para una barra de directriz recta sometida a esfuerzos de tracción-compresión) de la deformación bidimensional de Green se dene como: CAPÍTULO 11. 87 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES ds2 − dx2 2dx2 !   2  2 ∂u ∂w 1 1+ + −1 = 2 ∂x ∂x  2  2 ! ∂u 1 ∂u ∂w = + + ∂x 2 ∂x ∂x ε(x) = siendo u w el desplazamiento horizontal y (11.6) el vertical del sólido dx. Por consi- guiente, la ecuación de compatibilidad quedará como (11.7). ∂u 1 + ∂x 2 ε(x) =  ∂u ∂x 2 1 + 2  2 ! ∂w ∂x (11.7) Combinando las tres ecuaciones diferenciales (11.1), (11.2) y (11.7) se obtiene: ∂ EA ∂x ∂u 1 + ∂x 2  ∂u ∂x 2  1 + 2 ∂w ∂x 2 ! = ρA ∂2u − n(x, t) ∂t2 (11.8) Aplicando el principio de los trabajos virtuales, Z ∂W = σ∂εdΩ (11.9) Ω siendo Ω el dominio de integración, se obtiene la formulación débil de este problema. Sustituyendo en ella la expresión (11.8) y aplicando las funciones de forma, ψi , sobre las deformaciones virtuales se obtiene Z 0 l ∂ψi EA ∂x ∂u 1 + ∂x 2  ∂u ∂x 2 1 + 2  ∂w ∂x 2 ! Z dx + ψi ρA 0 Z l ∂2u dx = ∂t2 l ψi n(x, t)dx + (ψi N )0l (11.10) 0 donde habiendo integrado por partes y aplicado el teorema de Green al dominio Ω, los trabajos interno y externo equivalen a los términos izquierdo y derecho de la ecuación (11.10) respectivamente. Análogamente, la descripción general del movimiento en dirección transversal a una barra también puede puede describirse partiendo de las ecuaciones de equilibrio, CAPÍTULO 11. 88 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES z 2' b' w +dw | |dx dx w2 a' w 1' u1 w1 x+u x+u+dx+du | dx dx 1 a x b 2 dx u2 L x Figura 11.3: Deformación plana comportamiento y compatibilidad. No obstante, se puede obtener considerando los efectos de exión en la expresión de la deformación de Green, como se observa en la gura 11.3, desarrollada en la ecuación (11.6) y despreciando los efectos de rotación y alargamiento del eje debido a la exión. Haciendo esto la deformación total queda expresada como ∂2w 1 ∂u −z 2 + ε= ∂x ∂x 2  ∂u ∂x 2  + ∂w ∂x 2 ! (11.11) Insertando esta nueva deformación en la ecuación (11.9) y considerando que σ = −z MIyy obtenemos la expresión generalizada ante la combinación de cargas axiales y ectoras Z 0 l ∂ψi EA ∂x ∂u 1 + ∂x 2  ∂u ∂x  2 !  2  Z l 2 1 ∂w ∂ ψi ∂ w dx + + EI dx y 2 2 ∂x ∂x2 0 ∂x Z l Z l ∂2u + ψi ρA 2 dx = ψi n(x, t)dx + (ψi N )0l (11.12) ∂t 0 0 2 Empleando la aproximación de Galerkin para el campo de desplazamientos, CAPÍTULO 11. FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES 89 (u, w) u= X w= X ψj Uj ψj Wj (11.13) y empleando funciones de soporte local, se llega a la discretización clásica del método de los elementos nitos. Añadiendo los efectos dinámicos relativos a la aceleración como esfuerzos proporcionales a la masa y los relativos a la velocidad como proporcionales al amortiguamiento, la expresión (11.12) queda representada como Mu¨ + Cu˙ = f − q donde (11.14) M es la matriz de masa, C es la matriz de amortiguamiento, f es el vector de fuerzas externas, q es el vector de fuerzas internas, y u es el vector (incógnita) de desplazamientos nodales. Para resolver este sistema se deben utilizar métodos de integración numérica. Estos métodos utilizan, generalmente, la matriz tangente, o jacobiano, al vector de fuerzas internas, q, y se denota como Kt . 11.1. Formulación del elemento co-rotacional El enfoque co-rotacional, entendido como una vía alternativa para formular elementos nitos no lineales, ha suscitado un creciente interés durante la última década, ver [RB86, RnO88, Cri90, nNOR91, PE95, PE97, HL00]. La idea subyacente en este contexto es la descomposición del movimiento del elemento en una parte rígida y otra deformable mediante la referencia a un sistema de coordenadas locales (xl , zl ) solidarias al elemento y, por tanto, a sus movimientos de rotación y traslación, ver gura 11.4. Así, en el movimiento del elemento desde su estado original hasta su conguración deformada se distingue: en primer lugar una traslación y rotación rígida del elemento y; en segundo lugar, una deformación en el sistema de coordenadas local. CAPÍTULO 11. 90 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES z xl ¹ u 2' µ2 µ l2 zl µ l1 µ1 w2 ® ¯ u1 2 1' u2 w1 ¯0 1 x Figura 11.4: Deformación del elemento corrotacional Asumiendo que la longitud del elemento se hubiera elegido convenientemente, estas deformaciones serán siempre pequeñas respecto a los ejes locales y, consecuentemente, podrán expresarse como una no linealidad de bajo orden. No obstante, ha de hacerse énfasis en que el bajo orden de la no linealidad es sólo aparente, pues la no linealidad geométrica permanece realmente incluida en el movimiento del sistema coordenado local. La principal ventaja del enfoque co-rotacional estriba en la separación articial de las no linealidades material y geométrica en el elemento. Mientras la parte material tiene lugar en el sistema local donde se asume linealidad geométrica, la no linealidad geométrica se hace presente en las rotaciones rígidas y traslaciones del elemento sin deformar. Incluso incluyendo en la denición de la deformación términos que reejen una no linealidad de bajo orden, las expresiones del vector de CAPÍTULO 11. 91 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES esfuerzos internos y de la matriz de rigidez tangente resultantes son muy simples. El objeto de esta sección es denir el marco co-rotacional, siguiendo el enfoque de Criseld [Cri91], que dene las relaciones entre las expresiones locales y globales del vector de fuerzas internas y la matriz de rigidez tangente representadas como y Kt en la expresión (11.14). las coordenadas de los nodos son q (x1 , z1 ) y (x2 , z2 ), 1 y 2 en el sistema de coordenadas global (x, z) siendo el vector de desplazamientos globales: pg = (u1 w1 θ1 u2 w2 θ2 )T (11.15) El vector de desplazamientos locales se dene según 11.4: pl = (ul θl1 θl2 )T calculando las componentes de pl (11.16) como sigue: ul = ln − lo (11.17) θl1 = θ1 − α (11.18) θl2 = θ2 − α (11.19) denotando lo y ln en la (11.17) las longitudes inicial y actual del elemento, y α con: lo = ((x2 − x1 )2 + (z2 − z1 )2 )1/2 (11.20) ln = ((x2 + u2 − x1 − u1 )2 + (z2 + w2 − z1 − w1 )2 )1/2 (11.21) la rotación rígida, calculándose como: sin α = cos βo sin β − sin βo cos β (11.22) cos α = cos βo cos β − sin βo cos β (11.23) CAPÍTULO 11. 1 (x2 − x1 ) lo 1 so = sin βo = (z2 − z1 ) lo 1 c = cos β = (x2 + u2 − x1 − u1 ) ln 1 s = sin β = (z2 + w2 − z1 − w1 ) ln co = cos βo = siendo α, 92 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES tal que (11.24) (11.25) (11.26) (11.27) |α| < π : α= sin−1 (sin α) si sin α ≥ 0 y cos α ≥ 0 (11.28) α= cos−1 (cos α) si sin α ≥ 0 y cos α < 0 (11.29) α= sin−1 (sin α) si sin α < 0 y cos α ≥ 0 (11.30) si sin α < 0 y cos α < 0 (11.31) α = − cos−1 (cos α) Para aplicar el principio de los trabajos virtuales es preciso obtener los desplazamientos locales virtuales derivando las ecuaciones (11.17), (11.18) y (11.19): δul = cos β (δu2 − δu1 ) + sin β (δw2 − δw1 ) = = (− cos β − sin β 0 cos β sin β 0) δ pg = rt δ pg δθl1 = δθ1 δα = δθ1 δβ con (α = β − β0 ) (11.33) δθl2 = δθ2 − δα = δθ2 − δβ Asimismo, derivando la ecuación (11.27) es posible obtener δβ = tomando δln = δul (11.32) (11.34) δβ 1 ((δw2 − δw1 )ln − (z2 + w2 − z1 − w1 )δln ) cln2 (11.35) de la ecuación (11.32) y simplicando: 1 ((δw2 − δw1 ) − sc(δu2 − δu1 ) − s2 (δw2 − δw1 )) cln 1 (s − c 0 − s c 0)δ pg = zt δ pg δβ = ln δβ = (11.36) CAPÍTULO 11. 93 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES Aplicando (11.36) en (11.33) y (11.34):  δθl =  001000   − 1  ln 000001  z  pg = AT δpg t z t (11.37) Quedando nalmente la relación de desplazamientos locales y globales como sigue:  δul   δ pl =  δθl1  δθl2   −c −s 0 c s 0        =   −s/ln c/ln 1 s/ln −c/ln 0  δ pg    −s/ln c/ln 0 s/ln −c/ln 1  δ pl =  (11.38)  r  δ pg = Bδ pg t A t (11.39) La relación entre los vectores de esfuerzos internos local y global, ql y qg res- pectivamente, se obtiene igualando los trabajos virtuales,W , en ambos sistemas de referencia según reeja la ecuación (11.40). Dependiendo el vector de esfuerzos internos local de esta relación, qtl = (N, M 1 , M 2 ), de la denición propia del elemento. Wint = δ ptgv qg = N δuv + M 1 δθl1v + M 2 δθl2v = δ ptlv ql = δ ptgv Bt ql Así, para cualquier vector de desplazamientos virtuales de fuerzas internas δ qg δ pg (11.40) arbitrario, el vector queda como sigue: qg = Bt ql Siendo preciso conocer las tensiones (11.43) para el cálculo de δ ql (11.41) resultantes de las ecuaciones (11.42) y δ qg . N= EAul lo (11.42) CAPÍTULO 11. FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES   M1    94  2 1 θ  = 2EI    l1  lo M2 1 2 θl2 (11.43) No obstante, la ecuación (11.42) asume que la deformación axil del elemento es igual al cociente de la deformación relativa entre sus extremos y la longitud del elemento recto. Tal aproximación no permite considerar cualquier otro tipo de deformación sobre un elemento inicialmente recto, por ejemplo la resultante de su exión (ver gura 11.4). Este efecto puede ser tenido en cuenta incluyendo los términos de segundo orden en la deformación de Green relativa al sistema co-rotacional para un elemento inicialmente recto, quedando la deformación local como sigue: dul 1 + εxl = dxl 2  dul dxl 2 1 + θl2 2 (11.44) Deniendo el cambio de base isoparamétrico (11.45), con zamiento local ul (ξ) ξ ∈ (−1 1), el despla- puede expresarse como: 1 xl = (1 + ξ)lo 2 (11.45) 1 ul (ξ) = (1 + ξ)ul 2 (11.46) Diferenciando la ecuación (11.46) se obtiene la deformación local εxl = dul dul dξ ul = = dxl dξ dxl lo (11.47) Asimismo, deniendo el desplazamiento transversal local mediante un polinomio cúbico de ξ wl , ver gura 11.4, según muestra (11.48), es posible determinar el giro (11.49) nuevamente mediante diferenciación.  wl (ξ) = 2 t   lo  (ξ − 1)(ξ − 1)   θl1  8 (ξ 2 − 1)(ξ + 1) θl2 (11.48) CAPÍTULO 11. 95 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES  θl (ξ) = 3ξ 2 − 2ξ − 1 t dwl dwl dξ lo  θl = st θl = =  dxl dξ dxl 4 2 3ξ + 2ξ − 1 (11.49) Con la ayuda de las ecuaciones (11.47) y (11.49) es posible expresar (11.44) como (11.50): ul 1 + εxl (ξ) = xl 2  ul xl 2 1 + θtl sst θl 2 (11.50) Asumiendo deformación constante, el último término de (11.50) puede ser modicado por su valor medio, quedando: ul 1 εxl (ξ) = + lo 2  ul lo 2 1 t θ + 2lo l Z sst dxl θl (11.51) Realizando el cambio de variable (11.45) e integrando en (11.51), se obtiene: ul 1 εxl (ξ) = + lo 2  ul lo 2  2 2 2 2 2  1  (3ξ − 2ξ − 1) (3ξ − 1) − 4ξ  lo dξθl 2 2 2 2 2 2 2 −1 4 (3ξ − 1) − 4ξ (3ξ − 2ξ − 1)    2 4 −1 ul 1 u l 1  θl = + + θtl  (11.52) lo 2 lo 60 −1 4 + 1 t θ 2lo l Z 1 Derivando (11.52) para su inclusión en el trabajo virtual, con la ayuda de (11.32) y (11.37), la variación de la deformación queda:  4 −1  δul ul δul 1  δθl + + 2θtl  2 lo xl 60 −1 4     4 −1 1 ul 1  At δ pg = 1+ rt δpg + θtl  lo lo 30 −1 4 δεxl (ξ) = (11.53) Con lo que nalmente, con la inclusión de los términos de segundo orden en la deformación de Green, la ecuación (11.39) queda modicada como sigue: CAPÍTULO 11.  δ pl 96 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES l δε  o xl  =  δθl 1  δθl 2         =   1+ ul lo  lo t  rt δpg + 30 θl 4 −1 At −1    At    δ pg = Bδ pg 4  (11.54) Derivando ahora (11.54) y deniendo la matriz de rigidez tangente δ qg = Kg δ pg , como se llega a la ecuación (11.55): δ qg = Bt δ pl + N δ B1 + M 1 δ B2 + M 2 δ B3 = Kgt1 δ pg + Kgtσ δ pg donde Kg (11.55) B2 , por ejemplo, corresponde a la segunda la de B (ver (11.38) y (11.39)). Asumiendo un comportamiento lineal del material y derivando las ecuaciones (11.42) y (11.43) se obtiene:  δN   1 0 0     EA      2 2 =  δM 1   0 4r 2r  δ pl = Cl δ pl lo     δM 2 0 2r2 4r2 donde r (11.56) es el radio de giro de la sección. Así, empleando la ecuación (11.56) se obtiene el término de la matriz de rigidez tangente estándar Kgt1 de la ecuación (11.55): Kgt1 = Bt Cl B (11.57) A su vez, la matriz de rigidez geométrica se obtiene de los tres últimos términos de (11.55). De forma que derivando la primera columna de B se obtienen los siguientes términos: δ B1 =  1+ ul lo   δ rt +  δul t lo t  4 −1  t r + δθl A + lo 30 −1 4    lo t  4 −1   δ B2  θ 30 l −1 4 δ B3 (11.58) CAPÍTULO 11. 97 FORMULACIÓN DEL PROBLEMA DINÁMICO EN CABLES con ayuda de las ecuaciones (11.32) y (11.36) y observando que gura 11.4) y que δβ = δα (ver δ B2 = δ B3 1 ln r = δβ z = zzt δpg δ B2 = (11.59) 1 1 1 δ z + 2 zδul = 2 (rzt + zrt )δ pg ln ln ln (11.60) con (11.60),   lo lo t  4 −1  (δ B2 δ B3 )t = (4θl1 − θl2 − θl1 + 4θl2 )(δ B2 δ B3 )t = θl 30 30 −1 4 lo 1 (θl1 + θl2 ) 2 (rzt + zrt )δ pg(11.61) 10 ln y por otra parte, con (11.37):  4 lo t  δθ 30 l −1 −1   4 −1   At = lo δpt A   At = 30 4 −1 4   0 0 0  lo t   B  0 4 −1   Bδpg 30   0 −1 4 (11.62) Por último, sustituyendo en (11.55) y simplicando, la expresión resultante de la matriz de rigidez tangente es:  Kgtσ 0 0 0   N lo t  N (1 + ul /lo ) t rrt  = B  0 4 −1  zz + N + B + 30 ln lo   0 −1 4 1 1 (M + M + N lo (θl1 + θl2 ))(rzt + zrt ) 1 2 2 ln 10 (11.63) Capítulo 12 Formulación del contacto catenaria-pantógrafo Los algoritmos matemáticos que abordan el problema del contacto son generalizables a un número indeterminado n de cuerpos. Sin embargo, para simplicar el proceso, este epígrafe se limitará al caso en que dos cuerpos interaccionan. La terminología seguida es la misma que recoge Ted Belytschko en el capítulo de [BLM00] en el que se trata en profundidad el método del ¡a penalty. ¡b na nb ¡c ­b ­a Figura 12.1: Problema de contacto generalizado 98 Contact-Impact CAPÍTULO 12. 99 FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO En la gura 12.1 se muestra el esquema que se empleará para analizar el problema del contacto. Para los cuerpos en contacto Ωa y Γa y na Ωb , Γb . cuerpos: y nb a y b, cuyos dominios se denotan por representan los vectores normales a sus respectivas supercies Asimismo, Γc hace referencia a la supercie de contacto común a ambos Γc = Γa ∩ Γa . Desde el punto de vista matemático, la formulación del contacto es independiente del cuerpo de referencia tomado, no obstante se designará un cuerpo maestro o master y otro esclavo o la formulación subsiguiente al cuerpo maestro slave, reriendo por comodidad toda a, tomando el cuerpo b como esclavo. En general, la supercie de contacto será una función temporal y su determinación es una de las dicultades mayores que presenta el análisis matemático del problema de contacto. Aunque las supercies Γa y Γb que denen el contacto no siempre sean numéricamente coincidentes, se hará referencia a la supercie de contacto Γc como una entidad única y solidaria al cuerpo maestro. Habiendo denido las supercies de contacto, los desplazamientos en los puntos de las mismas puede descomponerse como ua = uaN · eaN + uaT · eaT siendo uN de contacto y , ub = ubN · ebN + ubT · ebT el módulo de la desplazamiento en la dirección normal a la supercie uT el módulo de el desplazamiento en la dirección tangencial a la supercie de contacto de los cuerpos referidos por el superíndice a o b, lo cual reeja con mayor claridad la gura 12.2. Para obtener el desplazamiento pos a y b, (12.1) u de los cuer- se deben sumar vectorialmente los desplazamientos normal y tangencial multiplicados por los vectores unitarios, eT y eN , para cada una de estas direcciones. Además de los principios generales de la mecánica (conservación de masas, momentos y energías), los cuerpos de la gura 12.2 han de satisfacer una condición adicional de impenetrabilidad en virtud del contacto entre las supercies de los mismos. Es decir, dado que teóricamente ningún subdominio puede pertenecer simultáneamente a ambos cuerpos, se ha de cumplir que Ω a ∩ Ω b = ∅. En general esta condición de impenetrabilidad es altamente no lineal y no puede ser expresa- CAPÍTULO 12. 100 FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO ­b uan ubn ¡c ua ub uaT ­a ubT Figura 12.2: Sistemas de referencia locales da como una ecuación algebraica o diferencial en términos del desplazamiento. Por otra parte, tampoco la predicción de la región de contacto es trivial cuando existe arbitrariedad en los movimientos, lo cual incrementa notablemente la complejidad del problema. Sin embargo, es posible expresar la condición de impenetrabilidad de forma incremental para cada etapa del proceso aplicándola a aquellas partes de los cuerpos que estén en contacto, Γc , aplicando que γN = (ua + ub ) · na ≡ uaN − ubN ≤ 0 Suponiendo que los cuerpos interpenetración nula, γN = 0, a y b en Γc (12.2) están inicialmente en contacto con tasa de la ecuación (12.2) condiciona a éstos a permanecer juntos o separarse, en modo alguno solaparse. Esto es, en caso de permanencia del contacto los cuerpos se desplazarán de forma solidaria satisfaciendo la igualdad cinemática. Pero en caso de pérdida de contacto el desplazamiento del cuerpo mayor que el de a, b será produciéndose como resultado la pérdida de contacto prescrita. Cuando la expresión anterior se hace cumplir para todos los puntos de la región de contacto y de manera continua en el tiempo, entonces la condición de impenetrabilidad se cumple de manera exacta. Sin embargo, cuando esta expresión sólo CAPÍTULO 12. FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO 101 se evalúa en instantes discretos de tiempo, como sucede en la gran mayoría de métodos numéricos, entonces no puede decirse que la condición de impenetrabilidad se cumpla de manera estricta ya que el paso de tiempo empleado en la resolución puede enmascarar interpenetraciones de puntos cercanos entre sí en los distintos instantes considerados. Conviene advertir que la expresión (12.2) puede introducir discontinuidades en la evolución temporal de los desplazamientos, ya que antes del contacto los desplazamientos en la región de contacto de ambos sólidos serán distintas, igualándose a partir de la aparición del contacto. Este fenómeno introduce complicaciones a la hora de proceder a la integración de las ecuaciones que rijan el movimiento de los cuerpos involucrados en el contacto. Igualmente es importante resaltar que la condición de impenetrabilidad sólo es recomendable para puntos que estén en contacto o en disposición de estarlo próximamente, dado el carácter no integrable de la tasa de interpenetración en caso contrario. Por último, considerando la hipótesis de ausencia de fricción, no hay restricción alguna que regule los desplazamientos tangenciales en Γc . Además de las condiciones cinemáticas enunciadas en el anterior epígrafe, otras cinéticas han de satisfacerse en la región de contacto. Concretamente, la suma de los esfuerzos existentes en la región de contacto ha de ser nula, ecuación (12.3). ta + tb = 0 Descomponiendo los esfuerzos en la supercie de contacto (12.3) Γc de manera análoga a la efectuada para desplazamientos en la ecuación (12.1), se obtiene (12.4): ta = taN · eaN + taT · eaT , tb = tbN · ebN + tbT · ebT (12.4) Por tanto, con las ecuaciones (12.3), (12.4) y adoptando la hipótesis de ausencia de adhesión (por la que las tensiones normales en la región de contacto sólo podrían ser de compresión), aplicando la condición de equilibrio sobre las componentes normal y tangencial se obtienen las ecuaciones (12.5) y (12.6) respectivamente. CAPÍTULO 12. FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO 102 taN + tbN = 0 en Γc (12.5) taT + tbT = 0 en Γc (12.6) Las condiciones cinemáticas y cinéticas anteriormente enunciadas ((12.2), (12.5) y (12.6)) pueden ser combinadas en una única ecuación, la condición unitaria de contacto (12.7). tN · γn = 0 (12.7) Ésta expresa que las fuerzas de contacto no realizan trabajo ya que: cuando tN > 0 cero, la tasa de interpenetración es nula; y cuando ésta toma valores distintos de tN = 0 por no existir contacto. Tal como se ha anticipado con la condición de impenetrabilidad, en principio no es posible la interferencia entre las supercies de dos cuerpos en contacto. Sin embargo, esta condición puede ser demasiado brusca para ser implementada en métodos numéricos de cálculo, fundamentalmente debido a las discontinuidades inducidas en la evolución temporal de los desplazamientos y a las dicultades de integración que presenta γN . Por este motivo es común relajar la condición permitiendo un cierto solapamiento entre las supercies en contacto, ver gura 12.3. Deniendo la interpenetración como la mínima distancia entre el punto y la supercie del cuerpo punto situado en b b, Γb , la distancia que separa el punto P de cualquier otro viene dada por la siguiente expresión p lP b = xb (eb , t) − xP (ea , t) = (xb − xP )2 + (y b − y P )2 + (z b − z P )2 donde ea y eb P ∈ Γa (12.8) hacen referencia al sistema local de coordenadas situado en la supercie de los cuerpos a y b respectivamente. La interpenetración será por tanto, según reeja la ecuación (12.9), la mínima distancia entre P y la supercie de b CAPÍTULO 12. 103 FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO ­b eb P gN(xP) > 0 ea ­a Figura 12.3: Penetración cuando P esté dentro del cuerpo b. Si P estuviera fuera del cuerpo b no existiría interpenetración ya que no se daría contacto.   xb (eb , t) − xP (ea , t) gN (e, t) =  0 si  b  x (eb , t) − xP (ea , t) · na ≤ 0 (12.9) en otro caso Es importante observar que el punto de mínima distancia es la proyección ortogonal desde el punto P a la supercie del cuerpo b. No obstante, la conjetura es válida cuando la geometría de los sólidos en contacto es suave, pues en caso de producirse interferencia con algún cuerpo de contorno anguloso la máxima interpenetración ya no coincidirá con la proyección ortogonal desde el punto P a la supercie de b, ver CAPÍTULO 12. 104 FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO por ejemplo la gura 12.4. ­b eb P ea ­a Figura 12.4: Penetración en arista Habiendo relajado la condición de impenetrabilidad de forma que se admita un cierto grado de solapamiento entre los sólidos involucrados en el contacto para evitar restricciones demasiado bruscas en los métodos numéricos, es preciso asociar una fuerza normal que dependa de la interpenetración de los cuerpos, gN , de forma que cuanto mayor sea la interpenetración la fuerza asociada a ésta aumente también. Según este tipo de penalización y de acuerdo con la deducción hecha por Belytschko en [BLM00] para el método del penalty, las fuerzas normales en la supercie de contacto pueden expresarse como (12.10): taN + p = 0 , tbN − p = 0 (12.10) Entre las diversas posibilidades de denición de la fuerza de penalización p in- troducida por la interferencia entre cuerpos, la más general es la siguiente p = (β1 gN + β2 γN ) · H(gN + γN ) siendo H la función escalón de Heaviside (12.11) CAPÍTULO 12. FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO   1 H(x) =  0 si 105 x≥0 (12.12) en otro caso Las ecuaciones que describen la dinámica de un sistema discreto en el que se haya implementado el método del penalty pueden expresarse según la forma deducida en [BLM00]: Mu¨ + q + fc = f En esta expresión (12.13), discreto; q M hace referencia a la matriz de masas del sistema a las fuerzas internas del sistema; fc sobre el sistema y; nalmente, penalty. (12.13) f a las fuerzas externas que actúan introduce el efecto de las fuerzas de contacto o De manera general, este último término puede ser expresado del siguiente modo fc = β1 H(gN )GT Gx + β2 H(γN )GT Gx˙ (12.14) Sin embargo, por la condición de impenetrabilidad (12.2), la penalización (12.14) puede reducir su dependencia únicamente a la interpenetración fc = βH(gN )GT Gx gN según (12.15) Por tanto, en el caso particular del contacto dinámico catenaria-pantógrafo, la interpenetración entre ambas supercies será equivalente al de un resorte de constante β intercalado entre la catenaria, xc , y el pantógrafo, xp , según el vector G que dene el contacto como se expresa en (12.16):  gn = Gx =  1 −1  xp    = xp − xc (12.16) xc Por tanto puede decirse que este método presenta dos importantes ventajas: por un lado su implementación es muy simple, dado que se reduce a insertar un resorte CAPÍTULO 12. FORMULACIÓN DEL CONTACTO CATENARIA-PANTÓGRAFO de característica β 106 entre los sólidos que puedan estar en contacto y; por otro lado, el segundo aspecto favorable radica en la no introducción de variables adicionales. Es importante observar que las fuerzas de contacto no se calculan de manera exacta al emplear el método del mático, el parámetro β penalty. Desde un punto de vista puramente mate- ha de tender a innito para que la aproximación sea lo más precisa posible, ya que el efecto de la penalización consiste en relajar la condición de impenetrabilidad. Cuanto mayor sea β , más cercano a la realidad será el compor- tamiento del sistema. Sin embargo, valores excesivamente elevados pueden originar serios problemas de condicionamiento de la matriz de rigidez, lo cual afecta de forma signicativa a la convergencia del método numérico empleado en la resolución del sistema de ecuaciones diferenciales obtenido. Esto constituye el principal inconveniente del método, ya que la solución obtenida depende enormemente del valor de penalización elegido. Para solventar este problema se han desarrollado numerosos algoritmos de actualización del penalty, de forma que en cada paso de tiempo se em- plee el valor más apropiado logrando un compromiso entre exactitud de la solución y facilidad de resolución numérica. En [Cha02] puede consultarse por ejemplo una propuesta de penalty adaptativo en función del valor que tome la interpenetración gN . taN = −p = −φ(gna )βna   ga    n a a )2 (gn φ(gna ) = + g2n + β4 4β     0 (12.17) si gna ≤ −β si kgna k < β en otro caso (12.18) Capítulo 13 Integración temporal Un modelo de elementos nitos de dinámica estructural discretiza las ecuaciones diferenciales que representan el comportamiento de un sistema físico. Éste tipo de modelos simula exactamente el comportamiento de los modos de baja frecuencia del sistema. Sin embargo, su contenido de alta frecuencia suele ser dominado por la inuencia numérica resultante de la discretización de los elementos nitos. Debido a esto, el sistema de ecuaciones es sti, concepto denido y ampliamente discutido por Hairer y Wanner en [HW91]. Por tanto, la discretización en el tiempo y la denición de un método de integración es un aspecto de enorme relevancia en la simulación dinámica de estructuras exibles. El formato tradicional del sistema de ecuaciones de partida que se encuentra en la bibliografía sobre algoritmos de integración temporal es de primer orden y se rige por la expresión general 13.1: y˙ = f (y, t) siendo (13.1) f , en general, una función no lineal de las variables (y(t),t). No obstante, el sistema de ecuaciones diferenciales que gobierna dinámica de un sistema se expresa en forma matricial compacta: ˙ t) M · u¨ = q (u, u, 107 (13.2) CAPÍTULO 13. 108 INTEGRACIÓN TEMPORAL donde el vector q y la matriz M son en general funciones no lineales de (u, u˙ , t). Es sencillo comprobar que el sistema 13.2 puede plantearse como un caso particular del 13.1 sin más que redenir el vector incógnita, mediante el cambio: y= nuo u˙ = f (y, t) y˙ = −1 M · q (u, u˙ , t)   tal que u˙ un método numérico de integración temporal calcula secuencialmente valores de la incógnita {yn /n = 0, 1, . . . , N } El cálculo de espaciados mediante intervalos hn+1 = tn+1 − tn . y en cada instante se realiza a partir de su valor en k instantes previos de tiempo. El número k es lo que se conoce como el número de pasos del integrador. En consecuencia, el valor de yn+k se calcula a partir de los k valores anteriores de y con la expresión general [Lam91]: k X αj yn+j = hφf (yn+k ,yn+k−1 ,...,yn ,tn+k ;hn+k ) (13.3) j=0 donde se ha usado el subíndice zar que la dependencia de φ con f en la función del segundo miembro para enfati- yn+k , . . . , yn , tn+k ; hn+k
es a través de la función f (y, t). Pueden clasicarse los métodos dados por la expresión 13.3 en base a diferentes criterios: Métodos de un paso / Métodos multipaso. En los de un paso, el cálculo de yn+1 k = 1. Es decir, se realiza haciendo uso exclusivamente de información del paso anterior. En los multipaso se requiere, sin embargo, un número k > 1 de valores iniciales para arrancar, utilizando habitualmente información de los valores anteriores. En cambio, suelen tener una estructura más simple que los de un paso para una precisión similar. Entre los métodos de un paso destacan la familia Newmark, cuya formulación original fue propuesta en [nN59], y otros algoritmos desarrollados a partir de la formulación de Newmark como el Hilber-Hughes-Taylor [HHT77] y el α-generalizado de Chung y Hulbert [CH93]. Otros métodos como los Runge-Kutta son también ampliamente re- CAPÍTULO 13. 109 INTEGRACIÓN TEMPORAL conocidos, aunque esta proyecto no profundizará en ellos. Tampoco en otros multipaso como por ejemplo los métodos lineales, en los que la expresión 13.3 es lineal en  yn+j , f yn+j , tn+j ; j = 0, 1, . . . , k
. Métodos explícitos / Métodos implícitos. En un método explícito, la expresión 13.3 permite despejar yn+1 conocidos los valores yn+j ; j = 0, 1, . . . , k − 1  ; si esto no es posible, el método es implícito. Es preciso destacar que la estabilidad de un método explicito sólo puede ser garantizada si el paso de tiempo elegido es sucientemente pequeño respecto a las frecuencias naturales del sistema. Mientras que los métodos implícitos comúnmente empleados son incondicionalmente estables (o A-estables), lo cual signica que la solución numérica será estable cualquiera que sea el contenido en frecuencias del sistema mecánico. Esta propiedad se antoja muy deseable en la simulación de sistemas sti, siendo la mayor complejidad computacional que requieren los métodos implícitos el precio que hay que pagar por ello. Esta clasicación no es internamente excluyente entre sí, puesto que tanto los métodos de un paso como los multipaso pueden ser explícitos o implícitos, e incluso pueden construirse algoritmos que combinen un método explícito con uno implícito (algoritmos predictor-corrector). De acuerdo con Hughes en [Hug87], puede decirse que un método de integración temporal para dinámica estructural debería cobinar las siguientes propiedades: estabilidad incondicional para sistemas lineales, no más de un sistema de ecuaciones implícitas a resolver en cada paso, precisión de segundo orden, control de la disipación numérica en los modos de alta frecuencia e inicialización autónoma. En los siguientes epígrafes se presenta el algoritmo original de Newmark [nN59] y la extensión α-generalizado propuesta por Chung y Hulbert en [CH93] así como las adaptaciones particulares desarrolladas en esta proyecto para al cálculo dinámico de la interacción catenaria-pantógrafo. CAPÍTULO 13. 110 INTEGRACIÓN TEMPORAL 13.1. La familia β -Newmark De amplia utilización en la dinámica estructural, la familia β -Newmark está especialmente diseñada para la resolución de sistemas de segundo orden, por lo que se aplica a las ecuaciones con el formato dado en 13.2. El término `familia' viene del hecho de que su planteamiento más general tiene dos parámetros (γ y β) cuya variación genera todos los distintos métodos de la familia. La formulación clásica proporciona la posición y la velocidad en el instante (n+1) a partir de la posición y velocidad en (n) y de la aceleración en (n + 1). Partiendo del desarrollo en serie de Taylor de desplazamientos y velocidades respecto al paso de tiempo h se obtienen las expresiones 13.4: un+1 = un + hu˙ n + h2
¨ n + h2 β u ¨ n+1 −β u 1 2 (13.4) u˙ n+1 = u˙ n + h (1 − γ) u¨ n + hγ u¨ n+1 siendo γ y β los parámetros numéricos que dan lugar a los distintos métodos. Mediante la combinación de estos parámetros puede observarse que en régimen lineal los métodos de la familia en los que β≥  1 16 + γ 2 +γ 4  y γ≥ 1 son incondicionalmente 2 estables. Los métodos más representativos obtenidos para distintos valores de γ y β son presentados por Geradin y Rixen en [GR97] con el estudio de un sistema `patrón'. Como conclusión general del análisis de estos métodos puede decirse que todos salvo la regla trapezoidal y la regla trapezoidal modicada tienen un interés limitado, ya que son inestables o condicionalmente estables incluso en el régimen lineal. La regla trapezoidal se obtiene para γ= 1 2 y β= 1 , y en el contexto de 4 la familia de Newmark se le llama también método de aceleración media constante, ya que las expresiones 13.4 se pueden interpretar como actualizaciones en posición y velocidad suponiendo una aceleración media constante entre tn y tn+1 . En el régimen lineal, éste es además el método absolutamente estable más preciso. No obstante, puede introducirse amortiguamiento numérico en la formulación según: CAPÍTULO 13. 1 γ = +α 2 donde α 111 INTEGRACIÓN TEMPORAL y 1 β= 4  1 γ+ 2 2 α≥0 (13.5) es el parámetro de amortiguamiento numérico. Éste es el llamado mé- todo de aceleración constante modicado o regla trapezoidal modicada. Permite aumentar el amortiguamiento numérico en el sistema manteniendo la condición de estabilidad en el algoritmo de integración, lo cual en ciertos casos puede ser muy útil aún a costa de una degradación de la precisión. Géradin y Rixen presentan una sistematización del algoritmo de Newmark para dinámica estructural en [GR97]. La resolución numérica sigue un algoritmo predictor-corrector partiendo de los valores un , u˙ n y u¨ n del instante tn , para cuya predicción inicial se asume aceleración nula: u¨ 0n+1 = 0 u˙ 0n+1 = u˙ n+1 + h (1 − γ) u¨ n u0n+1 = un + hu˙ n + h2 1 2 (13.6) ¨n −β u
llegando a satisfacer la formulación de Newmark 13.4 mediante las correcciones iterativas: u ∆¨ un+1 = 1 ∆ n+1 βh2 ∆u˙ n+1 = γ ∆ n+1 βh u En el caso de la interacción dinámica catenaria-pantógrafo, la corrección (13.7) ∆un+1 en cada paso de integración se calcula a través de la ecuación 13.8: (i) Kt∗(i) ∆u(i+1) n+1 = Rn+1 cuyos cálculos previos del residuo (13.8) R y la matriz de rigidez tangente consistente K∗t se expresan en 13.9: R = Mu¨ n+1 + Cu˙ n+1 + q − f K∗t = Kt + βhγ C + βh1 2 M (13.9) CAPÍTULO 13. donde estático, Kt C 112 INTEGRACIÓN TEMPORAL es la matriz de rigidez tangente calculada para un problema cuasi- es la matriz de amortiguamiento y M la matriz de masas. 13.2. El método α-Generalizado La integración de ecuaciones algebraicas-diferenciales (DAE) puede conducir a inestabilidad numérica cuando se usa un método de integración de la familia de Newmark debido a las restricciones algebraicas, que se maniesta a través de oscilaciones crecientes en la respuesta en aceleraciones. Introduciendo una pequeña disipación en el algoritmo para las altas frecuencias se logra controlar esta inestabilidad, manteniendo la estabilidad de la integración de dinámica lineal con restricciones. Destaca en este aspecto el método α-generalizado descrito por Chung y Hulbert en [CH93]. Éste incluye como casos particulares algunos de los algoritmos de integración temporal más importantes en dinámica estructural, como el algoritmo Hilber-Hughes-Taylor [HHT77], constituyendo así un marco general para investigaciones teóricas. El método se fundamenta en las fórmulas de Newmark 13.4, aunque para el cálculo del residuo este algoritmo promedia la diferente contribución de dos instantes consecutivos según los parámetros numéricos αm y αf como reeja 13.10: R = (1 − αm ) M u¨ n+1 + αm M u¨ n + (1 − αf ) C u˙ n+1 +
αf C u˙ n + (1 − αf ) q n+1 − f n+1 + αf (q n − f n ) En particular, el algoritmo Hilber-Hughes-Taylor se obtiene para   αf ∈ 0, 31 . No obstante, estos parámetros del método calculados en función del radio espectral αm = Deniendo en 13.12: 2ρu∞ − 1 ρu∞ + 1 α-generalizado (13.10) αm = 0 y pueden ser ρq∞ : y ρu αf = u ∞ ρ∞ + 1 (13.11) αf m = αf −αm , los parámetros de Newmark quedan como se muestra CAPÍTULO 13. 113 INTEGRACIÓN TEMPORAL 1 γ = + αf m 2 y 1 β= 4 Demostrándose que incluso para valores de  1 γ+ 2 αf m > 0 2 (13.12) el método presenta una precisión de segundo orden. La solución numérica se obtiene mediante un algoritmo predictor-corrector como el empleado en la familia Newmark, quedando la matriz tangente aumentada en cada paso corrector como reeja la ecuación 13.13, análoga a la 13.9: Kt∗ = (1 − αf ) Kt + (1 − αf ) 1 γ C + (1 − αm ) 2 M βh βh (13.13) Capítulo 14 Validación con la norma EN50318 La norma EN50318 [CEN99] fue aprobada el 1 de Abril de 2002. Esta versión de la norma europea ha sido realizada por CENELEC a petición de la comisión europea siguiendo la normativa de interoperatibilidad 96/48/EC. Para obtener el certicado de comprobación EC para los elementos constituyentes de una línea de contacto de transmisión de energía se necesita una simulación con un programa que debe ser validado previamente por la norma EN50318. Como especica la norma EN50318, capítulo 11,el primer paso de la validación de un código debe ser la comparación con un modelo de referencia, ver Figura 14.1, para tener conanza en la precisión de la simulación. Si los resultados están dentro del rango señalado en la norma EN50318, Tabla 14.1, el método de simulación puede usarse para validar el código mediante la comparación con resultados medidos en líneas de alta velocidad. Para ello, las medidas deben haberse realizado de acuerdo a la norma EN50317. La tabla 14.1 muestra el rango de resultados admisibles en la simulación de la catenaria EN50318, así como los resultados obtenidos con el modelo presentado en este proyecto. Además, de acuerdo con la norma EN50318, sección 3.17 y apéndice A.3, el máximo desplazamiento se debe calcular en los apoyos de los vanos 5 y 6. Los resultados obtenidos se muestran a continuación en las guras 14.2 y 14.3 a 250 km/h y en las guras 14.4 y 14.5 a 300 km/h. 114 CAPÍTULO 14. 115 VALIDACIÓN CON LA NORMA EN50318 1.8 1.6 1.4 1.2 Height [m] 1 0.8 0.6 0.4 0.2 0 −0.2 0 60 120 180 240 300 Length [m] 360 420 480 540 600 Figura 14.1: Catenaria de referencia EN50318 (10 vanos) En cada gura se representan conjuntamente los cuatro vanos centrales con la fuerza y los desplazamientos superpuestos en cada gráca. Por ello, se ha aplicado la misma escala en el eje de abscisas mientras la misma varía en el eje de ordendas para poder observar los resultados sin problemas. Para terminar, remarcar que todas las unidades usadas siguen el sistema internacional: longitudes y desplazamientos en [m], fuerzas de contacto en [N] y tiempo en [s]. Los resultados de la simulación se encuentran dentro del rango admisible cumpliendo los requisitos de la norma EN50318. CAPÍTULO 14. 116 VALIDACIÓN CON LA NORMA EN50318 Velocidad Rango de frecuencias 20Hz Fuerza media de contacto Desviación típica (σ) Máx. estad. de fuerza Mín. estad. de fuerza Máx. fuerza de contacto Mín. fuerza de contacto Máx. despl. del apoyo 5 Máx. despl. del apoyo 6 Máx. despl. del apoyo 7 % de pérdida de contacto 250 km/h Rango Admis. Simulación 110 N - 120 N 116.07 N 26 N - 31 N 27.38 N 190 N - 210 N 198.20 N 20 N - 40 N 33.91 N 175 N - 210 N 177.57 N 50 N - 75 N 60.05 N 48 mm - 55 mm 53.4 mm 48 mm - 55 mm 51.9 mm 48 mm - 55 mm 53.0 mm 0% 0% 300 km/h Rango Admis. Simulación 110 N - 120 N 115.35 N 32 N - 40 N 33.64 N 210 N - 230 N 216.27 N -5 N - 20 N 14.44 N 190 N - 225 N 210.35 N 30 N - 55 N 40.87 N 55 mm - 65 mm 62.7 mm 55 mm - 65 mm 63.1 mm 55 mm - 65 mm 62.0 mm 0% 0% Tabla 14.1: Validación con el modelo de referencia Time [s] 2.592 3.024 3.456 3.888 4.32 4.752 5.184 5.616 0.08 1.9 0.06 1.6 0.04 1.3 0.02 1 0 0.7 −0.02 0.4 −0.04 Mast 5 Mast 6 Mast 7 Pantograph 0.1 −0.2 180 210 240 270 300 Distance [m] 330 360 390 Uplift [m] Height [m] 2.16 2.2 −0.06 −0.08 420 Figura 14.2: Geometría y desplazamiento en los vanos centrales a 250 km/h Height [m] 180 2.5 117 VALIDACIÓN CON LA NORMA EN50318 210 240 270 300 330 360 390 420 180 2.2 160 1.9 140 1.6 120 1.3 100 1 80 0.7 60 0.4 40 0.1 20 Contact force [N] CAPÍTULO 14. Contact force −0.2 180 210 240 270 300 Distance [m] 330 360 390 0 420 Figura 14.3: Geometría y fuerza de contacto en los vanos centrales a 250 km/h Time [m] 2.16 2.52 2.88 3.24 3.6 3.96 4.32 4.68 0.08 1.9 0.06 1.6 0.04 1.3 0.02 1 0 0.7 −0.02 0.4 −0.04 Mast 5 Mast 6 Mast 7 Pantograph 0.1 −0.2 180 210 240 270 300 Distance [m] 330 360 390 Uplift [m] Height [m] 1.8 2.2 −0.06 −0.08 420 Figura 14.4: Geometría y desplazamiento en los vanos centrales a 300 km/h Height [m] 2.3 118 VALIDACIÓN CON LA NORMA EN50318 200 250 300 350 400 220 2 190 1.7 160 1.4 130 1.1 100 0.8 70 0.5 40 0.2 10 Contact force [N] CAPÍTULO 14. Contact force −0.1 180 210 240 270 300 Distance [m] 330 360 390 −20 420 Figura 14.5: Geometría y fuerza de contacto en los vanos centrales a 300 km/h Capítulo 15 Conclusiones En esta segunda parte se propone un nuevo método para simulación de la interacción dinámica pantógrafo-catenaria. Las aportaciones más relevantes que se han llevado a cabo se pueden separar en tres aspectos: Se ha implementado el efecto del pandeo en las péndolas mediante la denición de una rigidez variable. Esto conere un tratamiento más realista del sistema de cables abordado. La simulación de la interacción dinámica catenaria-pantógrafo se ha desarrollado con un modelo de elementos nitos empleando una formulación corotacional para los hilos sustentador y de contacto. Para la integración numérica de la evolución dinámica del sistema se ha implementado el método α-generalizado que ha permitido incrementar la estabilidad de la solución ltrando las vibraciones de alta y muy alta frecuencia. Frente a todo lo publicado anteriormente al respecto en la literatura cientíca, este modelo sí ha sido validado mediante lo propuesto en la norma de validación europea EN50318 [CEN99]. Los resultados obtenidos se encuentran dentro de los margenes de aceptación establecidos en dicha norma, por lo que es apto para la certicación ocial de catenarias de alta velocidad. 119 Parte III Reducción Dinámica mediante Física Multicuerpo 120 121 Dado al elevado tiempo de computación del método presentado en la parte anterior, se ha desarrollado un método reducido para obtener soluciones aproximadas de la fuerza de contacto de la interacción catenaria pantógrafo. Este método se apoya en el uso de la física multicuerpo para aplicar una jerarquía variable de modelos a cada cuerpo (vano a vano en este caso) y que depende de si la respuesta dinámica de los cuerpos se comporta de manera lineal o de manera no lineal. Esta parte está estructurada de la siguiente manera: En primer lugar se expone una revisión de las técnicas multicuerpo, 16. A continuación, los capítulos 17, 18 y 19 presentan la formulación que permitirá reducir el tamaño del problema de interacción pantógrafo-catenaria sin perder precisión. El capítulo,20 muestra el resultado de la implementación del método propuesto y presenta los resultados obtenidos. Por último, el capítulo 15 presenta brevemente las conclusiones del trabajo. Las referencias empleadas en el desarrollo del trabajo serán presentadas en orden alfabético al nal del documento. Capítulo 16 Estado del arte Para poder simular un mecanismo compuesto por cuerpos rígidos conectados por uniones cinemáticas deben obtenerse las ecuaciones que denan su movimiento. Existen dos métodos fundamentales para ello: la formulación Newton-Euleriana y la formulación Lagrangiana. Las ecuaciones de Newton-Euler provienen de la aplicación de las leyes del movimiento a cada cuerpo, aplicando el principio de acción y reacción a fuerzas y pares. De acuerdo con el principio de D'Alembert, las reacciones pueden aparecer como fuerzas aplicadas en cada cuerpo y que éstos tan sólo estén relacionados por ligaduras algebraicas. La formulación lagrangiana [DJB94] describe un sistema dinámico en términos de trabajo y energía usando coordenadas generalizadas, e.g. coordenadas relativas o cartesianas. Si las coordenadas generalizadas son independientes, las reacciones y los momentos se eliminan automáticamente, con lo que se deriva un sistema compacto de ecuaciones del movimiento. Para modelar uniones cinemáticas se deben utilizar restricciones denominadas multiplicadores de Lagrange como se estudiará más adelante. Para modelar mecanismos compuestos por cuerpos exibles, debe aplicarse sobre dichos cuerpos alguna técnica de discretización. Como se ha desarrollado en la parte anterior, uno de los métodos más utilizados es la formulación por elementos nitos. Aplicando una formulación lagrangiana se pueden encontrar multiplicado- 122 CAPÍTULO 16. ESTADO DEL ARTE 123 res de Lagrange para enlazar los diferentes cuerpos y modelar mecanismos exibles [GC01]. La formulación moderna de la mecánica multicuerpo permite obtener una representación detallada y able de sistemas mecánicos complejos. No obstante, la consecución de alta precisión tan sólo puede hacerse a costa de algoritmos más sosticados que requieren un mayor esfuerzo computacional. Por ello, Eberhard y Schiehlen presentan en [ES98] un modelado jerárquico de diferentes modelos para conseguir la precisión adecuada en cada cuerpo. Dependiendo de si el comportamiento de un cuerpo exible es lineal o no lineal, los modelos disponibles para analizarlos pueden ser diferentes. Existen técnicas muy maduras de bases reducidas para estructuras con comportamiento lineal. No obstante, para cuerpos no lineales es preferible usar modelos completos de elementos nitos. Una reducción lineal transforma un modelo original en un modelo de orden reducido minimizando la pérdida de precisión. En sistemas dinámicos lineales suelen utilizarse técnicas de descomposición modal que serán presentadas de forma más detallada en este capitulo. De Fonseca [DF00] ha realizado interesantes estudios en este campo y las ha aplicado a la dinámica estructural. Capítulo 17 Frecuencias naturales y modos de vibración Discretizando las ecuaciones del movimiento por elementos nitos en un cuerpo cuyo comportamiento es no lineal se obtiene que Mu ¨ + C u˙ = f − q donde u (17.1) representa los desplazamientos en cada una de los grados de libertad del problema. En problemas con una discretización espacial de tamaño medio, el tamaño de las matrices M ,C y la matriz tangente de q , Kt , puede llegar a ser relativamente grande y los tiempos de ejecución demasiado elevados al realizar análisis dinámicos. Por ello, cuando se puede asumir un comportamiento lineal, se utilizan modelos reducidos que aproximan el fenómeno físico con un menor número de ecuaciones. El método más conocido de obtención de bases reducidas es la superposición modal, ya que, además de reducir el sistema, aporta información muy útil para evitar problemas de resonancia. No obstante, este método tan sólo es válido si el sistema se comporta de una manera lineal y los coecientes de las matrices de masa y rigidez Kt M son constantes. Por ello, si un problema de naturaleza no lineal se abordara con esta técnica, se habría de linealizar la ecuación (17.1) en el entorno del punto de equilibrio u0 . Expandiendo la ecuación mediante una aproximación de 124 CAPÍTULO 17. 125 FRECUENCIAS NATURALES Y MODOS DE VIBRACIÓN Taylor de primer orden Mu ¨ + C u˙ + Kt (u0 ) (u − u0 ) + q (u0 ) = f (17.2) o, lo que es lo mismo, Mu ¨ + C u˙ + Kt (u0 ) u = f − q (u0 ) + Kt (u0 ) u0 M, C y Kt serán válidas siempre que la posición de equilibrio dinámico, demasiado alejada de la posición de equilibrio estático, (17.3) u, no esté u0 . Al analizar las frecuencias naturales y los modos de vibración de un sistema mecánico, realmente se estudia su respuesta libre sin amortiguamiento; esto es, sin que existan fuerzas exteriores aplicadas ni fuerzas dependientes de la velocidad, por lo que las ecuaciones del movimiento quedarán reducidas a Mu ¨ + Kt (u0 ) u = 0 Se asumirá que dichas respuestas son del tipo compleja, u = χϕ eiωt ϕ es un vector de constantes reales y eiωt (17.4) donde χ es una constante indica la notación de Euler para la exponencial compleja. Dicho cambio, desarrollando la parte compleja, quedará expresado como u = αϕ ei(ωt−δ) donde α es una constante real y δ (17.5) un desfase en tiempo. Derivando dos veces la ecuación (17.5) con respecto al tiempo se obtiene que u ¨ = −ω 2 αϕ ei(ωt−δ) = −ω 2 u (17.6) y sustituyendo (17.5) y (17.6) en la ecuación de movimiento (17.4) resulta
−ω 2 M + Kt (u0 ) ϕα ei(ωt−δ) = 0 Eliminando el escalar α ei(ωt−δ) (17.7) se deduce que
−ω 2 M + Kt (u0 ) ϕ = 0 (17.8) CAPÍTULO 17. 126 FRECUENCIAS NATURALES Y MODOS DE VIBRACIÓN lo cual dene un problema de autovalores generalizado. Dicho problema trata de encontrar los vectores ϕ que, tomando λ = ω2, veriquen que Kϕ = λM ϕ (17.9) Para que existan soluciones diferentes a la trivial, el sistema debe ser indeterminado. Por eso, un método para hallar los autovalores del sistema consiste en buscar aquéllos que cumplan det
−ω 2 M + Kt (u0 ) = 0 (17.10) La expresión (17.10) se denomina ecuación característica. Si los problemas son pequeños, se buscan las raíces de dicha ecuación característica con lo que se obtienen los autovalores λi y sus autovectores asociados ϕi . En caso de que los sistemas sean mayores, se utilizan algoritmos que obtienen todos estos autovalores o parte de ellos. De esta forma se obtiene un sistema reducido con tan solo la información de las frecuencias útiles. Los valores tovectores, ϕi , ωi = √ λi se denominan frecuencias propias del sistema y los au- son los modos de vibración asociados a las frecuencias ωi . Éstas representan las frecuencias de vibración de las posibles soluciones armónicas del sistema. 17.1. Frecuencias propias en catenarias ferroviarias Pese a la no linealidad de las estructuras de cables, se realizan a menudo análisis modales para su estudio. Al conocer las frecuencias naturales de vibración se pueden evitar problemas de resonancia muy peligrosos para la estabilidad de cualquier estructura. Para mostrar la forma de los diferentes modos de vibración, se ha realizado un análisis modal sobre la catenaria denida por la norma EN50318 representada en la gura 17.1. La tabla 17.1 presenta la sensibilidad que presentan las frecuencias naturales de dicha catenaria al mallado. El número de elementos se dene como el CAPÍTULO 17. FRECUENCIAS NATURALES Y MODOS DE VIBRACIÓN 127 Figura 17.1: Catenaria denida por la norma EN50318 número de divisiones introducidas entre péndolas en el hilo de contacto. Se observa que el mallado no afecta en las frecuencias más bajas y tan sólo de forma leve en a partir de la décima. Por ello, para recoger el comportamiento lineal a baja frecuencia no será necesario un modelo con una malla muy na. Los cuatro primeros modos se muestran en las guras 17.2,17.3,17.4 y 17.5. Puede comprobarse cómo gran parte de la información dinámica se recoge en relativamente pocos modos de vibración. Figura 17.2: Modo de vibración 1 (1.0182 Hz) CAPÍTULO 17. 128 FRECUENCIAS NATURALES Y MODOS DE VIBRACIÓN 2 4 6 8 10 12 14 16 18 20 Modo 1 1.07 1.07 1.07 1.08 1.09 1.10 1.11 1.13 1.15 1.16 Modo 2 2.14 2.13 2.14 2.16 2.17 2.19 2.22 2.25 2.28 2.32 Modo 3 3.21 3.20 3.21 3.23 3.26 3.29 3.32 3.36 3.41 3.46 Modo 4 4.29 4.27 4.28 4.30 4.34 4.37 4.42 4.47 4.52 4.59 Modo 5 5.37 5.34 5.35 5.38 5.42 5.46 5.51 5.57 5.64 5.71 Modo 6 6.45 6.42 6.42 6.45 6.49 6.54 6.60 6.67 6.74 6.82 Modo 7 7.49 7.48 7.48 7.51 7.55 7.60 7.67 7.73 7.81 7.89 Modo 8 8.46 8.51 8.51 8.54 8.57 8.62 8.67 8.72 8.78 8.83 Modo 9 9.32 9.65 9.64 9.66 9.67 9.71 9.77 9.84 9.91 9.99 Modo 10 10.06 9.88 9.70 9.67 9.72 9.77 9.83 9.90 9.98 10.06 Modo 11 10.72 9.90 9.74 9.71 9.73 9.78 9.86 9.94 10.03 10.13 Modo 12 11.32 9.94 9.82 9.79 9.82 9.87 9.93 10.01 10.09 10.19 N. de elem. Tabla 17.1: Sensibilidad del mallado de las frecuencias naturales Figura 17.3: Modo de vibración 3 (3.0555 Hz) CAPÍTULO 17. FRECUENCIAS NATURALES Y MODOS DE VIBRACIÓN Figura 17.4: Modo de vibración 5 (5.0938 Hz) Figura 17.5: Modo de vibración 7 (7.1341 Hz) 129 Capítulo 18 El método de la superposición modal La linealidad de la ecuación (17.4) implica también la linealidad de sus soluciones. Por lo tanto si u1 y u2 son soluciones de (17.4), νu1 + µu2 también lo será. Por ello, la solución más general al sistema lineal anterior vendrá dada como una combinación lineal de sus soluciones del tipo u = α1 ϕ1 ei(ωt−δ) + α2 ϕ2 ei(ωt−δ) + · · · + αn ϕn ei(ωt−δ) n P αk ϕk ei(ωk t−δk ) = (18.1) k=1 Dicha solución general se puede interpretar como la suma de brando a su frecuencia natural ωk con una amplitud n modos ϕk vi- αk . El sistema completo tendrá tantos modos de vibración como grados de libertad tenga el sistema original. No obstante, cada uno de estos modos contendrá la información relativa a la frecuencia propia con la que esté relacionado. Por ello, es posible reducir el tamaño original del sistema tomando tan sólo las frecuencias que interesen en el estudio concreto que se realice siempre que dicho sistema se comporte de manera lineal. Denominando ϕki a la componente i del vector propio ϕk y siendo j= √ −1, se puede reescribir la ecuación (18.1) como ui = n X αk ϕki ej(ωk t−δk ) k=1 130 (18.2) CAPÍTULO 18. 131 EL MÉTODO DE LA SUPERPOSICIÓN MODAL Un sistema de coordenadas normales asociadas al problema inicial se dene como ξ (t) = αk ej(ωk t−δk ) y, sustituyendo en (18.2), se obtiene lo que se puede interpretar como un cambio de coordenadas: de grados de libertad asociados a los desplazamientos y giros, a la base de coordenadas normales. Esto queda expresado de la siguiente manera ui (t) = n X ϕki ξk (t) (18.3) k=1 donde las componentes ϕki de la ecuación (18.3) constituyen la matriz de cambio de base       φ = (ϕ1 , ϕ2 , · · · , ϕn ) =       ϕ11   ϕ12   .  ..   ϕ1n   ϕ21    ϕ22   .   ..   ϕ2n     · · ·   ϕn1     · · ·   ϕn2    .  .   . .   ..     ··· ϕnn                   (18.4) con la que, expresando matricialmente la ecuación (18.3), se puede escribir u = φξ (18.5) o lo que es lo mismo, u (t) = ξ1 (t) ϕ1 + ξ2 (t) ϕ2 + · · · + ξn (t) ϕn Al aplicar este cambio de coordenadas al sistema de ecuaciones (17.3) se obtiene M φξ¨ + Kφξ = 0 y, premultiplicando por la transpuesta de la matriz de cambio de base, φT , φT M φ ξ¨ + φT Kφ ξ = 0 se obtiene el sistema de ecuaciones del movimiento denido en coordenadas normales ˜ ξ¨ + K ˜ξ=0 M CAPÍTULO 18. 132 EL MÉTODO DE LA SUPERPOSICIÓN MODAL en el que las matrices características del sistema ˜ M y ˜ K son matrices diagonales por las propiedades de ortogonalidad de los modos respecto a las matrices Este cambio de base es aplicable a las fuerzas, Con esto se le asigna a f˜ la aportación de f f, M y K. que actúan sobre el sistema. a cada modo y frecuencia propia. Un sistema con oscilaciones forzadas no amortiguado en el sistema de coordenadas normales quedaría denido por: ˜ ξ¨ + K ˜ ξ = φT f = f˜ M Este cambio, sin embargo, no es aplicable a la matriz de amortiguamiento que no es posible lograr una diagonalización simultanea de M ,C y Kt . C ya Por ello, para añadir el amortiguamiento, se utiliza la hipótesis de Rayleigh, que considera la matriz de amortiguamiento proporcional a la matriz de masas y a la matriz de rigidez. Entonces, siendo ˜ = αK ˜t + β M ˜, C el sistema con oscilaciones forzadas y amor- tiguado en coordenadas normales quedaría expresado como: ˜ ξ¨ + C ˜ ξ˙ + K ˜ ξ = f˜ M (18.6) Utilizando la ecuación (18.6) se puede estudiar el sistema reducido y obtener información del sistema completo usando el cambio de variables (18.4). De este modo se pasa del sistema de ecuaciones (17.1) de tamaño tamaño m, n a un sistema (18.6) de siendo éste el número de modos utilizados y mucho menor que n. 18.1. Condiciones iniciales Para utilizar la superposición modal en el resolución de problemas estructurales es preciso evaluar los 2n coecientes α0 y δ0 de la ecuación (17.5) o, lo que es lo mismo, las condiciones iniciales en coordenadas modales ξ (t = 0) ξ˙ (t = 0). Para (u0 , u˙ 0 ), para poder   ˙ reducidas, ξ0 , ξ0 . La ello se utilizarán las condiciones iniciales del modelo completo, incorporar unas condiciones iniciales al modelo en bases y CAPÍTULO 18. EL MÉTODO DE LA SUPERPOSICIÓN MODAL 133 expresión (17.5) expresada sin notación de Euler es X u= αk ϕk cos (ωk t − δk ) (18.7) k que, desarrollando el coseno de la suma, queda expresada como u= X αk ϕk cosδk cosωk t + X k αk ϕk senδk senωk t (18.8) k Por otro lado, la derivada de la expresión (18.7) es u˙ = X −αk ϕk ωk sen (ωk t − δk ) (18.9) k Estas expresiones, particularizadas para u0 = X u˙ 0 = X t = 0, se transforman en αk ϕk cosδk k αk ϕk ωk senδk (18.10) k y en coordenadas modales u0 = φξ 0 u˙ 0 = φξ˙0 (18.11) Premultiplicando ambas ecuaciones (18.11) por φT M resulta ˜ −1 φT M u0 ξ0 = M ˜ −1 φT M u˙ 0 ξ˙0 = M (18.12) Las expresiones (18.12) permiten estudiar mediante superposición modal un sistema que ya se encuentre en movimiento; es decir, que en lugar de partir de condiciones de reposo, tenga unas condiciones iniciales (u0 , u˙ 0 ). Capítulo 19 Mecánica multicuerpo Un sistema multicuerpo es el resultado de describir un sistema mecánico como la composición de cuerpos sólidos (rígidos o exibles) unidos por enlaces que restringen su movimiento relativo. El estudio de la dinámica multicuerpo es el análisis de cómo se mueven dichos sistemas bajo la inuencia de fuerzas. 19.1. Acoplamiento de modelos físicos Para aplicar este tipo de formulaciones es necesario tener un modelo que describa el comportamiento de cada uno de los cuerpos que se desea unir. Uno de los modelos más utilizados es la discretización con elementos nitos. Gracias a la mecánica multicuerpo es posible analizar sistemáticamente el comportamiento de varios modelos conectados. Para ello, al sistema inicial de ecuaciones de cada cuerpo, denido en la expresión (17.1), se le añadirán las ecuaciones de ligadura. Dichas ecuaciones de ligadura afectarán a los puntos de contacto de cada uno de los N cuerpos imponiendo una restricción de movimiento y aplicando el principio de acción y reacción en dichos grados de libertad; es decir,   ∆fiA = λAB  ij  
A B AB L u ,u ,λ = ∆fjB = −λAB ij     A ui − uB j = 0 134 (19.1) CAPÍTULO 19. donde L uA , uB , λAB ij modelo A y 135 MECÁNICA MULTICUERPO λAB ij
son las ecuaciones de ligadura, uA son las incógnitas del son las incógnitas relativas a las ligaduras de los modelos en los grados de libertad i y j A y B respectivamente. En el método de los elementos nitos se representa la rigidez de un cuerpo exible como  A K11 A K11 A Knn ··· ···   . . . .. . .  .. . . .   A A A KA =   Ki1 · · · Kii · · · Kin  . . . ..  . . . . . .  .  A A A Kn1 · · · Kni · · · Knn n grados de   A  u1   .   ..      A  A u =  ui    .   .   .    A un y que está vinculada a cada uno de los             (19.2) libertad del cuerpo (19.3) . Por lo tanto, el sistema resultante será del tipo
M u¨A + C u˙ A + L uA , uB , λAB = f A − qB ij (19.4) Para enlazar dicha rigidez se añaden las ecuaciones de ligadura (19.1). Éstas se representan matricialmente como dos vectores relacionados con los modelos a conectar. La columna de la matriz asociada al modelo modelo B será:  LAB A que lo conecta con el   0   .   ..       = δ ki =   1   .   .   .    0 (19.5) CAPÍTULO 19. donde δ ki 136 MECÁNICA MULTICUERPO es la función delta de Kronecker e i el nodo de conexión. L es una matriz con tantas columnas como ligaduras a otros modelos y tantas las como grados de libertad tenga el modelo. Cada columna de L aporta la información sobre una conexión a un modelo colindante. Para ello, dicha columna valdrá k, que coincida con el grado de libertad conectado, 1 en aquella la, i. A i ¢fAB l ¢f BA ¢fAC j B ¢fCA k C Figura 19.1: Sistema multicuerpo Utilizando esta nomenclatura, las matrices del sistema (19.4) aplicado a los cuerpos A, B y C se ensamblan de la siguiente manera. La matriz tangente de las fuerzas internas del sistema vendrá dada por  KA      Kt =   
T   LAB 
T LAC LAB K B BA −L KC − LBA LAC −LCA
T − LCA
T             (19.6) CAPÍTULO 19. 137 MECÁNICA MULTICUERPO las incógnitas del problema serán       u=      A u u B uC λAB ij λAC lk             (19.7) y el vector de fuerzas externas estará denido como  A AB  f + ∆f + ∆f   f B + ∆f BA   f = f C + ∆f CA    0   0 AC             Gracias a esta formulación sistemática, los diferentes cuerpos (19.8) A,B ,..., N se pue- den estudiar conjuntamente utilizando distintos modelos para cada uno. En este caso, todos los cuerpos se han modelado utilizando elementos nitos, pero es posible aplicar esta sistematización a cuerpos modelados por diferentes técnicas, como se propondrá más adelante. 19.2. Aplicación a catenarias con modelos FEM La mecánica multicuerpo puede ser de utilidad en estructuras compuestas por partes que se repiten regularmente, como es el caso de las catenarias ferroviarias. Dichas catenarias están compuestas por la repetición de vanos idénticos a lo largo de cada uno de los cantones que suelen medir cientos de metros. Al utilizar esta técnica se reduce el tamaño inicial del modelo ya que sólo es preciso denir un cuerpo que se repetirá tantas veces como se desee. Con esto no se reduce, no obstante, el tamaño del problema, aunque se simplica el cálculo del equilibrio inicial. Para analizar y vericar los diferentes modelos, se estudiará la catenaria de- CAPÍTULO 19. 138 MECÁNICA MULTICUERPO nida por la norma de validación de modelos numéricos de interacción catenariapantógrafo EN50318 y que se muestra en la gura 19.2. Figura 19.2: Catenaria ferroviaria EN50318 de 3 vanos Cada vano de la catenaria constituye un cuerpo diferente enlazado por sus extremos mediante multiplicadores de Lagrange a los cuerpos colindantes como se observa en la gura 19.3. Cada uno de estos cuerpos se puede analizar independientemente siempre que se obtengan los multiplicadores de los modelos colindantes. 2 ¸AB µ ¸u AB2 ¸u ¸µ ¸v 1 ¸AB µ ¸u ¸v AB1 AB1 ¸u 2 ¸BC µ BA2 ¸u BC2 BA1 1 ¸BC µ BA1 ¸µ ¸u CB2 ¸µ CB2 BA2 1 ¸CB v 1 ¸BC u BA1 1 ¸CB u 1 ¸CB µ 1 ¸BC v Figura 19.3: Descomposición de la catenaria por vanos Por ejemplo, es posible usar un paso de tiempo para la integración temporal más grueso en la zona alejada del punto de contacto del pantógrafo y uno más no en la zona cuyo cálculo sea más crítico. Para ello habrá que extrapolar el valor de los multiplicadores (fuerzas externas) en el vano de interés en los pasos de tiempo en los que sólo se resuelva la zona de crítica. En general, el uso de técnicas multicuerpo ofrecerá versatilidad a la resolución de cualquier problema. En la gura 19.4 están representados los multiplicadores necesarios en uno de los nodos de unión. Dichos multiplicadores representan las fuerzas que aparecen en la unión de los cuerpos. Estos multiplicadores serán las incógnitas de las ecuaciones algebraicas que ligan los cuerpos. Gracias a ellos se tiene una gran versatilidad en el tratamiento interno CAPÍTULO 19. 139 MECÁNICA MULTICUERPO de cada cuerpo. ¸v ¸AB µ ¸u ¸v AB AB ¸u BA BA ¸µ BA Figura 19.4: Ampliación de la ligadura en el hilo de contacto entre los vanos A y B Capítulo 20 Modelo multicuerpo jerárquico para la reducción del sistema Además de su uso en estructuras compuestas, la mecánica multicuerpo se usa habitualmente para resolver problemas compuestos por cuerpos de comportamiento físico muy diferente. También permite abordar el análisis dinámico de mecanismos exibles utilizando técnicas de elementos nitos. Éste es el caso de la interacción catenaria-pantógrafo, donde la catenaria se comporta de manera no lineal y el pantógrafo se puede modelar correctamente de manera lineal. Cuando se aplica una fuerza puntual sobre una catenaria, como la aplicada por el pantógrafo, esta se desplaza de forma no lineal, no sólo debido al aumento de la no linealidad geométrica, sino también al posible pandeo de las péndolas y pérdidas de contacto. Dichos desplazamientos son mucho mayores cerca de la zona de contacto y, si la longitud de la catenaria es sucientemente grande, en la zona alejada del punto de contacto los desplazamientos serán pequeños. Pese al comportamiento dinámico no lineal de la catenaria, en los lugares donde los desplazamientos sean pequeños este comportamiento se podrá considerar lineal. Gracias a ello, en estas zonas es posible aplicar técnicas de reducción de variables mediante superposición modal, como se ha detallado en el capítulo 18. Hay que tener en cuenta que el pantógrafo, al desplazarse, modica la zona que 140 CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA141 Analisis Modal ¸ ¸ ¸ Analisis Modal + + + + ¸ ¸ ¸ Zona FEM Figura 20.1: Descomposición de la catenaria por vanos Analisis Modal ¸ ¸ ¸ ¸ ¸ + + ¸ Zona FEM Figura 20.2: Paso de vano modal a vano FEM Analisis Modal ¸ + Analisis Modal ¸ ¸ + + + ¸ ¸ ¸ Zona FEM Figura 20.3: Paso de vano FEM a vano Modal se comporta de manera lineal. Por ello, para analizar la interacción dinámica que se produce entre el pantógrafo y la catenaria, no siempre se podrá utilizar el mismo tipo de modelo. Un cuerpo modelado con superposición modal deberá transformarse en un modelo completo a medida que el pantógrafo se acerque para poder capturar el comportamiento no lineal de la estructura. Asimismo, cuando el pantógrafo se aleja del cuerpo, se debe volver al modelo reducido para minimizar el número de variables del sistema, y de este modo se conseguirá reducir signicativamente el CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA142 tiempo de cálculo. 20.1. Formulación En la resolución combinada se aprovecharán las propiedades de la formulación multicuerpo. Como se ha explicado anteriormente, este tipo de formulaciones permite que, siempre que los cuerpos estén ligados correctamente, puedan aplicarse formulaciones diferentes para cada cuerpo con mayor libertad. Para resolver la interacción dinámica pantógrafo-catenaria se propone: utilizar superposición modal para analizar aquellos cuerpos que se encuentren más alejados del punto de contacto (comportamiento lineal); y un análisis completo por elementos nitos en aquéllos que estén más cerca (comportamiento no lineal). La dinámica de los cuerpos que se comportan linealmente quedará denida por el sistema de ecuaciones (18.6), mientras que los que se analizan utilizando un modelo completo lo estarán por el sistema (17.1). Para ligarlos, se utilizarán los mismos conceptos que se han aplicado anteriormente: unicidad de desplazamientos y la aparición de fuerzas de acción y reacción. No obstante, cuando uno de los cuerpos esté modelado utilizando superposición modal la ligadura deberá relacionar un desplazamiento con la amplitud de los modos de vibración del modelo colindante. Para ello se realizará un cambio de base en el nodo de contacto. Por lo tanto, si a las ligaduras (19.1) se les aplica el cambio de base (18) se obtiene que la ligadura entre un cuerpo modelado utilizando superposición modal y uno modelado por elementos nitos es donde  A  ∆f˜ = ϕi λAB  ij 
 A B AB B AB L ξ ,u ,λ = ∆fj = −λij     ϕ i ξ A − uB j = 0
A L ξ A , uB , λAB son las ecuaciones de ligadura, ξ ij los modos del modelo A, u B (20.1) son las amplitudes de son los desplazamientos de los grados de libertad del CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA143 modelo B y λAB ij son las incógnitas relativas a las ligaduras de los modelos en los grados de libertad i y j, A y B respectivamente. Se observa cómo en el modelo A, al que se le ha aplicado la superposición modal, la fuerza de contacto se reparte en todos los modos del modelo. Por ello, A ∆f˜ es un vector en lugar de ser un escalar. Esta ligadura, expresada de forma vectorial, es diferente si el modelo al que va asociado es modal o completo. Mientras que si es completo coincide con la expresión (19.5), si es modal la ligadura será LAB = ϕi = ΦAB (20.2) Para diferenciar las ligaduras de modelos con superposición modal y modelos completos, las ligaduras asociadas a modelos con superposición modal se designarán con la letra Φ. Por lo tanto, un sistema multicuerpo para análisis combinado se empleará como matriz tangente al vector de esfuerzos internos  ˜A K        
T   ΦAB 
T ΦAC ΦAB K B BA −L KC − LBA ΦAC −LCA
T − LCA
T             (20.3) , como incógnitas,             ξ A uB uC λAB ij λAC lk             (20.4) CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA144 y como matriz de fuerzas externas,   A AB AC ˜ ˜ ˜  f + ∆f + ∆f    B BA   f + ∆f     C CA   f + ∆f       0     0 (20.5) 20.2. Resultados y vericación El modelo multicuerpo FEM-Modal se ha implementado en el lenguaje de propósito general MATLAB R . Para realizar la vericación se ha utilizado un cantón de 10 vanos de la catenaria denida por la norma europea EN50318 y representada en la gura 19.2. Se han simulado 6 segundos de desplazamiento del pantógrafo bajo la catenaria a 300 km/h lo que supone una distancia de 500 metros a partir del punto de carga del pantógrafo. Se han modelado utilizando elementos nitos aquellos vanos que tuvieran algún punto a menos de 20 metros del pantógrafo en cada instante y por superposición modal el resto de vanos. Por lo tanto, dependiendo de la posición del pantógrafo en cada instante, cada vano ha sido modelado por superposición modal o usando la malla completa de elementos nitos. Para la resolución numérica se ha utilizado un integrador α-generalizado con radio espectral, ρ = 0,9. En la gura 20.2 se presentan superpuestas las fuerzas de contacto resultantes utilizando el modelo completo de elementos nitos y el modelo multicuerpo FEMModal. Dichas fuerzas han sido ltradas a 20Hz tal y como dene la norma EN50318. Los resultados muestran que el modelo responde a la perfección en el cálculo de la fuerza de contacto y además reduce el tiempo de ejecución en un 60 % utilizando estos parámetros. Mientras el modelo completo necesitó 3h 23m 32s para completar la simulación, el modelo multicuerpo tan solo requirió 1h 6m 55s. La reducción en el tiempo de ejecución depende criticamente de las variables utilizadas como la distancia del pantógrafo en la que el cálculo completo, el número de modos utilizados CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA145 250 Completo Modal+FEM Fuerza de contacto [N] 200 150 100 50 0 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.4: Fuerza de contacto: FEM vs. multicuerpo FEM-Modal y el tiempo total de simulación. Un análisis de sensibilidad a estas variables se llevará a cabo en la sección siguiente. La gura 20.5 representa el desplazamiento vertical de un poste mientras el pantógrafo avanza. Como ya se ha comentado anteriormente, el máximo de este valor es crítico para cumplir los requisitos denidos por la norma europea EN50318. Se observa cómo en la zona de interés, dónde los desplazamientos son máximos, los resultados usando el modelo multicuerpo FEM-Modal son muy precisos. Conforme el pantógrafo se aleja del poste, situado a 300 metros, dicha precisión disminuye. No obstante, pese a que los resultados son puntualmente diferentes, el rango de desplazamientos es similar, y por ello, al alejarse del cuerpo reducido por superposición modal el efecto sobre la fuerza de contacto y el desplazamiento del pantógrafo es mínimo. A la vista de los resultados se puede concluir que el método garantiza la resolu- CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA146 Completo Modal+FEM 0.06 0.05 Desplazamiento [m] 0.04 0.03 0.02 0.01 0 −0.01 −0.02 180 195 210 225 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 375 Figura 20.5: Desplazamientos: FEM vs. multicuerpo (FEM+Modal) ción de este problema ya que ofrece una solución precisa y rápida. Dependiendo de la precisión necesaria en los cálculos se podrán, además, ajustar los parámetros consiguiendo simulaciones más o menos rápidas. Para determinar el efecto que tienen los distintos parámetros sobre la solución, se ha realizado un análisis de sensibilidad tal y como se expone en la sección 20.3. 20.3. Análisis de sensibilidad El tiempo y la precisión del cálculo multicuerpo FEM-modal depende de dos variables fundamentalmente: qué zona se considera lineal y cuántos modos de vibración se utilizan para modelar dicha zona lineal. Cuanto más alejada del punto de contacto esté la zona considerada lineal, mayor será el número de cuerpos modelados utilizando la malla completa de elementos CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA147 nitos y, por tanto, el número de incógnitas a resolver en cada iteración será mucho mayor. Esto afectará negativamente al tiempo de cálculo, ya que el proceso será computacionalmente más costoso; pero, por otro lado, la zona no lineal será más amplia y, consecuentemente, la precisión de los cálculos será mejor. Por otro lado, aumentar el número de modos mejora la respuesta en la zona lineal, ya que ésta se aproxima más a la real y, por lo tanto, será mayor la precisión en el cálculo de la fuerza de contacto y de los desplazamientos. Es de esperar que el número de modos de vibración también afecte al tiempo de ejecución, ya que a cada modo está asociada una incógnita en cada vano modelado por superposición modal. No obstante, teniendo en cuenta el tamaño del problema, el efecto del número de modos en el tiempo de computación debe ser menor que el del tamaño de la zona considerada no lineal y modelada por elementos nitos. En las guras 20.6, 20.7, 20.8 y 20.9 se presentan fuerzas de contacto calculadas utilizando diferentes número de modos de vibración y en las guras 20.11, 20.12, 20.13 y 20.14 se presentan desplazamientos. En dichas grácas se representa la fuerza de contacto entre el pantógrafo y la catenaria y el desplazamiento del quinto poste durante el paso del pantógrafo por el quinto y el sexto vano a 300km/h. Estos vanos están comprendidos entre el cuarto poste, situado a 240 metros del inicio del cantón, y el sexto poste, situado a 360 metros. Además, los resultados se presentan ltrados a 20Hz como dene la norma europea EN50318. Cada una de las grácas conserva la zona modelada por elementos nitos y utiliza 10, 20 y 30 modos para cada cuerpo modelado utilizando superposición modal. Estos resultados son contrastados con los obtenidos utilizando el modelo completo de elementos nitos. En las guras 20.10 y 20.15 se observa el efecto de variar el tamaño de la zona modelada usando elementos nitos. Ambas grácas se han realizado utilizando 30 modos en los cuerpos reducidos. Los resultados obtenidos utilizando 10 modos de vibración y 15 metros de zona CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA148 de análisis por elementos nitos son aceptablemente buenos, aunque la forma de la fuerza de contacto varía signicativamente respecto al modelo completo en los primeros y últimos metros de cada vano. Sin embargo, al aumentar el número de modos o el tamaño de la zona no lineal se observa una clara convergencia hacia los resultados obtenidos utilizando el modelo completo. 250 Completo 10 modos 20 modos 30 modos Fuerza de contacto [N] 200 150 100 50 0 240 260 280 300 320 Posición del pantógrafo [m] 340 360 Figura 20.6: Fuerza de contacto con 15 metros de análisis FEM En la gura 20.10 se pone de maniesto el efecto que tiene sobre la fuerza de contacto la distancia modelada utilizando elementos nitos. Se puede observar cómo la diferencia es más acusada en los máximos y mínimos, mientras que en las zonas ascendentes y descendentes los resultados son exactamente los mismos en el modelo completo y los reducidos. Respecto a los desplazamientos se observa cómo el análisis modal afecta más a la zona que queda detrás del pantógrafo. La zona que se encuentra inmediatamente delante del pantógrafo responde de manera idéntica utilizando ambos métodos de cálculo. CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA149 250 Completo 10 modos 20 modos 30 modos Fuerza de contacto [N] 200 150 100 50 0 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.7: Fuerza de contacto con 20 metros de análisis FEM 250 Completo 10 modos 20 modos 30 modos Fuerza de contacto [N] 200 150 100 50 0 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.8: Fuerza de contacto con 30 metros de análisis FEM CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA150 De la gura 20.15 se deduce que el tamaño de la zona considerada no lineal y, por lo tanto, modelada por elementos nitos es especialmente importante para recoger con detalle el desplazamiento de los puntos que van quedando detrás del pantógrafo a medida que el tren se desplaza. En las tablas 20.1 y 20.2 se han recogido las principales variables de la validación de la norma EN50318. Se puede comprobar que todos los valores menos uno están dentro del rango delimitado por la norma y solo hay pequeñas variaciones dentro de dichos márgenes. Además se presenta de forma numérica el tiempo necesario para completar cada uno de los cálculos. Esta comparativa de tiempos se presenta también de forma gráca en la gura 20.16. A la vista de los resultados se puede concluir que el método presenta una precisión más que aceptable con unos tiempos de ejecución mucho más reducidos. Gracias a ello, este modelo permite un mejor proceso de diseño y optimización de estructuras ferroviarias. Modos Dist. FEM max F min F media F std F 205.72 39.91 115.82 34.41 15 m. 192.96 48.21 115.13 33.77 20 m. 188.54 31.08 115.69 33.52 30 m. 206.17 41.01 115.74 34.79 50 m. 214.26 39.10 116.20 36.79 15 m. 197.72 40.65 115.26 34.10 20 m. 197.86 30.25 115.63 35.26 30 m. 205.30 42.59 115.92 33.83 50 m. 212.30 40.59 116.12 35.52 15 m. 201.67 41.41 115.43 34.46 20 m. 206.88 28.07 115.76 35.95 30 m. 207.17 41.61 115.98 34.21 50 m. 213.37 40.30 116.15 35.55 Completo 10 20 30 Tabla 20.1: Comparativa de resultados en fuerzas CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA151 Modos Dist. FEM desp 4 desp 5 desp 6 Tiempo 58.38 59.05 57.47 2h 22m 28.4s 15 m. 56.63 59.15 57.78 0h 41m 20.7s 20 m. 57.32 58.48 55.28 0h 46m 18.5s 30 m. 57.36 58.67 56.74 0h 51m 36.4s 50 m. 57.84 58.52 56.42 1h 6m 22.5s 15 m. 57.03 58.23 56.73 0h 53m 28.5s 20 m. 57.18 58.31 55.40 0h 56m 16.0s 30 m. 57.37 58.05 56.31 1h 0m 51.0s 50 m. 57.69 58.34 56.21 1h 9m 48.9s 15 m. 57.11 57.90 54.37 0h 58m 29.4s 20 m. 57.19 58.66 55.12 1h 1m 39.7s 30 m. 57.42 58.12 56.10 1h 18m 20.7s 50 m. 57.65 58.05 56.24 1h 33m 25.8s Completo 10 20 30 Tabla 20.2: Comparativa de resultados en desplazamientos 250 Completo 10 modos 20 modos 30 modos Fuerza de contacto [N] 200 150 100 50 0 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.9: Fuerza de contacto con 50 metros de análisis FEM CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA152 250 Completo 15 metros 20 metros 30 metros 50 metros Fuerza de contacto [N] 200 150 100 50 0 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.10: Fuerza de contacto con análisis modal de 30 modos de vibración Completo 10 modos 20 modos 30 modos 0.06 0.05 Desplazamiento [m] 0.04 0.03 0.02 0.01 0 −0.01 −0.02 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.11: Desplazamiento con 15 metros de análisis FEM CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA153 Completo 10 modos 20 modos 30 modos 0.06 0.05 Desplazamiento [m] 0.04 0.03 0.02 0.01 0 −0.01 −0.02 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.12: Desplazamiento con 20 metros de análisis FEM Completo 10 modos 20 modos 30 modos 0.06 0.05 Desplazamiento [m] 0.04 0.03 0.02 0.01 0 −0.01 −0.02 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.13: Desplazamiento con 30 metros de análisis FEM CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA154 Completo 10 modos 20 modos 30 modos 0.06 0.05 Desplazamiento [m] 0.04 0.03 0.02 0.01 0 −0.01 −0.02 240 255 270 285 300 315 Posición del pantógrafo [m] 330 345 360 Figura 20.14: Desplazamiento con 50 metros de análisis FEM Completo 15 metros 20 metros 30 metros 50 metros 0.06 0.05 0.04 0.03 0.02 0.01 0 −0.01 −0.02 240 255 270 285 300 315 330 345 360 Figura 20.15: Desplazamiento con análisis modal de 30 modos de vibración CAPÍTULO 20. MODELO MULTICUERPO JERÁRQUICO PARA LA REDUCCIÓN DEL SISTEMA155 1.4 10 modos 20 modos 30 modos 1.3 1.2 Tiempo [h] 1.1 1 0.9 0.8 0.7 15 20 25 30 35 Distancia de modelado FEM [m] 40 Figura 20.16: Análisis de tiempos 45 50 Capítulo 21 Conclusiones En esta cuarta parte se propone una estrategia apoyada en una formulación multicuerpo jerárquica para simulación de la interacción dinámica pantógrafo-catenaria. La jerarquización se modica de forma variable, para ello se establece un modelo de asignación dinámica de modelización de cuerpos independientes (en este caso vano a vano) contemplando dos posibilidades: elementos nitos en bases completas o en bases reducidas por descomposición modal. Las particularidades más relevantes que se han desarrollado en este modelo son las siguientes: Se ha implementado un algoritmo general para el acoplamiento de cuerpos rígidos o exibles asociado a un código de elementos nitos. Se ha desarrollado un modelo jerárquico con asignación dinámica que no ha sido publicado anteriormente en el campo de la dinámica estructural. Respecto al modelo presentado en la parte anterior cabe destacar las siguientes aportaciones: Reducción considerablemente el tiempo de computación. Los resultados obtenidos son similares a los obtenidos con el modelo completo. Si bien la precisión se reduce, la diferencia en los tiempos de cálculo lo compensa. 156 CAPÍTULO 21. CONCLUSIONES 157 La reducción de tiempo es sucientemente alta como para utilizar el método propuesto en aplicaciones de diseño u optimización con algorítmos metaheurísticos. El modelo propuesto ha sido comparado con los resultados publicados en la norma europea EN50318 [CEN99]. Los resultados obtenidos demuestran que la validación es altamente satisfactoria a pesar de no cumplir con todos los requisitos para la validación. Por lo tanto, la utilización complementaria de las herramientas presentadas en las partes I, II y III permite el diseño, optimización, cálculo y validación de cualquier catenaria ferroviaria existente en el mundo. Parte IV Conclusiones y Aportaciones Originales 158 159 De la misma forma en que se ha estructurado el presente proyecto nal de carrera, también las conclusiones se expondrán en tres bloques relativos a cada parte. En la primera parte se propone un nuevo método para el cálculo de la posición de equilibrio estático de estructuras tridimensionales de cables. Este método se basa en las ecuaciones analíticas de la catenaria y supone una generalización de la aplicación previa para el cálculo de equilibrio inicial de catenarias realizado por el equipo de investigación en mecánica computacional del ICAI coordinado por el director de este proyecto. La validación del método se ha llevado a cabo mediante la comparación con problemas de estructuras de cables publicados en la literatura cientíca. La precisión del método ha quedado altamente contrastada ya que las diferencias entre los resultados obtenidos y los publicados por otros autores son del de los casos y de menos del 1% 7% en el peor en valor medio. Las aportaciones más relevantes de este método respecto al anteriormente citado son las siguientes: El método anterior consideraba sólo cables estructuras de cables bidimensionales en el plano vertical. En el nuevo método, aprovechando que las catenarias están contenidas en un plano vertical, se incorpora una formulación semianalítica completamente tridimensional. Además se tiene en cuenta la elasticidad y la deformación térmica de los cables . A pesar de que esta suele ser despreciable algunas aplicaciones requieren la consideración de deformación elástica en estructuras de cables. Las restricciones propias de las estructuras de cables son incorporadas en la formulación de una forma natural. Además pueden emplearse para incluir modelos de dispositivos mecánicos como poleas, pretensores de cables, soportes, muelles, etc. Respecto a otros métodos de simulación de estructuras de cables las principales ventajas que éste presenta son las siguientes: 160 Debido a que este método está basado en las ecuaciones exactas de la catenaria, las ventajas propias del método propuesto en [LGCT06] le son inherentes. Desde un punto de vista numérico, el método presenta una alta eciencia ya que el tamaño del sistema de ecuaciones a resolver, en lugar de depender de la discretización espacial, depende de la topología de la estructura. Por otra parte la mayoría de algoritmos de resolución de sistemas de ecuaciones no lineales algebraicos requieren el cálculo de la matriz jacobiana. A pesar de que su cálculo es laborioso, las expresiones analíticas de la matriz jacobiana son intrínsecas al método lo cual le reporta una mayor rapidez y precisión. El problema de equilibrio inicial y el cálculo del equilibrio estático bajo cargas se resuelven siguiendo el mismo algoritmo. Gracias al tratamiento de los parámetros conocidos y desconocidos ambos tipos de problemas se pueden resolver fácilmente. El método tratado en esta parte es perfectamente válido para la mayoría de las aplicaciones ingenieriles tales como catenarias ferroviarias, sistemas de transporte de energía eléctrica, redes de metro, funiculares, etc. Además, a través de la implementación de dicho modelo en MATLAB es posible obtener los datos necesarios para resolver cálculos dinámicos por elementos nitos. Aprovechando los resultados obtenidos en la primera parte, en la segunda parte se propone un nuevo método para simulación de la interacción dinámica pantógrafo-catenaria. Las aportaciones más relevantes que se han llevado a cabo se pueden separar en tres aspectos: Se ha creado una nueva herramienta general de elementos nitos en la que ha sido posible implementar todos los avances presentados en este proyecto nal de carrera. Además, dicha herramienta es sucientemente versátil como para resolver otros problemas físicos de transmisión de calor o de electromagnetismo. Se ha implementado el efecto del pandeo en las péndolas mediante la denición 161 de una rigidez variable. Esto conere un tratamiento más realista del sistema de cables abordado. La simulación de la interacción dinámica catenaria-pantógrafo se ha desarrollado con un modelo de elementos nitos empleando una formulación corotacional para los hilos sustentador y de contacto. Para la integración numérica de la evolución dinámica del sistema se ha implementado el método α-generalizado que ha permitido incrementar la estabilidad de la solución ltrando las vibraciones de alta y muy alta frecuencia. Frente a todo lo publicado anteriormente al respecto en la literatura cientíca, este modelo sí ha sido validado mediante lo propuesto en la norma de validación europea EN50318 [CEN99]. Los resultados obtenidos se encuentran dentro de los margenes de aceptación establecidos en dicha norma, por lo que es apto para la certicación ocial de catenarias de alta velocidad. Por último, en la tercera parte de este proyecto se propone una estrategia apoyada en una formulación multicuerpo jerárquica para simulación de la interacción dinámica pantógrafo-catenaria. La jerarquía de modelos se modica según la respuesta de cada cuerpo (en este caso vano a vano) sea lineal o no lineal contemplando dos posibilidades: bases reducidas por descomposición modal o elementos nitos en bases completas. Las particularidades más relevantes que se han desarrollado en este modelo son las siguientes: Se ha implementado un algoritmo general para el acoplamiento de cuerpos rígidos o exibles asociado a un código de elementos nitos. Se ha desarrollado un modelo jerárquico con asignación dinámica que no ha sido publicado anteriormente en el campo de la dinámica estructural. Respecto al modelo presentado en la parte anterior cabe destacar las siguientes aportaciones: 162 Se consigue reducir considerablemente el tiempo de computación. Los resultados obtenidos son similares a los obtenidos con el modelo completo. Si bien la precisión se reduce, la diferencia en los tiempos de cálculo lo compensa. La reducción de tiempo es sucientemente alta como para utilizar el método propuesto en aplicaciones de diseño u optimización con algorítmos metaheurísticos. El modelo propuesto ha sido comparado con los resultados publicados en la norma europea EN50318 [CEN99]. Los resultados obtenidos demuestran que la validación es altamente satisfactoria a pesar de no cumplir con todos los requisitos para la validación. Por lo tanto, la utilización complementaria de las herramientas presentadas en las partes I, II y III permite el diseño, optimización, cálculo y validación de cualquier catenaria ferroviaria existente en el mundo. Parte V Bibliografía 163 Bibliografía [AAB74] J. H. Argyris, T. Angelopoulus, and B. Bichat. A general method for the shape nding of lightweight tension structures. in Applied Mechanics and Engineering, [AGR06] A. Andreu, L. Gil, and P. Roca. Computer Methods 3(1):135149, 1974. A new deformable catenary element Computers and Structures, for the analysis of cable net structures. 84(29):18821890, 2006. [Bar88] MR Barnes. Form-nding and analysis of prestressed nets and membranes. Computers and Structures, [BL99] D. Bruno and A. Leonardi. transport systems. 30(3):685695, 1988. Nonlinear structural models in cableway Simulation Practice and Theory, 7(3):207218, 1999. [BLM00] T. Belytschko, W.K. Liu, and B. Moran. continua and structures. J Wiley & Sons, 2000. [CCH84] F. Chaplin, G. Calderbank, and J. Howes. cable net structures. Nonlinear nite elements for The technology of suspended Longman Group Limited, New York, 1984. [Cel06] P. Cella. Methodology for Exact Solution of Catenary. tural Engineering, Journal of Struc- 125(12):14511453, 2006. [CEN99] CENELEC. Validation of simulation of the dynamic interaction between pantographs and overhead contact line. Railway Applications. Current Collection Systems prEN 50318, CENELEC, 1999. 164 165 BIBLIOGRAFÍA [CH93] J. Chung and G.M. Hulbert. A time integration algorithm for struc- tural dynamics with improved numerical dissipation: The generalized-α method. Journal of Applied Mechanics - ASME, [Cha02] D. Chamoret. 60:371375, 1993. Modélisation du contact: nouvelles approches numériques. PhD thesis, École Central de Lyon, 2002. [Chr96] P.M. Christou. An Integrated Technique for the Analysis of Frame and Cable Structures. PhD thesis, University of Florida, 1996. [CI99] C.S. Chen and D.E. Ingber. Osteoarthritis Cartilage, skeleton to cytoskeleton. [Cri90] M.A. Criseld. Tensegrity and mechanoregulation: from A consistent co-rotational formulation for non-linear, three-dimensional, beam-elements. chanics and Engineering, [Cri91] M.A. Criseld. 7(1):8194, 1999. Computer Methods in Applied Me- 81:131150, 1990. Non-linear Finite Element Analysis of Solids and Struc- tures, Vol.1: Essentials. J Wiley & Sons, 1991. [CV01] T.F. Coleman and A. Verma. A Preconditioned Conjugate Gradient Approach to Linear Equality Constrained Minimization. Optimization and Applications, [DF00] P. De Fonseca. Computational 20(1):6172, 2001. Simulation and Optimisation of the Dynamic Behaviour of Mechatronic Systems. PhD thesis, Katholieke Universiteit Leuven, 2000. [DJB94] J.G. De Jalón and E. Bayo. Kinematic and Dynamic Simulation of Multibody Systems - The Real-Time Challenge. [Dre32] C.S. Drewry. A memoir on suspension bridges. London, 1832. Springer Verlag, 1994. A. and R. Spottiswoode, 166 BIBLIOGRAFÍA [EOSW88] S.D. Eppinger, D.N. O'Connor, W.P. Seering, and D.N. Wormley. Modeling and experimental evaluation of asymmetric pantograph dynamics. Transactions of the ASME. Journal of Dynamic Systems, Measurement and Control, 110(2):p168  174, 1988. [ES98] P. Eberhard and W. Schiehlen. Hierarchical modeling in multibody dynamics. Archive of Applied Mechanics (Ingenieur Archiv), 68(3):237 246, 1998. [Fle87] R. Fletcher. Practical Methods of Optimization [M]. John Wiley & Sons, 1987. [GB88] L. Grundig and J. Bahndorf. The design of wide-span roof structures using microcomputers [J]. Computers and Structures, 30(3):495501, 1988. [GC01] M. Géradin and A. Cardona. Element Approach. Flexible Multibody Dynamics: A Finite Wiley, 2001. [GR97] M. Géradin and D. Rixen. Mechanical Vibrations: Theory and Applica- tion to Structural Dynamics. John Wiley & Sons, 1997. [HA82] RB Haber and JF Abel. Initial Equilibrium Solution Methods for Cable Reinforced Membranes. Engineering, Computer Methods in Applied Mechanics and 30(3):263306, 1982. [HHT77] H. Hilber, T. Hughes, and R. Taylor. Improved numerical dissipation for time integration algorithms in structural dynamics. gineering and Structural Dynamics, Earthquake En- 5:283292, 1977. [HL00] K.M. Hsiao and W.Y. Lin. A co-rotational nite element formulation for buckling and postbuckling analyses of spatial beams. in Applied Mechanics and Engineering, Computer Methods 188:567594, 2000. 167 BIBLIOGRAFÍA [HL06] Y. Huang and W. Lan. Static analysis of cable structure. hematics and Mechanics, [Hug87] T.J.R. Hughes. Applied Mat- 27(10):14251430, 2006. The Finite Element Method, Linear Static and Dynamics Finite Element Analysis. Prentice Hall, 1987. Solving Ordinary Dierential Equations II - [HW91] E. Hairer and G. Wanner. Sti and Dierential-Algebraic Problems. Springer-Verlag, 1991. [Ing93] D.E. Ingber. Cellular tensegrity: dening new rules of biological design that govern the cytoskeleton. [Irv81] H.M. Irvine. Journal of Cell Science, 104:613627, 1993. Cable Structures. [KS72] K. Kanzaki and T. Sakai. MIT Press, Cambridge, MA, 1981. Studies on the cable crane hung at three supports (I)(1) A way of calculation supposing static balanced state. Journal of the Japanese Forestry Society, [Lam91] J.D. Lambert. 54:103112, 1972. Numerical Methods for Ordinary Dierential Systems. The initial value problem. John Wiley & Sons, 1991. [LGCM07] O. Lopez-Garcia, A. Carnicero, and JL Maroño. Inuence of stiness and contact modelling on catenarypantograph system dynamics. nal of Sound and Vibration, Jour- 299(4-5):806821, 2007. [LGCT06] O. Lopez-Garcia, A. Carnicero, and V. Torres. Computation of the initial equilibrium of railway overheads based on the catenary equation. Engineering structures, 28(10):13871394, 2006. [LS71] K. Linkwitz and H.J. Schek. Einige Bemerkungen zur Berechnung von vorgespannten Seilnetzkonstruktionen. (Ingenieur Archiv), 40(3):145158, 1971. [LYP98] CY Lai, Z. You, and S. Pellegrino. Reectors. Archive of Applied Mechanics Shape of Deployable Membrane Journal of Aerospace Engineering, 11(3):7380, 1998. 168 BIBLIOGRAFÍA [MC02] Jesús Montesinos and Antonio Carmona. Tecnología de la catenaria. Mantenimientos de Infraestructura RENFE, Madrid, 2002. [Mol84] M. Mollaert. Formnding of mixed structures . In Space Structures, Third Conference on 1984. [MS83] J.J. Moré and DC Sorensen. Computing a Trust Region Step. Journal on Scientic and Statistical Computing, SIAM 4:553, 1983. [MT90] E. Moncrie and BHV Topping. Computer methods for the generation of membrane cutting patterns. Computers and Structures, 37(4):441450, 1990. [nN59] N. Newmark. A method of computation for structural dynamics. of Engineering Mechanics Division - ASCE, Journal 85:6794, 1959. [nNOR91] B. Nour-Omid and C.C. Rankin. Finite rotation analysis and consistent linearization using projectors. and Engineering, Computer Methods in Applied Mechanics 93:353384, 1991. [nNW99] J. Nocedal and S.J. Wright. operations research. Numerical Optimization Springer series in Springer, 1999. [OC97] José Luis de la Fuente O' Connor. Técnicas de Cálculo para Sistemas de Ecuaciones, Programación Lineal y Programación Entera. Reverté, Barcelona, 1997. [OT71] J.R. Ockendon and A.B. Tayler. The Dynamics of a Current Collection System for an Electric Locomotive. Proceedings of the Royal Society of London, Series A (Mathematical and Physical Sciences), 322(1551):p447  468, 1971. [OTS67] F. Otto, R. Trostel, and F.K. Schleyer. Tensile Structures; Design, Struc- ture, and Calculation of Buildings of Cables, Nets, and Membranes. Press, 1967. MIT 169 BIBLIOGRAFÍA [PE95] C. Pacoste and A. Eriksson. Element behaviour in post-critical plane frame analysis. Computer Methods in Applied Mechanics and Engineering, 125:319343, 1995. [PE97] C. Pacoste and A. Eriksson. Beam elements in instability problems. Computer Methods in Applied Mechanics and Engineering, 144:163197, 1997. + [PEM 97] G. Poetsch, J. Evans, R. Meisinger, W. Kortüm, W. Baldauf, A. Veitl, and J. Wallaschek. Pantograph-catenary dynamics and control. System Dynamics, 28:159195, 1997. [PG79] A.H. Peyrot and AM Goulois. Analysis of cable structures. and Structures, Vehicle Computers 10(5):805813, 1979. [PHJ03] T.J. Park, C.S. Han, and J.H. Jang. Dynamic sensitivity analysis for the pantograph of a high-speed rail vehicle. Journal of Sound and Vibration, 266(2):235260, September 2003. [Pug68] Alfred Grenville Pugsley. The theory of Suspension Bridges. Edward Arnold Ltd, London, 1968. [RB86] C.C. Rankin and F.A. Brogan. An element-independent co-rotational procedure for the treatment of large rotations. Technology - ASME, Journal of Pressure Vessel 108:165174, 1986. [RnO88] C.C. Rankin and B. Nour Omid. The use of projectors to improve nite element performance. [Ron78] Colin A. Ronan. Computers and Structures, 30:257267, 1988. The shorter Science and civilisation in China : an abridgement of Joseph Needham's original text. Press, Cambridge, 1978. Cambridge University 170 BIBLIOGRAFÍA [Sch74] HJ Schek. The force density method for form nding and computation of general networks. gineering, Computer Methods in Applied Mechanics and En- 3(1):115134, 1974. [SE64] A. Siev and J Eidelman. Stress analysis of prestressed suspended roofs. Journal of Structural Division, ASCE, 90:103121, 1964. [Sor97] DC Sorensen. Minimization of a large scale quadratic function subject to an ellipsoidal constraint. SIAM Journal on Optimization, 7:141161, 1997. [Tim83] S. Timoshenko. History of Strength of Materials: with a brief account of the history of theory of elasticity and theory of structures. Courier Dover Publications, 1983. [TWK59] S.P. Timoshenko and S. Woinowsky-Krieger. Theory of plates and shells. McGraw-Hill, second edition, 1959. [Vin83] T. Vinayagalingam. Computer evaluation of controlled pantographs for current collection from simple catenary overhead equipment at high speed. ASME Journal of Dynamic Systems, Measurement, and Control, 105:287294, 1983. [WB98a] TX Wu and MJ Brennan. Basic analytical study of pantograph-catenary system dynamics. Vehicle System Dynamics, [WB98b] T.X. Wu and M.J. Brennan. catenary system dynamics. 30(6):443456, 1998. Basic analytical study of pantograph- Vehicle System Dynamics, 30:443456, 1998. [WB99] T.X. Wu and M.J. Brennan. Dynamic stiness of a railway overhead wire system and its eect on pantograph-catenary system dynamics. of Sound and Vibration, [Wri00] Kenneth R. Wright. Journal 219(3):483502, 1999. Machu Picchu: A Civil Engineering Marvel. Publications, Reston, VA, 2000. ASCE