double integration for below defined problem

To utilize in , we substitute the expression for into the integral for , yielding an expression involving and . Here's the updated integral expression:

Substituting the expression for gives:

Now, you can address this integral either analytically or numerically, depending on the specifics of and the strategies employed during the evaluation of 's integral. If is simple enough, you might be able to solve it analytically.

below my code:

% Constants and parameters
m_star = 0.5 * 9.1e-31;   % Effective mass of electrons (kg)
h_bar = 1.054e-34;        % Reduced Planck's constant (Joule seconds)
V_0 = 0.1;                % Barrier height (eV)
V_b = 0.5;                % Classical turning point (eV)
L_w = 2e-9;               % Width of the barrier (m)
E_1 = 0.5;                % Lower limit of the energy integral (eV)
F = 1e3;                  % Electric field (V/m)
k = 8.6173e-5;            % Boltzmann constant (in eV/K)
T = 300;                  % Kelvin

% Define the function for T(E, F)
TEF = @(E, F) exp(-2 * integral(@(z) sqrt(2 * m_star / h_bar^2 * (V_0 - E - F * z)), 0, (V_b - 0.5 * F * L_w - E) / F));

% Define the Fermi-Dirac distribution function with E_F = E_1 + k * T
f = @(E) 1 ./ (1 + exp((E_1 + k * T - E) / T));

% Define the integrand
integrand = @(E, F) f(E) .* TEF(E, F);

% Numerically evaluate TEF
TEF_value = TEF(E_1, F)

TEF_value = 0.9013 - 0.4332i

% Numerically evaluate the integral
nF = integral(integrand, E_1, Inf) * m_star / (pi * h_bar^2 * L_w);

Not enough input arguments.

Error in solution>@(E,F)f(E).*TEF(E,F) (line 19)
integrand = @(E, F) f(E) .* TEF(E, F);

Error in integralCalc/iterateScalarValued (line 314)
fx = FUN(t);

Error in integralCalc/vadapt (line 132)
[q,errbnd] = iterateScalarValued(u,tinterval,pathlen);

Error in integralCalc (line 83)
[q,errbnd] = vadapt(@AToInfInvTransform,interval);

Error in integral (line 87)
Q = integralCalc(fun,a,b,opstruct);

% Display the results
disp(['The value of TEF at E_1 is: ', num2str(TEF_value)]);
disp(['The value of n(F) is: ', num2str(nF)]);

it show error

Answer 1

[Torsten]

% Constants and parameters
m_star = 0.5 * 9.1e-31;   % Effective mass of electrons (kg)
h_bar = 1.054e-34;        % Reduced Planck's constant (Joule seconds)
V_0 = 0.1;                % Barrier height (eV)
V_b = 0.5;                % Classical turning point (eV)
L_w = 2e-9;               % Width of the barrier (m)
E_1 = 0.5;                % Lower limit of the energy integral (eV)
F = 1e3;                  % Electric field (V/m)
k = 8.6173e-5;            % Boltzmann constant (in eV/K)
T = 300;                  % Kelvin
% Define the function for T(E, F)
TEF = @(E, F) exp(-2 * integral(@(z) sqrt(2 * m_star / h_bar^2 * (V_0 - E - F * z)), 0, (V_b - 0.5 * F * L_w - E) / F));
% Define the Fermi-Dirac distribution function with E_F = E_1 + k * T
f = @(E) 1 ./ (1 + exp((E_1 + k * T - E) / T));
% Define the integrand
integrand = @(E, F) f(E) .* TEF(E, F);
% Numerically evaluate TEF
TEF_value = TEF(E_1, F)
% Numerically evaluate the integral
nF = integral(@(E)integrand(E,F), E_1, Inf, 'ArrayValued',1) * m_star / (pi * h_bar^2 * L_w);
% Display the results
disp(['The value of TEF at E_1 is: ', num2str(TEF_value)]);
disp(['The value of n(F) is: ', num2str(nF)]);