Trabajo grupo 4

De MateWiki
Saltar a: navegación, buscar
Trabajo realizado por estudiantes
Título Visualización de Campos Escalares y Vectoriales en Elasticidad. Grupo 4-C
Asignatura Teoría de Campos
Curso 2019-20
Autores Ana Regaliza Rodríguez

Andrea Palomar Expósito
Bertha Alicia Rodríguez
Marcos Nieto Horna

Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura


Consideramos una placa plana que ocupa un cuarto de anillo circular centrado en el origen de coordenadas y comprendido entre los radios [1,3], en el primer cuadrante [math]x,y\geq 0[/math]. En este campo tendremos definidas dos cantidades físicas, la temperatura [math]T(x,y)[/math] y los desplazamientos [math]\overrightarrow { u } (x,y)[/math].

1 Introducción

Durante este trabajo expondremos los cambios que sufre una placa cuando está sometida a distintas situaciones como temperaturas, desplazamientos o rotaciones.

2 Mallado de la placa


Con lo definido anteriormente, comenzaremos dibujando con ayuda del programa Matlab el mallado que representa los puntos interiores del sólido; consideramos como paso de muestreo h=1/10. Como podemos observar en el mallado el aro no está completamente cerrado en la coordenada [math]\theta =\pi [/math]. Esto se debe al muestreo que tenemos, ya que en el intervalo final correspondiente, éste será mayor que [math]\theta [/math], por lo que se omite.

Código Matlab
Mallado de la Placa

3 Temperatura

Una de las funciones definidas al principio es la temperatura, cuya ecuación es [math]T(x,y)=log(y+2)[/math]. Como ya conocemos la representación de la placa, ahora podremos observar como cambian las curvas de nivel y la barra de colores a lo largo del sólido debido a la temperatura que la afecta.

Código
Temperatura en la Placa










Al tener las curvas de nivel en nuestro sólido y la función en la que varía la temperatura podemos obtener el gradiente de la misma [math]\nabla T[/math]. Aprovechamos para verificar que nuestras ecuaciones y los diagramas están correctos al utilizar una de las propiedades del gradiente. Como podemos observar el [math]\nabla T[/math] es ortogonal a las curvas de nivel del anillo.

Recordamos la fórmula general del cálculo de [math]\nabla T[/math] en coordenadas cartesianas [math]\nabla T=\frac { \partial T(x,y) }{ \partial x } \overrightarrow { i } \quad +\frac { \partial T(x,y) }{ \partial y } \overrightarrow { j } \quad +\frac { \partial T(x,y) }{ \partial z } \overrightarrow { k } [/math].

Particularizándolo a nuestra temperatura [math]T(x,y)=log(y+2)[/math], obtenemos:
[math]\nabla T=0\overrightarrow { i } \quad +\frac { 1 }{ y+2 } \overrightarrow { j } \quad +0\overrightarrow { k } \longrightarrow [/math] [math] \nabla T=\frac { 1 }{ y+2 } \overrightarrow { j } [/math]

Código
Gradiente de Temperatura en la Placa

4 Campo de Desplazamiento

Queremos considerar un campo de desplazamiento [math]\overrightarrow { u } (\rho ,\theta )=f\left( \rho \right) { \overrightarrow { g } }_{ \theta }[/math] con las siguientes características:

  • Los puntos situados en [math]\rho =1[/math] no sufren desplazamiento.
  • El [math]\nabla x\overrightarrow { u } =\frac { 3\rho -2 }{ 10 } { \overrightarrow { g } }_{ z }[/math].

Para poder calcular el campo de desplazamiento [math]\overrightarrow { u } [/math] debemos recordar el cálculo general del rotacional en las coordenadas cilíndricas: [math]\nabla x \overrightarrow { u } =\frac { 1 }{ \rho } \left| \begin{matrix} { \overrightarrow { g } }_{ \rho } & { \overrightarrow { g } }_{ \theta } & { \overrightarrow { g } }_{ z } \\ \partial \rho & \partial \theta & \partial z \\ { \overrightarrow { u } }_{ \rho } & { \overrightarrow { u } }_{ \theta } & { \overrightarrow { u } }_{ z } \end{matrix} \right| \quad [/math]

