how to use vpasolve to solve the equation sin(x)==0 in the range of (0,5] instead of [0,5]?

I want to solve sin(x)=0 in the open interval of (0,5], but it seems vpasolve supports close interval only. So it is not able to give the answer I want.

Of course we can use solve and assume instead, but in my real case I have to use vpasolve because of accuracy reason.

[Dyuman Joshi]

"Of course we can use solve and assume instead, but in my real case I have to use vpasolve because of accuracy reason."

Can you elaborate on this?

Are you saying that there accuracy issues when using symbolic variables with solve?

[Sam Chak]

Hi @Kai Wang

You mentioned using the vpasolve() due to accuracy reasons. Could you please specify the answer you expect from solving the equation, in the left-open interval , or ?

Answer 1

[Walter Roberson]

vpasolve() does not support open or semi-open intervals.

The trick is to use

syms x
sol = vpasolve(sin(x)==0, x, [1e-63 5])

sol = 

3.1415926535897932384626433832795

I do not understand at the moment why left bounds smaller than 1e-63 appear to be treated as 0, even if you state them symbolically such as sym(10)^-80

[Dyuman Joshi]

"I do not understand at the moment why left bounds smaller than 1e-63 appear to be treated as 0, even if you state them symbolically such as sym(10)^-80"

Something one would not expect while using symbolic numbers.

How about this?

syms x
eqn = [sin(x)==0; x>0; x<=5];
x0 = solve(eqn,x)

x0 = 

π

[Walter Roberson]

vpasolve() uses symbolic numbers.

syms x
sol = vpasolve(sin(x)==sym(0), x, [sym(1e-80) sym(5)])

sol = 

0

sol = vpasolve(sin(x)==sym(0), x, [sym(1e-63) sym(5)])

sol = 

3.1415926535897932384626433832795

Notice that is more than 16 digits of output as would be expected for symbolic floating point numbers

[Dyuman Joshi]

What I meant was - vpasolve() should output pi as the output when provided with sym(1e-80) as the starting point of searching range for the given equation, but it does not.

[Walter Roberson]

syms x
sol = vpasolve(sin(x)==0, x, [1e-63 5])

sol = 

3.1415926535897932384626433832795

vpa(sol - sym(pi))

ans = 

vpa(sym(pi))

ans = 

3.1415926535897932384626433832795

Remember that π is a transcendental number -- irrational, infinite length. And that vpasolve() has to work at finite precision, because it is a numeric solver. How many digits? Well, the default is 32 digits as documented in digits and then there are also several (exact count is poorly documented) "guard" digits -- if I recall correctly it is 7 guard digits. 32+7 = 39, and error difference between the numeric solution and π is on the order of 10^-39 ... and if you vpa(pi) then you get the same result as the solution out of the vpasolve()

All of which is to say that No, vpasolve() should not have output pi it should have output the same finite approximation to π that it output.

Why it did output 0 for a left bound of 1e-80 but succeeds in finding the numeric approximation to pi with a bound of 1e-63 is something that I am not clear on at the moment.

[Walter Roberson]

Several years ago, I think I chased part way into the logic of vpasolve() . If I recall correctly, then (at least at the time) if it did a Newton-Raphson projection based upon its then-current point, and the projected value was beyond the user-supplied bounds, then vpasolve() would stop searching -- stopping instead of trying to recover by moving to the boundary or to a some fraction of the way to the boundary.

If that is what is happening, then it might be the case that with the way it chooses the starting position when given an interval, that within a small number of steps it tries to project lower than the boundary, and stops.

If you have a Newton-Raphson type search and that is leading you towards a root based upon derivative information... but you encounter a boundary in the way... then Newton-Raphson type of algorithms do not have any way of recovering from that kind of situation.

Generally speaking, if you have a function with several close 0s, then Newton-Raphson type of algorithms do not necessarily end up at the "closest" root to your current guess. N-R algorithms can "bounce around" a fair bit while searching out a root.

[Dyuman Joshi]

I am saying that for the given case it should not have given 0 as output, but it does, which is opposite of what is expected.

[Paul]

I've always found it odd that vpasolve states:

"If init_param is a matrix with two columns, then the two entries of the rows specify the bounds of an initial guess for the corresponding variables."

Taken literally, it doesn't claim to bound the solution within the bounds, only the inital guess.

Why wouldn't the doc be clearer about bounding the solution?

[David Goodmanson]

syms x
sol = vpasolve(sin(x)/x==0, x, [0 5])
sol =

3.1415926535897932384626433832795

Ok, admittedly it's a method particular to this problem. On the other hand it doesn't hurt to keep in mind that zeros can be removable, just as poles can be removable.

[Walter Roberson]

That suggests,

syms x
sol = pi * vpasolve(sinc(x) == 0, x, [0 5/pi])

sol =