Visualización de campos en un anillo plano (Grupo 9-C)
Contenido
1 Introducción
Dada una placa plana con forma de anillo circular centrado en el origen, y comprendido entre los radios 1 y 2 (Figura 1), se estudiarán las variaciones de temperatura en este objeto y las deformaciones causadas por el campo:
- [math] \vec u(\rho,\theta)=\frac{\sin(\pi \theta/2)}{30\rho}\vec g_{\rho} [/math]
En esta placa circular consideramos dos magnitudes físicas: la temperatura [math]T(\rho,\theta,t)[/math], que depende de las variables espaciales [math]ρ[/math] y [math]θ[/math] y del tiempo [math]t[/math], y los desplazamientos [math] \vec u (ρ,θ,t)[/math], que dependen de los mismos parámetros que la temperatura. De esta forma, si definimos [math] \vec r_0(\rho,\theta)[/math] como el vector de posición de los puntos de la placa en reposo, la posición de cada punto [math](\rho,\theta)[/math] de la placa para un instante de tiempo [math]t[/math] viene dada por:
- [math] \vec r (\rho,\theta,t)= \vec r_{0}(\rho,\theta)+\vec u(\rho,\theta,t) [/math]
2 Variaciones de la temperatura
Supongamos que tenemos un foco de calor en el centro de la placa circular, esta temperatura varía dependiendo solo del radio [math]ρ[/math], rigiéndose según la función [math]T(\rho,\theta)=-\log(\rho+0.1)[/math], como muestra la figura 2. Al ser una función logarítmica disminuye cuando aumenta el radio, por lo que la temperatura es más alta cuanto más nos acercamos al centro, como se observa en la figura 3, a la izquierda, en la que el centro queda el color rojo (más caliente, mayor temperatura) variando gradualmente hasta el azul (más frío, menor temperatura). De esta manera se define para cada radio una superficie equipotencial, es decir la variación de temperatura a lo largo de esa circunferencia es nula. En la misma figura, a la derecha, se puede observar como varía la temperatura en función del radio, estando ésta representada en el eje Z de la gráfica.
Dadas estas superficies equipotenciales se puede definir el gradiente de temperatura ([math]\nabla T[/math]), es decir, la dirección en la que el campo de temperaturas varía más rápidamente, esto conlleva que el [math]\operatorname{grad}(T)[/math] sea perpendicular a las superficies equipotenciales (variación nula), tal como se muestra en la figura 4.
3 Definición del campo de deformaciones [math] \vec u [/math]
Consideramos el campo de deformaciones:
[math] \vec u(\rho,\theta)=\frac{\sin(\pi \theta/2)}{30\rho}\vec g_{\rho} [/math]
El campo provoca deformaciones en el plano OXY, como se puede ver en la figura 5, que representa el campo de deformaciones en 3 dimensiones a la derecha. A la izquierda aparece representado el campo en el plano OXY. Los vectores representan la dirección, módulo y sentido del campo en cada punto. Podemos ver que el valor del campo no es uniforme, sino que varía según el punto del anillo en el que nos encontremos. Hay tres zonas de valor máximo y tres de valor mínimo, como muestra el módulo de los vectores. Tomaremos el valor del campo vectorial [math] \vec u [/math] en cada punto como el desplazamiento en dicho punto. Así tenemos que el vector posición será:
- [math] \vec r (\rho,\theta,t)= \vec r_{0}(\rho,\theta)+\vec u(\rho,\theta,t) [/math]
En la figura 6 hemos representado la posición inicial y final del anillo según el desplazamiento descrito. Apenas se aprecia la deformación, pero esto es debido a que esta tiene un valor muy pequeño, al encontrarse dividido entre [math] 30\rho [/math] .
Para ejemplificar esto hemos variado el denominador y hemos representado las posiciones inicial y final para el campo:
- [math] \vec u(\rho,\theta)=\frac{\sin(\pi \theta/2)}{\rho}\vec g_{\rho} [/math]
Si aplicamos el valor de este campo como desplazamiento en cada punto obtenemos una deformación clara del cuerpo, como se aprecia en la figura 7.
4 Divergencia y rotacional
La divergencia es una medida del cambio de volumen local debido al desplazamiento.
En este caso, al considerar el campo [math] \vec u(\rho,\theta)=\frac{\sin(\pi \theta/2)}{30\rho}\vec g_{\rho} [/math] , la divergencia la obtenemos de la expresión en coordenadas cilíndricas:
- [math] \nabla \cdot \vec u = \frac{1}{\sqrt{|g|}} \frac{\partial}{\partial x^i} \left(\sqrt{|g|} u^i \right)[/math]
Dado que [math] {u^\theta}[/math] y [math]{u^z}[/math] son cero, no se considerarán en el cálculo de la divergencia.
- [math] \nabla \cdot \vec u = \frac{1}{\rho}[\frac{\partial}{\partial \rho}\left(\rho\frac{\sin(\pi \theta/2)}{30\rho}\right)]=0[/math]
A raíz de este resultado, diremos que el campo vectorial [math] \vec u[/math] es adivergente [math] ( \nabla \cdot \vec u = 0)[/math] , esto es, el anillo no actúa como fuente [math] ( \nabla \cdot \vec u \gt 0)[/math] o sumidero [math] ( \nabla \cdot \vec u \lt 0)[/math] en ningún punto de la superficie. Se trata de un sólido incompresible.
En el anillo de estudio, donde [math]\vec u[/math] simboliza los desplazamientos, el rotacional representa el efecto de giro del sólido alrededor del vector [math] \vec n[/math] [math]\bot[/math] al anillo debido al desplazamiento:
- [math] \vec u=\vec u(\rho,\theta)=\frac{\sin(\pi \theta/2)}{30\rho}\vec g^{\rho} [/math]
- [math] \nabla\times\vec u = \begin{vmatrix} \vec g_{\rho} & \vec g_{\theta} & \vec g_{z} \\ \frac{\partial}{\partial \rho} & \frac{\partial}{\partial \theta} & \frac{\partial}{\partial z} \\ \frac{\sin(\pi \theta/2)}{30\rho} & 0 & 0 \end{vmatrix} =-\frac{\pi}{2}\frac{\cos(\pi \theta/2)}{30\rho}\vec g_{z} [/math]
- [math]|\nabla\times\vec u|=\frac{\pi}{2}\frac{\cos(\pi \theta/2)}{30\rho}[/math]
La representación del rotacional tiene forma de función cosenoidal. El valor del rotacional varía además en función del valor de [math]\rho[/math], haciéndose más pequeño a medida que [math]\rho[/math] crece (son inversamente proporcionales). Todo esto aparece representado en la parte derecha de la figura 8. En la parte izquierda de la figura tenemos representada la misma gráfica vista desde arriba (dos dimensiones). El color refleja el valor del rotacional, siendo el rojo los valores más altos y el azul los más bajos.
5 Deformaciones elásticas
Un material elástico lineal es aquel cuya relación tensión-deformación viene representada por una función lineal (recta). Un material será homogéneo si tiene las mismas propiedades en todos los puntos que lo componen y será isótropo si estas propiedades se extienden de igual forma en todas las direcciones (las propiedades no dependen de la posición). Considerando un medio elástico lineal, isótropo y homogéneo, se puede definir el tensor de tensiones [math]\sigma_{ij}[/math] mediante la siguiente fórmula: [math] \sigma_{ij}=\lambda \nabla \cdot \vec u \delta_{ij} + 2\mu \epsilon_{ij} [/math]
Siendo:
- [math] \vec u [/math] el campo vectorial con el que estamos trabajando
- [math] \epsilon_{ij} [/math] la parte simétrica del tensor gradiente de dicho [math] \vec u [/math]
- [math] \nabla \cdot \vec u [/math] la divergencia de dicho [math] \vec u [/math]
- [math] \lambda , \mu [/math] los coeficientes de Lamé
5.1 Coeficientes de Lamé
El coficiente [math] \mu [/math] es el módulo de elasticidad transversal. El coeficiente [math]\lambda[/math] no tiene interpretación física directa, pero junto al parámetro [math] \mu [/math] forma una parametrización del módulo de elasticidad para medios isótropos homogéneos. En otras palabras, los parámetros [math]\lambda[/math] y [math] \mu [/math] son los que determinan el módulo elástico de un material (la relación tensión-deformación propia de cada material). Existen otros coeficientes que también afectan al módulo de elasticidad como el módulo de Young o el módulo de compresibilidad. Existe una tabla de conversión entre las diferentes constantes que afectan al módulo de elasticidad.
5.2 Tensiones normales
Al aplicar el tensor de tensiones a nuestro campo vectorial tenemos que la divergencia obtenida es igual a 0 (ver ap. Divergencia y rotacional) y por tanto queda que [math] \sigma_{ij}=2\mu \epsilon_{ij}, [/math] Si tomamos [math]\lambda=\mu=1[/math] tenemos que [math]\sigma_{ij}=2·1·(\nabla \vec u + \nabla \vec u^t)/2=(\nabla \vec u + \nabla \vec u^t)[/math]y podemos ver que únicamente depende del parámetro μ, módulo de elasticidad transversal, por lo que las tensiones que analizaremos ahora serán únicamente las normales.
Matriz de tensiones del tensor de tensiones [math]\sigma_{ij}[/math] → [math](\sigma^i_{ j})= \begin{bmatrix} \frac{-\sin(\pi \theta/2)}{30\rho^2} & \frac{\pi·\cos(\pi \theta/2)}{60\rho} & 0 \\ \frac{\pi·\cos(\pi \theta/2)}{60\rho} & \frac{2·\sin(\pi \theta/2)}{30\rho^2} & 0 \\ 0 & 0 & 0 \end{bmatrix} [/math]
Expresión del tensor de deformaciones en forma diádica
[math]\sigma^i_{j}= \frac{-\sin(\pi \theta/2)}{30\rho^2}\vec g_{\rho}\otimes\vec g^ρ+\frac{\pi·\cos(\pi \theta/2)}{60\rho}\vec g_{\rho}\otimes\vec g^θ+\frac{\pi·\cos(\pi \theta/2)}{60\rho}\vec g_{\theta}\otimes\vec g^ρ+\frac{2·\sin(\pi \theta/2)}{30\rho^2}\vec g_{\theta}\otimes\vec g^θ[/math]
5.2.1 Tensiones normales según la dirección de [math] \vec g_ρ [/math]:
- [math] \vec g_{\rho} \cdot \sigma \cdot \vec g_\rho = -\frac{2·\sin(\pi \theta/2)}{30\rho^2} \vec g_ρ [/math]
5.2.2 Tensiones normales según la dirección de [math]\vec g_\theta/\rho[/math]:
- [math] \vec g_\theta/\rho \cdot \sigma \cdot \vec g_\theta/\rho = \frac{2·\sin(\pi \theta/2)}{30\rho^2} \vec g_θ [/math]
5.3 Tensiones tangenciales
En las figuras 12 y 14 que se definen a continuación vemos que las deformaciones debidas a las tensiones tangenciales son inapreciables. Podemos reducir en nuestro programa de Matlab el valor del denominador, que es inversamente proporcional a la deformación, y ver como se deforma la placa. La deformación es mayor en los puntos cuya tensión es mayor. Gráficamente, el módulo de la tensión viene determinado por la longitud de cada flecha de tensión. También apreciamos que la deformación es mayor en las regiones de la placa con mayor densidad de flechas de campo.
5.3.1 Respecto al plano ortogonal a [math] \vec g_ρ [/math]:
- [math]|\sigma \cdot \vec g_\rho-(\vec g_\rho \cdot \sigma \cdot \vec g_\rho) \vec g_\rho| = \frac{\pi·\cos(\pi \theta/2)}{60\rho} \vec g_θ [/math]
5.3.2 Respecto al plano ortogonal a [math]\vec g_\theta/\rho[/math]:
- [math]|\sigma \cdot \vec g_\theta/\rho-(\vec g_\theta/\rho \cdot \sigma \cdot \vec g_\theta/\rho) \vec g_\theta/\rho| = \frac{\pi·\cos(\pi \theta/2)}{60\rho^2} \vec g_ρ [/math]
Para ver mejor las deformaciones debidas a las tensiones tangenciales podemos hacer algo similar a lo que hicimos para visualizar mejor el campo de deformaciones debidas al campo [math] \vec u [/math] en el apartado 3.
Como ya se ha explicado, la deformación es inversamente proporcional al denominador de la expresión de las tensiones tangenciales, y un valor del denominador tan elevado como el que tenemos (tanto en la dirección de [math] \vec g_ρ [/math] como de [math] \vec g_θ [/math] tenemos que [math] ρ [/math] aparece multiplicado por 60) hace que las deformaciones sean inapreciables.
Si reducimos este valor podremos apreciar una mayor deformación, tal y como se muestra a continuación:
Nueva tensión tangencial respecto al plano ortogonal a [math] \vec g_ρ [/math]:
- [math] \frac{\pi·\cos(\pi \theta/2)}{3\rho} \vec g_θ [/math]
Nueva tensión tangencial respecto al plano ortogonal a [math]\vec g_\theta/\rho[/math]:
- [math] \frac{\pi·\cos(\pi \theta/2)}{3\rho^2} \vec g_ρ [/math]
6 ANEXO: Programas Matlab
Visualización de la malla:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
mesh(xx,yy,0*xx) % dibujar el mallado
axis([-3,3,-3,3])
view(2) % vista 2D
Función logarítmica:
fplot('-log(u+0.1)',[0 50])
title('Funcion de temperatura')
grid on
Temperaturas en 2D y en 3D:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
T=-log(uu+0.1);
subplot(1,2,1);
surf(xx,yy,T)
view(2) % dibujar el mallado
title('temperaturas 2D');
subplot(1,2,2);
plot3(xx,yy,T)
title('temperaturas 3D');
Gradiente de temperaturas:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
T=-log(uu+0.1);
gradx=-xx./(uu.*log(uu+0.1));
grady=-yy./(uu.*log(uu+0.1));
mesh(xx,yy,0.*xx)
view(2)
axis([-3,3,-3,3])
hold on
quiver(xx,yy,gradx,grady)
Campo de vectores dado y campo de deformaciones:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
ux=(xx.*sin(vv.*pi./2))./(30.*uu);
uy=(yy.*sin(vv.*pi./2))./(30.*uu);
subplot(1,2,1);
mesh(xx,yy,0.*xx)
view(2)
axis([-3,3,-3,3])
hold on
quiver(xx,yy,ux,uy)
title('campo con vectores')
subplot(1,2,2)
plot3(ux,uy,0.*xx)
title('campo de deformaciones')
Deformación de la placa debida al campo dado:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
X=xx+(xx.*sin(vv.*pi./2))./(30.*uu);
Y=yy+(yy.*sin(vv.*pi./2))./(30.*uu);
subplot(1,2,1);
mesh(xx,yy,0.*xx)
view(2);
axis([-3,3,-3,3]);
title('inicial')
subplot(1,2,2);
mesh(X,Y,0.*X)
view(2);
title('final')
axis([-3,3,-3,3]);
Deformación de la placa debida al campo dado [math] \rightarrow [/math] VARIACIÓN DE LOS DATOS PARA MEJORAR LA VISUALIZACIÓN DE LAS DEFORMACIONES:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
X=xx+(xx.*sin(vv.*pi./2))./(1.*uu);
Y=yy+(yy.*sin(vv.*pi./2))./(1.*uu);
subplot(1,2,1);
mesh(xx,yy,0.*xx)
view(2);
axis([-3,3,-3,3]);
title('inicial')
subplot(1,2,2);
mesh(X,Y,0.*X)
view(2);
title('final')
axis([-3,3,-3,3]);
Rotacional:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
zz=abs((pi.*cos((pi.*vv)./2))./(60.*uu));
subplot(1,2,1);
surf(xx,yy,zz)
view(2)
axis equal
title('rotacional 2D')
subplot(1,2,2);
surf(xx,yy,zz)
title('rotacional 3D')
Tensiones normales según la dirección de [math] \vec g_ρ [/math]:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=(-2.*xx.*sin((pi.*vv)./2))./(30.*uu.^3);
ty=(-2.*yy.*sin((pi.*vv)./2))./(30.*uu.^3);
mesh(xx,yy,0.*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
Tensiones normales según la dirección de [math]\vec g_\theta/\rho[/math]:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=(-2.*yy.*sin((pi.*vv)./2))./(30.*uu.^3);
ty=(2.*xx.*sin((pi.*vv)./2))./(30.*uu.^3);
mesh(xx,yy,0.*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
Tensiones tangenciales respecto al plano ortogonal a [math] \vec g_ρ [/math]:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=abs(-pi.*yy.*cos((pi.*vv)./2))./(60.*uu);
ty=abs(pi.*xx.*cos((pi.*vv)./2))./(60.*uu);
mesh(xx,yy,0.*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
Deformaciones debidas a las tensiones tangenciales respecto al plano ortogonal a [math] \vec g_ρ [/math]:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=abs(-pi.*yy.*cos((pi.*vv)./2))./(60.*uu);
ty=abs(pi.*xx.*cos((pi.*vv)./2))./(60.*uu);
X=xx+abs(-pi.*yy.*cos((pi.*vv)./2))./(60.*uu); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
Y=yy+abs(pi.*xx.*cos((pi.*vv)./2))./(60.*uu); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
subplot(1,2,1);
mesh(xx,yy,0*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
axis([-3,3,-3,3]);
title('Tensiones')
subplot(1,2,2);
mesh(X,Y,0.*X)
view(2);
axis([-3,3,-3,3]);
title('Deformaciones')
Tensiones tangenciales respecto al plano ortogonal a [math]\vec g_\theta/\rho[/math]:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=abs(pi.*xx.*cos((pi.*vv)./2))./(60.*(uu.^2));
ty=abs(pi.*yy.*cos((pi.*vv)./2))./(60.*(uu.^2));
mesh(xx,yy,0.*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
Deformaciones debidas a las tensiones tangenciales respecto al plano ortogonal a [math]\vec g_\theta/\rho[/math]:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=abs(pi.*xx.*cos((pi.*vv)./2))./(60.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
ty=abs(pi.*yy.*cos((pi.*vv)./2))./(60.*(uu.^2 )); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
X=xx+abs(pi.*xx.*cos((pi.*vv)./2))./(60.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
Y=yy+abs(pi.*yy.*cos((pi.*vv)./2))./(60.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
subplot(1,2,1);
mesh(xx,yy,0.*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
axis([-3,3,-3,3]);
title('Tensiones')
subplot(1,2,2);
mesh(X,Y,0.*X)
view(2);
axis([-3,3,-3,3]);
title('Deformaciones')
Deformaciones debidas a las tensiones tangenciales respecto al plano ortogonal a [math] \vec g_ρ \rightarrow [/math] VARIACIÓN DE LOS DATOS PARA MEJORAR LA VISUALIZACIÓN DE LAS DEFORMACIONES:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=abs(-pi.*yy.*cos((pi.*vv)./2))./(3.*uu);
ty=abs(pi.*xx.*cos((pi.*vv)./2))./(3.*uu);
X=xx+abs(-pi.*yy.*cos((pi.*vv)./2))./(3.*uu); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
Y=yy+abs(pi.*xx.*cos((pi.*vv)./2))./(3.*uu); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
subplot(1,2,1);
mesh(xx,yy,0*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
axis([-3,3,-3,3]);
title('Tensiones')
subplot(1,2,2);
mesh(X,Y,0.*X)
view(2);
axis([-3,3,-3,3]);
title('Deformaciones')
Deformaciones debidas a las tensiones tangenciales respecto al plano ortogonal a [math]\vec g_\theta/\rho \rightarrow [/math] VARIACIÓN DE LOS DATOS PARA MEJORAR LA VISUALIZACIÓN DE LAS DEFORMACIONES:
h=0.1;
u=1:h:2; % definir intervalo rho [1,2]
v=0:h:2*pi+h; % definir intervalo theta [0,2*pi]
[uu,vv]=meshgrid(u,v); % definir mallado
xx=uu.*cos(vv); % parametrizacion
yy=uu.*sin(vv);
tx=abs(pi.*xx.*cos((pi.*vv)./2))./(3.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
ty=abs(pi.*yy.*cos((pi.*vv)./2))./(3.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
X=xx+abs(pi.*xx.*cos((pi.*vv)./2))./(3.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
Y=yy+abs(pi.*yy.*cos((pi.*vv)./2))./(3.*(uu.^2)); %Para ver las deformaciones, cambiar el 60 por 3 aprox.
subplot(1,2,1);
mesh(xx,yy,0.*xx)
view(2)
hold on
quiver(xx,yy,tx,ty)
axis([-3,3,-3,3]);
title('Tensiones')
subplot(1,2,2);
mesh(X,Y,0.*X)
view(2);
axis([-3,3,-3,3]);
title('Deformaciones')