Particularizamos el cálculo del rotacional a nuestro campo. una vez calculado lo igualamos a [math]\frac { 3\rho -2 }{ 10 } { \overrightarrow { g } }_{ z }[/math].
Nota: Como sabemos que nuestro rotacional solo depende de [math]{ \overrightarrow { g } }_{ z }[/math], no tenemos en cuenta [math]{ \overrightarrow { g } }_{ \rho }[/math].

[math]\nabla x\overrightarrow { u } =\frac { 1 }{ \rho } \left| \begin{matrix} { \overrightarrow { g } }_{ \rho } & { \overrightarrow { g } }_{ \theta } & { \overrightarrow { g } }_{ z } \\ \partial \rho & \partial \theta & \partial z \\ 0 & { \rho }^{ 2 }\cdot f\left( \rho \right) & 0 \end{matrix} \right| =\frac { 1 }{ \rho } \frac { \partial }{ \partial \rho } \left( { \rho }^{ 2 }f\left( \rho \right) \right) { \overrightarrow { g } }_{ z }+\frac { 1 }{ \rho } \frac { \partial }{ \partial z } \left( { \rho }^{ 2 }f\left( \rho \right) \right) { \overrightarrow { g } }_{ \rho }[/math][math]\qquad \longrightarrow \qquad\frac { 1 }{ \rho } \frac { \partial }{ \partial \rho } \left( { \rho }^{ 2 }f\left( \rho \right) \right) { \overrightarrow { g } }_{ z }=\frac { 3\rho -2 }{ 10 } { \overrightarrow { g } }_{ z }[/math].

Calculamos nuestra función [math]f\left( \rho \right)[/math] .
[math]\int { \frac { \partial }{ \partial \rho } { \rho }^{ 2 }f\left( \rho \right) } =\int { \frac { 3{ \rho }^{ 2 }-2\rho }{ 10 } } \qquad \longrightarrow \qquad { \rho }^{ 2 }f\left( \rho \right) =\frac { 3{ \rho }^{ 3 } }{ 3\cdot 10 } -\frac { { 2\rho }^{ 2 } }{ 2\cdot 10 } =\frac { { \rho }^{ 3 }-{ \rho }^{ 2 } }{ 10 } +C[/math].

[math]{ \rho }^{ 2 }f\left( \rho \right) =\frac { { \rho }^{ 3 }-{ \rho }^{ 2 } }{ 10 }+C \qquad \longrightarrow \qquad f\left( \rho \right) =\frac { \rho -1 }{ 10 } +\frac { C }{ { \rho }^{ 2 } } [/math].
[math]\overrightarrow { u } =\left( \frac { \rho -1 }{ 10 } +\frac { C }{ { \rho }^{ 2 } } \right) { \overrightarrow { g } }_{ \theta }[/math].
Imponemos la primera condición la cual nos dice que los puntos situados en [math]\rho =1[/math] no tienen desplazamiento. Esto quiere decir que nuestro campo [math]\overrightarrow { u } [/math], en [math]\rho =1[/math], es cero.

[math]\rho =1\longrightarrow \overrightarrow { u } =0\qquad \frac { C }{ { \rho }^{ 2 } } { \overrightarrow { g } }_{ { \theta } }=0\longrightarrow C=0[/math].
[math]\overrightarrow { u } =\left( \frac { \rho -1 }{ 10 } \right) { \overrightarrow { g } }_{ \theta }[/math].

Una vez obtenido el campo de desplazamientos, a continuación podemos ver el sólido antes y después de aplicar [math]\overrightarrow { u } [/math].

Código
Campo de vectores en la malla



Código
Sólido antes y después del desplazamiento













5 Divergencia

