Solving difficult trigonmetric equations for theta.

First, thank you for taking the time to read my problem. I apprecaiate any help.

Trying to find roots of the below equation. It doesn't seem to work on octave. When I use an online tool like wolframalpha, I can get the asnwer that I am looking for. I can always approximate the answer graphically but I am interested in the code to cruch answers directly. The below is the input and output. There are roots to my problem. I'm not sure why Matlab can't solve for them. The larger picture problem that I am working on is in finding vertical tangent lines to polar equations which I have done on paper but would like to code.

solve(-1*sin(x)*(4 + 2*cos(3*x + 21*3.14/6)) - cos(x)*6*sin(3*x + 21*3.14/6) == 0,x)

ans = {}(0x0)

Interestingly enough, when I do a simple example like the below, I get the obvious solution. So, I know that "solve" at least works with simple examples:

solve(sin(x) == 0,x)
ans = (sym 2×1 matrix)
  ⎡0⎤
  ⎢ ⎥
  ⎣π⎦

Answer 1

[Turlough Hughes]

polyxpoly works quite well for this situation. I use it here to find the intersection of two lines; your equation (over a specified range) and y=0

range=[0 100];
x=linspace(range(1),range(2),10000); %generate data in a given range
y=-1*sin(x).*(4 + 2*cos(3*x + 21*3.14/6)) - cos(x).*6.*sin(3*x + 21*3.14/6);

[xroots,y0] = polyxpoly(x,y,range,[0 0]); % find intersection
plot(x,y), hold on, plot(xroots,y0,'or');

xroots

[Gabriel Schneck]

Unfortunately, I do not have the ability to use polyxpoly. I decided to d0 what I did below. However, It would be nice if I can get more precise and If I can get exactly 6 solutions.

x = 0 : .000001 : 2*pi;

y1=-1*sin(x).*(4 + 2*cos(3*x + 21*pi/6)) - cos(x).*6.*sin(3*x + 21*pi/6);

y2=0*x;

s = y1 - y2;

plot(x,y1);

ix = find(s > -.00001 & s < .00001);

x_sol = x(ix)

[Turlough Hughes]

Thankfully there are many ways to do it...

First generate the data as you have done:

x = 0 : .000001 : 2*pi;
y=-1*sin(x).*(4 + 2*cos(3*x + 21*pi/6)) - cos(x).*6.*sin(3*x + 21*pi/6);

This provides initial approximations, but not the precise roots. To get initial values (one for each root) I look for indices of y values where the next one changes sign:

idx=find(y(1:end-1)<0 & y(2:end)>=0 | y(1:end-1)>=0 & y(2:end)<0);

Precise roots can be found with fzero(). To do this you would make an m-file as follows:

function  y = f(x)
y=-1*sin(x).*(4 + 2*cos(3*x + 21*pi/6)) - cos(x).*6.*sin(3*x + 21*pi/6);
end

I have called it f.m, and saved it in my working directory (current folder). Your roots are then:

xroots=nan(size(idx));
for c=1:length(idx)
    xroots(c)=fzero(@f,x(idx(c)));
end

edit: modified the logic statements to >=0 as opposed to >0

[Gabriel Schneck]

That worked very well. I appreciate your help and I learned a few coding tips as well. Thanks again,

[Turlough Hughes]

No problem, good to know. The best way to show appreciation is by accepting an answer :P

Answer 2

[Walter Roberson]

syms x
sol = solve(-1*sin(x)*(4 + 2*cos(3*x + 21*sym(pi)/6)) - cos(x)*6*sin(3*x + 21*sym(pi)/6) == 0,x,'returnconditions',true);

>> sol.x
 
ans =
 
2*pi*k - log(z1)*1i
>> sol.conditions
 
ans =
 
in(k, 'integer') & (z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 1) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 2) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 3) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 4) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 5) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 6) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 7) | z1 == root(z^8 + z^6/2 + z^5*1i - z^3*1i + z^2/2 + 1, z, 8))

What this tells you is that there are an infinite number of solutions spaced 2*pi apart, with each solution being sqrt(-1) times the log() of a root of a polynomial of degree 8.

>> vpa(simplify(subs(sol.conditions,sol.parameters(1),0),'steps',10),16)
 
ans =
 
z1 in {- 0.925847879709287 + 0.377896419191579i, - 0.8236278734749009 - 0.5671306075633836i, - 0.5339703343863142 - 0.8455032122915724i, 0.768877497083311i, 1.300597304243443i, 0.5339703343863142 - 0.8455032122915724i, 0.8236278734749009 - 0.5671306075633836i, 0.925847879709287 + 0.377896419191579i}

when you take 1i*log() of those values, 6 out of the 8 of them are real valued, so there are 6 real solutions for every period of 2*pi

It is possible to extract the values mechanically using children()

By the way:

>> simplify(rewrite(6*sin(3*x + (7*pi)/2)*(2*sin(x/2)^2 - 1) + sin(x)*(4*sin((3*x)/2 + (7*pi)/4)^2 - 6) == 0,'sin'),'steps',10)
 
ans =
 
16*sin(x)^4 + 3 == 2*sin(x)*(9*sin(x) + 1)