Nguyên tử hidro trong trường hợp đối xứng cầu

Hidro là loại nguyên tử có cấu trúc đơn giản nhất trong tất cả các nguyên tố: chỉ một electron bao quanh hạt nhân cấu thành từ một proton. Hạt nhân proton này tạo ra xung quanh nó một điện trường, có xu hướng hút electron vào gần nó. Electron lúc này không còn chuyển động tự do, mà rơi vào hố thế của trường tĩnh điện Coulomb:

\[U(r)=-\frac{k_ee^2}{r},\qquad k_e=\frac{1}{4\pi\varepsilon_0},\label{eq:1}\tag{1}\]

với \(r\) –  khoảng cách đến hạt nhân nguyên tử, \(e\) – điện tích nguyên tố, \(\varepsilon_0\) – hằng số điện. Nếu xem electron như một hạt mang điện tích thông thường, ta sẽ có mô hình một điện tích điểm chuyển động tròn với lực hút Coulomb đóng vai trò lực hướng tâm. Chuyển động có gia tốc ấy sẽ bức xạ sóng điện từ, khiến electron dần mất năng lượng và rơi vào hạt nhân.

Thuyết nguyên tử Bohr đã tránh được hạn chế đó, nhưng bằng cách “ép” các electron tuân theo tiên đề về quỹ đạo bền. Điều đó không khác gì nói rằng: hạt không bị rơi xuống hạt nhân vì Bohr cấm nó rơi xuống!

Lý thuyết sóng-hạt của Luis de Broglie giúp ta tránh được những khó khăn ấy. Bản thân electron không còn là chất điểm, mà trở thành một sóng. Nó lấp đầy lòng nguyên tử, nên không còn rơi vào hạt nhân nữa.

Phương trình Schrodinger cho nguyên tử hidro đối xứng

Phương trình Schrodinger viết theo dạng ba chiều:

\[-\frac{\hbar^2}{2m}\left(\frac{\partial^2}{\partial x^2}+\frac{\partial^2}{\partial y^2}+\frac{\partial^2}{\partial z^2}\right)\Psi(x,y,z)+U(r)\Psi(x,y,z)=E\Psi(x,y,z).\label{eq:2}\tag{2}\]

Sẽ thuận tiện hơn nếu khai triển đạo hàm trong hệ toạ độ cầu:

\[\frac{\partial^2}{\partial x^2}+\frac{\partial^2}{\partial y^2}+\frac{\partial^2}{\partial z^2}=\frac{1}{r^2}\frac{\partial}{\partial r}\left(r^2\frac{\partial}{\partial r}\right)+\frac{1}{r^2\sin\vartheta}\frac{\partial}{\partial\vartheta}\left(\sin\vartheta\frac{\partial}{\partial\vartheta}\right)+\frac{1}{r^2\sin^2\vartheta}\frac{\partial^2}{\partial\varphi^2}.\]

Ta sẽ đi tìm những hàm sóng \(\Psi(r)\) chỉ phụ thuộc vào khoảng cách \(r\) mà không phụ thuộc vào góc \(\vartheta\) cũng như góc cực \(\varphi\):

\[\Psi(x,y,z)\equiv R(r).\]

Lúc này đạo hàm chỉ giữ lại các thành phần theo \(r\), vì đạo hàm theo \(\vartheta\) và \(\varphi\) bằng không. Phương trình \eqref{eq:2} viết thành:

\[-\frac{\hbar^2}{2m}\frac{1}{r^2}\frac{\partial}{\partial r}\left(r^2\frac{\partial}{\partial r}\right)R(r)+U(r)R(r)=ER(r).\label{eq:3}\tag{3}\]

Hàm thế \(U(r)\) có dạng một hố hyperbol như phương trình \eqref{eq:1}.

Giải bài toán nguyên tử hidro đối xứng

Thế năng của trường Coulomb do hạt nhân proton tạo ra tiến đến sâu vô cùng tại \(r=0\). Biểu hiện của hàm sóng \(\psi(r)\) tại vị trí đặc biệt này còn chưa rõ. Vì vậy ta sẽ đi tìm \(R(r)\) dưới dạng:

\[R(r)=\frac{u(r)}{r}.\]

Qua hệ thức này, ta biết chắc điều kiện biên của hàm \(u(r)\) tại \(r=0\) phải bằng \(0\). Bởi nếu \(u(0)>0\), \(\psi(r=0)\) sẽ tiến đến \(\infty\). Nhờ biến đổi toán học:

\[\frac{1}{r^2}\frac{\partial}{\partial r}\left(r^2\frac{\partial}{\partial r}\right)\frac{u(r)}{r}=\frac{1}{r}\frac{\partial^2u(r)}{\partial r^2},\]

phương trình Schrodinger \eqref{eq:3} có thể viết lại theo hàm \(u(r)\):

\[-\frac{\hbar^2}{2m}\frac{\partial^2u(r)}{\partial r^2}+U(r)u(r)=Eu(r).\label{eq:4}\tag{4}\]

