-
Notifications
You must be signed in to change notification settings - Fork 2
/
simple_simulation.jl
136 lines (99 loc) · 3.2 KB
/
simple_simulation.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
# NOTE: For improved performence run this script from the terminal as:
# >julia --optimize=3 --math-mode=fast --check-bounds=no simple_simulation.jl
# using Profile
# Profile.clear()
# using ProfileView
using Random
using IterativeSolvers
using SparseArrays
using Statistics
using Langevin
using Langevin.Geometries: Geometry
using Langevin.Lattices: Lattice
using Langevin.HolsteinModels: HolsteinModel
using Langevin.HolsteinModels: assign_μ!, assign_ω!, assign_λ!
using Langevin.HolsteinModels: assign_tij!, assign_ωij!
using Langevin.HolsteinModels: setup_checkerboard!
using Langevin.InitializePhonons: init_phonons_single_site!
using Langevin.FourierAcceleration: FourierAccelerator
using Langevin.LangevinSimulationParameters: SimulationParameters
using Langevin.RunSimulation: run_simulation!
#############################
## DEFINING HOLSTEIN MODEL ##
#############################
# NOTE: Currently Defining Cubic Lattice Geometry
# number of dimensions
ndim = 2
# number of orbitals per unit cell
norbits = 1
# lattice vectors
lvecs = [[1.0,0.0],
[0.0,1.0]]
# basis vectors
bvecs = [[0.0,0.0]]
# defining square lattice geometry
geom = Geometry(ndim, norbits, lvecs, bvecs)
# defining lattice size
L = 4
# constructing finite square lattice
lattice = Lattice(geom,L)
# discretization
Δτ = 0.1
# setting temperature
β = 2.0
# constructing holstein model
holstein = HolsteinModel(geom,lattice,β,Δτ)
# defining nearest neighbor hopping on square lattice
t = 1.0
assign_tij!(holstein, t, 0.0, 1, 1, [1,0,0]) # x direction hopping
assign_tij!(holstein, t, 0.0, 1, 1, [0,1,0]) # y direction hopping
# hamiltonian parameter values
ω = 1.0
g = 0.0
λ = sqrt(2.0*ω)*g
μ = -0.75 # for half-filling
assign_ω!(holstein, ω, 0.0)
assign_λ!(holstein, λ, 0.0)
assign_μ!(holstein, μ, 0.0)
# organize electron hoppings for checkerboard decomposition
setup_checkerboard!(holstein)
# intialize phonon field
init_phonons_single_site!(holstein)
###################################
## DEFINING FOURIER ACCELERATION ##
###################################
# langevin time step
Δt = 1e-3
# mass for fourier acceleration: increasing the mass reduces the
# amount of acceleration
mass = 0.5
# defining FourierAccelerator type
fa = FourierAccelerator(holstein,mass,Δt)
####################################
## DEFINING SIMULATION PARAMETERS ##
####################################
# tolerace of IterativeSolvers
tol = 1e-4
# number of thermalization steps
burnin = 250000
# total number of steps
nsteps = 1000000
# measurement frequency
meas_freq = 100
# number of bins
num_bins = 1
# euler or runge-kutta updates
euler = false
# filepath to where to write data
filepath = "."
# name of folder for data to get dumped into
foldername = "test"
sim_params = SimulationParameters(Δt,euler,tol,burnin,nsteps,meas_freq,num_bins,filepath,foldername)
########################
## RUNNING SIMULATION ##
########################
simulation_time, measurement_time, write_time, iters = run_simulation!(holstein, sim_params, fa)
println("Simulation Time (min) = ", simulation_time)
println("Measurement Time (min) = ", measurement_time)
println("Write Time (min) = ", write_time)
println("Average Iterative Solver Iterations = ", iters)