-- Matrix mechanics 1 N = 10 x = zero(N,N) E = zero(N,N) for(n, 1, N - 1, x[n, n + 1] = exp(-i omega0 t) sqrt(n hbar / (2 m omega0))) x = x + conj(transpose(x)) v = d(x,t) x2 = dot(x,x) v2 = dot(v,v) for(n, 1, N, E[n,n] = (n - 1/2) hbar omega0) "Verify equation (1)" -- row N and column N are incomplete due to finite N for(n, 1, N - 1, for(m, 1, N - 1, check(1/2 m v2[n,m] + 1/2 m omega0^2 x2[n,m] == E[n,m]) )) "ok" Run