×
INTELLIGENT WORK FORUMS
FOR ENGINEERING PROFESSIONALS

Log In

Come Join Us!

Are you an
Engineering professional?
Join Eng-Tips Forums!
  • Talk With Other Members
  • Be Notified Of Responses
    To Your Posts
  • Keyword Search
  • One-Click Access To Your
    Favorite Forums
  • Automated Signatures
    On Your Posts
  • Best Of All, It's Free!
  • Students Click Here

*Eng-Tips's functionality depends on members receiving e-mail. By joining you are opting in to receive e-mail.

Posting Guidelines

Promoting, selling, recruiting, coursework and thesis posting is forbidden.

Students Click Here

Jobs

Gradual changes in arrows acording to a equation

Gradual changes in arrows acording to a equation

Gradual changes in arrows acording to a equation

(OP)
Hi, I'm trying to simulate Electromagnetic Fields in a waveguide. The problem is that, the arrows of the representation doesn't change in a gradual way.  
Perhaps the problem is that quiver3.m doesn't take N-dimensional arrays like parameters.
I attach the code of two programs, the first one is mathematically wrong, but shows the result request, and the second one is mathematically correct. I don't know the reason why the second one doesn't represent the function as the first one.
Can anybody help me? thanks

FIRST ONE
%function TE(m,n,Reldim,epsilonr,mur,Relfrec)
close all


m=2;n=2;Reldim=2;epsilonr=2;mur=1;Relfrec=3;%Reldim es la relación entre a y b(anchura y altura respectivamente)
%Relfrec es la relación entre la frecuencia de corte y la frecuencia de trabajo.Tiene que ser >1 para que el modo se propague.
%Reldim=a/b;
%Relfrec=ftrabajo/fcorte;

Eruptura=10;%La verdad es que lo pongo por ponerlo.
c=3*10^8;%Velocidad de la luz en espacio libre.
mu0=4*pi*10^(-7);
epsilon0=1/(mu0*c^2);

%Normalmente, Reldim será siempre>=1 porque la guía es siempre más ancha que alta.

if Reldim>=1
   anchura=1;altura=anchura/Reldim;
else
   altura=1;anchura=Reldim*altura;
end

%Este if lo he hecho con el fin de que al hacer el dibujo de la guía, sus dimensiones esten normalizadas.
%La anchura=1 y altura<1 ó altura=1 y anchura<1;

mu=mur*mu0;
epsilon=epsilonr*epsilon0;
fcorte=(sqrt((m*pi)^2+(n*pi*Reldim)^2))/(anchura*2*pi*sqrt(mu*epsilon));%Ecuación que relaciona la anchura de la guía con la frecuencia de corte.
f=Relfrec*fcorte;
w=2*pi*f;
K=w*sqrt(mu*epsilon);
Kx=m*pi/anchura;
Kz=n*pi/altura;
Kc=((Kx)^2+(Kz)^2)^(1/2);
fcorte1=Kc/(2*pi*sqrt(mu*epsilon));%Frecuencia de corte por debajo de la cual los modos se devanecen rapidamenete.
BETA=K*sqrt(1-(fcorte/f)^2);%Constante de fase.
Zte=(i*w*mu)/(i*BETA);%Impedancia de la onda.
Vp=1/sqrt(mu*epsilon);%Velocidad de propagación en ese medio.
Lambda=(2*pi)/K;%Longitud de onda.
LambdaG=(2*pi)/BETA;%Longitud de onda en la guía.
LambdaC=Vp/fcorte;%Longitud de onda de corte.
longitud=LambdaG;%Los cortes longitudinales de la guía los vamos a representar para esa longitud.

B=Eruptura*Kc^2/(Zte*BETA*Kx);%Aqui no sé si tengo que dejar la 'i'y si hay un '-'.
Cte=-Zte*i*BETA*B*Kx/(Kc^2);
Cte1=Zte*i*BETA*B*Kz/(Kc^2);

