Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I've never used Mathematica before and am trying to numerically solve equation (12) from this paper. Ideally I'd be able to find the smallest value of $x_{n\nu}$ for $\exp(-kr\pi)$ close to 1, and close to 0 for some range of $\nu$ (must be larger than 2) and then plot it.

I've tried a few methods I found after googling, but none seem to work.

I tried to use this by writing:

f= (2*BesselJ[a, x] + 
      x*(BesselJ[a + 1, x] + BesselJ[a - 1, x]))*(2*
       BesselY[a, x*exp[-Pi]] + 
      x*exp[-Pi]*(BesselY[a + 1, x*exp[-Pi]] + 
         BesselY[a - 1, x*exp[-Pi]])) - (2*BesselY[a, x] + 
      x*(BesselY[a + 1, x] + BesselY[a - 1, x]))*(2*
       BesselJ[a, x*exp[-Pi]] + 
      x*exp[-Pi]*(BesselY[a + 1, x*exp[-Pi]] + 
         BesselY[a - 1, x*exp[-Pi]]));
 sol[_a] = NSolve[f ==0 && x>0 && x < 10 && a > 2 && a < 12, x, Reals);

and was told that NSolve can't solve it.

I've also tried using

sol[a_] := x /. FindRoot[f, {{x, 0}, {a, 2}}]; 

which throws up an error:

SetDelayed::write : Tag Plus in *My definition of f here*[a_, x_] is Protected.
share|improve this question
One note, Exp is not exp. _a is an element with Head a when you are probably looking for a_ which is named pattern. You can't go on deep water when you are not familiar with basics of syntax. – Kuba Mar 20 '14 at 13:08
I mean, please take some time to parse documentation examples for NSolve. Also functions construction tutorial will explain you a lot. – Kuba Mar 20 '14 at 13:39

Well, first in Mathematica one should write Exp instead of exp. Besides that, do you know from some source that your equation has analytic solution? Then it would be helpful, if you write this along with the reason. If you have no such information, I would bet that it has not.

In this case you might want to solve it numerically. It may be done, for instance as follows.

First let us make sure that there are solutions in the interval you gave:

 f[x_, a_] := (2*BesselJ[a, x] +  x*(BesselJ[a + 1, x] + BesselJ[a - 1, x]))*(2*
   BesselY[a, x*Exp[-Pi]] + x*Exp[-Pi]*(BesselY[a + 1, x*Exp[-Pi]] + 
     BesselY[a - 1, x*Exp[-Pi]])) - (2*BesselY[a, x] + x*(BesselY[a + 1, x] + 
  BesselY[a - 1, x]))*(2* BesselJ[a, x*Exp[-Pi]] + x*Exp[-Pi]*(
  BesselY[a + 1, x*Exp[-Pi]] + BesselY[a - 1, x*Exp[-Pi]]));

Manipulate[Plot[f[x, a], {x, 0, 10}, PlotRange -> {-200, 1000}], 
           {a, 2, 12}]

You should see this after evaluating: enter image description here Thus there are 2 solutions, the larger disappears at a > 3.74, while the smaller at a > 6.89

Let us now solve the equation for the smallest root:

lst1 = Table[{a, FindRoot[f[x, a] == 0, {x, 5.7}][[1, 2]]}, {a, 2, 6.7, 0.1}];

and with respect to the larger one:

lst2 = Table[{a, FindRoot[f[x, a] == 0, {x, 9}][[1, 2]]}, {a, 2, 3.6, 0.1}];

Now let us plot the solutions:

 ListPlot[{lst1, lst2}, AxesLabel -> {Style["a", 16, Italic], Style["x", 16, Italic]}]

This should appear on the screen: enter image description here Have fun.

share|improve this answer

Your Answer


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.