Ver el código de «Carga crítica de una columna (Grupo 7B)»

{{ TrabajoED | Carga crítica de una columna. Grupo 7B | [[:Categoría:Ecuaciones Diferenciales|Ecuaciones Diferenciales]]|[[:Categoría:ED13/14|Curso 2013-14]] | María Ramírez
Ignacio Posada 

Antonio López-Mateos

Pablo Bueno }}

= Introducción =
En este proyecto numérico hemos estudiado la estabilidad y deflexión de una columna vertical, sometida a una carga constante P. 

La base de la teoría de las columnas es la formula de Euler, que fue publicada en 1757 por Leonardo Euler, un matemático suizo. Esta formula es válida solamente para columnas largas. También tenemos que suponer que esta columna es recta, homogénea, y de sección transversal constante en toda su longitud.
Se aplica la ley de Hooke y los esfuerzos son inferiores al limite de proporcionalidad del material, para evitar fracturas en la estructura de la columna.

La carga critica de pandeo es la carga axial máxima que una columna puede soportar cuando esta a punto de pandearse.

Un esquema representativo sería:

[[Archivo:Columna1.JPG|marco|centro|Esquema de la deflexión de la columna]]

= Planteamiento matemático (LISTO)=

Examinando los datos iniciales del problema, nos encontramos ante una columna vertical larga y esbelta de sección transversal uniforme, y de longitud total igual a L. Su eje de simetría, por tanto, ocupa el intervalo espacial '''<math>x  \epsilon  [0, L]</math>''', donde '''<math>x=0</math>''' corresponde al
extremo superior de la columna, donde se aplica la carga, y '''<math>x=L</math>''' al extremo inferior de la columna. Sea '''<math>y(x)</math>''' la '''función elástica''' o función que expresa la curvatura del eje de la columna al aplicarle una fuerza vertical de compresión o carga '''<math>P</math>''', en su extremo superior (ver figura Apartado 1).  Al comparar los momentos flectores en cualquier punto de la columna se obtiene:

\begin{matrix}
y''(x)=\frac{M(x)}{E I(x)}
\end{matrix}

Siendo:

·     '''<math>E</math>''': módulo de elasticidad de Young.

·      '''<math>I(x)</math>''': momento de inercia de una sección trasversal respecto a una recta vertical por el centro.

·      '''<math>M(x)</math>''': momento flector. 

En el caso de la columna, el momento depende de la deflexión, de manera que: ''<math>M(x)=-Py(x)</math>''.