vectorX=linspace(0,anchura,20)';%Lo pongo en forma de columna. Esto se puede hacer siempre que no tenga parte compleja ya que sino me devolvería el traspuesto conjugado(le cambiaría el signo a la parte compleja).
vectorY=linspace(0,longitud,length(vectorX))';
vectorZ=linspace(0,altura,length(vectorX))';
uso=ones(1,length(vectorX));
for cont=1:length(vectorX);
   x(:,:,cont)=vectorX*uso;
   y(:,:,cont)=(vectorY*uso)';
   z(:,:,cont)=vectorZ(cont)*ones(length(vectorX),length(vectorX));
end

%Ya he creado los ejes X,Y y Z.

%Ex=REAL[Cte*sin(Kz*z)*cos(Kx*x)*imag(exp(-i*w*t).*exp(i*BETAte10*y))];
%Ez=REAL[Cte*sin(Kx*x)*cos(Kz*z)*imag(exp(-i*w*t).*exp(i*BETAte10*y))];


miembro1=sin(Kx*vectorX);%Este es el sin(Kx*x).
miembro2=(exp(i*BETA*vectorY));%exp(i*BETAte10*y)
miembro2=transpose(miembro2);
miembro3=cos(Kz*vectorZ);%Este es el cos(Kz*z).

mem1=miembro1*miembro2;

miembro4=sin(Kz*vectorZ);%Este es el sin(Kz*z).
miembro5=cos(Kx*vectorX);%Este es el cos(Kx*x).

mem2=miembro5*miembro2;%cos(Kx*x)*exp(i*BETAte10*y)

t=linspace(0,LambdaG/Vp,length(vectorX));

for cont=1:length(vectorX);
   Ex(:,:,cont)=real(Cte1*mem2.*miembro4(cont)*exp(-i*w*t(cont)));
   Ey(:,:,cont)=zeros(length(vectorX),length(vectorX));
   Ez(:,:,cont)=real(Cte*mem1.*miembro3(cont)*exp(-i*w*t(cont)));
end

M=moviein(20);
   set(gca,'NextPlot','replacechildren');%Este de aqui no sé si sirve.

for j=1:length(vectorX);
   
   
   %figure(1);
   %set(gca,'Units','normalized');
   %scrnsize=get(0,'ScreenSize');
   %set(gca,'RendererMode','manual','Position',[1 1 1024 700]);
   
   
   
   set(gca,'NextPlot','replacechildren'); %Esto aqui es importantísimo.
   set(gca,'Xlim',[0 anchura],'Ylim',[vectorY(j)-longitud/1.25 vectorY(j)+longitud/1.25],'Zlim',[0,altura]);
   
   %En modos muy superiores, lo que hay que hacer es aumentar la distancia desde la que se ven las flechas,
   %mediante la orden que viene dentro del set() => [vectorY(j)-longitud/2 vectorY(j)+longitud/2]. En vez de 2 ponemos 1.25,
   %y además, habría que poner cada vez más puntos con el fin de que se vean bien las 'circulaciones'.
   
    quiver3(x(:,j,:),y(:,j,:),z(:,j,:),Ex(:,j,:),Ey(:,j,:),Ez(:,j,:))%Esto es para un corte.
 
   
   %Como e cubo es de 15*15*15, solo puedo pedir que me enseñe 15 capas en el quiver3.
   xlabel('Eje x');
   ylabel('Eje y');
   zlabel('Eje z');
   view(19,19);
    grid on;
    box on;
   hold off;
   M(:,j)=getframe;
    close all;%He puesto este close all porque no sé porqué el replace children no funciona. Resulta que me hace una figura por cada frame que hago.
 end;
 
 %Para el campo magnético.

Hx=real(Ez./(-Zte));
for cont=1:length(vectorX);
     Hy(:,:,cont)=real(B*mem2.*miembro3(cont));
end
Hz=real(Ex./Zte);



