Python Calculating Spin-Loss of a Particle Using Integral Form

[amjad-sh]

TL;DR Summary

I'm trying to integrate an imaginary function by using the package quad.
The function is a complex function and also a very complicated one.
This typing error is showing up to me after running the program: """unsupported operand type(s) for *: 'complex' and 'function' """".
I will show you the details below.

The integral has the form:
$$\frac{s^2\nu^4}{(2\pi)^2}\int_{-1}^1 u(1-u^2)k_f^5[|r_1\chi_1|^2+|r_1\chi_2|^2-|r_1|^2\chi_1^*\chi_2\cos(2k_f\sqrt{u^2-\nu^2}a)-|r_1|^2\chi_2^*\chi_1cos(2k_f\sqrt{u^2-\nu^2}a)]\, du$$
$r_1,\chi_1$ and $\chi_2$ are also imaginary functions of u, because the form of these imaginary functions is so complicated, I also represented them interms of another imaginary functions like a,b,c,m1,m2 and f.
All the functions $r_1,r_2,\chi_1,\chi_2,b,a$ and $c$ are functions of u.
So I defined them as functions with one parameter u as below:

import numpy as np
from scipy.integrate import quad
import scipy
import matplotlib.pyplot as plt
import cmath as math

spinloss=[]
lnu=[]
a = 1.6*10**(-6)
nu = -0.01 \\defining some constants
kf = 4.2*1.6*10**(-19)

for i in range(0,101):
spinloss.append([])
for i in range(0,101):
lnu.append([])

def b(u):
q=u**2-nu**2
if q<0: \\defining the function b
q = -q
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u**2-nu**2)
return(kf*(u+v)*math.exp(-2j*kf*v*a)+kf*(u-v)*math.exp(2j*kf*v*a))def c(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q
v = complex(0, math.sqrt(q)) \\defining the function c
else:
v = math.sqrt(u ** 2 - nu ** 2)
return (kf**2*(u+v)**2*np.exp(-2j*kf*v*a)-kf**2*(u-v)**2*math.exp(2j*kf*v*a))

def a(u):
q = u ** 2 - nu ** 2
if q < 0:\\defining the function a
q = -q
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u ** 2 - nu ** 2)
return (((kf*b(u)*(u-v)+c(u))*math.exp(1j*kf*v*a)-2*kf**2*u*(u+v)*math.exp(-1j*kf*v*a))/((kf*b(u)*(u+v)+c(u))*math.exp(-1j*kf*v*a)-2*kf**2*u*(u-v)*math.exp(1j*kf*v*a)))def x1(u):
\\defining x1
return((a(u)*b(u)+2*kf*u)/c(u))

def x2(u):
\\defining x2
return((2*a(u)*kf*u+b(u))/c(u))def m1(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q
v = complex(0, math.sqrt(q)) \\defining m1
else:
v = math.sqrt(u ** 2 - nu ** 2)
return(kf*v*(1-a(u)))def m2(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q \\defining m2
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u ** 2 - nu ** 2)
return((-kf*v*(1-a(u))))def f(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q
v = complex(0, math.sqrt(q)) \\defining f
else:
v = math.sqrt(u ** 2 - nu ** 2)
return((a(u)+1)*math.exp(-1j*math.sqrt(kf*a))*math.cos(kf*v*a))def r1(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q \\defining r1
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u ** 2 - nu ** 2)
return((2*math.exp(-2j*kf*u*a))/(2*f(u)-((1/(kf*u))*(m1(u)*math.exp(1j*(kf*v-kf*u)*a)-m2(u)*math.exp(-1j*(kf*v+kf*u)*a))))) \\the error started from this linedef r0(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q \\defining r0
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u ** 2 - nu ** 2)
return(math.exp(-2j*kf*u*a)-(r1(u)/(kf*u))*(2*kf*v*math.exp(-1j*kf*(v+u)*a)-(2*kf*v*math.exp(1j*kf*(v-u)*a))))def t0(u):
q = u ** 2 - nu ** 2
if q < 0:
q = -q \\defining t0
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u ** 2 - nu ** 2)
return(((r1(u))/(u))*v*a(u)*(math.exp(-1j*kf*(v+u)*a)-math.exp(1j*kf*(v-u)*a)))
def ls(u):
q = u ** 2 - nu ** 2
if q < 0: \\defining ls
q = -q
v = complex(0, math.sqrt(q))
else:
v = math.sqrt(u ** 2 - nu ** 2)
return(u*(1-u**2)*kf**5*((abs(r1(u)*x1(u)))**2+(abs(r1(u)*x2(u)))**2-(abs(r1(u)))**2)*np.conjugate(x1(u))*x2(u)*math.cos(2*kf*v*a)-(abs(r1(u)))**2*np.conjugate(x2(u))*x1(u)*math.cos(2*kf*v*a))\\ it showed also here.(the error)
def jx(u):
return((((-kf)**2)/(8*(math.pi)**2))*u*(abs(r0(u)))**2) \\defining jx
def complex_quadrature(func, a, b):
def real_func(u):
return scipy.real(func(u))
\\for integration
real_integral = quad(real_func, a, b)

return (real_integral[0])

for i in range(0, 101):
nu = nu + 0.01
if nu == 1: \\ nu represents $\nu$ and it must be between zero and 1.
spinloss = 0
lnu = nu

if nu != 1:
spinloss = ((complex_quadrature(ls, -1, 1)) * (nu ** 4) /(2*math.pi)**2) / (complex_quadrature(jx, -1, 1)) \\ this gives us the result of the integration of ls from -1 to 1 divided by the integration of jx from -1 to 1.
lnu = nu
plt.plot(lnu,spinloss,color='b')
plt.show()

If some python expert can point out for me the reason why this error is appearing, I will be so appreciated.

 

[mfb]

It would help to know which multiplication exactly causes the error. The debugger might show that, alternatively simplify this equation step by step until you find a minimal expression that causes the error:

return((2*math.exp(-2j*kf*u*a))/(2*f(u)-((1/(kf*u))*(m1(u)*math.exp(1j*(kf*v-kf*u)*a)-m2(u)*math.exp(-1j*(kf*v+kf*u)*a))))) \\the error started from this line