Điện tích chuyển động trong từ trường đều

Kết quả mô phỏng chuyển động của điện tích trong từ trường đều

Dòng electron chuyển động theo đường tròn khép kín

Code chương trình Matlab

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91

function charged_particle_UNIFORM_mag_field
% Author: Tran Hai Cat
% Lecturer in Physics, HCM University of Technology and Education
% - Dai hoc Su pham Ky thuat Tp. Ho Chi Minh
% Created: 2019.08.16
clc
clear variables
close all

global m q Bx By Bz

%% CONST
q_element = 1.6e-19;

%% PARAMETR
B0 = 840e-6; % [Tesla]
Bx = B0;
By = 0;
Bz = 0;
m = 9.1095e-31;
q = -q_element;
v0 = 1.11e7; % [m/s]

% Initial position:
x0 = 0;
y0 = 0;
z0 = 0;
v0x = v0/10;
v0y = 0;
v0z = v0;

dt = 7e-11; % [second]
tmax = 5e3*dt;
xmin = 0; xmax = 50e-2;
ymin = -16e-2; ymax = 1e-2;
zmin = -10e-2; zmax = 10e-2;

%% CALCULATION
Nx = 20;
Ny = 10;
Nz = 10;
x_field = linspace(xmin,xmax,Nx);
y_field = linspace(ymin,ymax,Ny);
z_field = linspace(zmin,zmax,Nz);
[X_field,Y_field,Z_field] = meshgrid(x_field,y_field,z_field);
Bx_matrix = B0*ones(size(X_field));
By_matrix = zeros(size(Y_field));
Bz_matrix = zeros(size(Z_field));

[T,Y] = ode45(@lorentz,linspace(0,tmax,500),[x0 y0 z0 v0x v0y v0z]);
x = Y(:,1);
y = Y(:,2);
z = Y(:,3);

%% DISPLAY RESULT
figure('name','CALCULATION PROCESS',...
'color','k','numbertitle','off');
set(gcf,'Units','normalized');
set(gcf,'Position',[0 0.05 1 0.85]);
set(gca,'color','k','xcolor','w','ycolor','w','zcolor','w')
axis equal
axis([xmin xmax ymin ymax zmin zmax]);
rotate3d on

xlabel('x, [m]','fontsize',14);
ylabel('y, [m]','fontsize',14);
zlabel('z, [m]','fontsize',14);
view(30,30);
hold on
quiver3(X_field,Y_field,Z_field,Bx_matrix,By_matrix,Bz_matrix,...
'AutoScaleFactor',1,'color','g');
plot1 = plot3(x(1),y(1),z(1),'ro','markersize',5,'markerfacecolor','r');
plot2 = animatedline('linewidth',2,'color','b');

for i = 1:length(T)
set(plot1,'xdata',x(i),'ydata',y(i),'zdata',z(i));
addpoints(plot2,x(i),y(i),z(i));
drawnow limitrate
pause(0.01)
end

function dy = lorentz(t,y)
%% ODE
global m q Bx By Bz
dy = zeros(6,1);
dy(1) = y(4);
dy(2) = y(5);
dy(3) = y(6);
dy(4) = q*(y(5)*Bz-y(6)*By)/m;
dy(5) = q*(y(6)*Bx-y(4)*Bz)/m;
dy(6) = q*(y(4)*By-y(5)*Bx)/m;