How to remove the straight line coming at the x-axis at x = -0.11532 in my code
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
vetri veeran
am 21 Apr. 2017
Kommentiert: vetri veeran
am 22 Apr. 2017
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.
0 Kommentare
Akzeptierte Antwort
Walter Roberson
am 21 Apr. 2017
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 Kommentare
Walter Roberson
am 21 Apr. 2017
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?
Weitere Antworten (0)
Siehe auch
Kategorien
Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!