Để giải phương trình \eqref{eq:4}, cần viết lại theo thể thực hành:

\[u”(r)=-\frac{2m}{\hbar^2}\left[E-U(r)\right]u(r).\label{eq:5}\tag{5}\]

Phương trình \eqref{eq:5} trình bày trong Matlab theo code:

1
2
3
4
5
6
function dy = Schrodinger(r,y)
%% Schrodinger Equation
global E k
dy = zeros(2,1);
dy(1) = y(2);
dy(2) = -k<em>(E-U_hydro(r))</em>y(1);

Với hàm thế năng hyperbol:

1
2
3
4
5
6
7
8
function y = U_hydro(r)
%% Potential Energy
global epsilon0 q_element
lenr = length(r);
y = zeros(1,lenr);
for i = 1:lenr
y(i) = -1/4/pi/epsilon0<em>q_element</em>q_element/r(i);
end

Phạm vi giải phương trình của chúng ta trải dài từ tâm nguyên tử \(r=0\) đến “rìa” nguyên tử, một khoảng cách lựa chọn khéo léo để có thể xem là đủ xa vô cực. Nghiệm của phương trình (4) hoặc (5) phải tuân theo điều kiện biên:

\[\begin{cases}u(0)&=0,\\ u(+\infty)&=0.\end{cases}\]

Áp dụng phương pháp “bắn tên”, ta sẽ chọn một mức năng lượng \(E\) nào đó rồi bắn thử từ biên bên này và nghiệm xem mũi tên liệu có trúng đích, tức điều kiện biên của phía biên còn lại hay không. Để thuận tiện, ta sẽ bắn mũi tên không phải từ \(r=0\), mà từ ngoài rìa bắn vào với điều kiện ban đầu:

