Answer for Problem 6.2
Taking the derivative of both sides of :
gives:
g'(x,t) = (x/2)(g(x,t) + g(x,t)/t^2) = SnJn(x)t^(n-1)
Because the bounds on the summation are infinite:
SJn(x)t^(n-1) = SJn+1(x)t^(n)
and
SJn(x)t^(n-2) = SJn+2(x)t^(n)
substituting the original formula for g into the derivative we obtain
g'(x,t) = (x/2)(SJn(x)t^n + SJn(x)t^(n-2) =SnJn(x)t^(n-1)
combining or infinite sum formulas with this one leads to
g'(x,t) = (x/2)(SJn(x)t^n + SJn-2(x)t^n =SnJn-1(x)t^n
S(t^n)(x/2)(Jn(x)+Jn-2(x)) =Sn(t^n)Jn-1(x)
After shifting the summation to n=n+1, using the orthogonal polynomial method discussed in class we can eliminate the summation, divide by t^n to achieve the desired recursion formula.
(2n/x)Jn(x) = Jn+1(x) + Jn-1(x) (1)
Knowing that Jn-1(x) = (-1)^n Jn(x) and the summation goes from n=-infinity to n=+infinity it is obvious that all even terms will be counted twice and all odd terms will not be counted. This gives us
g(x,t) = J0(x) + 2*J2(x)t^2 + 2*J4(x)t^4 .....
Seeing as the J terms are independent of t, we can choose t=1 and achieve the desired normalization equation.
g(x,1) = e^0 = 1 = J0(x) + 2*J2(x) + 2*J4(x) + 2*J6(x) .....
These results were used to program the class BesselFunc below. The equation
J'n(x) = (Jn-1(x) - Jn+1(x))/2
was used to calculate the derivative, allowing for a Newton's Approximation of the zero's. This allowed for a calculation of the cut-off wavelength of a wave guide given by Tnm where n is the order of the Bessel Function and the solution is given by the mth zero of
Jn(2p/lc) = 0.
Below is a plot of all the functions that were used in this calculation as well as the code developed.
#include <functional> #include <iostream> #include <math.h>
using namespace std; const double pi = acos(-1);
class Function : public unary_function<double,double>
{
public:
Function() {};
virtual result_type operator()(argument_type x){return x;};
virtual Function& getDerivative() { return *this;};
};
class Approximation
{
public:
Approximation(double lowerBound, double upperBound, Function& function)
:a(lowerBound), b(upperBound), f(function)
{};
virtual double calculate(double precision, int maxTries = 3000)
{
cout << "Default approximation not implemented" << endl;
return -1;
}
protected: double a; double b; Function& f; };
class NewtonApproximation : public Approximation
{
public:
NewtonApproximation(double lower, double upper, Function& f, Function& derivative)
:Approximation(lower,upper,f), d(derivative)
{};
double calculate(double precision=0.0, int maxTries=3000)
{
double x0 = (a+b)/2;
double xn;
double p = b-a;
while ((p > precision) && (--maxTries))
{
//cerr << "x=" << x0 << "\tf=" << d(x0) <<"\td=" << f(x0) << endl;
xn = x0 - f(x0)/d(x0);
if (xn>b) xn = (x0+b)/2;
if (xn<a) xn = (x0+a)/2;
p = fabs(x0-xn);
x0=xn;
//cout << p << endl;
} if (!maxTries) cerr << "maxTries reached in Newton Approximation" << endl; return xn; }; protected: Function& d; };
//
//This implementation uses Miller's Algorithm to compute Bessel Function using it's
// the backward recursive algorithm Jn+1(x) + Jn-1(x) = (2n/x)Jn(x). This recursion
// will converge up to a constant coefficient given by the renormalization formula
// 1 = J0(x) + 2*J2(x) + 2*J4(x) ...
//
class BesselFunc : public Function
{
public:
BesselFunc(int n, int maxN, double jInit=1)
: _n(n), _maxN(maxN), _jInit(jInit), j(0), jprime(0), d(0), sum(0)
{
//cout << "BesselFunc Constructor" << endl;
};
virtual ~BesselFunc()
{
if (j)
{
free(j);
free(jprime);
}
if (d) delete d;
};
BesselFunc(const BesselFunc& bf)
{
_n = bf._n;
_maxN = bf._maxN;
_jInit = bf._jInit;
sum = bf.sum;
j = bf.j;
jprime = bf.jprime;
last_x = bf.last_x;
d = bf.d;
}
virtual double operator()(double x)
{
if (j && (last_x==x)) return j[_n]/sum;
last_x = x;
sum = 0;
if (!j)
{
j = (double*)malloc(sizeof(double)*(_maxN+2));
jprime = (double*)malloc(sizeof(double)*(_maxN+2));
memset(j, 0, sizeof(double)*(_maxN+2));
};
j[_maxN]=_jInit;
for (int i=_maxN; i > 0;i--)
{
j[i-1] = (2*i/x)*j[i] -j[i+1];
jprime[i] = (j[i-1]-j[i+1])/2;
if (!(i%2)) sum+= j[i];
if (sum > 1.0e+100) //for large _maxN parameters will run over maxDouble
{
sum /=1.0e+100;
for (int k=_maxN;k>=i;k--)
{
jprime[i] /=1.0e+100;
j[k] /= 1.0e+100;
}
j[i-1] /=1.0e+100;
};
//cout << "j[" << i << "]=" << j[i] << "\tsum=" << sum << endl;
};
sum = 2*sum + j[0]; //normalize using equality 1 = j0(x) + 2*j2(x) + 2*j4(x) + 2*j6(x) ...
jprime[0] = -j[1];
//cout << "BesselFunc(" << _n << ", " << x << " ) returned : " << (j[_n]/sum) << endl;
return j[_n]/sum;
};
friend class BesselDeriv; //recursive derivative structure (should not be applied in this case)
virtual Function& getDerivative();
protected: int _n; int _maxN; double _jInit; double sum; double *j; double *jprime; double last_x; Function * d; };
//class used to demonstrate structure of recursive getDerivative function
//was not necessary for this example
class BesselDeriv : public BesselFunc
{
public:
BesselDeriv(BesselFunc* parent) :_parent(parent) ,BesselFunc(*parent)
{
d = 0;
if (parent->j)
{
_n=parent->_n;
sum = parent->sum;
j = parent->jprime;
_maxN--;
//calc_jprime(); //pro-active recursion
}
};
virtual double operator()(double x)
{
if (_parent->j && (_parent->last_x==x))
{
j = _parent->jprime;
sum = _parent->sum;
}
else
{
last_x=x;
(*_parent)(x);
_n = _parent->_n;
sum = _parent->sum;
j = _parent->jprime;
//calc_jprime(); //recursion unneeded (only first derivative needed)
}
return j[_n]/sum;
}
protected:
void calc_jprime()
{
jprime[0] = -j[1];
for (int i = 1;i<_maxN;i++)
jprime[i] = (j[i+1]-j[i-1])/2;
}
BesselFunc* _parent;
};
Function& BesselFunc::getDerivative()
{
if (!d)
d = new BesselDeriv(this);
return *d;
};
int main(int argc, char* argv[])
{
double maxN = 500;
cout.precision(15);
cout.width(15);
cout << "# Welcome to the waveguide mode calculation using the MillerApproximation\n"; cout << "# of the Bessel Recursion:" << endl << "#\n"; cout << "# (2n/x)*J(n;x)-J(n+1;x) = J(n-1;x)\n" << endl; cout << "# The waveguide modes are defined by T(n;m;) where n is the order of the\n"; cout << "# Bessel Function and the mth zero gives the solution to the longest\n"; cout << "# wavelength (a cut-off wavelength lc) where J(2pi/lc) =0 where pi=" << pi << endl; cout << "# and the recursive relation was begun at " << maxN; cout << endl;
BesselFunc millers0(0,maxN); BesselFunc millers1(1,maxN); BesselFunc millers2(2,maxN); NewtonApproximation m0(2,4, millers0, millers0.getDerivative()), m1(3,4, millers1, millers1.getDerivative()), m2(4,6, millers2, millers2.getDerivative());
cout << "# For 0th order mode, lc = " << (2*pi/(double)m0.calculate()) <<endl; cout << "# For 1st order mode, lc = " << (2*pi/(double)m1.calculate()) << endl; cout << "# For 2nd order mode, lc = " << (2*pi/(double)m2.calculate()) << endl;
cout << "\n# The graphs of the Bessel functions and their derivatives are given by" << endl;
double x;
for (x=0.1;x<10;x+=0.1)
{
cout << x << "\t" << millers0(x) << "\t" << millers1(x) << "\t" << millers2(x)
<< "\t" << millers0.getDerivative()(x) << "\t" << millers1.getDerivative()(x) << "\t" << millers2.getDerivative()(x)<< endl;
};
for (int n=0;n<3;n++)
{
int m = 1;
double b;
BesselFunc millers(n, 500, 1);
NewtonApproximation zeroOfMiller(1,2*sqrt(n)+6, millers, millers.getDerivative());
b = zeroOfMiller.calculate();
cout << "#Bessel = " << b << "\tT(" << n << ";" << m << ";) = " << 2*pi/b << endl;
};
return 0;
}
millersdata.txt (output of millers.cpp)
# Welcome to the waveguide mode calculation using the MillerApproximation # of the Bessel Recursion: # # (2n/x)*J(n;x)-J(n+1;x) = J(n-1;x)
# The waveguide modes are defined by T(n;m;) where n is the order of the # Bessel Function and the mth zero gives the solution to the longest # wavelength (a cut-off wavelength lc) where J(2pi/lc) =0 where pi=3.14159265358979 # and the recursive relation was begun at 500 # For 0th order mode, lc = 2.61274057366553 # For 1st order mode, lc = 1.63978795764418 # For 2nd order mode, lc = 1.22345159707862
# The graphs of the Bessel functions and their derivatives are given by 0.1 0.99750156206604 0.049937526036242 0.00124895865879992 -0.049937526036242 0.49812630170362 0.0249583528602436 0.2 0.990024972239576 0.099500832639236 0.00498335415278356 -0.099500832639236 0.492520809043396 0.0496672911114004 0.3 0.977626246538296 0.148318816273104 0.011165861949064 -0.148318816273104 0.483230192294616 0.0738797366126776 0.4 0.960398226659564 0.196026577955319 0.0197346631170303 -0.196026577955319 0.470331781771267 0.0973532623701674 0.5 0.938469807240813 0.242268457674874 0.0306040234586826 -0.242268457674874 0.453932891891065 0.119852363840143 0.6 0.912004863497211 0.286700988063916 0.0436650967158417 -0.286700988063916 0.434169883390685 0.141150665677777 0.7 0.881200888607405 0.328995741540059 0.0587869443641917 -0.328995741540059 0.411206972121607 0.161033043356654 0.8 0.84628735275048 0.36884204609417 0.0758177624849447 -0.36884204609417 0.385234795132768 0.179297639881808 0.9 0.807523798122545 0.405949546078806 0.0945863042748011 -0.405949546078806 0.356468746923872 0.19575775880147 1 0.765197686557967 0.440050585744933 0.1149034849319 -0.440050585744933 0.325147100813033 0.210243615881133 1.1 0.719622018527511 0.470902394866293 0.136564153956658 -0.470902394866293 0.291528932285427 0.222603933126915 1.2 0.671132744264363 0.498289057567215 0.159349018347663 -0.498289057567215 0.25589186295835 0.23270736032111 1.3 0.620085989561509 0.52202324741466 0.183026698768738 -0.52202324741466 0.218529645396386 0.240443710847372 1.4 0.566855120374289 0.541947713930855 0.207355899526932 -0.541947713930855 0.179749610423678 0.245725000320952 1.5 0.511827671735918 0.5579365079101 0.232087672144215 -0.5579365079101 0.139869999795852 0.24848627838448 1.6 0.455402167639381 0.56989593526168 0.25696775143772 -0.56989593526168 0.0992172081008304 0.248686245964531 1.7 0.397984859446109 0.577765231529023 0.281738942352741 -0.577765231529023 0.0581229585466839 0.246307652290504 1.8 0.339986411042558 0.581516951731165 0.306143535325403 -0.581516951731165 0.0169214378585774 0.241357468036273 1.9 0.281818559374385 0.581157072713434 0.329925727692387 -0.581157072713434 -0.024053584159001 0.233866833037237 2 0.223890779141235 0.576724807756873 0.352834028615638 -0.576724807756873 -0.0644716247372012 0.223890779141236 2.1 0.16660698033199 0.568292135757039 0.374623625150904 -0.568292135757039 -0.104008322409457 0.211507730851416 2.2 0.110362266922174 0.555963049819064 0.395058687458793 -0.555963049819064 -0.14234821026831 0.196818788492888 2.3 0.0555397844456016 0.539872532604314 0.413914591732062 -0.539872532604314 -0.17918740364323 0.17994680066339 2.4 0.00250768329724337 0.520185268181931 0.430980040187699 -0.520185268181931 -0.214236178445228 0.161035234692182 2.5 -0.0483837764681984 0.497094102464274 0.446059058439617 -0.497094102464274 -0.247221417453908 0.14024685571258 2.6 -0.0968049543970388 0.470818266517578 0.458972851718253 -0.470818266517578 -0.277888903057646 0.117762226734307 2.7 -0.142449370046012 0.441601379118253 0.469561502726199 -0.441601379118253 -0.306005436386106 0.0937780437655126 2.8 -0.185036033364388 0.409709246852288 0.477685495401737 -0.409709246852288 -0.331360764383062 0.0685053215653339 2.9 -0.224311545791968 0.375427481813095 0.483227050490655 -0.375427481813095 -0.353769298141312 0.0421674469919543 3 -0.260051954901934 0.339058958525936 0.486091260585891 -0.339058958525936 -0.373071607743913 0.0149981181353421 3.1 -0.292064347650698 0.300921133101057 0.486207014167509 -0.300921133101057 -0.389135680909104 -0.0127608115231421 3.2 -0.320188169657123 0.261343248780504 0.483527700144938 -0.261343248780504 -0.401857934901031 -0.0408615638100822 3.3 -0.344296260398885 0.22066345298524 0.478031686450546 -0.22066345298524 -0.411163973424715 -0.0690527206211508 3.4 -0.364295596762001 0.179225851681506 0.469722568339357 -0.179225851681506 -0.417009082550679 -0.097081541459292 3.5 -0.380127739987264 0.137377527362327 0.458629184194307 -0.137377527362327 -0.419378462090786 -0.124696292177278 3.6 -0.391768983700798 0.0954655471778758 0.444805398799618 -0.0954655471778758 -0.418287191250208 -0.151648563266356 3.7 -0.399230203371191 0.053833987745461 0.428329656206576 -0.053833987745461 -0.413779929788883 0 3.8 -0.402556410178564 0.0128210029267309 0.409304306455791 -0.0128210029267309 -0.405930358317178 -0.202602316260527 3.9 -0.40182601488764 -0.0272440396207807 0.387854712518009 0.0272440396207807 -0.394840363702824 -0.226143892194118 4 -0.397149809863847 -0.0660433280235498 0.364128145852072 0.0660433280235498 -0.38063897785796 -0.248107400949586 4.1 -0.388669679835854 -0.103273257747339 0.338292480934712 0.103273257747339 -0.363481080385283 -0.268293980154516 4.2 -0.376557054367567 -0.138646942126047 0.310534700974212 0.138646942126047 -0.34354587767089 -0.286520609256624 4.3 -0.361011117236535 -0.171896560221541 0.2810592287614 0.171896560221541 -0.321035172998967 -0.302621782901262 4.4 -0.342256790003886 -0.202775521923087 0.250086098220664 0.202775521923087 -0.296171444112275 -0.316451021114298 4.5 -0.320542508985121 -0.231060431923371 0.217848983685846 0.231060431923371 -0.269195746335483 -0.327882202450413 4.6 -0.296137816574141 -0.256552836097445 0.184593105227426 0.256552836097445 -0.240365460900784 -0.336810707935456 4.7 -0.269330789419753 -0.279080735843335 0.150573029486419 0.279080735843335 -0.209951909453086 -0.343154365412024 4.8 -0.240425327291184 -0.298499858099558 0.116050386416368 0.298499858099558 -0.178237856853776 -0.346854185773044 4.9 -0.209738327585327 -0.31469467101519 0.0812915230893305 0.31469467101519 -0.145514925337329 -0.34787488452104 5 -0.177596771314339 -0.327579137591465 0.0465651162777529 0.327579137591465 -0.112080943796046 -0.346205184102566 5.1 -0.144334747060501 -0.337097202018232 0.0121397658768809 0.337097202018232 -0.0782372564686911 -0.341857894518969 5.2 -0.110290439790987 -0.343223005871922 -0.0217184086212903 0.343223005871922 -0.0442860155848486 -0.33486977178681 5.3 -0.0758031115855851 -0.345960833801186 -0.0547481464525984 0.345960833801186 -0.0105274825664934 -0.325301155894545 5.4 -0.0412101012449924 -0.345344790779586 -0.0866953768215211 0.345344790779586 0.0227426377882644 -0.313235391956801 5.5 -0.00684386941782036 -0.341438215429044 -0.117315481647286 0.341438215429044 0.055235806114733 -0.298778040284576 5.6 0.0269708846851132 -0.334332836291008 -0.146375469074759 0.334332836291008 0.086673176879936 -0.282055883050022 5.7 0.0599200097240361 -0.324147680222857 -0.173656037872407 0.324147680222857 0.116788023798222 -0.263215737109731 5.8 0.0917025675748147 -0.311027744303943 -0.198953513886519 0.311027744303943 0.145328040730667 -0.242423084343074 5.9 0.122033354592821 -0.295142444729017 -0.222081640941641 0.295142444729017 0.172057497767231 -0.21986053254541 5.99999999999999 0.150645257250996 -0.276683858127567 -0.242873209960184 0.276683858127567 0.19675923360559 -0.195726121474172 6.09999999999999 0.177291422242742 -0.255864772558384 -0.261181511606147 0.255864772558384 0.219236466924444 -0.170231490064566 6.19999999999999 0.201747222948903 -0.232916567073224 -0.276881599424136 0.232916567073224 0.23931441118652 -0.143599922097696 6.29999999999999 0.22381200613219 -0.208086940207244 -0.289871352229728 0.208086940207244 0.256841679180959 -0.116064288705743 6.39999999999999 0.243310604823405 -0.181637509024185 -0.300072326393463 0.181637509024185 0.271691465608434 -0.0878649070262278 6.49999999999999 0.260094605581605 -0.153841301409974 -0.307430390630828 0.153841301409974 0.283762498106217 -0.0592473350620266 6.59999999999999 0.274043360624145 -0.124980165160566 -0.311916137945529 0.124980165160566 0.292979749284837 -0.0304601233588903 6.69999999999999 0.285064737710575 -0.0953421180413887 -0.313525071454273 0.0953421180413887 0.299294904582424 -0.00175254447294877 6.79999999999999 0.293095603104272 -0.0652186634016882 -0.312277562928298 0.0652186634016882 0.302686583016285 0.0266276786360466 6.89999999999999 0.29810203540482 -0.0349020961046209 -0.308218585000362 0.0349020961046209 0.303160310202591 0.0544366241853392 6.99999999999999 0.300079270519556 -0.00468282348234844 -0.301417220085941 0.00468282348234844 0.300748245302748 0.0814363822564919 7.09999999999999 0.29905138050155 0.0251532742539077 -0.291965951134252 -0.0251532742539077 0.295508665817901 0.10739720415088 7.19999999999999 0.295070691400958 0.0543274202223642 -0.279979741339191 -0.0543274202223642 0.287525216370075 0.132099570594362 7.29999999999999 0.288216947635015 0.0825704304932552 -0.265594911883438 -0.0825704304932552 0.276905929759227 0.155336159776389 7.39999999999999 0.278596232657479 0.109625094853988 -0.248967828642887 -0.109625094853988 0.263782030650183 0.176913697189904 7.49999999999999 0.26633965788038 0.135248427579703 -0.230273410525792 -0.135248427579703 0.248306534203086 0.196654670386581 7.59999999999999 0.251601833849978 0.159213768396354 -0.209703473745674 -0.159213768396354 0.230652653797826 0.214398893066269 7.69999999999999 0.234559139586466 0.181312715324586 -0.187464927813847 -0.181312715324586 0.211012033700157 0.230004904367143 7.79999999999999 0.215407807746265 0.201356872755894 -0.163777840372959 -0.201356872755894 0.189592824059612 0.243351190800242 7.89999999999999 0.194361844841281 0.219179399921749 -0.138873389164889 -0.219179399921749 0.166617617003085 0.254337219963493 7.99999999999999 0.171650807137557 0.234636346853913 -0.112991720424078 -0.234636346853913 0.142321263780818 0.262884276959932 8.09999999999999 0.147517454044381 0.247607766981591 -0.0863797338020124 -0.247607766981591 0.116948593923197 0.268936096315422 8.19999999999999 0.122215301784141 0.25799859764868 -0.0592888145527554 -0.25799859764868 0.0907520581684481 0.272459284124962 8.29999999999999 0.096006100895014 0.265739302041863 -0.0319725341379384 -0.265739302041863 0.0639893175164762 0.273443527135342 8.39999999999999 0.0691572616569891 0.270786268276835 -0.00468434063869496 -0.270786268276835 0.036920801147842 0.271901587476524 8.49999999999999 0.0419392518429383 0.273121963674054 0.0223247396097803 -0.273121963674054 0.00980725611657903 0.26786908376587 8.59999999999999 0.0146229912787451 0.272754844545881 0.0488083679179715 -0.272754844545881 -0.0170926883196132 0.261404061309143 8.69999999999999 -0.0125227324496609 0.269719024092108 0.0745271058041686 -0.269719024092108 -0.0435249191269148 0.25258635609115 8.79999999999998 -0.039233803176538 0.264073703239679 0.0992505539128287 -0.264073703239679 -0.0692421785446833 0.241516759168581 8.89999999999998 -0.0652532468512403 0.255902371443977 0.122759397737527 -0.255902371443977 -0.0940063222943839 0.228315989929926 8.99999999999998 -0.0903336111828722 0.245311786573327 0.144847341532501 -0.245311786573327 -0.117590476357686 0.213123488454994 9.09999999999998 -0.114239232683195 0.232430745005859 0.165322912904263 -0.232430745005859 -0.139781072793729 0.196096038873053 9.19999999999998 -0.13674837076486 0.217408654960452 0.184011121843219 -0.217408654960452 -0.16037974630404 0.177406237168448 9.29999999999998 -0.1576551899434 0.200413927843705 0.200754959372154 -0.200413927843705 -0.179205074657777 0.157240818301307 9.39999999999998 -0.176771572751505 0.181632204007025 0.215416722540233 -0.181632204007025 -0.196094147645869 0.135798858785699 9.49999999999998 -0.193928747687419 0.161264430757534 0.22787915416269 -0.161264430757534 -0.210903950925055 0.113289871986441 9.59999999999998 -0.20897871836887 0.13952481174069 0.238046387481514 -0.13952481174069 -0.223512552925192 0.0899318143487076 9.69999999999998 -0.221795482031721 0.116638647900217 0.245844687784343 -0.116638647900217 -0.233820084908032 0.06594902155293 9.79999999999998 -0.232276027579366 0.0928400911128149 0.251222984949328 -0.0928400911128149 -0.241749506264347 0.0415700941843805 9.89999999999998 -0.240341105534759 0.0683698322836969 0.254153192864799 -0.0683698322836969 -0.247247149199779 0.0170257529170708 9.99999999999998 -0.245935764451347 0.0434727461688663 0.254630313685121 -0.0434727461688663 -0.250283039068234 -0.00745331656815796 #Bessel = 2.40482555769577 T(0;1;) = 2.61274057366553 #Bessel = 3.83170597020751 T(1;1;) = 1.63978795764418 #Bessel = 5.13562230184068 T(2;1;) = 1.22345159707862