import math 

def SQ(x):
	return x**2
def funx1(x1,x2,x3):
	return 8*(x2-x1)
def funx2(x1,x2,x3):
	return x1*x3+0.5*x1-(1+0.5)*x2
def funx3(x1,x2,x3):
	return 50*(1-(1+0.5)*x2+0.5*x2*x2)-0.5*x3
def rk4(f, dt, x1, x2, x3,i):
	k1 = dt * f(x1,x2,x3)
	k2 = dt * f(x1 + 0.5*k1,x2 + 0.5*k1,x3 + 0.5*k1)
	k3 = dt * f(x1 + 0.5*k2,x2 + 0.5*k2,x3 + 0.5*k2)
	k4 = dt * f(x1 + k3,x2 + k3,x3 + k3)
	if i==1 :
		return x1 + (k1 + 2 * k2 + 2 * k3 + k4) /6.0
	if i==2 :
		return x2 + (k1 + 2 * k2 + 2 * k3 + k4) /6.0
	if i==3 :
		return x3 + (k1 + 2 * k2 + 2 * k3 + k4) /6.0
def transpose(X):
	return map(list, zip(*X))

def Addvector(X,Y):
	#X=[]
	#Y=[]
	K=[0,0,0]
	for i in range(len(X)):
		K[i]=X[i] + Y[i]
		return K
def Subvector(X,Y):
	#X=[]
	#Y=[]
	K=[0,0,0]
	for i in range(len(X)):
		K[i]=X[i] - Y[i]
	return K
def ScalarProd(X,Y):
	if len(K)==len(L) and len(K)!=0:
            return sum([K[n]*L[n] for n in range(len(K))])
        else:
            return 0
	
def ScalarMul(X , p):
	res =[0,0,0]
	for i in range(len(X)):
		res[i]= X[i]*p
	return res
def MatrixMulti(X,Y):
	res=[[0,0,0],[0,0,0],[0,0,0]]
	for i in range(len(X)):
		for j in range(len(Y[0])):
			for k in range(len(Y)):
				res [i][j] += X[i][k]*Y[k][j]
	return res

def NormalizeVecto(X):
	res = [0,0,0]
	for i in range(len(X)):
		suma += X[i]*X[x]
	for i in range(len(X)):
		res [i] = X[x]/math.sqrt(suma)
	return res
	
def Modvecto(X):
	suma = 0
	for i in range(len(X)):
		suma += X[i]**2
	return math.sqrt(suma)
def SQ(x):
	return x*x

def Orthogonalization(X):
	out = [[0,0,0],[0,0,0],[0,0,0]]
	suma = [0,0,0]
	out = X
	for i in range(len(X)):
		for j in range(i-1):
			suma = Addvector(suma,ScalarMul(out[j],ScalarProd(out[j],X[i])/SQ(Modvecto(out[j]))))
		out[i] = Subvector(X[i],suma)
	for k0 in range(len(X)):
		if math.fabs(out[i][k0])< 1e-10 :
			out[i][k0] = 0
	return out
		

