``````

In [4]:

#script to plot the koulakov dynamics
import numpy as np
import matplotlib.pyplot as plt

def plotplane_nonlin_loop(X,Y,fx,fy,showplot):
#two distinct real eigenvalues. stable node

[nsx,nsy] = X.shape
U = np.zeros(X.shape)
V = np.zeros(X.shape)
for m in range(nsx):
for n in range(nsy):
U[m,n] = fx(X[m,n],Y[m,n])
V[m,n] = fy(X[m,n],Y[m,n])
#M = np.array([[-2,0],[1,-4]]) #stable

speed = np.sqrt(U*U + V*V)

fig0, ax0 = plt.subplots()
strm = ax0.streamplot(X, Y, U, V, color=speed, linewidth=2)
fig0.colorbar(strm.lines)

if showplot:
plt.show()
return fig0,ax0

``````
``````

In [2]:

def koulakov_line(S1,S2,I):
import math
import numpy as np
#koulakov_line(X,Y,L,n,a,dt,I_K)
#function to calculate the dynamics for a koulakov line

#S1 = S[0]
#S2 = S[1]
dt = 1.0e-3 #multiplicative factor for gradient
pi = math.pi

#parameters for number of stable points, length of line, etc
L = 0.7 #length of stable line in state space
n = 7.0 #number of stable points
theta = pi/4.0 #the rotation angle of line from horizontal in state space (make it diagonal
a = 0.2; #amplitude of oscillations in state space

#the rotation of the line
A = np.array([[math.cos(theta),-math.sin(theta)],
[math.sin(theta),math.cos(theta)]]) #the rotation matrix
R0 = np.array([L/2.0,1.0]) #point around which to roate
F = A.dot([S1,S2]-R0)+R0
#F = mvMul(A,(S-R0))+R0
S1r = F[0] #rotated S1
S2r = F[1] #rotated S2

#the main functions of the koulakov line, in rotated coordinates
u = a*math.sin((n-1)*pi/L*S1r+pi)
f = math.tanh( (pi/a)*(u-2.0*(S2r-0.5)))
g = math.tanh( (pi/a)*(-u-2.0*(S2r-0.5)))

#boundary conditions. don't use rotated coordinates
steep = 3.0
ramprange = 4.0*(n-1); #support of which the tanh will ramp up on bc in x direction
offset = (steep-1)/(ramprange)

bcS1 = (10.0*math.tanh(-pi/steep*ramprange*(S1+offset)) -
10.0*math.tanh(pi/steep*ramprange*(S1-0.65-offset)))
bcS2 = (10.0*math.tanh(- 2.0*pi/(steep*a)*(S2-L-offset))+
10.0*math.tanh( -2.0*pi/(steep*a)*(S2+offset)))

K1 = dt*(f+bcS1) + I[0]
K2 = dt*(g+bcS2) + I[1]

return K1,K2

``````
``````

In [3]:

Y, X = np.mgrid[0:1:500j, 0:1:500j]

ifact = 7.5e-4 #scaling factor, found by inspection
coh = -7.0
I_k = [0,0]
I_kbase = ifact*(1+abs(coh)/100.0);
if coh < 0:
I_k = [-I_kbase,I_kbase]
elif coh > 0:
I_k = [I_kbase,-I_kbase]

fx = lambda x,y: koulakov_line(x,y,I_k)[0]
fy = lambda x,y: koulakov_line(x,y,I_k)[1]
fplane = plotplane_nonlin_loop(X,Y,fx,fy,True)

``````
``````

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-3-074b12b34a41> in <module>()
12 fx = lambda x,y: koulakov_line(x,y,I_k)[0]
13 fy = lambda x,y: koulakov_line(x,y,I_k)[1]
---> 14 fplane = plotplane_nonlin_loop(X,Y,fx,fy,True)

<ipython-input-1-af93cc2e4195> in plotplane_nonlin_loop(X, Y, fx, fy, showplot)
18
19     fig0, ax0 = plt.subplots()
---> 20     strm = ax0.streamplot(X, Y, U, V, color=speed, linewidth=2, cmap=jet)
21     fig0.colorbar(strm.lines)
22

NameError: global name 'jet' is not defined

``````
``````

In [5]:

fplane = plotplane_nonlin_loop(X,Y,fx,fy,True)

``````
``````

``````
``````

In [ ]:

``````