Solution of the equation

We have such a function:

We set the array d = [10:10:10000], d(i) is substituted in and equated to the constant C_2, i.e. we get the equation F(d(i), t)=C_2 (C_2 = 6.81), and we need to find this 't' for the corresponding d(i), that is, as a result, we should have a graph d(t). I think that most likely it will be a discontinuous fast-oscillating function.

My code:

%% initial conditions
global d k0 h_bar ksi m E;
Ef = 2.77*10^3; 
Kb = physconst('boltzmann'); % 1.38*10^(-23)
T = 0.12:0.24:6.4;
m = 9.1093837*10^(-31);
Tc = 1.2;
%t = T./Tc;
t = 0.1:0.1:2;
nD = floor(375./(2.*pi.*t.*1.2) - 0.5);
D = 10^(-8); % толщина пленки
ksi = 10^(-9);
%d = D/ksi;
d = 1000;
E = Ef/(pi*Kb*Tc);
h_bar = (1.0545726*10^(-34));
k0 = (ksi/h_bar)*sqrt(2.*m.*pi.*Kb.*Tc);
C_2 = 6.81;
%% calculation
F = f_calc(t,nD);
plot(t,F)
grid on

function F = f_calc(t,nD) 
global d k0 h_bar ksi m;
F = zeros(1,numel(t));
for i = 1:numel(t) 
    for k = 0:nD(i)
        F(i)  = F(i) + 1/(2*k+1).*(k0.*real(((f_p1(k,t(i))-f_p2(k,t(i)))./2))+(f_arg_2(k,t(i))-f_arg_1(k,t(i)))./d); 
    end
end
F = -F;
%F = -(1/d).*F;
%F = F - C_2;
end
function p1 = f_p1(n,t)
p1 = ((1+1i)./sqrt(2)).*sqrt(t.*(2.*n+1));
end
function p2 = f_p2(n,t)
global E;
p2 = sqrt(3601+1i.*t.*(2.*n+1));
end 
function n_lg = f_lg(n,t)
global d k0;
arg_of_lg = (1+exp(-1i*d*k0.*f_p1(n,t)))/(1+exp(-1i*d*k0.*f_p2(n,t)));
n_lg = log(abs(arg_of_lg));
end
function arg_1 = f_arg_1(n,t)
global d k0;
arg_1 = angle(1+exp(-1i*d*k0.*f_p1(n,t)));
end
function arg_2 = f_arg_2(n,t)
global d k0;
arg_2 = angle(1+exp(-1i*d*k0.*f_p2(n,t)));
end

[Walter Roberson]

What are you asking of us?

[Torsten]

Since t is unknown, which nD(t) should be used ? Should it be changed according to

nD = floor(375./(2.*pi.*t.*1.2) - 0.5);

during the iteration process ?

[Dmitry]

Torsten, yes, you're right.

[Dmitry]

Walter Roberson, How can I convert my code to solve this problem

[Walter Roberson]

What problem? You showed an equation, you showed code. Are you getting an error message? If not then at what point in the code do you start seeing unexpected results? As outsiders how would we know if the code was producing correct answers or not?

[Dmitry]

Walter Roberson, I don't understand how to take into account in this problem that the upper limit of the sum depends on t

Answer 1

[Torsten]

This code doesn't work since there does not seem to be a solution for some or all of the d values you supplied.

I suggest you fix a value for d first and plot f_calc(t) for a number of t values to see if the curve passes the t-axis somewhere.

%% initial conditions
global d k0 h_bar ksi m E C_2
Ef = 2.77*10^3; 
Kb = physconst('boltzmann'); % 1.38*10^(-23)
T = 0.12:0.24:6.4;
m = 9.1093837*10^(-31);
Tc = 1.2;
D = 10^(-8); % толщина пленки
ksi = 10^(-9);
dd = [10:10:10000];
E = Ef/(pi*Kb*Tc);
h_bar = (1.0545726*10^(-34));
k0 = (ksi/h_bar)*sqrt(2.*m.*pi.*Kb.*Tc);
C_2 = 6.81;
t0 = 1.0;
for i = 1:numel(dd)
  d = dd(i);
  t(i) = fsolve(@f_calc,t0);
  t0 = t(i);
end
plot(t,F)

function F = f_calc(t) 
global d k0 h_bar ksi m C_2
nD = floor(375/(2*pi*t*1.2) - 0.5);
F = 0;
for k = 0:nD
    F  = F + 1/(2*k+1).*(k0.*real(((f_p1(k,t)-f_p2(k,t))./2))+(f_arg_2(k,t)-f_arg_1(k,t))./d); 
