Julia – 15 – stringhe – 1

Continuo da qui, copio qui.

Strings are finite sequences of characters. Of course, the real trouble comes when one asks what a character is. The characters that English speakers are familiar with are the letters A, B, C, etc., together with numerals and common punctuation symbols. These characters are standardized together with a mapping to integer values between 0 and 127 by the ASCII standard. There are, of course, many other characters used in non-English languages, including variants of the ASCII characters with accents and other modifications, related scripts such as Cyrillic and Greek, and scripts completely unrelated to ASCII and English, including Arabic, Chinese, Hebrew, Hindi, Japanese, and Korean. The Unicode standard tackles the complexities of what exactly a character is, and is generally accepted as the definitive standard addressing this problem. Depending on your needs, you can either ignore these complexities entirely and just pretend that only ASCII characters exist, or you can write code that can handle any of the characters or encodings that one may encounter when handling non-ASCII text. Julia makes dealing with plain ASCII text simple and efficient, and handling Unicode is as simple and efficient as possible. In particular, you can write C-style string code to process ASCII strings, and they will work as expected, both in terms of performance and semantics. If such code encounters non-ASCII text, it will gracefully fail with a clear error message, rather than silently introducing corrupt results. When this happens, modifying the code to handle non-ASCII data is straightforward.

There are a few noteworthy high-level features about Julia’s strings:

  • The built-in concrete type used for strings (and string literals) in Julia is String. This supports the full range of Unicode characters via the UTF-8 encoding. (A transcode() function is provided to convert to/from other Unicode encodings.)
  • All string types are subtypes of the abstract type AbstractString, and external packages define additional AbstractString subtypes (e.g. for other encodings). If you define a function expecting a string argument, you should declare the type as AbstractString in order to accept any string type.
  • Like C and Java, but unlike most dynamic languages, Julia has a first-class type representing a single character, called Char. This is just a special kind of 32-bit primitive type whose numeric value represents a Unicode code point.
  • As in Java, strings are immutable: the value of an AbstractString object cannot be changed. To construct a different string value, you construct a new string from parts of other strings.
  • Conceptually, a string is a partial function from indices to characters: for some index values, no character value is returned, and instead an exception is thrown. This allows for efficient indexing into strings by the byte index of an encoded representation rather than by a character index, which cannot be implemented both efficiently and simply for variable-width encodings of Unicode strings.

A Char value represents a single character: it is just a 32-bit primitive type with a special literal representation and appropriate arithmetic behaviors, whose numeric value is interpreted as a Unicode code point. Here is how Char values are input and shown:

You can convert a Char to its integer value, i.e. code point, easily:

On 32-bit architectures, typeof(ans) will be Int32. You can convert an integer value back to a Char just as easily:

Not all integer values are valid Unicode code points, but for performance, the Char() conversion does not check that every character value is valid. If you want to check that each converted value is a valid code point, use the isvalid() function:

As of this writing, the valid Unicode code points are U+00 through U+d7ff and U+e000 through U+10ffff. These have not all been assigned intelligible meanings yet, nor are they necessarily interpretable by applications, but all of these values are considered to be valid Unicode characters.

You can input any Unicode character in single quotes using \u followed by up to four hexadecimal digits or \U followed by up to eight hexadecimal digits (the longest valid value only requires six):

Julia uses your system’s locale and language settings to determine which characters can be printed as-is and which must be output using the generic, escaped \u or \U input forms. In addition to these Unicode escape forms, all of C’s traditional escaped input forms can also be used:

You can do comparisons and a limited amount of arithmetic with Char values:


SICP – cap. 2 – Sequenze come interfacce convenzionali – 55 – esercizi

Continuo da qui, copio qui.

Exercise 2.44: Define the procedure up-split used by corner-split. It is similar to right-split, except that it switches the roles of below and beside.

Da vedersi in stretto collegamento con lo spiegone del post precedente 😊
Bill the Lizard 🚀 al solito spiega bene, lo seguo 😊 Inoltre, contrariamente alla mia predilizione per il (vecchio) terminale, per via della grafica uso Dr.Racket, quando ci vuole ci vuole 😉

I’m going to be using the PLT Scheme SICP Picture Language package to run the examples and exercises. You can load the picture package by putting the following (require...) expression at the beginning of your Scheme file.