SECOND ONE
%function TE(m,n,Reldim,epsilonr,mur,Relfrec)
close all


m=2;n=2;Reldim=2;epsilonr=2;mur=1;Relfrec=3;%Reldim es la relación entre a y b(anchura y altura respectivamente)
%Relfrec es la relación entre la frecuencia de corte y la frecuencia de trabajo.Tiene que ser >1 para que el modo se propague.
%Reldim=a/b;
%Relfrec=ftrabajo/fcorte;

Eruptura=10;%La verdad es que lo pongo por ponerlo.
c=3*10^8;%Velocidad de la luz en espacio libre.
mu0=4*pi*10^(-7);
epsilon0=1/(mu0*c^2);

%Normalmente, Reldim será siempre>=1 porque la guía es siempre más ancha que alta.

if Reldim>=1
   anchura=1;altura=anchura/Reldim;
else
   altura=1;anchura=Reldim*altura;
end

%Este if lo he hecho con el fin de que al hacer el dibujo de la guía, sus dimensiones esten normalizadas.
%La anchura=1 y altura<1 ó altura=1 y anchura<1;

mu=mur*mu0;
epsilon=epsilonr*epsilon0;
fcorte=(sqrt((m*pi)^2+(n*pi*Reldim)^2))/(anchura*2*pi*sqrt(mu*epsilon));%Ecuación que relaciona la anchura de la guía con la frecuencia de corte.
f=Relfrec*fcorte;
w=2*pi*f;
K=w*sqrt(mu*epsilon);
Kx=m*pi/anchura;
Kz=n*pi/altura;
Kc=((Kx)^2+(Kz)^2)^(1/2);
fcorte1=Kc/(2*pi*sqrt(mu*epsilon));%Frecuencia de corte por debajo de la cual los modos se devanecen rapidamenete.
BETA=K*sqrt(1-(fcorte/f)^2);%Constante de fase.
Zte=(i*w*mu)/(i*BETA);%Impedancia de la onda.
Vp=1/sqrt(mu*epsilon);%Velocidad de propagación en ese medio.
Lambda=(2*pi)/K;%Longitud de onda.
LambdaG=(2*pi)/BETA;%Longitud de onda en la guía.
LambdaC=Vp/fcorte;%Longitud de onda de corte.
longitud=LambdaG;%Los cortes longitudinales de la guía los vamos a representar para esa longitud.

B=Eruptura*Kc^2/(Zte*BETA*Kx);%Aqui no sé si tengo que dejar la 'i'y si hay un '-'.
Cte=-Zte*i*BETA*B*Kx/(Kc^2);
Cte1=Zte*i*BETA*B*Kz/(Kc^2);

vectorX=linspace(0,anchura,20)';%Lo pongo en forma de columna. Esto se puede hacer siempre que no tenga parte compleja ya que sino me devolvería el traspuesto conjugado(le cambiaría el signo a la parte compleja).
vectorY=linspace(0,longitud,length(vectorX))';
vectorZ=linspace(0,altura,length(vectorX))';
uso=ones(1,length(vectorX));
for cont=1:length(vectorX);
   x1(:,:,cont)=vectorX*uso;
   y1(:,:,cont)=(vectorY*uso)';
   z1(:,:,cont)=vectorZ(cont)*ones(length(vectorX),length(vectorX));
end

%Ya he creado los ejes X,Y y Z.

%Ex=REAL[Cte*sin(Kz*z)*cos(Kx*x)*imag(exp(-i*w*t).*exp(i*BETAte10*y))];
%Ez=REAL[Cte*sin(Kx*x)*cos(Kz*z)*imag(exp(-i*w*t).*exp(i*BETAte10*y))];


miembro1=sin(Kx*vectorX);%Este es el sin(Kx*x).
miembro2=(exp(i*BETA*vectorY));%exp(i*BETAte10*y)
miembro2=transpose(miembro2);
miembro3=cos(Kz*vectorZ);%Este es el cos(Kz*z).

