How to remove the straight line coming at the x-axis at x = -0.11532 in my code

3 次查看（过去 30 天）
vetri veeran
vetri veeran 2017-4-21
评论： vetri veeran 2017-4-22
Hi all,
In my code, I want to remove the straight line coming at the x-axis at x = -0.11532.
I attached my code below,
B = 1e-4;
sigma_on = 0.45;
x_on = 0.06;
sigma_p = 4e-5;
A = 1e-10;
sigma_off = 0.013;
x_off = 0.4;
G_m = 0.025;
a = 7.2e-6;
b = 4.7;
beta = 500;
rho = 1e-3;
v_m = 1;
k=50;
t = -1:0.001:1;
for x = 1:length(t)
G(x) = G_m*t(x)+a*exp(b*sqrt(v_m/(1+exp(-v_m/rho))-v_m/(1+exp(v_m/rho))))*(1-t(x));
f1(x) = A*sinh(v_m/sigma_off)*exp(-(x_off^2/t(x)^2))*exp(1/(1+beta*G(x)*v_m^2))*(1/(1+exp(k*v_m)));
f2(x) = B*sinh(v_m/sigma_on)*exp(-(t(x)^2/x_on^2))*exp(G(x)*v_m^2/sigma_p)*(1/(1+exp(-k*v_m)));
f(x) = f1(x) + f2(x);
end
semilogy(t,f);
xlabel('x','fontsize', 20);
ylabel('y','fontsize', 20);
Could someone help me.

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

采纳的回答

Walter Roberson
Walter Roberson 2017-4-21
You cannot do much about it. Your expression for f1 includes
exp(1/(1+beta*G(x)*v_m^2))
That value can be arbitrarily high if (1+beta*G(x)*v_m^2) approaches 0; with your beta = 500 and v_m = 1, that is the condition that G(x) approximately equal -1/500 .
You can go to the line above,
G(x) = G_m*t(x)+a*exp(b*sqrt(v_m/(1+exp(-v_m/rho))-v_m/(1+exp(v_m/rho))))*(1-t(x));
and construct
-1/500 == G_m*T+a*exp(b*sqrt(v_m/(1+exp(-v_m/rho))-v_m/(1+exp(v_m/rho))))*(1-T);
and solve for T. You get a result of T about -0.115316250006499 . You are incrementing by 0.001 so when your T becomes -0.115 you have a near singularity .
If your equations are correct then the only way to avoid the singularity is not to calculate near it.
  3 个评论
显示 1更早的评论
Walter Roberson
Walter Roberson 2017-4-21
Considering that you have 1/(1+beta*G(x)*v_m^2) and the denominator switches between positive and negative, what kind of output were you hoping for?

请先登录，再进行评论。

更多回答（0 个）

类别

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

标签

Community Treasure Hunt

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

Start Hunting!

Translated by