--- a/lattice2d.jl	Wed Jan 05 21:18:21 2022 +0100
+++ b/lattice2d.jl	Thu Jan 06 20:07:01 2022 +0100
@@ -142,7 +142,7 @@
 u0 = initial_state(n, initial_corner_displacements)
 
 # Essential part: the coupling matrix!
-coupling = coupling_matrix(n, 1.; damping=0.0001)
+coupling = coupling_matrix(n, 1.; damping=0.0000)
 
 # ODE formalities
 problem = ODEProblem(f2d, u0, (0, 40), coupling)

--- a/linear_chain.jl	Wed Jan 05 21:18:21 2022 +0100
+++ b/linear_chain.jl	Thu Jan 06 20:07:01 2022 +0100
@@ -4,7 +4,8 @@
 using Plots
 theme(:juno)
 
-
+# The coupling matrix has 4nx4n elements for a linear chain of two degrees of freedom per oscillator.
+# We use Born-von Karman boundary conditions, i.e. the chain is a ring.
 function n_coupled_masses(n::Int, ks::Vector, masses::Vector)
     # Prepares coupling matrix.
     if length(ks) != (n+1)
@@ -86,8 +87,8 @@
 N = 10
 k1, k2 = 1, 2
 mass = 1
-initial_displacement = [(1, .1, .1), (N, -.1, .1)]
-timespan = (0, 40)
+initial_displacement = [(1, .4, .0), (N, -.0, .0)]
+timespan = (0, 20)
 coupling_H = n_coupled_masses(N, [k1, k2, k1, k2, k1, k2, k1, k2, k1, k2, k1], mass * ones(N))
 
 E = eigenmodes(N, coupling_H)
@@ -108,5 +109,7 @@
 plot(solution, vars=[(1,2), (3,4), (5,6), (7,8), (9,10)])
 savefig("phasefig_$(N).svg")
 
+plot(solution, vars=[(0, 2i-1) for i in 1:N])
+
 # Animation
 plot_chain_animation(N, solution)