\[u(+\infty)\approx 0,\qquad u'(+\infty)\approx 0.\]

Code Matlab tương ứng có thể viết:

1
2
3
u0 = 1e-10;
u10 = 0;
[r,U] = ode45(@Schrodinger,linspace(rmax_cal,rmin_cal,1000),[u0 u10]);

Ở đây \(\mathrm{rmax\_cal}\) tương ứng với một điểm nằm cách đủ xa hạt nhân nguyên tử. \(\mathrm{rmin\_cal}\) ta chọn sát hạt nhân, nhưng không được trùng với vị trí \(r=0\) của hạt nhân, bởi lực hút sẽ tiến về vô cùng, tâm hố thế sâu vô cực. Trên thực tế, hạt nhân cũng có kích thước, tuy rất nhỏ có thể xem gần như chất điểm, nhưng cũng vừa đủ để tránh lực điện từ vô hạn. Bên trong phạm vi rất nhỏ của hạt nhân có hố thế đặc biệt khác và diễn biến khác, ta không xét ở đây. Có thể chọn \(\mathrm{rmin\_cal}\) tầm một phần nghìn Angstrom. Nhớ rằng bán kính Bohr của hidro bằng \(0.529\,\mathrm{A}\).

Mục đích của chúng ta là phải bắn sao cho hàm \(u(r)\) rơi trúng gốc toạ độ: \(u(0)\approx 0\). Kết quả tính toán cho thấy, không phải mọi giá trị của năng lượng \(E\) đều thoả điều kiện biên này, mà chỉ tại những mức \(E_n\) rời rạc nhất định. Hoá ra những mức rời rạc này đúng theo công thức theo thuyết nguyên tử Bohr đã đề cập trước đó:

\[E_n=-\frac{me^4}{32\pi^2\epsilon_0^2\hbar^2}\cdot\frac{1}{n^2},\qquad n=1,2,3,\ldots\]

Theo đó ba mức năng lượng đầu tiên: \(E_1=-13.566\,\mathrm{eV}\), \(E_2=-3.391\,\mathrm{eV}\), \(E_1=-1.507\,\mathrm{eV}\). Những hình dưới đây thể hiện ba nghiệm đầu tiên ấy.

Hình 1: Nghiệm bậc 1

Ta thấy rằng nghiệm \(u(r)\) của phương trình \eqref{eq:5}, biểu diễn theo đường màu vàng có đi qua gốc toạ độ, đảm bảo đúng điều kiện biên. Vùng tô màu cam thể hiện mật độ của electron theo khoảng cách tính từ tâm nguyên tử:

\[P(r)=4\pi r^2R(r)^2.\]

Lý thuyết bán cổ điển Bohr cho rằng electron như một chất điểm chuyển động trên các quỹ đạo tròn có bán kính \(1\,r_{Bohr}\), \(4\,r_{Bohr}\), \(9\,r_{Bohr}\ldots\). Trong khi đó theo kết quả của cơ học lượng tử, electron hình thành đám mây bao quanh hạt nhân. Nhưng quan sát kĩ cho thấy, khi nguyên tử nằm ở trạng thái dừng có mức năng lượng thấp nhất, ta thấy mật độ \(P(r)\) của electron quanh hạt nhân proton đậm đặc nhất tại vị trí đúng bằng \(1\) bán kính Bohr! Một kết quả trùng hợp khá thú vị. Lúc này đám mây electron hình thành một vỏ cầu nhỏ bao quanh hạt nhân.

Hình 2: Nghiệm bậc 2

Khi chuyển sang các mức năng lượng cao hơn, phần không gian sóng có thể dao động trở nên rộng hơn, nên đám mây electron bị tản ra xa hơn. Với nghiệm bậc 2, electron hình thành hai lớp vỏ cầu, với lớp vỏ ngoài dày hơn và tập trung mật độ electron cao hơn. Có thể hình dung rằng electron đã “chuyển lên” mức quỹ đạo cao hơn, cách xa hạt nhân hơn. Có điều ta thấy, bán kính quỹ đạo này không còn khớp với lý thuyết Bohr nữa.

Hình 3: Nghiệm bậc 3

Tình hình tương tự khi năng lượng tiếp tục dâng cao lên mức 3, tương ứng với nghiệm bậc 3. Một cách hình tượng, có thể xem electron như nước rót vào trong cốc, một chiếc cốc hyperbol tròn xoay. Nước càng dâng cao, năng lượng càng đầy, mặt nước càng rộng, bán kính nguyên tử càng lớn.

Video dưới đưa ra sự so sánh thú vị về sóng dừng hình thành trong nguyên tử hidro với sóng dừng cơ học hình thành trên một đĩa tròn rung động.

[Video]

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
92
93
94
95
96
97
98
99
100
function R_hidro_symmetric_2
% Created by Tran Hai Cat
% 2019.04.27
clc;
clear variables
close all

global E k epsilon0 q_element

%% CONSTANTS
h = 1.054e-34;
m = 9.1095e-31;
q_element = 1.6e-19;
Angstrom = 1e-10;
k = 2*m/h/h;
epsilon0 = 8.854e-12;
r_bohr = 0.529e-10;

%% INPUT PARAMETRS FOR USER
n = 2;
E = -m*q_element^4/32/pi^2/epsilon0^2/h^2/n^2;

% E = -1.5*q_element;

%% DATA PROCESSING
Emin = -20*q_element;
Emax = 20*q_element;
rmin = -5*r_bohr*n^2;
rmax = 5*r_bohr*n^2;
rmin_cal = 0.001*r_bohr;
rmax_cal = 2*rmax;

%% CALCULATION
u0 = 1e-10;
u10 = 0;

[r,U] = ode45(@Schrodinger,linspace(rmax_cal,rmin_cal,1000),[u0 u10]);
u = U(:,1);

r = [-r;flipud(r)];
u = [-u;flipud(u)];
R = u./r;
P = 4*pi*r.^2.*R.^2;
norm = abs(trapz(r,P));
R = R/sqrt(norm);
u = R.*r;
P = P/norm;

u_draw = u*3e-4;
R_draw = R*2e-14;
P_draw = P*3e-9;

%% FIGURE
r_plot = [linspace(-rmax,-0.5*Angstrom,1000),...
linspace(0.5*Angstrom,rmax,1000)];
U_plot = U_hydro(abs(r_plot));

figure('name','R-function hydrogen atom',...
'color','black','numbertitle','off')
hold on
line([0,0],[Emin/q_element,Emax/q_element],'color','w');
patch([r/r_bohr,fliplr(r/r_bohr)],...
[E/q_element+P_draw,fliplr(E/q_element+zeros(size(P_draw)))],...
[0.6 0.3 0.0]);

plot(r/r_bohr,E/q_element+u_draw,'linewidth',1.2,'color','y');
plot(r/r_bohr,E/q_element+R_draw,'linewidth',1.2,'color','m');
plot(r_plot/r_bohr,U_plot/q_element,'linewidth',2,'color','w');
line([rmin/r_bohr,rmax/r_bohr],[E/q_element,E/q_element],...
'linewidth',1.5,'color','b');

title(sprintf('E = %0.3f [eV]',E/q_element),'color','w');
axis([rmin/r_bohr rmax/r_bohr Emin/q_element Emax/q_element]);
set(gca,'color','k','xcolor','w','ycolor','w')
xlabel('r/r_{bohr}');
ylabel('E, U [eV]');
legend(...
'\color{white} r=0',...
'\color{white} 4\pi r^2R^2(r)',...
'\color{white} u(r)',...
'\color{white} R(r)',...
'\color{white} U(r)',...
'\color{white} E'...
);

function dy = Schrodinger(r,y)
%% Schrodinger Equation
global E k
dy = zeros(2,1);
dy(1) = y(2);
dy(2) = -k*(E-U_hydro(r))*y(1);

function y = U_hydro(r)
%% Potential Energy
global epsilon0 q_element
lenr = length(r);
y = zeros(1,lenr);
for i = 1:lenr
y(i) = -1/4/pi/epsilon0*q_element*q_element/r(i);
end