Matlab Pdetool Lezione 3 Magnetostatica · Equazioni integrali magnetostatica. Formulazione...

Matlab Pdetool

Lezione 3Magnetostatica

Ing. Flavio Calvano

• Modello magnetostatica;

• Utilizzo del pdetool in modalità grafica;

• Esercitazione di laboratorio

Calcolo auto-induttanza;

Calcolo mutua-induttanza;

Confronto dati numerici-sperimentali.

Magnetostatica

Equazioni integrali magnetostatica

Formulazione differenziale

coulombdigauge0

''

0

S0S0

A

AABAA

JAJB

ABB

Il potenziale vettore A è univocamente

determinato dalla gauge di coulomb

Caso 2d assial-simmetrico

JeB

z

Be

Je)z,(

0z

0B

Caso 2d assial-simmetrico

L’incognita nel caso assial-simmetrico

diventa uno scalare

J)A

z(

z

1)A(

10

J

eAeJ

0con)(2

AJ

AAAA

Pdetool

Creazione geometria

xmin=0;

xmax=2.2*10^-2;

dh=0.04;

dr=0.0025

pderect([xmax-dr xmax -dh dh],'R1');

pderect([xmin 10*xmax 10*(-dh) 10*dh],'R2');

Condizioni al contornopdesetbd(8,...

'dir',...

1,...

'1',...

'0')

pdesetbd(5,...

'dir',...

1,...

'1',...

'0')

pdesetbd(4,...

'dir',...

1,...

'1',...

'0')

pdesetbd(3,...

'dir',...

1,...

'1',...

'0')

Mesh

R2=box aria ->t(4,:)=1R1=sonda ->t(4,:)=2

Dall’indice della

quarta riga della

matrice t è possibile

riconoscere i domini

e imporre il termine

noto J

0 0.05 0.1 0.15 0.2 0.25-0.25

-0.2

-0.15

-0.1

-0.05

0

0.05

0.1

0.15

0.2

0.25

set(findobj(get(pde_fig,'Children'),'Tag','PDEEval'),'String','R2+R1');

gd=get(findobj(get(pde_fig,'Children'),'flat‘,'Tag','PDEMeshMenu'),'UserData');

dl=decsg(gd);

[p,e,t]=initmesh(dl,'Hmax',5e-1,'init','off');

[p,e,t]=refinemesh(dl,p,e,t,'regular');

h=findobj(get(pde_fig,'Children'),'flat','Tag','PDEBoundMenu');

bl=get(findobj(get(h,'Children'),'flat','Tag','PDEBoundMode'),'UserData');

Definizione del forzamento

Nspire=500;

dr=0.002;

dh=0.08;

I=1;

J0=Nspire*I/(dh*dr);

%nell’area occupata dalla bobina

JAA 0)()(

zz

Calcolo Soluzione

for n=1:size(t,2)

ni=t(1,n); % indici dei nodi del triangolo n-mo

nj=t(2,n);

nk=t(3,n);

xi=p(1,ni); % coordinate dei vertici del tr. n-mo

xj=p(1,nj);

xk=p(1,nk);

rc(n)=(xi+xj+xk)/3; % baricentro del tr. n-mo

if t(4,n)==1

J(n)=0;

else

J(n)=J0;

end

end J definita solo in R2 (t(4,:)=2)

corrispondente all’area della

bobina, nulla all’esterno

A = assempde(bl,p,e,t,c,'0',f);

-div∙c(grad(A))+aA=f

Con:

c=rcf= rc.*(4e-7*pi).*J;

Calcolo del Campo magnetico

Br=-d(Aphi)/dz

Bz=(d(Aphi*r)/dr)/r

r=p(1,:);

z=p(2,:);

u=A.*r';

[g,dAdz]=pdegrad(p,t,A);

Br=-dAdz;

[drAdr,g2]=pdegrad(p,t,u);

Bz=(1./rc).*drAdr;

Plot linee di Campo

figure

pdegplot(dl)

hold on

xm=min(p(1,t)); xma=max(p(1,t));

ym=min(p(2,t)); yma=max(p(2,t));

na=100;

x=linspace(xm,xma,na);

y=linspace(ym,yma,na);

flowdata1=pdeprtni(p,t,Br);

flowdata2=pdeprtni(p,t,Bz);

uu=tri2grid(p,t,flowdata1,x,y);

vv=tri2grid(p,t,flowdata2,x,y);

[msg,x,y]=xyzchk(x,y,uu,vv);

scale=2;

oks=find(~isnan(uu));

quiver(x(oks),y(oks),-uu(oks),-vv(oks),scale,'r-');

0 0.005 0.01 0.015 0.02 0.025

-0.1

-0.05

0

0.05

0.1

Campo Magnetico

Plot linee di flusso

-0.05 0 0.05 0.1 0.15

-0.08

-0.06

-0.04

-0.02

0

0.02

0.04

0.06

0.08

rA2dlndVdV

SS

AAnBΦ pdecont(p,t,A.*p(1,:)*2*pi,10)

S

Esempio 2

pderect([xmax-dr xmax ymax-dh ymax],'R1');

pderect([xmax-dr xmax ymax-2*dh-h ymax-1*dh-h],'R2');

pderect([xmax-dr xmax ymax-3*dh-2*h ymax-2*dh-2*h],'R3');








pderect([0 10*xmax -5*xmax 5*xmax],'R11');

xmax=1*10^-2;

dh=0.001;

dr=0.001;

ymax=0.02;

h=0.003;

Condizioni al contornopdetool('changemode',0)

pdesetbd(5,...

'dir',...

1,...

'1',...

'0')

pdesetbd(6,...

'dir',...

1,...

'1',...

'0')

pdesetbd(3,...

'dir',...

1,...

'1',...

'0')

pdesetbd(2,...

'dir',...

1,...

'1',...

'0')

Assemblaggio termine notofor n=1:size(t,2)


nj=t(2,n);

nk=t(3,n);


nj=t(2,n);

nk=t(3,n);


xj=p(1,nj);

xk=p(1,nk);


I1=1;

if t(4,n)==1

J(n)=0;

else

J(n)=I1/(dh*dr);

end

end

J uguale a zero in R11 (t(4,:)=1)

corrispondente al box aria ,

mentre è diversa da zero negli

avvolgimenti che hanno tutte la

stessa corrente

A = assempde(bl,p,e,t,c,'0',f);

Con:c=rcf= rc.*(4*pi*1e-7).*J;

Campo magnetico e flusso

5 6 7 8 9 10 11 12 13

x 10-3

0.015

0.016

0.017

0.018

0.019

0.02

0.021

Campo Magnetico

Campo e Flusso tra le spire

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1

-0.04

-0.03

-0.02

-0.01

0

0.01

0.02

0.03

0.04

Esercitazione laboratorio

Definizione Geometria

xmin=0;

xmax=0.018;

dh=0.044;

dr=0.002;

xmax2=0.025;

pderect([xmax-dr/2 xmax+dr/2 -dh/2 dh/2],'R1');

pderect([xmax2-dr/2 xmax2+dr/2 -dh/2 dh/2],'R2');

pderect([20*xmin 20*xmax 10*(-xmax) 10*xmax],'R3');

set(findobj(get(pde_fig,'Children'),'Tag','PDEEval'),'String',‘R3+R1+R2');

t(4,:)=1 t(4,:)=2 t(4,:)=3

Nspire=500;

I1=0.0021; %dati laboratorio

I2=0.0013;

J1=zeros(1,size(t,2));

J2=zeros(1,size(t,2));

for n=1:size(t,2)


nj=t(2,n);

nk=t(3,n);


xj=p(1,nj);

xk=p(1,nk);


if t(4,n)==2

J1(n)=Nspire*I1/(dh*dr);

elseif t(4,n)==3

J2(n)= Nspire*I2/(dh*dr);

end %endif

end %endfor

A1 = assempde(bl,p,e,t,rc.*(1/(4*pi*1e-7)),'0',rc.*J1,'0');

A2= assempde(bl,p,e,t,rc.*(1/(4*pi*1e-7)),'0',rc.*J2,'0');

Definizione forzamento

Calcolo flusso e induttanza

21222

211m

V V

22

V

2211

V

12111

V

22211211

V

2121

V

m

Vsol

m

MiiiL2

1iL

2

1W

.dVJA2

1dVJA

2

1dVJA

2

1dVJA

2

1

dV)JAJAJAJA(2

1

dV)JJ()AA(2

1JdVA

2

1W

AdVJ2

1W

1 221

Calcolo AUTO-INDUTTANZA11

S

11121

1

V

1121

1 dzdJAi

2dVJA

i

1L

11

A1 = assempde(bl,p,e,t,rc.*(1/(4*pi*1e-7)),'0',rc.*J1,'0');

AT1=pdeintrp(p,t,A1);

n1=find(t(4,:)==2);

for n=n1


nj=t(2,n);

nk=t(3,n);


xj=p(1,nj);

xk=p(1,nk);


D=[1 z1i z2i; 1 z1j z2j; 1 z1k z2k];

area=det(D)/2;

aa=(J1(n)*AT1(n));

ce(n)=(rc)*(aa)*area;

end

Lenergia=(1/I1^2)*sum(ce)*2*pi

Calcolo MUTUA-INDUTTANZA

222

S

21

21

2

V

21

21

21 dzdJAii

2dVJA

ii

1M

22

AT1=pdeintrp(p,t,A1); %A1 calcolato imponendo J1 in assempde

n2=find(t(4,:)==3);% trova i triangoli del dominio della seconda bobina

for n=n2


nj=t(2,n);

nk=t(3,n);


xj=p(1,nj);

xk=p(1,nk);


D=[1 z1i z2i; 1 z1j z2j; 1 z1k z2k];

area=det(D)/2;

aa=(J2(n)*AT1(n));

ce(n)=(rc)*(aa)*area;

end

Menergia=(1/(I1*I2))*sum(ce)*2*pi

Confronto dati numerici-sperimentali

M21energia = 14.6 mHM21misurata=14.8 mH

L11energia = 10 .6 mHL11misurata=10.6 mH

L22energia= 16.4 mHL22misurata= 16 .0 mH

M12energia = 14.6 mHM12misurata=14.6 mH

Matlab Pdetool Lezione 3 Magnetostatica · Equazioni integrali magnetostatica. Formulazione...

Documents

Transcript of Matlab Pdetool Lezione 3 Magnetostatica · Equazioni integrali magnetostatica. Formulazione...