Al tener el vector [math]\overrightarrow { u } [/math] podremos aprovechar para calcular el cambio de densidad que afecta a nuestra placa en cada punto, es decir su divergencia. Ésta la calcularemos a continuación:
Nuestro campo anteriormente calculado: [math]\overrightarrow { u } =\left( \frac { \rho -1 }{ 10 } \right) { \overrightarrow { g } }_{ \theta }[/math].
Recordamos la fórmula general:
[math]\nabla \cdot \overrightarrow { u } =\frac { 1 }{ \sqrt { g } } \left( \frac { \partial { u }_{ \rho } }{ \partial \rho } +\frac { \partial { u }_{ \theta } }{ \partial \theta } +\frac { \partial { u }_{ z } }{ \partial z } \right) [/math].

Ahora calculamos la divergencia particularizada a nuestro campo:

[math]\nabla \cdot \overrightarrow { u } =\frac { 1 }{ \rho } \left( \frac { \partial \left( \rho \cdot 0 \right) }{ \partial \rho } +\frac { \partial \left( \rho \cdot \frac { \rho -1 }{ 10 } \right) }{ \partial \theta } +\frac { \partial \left( \rho \cdot 0 \right) }{ \partial z } \right) =\frac { 1 }{ \rho } \cdot \left( 0 \right) [/math].

[math]\nabla \cdot \overrightarrow { u } =0[/math]
Al realizar el cálculo de la divergencia podemos observas que nuestro resultado a dado cero, por lo que no afecta a nuestra placa.

Código
Divergencia en la Placa














6 Rotacional

Con una de las características que nos des daban para calcular en vector de desplazamiento [math]\overrightarrow { u } [/math], nos dan el rotacional [math]\nabla x\overrightarrow { u } =\frac { 3\rho -2 }{ 10 } { \overrightarrow { g } }_{ z }[/math]. Recordamos que el rotacional que es el vector perpendicular al plano que gira sobre si mismo. Lo aplicaremos a nuestra placa para observar como lo afecta.
Tenemos:
[math]\nabla x\overrightarrow { u } =\frac { 3\rho -2 }{ 10 } { \overrightarrow { g } }_{ z }[/math] que ya nos lo daban

Y su módulo será: [math]\left| \nabla x\overrightarrow { u } \right| =\frac { 3\rho -2 }{ 10 } [/math]

Código
Rotacional en la Placa













7 Tensiones

Definimos [math]\epsilon (\overrightarrow { u } )=\frac { \nabla \overrightarrow { u } +{ \nabla }\overrightarrow { u } ^{ t } }{ 2 } [/math]. Recordamos que el vector [math]\overrightarrow { u } =\left( \frac { \rho -1 }{ 10 } \right) { \overrightarrow { g } }_{ \theta }[/math].
[math]\nabla \cdot \overrightarrow { F } =\frac { 1 }{ \sqrt { g } } \frac { \partial }{ { \partial x }^{ j } } (\sqrt { g } { F }^{ j })\qquad [/math][math]g={ \left| { \overrightarrow { g } }_{ 1 } \right| }^{ 2 }{ \left| { \overrightarrow { g } }_{ 2 } \right| }^{ 2 }{ \left| { \overrightarrow { g } }_{ 3 } \right| }^{ 2 }[/math].

Aplicado a nuestro campo tendríamos: [math]\nabla \overrightarrow { u } =\left( \begin{matrix} \frac { 1 }{ { \left| { \overrightarrow { g } }_{ \rho } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 1 } }{ \partial \rho } & \frac { 1 }{ { \left| { \overrightarrow { g } }_{ \theta } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 1 } }{ \partial \theta } & \frac { 1 }{ { \left| { \overrightarrow { g } }_{ z } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 1 } }{ \partial z } \\ \frac { 1 }{ { \left| { \overrightarrow { g } }_{ \rho } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 2 } }{ \partial \rho } & \frac { 1 }{ { \left| { \overrightarrow { g } }_{ \theta } \right| }^{ 2 } } \frac { { \partial \overrightarrow { u } }_{ 2 } }{ \partial \theta } & \frac { 1 }{ { \left| { \overrightarrow { g } }_{ z } \right| }^{ 2 } } \frac { { \partial }\overrightarrow { u } _{ 2 } }{ \partial z } \\ \frac { 1 }{ { \left| { \overrightarrow { g } }_{ \rho } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 3 } }{ \partial \rho } & \frac { 1 }{ { \left| { \overrightarrow { g } }_{ \theta } \right| }^{ 2 } } \frac { { \partial \overrightarrow { u } }_{ 3 } }{ \partial \theta } & \frac { 1 }{ { \left| { \overrightarrow { g } }_{ z } \right| }^{ 2 } } \frac { { \partial \overrightarrow { u } }_{ 3 } }{ \partial z } \end{matrix} \right) [/math].

