Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I'm trying to find the roots of the following equation:

enter image description here

I need to find λs for different values of ξ. I know that for all values of ξ, I must have n<λ<n+1/2. I use FindRoot to do this. When I put ξ^2=2, I can find all the correct roots. For example, for the first 5 roots, I have

FindRoot[x - Pi*Cot[Pi*x] == 0, {x, 0.5}]
(* Out[173]= {x -> 0.454288} *)

FindRoot[x - Pi*Cot[Pi*x] == 0, {x, 1.5}]
(* Out[174]= {x -> 1.36917} *)

FindRoot[x - Pi*Cot[Pi*x] == 0, {x, 2.5}]
(* Out[175]= {x -> 2.29891} *)

FindRoot[x - Pi*Cot[Pi*x] == 0, {x, 3.5}]
(* Out[176]= {x -> 3.24485} *)

FindRoot[x - Pi*Cot[Pi*x] == 0, {x, 4.5}]
(* Out[177]= {x -> 4.20427} *)

Comparing with the following plot, these roots are correct (The roots are where the red line is coinciding with the green curve):

enter image description here

But when I put some different value for ξ^2, for example ξ^2=0.5 (which is a value which I need to use in my calculations), FindRoot doesn't give me the correct roots anymore:

FindRoot[x - 0.25*Pi*Cot[Pi*x] == 0, {x, 0.5}]
(* Out[191]= {x -> 0.362393} *)

FindRoot[x - 0.25*Pi*Cot[Pi*x] == 0, {x, 1.5}]
(* Out[192]= {x -> 1.18617} *)

FindRoot[x - 0.25*Pi*Cot[Pi*x] == 0, {x, 2.5}]
(* Out[193]= {x -> 2.11326} *)

FindRoot[x - 0.25*Pi*Cot[Pi*x] == 0, {x, 3.5}]
(* Out[194]= {x -> 2.11326} *)

FindRoot[x - 0.25*Pi*Cot[Pi*x] == 0, {x, 4.5}]
(* Out[195]= {x -> 3.07949} *)

FindRoot[x - 0.25*Pi*Cot[Pi*x] == 0, {x, 5.5}]
(* Out[196]= {x -> 5.04912} *)

Which are not the correct roots (specially the 4th and 5th ones) compared to this plot:

enter image description here

Why doesn't FindRoot work correctly? What am I doing wrong?

share|improve this question
1  
It's finding a root, which is all it is asked to do. Starting a bit closer might give the root that you want. Try changing your .5 to x.1. –  wxffles Jun 26 '13 at 22:01
    
Indeed, put x - 0.25*Pi*Cot[Pi*x] /. in front of each FindRoot and you'll see it has found a correct root. –  Sjoerd C. de Vries Jun 26 '13 at 22:02
    
@Sjoerd C. de Vries, Thanks a lot for your answer. But I didn't understand what I should do. Where should I put x - 0.25*Pi*Cot[Pi*x]? –  ZKT Jun 26 '13 at 22:06
    
x - 0.25*Pi*Cot[Pi*x] is the expression that you want to set to zero for a certain x. FindRoot finds this x. If you now write x - 0.25*Pi*Cot[Pi*x] /. FindRoot[...] the found value of x is filled in in the expression and you see that it is indeed (close to) zero. –  Sjoerd C. de Vries Jun 27 '13 at 13:30
add comment

2 Answers 2

up vote 3 down vote accepted

You can use NSolve with a condition, instead. See:

    eq = x - (Pi/2)*0.5*Cot[Pi*x] == 0
    NSolve[{eq, 3.5 < x < 4.5}, x][[1]]

The output being

    {x -> 4.06081}

which correctly falls between 3.5 and 4.5.

share|improve this answer
    
This looks only good at the first glance because it doesn't find the root at 17.0059 when you use 16.5<x<17.5 and xi is 0.1. On the other hand I have to admit, that NSolve finds most of the roots. I wonder whether one can instruct FindRoot to use the same method, because besically NSolve only knows the interval boundaries. The same is possible with FindRoot but it doesn't seem to work reliable with this. –  halirutan Jun 27 '13 at 20:31
    
But if you change the lower bound a bit, it gives the corretc result: NSolve[{x - (Pi)*0.1*Cot[Pi*x] == 0, 16.4 < x < 17.5}, x], the output being {x -> 17.0059}. –  fcpenha Jun 27 '13 at 21:27
    
The problem is, that you maybe don't want to fix all those things when the OP wants to extract all roots in the Range[0.001,1.,0.001] as stated in one of his/her comments. –  halirutan Jun 27 '13 at 21:29
    
Indeed. In this case I would subtract a small number of the lower bound for each computation. But probably this fix is specific to this function. I wouldn't know what to do in a general situation. –  fcpenha Jun 27 '13 at 21:35
    
