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

1 回表示 (過去 30 日間)
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 件のコメント
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?
vetri veeran
vetri veeran 2017 年 4 月 22 日
Thank you very much Walter Robertson

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeMATLAB についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by