[math]\frac { 1 }{ { \left| { \overrightarrow { g } }_{ \rho } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 2 } }{ \partial \rho }=\frac { \partial { \overrightarrow { u } }_{ 2 } }{ \partial \rho } =\frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }+\frac { \rho -1 }{ 10 } \frac { { \partial \overrightarrow { g } }_{ \theta } }{ \partial \rho } =\frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }+\frac { \rho -1 }{ 10 } \frac { 1 }{ \rho } { \overrightarrow { g } }_{ \theta }[/math].
Calculamos las tensiones normales en la direccion que marca el eje g(rho), g(tetha), g(z)
[math]\frac { \partial { \overrightarrow { g } }_{ \theta } }{ \partial \rho } ={ \Gamma }_{ \theta \rho }^{ \rho }{ \overrightarrow { g } }_{ \rho }+{ \Gamma }_{ \theta \rho }^{ \theta }{ \overrightarrow { g } }_{ \theta }+{ \Gamma }_{ \theta \rho }^{ z }{ \overrightarrow { g } }_{ z }=\frac { 1 }{ \rho } { \overrightarrow { g } }_{ \theta }[/math].

[math]\frac { 1 }{ { \left| { \overrightarrow { g } }_{ \theta } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 1 } }{ \partial \theta }=\frac { { \partial \overrightarrow { u } }_{ 1 } }{ \partial \theta } =\frac { \rho -1 }{ 10 } \frac { { \partial \overrightarrow { g } }_{ \theta } }{ \partial \theta } =\frac { \rho -1 }{ 10 } \left( -\rho \right) { \overrightarrow { g } }_{ \rho }=\frac { -{ \rho }^{ 2 }+\rho }{ 10 } { \overrightarrow { g } }_{ \rho }[/math]
[math]\frac { { \partial \overrightarrow { g } }_{ \theta } }{ \partial \theta } ={ \Gamma }_{ \theta \theta }^{ \rho }{ \overrightarrow { g } }_{ \rho }+{ \Gamma }_{ \theta \theta }^{ \theta }{ \overrightarrow { g } }_{ \theta }+{ \Gamma }_{ \theta \theta }^{ z }{ \overrightarrow { g } }_{ z }=-\rho { \overrightarrow { g } }_{ \rho }[/math].

Para poder resolver la divergencia nos apoyamos en los símbolos de Christoffel, donde han establecido una matrices para las coordenadas cilíndricas que son las siguientes:
[math]{ \Gamma }^{ \rho }=\left( \begin{matrix} 0 & 0 & 0 \\ 0 & -\rho & 0 \\ 0 & 0 & 0 \end{matrix} \right) \qquad { \Gamma }^{ \theta }=\left( \begin{matrix} 0 & \frac { 1 }{ \rho } & 0 \\ \frac { 1 }{ \rho } & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) \qquad { \Gamma }^{ z }=\left( \begin{matrix} 0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) [/math].

Donde podemos ahora sustituir en la gráfica y finalmente tenemos: [math]\nabla \overrightarrow { u } =\left( \begin{matrix} 0 & \frac { -{ \rho }^{ 2 }+\rho }{ 10{ \rho }^{ 2 } } & 0 \\ \frac { 1 }{ 10 } (2-\frac { 1 }{ \rho } ) & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) [/math].


Siguiendo definiendo [math]\epsilon ({ \overrightarrow { u } })[/math], teniendo [math]\nabla \overrightarrow { u } [/math]. Ahora encontramos tu traspuesta [math]\nabla \overrightarrow { u } ^{ t }[/math]:
[math]\nabla \overrightarrow { u } ^{ t }=\left( \begin{matrix} 0 & \frac { 1 }{ 10 } (2-\frac { 1 }{ \rho } ) & 0 \\ \frac { -\rho ^{ 2 }+\rho }{ 10\rho ^{ 2 } } & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) [/math].