Por tanto, sustituyendo la fórmula del momento flector en la ecuación inicial: 
\[\begin{array}{crl}
 \\
y''(x)=\frac{M(x)}{E I(x)} \Longrightarrow   y''(x)=\frac{-Py(x)}{E I(x)} \Longrightarrow  y''(x)+\frac{Py(x)}{E I(x)}=0 \\
 \\
\end{array}\]
Para poder resolver este problema de contorno faltarían las condiciones de contorno que serían:
<math>
\left\{\begin{matrix} y(0)=0 \\ y(L)=0\\ \end{matrix}\right.
</math> 
ya que en los extremo de la viga: ''<math>x=0</math>'' y ''<math>x=L</math>'' tiene deflexión cero. 


Por tanto, el problema de contorno planteado para resolver el ejercicio sería:

<math>
\left\{\begin{matrix} y''(x)+\frac{Py(x)}{E I(x)}=0 \\ y(0)=0 \\ y(L)=0\\ \end{matrix}\right.
</math>

= Estabilidad de la columna =
Vamos a estudiar ahora la estabilidad de la columna en función de la carga P aplicada. El problema que tenemos que resolver para determinar los valores de P que proporcionan estabilidad a la misma, o lo que es lo mismo, los puntos donde la columna no flecta,es el problema de contorno ya planteado en el Apartado 2:

<math>
\left\{\begin{matrix} y''(x)+\frac{Py(x)}{E I(x)}=0 \\ y(0)=0 \\ y(L)=0\\ \end{matrix}\right.
</math>


Con el problema de contorno ya planteado, nos disponemos ahora a resolverlo por los procedimientos conocidos.
(Nota: para facilitar los cálculos, sustituimos las variables que dependen de '''x''' por variables independientes, como una herramienta de cálculo, esto es '''<math>y(x)=y</math>''', '''<math>I(x)=I</math>''').

Por tanto, el problema quedaría

\[\begin{array}{crl}
 \\
E Iy''+Py=0 \\
 \\
\end{array}\]


Resolución:

\[\begin{array}{crl}
 \\
E Iy''+Py=0 \Longrightarrow   E Ir^2+P=0 \Longrightarrow  r^2=\frac{-P}{E I}=0 \Longrightarrow r=\; \pm\;\sqrt{\frac{P}{E I}} i \\
 \\
\end{array}\]
La solución general de la ecuación diferencial homogénea será
\[\begin{array}{crl}
 \\
y\;= y(x)\; =\; A\; cos\; (\sqrt{\frac{P}{E I}} x) + \; B\; sin\; (\sqrt{\frac{P}{E I}} x) = 0\\
 \\
\end{array}\]

Teniendo en cuenta las condiciones de contorno iniciales '''<math>y(0)=0</math>''' y '''<math>y(L)=0</math>''' y
que, por definición, se dice que una columna es estable si la única solución posible de la función elástica '''<math>y(x)</math>''' es la trivial, '''<math>y(x)=0</math>''', sustituimos dichos valores en la solución general obtenida.

Para '''<math>y(0)=0</math>'''
\[\begin{array}{crl}
 \\
y(0)\; =\; A\; cos\; (0) + \; B\; sin\; (0) = 0 \Longrightarrow y(0)\; =\; A\; 1 + \; B\; 0 = 0\\
 \\
\end{array}\]

Sabiendo que la solución de '''<math>y(0)</math>''' tiene que ser 0, se llega a la conclusión de que, para que se cumpla la ecuación, tiene que ocurrir que '''<math>A=0</math>''' necesariamente. 

Por otro lado, para '''<math>y(L)=0</math>'''
\[\begin{array}{crl}
 \\
y(L)\; =\; 0\; cos\; (\sqrt{\frac{PL^2}{E I}}) + \; B\; sin\; (\sqrt{\frac{PL^2}{E I}}) = 0 \Longrightarrow y(L)\; =\; 0\ + \; B\; sin\; (\sqrt{\frac{PL^2}{E I}}) = 0\\
 \\
\end{array}\]

Como se puede comprobar, '''<math>y(L)=0</math>''', implica que '''<math>\; B\; sin\; (\sqrt{\frac{PL^2}{E I}}) = 0</math>'''.
Si '''<math>B=0</math>''', se tiene '''<math>y=0</math>''', pero, si '''<math>B ≠0</math>''', entonces '''<math>sin\; (\sqrt{\frac{P}{E I}} x) = 0</math>'''. La última condición indica que, para que la igualdad planteada se cumpla, el argumento de la función seno ha de ser un múltiplo entero de π. 
\[\begin{array}{crl}
 \\
\sqrt{\frac{P}{E I}} L= nπ \Longrightarrow {\frac{P}{E I}}= {\frac{n^2π^2}{L^2}} \\
 \\
\end{array}\]

Por lo tanto, para todo real B distinto de cero
\[\begin{array}{crl}
 \\
y = \; B\; sin\; (\frac{nπx}{L}) \\
 \\
\end{array}\]
Por lo tanto, para todo real B distinto de cero,  es una solución del problema para cada n. Puesto que la ecuación diferencial es homogénea, y ya hemos supuesto que <math>sin\; (\frac{nπx}{L})= 0</math>, no necesitamos escribir B si así lo deseamos, por lo que la ecuación nos quedaría
\[\begin{array}{crl}
 \\
y = sin\; (\frac{nπx}{L}) \\
 \\
\end{array}\]


Volviendo a la ecuación anteriormente obtenida, buscamos hallar el valor de P donde la columna se desvía. Esto sólo ocurre cuando la fuerza de compresión toma algunos de los valores de dicha ecuación. Esas fuerzas se llaman cargas críticas, que corresponden a la ecuación
\[\begin{array}{crl}
 \\
\sqrt{\frac{P}{E I}} L= nπ \Longrightarrow {\frac{P}{E I}}= {\frac{n^2π^2}{L^2}} \Longrightarrow P_{cr}\; = {\frac{n^2π^2EI}{L^2}} \\
 \\
\end{array}\]
Suponemos, según los datos proporcionados, que el módulo de elasticidad de Young es '''<math>E=1</math>''', que la sección de la columna circular es de radio constante '''<math>R=1\ m</math>''', que su densidad es '''<math>\rho=1\ Kg/m^3</math>''', que la longitud total de la columna es '''<math>L=10\ m.</math>''' y que el momento de inercia de una sección circular trasversal respecto a una recta vertical por el centro es '''<math>I = {\frac{πr^4}{2}}</math>'''. Por tanto, la nueva ecuación de P_{cr} quedaría expresada así
\[\begin{array}{crl}
 \\
P_{cr}\; = {\frac{n^2π^3r^4}{2L^2}} \\
 \\
\end{array}\]

El número natural n representa el número de curvas de desviación que posee la columna, por lo que en nuestro caso '''<math>n=1</math>'''.
Según este último dato, sustituyéndolo en las ecuaciones anteriormente halladas,obtenemos que 
\[\begin{array}{crl}
 \\
P_{cr}\; = {\frac{π^3r^4}{2L^2}} \\
 \\
\end{array}\]
\[\begin{array}{crl}
 \\
y(x) = sin\; (\frac{nπx}{L})\\
 \\
\end{array}\]
donde '''<math>P_{cr}</math>''' corresponde a la mínima carga crítica de la columna y '''y(x)''' a la curva de deflexión correspondiente a dicha carga crítica mínima, conocida también como '''Carga de Euler'''.

= Dibujar la función que a cada radio R 2 [1; 5] le asocia el valor más pequeño de Pcr para el cual la columna deja de ser estable. El valor Pcr se conoce como carga crítica (obtenido por primera vez por L. Euler)(3) =

== Proceso de obtencion de la ecuacion ==

== Resolucion numerica ==

= Punto 5 de despejar en funcion del radio =
Ahora vamos a trabajar con una columna de masa constante '''<math>M=31 Kg</math>'''. La columna tendrá una seccion circular con radio variable. Este radio cambiará en funcion de un polinomio de grado 2 respecto al punto superior de la columna, donde se aplica la carga (como hemos mencionado antes, '''<math>x=0</math>''' corresponde con la parte superior de la columna).

Este polinomio de grado 2, es simétrica respecto al punto '''<math>x=L/2</math>''', lo que quiere decir que tiene su extremo absoluto en dicho punto. Según la formula del vertice para parabolas o derivando e igualando a cero, obtenemos: '''<math>-b= L·a</math>'''.

\[\begin{array}{crl}
R(x) = a x^2 + b x + C  \Longrightarrow  R(x) = a x^2 - a L x + C\\
\end{array}\]

Vamos a trabajar con que el radio de la columna en la base es '''<math>Ro = R_{0}</math>'''. Por lo tanto:

\[\begin{array}{crl}
R(0) = a0^2 - aL0 + C  \Longrightarrow R(0) = R_{0} = C \Longrightarrow R(x) = ax^2 - aLx + R_{0}
\end{array}\]

== Proceso de obtencion de la ecuacion ==

Para obtener la ecuacion del radio igualaremos la masa a '''<math>31Kg</math>''':

\[\begin{array}{crl}
M = p V = 1 \int_0^L \! \pi R(x)^2 \, \mathrm{d}x =31  \;\;\; \Longrightarrow \;\;\; \pi ( \frac{L^5}{30}a^2 - \frac{L^3 R_{0}}{3}a + R_{0}^2 L) = 31
\end{array}\]

Despejando obtenemos el valor de '''a''' en funcion del '''<math>R_{0}</math>''':

\[\begin{array}{crl}
a = \frac{1}{20}R_{0} \; \pm\; \frac{3}{200}\sqrt{\frac{-20 R_{0}^2}{9}+\frac{4·31}{3\pi}}
\end{array}\]

Para saber si tenemos que coger el '''<math>\; \pm\;</math>''' de la raiz, hacemos la siguiente comprobación.

Si <math>R_{0} = 1 \Longrightarrow a=0</math>, para que <math>R(x)=1</math> sea constante. Sustituyendo:

\[\begin{array}{crl}
a = \frac{1}{20}R_{0} \; \pm\; \frac{3}{200}\sqrt{\frac{-20 \frac{1}{20}R_{0}^2}{9}+\frac{4·31}{3\pi}} \;\;\; \Longrightarrow \;\;\; a=0 \;\;\; \Longrightarrow \;\;\; a = \frac{1}{20}R_{0} \; -\; \frac{3}{200}\sqrt{\frac{-20 R_{0}^2}{9}+\frac{4·31}{3\pi}} \;\;\; \Longrightarrow \;\;\; \frac{1}{20}R_{0} = \frac{3}{200}\sqrt{\frac{-20 R_{0}^2}{9}+\frac{4·31}{3\pi}}
\end{array}\]

Por lo tanto:

 \[\begin{array}{crl}
R(x)= \left(\frac{1}{20}R_{0} \; -\; \frac{3}{200}\sqrt{\frac{-20 R_{0}^2}{9}+\frac{4·31}{3\pi}}\right) x^2 - \left(\frac{1}{20}R_{0} \; -\; \frac{3}{200}\sqrt{\frac{-20 R_{0}^2}{9}+\frac{4·31}{3\pi}}\right)10x + R_{0}
\end{array}\]

Vamos a trabajar con un intervalo '''<math>\; R_{0}\in (0.5,1.5)\;</math>''', para ver el radio máximo y mínimo que alcanza para cada valor. Lo estudiaremos en los extremos del intervalo para ver el cambio mas significativamente:

Para '''<math>\; R_{0}=0.5</math>''':
\[\begin{array}{crl}
a = -0.028247  \;\;\; \Longrightarrow \;\;\; R(x)=-0.028247x^2 +0.28247x + 0.5 \;\;\; \Longrightarrow \;\;\; R(L/2)=R(5)=1.206175
\end{array}\]

Para '''<math>\; R_{0}=1.5</math>''':
\[\begin{array}{crl}
a = 0.032156  \;\;\; \Longrightarrow \;\;\; R(x)=0.032156x^2 +0.32156x + 1.5 \;\;\; \Longrightarrow \;\;\; R(L/2)=R(5)=0.793825
\end{array}\]


Lógicamente, para '''<math>\; R_{0}=0.5</math>''' la función es cóncava ya que, a parte de '''<math>\; a<0</math>''' si '''<math>\; R_{0}<1</math>''' habrá que compensar las pérdidas iniciales de masa. Y obviamente viceversa para caso contrario. Adjuntamos un dibujo comparativo para ver lo que sucede:

 [DIBUJO]

== Comparacion de valores de Pcr segun el radio ==

= Resolucion por Fourier =








[[Categoría:Ecuaciones Diferenciales]]
[[Categoría:ED13/14]]
[[Categoría:Trabajos 2013-14]]