MATLAB Answers

0

how to solve transcedental equation in matlab

Asked by alburary daniel on 13 Jun 2018
Latest activity Edited by Walter Roberson
on 27 Jul 2018
I am practicing solving the next transcendental equation in matlab
(a/c)*sqrt((b*m1)^2-(p*c)^2)-atan( sqrt(( (p*c)^2-(b*m2)^2 ) /( (b*m1)^2-(p*c)^2 ) ) )-atan( sqrt(( (p*c)^2-(b*m3)^2 ) /( (b*m1)^2-(p*c)^2 ) ) ) == r*pi
here
a=1x10^-6;
c= 3x10^8;
m1=2.2;
m2=1.5;
m3=1;
I was trying to plot "p" vs "b" where "b" runs from 0 to 3x10^15 and r is a parameter that takes values of 0, 1 and 2. I already tried all day but I cannot find solution, I tried with fzero(fun,xo) without success, can you give any suggestion?

  4 Comments

Show 1 older comment
You ask about plotting p vs b where b runs from 0 to 3E15, but the axis that runs to about 3E15 on the plot is the left hand axis. That suggest that the plot is b vs p instead of p vs p.
It is possible to demonstrate that if a root exists at all, that it is in the range p = b/200000000 to 11*b/1500000000 so you could use those as the bounds on your fzero.
without success today, Can you help me?
Please do not close questions that have an answer. You were informed about that before https://www.mathworks.com/matlabcentral/answers/405170-how-to-solve-transcendental-equations-in-matlab#comment_578186
I spent several hours working on your question, and it is very frustrating that my answer just disappeared.

Sign in to comment.

1 Answer

Answer by Walter Roberson
on 15 Jun 2018
 Accepted Answer

I got code to work... but the trend is exactly opposite of what you are looking for.
Q = @(v) sym(v, 'r');
arctan = @atan;
a = Q(1*10^-6);
c = Q(3*10^8);
m1 = Q(2.2);
m2 = Q(1.5);
m3 = Q(1);
syms b p r
eqn = (a/c)*sqrt((b*m1)^2-(p*c)^2)-atan( sqrt(( (p*c)^2-(b*m2)^2 ) /( (b*m1)^2-(p*c)^2 ) ))-atan( sqrt(( (p*c)^2-(b*m3)^2 ) /( (b*m1)^2-(p*c)^2 ) ) ) - r*pi;
B0min = Q(3000000000000000)*sqrt(Q(259))*arctan(5*sqrt(Q(259))*sqrt(Q(5))*(1/Q(259)))*(1/Q(259));
B0max = Q(3*10^15);
B0 = linspace(B0min, B0max, 100);
nB0 = length(B0);
R = 0 : 2;
nR = length(R);
ps = zeros(nB0, nR, 'sym');
for ridx = 1 : nR
this_r = R(ridx);
eqnr = subs(eqn, r, this_r);
for K = 1 : nB0
this_b = B0(K);
this_eqn = subs(eqnr, b, this_b);
sols = vpasolve(this_eqn, p, [this_b/200000000, 11*this_b/1500000000]);
if isempty(sols)
fprintf('No symbolic solution for eqn, r = %d, b = %g\n', this_r, this_b);
sols = nan;
else
% fprintf('symbolic okay for eqn, r = %d, b = %g\n', this_r, this_b);
end
ps(K, ridx) = sols;
end
end
plot(ps, B0)
xlabel('p')
ylabel('b')
legend( sprintfc('r = %d', R) )

  13 Comments

It is not clear how the variable names there relate to your m1, m2, m3 ?
I do not know. I notice that they are using different c values in the diagram, but you set up your equation with only one c value.

Sign in to comment.