if __name__ == '__main__':
	raw_input('rikitake system')
	dimension = 3
	xc1 = 0.0
	xc2 = 0.0
	xc3 = 0.0	
	N0 = 1e6
	r=8
	g = 50
	f = 0.5
	m = 0.5
	T = 1e4
	dT = 1e-3
        g = 50
	dy = [0.0,0.0,0.0]
	dummy = [0,0,0]
	lyapunov = []
	sumly = [0,0,0]	
	v = [[1,0,0],[0,1,0],[0,0,1]]
	v1 = [[0,0,0],[0,0,0],[0,0,0]]
	unit = [[0,0,0],[0,0,0],[0,0,0]]
	J = [[0,0,0],[0,0,0],[0,0,0]]
	unit = v
	for n in range(int(T/dT)):
		xc1 = rk4(funx1, dT, dy[0],dy[1],dy[2],0)
		xc2 = rk4(funx2, dT, dy[0],dy[1],dy[2],1)	
		xc3 = rk4(funx3, dT, dy[0],dy[1],dy[2],2)
		dy[0] = xc1
		dy[1] = xc2
		dy[2] = xc3
		if n > int(N0):
			J[0][0] = -r
			J[0][1] = r
			J[0][2] = 0
			J[1][0] = xc3+m
			J[1][1] = -1-m
			J[1][2] = xc1
			J[2][0] = -g*(1+m)*xc2+g*2*xc1
			J[2][1] = -g*xc1*(1+m)
			J[2][2] = -f
			J[0] = Addvector(unit[0],ScalarMul(J[0], dT))
			J[1] = Addvector(unit[1],ScalarMul(J[1], dT))
			J[2] = Addvector(unit[2],ScalarMul(J[2], dT))
			v1 = Orthogonalization(MatrixMulti(J,transpose(v)))
			for k in range(3):
				lyapunov[k] = math.log(Modvecto(v1[k]))/dT
				v[k]=NormalizeVecto(v1[k])
				sumly[k] += lyapunov[k]
	print("Time Averaged Lyapunov Exponents.")
	print(sumly/(int(T/dT - N0)))
	#print("l1=" sumly[0]/(int(T/dT - N0)))
	#print("l2=" sumly[1]/(int(T/dT - N0)))
	#print("l3=" sumly[2]/(int(T/dT - N0)))
	#print("ltot=" sumly[0]/(int(T/dT - N0))+sumly[1]/(int(T/dT - N0))+sumly[2]/(int(T/dT - N0)))

Error:
rikitake system
Traceback (most recent call last):
  File "lyapunov.py", line 112, in <module>
    xc1 = rk4(funx1, dT, dy[0],dy[1],dy[2],0)
  File "lyapunov.py", line 12, in rk4
    k1 = dt * f(x1,x2,x3)
  File "lyapunov.py", line 6, in funx1
    return 8*(x2-x1)
TypeError: unsupported operand type(s) for -: 'float' and 'NoneType'

def rk4(f, dt, x1, x2, x3,i):
    print(f'data types f: {type(f)}, dt: {type(dt)}, x1: {type(x1)}, x2: {type(x2)}, x3: {type(x3)}, i: {type(i)}')

(Mar-19-2018, 10:58 AM)Larz60+ Wrote: The best way to solve this is to break the equations into smaller pieces.
Then print out the results.
one of the operators is attempting to act on non numeric data.
You have to isolate where that is happening in order to correct.
print is a wonderful debugging tool, use it.

for example:
(f-strings require python 3.6 or newer)

def rk4(f, dt, x1, x2, x3,i):
    print(f'data types f: {type(f)}, dt: {type(dt)}, x1: {type(x1)}, x2: {type(x2)}, x3: {type(x3)}, i: {type(i)}')

def rk4(f, dt, x1, x2, x3,i):
    print('data types f: {}, dt: {}, x1: {}, x2: {}, x3: {}, i: {}'.format(type(f), type(dt), type(x1), type(x2), type(x3), type(i)))

>>> def rk4(f, dt, x1, x2, x3,i):
...     print('data types f: {}, dt: {}, x1: {}, x2: {}, x3: {}, i: {}'.format(type(f), type(dt), type(x1), type(x2), type(x3), type(i)))
...
>>> z = [3, 4]
>>> def zingo():
...     print('zingo')
...
>>> y = (2, 3)
>>>
>>> rk4('tz', 0, 123.45, z, y, zingo)
data types f: <class 'str'>, dt: <class 'int'>, x1: <class 'float'>, x2: <class 'list'>, x3: <class 'tuple'>, i: <class 'function'>
>>>

Error:
rikitake system
method test
method test
Traceback (most recent call last):
  File "lyapunov.py", line 113, in <module>
    xc1 = rk4(fun1 , dT, dy[0],dy[1],dy[2],0)
  File "lyapunov.py", line 12, in rk4
    k1 = dt * f(x1,x2,x3)
  File "lyapunov.py", line 6, in fun1
    return 8*x2 - 8*x1
TypeError: unsupported operand type(s) for *: 'int' and 'NoneType'

return 8*x2 - 8*x1