[math]\epsilon \left( \overrightarrow { u } \right) =\frac { \nabla \overrightarrow { u } +{ \nabla \overrightarrow { u } }^{ t } }{ 2 } =\frac { 1 }{ 20 } \left( \begin{matrix} 0 & 1 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) [/math]

Tenemos la tensión definida por [math]\sigma =\lambda \nabla \cdot \overrightarrow { u } 1+2\mu \epsilon [/math], donde sabemos que [math]\nabla \cdot \overrightarrow { u } =0[/math] y [math]\mu =\lambda =1[/math]

Por lo que finalmente la tensión será [math]\sigma =\left( \begin{matrix} 0 & \frac { 1 }{ 10 } & 0 \\ \frac { 1 }{ 10 } & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) [/math].

[math]\sigma =\frac { 1 }{ 10 } { \overrightarrow { g } }_{ \rho }\otimes { \overrightarrow { g } }_{ \theta }+\frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }\otimes { \overrightarrow { g } }_{ \rho }[/math].
Calculamos las tensiones normales en la direccion que marca el eje [math]{ \overrightarrow { g } }_{ \rho }[/math],[math]{ \overrightarrow { g } }_{ \theta },{ \overrightarrow { g } }_{ z }[/math]
[math]{ \overrightarrow { g } }_{ \rho }\cdot \sigma \cdot { \overrightarrow { g } }_{ \rho }=0\quad\quad,\quad \quad { \overrightarrow { g } }_{ z }\cdot \sigma \cdot { \overrightarrow { g } }_{ z }=0\quad \quad y\quad \quad \frac { { \overrightarrow { g } }_{ \theta } }{ \rho } \cdot \sigma \cdot \frac { { \overrightarrow { g } }_{ \theta } }{ \rho } =0[/math].

A partir de la tensión obtenida anteriormente, calcularemos las tensiones tangenciales respecto al plano ortogonal a [math]{ \overrightarrow { g } }_{ \rho }[/math], es decir [math]\left| \sigma \cdot { \overrightarrow { g } }_{ \rho }-\left( { \overrightarrow { g } }_{ \rho }\cdot \sigma \cdot { \overrightarrow { g } }_{ \rho } \right) { \overrightarrow { g } }_{ \rho } \right| [/math].

Tenemos la tensión [math]\frac { 1 }{ { \left| { \overrightarrow { g } }_{ \rho } \right| }^{ 2 } } \frac { \partial { \overrightarrow { u } }_{ 2 } }{ \partial \rho }=\frac { \partial { \overrightarrow { u } }_{ 2 } }{ \partial \rho } =\frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }+\frac { \rho -1 }{ 10 } \frac { { \partial \overrightarrow { g } }_{ \theta } }{ \partial \rho } =\frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }+\frac { \rho -1 }{ 10 } \frac { 1 }{ \rho } { \overrightarrow { g } }_{ \theta }[/math].

Por lo calculado anteriormente sabemos que [math]\left( { \overrightarrow { g } }_{ \rho }\cdot \sigma \cdot { \overrightarrow { g } }_{ \rho } \right) { \overrightarrow { g } }_{ \rho }=0[/math].
Entonces nos quedaría [math]\sigma \cdot { \overrightarrow { g } }_{ \rho }=\frac { 1 }{ 10 } \left( { \overrightarrow { g } }_{ \rho }\otimes { \overrightarrow { g } }_{ \theta } \right) { \overrightarrow { g } }_{ \rho }\quad +\quad \frac { 1 }{ 10 } \left( { \overrightarrow { g } }_{ \theta }\otimes { \overrightarrow { g } }_{ \rho } \right) { \overrightarrow { g } }_{ \rho }[/math].