@fcpenha Thanks a lot for your answer, it indeed solved my problem :). –  ZKT Jun 27 '13 at 22:33
show 1 more comment

I see two potential pitfalls in your approach. First, you need to understand that FindRoot is a numerical procedure which starts at a certain point and tries to find a root by approximating the gradient and moving towards it$^{1}$. Therefore, different starting points might lead to different roots, although you probably would expect they give the same root.

Therefore, it can happen if you choose the start point badly, you end up with a completely unexpected root since maybe the gradient was too large and introduced a very big first step:

FindRoot[x - 0.25 Pi*Cot[Pi*x] == 0, {x, 4.99}]
(* {x -> 2.11326} *)

Second, when you are completely out of luck, then it might happen you choose a starting point which is a singularity. Although you get a warning, you might over-read this and use this wrong root:

FindRoot[x - 0.25 Pi*Cot[Pi*x] == 0, {x, 5.0}]
(* {x -> 5.} *)

This is why I suggest for your example a different approach. You could extract all roots in a certain interval and delete those, which are singularities. This is in your case easy, because you used Cot and the singularities are the roots where Sin[Pi x] is zero. After you collected all roots, you can use Nearest to create a function which indeed gives always the nearest root to your input:

With[{roots = 
   Quiet[Select[
     Union[x /. (FindRoot[x - 0.25 Pi*Cot[Pi*x] == 0, {x, #1}] & /@ 
         Range[0.1, 10, .1])], Abs[Sin[Pi #]] > 10.0^-5 &]]},
 giveRootAt[x_?NumericQ] := First[Nearest[roots, x]]
]

A quick check reveals that you get the correct roots now:

{#,giveRootAt[#]}&/@Range[.5,5.5]
(* 
  {{0.5,0.362393},{1.5,1.18617},{2.5,2.11326},
  {3.5,3.07949},{4.5,4.06081},{5.5,5.04912}}
*)

Footnote 1: This is not entirely correct. FindRoot has several modes and options where you can influence the behavior. See the reference page for more information.

Update regarding your comment

Beware that for very small $\xi$ it's getting harder to extract the roots reliably. Therefore, I will work with a minimal $\xi$ of 0.1. What you can do is to use the above method and create several giveRootAt functions. I will give for illustration purposes a slightly different version which defines several giveAllRoots function but the definition of giveRootAt is equivalent. The main difference in the code below is, that the whole definition is now a Function which takes $\xi$ as parameter. With this, we can map it over the list of $\xi$ values you want:

Function[xi,
  With[{roots = 
    Quiet[x /. (FindRoot[x - xi*Pi*Cot[Pi*x] == 0, {x, #}] & /@ Range[0.1, 20, xi/2])]},
    giveAllRoots[xi] = Select[Union[roots], Chop[Sin[Pi #]] != 0 &]
  ]
] /@ Range[0.1, 1.0, 0.001];

Now you have a function for which gives you the roots for different values of $\xi$ ranging from 0.1 to 1.0 in 0.001 steps.

To test this you can use for instance a Manipulate showing you the roots

Manipulate[
 Plot[x - \[Xi]*Pi*Cot[Pi*x], {x, 0, xmax}, 
  Epilog -> {Red, Point[{#, 0}] & /@ giveAllRoots[\[Xi]]}, 
  PlotRange -> {Automatic, {-3, 3}}],
 {xmax, 2, 20},
 {\[Xi], 0.1, 1.0, 0.001}
]

enter image description here

share|improve this answer
    
Thanks a lot for your answer. It solves the problem for only one value of ξ. What I did before was that I made a table with respect to ξ, and I used Map to put FindRoot in each element of the table: table=Table[{ξ,x-ξ*Pi*Cot[Pi*x]},{ξ,0.001,1.,0.001}]; final =table /. {{a_,c_} :> {a, FindRoot[c, {x,0.5}]}; This gives me all λ_0s for different values of ξ. Then I used a For loop to find other λs, putting FindRoot[c, {x,0.5+i}]}, i changing form 0 to 20. How can I do all these things using your method? –  ZKT Jun 26 '13 at 23:12
    
@ZKT Please see my update in the answer. –  halirutan Jun 27 '13 at 1:37
    
Dear @halirutan , thanks a lot for the time you put to answer my question. Your method was useful, but since I need to consider very small values of xi in my code (as they are physical parameters and I can't just ignore them), so NSolve works better in my case. Thanks again by the way :). –  ZKT Jun 27 '13 at 22:39
    
@ZKT Please, then don't accept my answer as final solution. fcpenha's solution is nice and small and if it solves your problem, you really should keep this as accepted answer. If you like my answer, upvote it, but accept the one the solves your problem. –  halirutan Jun 27 '13 at 22:43
    
Did as you said. I learnt things from your answer by the way :)... –  ZKT Jun 27 '13 at 23:06
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.