# How can I translate a MathWorld statement into Mathematica syntax?

The Wolfram Mathworld has the following statement here:

A derivative identity for expressing higher order modified Bessel functions in terms of $I_0(x)$ is $I_n(x)=T_n(d/(dx))I_0(x)$

However, I'm not sure how to translate that into Mathematica syntax. My simple minded attempts have failed:

BesselI[l, x] /. {l -> 0.5, x -> 0.100}
D[ChebyshevT[l, x], x]*BesselI[0, x] /. {l-> 0.5, x -> 0.100}
ChebyshevT[l, D[BesselI[0, x], x]]*BesselI[0, x] /. {l -> 0.5, x -> 0.100}
ChebyshevT[l, D[BesselI[0, x], x]] /. {l -> 0.5, x -> 0.100}


Apologies for being obtuse if I've missed something obvious

You can make use of my DifferentialOperator paclet to do this. Install with:

PacletInstall["https://github.com/carlwoll/DifferentialOperator/releases/download/0.1/DifferentialOperator-0.0.1.paclet"]


<<DifferentialOperator


Then (I had to use an image since the partial derivative in the input doesn't translate well to MSE):

Waiting for a more elegant solution, you can try this:

Creates a differential operator from a polynomial in x:

DiffOp[poly_] /; PolynomialQ[poly, x] := Block[{Dop},
Dop[c_, {ord_}] := c*(D[#, {x, ord - 1}] &);
Return[Total[MapIndexed[Dop, CoefficientList[poly, x]]]];
];


For instance:

op = DiffOp[ChebyshevT[4, x]]

FullSimplify[op[BesselI[0, x]]]


BesselI[4, x]

Update: maybe slightly better (and without Return[], see Carl Woll remark):

DiffOp[poly_, x_, f_] /; PolynomialQ[poly, x] :=
CoefficientList[poly, x].Table[D[f, {x, i}], {i, 0, Exponent[poly, x] }]


For instance:

DiffOp[ChebyshevT[4, x], x, f[x]]


prints: $$8 f^{(4)}(x)-8 f''(x)+f(x)$$

and

FullSimplify[DiffOp[ChebyshevT[4, x], x, BesselI[0,x]]]


returns

Bessel[4,x]

• It seems this is peculiar to the ChebyshevT or some such Std function. If I try to use (1/2*x^2*((Sqrt[1 - 1/x^2 + 1])^n + (Sqrt[1 - 1/x^2])^n)) inplace of ChebyshevT, the function is not evaluated. I expect I've forgotten some other MMA fu? – Hedgehog Dec 18 '17 at 12:17
• @Hedgehog to define the differential operator the function expects a polynomial x^n -> d^n. In your example, you have a 1/x^2 term, thus not a polynomial. The function does not accept it due to the /; PolynomialQ[poly, x] test. If you want to manage negative x powers this is certainly possible, but you will get integro-differential operators – Picaud Vincent Dec 18 '17 at 12:27
• If you want to compute ChebT, the formula I have is: 1/2 ((x - I Sqrt[1 - x^2])^n + (x + I Sqrt[1 - x^2])^n) and there is no 1/x^2 term. Do not forget to FullSimplify the expression to let Mathematica recognize a polynomial before calling DiffOp[] (I will be back in one hour, sorry) – Picaud Vincent Dec 18 '17 at 12:34

Assuming the idea is to translate the derivative operator (which must make this Q&A a duplicate, no?), this was what I came up with...

ClearAll[CircleDot];
(* assumes op is a polynomial with constant coefficients and the variable in f is x *)
CircleDot[op_, f_] /; PolynomialQ[op, Derivative[1]] && FreeQ[op, x] :=
With[{c = CoefficientList[op, Derivative[1]]},
c.NestList[D[#, x] &, f, Length@c - 1]
];

ChebyshevT[3, Derivative[1]] \[CircleDot] BesselI[0, x] // Expand
(*  BesselI[3, x]  *)


• Interesting: Nested D[#, x] simplifies by Expand, but direct D[#, {x, k}] gives derivatives that do not simplify to the desired form. – Michael E2 Dec 17 '17 at 19:46

Operate[] combined with Through[] are suitable for this task:

k = 5; (* example *)
FullSimplify[Through[Operate[Composition[Through, ChebyshevT[k, δ]], BesselI[0, x]]] /.
{(c_ /; FreeQ[c, δ])[BesselI][0, x] :> c BesselI[0, x],
(c_. δ^j_.)[BesselI][0, x] :> c Derivative[0, j][BesselI][0, x]}]
BesselI[5, x]


Like Picaud Vincent, I create a differential operator to represent $T_n(d/(dx))$. Only I do it somewhat differently.

chebDOp[n_] :=
Module[{rules, opExpr},
rules = CoefficientRules[ChebyshevT[n, \[FormalU]]];
opExpr =
Plus @@
(#[[2]] D[BesselI[0, \[FormalU]], {\[FormalU], #[[1, 1]]}] & /@ rules);
Function[\[FormalU], \[FormalV]] /. \[FormalV] -> opExpr]


Then

And @@ Table[FullSimplify[chebDOp[n][x] == BesselI[n, x]], {n, 4}]


gives

True`

• Thanks for the contribution. Your answer gives the same as Picaud, which, ironically is what I wanted. But his seems simpler and possibly more general? in that I can pass in any polynomial and looks like it should work - I acknowledge this is not I requirement of the question. Thanks again. – Hedgehog Dec 18 '17 at 12:20