The painter primitive in the PLT package behaves slightly differently than the one presented in the text, in that it doesn’t display itself. To display a painter you must call the paint procedure, which takes a painter as its argument.

The einstein [... is] built in to the PLT Scheme Picture Language package. The package also includes procedures for flipping, rotating and combining painter values, similar to those discussed in the lecture.

Now that we can work with painters, we can start to combine them in ways that form patterns. For example, the text shows us how we can combine an image beside a flipped representation of itself, then draw the resulting painter below itself.

We can also define recursive operations, like right-split, which makes painters split and branch towards the right as many times as we specify.

E a questo punto si arriva all’esercizio

Once we have both right-split and up-split defined, we can create balanced patterns by branching in both directions.

By combining four copies of the corner-split pattern, we can create the square-limit pattern from the Escher drawing that we saw in the lecture.

Secondo me Bill non solo rockz ma è anche oltre 👽

La soluzione di sicp-ex è molto simile a quella di Bill; come pure Drewiki.


SciPy – 31 – algebra lineare – 3

Continuo da qui, copio qui.

Risoluzione con il metodo dei minimi quadtrati e pseudo inversi
Linear least-squares problems occur in many branches of applied mathematics. In this problem a set of linear scaling coefficients is sought that allow a model to fit data. In particular it is assumed that data yi is related to data xi through a set of coefficients cj and model functions fj(xi) via the model

where ϵi represents uncertainty in the data. The strategy of least squares is to pick the coefficients cj to minimize

Theoretically, a global minimum will occur when



When AHA is invertible, then

where A is called the pseudo-inverse of A. Notice that using this definition of A the model can be written

The command linalg.lstsq will solve the linear least squares problem for c given A and y. In addition linalg.pinv or linalg.pinv2 (uses a different method based on singular value decomposition) will find A given A.

The following example and figure demonstrate the use of linalg.lstsq and linalg.pinv for solving a data-fitting problem. The data shown below were generated using the model:

where xi=0.1i for i=1...10, c1=5, and c2=4. Noise is added to yi and the coefficients c1 and c2 are estimated using linear least squares.


cit. & loll – 53

🚀 48 anni fa –OK da festeggiare 🚀 intanto ecco qualcosa di completamente differente 😜

yes it still works while vacuum packed
::: whitequark

Remember, only you can prevent OOP
::: reverentgeek

The computer is a very alluring machine
::: CodeWisdom

Another matplotlib optical illusion added to my collection
::: jakevdp

I like to read AI as Al
::: johnregehr ::: johnregehr

Software is only as shitty
::: rightfold

Google’s DeepMind AI just taught itself to walk
::: SarahJReports

Fools ignore complexity
::: CodeWisdom

Choice of programming language is still one of the most important project germination decisions
::: spakhm

Se vai a Milano e urli in piazza “HEY SOCIAL MEDIA DIGITAL STRATEGIST”
::: Tremenoventi

Over and over again
::: Nonicoc

Quello che ha inventato i programmi user friendly
::: Facebook

A language that doesn’t have everything is actually easier to program in than some that do
::: CodeWisdom

Moving out of my old 4th floor office
::: mattmight

Trust nothing, stay frosty
::: amritamaz

Enough with complaining about Windows. You got reliable alternatives
::: Donearm

La tastiera dell’indignato
::: Beppe Beppetti

Intel was born nearly 50 yrs ago #otd, & was almost called “Moore Noyce,” a homophone for “more noise”

the cycle: some langs are good → most langs are bad →
::: lexi_lambda

L’IA conquisterà il mondo
::: gignico

Wifi vs Cellular
::: xkcd

The only way to learn a new programming language
::: CodeWisdom

Write code as if you had to support it for the rest of your life
::: CodeWisdom

Julia – 14 – numeri complessi e razionali – 2

Continuo da qui, copio qui.

Numeri razionali
Julia has a rational number type to represent exact ratios of integers. Forse non tutti sanni che Julia viene dalla prestigiosa amiglia del Lisp 😁 Rationals are constructed using the // operator:

If the numerator and denominator of a rational have common factors, they are reduced to lowest terms such that the denominator is non-negative:

Si possono ovviamente usare spazi e variabili

This normalized form for a ratio of integers is unique, so equality of rational values can be tested by checking for equality of the numerator and denominator. The standardized numerator and denominator of a rational value can be extracted using the numerator() and denominator() functions:

