I am using Julia to solve a quantum system with a self-consistence iteration algorithm. In my current code, many global variables and arrays are used which, from many documents, are harmful to the performance. However, I don’t know how to avoid this. I post the structure of my code here and wonder how to improve it.
Part 1:
I will define some parameters that describe the system and calculate something from them. These things are global. They will be useful and stay unchanged in Part 2 but will be changed in Part 3.
# define some parameters that will be constant
const m = 1.1 # electron mass
const ManyOtherParameters = ...
# define some parameters that will be changed
B = 2.2 # magnetic field
d = 3.3 # some meaningful distance
ManyOtherThings = ...
# calculate some useful things from magnetic field and other things
# This step consumes a lot of time
function CalculateArray(Mag = B, Dis = d)
...
# calculation
...
return Array
end
A = CalculateArray(B, d)
# an interface to change magnetic field and other parameters that will be useful in Part 3
function ChangeParameters(; Mag = B, Dis = d, OtherThings...)
global B = Mag
global d = Dis
...
global A = CalculateArray(B, d)
end
Part 2:
I will solve the system using the global variables in Part 1.
# one step of iteration
function Solve_OneStep(initialSolution = 0.0)
...
# use B, d, A and initialSolution to find out renewedSolution
...
return renewedSolution
end
# Solve the system
function Solve(initialSolution = 0.0; ControlParameters...)
solution = initialSolution
for i in 1:IterationTimes # for example, 1000
...
# iteration
...
end
return solution
end
Part 3:
I will vary the system parameters and see how the result changes. For example, I will change the value of these global variables, solve it, and record the solution.
using CairoMakie
# set the parameters first
ChangeParameters(Mag = 0.5, Dis = 2.0, OtherThings...)
# the varying range for magnetic field
B_range = 0. : 0.1 : 20.
results = similar(B_range, Float64)
# solve the system for each magnetic fields
for i in eachindex(B_range)
ChangeParameters(Mag = B_range[i])
results[i] = Solve()
end
# plot the results
lines(B_range, results)