Campos en Elasticidad
| Trabajo realizado por estudiantes | |
|---|---|
| Título | Campos escalares y vectoriales en elasticidad. Grupo 16-A |
| Asignatura | Teoría de Campos |
| Curso | 2014-15 |
| Autores | Araceli Martín, Juan Carlos Durán, Francisco Javier Alcaraz, Álvaro Llera, Clara Callejo, Manuel Escudero |
| Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura | |
Para este análisis y representación de campos escalares en elasticidad nos situaremos en el contexto de una placa plana que ocupa la región comprendida entre las parábolas :
- [math]P1: 18y -81x^{2}-1=0 [/math]
- [math]P2: 2y +x^{2}-1=0 [/math]
Para representar la placa se utilizará un sistema de coordenadas curvilíneas tal que:
- [math]x=uv[/math]
- [math] y = \dfrac{1}{2}(u^2-v^2)[/math]
En nuestro análisis, el dominio de [math]u[/math] y [math]v[/math] comprenderá: [math](u,v) \in [1/3,1]*[-1,1][/math]
Contenido
1 º Situación inicial de la placa
Primero comenzaremos el estudio de la placa en cuestión planteando su representación gráfica mediante un mallado, utilizando, en el código, la conversión a coordenadas curvilíneas de [math]u[/math] y [math]v[/math]. El intervalo en el que representaremos comprende::
[math](x,y) \in [-1,1]*[-1,1][/math]
Para discretizar los vectores [math]u[/math] y [math]v[/math] utilizaremos un paso [math] h = \dfrac{1}{20}[/math].
1 h=1/20; % Paso de muestreo.
2 u=1/3:h:1; % Intervalo [1/3,1].
3 v=-1:h:1; % Intervalo [-1,1].
4 [uu,vv]=meshgrid(u,v); % Malla de datos.
5 xx=uu.*vv; % Parametrización X.
6 yy=(1/2).*((uu.^2)-(vv.^2)); % Parametrización Y.
7 plot(xx,yy); % Imagen.
8 mesh(xx,yy,0*xx) % Mallado.
9 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
10 view(2) % Ver imagen desde arriba.
Los resultados de la representación se ven en la imagen:
1.1 Líneas Coordenadas
Cuando se hace uso de cambio de coordenadas a ciertas coordenadas curvilíneas, es útil para el entendimiento de la transformación la representación de las Líneas Coordenadas. Las líneas coordenadas son aquellas que se obtienen variando una de las coordenadas de la transformación ([math]u[/math] o [math]v[/math]) y manteniendo fija la restante. En este estudio hemos representado varias líneas coordenadas a base de dar un valor concreto a [math]u[/math] o a [math]v[/math], dentro de sus respectivos intervalos:
1 xx11=uu.*0.5 ;xx12=uu.*-0.5;xx13=uu.*1;xx14=uu.*-1;
2 xx15=uu.*0.75;xx16=uu.*-0.75;xx17=uu.*0; % Parametrización X fijando v (0.5/-0.5/1/-1/0.75/-0.75/0) y variando u.
3 yy11=(1/2).*((uu.^2)-(0.5.^2));yy12=(1/2).*((uu.^2)-((-0.5).^2));
4 yy13=(1/2).*((uu.^2)-(1.^2));yy14=(1/2).*((uu.^2)-((-1).^2));
5 yy15=(1/2).*((uu.^2)-(0.75.^2));yy16=(1/2).*((uu.^2)-((-0.75).^2));
6 yy17=(1/2).*((uu.^2)-(0.^2)); % Parametrización Y fijando v (0.5/-0.5/1/-1/0.75/-0.75/0) y variando u.
7 xx21=vv.*0.5;xx22=vv.*0.4;xx23=vv.*1;xx24=vv.*0.9;xx25=vv.*0.75;
8 xx26=vv.*0.65;xx27=vv.*(1/3); % Parametrización X fijando u (0.5/0.4/1/0.9/0.75/0.65/0.333) y variando v.
9 yy21=(1/2).*((0.5.^2)-(vv.^2));yy22=(1/2).*((0.4.^2)-(vv.^2));
10 yy23=(1/2).*((1.^2)-(vv.^2));yy24=(1/2).*((0.9.^2)-(vv.^2));
11 yy25=(1/2).*((0.75.^2)-(vv.^2));yy26=(1/2).*((0.65.^2)-(vv.^2));
12 yy27=(1/2).*(((1/3).^2)-(vv.^2)); % Parametrización X fijando u (0.5/0.4/1/0.9/0.75/0.65/0.333) y variando v.
13 subplot(1,2,1); % Muestra varias imágenes. 1ª Imagen.
14 hold on % Inicio superposición de gráficos.
15 mesh(xx11,yy11,0*xx);mesh(xx12,yy12,0*xx);mesh(xx13,yy13,0*xx);
16 mesh(xx14,yy14,0*xx);mesh(xx15,yy15,0*xx);mesh(xx16,yy16,0*xx);
17 mesh(xx17,yy17,0*xx); % Mallados.
18 axis([-1,1,-1,1]) % Selecciona la regíon a dibujar.
19 view(2) % Ver imagen desde arriba.
20 hold off % Fin superposición de gráficos.
21 subplot(1,2,2); % Muestra varias imágenes. 2º Imagen.
22 hold on % Inicio superposición de gráficos.
23 mesh(xx21,yy21,0*xx);mesh(xx22,yy22,0*xx);mesh(xx23,yy23,0*xx);
24 mesh(xx24,yy24,0*xx);mesh(xx25,yy25,0*xx);mesh(xx26,yy26,0*xx);
25 mesh(xx27,yy27,0*xx); % Mallados.
26 axis([-1,1,-1,1]) % Selecciona la regíon a dibujar.
27 view(2) % Ver imagen desde arriba.
28 hold off % Fin superposición de gráficos.
A continuación se muestra las gráficas resultantes:
1.2 Base Natural
Cuando se realiza una transformación a coordenadas curvilíneas (de [math]x[/math] e [math]y[/math] a [math]u[/math] y [math]v[/math]), el vector de posición [math] \vec{r_o}[/math] se expresa como:
[math] \vec{r_o}=x\hat{e_1} +y \hat{e_2} =uv\hat{e_1} +\dfrac{1}{2}(u^2-v^2) \hat{e_2} [/math]
La base natural [math]\vec{g_u}, \vec{g_v}[/math] es la que tiene por vectores el resultado de derivar el vector posición [math] \vec{r_o}[/math] según las nuevas coordenadas [math]u[/math] y [math]v[/math]. Dado que es un problema sobre una placa plana, estamos en una situación de dos dimensiones en la que para cualquier base, sólo se requieren dos vectores, ([math]\vec{g_u}, \vec{g_v}[/math]). No obstante, para cuestiones que trataremos posteriormente será necesario considerar una tercera coordenada (por ejemplo, para el cálculo del rotacional), por ello, incluiremos en nuestra base natural el vector [math]\vec{g_w}[/math].
- [math] \vec{g_u}=v\hat{e_1} +u \hat{e_2}[/math]
- [math] \vec{g_v}=u\hat{e_1} -v \hat{e_2}[/math]
El código de Matlab utilizado para hallar los vectores [math]\vec{g_u}, \vec{g_v}[/math] es :
1 plot(xx,yy); % Imagen.
2 hold on % Inicio superposición de gráficos.
3 mesh(xx,yy,0*xx) % Mallado completo.
4 quiver(xx,yy,vv,uu); % Representación del primer vector de la base natural en cada punto.
5 quiver(xx,yy,uu,-vv); % Representación del segundo vector de la base natural en cada punto.
6 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
7 view(2) % Ver imagen desde arriba.
8 hold off % Fin superposición de gráficos.
La imagen obtenida es:
Resulta coherente porque los vectores dibujados son tangentes a las líneas coordenadas (Recordamos que [math]\vec{g_u}[/math] y [math]\vec{g_v}[/math] se obtienen al derivar el vector de posición [math] \vec{r_0}[/math] con respecto a [math]u[/math] y [math]v[/math]).
2 º Acción de la temperatura en la placa
2.1 Influencia de un foco de calor
La temperatura proviene de un foco de calor dada por el campo escalar [math]T(x,y)=(8-y^2+2y)e^{-x^2} [/math]
1 T=(8-yy.^2+yy.*2).*exp(-(xx.^2)); % Función Temperatura.
2 subplot(1,2,1); % Muestra varias imágenes. 1ª Imagen.
3 contour(xx,yy,T,20); % Define 20 líneas de nivel.
4 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
5 view(2) % Ver imagen desde arriba.
6 subplot(1,2,2); % Muestra varias imágenes. 2ª Imagen.
7 surf(xx,yy,T); colorbar; % Visualización de superficie en 3D más leyenda en color.
8 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
9 view(2) % Ver imagen desde arriba.
10 max(max(T)) % Valor máximo de la temperatura en toda la región.
Al ejecutar el código se obtienen estas gráficas:
El valor de temperatura máximo obtenido en la placa es 8.7332 y se da en el punto con el [math]|x|[/math] mínimo y la [math]y[/math] máxima.
2.2 Gradiente de T y curvas de nivel
El gradiente de una función escalar expresa la dirección en la cual el campo crece más rápido. Por otra parte, las curvas de nivel expresan los puntos que se encuentran a la misma altura, siendo en nuestro caso aquellos que tienen la misma temperatura. Por tanto, si tenemos dos curvas de nivel y estamos en un punto de la menor, el gradiente será el vector de mínima distancia a la otra curva de nivel, siendo en consecuencia perpendicular a ambas. [math]\nabla T(x,y)= \frac{\partial T}{\partial x} \hat{e_1}+\frac{\partial T}{\partial y}\hat{e_2}=[-2xe^{-x^2}(8-y^2+2y)]\hat{e_1}+[(-2y+2)e^{-x^2}]\hat{e_2}[/math]
1 h=1/20; % Paso de muestreo.
2 u=1/3:h:1; % Intervalo [1/3,1].
3 v=-1:h:1; % Intervalo [-1,1].
4 [uu,vv]=meshgrid(u,v); % Matrices.
5 xx=uu.*vv; % Parametrización X.
6 yy=(1/2).*((uu.^2)-(vv.^2)); % Parametrización Y.
7 f=(8-yy.^2+yy.*2).*exp(-(xx.^2)); % Función Temperatura.
8 fx =(-2.*xx).*(exp(-(xx.^2))).*(8-yy.^2+2.*yy); % Derivada con respecto a x de la función Temperatura.
9 fy=(exp(-(xx.^2))).*(-2.*yy+2); % Derivada con respecto a y de la función Temperatura.
10 hold on % Fin superposición de gráficos.
11 quiver(xx,yy,fx,fy) % Representación de los vectores gradiente.
12 contour(xx,yy,f,20);colorbar; % Visualización de superficie en 3D más leyenda en color.
13 view(2) % Ver imagen desde arriba.
14 hold off % Fin superposición de gráficos.
La imagen obtenida es la siguiente:
Efectivamente en la imagen se aprecia la perpendicularidad de los vectores del gradiente respecto a las curvas de nivel.
3 º Acción de una fuerza sobre el sólido
3.1 Campo de desplazamientos
Una fuerza determinada aplicada sobre nuestro sólido ha provocado un desplazamiento del mismo que viene dado por [math] \vec{u}(x,y) [/math] . Este vector será [math]\vec u(u,v)=\vec a(\vec b\cdot\vec r_{o})[/math] siendo :: [math] \vec{a}= \frac{\vec{g}_u}{|\vec{g_u}|}=\frac{v\vec {e_1}+ u\vec {e_2}}{\sqrt{v^2+u^2}}\qquad \vec{b}=-4 \frac{\vec{g}_u}{|\vec{g_u}|}=-4\frac{v\vec {e_1}+ u\vec {e_2}}{\sqrt{v^2+u^2}}[/math].
Como hemos hallado anteriormente [math]\vec{g_u}=v \vec{e_1} +u \vec{e_2} [/math]. Tomaremos:: [math]\vec{r_o}= x\vec {e_1}+y\vec {e_2}= uv\vec {e_1}+ \frac{1}{2}(u^2-v^2)\vec {e_2}[/math].
Con todo esto: [math]\vec{u}= \frac{\vec{g}_u}{|\vec{g_u}|}(-4 \frac{\vec{g}_u}{|\vec{g_u}|} \vec{r_o})= \frac{\vec{g_u}}{|\vec{g_u}|^2}(-4 \vec{g_u}\cdot\vec{r_o})=\frac{-4uv^2-2u(u^2-v^2)}{u^2+v^2} \vec{g_u}=\frac{-4uv^2-2u(u^2-v^2)}{u^2+v^2}(v \vec{e_1} +u \vec{e_2})=-2u\vec{g_u}=-2uv \hat{e_1} -2u^2 \hat{e_2}[/math] La representación del campo de desplazamiento [math]\vec{u}[/math] y la placa deformada será la siguiente:
1 Ux=-2.*uu.*vv; % Componente x del campo de desplazamientos u.
2 Uy=-2.*uu.^2; % Componente y del campo de desplazamientos u.
3 subplot(1,3,1) % Muestra varias imágenes. 1ª Imagen.
4 quiver(xx,yy,Ux,Uy); % Representación del campo vectorial de desplazamientos u.
5 axis([-1,1,-1,1]); % Selecciona la región a dibujar.
6 xd=xx+Ux; % Componente x final del sólido deformado.
7 yd=yy+Uy; % Componente y final del sólido deformado.
8 subplot(1,3,2), mesh(xx,yy,0*xx); % Muestra varias imágenes (2ª Imagen) y mallado completo.
9 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
10 view(2) % Ver imagen desde arriba.
11 subplot(1,3,3), mesh(xd,yd,0*xd); % Muestra varias imágenes (3ª Imagen) y mallado completo.
12 axis([-1,1,-2,1]); % Selecciona la región a dibujar.
13 view(2) % Ver imagen desde arriba.
Las imágenes obtenidas fueron las siguientes:
3.2 Divergencia de un campo
La divergencia de un campo vectorial [math]\vec{u}[/math] se halla por la expresión::
[math]\nabla\cdot\vec{u}= \frac{1}{ \sqrt{g} } \frac{\partial [\sqrt{g} u^{i}] }{\partial u^i}[/math]
La divergencia controla la diferencia entre el flujo saliente y el flujo entrante de un campo vectorial (sobre la superficie que rodea a un volumen de control), con lo cual, si la divergencia es positiva, el campo tiene "fuentes" , y si la divergencia es negativa se dice que tiene "sumideros" .
Por su parte, el campo [math]\vec{u}[/math] , en este caso, recordamos que estaba definido por:
[math]\vec{u}=-2uv \hat{e_1} -2u^2 \hat{e_2}[/math]
En el cálculo de la divergencia necesito las componentes contravariantes [math]u^u,u^v,u^w[/math] del campo. Éstas se calculan por la expresión [math]u^i=\vec{u}{g^i}[/math]. Recordando que [math]{g^i}[/math] son las componentes contravariantes de la base natural [math]\{\vec{g_u},\vec{g_v}\}=\{ v\hat{e_1} +u \hat{e_2} , u\hat{e_1} -v \hat{e_2}\}[/math], éstas se hallan a partir de [math]{g^i}=G^{ij}g_j[/math] siendo [math]G^{ij}[/math] la matriz inversa de la matriz de Gram de la base natural [math]G_{ij}[/math]:
[math]G_{ij}=\begin{bmatrix} \vec{g_u}\cdot\vec{g_u} & \vec{g_u}\cdot\vec{g_v} \\ \vec{g_v}\cdot\vec{g_u} & \vec{g_v}\cdot\vec{g_v}\end{bmatrix}=\begin{bmatrix} u^2+v^2 & 0 \\ 0 & u^2+v^2\end{bmatrix}\qquad G^{ij}=G^{-1}_{ij}[/math]
Operando la anterior expresión se calculan las [math]g^{i}[/math]:
- [math]g^{u}=\frac{1}{ u^2+v^2 }(v\hat{e_1} +u \hat{e_2}) [/math]
- [math]g^{v}=\frac{1}{ u^2+v^2 }(u\hat{e_1} -v \hat{e_2}) [/math]
De esta manera, las componentes contravariantes del campo quedan:
- [math]u^u=-4uv^2 -2u(u^2-v^2)[/math]
- [math]u^v=0[/math]
Y ya se pueden sustituir todos los términos en la expresión de la divergencia. El resultado final es [math]\nabla\cdot\vec u = -2\frac{3u^2+v^2}{u^2+v^2}[/math]. Dado que el resultado es negativo se puede concluir que existen sumideros en el flujo que atraviesa el campo vectorial [math]\vec{u}[/math]. Ésa es la expresión que se utilizará en el código de Matlab.
1 Divergencia=(-2.*(3*uu.^2+vv.^2)./(uu.^2+vv.^2)); % Campo Divergencia.
2 surf(xx,yy,Divergencia);colorbar; % Visualización de superficie en 3D más leyenda en color.
3 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
Podemos corroborar con la imagen de la placa deformada, que el mayor cambio de área se produce en la zona superior de la misma, la cual se deforma en sentido descendente; mientras que la menor variación se produce en los picos inferiores, los cuales se trasladan al lado opuesto, "doblándose" la placa sobre el eje de ordenadas.
3.3 Cálculo del rotacional de un campo vectorial
La expresión del rotacional de un campo vectorial [math]\vec{u}[/math] se halla por la siguiente expresión: [math]\nabla \times\vec{u}= \frac{1}{ \sqrt{g} } \begin{bmatrix} \vec{g_u} & \vec{g_v} & \vec{g_w} \\ \frac{\partial}{\partial u} & \frac{\partial}{\partial v} & \frac{\partial}{\partial w} \\ u_{u} &u_{v}&u_{w}\end{bmatrix}[/math]
De nuevo necesitaremos definir una tercera componente [math]\vec{g_w}[/math] para el cálculo de ese determinante :
- [math] \vec{g_u}=v\hat{e_1} +u \hat{e_2}[/math]
- [math] \vec{g_v}=u\hat{e_1} -v \hat{e_2}[/math]
- [math] \vec{g_w}=\hat{e_3}[/math] (suponiendo una tercera componente en la transformación [math] z=w[/math]).
El término [math]g[/math] es el determinante de la matriz de Gram [math]G[/math]de la base natural:
[math]G=\begin{bmatrix} \vec{g_u}\cdot\vec{g_u} & \vec{g_u}\cdot\vec{g_v} & \vec{g_u}\cdot\vec{g_w} \\ \vec{g_v}\cdot\vec{g_u} & \vec{g_v}\cdot\vec{g_v} & \vec{g_v}\cdot\vec{g_w} \\ \vec{g_w}\cdot\vec{g_u} &\vec{g_w}\cdot\vec{g_v}&\vec{g_w}\cdot\vec{g_w}\end{bmatrix}=\begin{bmatrix} u^2+v^2 & 0 & 0 \\ 0 & u^2+v^2 & 0 \\ 0 &0&1\end{bmatrix}[/math]
De modo que [math]g=(u^2+v^2)^2[/math]
Por su parte, el campo [math]\vec{u}[/math] había quedado definido por la expresión:
[math] \vec {u}(u,v)= \vec{a} (\vec{b}\cdot{r_{o}})=\frac{-4uv^2 -2u(u^2-v^2)}{ u^2+v^2 }(v\hat{e_1} +u \hat{e_2})[/math]
Volviendo al Rotacional, necesitamos, por último hallar las componentes covariantes [math]u_{u},u_{v},u_{w}[/math] del campo. Éstas se calculan por la expresión:
[math]u_{i}=\vec{u}{g_{i}}[/math]
De esta manera:
- [math]u_{u}=-4uv^2 -2u(u^2-v^2)[/math]
- [math]u_{v}=0[/math]
- [math]u_{w}=0[/math]
Y ahora se tienen todos los términos para sustituir en la expresión del rotacional:
[math]\nabla \times\vec{u}= \frac{1}{ \sqrt{g} } \begin{bmatrix} v\hat{e_1} +u \hat{e_2} & u\hat{e_1} -v \hat{e_2} & \hat{e_3} \\ \frac{\partial}{\partial u} & \frac{\partial}{\partial v} & \frac{\partial}{\partial w} \\ -4uv^2 -2u(u^2-v^2) &0&0\end{bmatrix}=\vec{0}[/math]
El significado físico de esta situación se traduce en que el campo [math]\vec{u}[/math] no tiene tendencia a rotación.
4 º Tensiones sobre la placa
4.1 Tensor de tensiones
En un medio elástico, lineal, isótropo y homogéneo, los desplazamientos permiten escribir el tensor de tensiones [math]\sigma ^ij=λ\nabla·\vec{u}1+2μ\epsilon[/math] siendo:
- [math]\epsilon (\vec{u})[/math] la parte simétrica del tensor gradiente de [math]\vec{u}[/math], [math]\nabla\vec{u}[/math]::
[math]\epsilon (\vec{u})=\frac{\nabla\vec{u}+ \nabla\vec{u}^t}{2}[/math]
- [math]λ[/math] y [math]μ[/math] los coeficientes de Lamé que dependen de las propiedades elásticas de cada material.
Se utilizarán estas expresiones para dibujar las tensiones normales en la dirección que marca [math]\vec{g_u}[/math] y la dirección que marca [math]\vec{g_v}[/math], es decir:
[math]\frac{\vec{g}_u}{|\vec{g_u}|}·σ·\frac{\vec{g}_u}{|\vec{g_u}|}\qquad\frac{\vec{g}_v}{|\vec{g_v}|}·σ·\frac{\vec{g}_v}{|\vec{g_v}|}[/math]
Bajo estas instrucciones se empieza a definir [math]\sigma ^ij[/math]. Para ello nos serviremos de la fórmula de las derivadas parciales covariantes [math]u^{i}, _j = \frac{\partial u^i}{\partial x^j}+\Gamma^i_{jk}u^k[/math]
1. Tomaremos [math]λ=μ=1[/math].
2. Se calcula el tensor gradiente de [math]\vec{u}[/math], [math]\nabla\vec{u}[/math]:
[math]\nabla\vec{u}= \begin{bmatrix} \frac{\partial u^u}{\partial u} + \Gamma ^1_{11} u^u & \frac{\partial u^u}{\partial v} + \Gamma ^1_{12} u^u & \frac{\partial u^u}{\partial w} + \Gamma ^1_{13} u^u \\ \frac{\partial u^v }{\partial u} + \Gamma ^{2}_{11}u^u & \frac{\partial u^v }{\partial v} + \Gamma ^{2}_{12}u^u & \frac{\partial u^v }{\partial w} + \Gamma ^{2}_{13}u^u \\ \frac{\partial u^w }{\partial u} + \Gamma ^{3}_{11}u^u &\frac{\partial u^w }{\partial v} + \Gamma ^{3}_{12}u^u&\frac{\partial u^w }{\partial w} + \Gamma ^{3}_{13}u^u\end{bmatrix}=\begin{bmatrix}-2- \frac{2u^2}{u^2+v^2} & \frac{-2uv}{u^2+v^2} & 0 \\\frac{2uv}{u^2+v^2} & \frac{-2u^2}{u^2+v^2} & 0\\0&0&0 \end{bmatrix}\qquad\nabla \vec{u}^t=\begin{bmatrix}-2- \frac{2u^2}{u^2+v^2} & \frac{2uv}{u^2+v^2} & 0 \\\frac{-2uv}{u^2+v^2} & \frac{-2u^2}{u^2+v^2} & 0\\0&0&0 \end{bmatrix}[/math]:
3. Parte simétrica del tensor gradiente [math]\nabla\vec{u}[/math], [math]\epsilon (\vec{u})[/math] :
Por lo tanto [math]\epsilon (\vec{u})=\begin{bmatrix}\frac{-4u^2-2v^2}{(u^2+v^2)} & 0 & 0 \\\ 0 & \frac{-2u^2}{(u^2+v^2)} & 0\\0&0&0 \end{bmatrix} [/math]
4. Cálculo de [math]σ[/math]: Se recuerda que la divergencia de [math]\vec{u}[/math] es [math]\nabla\cdot\vec u = -2\frac{3u^2+v^2}{u^2+v^2}[/math]
[math]\sigma ^ij=λ\nabla·\vec{u}1+2μЄ=\begin{bmatrix}\frac{-14u^2-6v^2}{u^2+v^2} & 0 & 0 \\\ 0 & \frac{-10u^2-2v^2}{u^2+v^2} & 0\\0 &0& \frac{-6u^2-2v^2}{u^2+v^2} \end{bmatrix}[/math]
Por último, para la representación de las tensiones normales en la dirección de [math]\vec{g_u}[/math] y la dirección de [math]\vec{g_v}[/math], se necesita definir:[math]\frac{\vec{g}_u}{|\vec{g_u}|}= \frac{v\vec {e_1}+ u\vec {e_2}}{\sqrt{v^2+u^2}}\qquad\frac{\vec{g}_v}{|\vec{g_v}|}=\frac{u\vec {e_1}-v\vec {e_2}}{\sqrt{v^2+u^2}}[/math]
Sabiendo esto, se puede proceder en Matlab al cálculo de las matrices que permitirán posteriormente la representación de las tensiones normales en la dirección que marca [math]\vec{g_u}[/math] y la dirección de [math]\vec{g_v}[/math]:
- [math]\frac{\vec{g}_u}{|\vec{g_u}|}·σ·\frac{\vec{g}_u}{|\vec{g_u}|}[/math]
- [math]\frac{\vec{g}_v}{|\vec{g_v}|}·σ·\frac{\vec{g}_v}{|\vec{g_v}|}[/math]
4.2 Tensión de Von Mises
La tensión de Von Mises es un escalar que indica el comportamiento plástico del sólido.: [math] \sigma _{vm}= \sqrt{ \frac{( \sigma _1- \sigma _2)^2+(\sigma_2-\sigma_3)^2+(\sigma_3-\sigma_1)^2 }{2} }[/math] [math] \sigma_1,\sigma_2 [/math] y [math] \sigma_3 [/math] son los autovalores del tensor [math]σ[/math] (tensiones principales).
1 sigma1=(-14*uu.^2-6*vv.^2)./(vv.^2+uu.^2); % Autovalor 1.
2 sigma2=(-10*uu.^2-2*vv.^2)./(vv.^2+uu.^2); % Autovalor 2.
3 sigma3=(-6*uu.^2-2*vv.^2)./(vv.^2+uu.^2); % Autovalor 3.
4 Mises=sqrt(((sigma1-sigma2).^2+(sigma2-sigma3).^2+(sigma3-sigma1).^2)./2); % Fórmula de Von Misses.
5 surf(xx,yy,Mises); % Visualización de superficie en 3D.
6 axis([-1,1,-1,1]) % Selecciona la región a dibujar.
7 max(max(Mises)) % Valor máximo de la tensión de Von Mises en toda la región.
Este fue el resultado:
El valor máximo de la tensión de Von Mises es 6.9282 y se alcanza en los puntos de la placa pertenecientes al eje de ordenadas.
5 º Masa de la placa
En último lugar procedemos al cálculo de la masa de la placa. Para ello partimos de una función de densidad [math]xye^{-1/x^2}[/math].
Para hallar la masa de la placa, utilizaremos un proceso numérico basado en el Método del Trapecio para la aproximación de integrales. El cálculo se basa en aplicar la función en cada punto y obtener una matriz de valores de la densidad en cada punto evaluado de la malla. Una vez hallada la matriz, multiplicándola por un vector fila y columna y sumando todos los elementos de la matriz obtendremos la masa total de la placa.
Como [math]x[/math] e [math]y[/math] pueden ser valores negativos, puede resultar que la función densidad sea negativa en algunos puntos de la placa. Para evitar este absurdo, los resultados finales se obtendrán a partir de convertir todos los valores de la matriz de densidades a su valor absoluto, multiplicándolos posteriormente por los pasos y sumándolos todos entre sí.
1 N1=200; N2=200; % Número de puntos.
2 a=1/3; b=1; c=-1; d=1; % Extremos de los intervalos.
3 h1=(b-a)/N1; h2=(d-c)/N2; % Pasos.
4 u=a:h1:b; v=c:h2:d; % Intervalos.
5 [uu,vv]=meshgrid(u,v); % Malla de datos.
6 xx=uu.*vv; % Parametrización X.
7 yy=(1/2).*((uu.^2)-(vv.^2)); % Parametrización Y.
8 d=(xx.*yy).*(exp(-1./(xx.^2))); % Función Densidad.
9 D=abs(d); % Función Densidad en valor absoluto.
10 w1=ones(N1+1,1); % Matriz de N1+1 elementos unidad.
11 w1(1)=1/2; w1(N1+1)=1/2; % Primer y último de la matriz w1.
12 w2=ones(N2+1,1); % Matriz de N2+1 elementos unidad.
13 w2(1)=1/2; w2(N1+1)=1/2; % Primer y último de la matriz w2.
14 result=h1*h2*w2'*D*w1 % Masa obtenida con la función densidad D.
El valor final de la masa obtenido es 0.0026