Direct comparison of the numerator and denominator is generally not necessary, since the standard arithmetic and comparison operations are defined for rational values:

Nota per me: gli spazi nel razionale sono una pessima idea.

Rationals can be easily converted to floating-point numbers:

Conversion from rational to floating-point respects the following identity for any integral values of a and b, with the exception of the case a == 0 and b == 0:

Constructing infinite rational values is acceptable:

Trying to construct a NaN rational value, however, is not:

As usual, the promotion system makes interactions with other numeric types effortless:

Non so voi ma per me Julia rockz! 🚀


SciPy – 30 – algebra lineare – 2

Continuo da qui, copio qui.

Risoluzione di sistemi lineari
Solving linear systems of equations is straightforward using the scipy command linalg.solve. This command expects an input matrix and a right-hand-side vector. The solution vector is then computed. An option for entering a symmetric matrix is offered which can speed up the processing when applicable. As an example, suppose it is desired to solve the following simultaneous equations:

We could find the solution vector using a matrix inverse:

However, it is better to use the linalg.solve command which can be faster and more numerically stable. In this case it however gives the same answer as shown in the following example:

Trovare il determinante
The determinant of a square matrix A is often denoted |A| and is a quantity often used in linear algebra. Suppose aij are the elements of the matrix A and let Mij=|Aij| be the determinant of the matrix left by removing the ith row and jth column from A. Then for any row i,

This is a recursive way to define the determinant where the base case is defined by accepting that the determinant of a 1×1 matrix is the only matrix element. In SciPy the determinant can be calculated with linalg.det. For example, the determinant of


In SciPy this is computed as shown in this example:

Calcolo di vettori e matrici normati
Matrix and vector norms can also be computed with SciPy. A wide range of norm definitions are available using different parameters to the order argument of linalg.norm. This function takes a rank-1 (vectors) or a rank-2 (matrices) array and an optional order argument (default is 2). Based on these inputs a vector or matrix norm of the requested order is computed.

For vector x, the order parameter can be any real number including inf or -inf. The computed norm is

For matrix A the only valid values for norm are ±2, ±1, ±inf, and 'fro' (or 'f') Thus,

where σi are the singular values of A.



Julia – 13 – numeri complessi e razionali – 1

Continuo da qui, copio qui.

Julia ships with predefined types representing both complex and rational numbers, and supports all standard Mathematical Operations and Elementary Functions on them. Conversion and Promotion [prossimamente] are defined so that operations on any combination of predefined numeric types, whether primitive or composite, behave as expected.

Numeri complessi
The global constant im is bound to the complex number i, representing the principal square root of -1. It was deemed harmful to co-opt the name i for a global constant, since it is such a popular index variable name. Since Julia allows numeric literals to be juxtaposed with identifiers as coefficients, this binding suffices to provide convenient syntax for complex numbers, similar to the traditional mathematical notation:

You can perform all the standard arithmetic operations with complex numbers:

The promotion mechanism ensures that combinations of operands of different types just work:

Note that 3/4im == 3/(4*im) == -(3/4*im), since a literal coefficient binds more tightly than division.

Standard functions to manipulate complex values are provided:

As usual, the absolute value (abs()) of a complex number is its distance from zero. abs2() gives the square of the absolute value, and is of particular use for complex numbers where it avoids taking a square root. angle() returns the phase angle in radians (also known as the argument or arg function). The full gamut of other Elementary Functions is also defined for complex numbers:

Note that mathematical functions typically return real values when applied to real numbers and complex values when applied to complex numbers. For example, sqrt() behaves differently when applied to -1 versus -1 + 0im even though -1 == -1 + 0im:

The literal numeric coefficient notation does not work when constructing a complex number from variables. Instead, the multiplication must be explicitly written out:

However, this is not recommended; Use the complex() function instead to construct a complex value directly from its real and imaginary parts:

This construction avoids the multiplication and addition operations.

Inf and NaN propagate through complex numbers in the real and imaginary parts of a complex number:


SciPy – 29 – algebra lineare – 1

Continuo da qui, copio qui.

Nota: salto tutto il capitolo sui segnali, troppo specifico e lontano dai miei interessi

When SciPy is built using the optimized ATLAS LAPACK and BLAS libraries, it has very fast linear algebra capabilities. If you dig deep enough, all of the raw lapack and blas libraries are available for your use for even more speed. In this section, some easier-to-use interfaces to these routines are described.

