Download 2008-09_AF1-PRÁCTICA-Maple 11

PRÁCTICAS DE AMPLIACIÓN DE FÍSICA I
Curso 2008-09
Basilio Carrascal Santa Olalla
Luis Seidel Gómez de Quero
Objetivos de las Prácticas
El objetivo de las Prácticas de Ampliación de Física I es doble:
1. Presentar una introducción al sistema de computación matemática Maple 11.
2. Utilizar Maple 11 para la resolución de problemas de examen de la asignatura.
Indicaciones a tener en cuenta
• Las prácticas son obligatorias para aprobar la asignatura.
• Para tener realizadas las prácticas es necesario asistir con aprovechamiento a las dos sesiones
programadas.
• La bonificación por prácticas (hasta 1 punto) se obtiene según la calificación de las memorias
entregadas en cada sesión, que deberán ser individuales y realizadas a lo largo de la misma sesión,
según las instrucciones que se den en la misma.
• Se recomienda realizar las prácticas en un documento nuevo de Maple y seguir los pasos que se
detallan en el guión.
Introducción al sistema de computación matemática Maple 11
Primeros pasos en Maple 11
Algunas ideas para empezar:
• El lenguaje de programación MAPLE requiere una sintaxis estricta, pero lógica y por tanto fácil
de entender y de familiarizarse con ella. Además, en caso de errores o sintaxis defectuosa, nos
avisa, sugiriéndonos con frecuencia la expresión correcta. MAPLE posee el comando Help (F1)
que es una ayuda completa, eficaz e inmediata.
• MAPLE distingue entre letras mayúsculas y minúsculas.
• Las expresiones que han de ser evaluadas, cuando se introducen por teclado, se ponen a la derecha
del símbolo > que se llama "prompt" y terminan o en dos puntos : , o en ;
• MAPLE realiza perfectamente todo tipo de cálculo simbólico, manipula con rapidez las fórmulas
y expresiones algebraicas y tiene una gran capacidad de representación en dos y tres dimensiones.
• Finalmente, MAPLE posee una serie de "PAQUETES" (Packages) que son ficheros que contienen
definiciones escritas en el propio lenguaje de MAPLE. Para utilizarlos, es necesario cargarlos
previamente en memoria con el comando with . Ejemplos: with(plots). with(linalg), etc.
Operaciones matemáticas elementales: ejemplos
DEFINICION DE FUNCIONES.
O f:=(5*x+7)/x^2;
f :=
5 x C7
x2
O Valor1:=subs(x=8,%);
Valor1 :=
47
64
O Valor2:=evalf(subs(x=-8.3,%));
Valor2 := 0.7343750000
RESOLUCIÓN DE SISTEMAS DE ECUACIONES.
O restart:
O a:=x+2*y+3*z+4*t=6:
O b:=5*x+5*y+4*z+3*t=1:
O c:=3*y+4*z-8*t=1:
O d:=x+y+z+t=1:
O solve({a,b,c,d});
6
40
33
x = 0, t =
,z=
, y=K
13
13
13
CALCULO DE LIMITES. limit( f(x), x=a, dirección).
O restart;
O Limit((x^2-1)/(3*x^2-3*x),x=1)=limit((x^2-1)/(3*x^2-3*x),x=
1);
x2 K1
2
lim
=
2
x/1
3
3 x K3 x
O Limit((x^2-1)/(3*x^2-3*x),x=infinity)=limit((x^2-1)/(3*x^2
-3*x),x=infinity);
x2 K1
1
lim
=
x/N
2
3
3 x K3 x
DESARROLLOS EN SERIES.
taylor( f(x), x=a,n). series( f(x), x, n).
O restart:
O f(x)=taylor(f(x),x=0,5);
1
f x = f 0 CD f 0 x C
D2 f
2
0 x2 C
1
D3 f
6
0 x3 C
1
D4 f
24
0 x4
CO x5
O ln(cos(a*x))=taylor(ln(cos(a*x)),x);
1
1 4 4
ln cos a x = K a2 x2 K
a x CO x6
2
12
O e^x=series(exp(x),x,8);
1 2
1 3
1 4
1 5
1
1
6
7
8
ex = 1 Cx C
x C
x C
x C
x C
x C
x CO x
2
6
24
120
720
5040
SUMAS a) INDEFINIDAS Y b) DEFINIDAS.
a) sum( f(n),n)
b) sum( f(n), n=a..b)
O restart:
O Sum((1+x)/(1+x^4),x=-1..15)=evalf(sum((1+x)/(1+x^4),x=-1.
.15));
15
>
1 Cx
= 2.270482282
4
x = K1 1 Cx
O Sum(((1/2)^n+(-1/2)^n),n=0..infinity)=sum(((1/2)^n+(-1/2)
^n),n=0..infinity);
N
>
n=0
1
2
n
1
2
C K
n
=
8
3
DERIVADAS DE UNA FUNCIÓN. diff( f(x,y,z,...), x,y,z..) . Derivada múltiple.
O restart;
O Diff(6*x^3*y^2/z^2,x,x,x,y,z,z)=diff(6*x^3*y^2/z^2,x$3,y,z,
z);
6 x3 y2
432 y
v6
=
2
3
2
vz vy vx
z
z4
INTEGRALES a) INDEFINIDAS y b) DEFINIDAS.
a) int( f(x),x) b) int( f(x),x=a..b). Ejemplos:
O restart:
O Int(sin(a*x)/(cos(a*x)^2),x)=int(sin(a*x)/(cos(a*x)^2),x)
;
sin a x
1
dx =
2
a cos a x
cos a x
O Int(1/sqrt(x),x=0..b)=int(1/sqrt(x),x=0..b);
b
1
0
dx = 2
b
x
O Int(1/(x*(1+ln(x))),x=1..exp(2))=int(1/(x*(1+ln(x))),x=1.
.exp(2));
e2
1
1
x 1 Cln x
dx = ln 3
Integrales dobles y triples. Ejemplos
O Int(Int(x^3-7*y*sin(3*x),x=0..6),y=-1..3)=evalf(int(int
(x^3-7*y*sin(3*x),x=0..6),y=-1..3));
3
6
x3 K7 y sin 3 x
dx dy = 1292.829623
K1 0
O Int(Int(Int(r^2*sin(theta),r=0..R),theta=0..Pi),phi=0..2*
pi)=int(int(int(r^2*sin(theta),r=0..R),theta=0..Pi),phi=
0..2*pi);
2π π R
r2 sin θ dr dθ dφ =
0
0 0
4 3
R π
3
Problema resuelto de Cálculo Vectorial
Dado el campo vectorial F en coordenadas cartesianas:
a) Comprobar que existe potencial escalar y en caso afirmativo determinarlo.
b) Comprobar que existe potencial vector y en caso afirmativo, determinarlo.
O restart;
O with(VectorCalculus):
O SetCoordinates( 'cartesian'[x,y,z]):
O X:=x/(x^2+y^2+z^2)^(3/2):
O Y:=y/(x^2+y^2+z^2)^(3/2):
O Z:=z/(x^2+y^2+z^2)^(3/2):
O F:=VectorField(<X,Y,Z>);
_
_
_
y
z
x
eC
eC
e
F :=
3/2
3/2
3/2
x
y
z
2
2
2
2
2
2
2
2
2
x Cy Cz
x Cy Cz
x Cy Cz
a) Comprobamos si el campo F es irrotacional, y en caso afirmativo, determinamos su potencial
escalar:
O RotacionalF=Curl(F);
_
RotacionalF = 0e
x
Como su rotacional es nulo, existe potencial escalar. Determinación:.
O Phi:=ScalarPotential(F);
1
Φ := K
x2 Cy2 Cz2
b) Comprobamos si el campo F es adivergente, y en caso afirmativo, determinamos su potencial
vector:
O divergenciaF:=simplify(Divergence(F));
divergenciaF := 0
Como su divergencia es nula, existe potencial vector. Determinación:.
O A:=VectorPotential(F):
4.1.2- El mismo ejercicio, pero en coordenadas esféricas.
Dado el campo vectorial G=K.ur / r2 en coordenadas esféricas:
a) Comprobar que existe potencial escalar y en caso afirmativo determinarlo.
b) Comprobar que existe potencial vector y en caso afirmativo, determinarlo.
c) Comparar este campo con el anterior.
d) ¿Son compatibles las expresiones de los potenciales escalar y vector en ambos problemas?.
O restart:
O with(VectorCalculus):
O SetCoordinates( 'spherical'[r,theta,phi]):
O G:=VectorField(<K/r^2,0,0>);
K _
G :=
e
r2 r
O RotacionalG=Curl(G);
_
RotacionalG = 0e
r
O Phi:=ScalarPotential(G);
K
r
Φ := K
O DivergenciaG:=Divergence(G);
DivergenciaG := 0
O A:=VectorPotential(G):
PRIMERA SESIÓN
Problema 1
Una distribución superficial y uniforme de carga con densidad σ0 se extiende
sobre un dominio circular de radio a. El potencial eléctrico en puntos del eje Oz
se puede expresar como
V z =
a
σ0
2 ε0
ρ
2
2
dρ
z Cρ
0
1. Determinar el potencial eléctrico V creado en puntos del eje Oz.
2. Determinar el campo eléctrico Ez en los mismos puntos.
3. Representar el potencial eléctrico y la componente z del campo eléctrico
2
K12
para z = K4 .. 4 en una misma gráfica. (σ = 1 nC/m ; ε = 8.85 10
unidades
0
0
SI, a = 1 m)
Se considera ahora una distribución de carga uniforme con la misma densidad
sobre una corona circular de radios a y b (b ! a). (b= 0.4 m).
4. Repetir los apartados 1., 2. y 3. en este caso.
Problema 2
Una distribución lineal y uniforme de carga λ se distribuye sobre una
circunferencia de radio R. El potencial eléctrico en un punto P de coordenadas P
(x,0,z) se puede expresar como
2π
λR
V x, 0, z =
4 π ε0
1
2
0
2
2
dφ
R Cx Cz K2 x R cos φ
1. Determinar el potencial eléctrico V creado en cualquier punto del espacio.
2. Obtener, a partir del resultado anterior, el potencial para los puntos del eje
Oz.
3. A partir del gradiente del potencial, obtener el vector campo en un punto
genérico P(x,0,z).
4. Representar el potencial eléctrico para algunos valores concretos de z (z= 0,
2; 0,4; 0,6; 0,8; 1) y para x de 0 a 4. (λ = 3 nC/m ; R= 1 m) en una sola
gráfica.
5. Representar el potencial eléctrico en el espacio para x=-2..2 y z=-2..2 .
6. Representar la componente z del campo eléctrico para x=-2..2 y z=-2..2 .
SEGUNDA SESIÓN
Problema 3
Dos hilos rectilineos infinitos de densidad de carga uniforme λ son paralelos
entre sí y al eje Oz. Se situan en el plano y=0 separados una distancia 2a de
manera que el hilo con densidad de carga λ pasa por el punto (-a,0,0) y el de
densidad de carga -λ por el punto (a,0,0). El potencial creado por esta
distribución de carga en el plano z=0 está dado en coordenadas polares por:
λ
V ρ, 4 =K
ln
2 π ε0
2
ρ C2 ρ a cos 4 C a2
2
ρ K2 ρ a cos 4 C a2
1. Determinar el campo eléctrico en los puntos del plano z=0.
2. Representar las curvas equipotenciales para V=-50V, -60V, -80V, entre x=-5.
.5 e y=-5..5. (λ = 3 nC/m ; a= 2 m)
3. Obtener y resolver la ecuación de las lineas de campo eléctrico.
4. Representar las líneas de campo eléctrico.
Problema 4
Un hilo no conductor rectilineo de longitud 2a cargado con una densidad
uniforme de carga λ está situado sobre el eje Oz de forma que su centro
coincide con el origen de coordenadas. El potencial que crea en cualquier punto
del espacio está dado en coordenadas cartesianas por:
V x, y, z =
λ
4 π ε0
a
1
Ka
2
2
x Cy C z Kζ
dζ
2
1. Determinar el potencial eléctrico.
2. Determinar el campo eléctrico.
3. Si se considera la solución de la ecuación de Laplace para el mismo
problema, dada para r >a en coordenadas esféricas en el plano x=0 por:
λ
a
a3
a5
V r O a, θ =
P cos θ C
P2 cos θ C
P4 cos θ C...
r 0
2 π ε0
3 r3
5 r5
donde Pi es el polinomio de Legendre de grado i, comparar ésta con la solución
exacta del apartado 1. Si λ=3 nC/m y a=0,2 m ¿cuántos términos debemos
tomar en el desarrollo en serie de la solución de la ecuación de Laplace para
que el error cometido sea inferior a 10K6 ?
Algunos problemas de Ampliación de Física I resueltos con
Maple
Problema de imágenes
Se considera una carga puntual q en el punto A(a,b,0); los planos cartesianos son conductores
y están a potencial nulo.
a) Determinación de la imagen.
b) Determinar la expresión cartesiana del potencial y del campo.
c) Hallar las expresiones de la densidad superficial de carga σ y, z y σ x, z .
d) Estudiar la posición de los puntos de máxima densidad superficial de carga inducida
respecto a la posición de la carga.
a.- Configuración de la imagen
El conjunto de imágenes que verifican las condiciones de contorno, y por consiguiente, que
resuelve el problema es:
Kq, A1 Ka, b.0 ,
q, A2 Ka, Kb, 0 ,
Kq, A3 a, Kb, 0 .
b.- Determinación de la expresión cartesiana del potencial y del campo eléctrico.
O restart:with(plots):
O V:=(x,y,z)->K[C]*(((x-a)^2+(y-b)^2+z^2)^(-1/2)+((x+a)^2+
(y+b)^2+z^2)^(-1/2)-((x-a)^2+(y+b)^2+z^2)^(-1/2)-((x+a)
^2+(y-b)^2+z^2)^(-1/2));
1
1
C
V := x, y, z /KC
x Ka 2 C y Kb 2 Cz2
x Ca 2 C y Cb 2 Cz2
1
1
K
K
2
2
2
2
x Ka C y Cb Cz
x Ca C y Kb 2 Cz2
Expresión cartesiana del campo eléctrico (Componentes).
O E[x]:=-diff(V(x,y,z),x);
1
2 x K2 a
1
2 x C2 a
Ex := KKC K
K
3/2
2
2
2
2
2
2
2 3/2
2
x Ka C y Kb Cz
x Ca C y Cb Cz
C
1
2
2 x K2 a
x Ka
2
2
2 3/2
C
C y Cb Cz
O E[y]:=-diff(V(x,y,z),y);
1
2 y K2 b
Ey := KKC K
2
x Ka 2 C y Kb 2 Cz2
C
1
2
2 y C2 b
x Ka
2
2 3/2
C y Cb 2 Cz
O E[z]:=-diff(V(x,y,z),z);
z
Ez := KKC K
x Ka 2 C y Kb 2 Cz2
C
3/2
2 x C2 a
1
2
3/2
1
2
K
x Ca
K
1
2
x Ca
2
C y Kb
2
Cz2
3/2
2 y C2 b
x Ca 2 C y Cb
2 y K2 b
2
2
Cz2
3/2
Cz2
3/2
C y Kb
2
Cz2
z
x Ca 2 C y Cb
2
3/2
z
C
x Ka
2
C y Cb
2
2 3/2
z
C
Cz
x Ca
2
C y Kb
2
Cz2
3/2
c.- Determinar la expresión de las densidades superficiales de carga inducida en los
semiplanos x = 0 e y = 0.
O sigma[y]:=subs(x=0,E[x]);
2a
2a
σy := KKC
K
3/2
3/2
a2 C y Kb 2 Cz2
a2 C y Cb 2 Cz2
y la densidad de carga inducida en el semiplano y=0:
O sigma[x]:=subs(y=0,E[y]);
2b
2b
σx := KKC
K
3/2
2
2
2 3/2
2
x Ka Cb Cz
x Ca Cb2 Cz2
d.- Cálculo de los puntos de máxima densidad inducida y representar ambas densidades en el
espacio E3 y el un plano.
Es evidente que la extremal de las densidades de carga inducida deben estar en el plano z=0,
que es donde se encuentra la carga puntual.
APLICACIÓN NUMÉRICA.
O a:=2:b:=3:
O densidad[1]:=(a^2+(y+b)^2)^(-3/2)-(a^2+(y-b)^2)^(-3/2);
1
1
densidad1 :=
K
3/2
2
2 3/2
4 C y C3
4 C y K3
O densidad[2]:=((x+a)^2+b^2)^(-3/2)-((x-a)^2+b^2)^(-3/2);
1
1
densidad2 :=
K
3/2
3/2
x C2 2 C9
x K2 2 C9
O plot3d({-5*densidad[1],-5*densidad[2]},x=0..10,y=0..10,
thickness=1,axes=boxed);
0
y
5
0,6
10
0
0,4
0,2
5
x
0,0 10
O diff(densidad[1],y):máximo[y]:=fsolve(%,y);
máximoy := 3.018771410
O diff(densidad[2],x):máximo[x]:=fsolve(%,x);
máximox := 2.273212314
O
O
O
O
O
a:=plot(-10^2*densidad[1],y=0..10,symbol=cross):
b:=plot(-10^2*densidad[2],x=0..10,symbol=diamond):
c:=textplot([4,1,`Densidad 2`]):
d:=textplot([5,8,`Densidad 1`]):
display([a,b,c,d],axes=boxed, title=`Densidades
inducidas`);
Densidades inducidas
12
10
Densidad 1
8
6
4
2
Densidad 2
0
0
2
4
6
8
10
Problema propuesto sobre la Ecuación de Laplace
Se considera el campo eléctrico plano que determinan los conductores que limitan un sector
circular de radio R=1 y apertura de ángulo α. Utilizando por conveniencia, coordenadas polares,
las condiciones de contorno son:
a) V(ρ,0)=V(ρ,α)=0
b) V R, φ = V1.
1.- Determinar el potencial en el sector circular.
2.- Determinar el campo eléctrico en el sector circular.
3.- Obtener la densidad de carga en los conductores OA y OB (a potencial nulo) y en AB (a
potencial V1).
4.- Representar la densidad de carga sobre OA y OB, en función de ρ, y la densidad de carga en
AB, en función de φ.
Inicialización
O restart:with(DEtools):addcoords(cylindrical,[r,theta,z],
[r*cos(theta),r*sin(theta),z]):
Warning, coordinates already exists, system redefined.
1.- Cálculo del potencial por separación de variables.
El potencial electrostático, en el sector circular, verifica la ecuación de Laplace. en
coordenadas cartesianas,
O Diff(V,x,x)+Diff(V,y,y)=0;
v2
v2
V
C
V=0
vx2
vy2
y en coordenadas polares
O Diff(V,rho,rho)+1/rho*Diff(V,rho)+1/rho^2*Diff(V,phi,phi)
=0;
v2
v
V
V
2
vρ
vφ
v2
VC
C
=0
2
2
ρ
ρ
vρ
El método de separación de variables determina soluciones al potencial de la forma:
O U:=(rho,phi)->R(rho)*Phi(phi);
U := ρ, φ /R ρ Φ φ
que sustituyendo en la ecuación de Laplace que:
O rho^2*Diff(R,rho,rho)/R+rho*Diff(R,rho)/R+Diff(Phi,phi,
phi)/Phi=0;
2
v2
v2
v
ρ
R
Φ
ρ
R
2
2
vρ
vρ
vφ
C
C
=0
R
R
Φ
De acuerdo con las técnicas del método se separación de variables, tenemos dos ecuaciones
diferenciales que serán:
O eq1:=rho^2*Diff(R,rho,rho)/R+rho*Diff(R,rho)/R=-K;
eq2:=Diff(Phi,phi,phi)/Phi=K;
2
v2
v
ρ
R
ρ
R
2
vρ
vρ
eq1 :=
C
= KK
R
R
v2
Φ
2
vφ
eq2 :=
=K
Φ
Resolvemos en primer lugar la eciación eq2
O dsolve(diff(Phi(phi),phi,phi)/Phi(phi)=K);
Φ φ = _C1 e
K φ
K φ
C_C2 eK
Como Φ debe verificar las condiciones de contorno: potencial nulo para φ = 0 y φ = α,
2
debe ser K ! 0Haciendo K = Kω
O Phi[1]:=dsolve(diff(Phi(phi),phi,phi)/Phi(phi)=-omega^2);
Φ1 := Φ φ = _C1 sin ω φ C_C2 cos ω φ
La verificación de las condiciones de contorno sobre la variable φ, se cumplen bien,
nπ
haciendo:ω =
, _C1 = 1 y _C2 = 0, siendo n un entero positivo. Es decir:
α
O Phi[1]:=sin(n*Pi*phi/alpha);
nπφ
Φ1 := sin
α
Resolución de la primera ecuación diferencial. Empezamos por el caso n = 0.
O dsolve(rho^2*diff(R(rho),rho,rho)/R(rho)+rho*diff(R(rho),
rho)/R(rho)=0);
R ρ = _C1 C_C2 ln ρ
Como esta solución presenta una singularidad para ρ = 0, debemos descartarla. Por
consiguiente n debe ser un número natural. Para cualquier valor de n, se tiene
O R1:=dsolve(rho^2*diff(R(rho),rho,rho)/R(rho)+rho*diff(R
(rho),rho)/R(rho)=(n*Pi/alpha)^2);
1
2
R1 := R ρ = K
ρ
2nπ
α
_C1 K_C2 α
ρ
nπ
α
πn
Para evitar singularidades en el origen, hacemos _C1 = 1 y _C2 = 0. O sea:
O simplify(subs(_C1=1,_C2=0,R1));
nπ
α
α
1 ρ
2
πn
R ρ =K
Entonces la función potencial es el producto de R ρ y Φ φ .
O U:=rho^(n*Pi/alpha)/(2*n*Pi)*sin(n*Pi*phi/alpha);
U :=
1
2
ρ
nπ
α
nπφ
α
sin
πn
En una solución de este tipo no es imponible la condición de contorno V R, φ = V1 para
todo φ. Utilizando la metodología al uso en la Teoría de Ecuaciones Diferenciales, se tiene:
O Sum(C[n]*R^(n*Pi/alpha)*sin(n*Pi*phi/alpha),n=0..
infinity)=V[1];
N
>C R
n=0
nπ
α
n
sin
nπφ
α
= V1
Comprobar que esta expresión es una combinación de términos linealmente independientes
(las funciones de φ son ortogonales).
O assume(n,integer):assume(m,integer):
O Int(sin(n*Pi*phi/alpha)*sin(p*Pi*phi/alpha),phi=0..alpha)
=int(sin(n*Pi*phi/alpha)*sin(m*Pi*phi/alpha),phi=0..
alpha);
α
sin
0
n~ π φ
α
sin
pπφ
α
dφ = 0
Determinación de los coeficientes Cn . Se utiliza el método de Fourier.
O V[1]*int(sin(m*Pi*phi/alpha),phi=0..alpha);
V1 α K1 C K1 m~
K
m~ π
O C[m]*R^(m*Pi/alpha)*int(sin(m*Pi*phi/alpha)^2,phi=0..
alpha);
1
C R
2 m~
m~ π
α
α
Se comprueba que los términos de m par se tiene que Cm = 0, cualquiera que sea m.
Ponemos entonces: m = 2 n C1
O V[1]*int(sin((2*n+1)*Pi*phi/alpha),phi=0..alpha)=C[2*n+1]
*R^((2*n+1)*Pi/alpha)*int(sin((2*n+1)*Pi*phi/alpha)^2,
phi=0..alpha);
2 V1 α
1
=
C
R
2 2 n~ C 1
2 n~ C1 π
2 n~ C 1 π
α
α
y de esta expresión obtenemos el valor de las constantes. En conclusión el potencial
eléctrico viene dado por:
O restart:with(plots):with(linalg):
O addcoords(z_cylindrical,[z,r,theta],[r*cos(theta),r*sin
(theta),z]):
O V(rho,phi)=4*V[1]/Pi*Sum(R^(-(2*n+1)*Pi/alpha)*rho^((2*
n+1)*Pi/alpha)/((2*n+1)*Pi)*sin((2*n+1)*Pi*phi/alpha),n=
0..N);
K
N
V ρ, φ =
4 V1
>
R
2 n C1 π
α
ρ
2 n C1 π
α
sin
2 n C1 π φ
α
2 n C1 π
π
n=0
π
y V1 = 1 :
4
O term:=4*R^(-(2*n+1)*Pi/alpha)*rho^((2*n+1)*Pi/alpha)/((2*
n+1)*Pi)*sin((2*n+1)*Pi*phi/alpha);
Representación de V ρ, φ , para N=200, R=1 , α =
K
4R
2 n C1 π
α
term :=
ρ
2 n C1 π
α
sin
2 n C1 π φ
α
2 n C1 π
O N:=200:R:=1:alpha:=Pi/4:
O potencial:=sum(subs({alpha=Pi/4,phi=theta,rho=r},term),n=
0..N):
O plot3d(potencial,r=0..R,theta=0..alpha,coords=
z_cylindrical,axes=BOXED,style=WIREFRAME,title="Potencial
en el sector circular",orientation=[-110,62]);
Potencial en el sector circular
1,0
0,75
0,5
0,4
0,25
0,0
0,0
0,0
0,5
1,0
Las líneas equipotenciales se visualizan por la sentencia:
O contourplot(potencial,r=0..R,theta=0..alpha,coords=
z_cylindrical,contours=50);
0,7
0,6
0,5
0,4
0,3
0,2
0,1
0,4 0,5 0,6 0,7 0,8 0,9 1,0
2.- Determinación del campo eléctrico
Obtención de las componentes del campo eléctrico en coordenadas polares (a partir del
término general). Reiniciamos otra vez el problema para controlar mejor el flujo del
programa.
O restart:with(linalg):with(plots):
O term:=4*R^(-(2*n+1)*Pi/alpha)*rho^((2*n+1)*Pi/alpha)/((2*
n+1)*Pi)*sin((2*n+1)*Pi*phi/alpha);
K
4R
2 n C1 π
α
ρ
2 n C1 π
α
term :=
sin
2 n C1 π φ
α
2 n C1 π
O E_component:=grad(-term,[rho,phi],coords=polar);
K
4R
2 n C1 π
α
ρ
2 n C1 π
α
E_component := K
K
4R
2 n C1 π
α
K
sin
αρ
ρ
2 n C1 π
α
cos
2 n C1 π φ
α
,
2 n C1 π φ
α
ρα
3 y 4.- Determinación y representación de las funciones densidad de carga.
La función densidad de carga en un conductor viene expresada por σ = ε0 E P . En
particular en este problema tenemos tres funciones densidad: σOA = ε0 E φ = 0 ,
σOB = ε0 E φ = α
y σAB = ε0 E ρ = R .
Comprobar que σOA = σOB para cualquier valor de ρ y ambas densidades son negativas.
Representar las funciones σOA y σAB.
O N:=400:alpha:=Pi/4:R:=1:
O densOA:=evalf(sum(subs(phi=0,E_component[2]),n=1..N)):
O densOB:=evalf(-sum(subs(phi=alpha,E_component[2]),n=1..N)
):
O diferencia:=densOA-densOB;
diferencia := 0
O
O
O
O
densAB:=evalf(sum(subs(rho=R,E_component[1]),n=1..N)):
den_OA:=plot(densOA,rho=0..R,y=-30..0,color=blue):
den_AB:=plot(densAB,phi=0..alpha,color=red):
display(den_OA,axes=boxed,title="Densidades en los
conductores OA y OB",thickness=2);
Densidades en los conductores OA y OB
0
K10
K20
K30
0,2
0,4
0,6
0,8
1,0
O display(den_AB,axes=boxed,title="Densidad en el conductor
AB");
Densidad en el conductor AB
0
K200
K400
K600
K800
K1.000
K1.200
0,0 0,1 0,2 0,3 0,4 0,5 0,6 0,7