Desarollándolo: [math]\frac { 1 }{ 10 } \left( { \overrightarrow { g } }_{ \theta }\cdot { \overrightarrow { g } }_{ \rho } \right) { \overrightarrow { g } }_{ \rho }\quad +\frac { 1 }{ 10 } \left( { \overrightarrow { g } }_{ \rho }\cdot { \overrightarrow { g } }_{ \rho } \right) { \overrightarrow { g } }_{ \theta }\quad =\quad \frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }[/math].

[math]\sigma \cdot { \overrightarrow { g } }_{ \rho }\quad =\quad \frac { 1 }{ 10 } { \overrightarrow { g } }_{ \theta }[/math].
[math]\left| \sigma \cdot { \overrightarrow { g } }_{ \rho } \right| =\sqrt { \frac { 1 }{ { 10 }^{ 2 } } { \overrightarrow { g } }_{ \theta }\cdot { \overrightarrow { g } }_{ \theta } } =\frac { 1 }{ 10 } \rho [/math].
Nota: Recordamos que [math]{ \overrightarrow { g } }_{ \theta }\cdot { \overrightarrow { g } }_{ \theta }={ \rho }^{ 2 }[/math]


En la siguiente gráfica observamos la tensión de Von Mises, la cual viene dada por la siguiente fórmula: [math]{ \sigma }_{ VM }=\sqrt { \frac { { \left( { \sigma }_{ 1 }-{ \sigma }_{ 2 } \right) }^{ 2 }+{ \left( { \sigma }_{ 2 }-{ \sigma }_{ 3 } \right) }^{ 2 }+{ \left( { \sigma }_{ 3 }-{ \sigma }_{ 1 } \right) }^{ 2 } }{ 2 } } [/math], donde sabemos que [math]{ \sigma }_{ 1 } ,{ \sigma }_{ 2 }y { \sigma }_{ 3 }[/math] son las tensiones principales, conocidas como los autovalores de [math]\sigma[/math] .

Código Matlab


Código Matlab
Tensiones de Von Mises
Tensiones de Von Mises















8 Masa

Al tener la densidad de la placa, que es: [math]d(x,y,z)=1+{ e }^{ -\left| x \right| /{ \left( y+1 \right) }^{ 2 } }[/math]. Podemos calcular la masa total mediante la siguiente integral:
[math]\int _{ S }^{ }{ f\cdot dS\quad =\quad \iint _{ D }^{ }{ f\left( s\left( u,v \right) \right) \left| \frac { \partial \overrightarrow { r } }{ \partial u } \left( u,v \right) x\frac { \partial \overrightarrow { r } }{ \partial v } \left( u,v \right) \right| dudv } } \quad =\quad \iint _{ D }^{ }{ f\left( s\left( u,v \right) \right) \left| { \overrightarrow { r } }_{ u }x{ \overrightarrow { r } }_{ v } \right| dudv } [/math].

Con la parametrización de nuestra placa: [math]\begin{cases} \rho =u \\ \theta =v \\ z=0 \end{cases}\qquad \begin{matrix} u\quad \epsilon \left( 1,3 \right) \\ v\quad \epsilon \left( 0,\frac { \pi }{ 2 } \right) \end{matrix}[/math].
Cambiando nuestra densidad de cartesianas a cilíndrica: [math]d\left( \rho ,\theta \right) =-1+{ e }^{ \frac { -\left| pcos\left( \theta \right) \right| }{ { \left( \rho sin\left( \theta \right) +1 \right) }^{ 2 } } }[/math].
Y sabiendo que [math]\left| { \overrightarrow { r } }_{ u }x{ \overrightarrow { r } }_{ v } \right| =\rho =u[/math].
Finalmente tenemos nuestra integral:
[math]\int _{ 1 }^{ 3 }{ \int _{ 0 }^{ \frac { \pi }{ 2 } }{ \left( (-1+{ e }^{ \frac { -\left| pcos\left( \theta \right) \right| }{ { \left( \rho sin\left( \theta \right) +1 \right) }^{ 2 } } })\rho \right) dudv } } [/math]
Al ser una integral complicada para realizarla a mano, nos apoyaremos de los métodos numéricos para aproximar la integral.

Código Matlab

Que nos da como resultado m=0.094601