All of these linear algebra routines expect an object that can be converted into a 2-dimensional array. The output of these routines is also a two-dimensional array.

scipy.linalg vs numpy.linalg
scipy.linalg contains all the functions in numpy.linalg. plus some other more advanced ones not contained in numpy.linalg.

Another advantage of using scipy.linalg over numpy.linalg is that it is always compiled with BLAS/LAPACK support, while for numpy this is optional. Therefore, the scipy version might be faster depending on how numpy was installed.

Therefore, unless you don’t want to add scipy as a dependency to your numpy program, use scipy.linalg instead of numpy.linalg.

numpy.matrix vs 2D numpy.ndarray
The classes that represent matrices, and basic operations such as matrix multiplications and transpose are a part of numpy. For convenience, we summarize the differences between numpy.matrix and numpy.ndarray here.

numpy.matrix is matrix class that has a more convenient interface than numpy.ndarray for matrix operations. This class supports for example MATLAB-like creation syntax via the [non mi è chiaro cosa vuol dire], has matrix multiplication as default for the * operator, and contains I and T members that serve as shortcuts for inverse and transpose:

Despite its convenience, the use of the numpy.matrix class is discouraged, since it adds nothing that cannot be accomplished with 2D numpy.ndarray objects, and may lead to a confusion of which class is being used. For example, the above code can be rewritten as:

scipy.linalg operations can be applied equally to numpy.matrix or to 2D numpy.ndarray objects.

Matrici inverse
The inverse of a matrix A is the matrix B such that AB=I where I is the identity matrix consisting of ones down the main diagonal. Usually B is denoted B=A§−1 . In SciPy, the matrix inverse of the Numpy array, A, is obtained using linalg.inv(A), or using A.I if A is a Matrix. For example, let


The following example demonstrates this computation in SciPy


Julia – 12 – operazioni matematiche e funzioni elementari – 3

Continuo da qui, copio qui.

Precedenze tra operatori
Julia applies the following order of operations, from highest precedence to lowest:

Category       Operators
Syntax	       . followed by ::
Exponentiation ^
Fractions      //
Multiplication * / % & \
Bitshifts      << >> >>>
Addition       + - | ⊻
Syntax	       : .. followed by |>
Comparisons    > < >= <= == === != !== <:
Control flow   && followed by || followed by ?
Assignments    = += -= *= /= //= \= ^= ÷= %= |= &= ⊻= <<= >>= >>>=

For a complete list of every Julia operator’s precedence, see the top of this file.
Uh! 😁 Fantastico il linguaggio con cui è scritto Julia: Scheme, cioè Lisp (una variante del) 🚀

You can also find the numerical precedence for any given operator via the built-in function Base.operator_precedence, where higher numbers take precedence:

Conversioni numeriche
Julia supports three forms of numerical conversion, which differ in their handling of inexact conversions.

  • The notation T(x) or convert(T, x) converts x to a value of type T.
    ° If T is a floating-point type, the result is the nearest representable value, which could be positive or negative infinity.
    ° If T is an integer type, an InexactError is raised if x is not representable by T.
  • x % T converts an integer x to a value of integer type T congruent to x modulo 2^n, where n is the number of bits in T. In other words, the binary representation is truncated to fit.
  • The Rounding functions take a type T as an optional argument. For example, round(Int, x) is a shorthand for Int(round(x)).

The following examples show the different forms.

See Conversion and Promotion [prossimamente] for how to define your own conversions and promotions.

Funzioni di arrotondamento

Function	Description                    Return type
round(x)        round x to the nearest integer typeof(x)
round(T, x)     round x to the nearest integer T
floor(x)        round x towards -Inf           typeof(x)
floor(T, x)     round x towards -Inf           T
ceil(x)	        round x towards +Inf               typeof(x)
ceil(T, x)      round x towards +Inf           T
trunc(x)        round x towards zero           typeof(x)
trunc(T, x)     round x towards zero           T

Funzioni di divisione

Function    Description
div(x,y)    truncated division; quotient rounded towards zero
fld(x,y)    floored division; quotient rounded towards -Inf
cld(x,y)    ceiling division; quotient rounded towards +Inf
rem(x,y)    remainder; satisfies x == div(x,y)*y + rem(x,y); sign matches x
mod(x,y)    modulus; satisfies x == fld(x,y)*y + mod(x,y); sign matches y
mod1(x,y)   mod() with offset 1; returns r∈(0,y] for y>0 or r∈[y,0) for
            y<0, where mod(r, y) == mod(x, y)
