Skip to content
Snippets Groups Projects
Commit 8e03403b authored by Georg Wolfgang Winkler's avatar Georg Wolfgang Winkler
Browse files

added skeletons for ex02 and Lanczos example

parent 29ee7cd7
No related branches found
No related tags found
No related merge requests found
%% Cell type:code id: tags:
``` python
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
```
%% Cell type:code id: tags:
``` python
def step(v1,v0,A,beta):
w=np.dot(A,v1)
alpha=np.real(np.dot(np.conj(w),v1))
w-=alpha*v1+beta*v0
beta=np.linalg.norm(w)
w/=beta
return w,alpha,beta
```
%% Cell type:code id: tags:
``` python
N=800
M=200
#mat=(np.random.rand(N,N)-0.5)+1j*(np.random.rand(N,N)-0.5)
#mat+=np.conj(mat.T)
#evs=np.linalg.eigvalsh(mat)
```
%% Cell type:code id: tags:
``` python
vs=np.zeros((2,N),dtype=complex)
vs[1]=np.random.rand(N)
vs[1]/=np.linalg.norm(vs[1])
alphas=np.zeros(M+1)
betas=np.zeros(M+1)
for i in xrange(M):
vs[i%2],alphas[i],betas[i+1]=step(vs[(i+1)%2],vs[i%2],mat,betas[i])
alphas[-1]=np.real(np.dot(np.conj(vs[i%2]),np.dot(mat,vs[i%2])))
```
%% Cell type:code id: tags:
``` python
plt.figure(figsize=(15,10))
t=np.diag(alphas[:],0)+np.diag(betas[1:],1)+np.diag(betas[1:],-1)
for i in xrange(2,M):
plt.plot(i*np.ones(i),np.linalg.eigvalsh(t[:i,:i]),'ok')
for i in xrange(20):
plt.axhline(evs[-(1+i)])
plt.ylim(1.01*evs[-1],1.01*evs[-20]);
```
%% Output
%% Cell type:code id: tags:
``` python
```
// ---- CQP FS16, Ex 3.1 Skeleton code ----
// ----------------------------------------
#include <valarray>
#include <vector>
#include <iostream>
#include <cmath>
#include <complex>
typedef std::valarray< std::complex<double> > state_vector_t;
typedef std::vector<double> param_vector_t;
typedef unsigned state_t;
// Implement the time evolution of state x by n*dt here. Use the split
// operator scheme we discussed in the exercise class.
void evolve(unsigned l, const param_vector_t &j, const param_vector_t &hx, double dt, unsigned n, state_vector_t& x)
{
// ........
}
void printMagnetization(unsigned l, const state_vector_t& v, std::ostream& outputStream)
{
std::vector<double> res(l,0.0);
for(int r=0; r<l; r++){
double m(0);
// Measure magnetization at site r here. To be completed.
outputStream << m << "\t\t";
}
outputStream << std::endl;
}
int main(){
// ======== input
unsigned l = 15; // length of chain
param_vector_t j(l-1,0.0); // Ising coupling
param_vector_t hx(l,0.4); // transverse magnetic field
double tmax = 20.; // maximum time reached in time evolutin
double tmeas = 0.5; // time after which a measurment is performed
double dt = 0.1; // time step used for the evolution
unsigned nmeas = tmax/tmeas; // number of measurements
unsigned nstep = tmeas/dt; // number of steps between measurements
// ======== input
// Starting configuration
state_t initialInd = 1<<7;
state_vector_t x(0., 1 << l);
x[initialInd] = 1;
printMagnetization(l,x,std::cout);
// Do time evolution
for(int n=1; n<=nmeas; ++n){
evolve(l, j, hx, dt, nstep, x);
printMagnetization(l,x,std::cout);
}
return 0;
}
# ---- CQP FS16, Ex3.1 skeleton ----
# ----------------------------------
import numpy as np
import matplotlib.pyplot as plt
def evolve(l, j, hx, v, dt, n):
""" Evolve the state v in time by n * dt
Parameters:
l length of chain
j Ising coupling for each neighbors, 1-d numpy array with length at least l-1
hx transverse field for each site, 1-d numpy array with length at least l
v the state at time t, complex 1-d numpy array with length at least 1 << l
dt timestep
n number of consecutive timesteps
Return:
v the state at time t+n*dt
"""
# Implement the time evolution of the state v here.
# Use the split-operator scheme we discussed in the
# exercise class.
return v
def magnetization(v, l):
""" Measure the magnetization per site for the state v
Parameters:
v the state to be measured, complex 1-d numpy array with length at least 1 << l
l length of chain
Return:
m the magnetization at each site, 1-d numpy array with length l
"""
# Implement the measurement of the magnetization here.
return m
l = 9 # length of chain
dim = 1 << l
j = 1.00*np.ones((l)) # Ising coupling
hx = 0.40*np.ones((l-1)) # transverse field strength Gamma
# Starting configuration
ind = int('000010000',base=2) # specify spin configuration: 0 = down, 1 = up
v = np.zeros((dim),dtype=complex)
v[ind] = 1. # start time evolution with a single basis state
# Parameters for time evolution
tend = 10. # duration of the time evolution
tmeas = 0.5 # time after which a measurment is performed
dt = 0.1 # time step used for the evolution
nmeas = int((tend)/tmeas) # number of measurements
nstep = int((tmeas)/dt) # number of steps between measurements
# measure magnetization
m=np.zeros((nmeas,l));
m[0,:] = magnetization(v, l)
# Do time evolution
for t in range(1,nmeas):
print(t*tmeas)
v=evolve(l,j,hx,v,dt,nstep)
m[t,:] = magnetization(v, l)
# plot magnetization
plt.pcolor(np.array(range(l+1)),np.linspace(0.,tend,nmeas),m,
vmin=-1.0, vmax=1.0)
plt.xlabel('site')
plt.ylabel('time')
plt.axis([0,l,0,tend])
plt.title("Magnetisation, J = "+str(j[0])+" , hx = "+str(hx[0]))
plt.colorbar()
plt.show()
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment