Sympy - speed up expressions evaluation

Question

I compute something in python with sympy but this require too much time. The code looks like this:

def compute_aux_lyupanov_val_k(y, A, B, C, D, k):
  res = []
  res.append(y)
  for i in range(k-1):  
    y = (C + D * y) / (A + B * y)
    res.append(y)
  return res


def compute_lyupanov_val(y, A, B, C, D, n):
  res = 0
  y_values = compute_aux_lyupanov_val_k(y, A, B, C, D, n) # values for y'_0. y'_1...

  for i in range(n):
    T = ( (A + B * y_values[i])**2 + (C + D * y_values[i])**2  ) / (1 + y_values[i]**2)
    res = res + math.log(T, 10) # default is base e = 2.71

  res = res * (1/(2*n))
  return res

n = 1000
y_inital = 0.01

beta_set_values = np.linspace(0, 6, 1000)
alpha_set_values = np.linspace(0, 5, 1000)
points = it.product(beta_set_values, alpha_set_values)

for point in points:
  alpha_val = point[1]
  beta_val = point[0]

  A_val = A.subs([  (beta, beta_val), (alpha, alpha_val) ])
  B_val = B.subs([  (beta, beta_val), (alpha, alpha_val) ])
  C_val = C.subs([  (beta, beta_val), (alpha, alpha_val) ])
  D_val = D.subs([  (beta, beta_val), (alpha, alpha_val) ])

  lyupanov_val = compute_lyupanov_val(y_inital, A_val, B_val, C_val, D_val, n)

  if lyupanov_val > 0:
    plt.scatter([beta_val], [alpha_val], color = 'r', s=1)


plt.show()

Maybe it seems long, but let me explain clearly: A, B, C, D are expression which can be evaluated depending on alpha, beta values. For a set of points (alphas and betas) I compute a value (lyupanov_val) and scatter the (beta, alpha) point only if computed value lyupanov_val is greater than 0. The computations which happened in auxiliary functions (compute_aux_lyupanov_val_k and compute_lyupanov_val) are only numerical (A, B, C, D and other parameters are scalars now - inside of these functions, not symbols)

I think the expensive computation part happen because of the substitution at every iteration. I know that sympy provide a special plot function which can speed up the computations sp.plot_implicit(expresion > 0, (beta, 0, 6), (alpha, 0, 5) but the problem with this approach is the fact that I don't use a direct expression (e.g alpha + beta / 2), I use auxiliary functions for numerical computations (e.g compute_aux_lyupanov_val_k)

If I used

sp.plot_implicit(compute_lyupanov_val(y_inital, A, B, C, D, n) > 0, (beta, 0, 6), (alpha, 0, 5))

then the python raise TypeError: can't convert expression to float

Any ideas?

Answer 1

subs does a symbolic substitution but the result is still a symbolic expression. May be what you want is instead

A_val = A.subs([  (beta, beta_val), (alpha, alpha_val) ]).evalf()

and so on for the 4 variables.

Even better would be to use

A_val = A.evalf(subs={beta: beta_val, alpha: alpha_val})

(subs does a replacement and a symbolic simplification; if what you want is instead just the numeric value a direct computation is better).