mod2pi(x)   modulus with respect to 2pi; 0 <= mod2pi(x)   < 2pi
divrem(x,y) returns (div(x,y),rem(x,y))
fldmod(x,y) returns (fld(x,y),mod(x,y))
gcd(x,y...) greatest positive common divisor of x, y,...
lcm(x,y...) least positive common multiple of x, y,...

Potenze, logaritmi e radici

Function       Description
sqrt(x), √x    square root of x
cbrt(x), ∛x    cube root of x
hypot(x,y)     hypotenuse of right-angled triangle with
               other sides of length x and y
exp(x)         natural exponential function at x
expm1(x)       accurate exp(x)-1 for x near zero
ldexp(x,n)     x*2^n computed efficiently for integer values of n
log(x)         natural logarithm of x
log(b,x)       base b logarithm of x
log2(x)        base 2 logarithm of x
log10(x)       base 10 logarithm of x
log1p(x)       accurate log(1+x) for x near zero
exponent(x)    binary exponent of x
significand(x) binary significand (a.k.a. mantissa) of 
               a floating-point number x

For an overview of why functions like hypot(), expm1(), and log1p() are necessary and useful, see John D. Cook’s excellent pair of blog posts on the subject: expm1, log1p, erfc, and hypot.

Nota: codici Unicode U+221A, U+221B. Io li cerco (se costretto a usarli) da FileFormat.Info.

Funzioni trigonometriche e iperboliche
All the standard trigonometric and hyperbolic functions are also defined:

sin    cos    tan    cot    sec    csc
sinh   cosh   tanh   coth   sech   csch
asin   acos   atan   acot   asec   acsc
asinh  acosh  atanh  acoth  asech  acsch
sinc   cosc   atan2

These are all single-argument functions, with the exception of atan2, which gives the angle in radians between the x-axis and the point specified by its arguments, interpreted as x and y coordinates.

Additionally, sinpi(x) and cospi(x) are provided for more accurate computations of sin(pi*x) and cos(pi*x) respectively.

In order to compute trigonometric functions with degrees instead of radians, suffix the function with d. For example, sind(x) computes the sine of x where x is specified in degrees. The complete list of trigonometric functions with degree variants is:

sind   cosd   tand   cotd   secd   cscd
asind  acosd  atand  acotd  asecd  acscd

Funzioni speciali

Function   Description
gamma(x)   gamma function at x
lgamma(x)  accurate log(gamma(x)) for large x
lfact(x)   accurate log(factorial(x)) 
           for large x; same as lgamma(x+1) 
           for x > 1, zero otherwise
beta(x,y)  beta function at x,y
lbeta(x,y) accurate log(beta(x,y)) for large x or y


SciPy – 28 – trasformate di Fourier – 2

Continuo da qui, copio qui.

Trasformazioni discrete di Fourier bi e n-dimensionali
The functions fft2 and ifft2 provide 2-dimensional FFT, and IFFT, respectively. Similar, fftn and ifftn provide n-dimensional FFT, and IFFT, respectively.

The example below demonstrates a 2-dimensional IFFT and plots the resulting (2-dimensional) time-domain signals. Codice troppo lungo, lo raccolgo nel file fft2.py

import numpy as np
from scipy.fftpack import ifftn
import matplotlib.pyplot as plt
import matplotlib.cm as cm

N = 30
f, ((ax1, ax2, ax3), (ax4, ax5, ax6)) = plt.subplots(2, 3, sharex='col', sharey='row')
xf = np.zeros((N,N))
xf[0, 5] = 1
xf[0, N-5] = 1
Z = ifftn(xf)
ax1.imshow(xf, cmap=cm.Reds)
ax4.imshow(np.real(Z), cmap=cm.gray)
xf = np.zeros((N, N))
xf[5, 0] = 1
xf[N-5, 0] = 1
Z = ifftn(xf)
ax2.imshow(xf, cmap=cm.Reds)
ax5.imshow(np.real(Z), cmap=cm.gray)
xf = np.zeros((N, N))
xf[5, 10] = 1
xf[N-5, N-10] = 1
Z = ifftn(xf)
ax3.imshow(xf, cmap=cm.Reds)
ax6.imshow(np.real(Z), cmap=cm.gray)

