How can I translate a MathWorld statement into Mathematica syntax?

Question

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

Not sure how you wish to use the identity. Related (defining differential operators): (5030), (15605), (20519), (71643), (124370). For fractional calculus: https://mathematica.stackexchange.com/search?q=fractional+derivative — Michael E2, Dec 17 '17 at 14:28

Picaud Vincent · Answer 1 · 2017-12-17T20:12:16.037

6

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]

edited Dec 17 '17 at 20:12

answered Dec 17 '17 at 14:56

Picaud Vincent

2,463
13
20

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

score 4 · Answer 2 · answered Dec 17 '17 at 19:33

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]  *)

Mathematica graphics

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

score 3 · Accepted Answer · answered Dec 25 '17 at 23:29

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"]

and load with:

<<DifferentialOperator`

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

score 2 · Answer 4 · answered Oct 01 '18 at 15:47

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]

m_goldberg · Answer 5 · 2017-12-17T16:36:58.617

1

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

edited Dec 17 '17 at 16:36

answered Dec 17 '17 at 16:29

m_goldberg

107,779
16
103
257

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

How can I translate a MathWorld statement into Mathematica syntax?

5 Answers5