mem1=miembro1*miembro2;

miembro4=sin(Kz*vectorZ);%Este es el sin(Kz*z).
miembro5=cos(Kx*vectorX);%Este es el cos(Kx*x).

mem2=miembro5*miembro2;

%t=vectorY./Vp;
t=linspace(0,LambdaG/Vp,length(vectorX));

for cont=1:length(vectorX);
   Ex1(:,:,cont)=real(Cte1*mem2.*miembro4(cont));
   Ey1(:,:,cont)=zeros(length(vectorX),length(vectorX));
   Ez1(:,:,cont)=real(Cte*mem1.*miembro3(cont));
end

for cont=1:length(vectorX);
   Ex(:,:,:,cont)=real(Ex1(:,:,:).*exp(-i*w*t(cont)));
   Ey(:,:,:,cont)=real(Ey1(:,:,:).*exp(-i*w*t(cont)));
   Ez(:,:,:,cont)=real(Ez1(:,:,:).*exp(-i*w*t(cont)));
   x(:,:,:,cont)=x1(:,:,:);
   y(:,:,:,cont)=y1(:,:,:);
   z(:,:,:,cont)=z1(:,:,:);
end

%Aqui lo que he creado, son cubos de 20*20*20*20. La primera dimensión es para la X, la segunda para la Y,
%la tercera para la Z y la cuarta para el tiempo t.

M=moviein(20);
  % set(gca,'NextPlot','replacechildren');%Este de aqui no sé si sirve.

for j=1:length(vectorX);
   
   
   %figure(1);
   %set(gca,'Units','normalized');
   %scrnsize=get(0,'ScreenSize');
   %set(gca,'RendererMode','manual','Position',[1 1 1024 700]);
   
   
   
   set(gca,'NextPlot','replacechildren'); %Esto aqui es importantísimo.
   set(gca,'Xlim',[0 anchura],'Ylim',[vectorY(5)-longitud/1.25 vectorY(5)+longitud/1.25],'Zlim',[0,altura]);
   
   %En modos muy superiores, lo que hay que hacer es aumentar la distancia desde la que se ven las flechas,
   %mediante la orden que viene dentro del set() => [vectorY(j)-longitud/2 vectorY(j)+longitud/2]. En vez de 2 ponemos 1.25,
   %y además, habría que poner cada vez más puntos con el fin de que se vean bien las 'circulaciones'.
   xx=x(:,5,:,j)
    quiver3(x(:,5,:,j),y(:,5,:,j),z(:,5,:,j),Ex(:,5,:,j),Ey(:,5,:,j),Ez(:,5,:,j))%Esto es para un corte.
   
keyboard
   %Como e cubo es de 15*15*15, solo puedo pedir que me enseñe 15 capas en el quiver3.
   xlabel('Eje x');
   ylabel('Eje y');
   zlabel('Eje z');
   view(19,19);
    grid on;
    box on;
   %hold off;
   M(:,j)=getframe;
    close all;%He puesto este close all porque no sé porqué el replace children no funciona. Resulta que me hace una figura por cada frame que hago.
 end;
 

%Para el campo magnético.

Hx=real(Ez./(-Zte));
Hz=real(Ex./Zte);

for cont=1:length(vectorX);
     Hy1(:,:,cont)=real(B*mem2.*miembro3(cont));
end
for cont=1:length(vectorX);
     Hy(:,:,:,cont)=real(Hy1.*exp(-i*w*t(cont)));
end

Red Flag This Post

Please let us know here why this post is inappropriate. Reasons such as off-topic, duplicates, flames, illegal, vulgar, or students posting their homework.

Red Flag Submitted

Thank you for helping keep Eng-Tips Forums free from inappropriate posts.
The Eng-Tips staff will check this out and take appropriate action.

Reply To This Thread

Posting in the Eng-Tips forums is a member-only feature.

Click Here to join Eng-Tips and talk with other members!


Resources