Convoluzione FFT

scipy.fftpack.convolve performs a convolution of two one-dimensional arrays in frequency domain.

Trasformazione discreta del coseno (DCT)
Scipy provides a DCT with the function dct and a corresponding IDCT with the function idct. There are 8 types of the DCT; however, only the first 3 types are implemented in scipy. “The” DCT generally refers to DCT type 2, and “the” Inverse DCT generally refers to DCT type 3. In addition, the DCT coefficients can be normalized differently (for most types, scipy provides None and ortho). Two parameters of the dct / idct function calls allow setting the DCT type and coefficient normalization.

For a single dimension array x, dct(x, norm=’ortho’) is equal to MATLAB dct(x).

DCT tipo I
Scipy uses the following definition of the unnormalized DCT-I (norm='None'):

Only None is supported as normalization mode for DCT-I. Note also that the DCT-I is only supported for input size > 1.

DCT tipo II
Scipy uses the following definition of the unnormalized DCT-II (norm='None'):

In case of the normalized DCT (norm='ortho'), the DCT coefficients y[k] are multiplied by a scaling factor f:

In this case, the DCT “base functions”

become orthonormal:

DCT tipo III
Scipy uses the following definition of the unnormalized DCT-III (norm='None'):

or, for norm='ortho':

The (unnormalized) DCT-III is the inverse of the (unnormalized) DCT-II, up to a factor 2N. The orthonormalized DCT-III is exactly the inverse of the orthonormalized DCT- II. The function idct performs the mappings between the DCT and IDCT types.

The example below shows the relation between DCT and IDCT for different types and normalizations.

The DCT exhibits the “energy compaction property”, meaning that for many signals only the first few DCT coefficients have significant magnitude. Zeroing out the other coefficients leads to a small reconstruction error, a fact which is exploited in lossy signal compression (e.g. JPEG compression).

The example below shows a signal x and two reconstructions (x20 and x15)from the signal’s DCT coefficients. The signal x20 is reconstructed from the first 20 DCT coefficients, x15 is reconstructed from the first 15 DCT coefficients. It can be seen that the relative error of using 20 coefficients is still very small (~0.1%), but provides a five-fold compression rate.

import numpy as np
from scipy.fftpack import dct, idct
import matplotlib.pyplot as plt

N = 100
t = np.linspace(0,20,N)
x = np.exp(-t/3)*np.cos(2*t)
y = dct(x, norm='ortho')
window = np.zeros(N)
window[:20] = 1
yr = idct(y*window, norm='ortho')
sum(abs(x-yr)**2) / sum(abs(x)**2)
plt.plot(t, x, '-bx')
plt.plot(t, yr, 'ro')
window = np.zeros(N)
window[:15] = 1
yr = idct(y*window, norm='ortho')
sum(abs(x-yr)**2) / sum(abs(x)**2)
plt.plot(t, yr, 'g+')
plt.legend(['x', '$x_{20}$', '$x_{15}$'])

Trasformazione discreta del seno (DST)

Scipy provides a DST with the function dst and a corresponding IDST with the function idst.

There are theoretically 8 types of the DST for different combinations of even/odd boundary conditions and boundary off sets, only the first 3 types are implemented in scipy.

DST tipo I
DST-I assumes the input is odd around n=-1 and n=N. Scipy uses the following definition of the unnormalized DST-I (norm='None'):

Only None is supported as normalization mode for DST-I. Note also that the DST-I is only supported for input size > 1. The (unnormalized) DST-I is its own inverse, up to a factor 2(N+1).

DST tipo II
DST-II assumes the input is odd around n=-1/2 and even around n=N. Scipy uses the following definition of the unnormalized DST-II (norm='None'):

DST tipo III
DST-III assumes the input is odd around n=-1 and even around n=N-1. Scipy uses the following definition of the unnormalized DST-III (norm='None'):

The example below shows the relation between DST and IDST for different types and normalizations.

Distruzione della cache
To accelerate repeat transforms on arrays of the same shape and dtype, scipy.fftpack keeps a cache of the prime factorization of length of the array and pre-computed trigonometric functions. These caches can be destroyed by calling the appropriate function in scipy.fftpack._fftpack. dst(type=1) and idst(type=1) share a cache (*dst1_cache). As do dst(type=2), dst(type=3), idst(type=3), and idst(type=3) (*dst2_cache).