Suppose that you want to find the solution(s) of a transcendental equation. For example, you might want to know what values of $x$ that solve the equation $$ \sin x = 2x\cos x.$$ The graphical approach is to plot both sides of the equation and see where they cross as shown below.
from pylab import *
#from scipy import *
x = linspace(0,3*pi,1000)
yL = sin(x)
yR = 2*x*cos(x)
figure()
plot(x,yL,label=r'$\sin x$')
plot(x,yR,label=r'$2x\cos x$')
legend(loc='upper left')
ylim(-15,15)
xlabel(r'$x$')
axhline(y=0, c="k") # draw horizontal line at y=0
show()
From the graph above, it is possible to find approximate solutions. There is a solution at $x=0$. You can also see that there are solutions between $n\pi$ and $(n+1)\pi$ for integer values of $n$.
Another way to find the solution graphically is to plot the difference of the two sides of the original equation. If we define the function
$$f(x) = \sin x - 2x\cos x,$$
the original equation has a solution where $f(x)=0$. These are known as the roots of $f(x)$. It is slightly easier to find the solutions by looking at a graph of this funciton as shown below.
x = linspace(0,3*pi,1000)
fx = sin(x) - 2*x*cos(x)
figure()
plot(x,fx)
#ylim(-15,15)
xlabel(r'$x$')
ylabel(r'$f(x)$')
axhline(y=0, c="k") # draw horizontal line at y=0
show()
There are several functions in the scipy.optimize
library for finding the roots of a function. For all of the methods, you must define the function whose roots you want to find. In general, you must know the approximate location of the root.
First, you can use the optimize.brentq
function if you know values that bracket the root. As mentioned above, the roots of the example function $f(x)$ are between $n\pi$ and $(n+1)\pi$ for integer values of $n$.
from scipy import optimize
# Define the function
def f(x):
return sin(x) - 2*x*cos(x)
# The function has roots between 0, pi, 2*pi, etc.
for n in arange(10):
print(optimize.brentq(f, n*pi, (n+1)*pi))
0.0 4.6042167772 7.78988375114 10.9499436485 14.1017251336 17.2497818346 20.3958423573 23.5407082923 26.6848024909 29.8283692131
Notice that for the range 0 to $\pi$, the function found the root $x=0$, not the root near $x=1$. The other root can be found by changing the lower limit of the range slightly as shown below.
print(optimize.brentq(f,0.1,pi))
1.16556118521
The function $f(x)$ must have opposite signs at the two limits given or the the root-finding function will fail as shown below.
print(optimize.brentq(f,4,8))
--------------------------------------------------------------------------- ValueError Traceback (most recent call last) <ipython-input-6-6e832f7132e2> in <module>() ----> 1 print(optimize.brentq(f,4,8)) /projects/sage/sage-7.3/local/lib/python2.7/site-packages/scipy/optimize/zeros.pyc in brentq(f, a, b, args, xtol, rtol, maxiter, full_output, disp) 436 if rtol < _rtol: 437 raise ValueError("rtol too small (%g < %g)" % (rtol, _rtol)) --> 438 r = _zeros._brentq(f,a,b,xtol,rtol,maxiter,args,full_output,disp) 439 return results_c(full_output, r) 440 ValueError: f(a) and f(b) must have different signs
Second, you can use the optimize.fsolve
function if you know approximate values of the roots. In the example below, initial guesses of $x\approx (n+\frac{1}{2})\pi$ for integers $n$ from 0 to 9 are used.
n = arange(10) # make a list of integers from 0 to 9
guesses = (n+0.5)*pi # make a list of guesses
print(optimize.fsolve(f,guesses))
[ 1.16556119 4.60421678 7.78988375 10.94994365 14.10172513 17.24978183 20.39584236 23.54070829 26.68480249 29.82836921]