end
F = F - C_2;
end
function p1 = f_p1(n,t)
p1 = ((1+1i)./sqrt(2)).*sqrt(t.*(2.*n+1));
end
function p2 = f_p2(n,t)
global E
p2 = sqrt(3601+1i.*t.*(2.*n+1));
end 
function n_lg = f_lg(n,t)
global d k0
arg_of_lg = (1+exp(-1i*d*k0.*f_p1(n,t)))/(1+exp(-1i*d*k0.*f_p2(n,t)));
n_lg = log(abs(arg_of_lg));
end
function arg_1 = f_arg_1(n,t)
global d k0
arg_1 = angle(1+exp(-1i*d*k0.*f_p1(n,t)));
end
function arg_2 = f_arg_2(n,t)
global d k0
arg_2 = angle(1+exp(-1i*d*k0.*f_p2(n,t)));
end

[Dmitry]

I think maybe I can try to solve the problem by setting the grid 't' 'd' and using the 'contour' function to build a graph?

[Dmitry]

Do you think maybe such a solution will work, but for some reason I get an error?

%% initial conditions

global k0 h_bar ksi m E C_2

Ef = 2.77*10^3;

Kb = physconst('boltzmann'); % 1.38*10^(-23)

m = 9.1093837*10^(-31);

Tc = 1.2;

ksi = 10^(-9);

E = Ef/(pi*Kb*Tc);

h_bar = (1.0545726*10^(-34));

k0 = (ksi/h_bar)*sqrt(2.*m.*pi.*Kb.*Tc);

C_2 = 6.81;

t = linspace(0.1, 2);

d = linspace(10, 1000);

[TT,DD] = meshgrid(t,d);

contour(TT, DD, @f_calc);

function F = f_calc(t, d)

global k0 h_bar ksi m C_2

nD = floor(375/(2*pi.*t*1.2) - 0.5);

F = 0;

for k = 0:nD

F = F + 1/(2*k+1).*(k0.*real(((f_p1(k,t)-f_p2(k,t))./2))+(f_arg_2(k,t, d)-f_arg_1(k,t, d))./d);

end

F = F - C_2;

end

function p1 = f_p1(n,t)

p1 = ((1+1i)./sqrt(2)).*sqrt(t.*(2.*n+1));

end

function p2 = f_p2(n,t)

global E

p2 = sqrt(3601+1i.*t.*(2.*n+1));

end

function arg_1 = f_arg_1(n,t,d)

global k0

arg_1 = angle(1+exp(-1i.*d*k0.*f_p1(n,t)));

end

function arg_2 = f_arg_2(n,t, d)

global k0

arg_2 = angle(1+exp(-1i.*d*k0.*f_p2(n,t)));

end

[Torsten]

As you can see, there seems to be no t value in the range you specified such that F(t,d(i)) = 6.81 for all d values.

%% initial conditions
global k0 h_bar ksi m E C_2
Ef = 2.77*10^3; 
Kb = physconst('boltzmann'); % 1.38*10^(-23)
m = 9.1093837*10^(-31);
Tc = 1.2;
ksi = 10^(-9);
E = Ef/(pi*Kb*Tc);
h_bar = (1.0545726*10^(-34));
k0 = (ksi/h_bar)*sqrt(2.*m.*pi.*Kb.*Tc);
C_2 = 6.81;
t = linspace(0.1, 2);
d = linspace(10, 1000);
%[TT,DD] = meshgrid(t,d);
%contour(TT, DD, @f_calc);
for i=1:numel(d)
  for j = 1:numel(t)
    F(i,j) = f_calc(t(j),d(i));
  end
end
plot(t,F)

function F = f_calc(t, d) 
global k0 h_bar ksi m C_2
nD = floor(375/(2*pi.*t*1.2) - 0.5);
F = 0;
for k = 0:nD
    F  = F + 1/(2*k+1).*(k0.*real(((f_p1(k,t)-f_p2(k,t))./2))+(f_arg_2(k,t, d)-f_arg_1(k,t, d))./d); 
end
F = F - C_2;
end
function p1 = f_p1(n,t)
p1 = ((1+1i)./sqrt(2)).*sqrt(t.*(2.*n+1));
end
function p2 = f_p2(n,t)
global E
p2 = sqrt(3601+1i.*t.*(2.*n+1));
end 
function arg_1 = f_arg_1(n,t,d)
global k0
arg_1 = angle(1+exp(-1i.*d*k0.*f_p1(n,t)));
end
function arg_2 = f_arg_2(n,t, d)
global k0
arg_2 = angle(1+exp(-1i.*d*k0.*f_p2(n,t)));
end

[Dmitry]

But you've got a graph

[Dmitry]

We also need to plot t(d)

Maybe build a plane F(t, d) and see if it intersects with the plane C_2?

[Torsten]

I've got a graph that never crosses F - C_2 = 0 (which would mean that there exists a t value or which F(t,d(i)) - C_2 = 0).

Simply spoken: The graphs must somewhere cross y=0 for that your problem has a solution.

It's the same as if you graph f(x) = x^2 - 2 to see where the zeros of x^2 - 2 = 0 are located.

[Dmitry]

Thank you so much!