# Source code for cromosim.ftl

# Authors:
#     Sylvain Faure <sylvain.faure@universite-paris-saclay.fr>
#     Bertrand Maury <bertrand.maury@universite-paris-saclay.fr>
import numpy as np
import scipy as sp
import sys
import random
import matplotlib
import matplotlib.pyplot as plt

from matplotlib.lines import Line2D
from matplotlib.collections import LineCollection

[docs]def update_positions_ftl_order_1(L, X, t, dt, Phi, periodicity=True,
"""
To update the positions according to the follow the leader model (order 1)

Parameters
----------
L: float
width of the domain
X: numpy array
coordinates of the individuals
t: float
time
dt: float
time step
Phi: function
speed as a function of the distance
periodicity: boolean
if true the domain is periodic
leader velocity as a function of the time

Returns
-------
Xnew: numpy array
new coordinates of the individuals
V: numpy array
new velocities of the individuals

"""
if (periodicity == True):
S = np.argsort(X)
W = X[S[1:]] - X[S[:-1]]
W = np.append(W,(X[S[0]]+L) - X[S[-1]])
V = Phi(W)
Xnew = X[S] + dt*V
Xnew1 = Xnew*(Xnew<L) + (Xnew-L)*(Xnew>=L)
Xnew[S] = Xnew1
V1 = V ; V[S] = V1
else:
# The last position in the array is the leader...
W = X[1:]-X[:-1]
VV = Phi(W)
V = np.append(VV,VL)
Xnew = X+dt*V
return Xnew, V

[docs]def update_positions_ftl_order_2(L, X, U, t, dt, tau, Phi, periodicity=True,
"""
To update the positions according to the follow the leader model (order 2)

Parameters
----------
L: float
width of the domain
X: numpy array
coordinates of the individuals
U: numpy array
velocities of the individuals
t: float
time
dt: float
time step
tau: float
the actual velocity of an individual relaxes toward the velocity that \
is associated to the current distance from the person in front of them \
with a characteristic time tau
Phi: function
speed as a function of the distance
periodicity: boolean
if true the domain is periodic
leader velocity as a function of the time

Returns
-------
Xnew: numpy array
new coordinates of the individuals
Unew: numpy array
new velocities of the individuals

"""
if (periodicity == True):
S = np.argsort(X)
W = X[S[1:]]-X[S[:-1]]
W = np.append(W,(X[S[0]]+L)-X[S[-1]])
V = Phi(W)
Unew = U[S]+dt*(V-U[S])/tau
Xnew = X[S]+dt*Unew
Xnew1 = Xnew*(Xnew<L)+(Xnew-L)*(Xnew>=L)
Xnew[S] = Xnew1
V1 = V ; V[S] = V1
Unew1 = Unew ; Unew[S] = Unew
else:
# The last position in the array is the leader...
Xnew = X+dt*U
W = X[1:]-X[:-1]
Unew = np.zeros(U.shape)
Unew[:-1] = U[:-1] + dt*(Phi(W)-U[:-1])/tau
return Xnew, Unew

[docs]def plot_people(L, X, t, data, tgrid, periodicity=True, V_leader=None,
speed_fct=None, savefig=False, filename='fig.png',
ifig=100, spheresize=300):
"""
To plot people and their individual paths

Parameters
----------
L: float
width of the domain
X: numpy array
positions of the individuals
data: numpy array
coordinates of indviduals at each time
tgrid: numpy array
temporal discretization
periodicity: boolean
if true the domain is periodic
leader velocity as a function of the time
speed_fct: function
speed as a function of the distance
savefig: boolean
writes the figure as a png file if true
filename: string
png filename used to write the figure
ifig: int
figure number
spheresize: int
size of the spheres used to represent the individuals
"""
N = X.shape[0]
XY = np.zeros((N,2))
XY[:,0] = X
itmax = np.where(tgrid<=t)[0].max()

fig = plt.figure(ifig,figsize=(16,9))
fig.clf()

if periodicity:
else:
ax1.set_xlim(0,L)
ax1.set_ylim(-0.5,0.5)
#ax1.set_xticks([])
ax1.set_yticks([])
#ax1.axis('off')
ax1.set_xlabel('x')
offsets = list(zip(XY))
colors = np.arange(N)/N
colors = np.ones(N)#np.arange(N)/N
scale = 20
sc = ax1.scatter(XY[:,0],XY[:,1], c=colors, s=spheresize, edgecolor='None',
marker='.',cmap=plt.get_cmap('Greys'),norm=plt.Normalize(-0.2, 1))
line = Line2D([ 0, 0 ],[ -0.05, 0.05 ], lw=3., alpha=1, c='k')
line = Line2D([ L, L ],[ -0.05, 0.05 ], lw=3., alpha=1, c='k')

if (periodicity==False):
ax3.set_xlabel('time')

wgrid = np.linspace(0,L,500)
if periodicity:
else:
ax2.plot(wgrid,speed_fct(wgrid))
S = np.argsort(X)
W = X[S[1:]] - X[S[:-1]]
if periodicity:
W = np.append(W,(X[S[0]]+L) - X[S[-1]])
ax2.plot(W,speed_fct(W),'ko')
ax2.set_xlabel('distance')
ax2.set_ylabel('speed')

if periodicity:
else:
for ip in np.arange(N):
ax4.plot(tgrid[:itmax], data[ip,0,:itmax])
ax4.set_ylim(0,max(L,data[:,0,:itmax].max()))
ax4.set_xlim(tgrid.min(),tgrid.max())
ax4.set_xlabel('time')
ax4.set_ylabel('x')

if periodicity:
else:
Y = np.linspace(0, 0.25*L*N, num=N)
for ip in np.arange(N-1):
ax5.plot(tgrid[:itmax],Y[ip]+(data[ip+1,0,:itmax]-data[ip,0,:itmax]))
if periodicity:
d0 = data[0,0,:itmax]
dN = data[N-1,0,:itmax]
ax5.plot(tgrid[:itmax],Y[-1]+(d0<dN)*(d0+L-dN)+(d0>=dN)*(d0-dN))
ax5.set_yticks([])
#ax5.set_ylim(0,1.1*Y.max())
ax5.set_xlim(tgrid.min(),tgrid.max())
ax5.set_xlabel('time')
ax5.set_ylabel('distances')

fig.set_tight_layout(True)

if (savefig):
fig.savefig(filename,dpi=150)