# Animation Nodes: the Si(x) function

I'm searching how to do the Si(x) function with animation nodes.

It is not the sinus function but the indefinite integral of sinx/x which can't be written as the combinations of other elementary functions.

• Python expression calculating an approximation? Commented Oct 28, 2020 at 14:51
• Hello, I feel like it would be more appropriate to post this on the Maths SE site rather than here, since it doesn't really involve using Animation Nodes but rather how to dissect an integral equation. Basically if you can write the equation in pseudo code or mathematical language, then translate it to python (again, arguably off-topic), then it should be as easy as using an expression node afterwards or a script node if the calculations are more involved. Commented Oct 28, 2020 at 14:53
• Please show what you have tried in order to accomplish what you are trying to do. Where exactly are you stuck in the process ? Commented Oct 28, 2020 at 14:55
• to Gorgious. I don't have a single knowledge about python. Commented Oct 28, 2020 at 14:57
• @lazare, I think we could (perhaps) help if you give the mathematical expression of si(x) approximation. Commented Oct 28, 2020 at 15:34

# Sine Integral

The function, denoted $$Si$$, by definition is the antiderivative of $$\frac{sin\ x}{x}$$ from $$0$$ to the input point whose value is zero at $$x = 0$$. The sine integral definitions are

$$Si(x) = \int_0^x{\frac{sin\ t}{t}}\ dt$$

To put it simply, we want to compute integral with animation nodes as the $$sin$$ function is already available. One technique we can use is Numerical Integration, which is a technique for approximating a definite integral and also quite easy to implement inside animation nodes. There are various methods in numerical integration but for this example I will use two methods that won't use too much node.

## Trapezoidal rule

The Trapezoidal rule works by approximating the region under graph of the given function as trapezium rather than rectangle. Let $$f(x)$$ be continuous on $$[a, b]$$. We partition the interval $$[a, b]$$ into $$N$$ subintervals, each of width

$$\Delta x = \frac{b-a}{N}$$

The Trapezoidal rule for approximating $$\int_a^b f(x)$$ is given by

$$\int_a^b f(x)\ dx \approx \Delta x \Biggl(\sum_{j=1}^{N-1} f(x_j) + \frac{f(x_N) + f(x_0)}{2}\Biggr)$$

Where $$x_j = a + j\Delta x$$. Because $$Si(0) = 0$$, we can simplify the equation into the following

$$\frac{b}{N} \Biggl(\sum_{j=1}^{N-1} f(x_j) + \frac{f(x_N)}{2}\Biggr)$$

Notice that $$\Delta x$$ is also substituted by $$\frac{b}{N}$$ because $$a$$ is always zero in this case.

### Implementation

The result for $$Si(2\pi)$$ with 2000 subintervals is $$1.41804$$. Let's compare with the result from WolframAlpha

To get more accurate approximation, try increase the subintervals.

## Monte Carlo integration

Another technique for numerical integration is Monte Carlo integration. This technique is using random numbers instead of subintervals or midpoint. While Monte Carlo integration is usually used for higher-dimensional integrals (eg. in Machine Learning), we can also use it for one-dimensional integral. Beside that, we can implement this technique in animation nodes with just a few nodes.

$$\int_a^b f(x)\ dx \approx \frac{1}{N} \sum_{j=1}^N f(\bar x_j)$$

where $$\bar x_j$$ is $$N$$ uniform random numbers on $$[a,b]$$. Monte Carlo integration use uniform sampling and fortunately random number generator in animation nodes use uniform sampling. Maybe… Or is it?.

### Implementation

Blend file for study and practice

## Python

### Scipy

Both techniques is not really accurate and a bit time consuming to make the node tree, another alternative is to use scipy. If your blender installation is using system-wide python and have scipy installed, you can use scipy.special.sici function or you can install scipy directly in blender.

from scipy.special import sici
# sici return two results, first is Si(x), second is Ci(x)
Result, _ = sici(Number) # Store in output variable


As you can see, using scipy is more accurate.

### Using Hypergeometric Series

If you don't have or don't want to install scipy, this function will do without any external dependency and also quite accurate. Adapted from mpmath then simplified and this for epsilon value.

EPS = abs(7.0/3 - 4.0/3 - 1)  # 2.2204460492503131e-16 on 64bit CPU

def si_(z, *, maxterms: int = 6000):
"""Sine Integral"""
oz = z
z = -0.25*z*z
s = t = 1.0
k = 0
while 1:
t *= 0.5 + k  # coeffs[0] + k
t /= 1.5 + k  # coeffs[1] + k
t /= 1.5 + k  # coeffs[2] + k
k += 1
t /= k
t *= z
s += t
if abs(t) < EPS:
break
if k > maxterms:
raise Exception("Hypergeometric series converges too slowly. Try increasing maxterms.")
return oz*s

Result = si_(Number)

• You may be have right. I don't have enought skill in python to say it. I'll train my self. Commented Nov 4, 2020 at 22:02
• I tried your tutorial for installing scipy. I enter the text that they told to enter and nothing happen. Commented Nov 8, 2020 at 22:28
• I succed at install scipy and sici fonction work like you said thank you ! Commented Nov 10, 2020 at 16:37