清空筛选条件
清空筛选条件

Script has no errors, but no plot is given

23 次查看（过去 30 天）
Sergio
Sergio 2024-7-10，15:45
回答： Mathieu NOE 2024-7-10，16:34
Hi , I wonder what is wrong here with the script. Nothinig happens when it should plot something... Can someone help? Thanks
% Constants
hbar = 1.0545718e-34; % Planck's constant
m = 9.10938356e-31; % Electron mass
l = 1; % angular momentum quantum number
E_ion = 1; % Ionization energy, example value
%Functions
function e =exp(r)
% Discretization parameters
r_min = 0.01; % Avoid division by zero
r_max = 10;
N = 100; % Number of grid points
r = linspace(r_min, r_max, N);
dr = r(2) - r(1);
% Initialize the R(r) vector
R = zeros(N, 1);
% Set up the finite difference matrix
A = zeros(N, N);
for i = 3:N-2
A(i, i-2) = -1 / (2 * dr^3);
A(i, i-1) = 2 / (dr^3) - 3 / (r(i) * dr^2);
A(i, i) = -2 / (dr^3) + 3 / (r(i)^2 * dr) + e^(r(i)^2) / (i * hbar^3 / (sqrt(2*m)^3));
A(i, i+1) = -2 / (dr^3) + 3 / (r(i) * dr^2);
A(i, i+2) = 1 / (2 * dr^3);
end
% Apply boundary conditions (example: R(0) = 0 and R(r_max) = 0)
A(1,1) = 1;
A(N,N) = 1;
% Solve the eigenvalue problem
[~, D] = eig(A);
% The eigenvalues are the diagonal elements of D
eigenvalues = diag(D);
% The solution R(r) corresponds to the eigenvector with eigenvalue closest to E_ion
[~, idx] = min(abs(eigenvalues - E_ion));
R = eigvecs(:, idx);
% Plot the solution
plot(r, R);
xlabel('r');
ylabel('R(r)');
title('Radial Wavefunction');
end
  2 个评论
Aquatris
Aquatris 2024-7-10，15:46
编辑：Aquatris 2024-7-10，15:48
You are not calling the function which you called exp() that produces the plot so why do you expect a plot?
Also never use built-in function/variable names as custom function names. exp() already exist in matlab.
Sergio
Sergio 2024-7-10，15:54
编辑：Sergio 2024-7-10，15:58
The problem lies in "A(i, i) = -2 / (dr^3) + 3 / (r(i)^2 * dr) + e^(r(i)^2) / (i * hbar^3 / (sqrt(2*m)^3));" I defined the function exp(r ) since it should be connected to the line here. Should I remove the function calling at the top, and rewrite exp(r(i)^2) instead? Did that, then I got a new error related to eigvecs which is unknown.

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

KSSV
KSSV 2024-7-10，16:26
Ran in:
% Constants
hbar = 1.0545718e-34; % Planck's constant
m = 9.10938356e-31; % Electron mass
l = 1; % angular momentum quantum number
E_ion = 1; % Ionization energy, example value
%Functions
% Discretization parameters
r_min = 0.01; % Avoid division by zero
r_max = 10;
N = 100; % Number of grid points
r = linspace(r_min, r_max, N);
dr = r(2) - r(1);
% Initialize the R(r) vector
R = zeros(N, 1);
% Set up the finite difference matrix
A = zeros(N, N);
for i = 3:N-2
A(i, i-2) = -1 / (2 * dr^3);
A(i, i-1) = 2 / (dr^3) - 3 / (r(i) * dr^2);
A(i, i) = -2 / (dr^3) + 3 / (r(i)^2 * dr) + exp(r(i)^2) / (i * hbar^3 / (sqrt(2*m)^3));
A(i, i+1) = -2 / (dr^3) + 3 / (r(i) * dr^2);
A(i, i+2) = 1 / (2 * dr^3);
end
% Apply boundary conditions (example: R(0) = 0 and R(r_max) = 0)
A(1,1) = 1;
A(N,N) = 1;
% Solve the eigenvalue problem
[eigvecs, D] = eig(A);
% The eigenvalues are the diagonal elements of D
eigenvalues = diag(D);
% The solution R(r) corresponds to the eigenvector with eigenvalue closest to E_ion
[~, idx] = min(abs(eigenvalues - E_ion));
R = eigvecs(:, idx);
% Plot the solution
plot(r, R);
xlabel('r');
ylabel('R(r)');
title('Radial Wavefunction');

更多回答（1 个）

Mathieu NOE
Mathieu NOE 2024-7-10，16:34
Ran in:
again...
let's fix it
I cannot check here if e^(r(i)^2) (probably wrong) must be replaced by exp(r(i)^2) - I let you double check that point , but that is actually what I did
then, you have to change one more line to define eigvecs
[~, D] = eig(A); => [eigvecs, D] = eig(A);
% Constants
hbar = 1.0545718e-34; % Planck's constant
m = 9.10938356e-31; % Electron mass
l = 1; % angular momentum quantum number
E_ion = 1; % Ionization energy, example value
% Discretization parameters
r_min = 0.01; % Avoid division by zero
r_max = 10;
N = 100; % Number of grid points
r = linspace(r_min, r_max, N);
dr = r(2) - r(1);
% Initialize the R(r) vector
R = zeros(N, 1);
% Set up the finite difference matrix
A = zeros(N, N);
for i = 3:N-2
A(i, i-2) = -1 / (2 * dr^3);
A(i, i-1) = 2 / (dr^3) - 3 / (r(i) * dr^2);
A(i, i) = -2 / (dr^3) + 3 / (r(i)^2 * dr) + exp(r(i)^2) / (i * hbar^3 / (sqrt(2*m)^3));
A(i, i+1) = -2 / (dr^3) + 3 / (r(i) * dr^2);
A(i, i+2) = 1 / (2 * dr^3);
end
% Apply boundary conditions (example: R(0) = 0 and R(r_max) = 0)
A(1,1) = 1;
A(N,N) = 1;
% Solve the eigenvalue problem
[eigvecs, D] = eig(A);
% The eigenvalues are the diagonal elements of D
eigenvalues = diag(D);
% The solution R(r) corresponds to the eigenvector with eigenvalue closest to E_ion
[~, idx] = min(abs(eigenvalues - E_ion));
R = eigvecs(:, idx);
% Plot the solution
plot(r, R);
xlabel('r');
ylabel('R(r)');
title('Radial Wavefunction');

类别

在 Help Center 和 File Exchange 中查找有关 Annotations 的更多信息

标签

产品


版本

R2024a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by