登录论坛

查看完整版本 : Matlab:用biot-savarts-law计算的磁场


poster
2019-12-10, 20:30
这是使用毕奥-萨伐尔定律计算磁场的代码。我希望获得一些优化该代码的技巧。遗憾的是,我使用德语:(我再也不会这样做了。:)

tic clear all; clc; clf skalierungsfaktor = 10^-6; % vom m-Bereich zum mm-Bereich wg. T = Vs / m^2 I = 1; % in A konstante = 10^-10; % mu0 / (4 pi) in Vs / (A mm) % Faktor 10^3 kleiner wg. mm statt m matrix = 30; B = zeros(200,200,200,5); zaehler = 0; % f眉r Zeitmessung xmax = matrix; ymax = matrix; zmax = 1; radius = 9; % Spulenradius in mm genauigkeit = 5; % 1 = 6.28 Elemente pro Kreis; 2 = 12.56 Elemente pro Kreis; 4 bis 5 scheint gut zu sein windungen = 10; leiterelemente = ceil(2 * 3.14152 * genauigkeit * windungen) % 2 * Pi * genauigkeit f眉r eine Umrundung leiter = repmat(leiterelemente+1,3); windungen = leiterelemente / genauigkeit / 2 / 3.1415927; spulenlaenge = 20; % Spulenlaenge in mm steigung = spulenlaenge / windungen for i = 1:leiterelemente+1; leiter(i,1) = i * steigung / (genauigkeit * 2 * 3.1415927) + matrix/2 - spulenlaenge/2; % x-Ausrichtung leiter(i,2) = radius * cos(i/genauigkeit) + matrix/2; % y-Ausrichtung leiter(i,3) = radius * sin(i/genauigkeit); % z-Ausrichtung end for x = 1:xmax zaehler = zaehler + 1; % f眉r Zeitmessung hhh = waitbar(0,num2str(zaehler*100/matrix)); % Wartebalken waitbar(zaehler/matrix) % Wartebalken for y = 1:ymax % wenn streamslice nicht genutzt wird, nur einen y-Wert berechnen for z = 1:zmax for i = 1:leiterelemente dl(1) = leiter(i+1,1)-leiter(i,1); dl(2) = leiter(i+1,2)-leiter(i,2); dl(3) = leiter(i+1,3)-leiter(i,3); vecs = [(leiter(i,1)+leiter(i+1,1))/2, ... (leiter(i,2)+leiter(i+1,2))/2, ... (leiter(i,3)+leiter(i+1,3))/2]; vecr = [xyz]; vecrminusvecs = vecr - vecs; einheitsvecr = vecrminusvecs./norm(vecrminusvecs); % ok r = sqrt(vecrminusvecs(1).^2 + vecrminusvecs(2).^2 + vecrminusvecs(3).^2); % ok vektorprodukt = [dl(2).*einheitsvecr(3) - dl(3).*einheitsvecr(2), ... dl(3).*einheitsvecr(1) - dl(1).*einheitsvecr(3), ... dl(1).*einheitsvecr(2) - dl(2).*einheitsvecr(1)]; dB = konstante * I * vektorprodukt / (r.^2); dB = dB / skalierungsfaktor; % nur hier wird der Wert ver盲ndert bzw. skaliert B(x,y,z,1) = B(x,y,z,1) + dB(1); B(x,y,z,2) = B(x,y,z,2) + dB(2); B(x,y,z,3) = B(x,y,z,3) + dB(3); B(x,y,z,4) = B(x,y,z,4) + sqrt(dB(1).^2 + dB(2).^2 + dB(3).^2); end; end; end; close(hhh) end; toc n = 1:leiterelemente; Lx = leiter(n,1); Ly = leiter(n,2); Lz = leiter(n,3); %subplot(2,1,2), line(Lx,Ly,Lz,'Color','k','LineWidth',2); hold on view(15,30); % view(0,0) = Blickwinkel, 2D-Perspektive grid on % Gitter anzeigen xlim([0 matrix]) ylim([0 matrix]) zlim([0 5]) xlabel('x-Achse'); ylabel('y-Achse'); zlabel('z-Achse'); daspect([1 1 1]) [X,Y]=meshgrid(1:matrix); U=(B(1:matrix,1:matrix,z,1))'; V=(B(1:matrix,1:matrix,z,2))'; streamslice(X,Y,U,V) % quiver, streamslice

回答:

不一定是速度优化,而是一些注意事项:


使用pi ,而不是3.14159或3.1415927
您通常可以向量化循环:
非向量化:

for i = 1:leiterelemente+1; leiter(i,1) = i * steigung / (genauigkeit * 2 * 3.1415297) + matrix/2 - spulenlaenge/2; % x-Ausrichtung leiter(i,2) = radius * cos(i/genauigkeit) + matrix/2; % y-Ausrichtung leiter(i,3) = radius * sin(i/genauigkeit); % z-Ausrichtung end 向量化:

ii = 1:leiterelemente+1; leiter(ii,1) = ii * steigung / (genauigkeit * 2 * 3.1415297) + matrix/2 - spulenlaenge/2; % x-Ausrichtung leiter(ii,2) = radius * cos(ii/genauigkeit) + matrix/2; % y-Ausrichtung leiter(ii,3) = radius * sin(ii/genauigkeit); % z-Ausrichtung 大多数matlab函数都将向量/矩阵作为参数,包括cos() , sin() , exp() , log()等。


对于少量元素(例如