Hello Julia users,

I use Mathematica to generate analytic expressions for some quite complicated integrals I calculate.

In Mathematica I generate tables of indexes and symbolic expression e.g.

```
{{-1, -1, -1, -1, -1, -1, 0}, {(256*Pi^4*(Log[b] - Log[a + b] - Log[b + d]))/(a^2*c^2*d^2)}}
{{-1, -1, -1, -1, -1, 0, -1}, {(256*Pi^4*(Log[b] - Log[a + b] + Log[a + b + c] + Log[a + b + d]))/(a^2*c^2*d^2)}}
{{-1, -1, -1, -1, 0, -1, -1}, {(256*Pi^4*(Log[b] - 2*Log[b + c] - Log[b + d] + Log[b + c + d]))/(a^2*c^2*d^2)}}
{{-1, -1, -1, -1, 0, -1, 0}, {(-256*(b^2 - c*d)*Pi^4)/(a^2*b*c^2*(b + c)*d^2*(b + d))}}
{{-1, -1, -1, -1, 0, 0, -1}, {(256*Pi^4)/(a^2*b*c^2*d^2)}}
```

and save it to the text file (in production the formulas are often much longer, I chose some of the shortest ones for this MWE).

Then, using SymPy.jl and some manual parsing I create function to compute values of integral with given indexes ({i,j,k,l,m,n,p} that are integers) and parameters (a,b,c,d which are floats).

```
f = open("I4O2s.txt")
lines = readlines(f)
using SymPy
@vars a b c d
const sympy_parsing_mathematica = SymPy.PyCall.pyimport("sympy.parsing.mathematica")
# create new array with pairs as arguments
A = []
for i in 1:length(lines)
push!(A, (tuple(parse.(Int, split(split(strip(lines[i], ['{', '}']), "}, {")[1], ","))...), SymPy.lambdify(sympy_parsing_mathematica."mathematica"(split(strip(lines[i], ['{', '}']), "}, {")[2]), (a,b,c,d))))
end
# create dict from array
i4dict = Dict(A)
A = 0
# use dict to create function
i4(i,j,k,l,m,n,p,a,b,c,d) = i4dict[(i,j,k,l,m,n,p)](a,b,c,d)
```

It may not be elegant, but it works (all sugestions of improvement are more than welcome!).

## Offtopic question: How to view Julia expression from SymPy.lambdify?

One, a bit off-topic question would be how to look–up the Julia expression for the givent formula, since after SymPy.lambdify i only get:

```
x = SymPy.lambdify(sympy_parsing_mathematica."mathematica"(split(strip(lines[6], ['{', '}']), "}, {")[2]), (a,b,c,d))
#79 (generic function with 1 method)
```

when creating dictionary i get expressions such as `##79#81(Box(##418))`

instead of formula

```
i4dict = Dict(A)
Dict{NTuple{7,Int64},getfield(SymPy, Symbol("##79#81"))} with 36 entries:
(-1, -1, 0, 0, -1, -1, -1) => ##79#81(Box(##418))
(-1, -1, -1, 0, -1, -1, 0) => ##79#81(Box(##410))
etc...
```

and trying to view it in the dictionary gets me:

```
i4dict(-1,-1,-1,-1,-1,-1,-1)
ERROR: MethodError: objects of type Dict{NTuple{7,Int64},getfield(SymPy, Symbol("##79#81"))} are not callable
Stacktrace: [1] top-level scope at none:0
```

My problem is following:

When I try to compute the expression e.g.

i4(-1, -1, -1, -1, 0, 0, -1, BigFloat(1.0), BigFloat(1.0), BigFloat(1.0), BigFloat(1.0))

i get different result than when I do it with higher precision in Mathematica (and different from the article with some high-precision values for these expressions, they agree with the ones calculated in Mathematica).

I managed to narrow it down to the following: the problem is with the precision of pi^4. After some tests I realised that the pi is converted to Float64 and then powered. Of course this is expected behaviour but it leads to precision loss. However, due to numerical instabilities in some of the more complicated expressions I would need extended precision: DoubleDouble, Quad or BigFloats (or ArbFloats, maybe I will even think about some port of quad double and/or double quad to Julia in the future but BigFloats are probably good enough for now).

**So my question is: is there a nice, Julian way to generate function that computes with arbitrary precision floating, from symbolic expression, with automatically adjusted precision of the constant (pi in this case, but I believe the question is more broad and applies to many mathematical constans in Base.MathConstants).**

In my case I can collect pi^4 using mathematica and manually multiply the function by BigFloat(pi)^4 or Float128(pi)^4 (or taking the inherited type from a,b,c,d arguments) but I am asking for general case, when in principle pi or/and e could be in various powers and therefore such manual manipulation would be infeasible. I believe that there is a better way to do it and therefore asking this question.

I would also like to add that it is amazing that such things are possible in Julia and thank all autors for this amazing language

## Why don't you use SymPy/Reduce from the start to generate expressions?

I would prefer to use open source tool for generation of formulas, but at the moment it is not really possible:

– Simplification of expressions is unfortunately much better in Mathematica,

– Mathematica has better support for some special functions that sometimes appear in these formulas,

– I know it better and my collaborators use it.

Hopefully latter I will be able to move to SymPy and FORM (https://github.com/HaraldHofstaetter/FORM.jl) but at the moment it is not possible and I would need the results from Mathematica for benchmarks anyway.

I am also looking forward to learn and try Reduce since there’s a great Julia package, but the same problems as above apply.