Jesús Andrés Espinoza Valverde.
The GNM is simple toy model in QCD that shares very interesting properties with several Condensed Matter theories, specially in the context of superconductors. In this sense, the calculations on the Gross-Neveu model can be ingeniously interpreted and used for the study of condensed matter systems. This fact is the main motivation for the present work.
This repository contains an algorithmic framework to approximate vacuum expectation values of observables of the Gross-Neveu model for the
The proposal is to generalized the algorithms implemented here to be compatible with
The Lagrangian of that describes the theory is the following:
$$
\mathcal{L}[\bar{\psi},\psi,\phi] = \sum^{N_f}{f=1}\bar{\psi}(\textbf{x})[\gamma{\mu}\partial_{\mu}+m^{(f)}+\sqrt{g}\phi(\textbf{x})]\psi^{(f)}(\textbf{x}) + \frac{1}{2}\phi^{2}(\textbf{x}).
$$
Here
We associate to the continuum space-time a square like lattice
On the lattice is easy to split the action on its fermionic and scalar parts:
$$
S_{F}[\hat{\psi},\psi,\phi] = \sum_{\textbf{n}, \textbf{m}\in \Lambda}\hat{\psi}(\textbf{n})M(\textbf{n}, \textbf{m})\psi(\textbf{m}), \hspace{0.5cm} S_S[\phi]=\frac{1}{2}\sum_{\textbf{n}\in \Lambda} \phi^2 (\textbf{n}),
$$
where
Using this result, we express the total partition function as: $$ Z_{GN} = \int \prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2} \mathrm{det}(M[\phi]), $$
Using the previous results we transform the expression for vacuum expectation values from the continuum to the lattice as: $$ \langle O \rangle = \frac{1}{Z_{GN}}\int \mathcal{D}[\bar{\psi},\psi,\phi]e^{-S[\bar{\psi},\psi,\phi]}O[\bar{\psi},\psi,\phi]. $$
$$ \langle O \rangle = \frac{\int^{+\infty}{-\infty} \dots \int^{+\infty}{-\infty}\prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2} \hspace{1mm} \mathrm{det}(M[\phi]) \hspace{1.5mm} O[M[\phi], \phi]}{\int^{+\infty}{-\infty} \dots \int^{+\infty}{-\infty} \prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2}\hspace{1mm} \mathrm{det}(M[\phi])} \hspace{3cm} (1) $$
We compute the required expectation values computing this expression numerically. The integrals contained
On the lattice, the fermionic part of the action is given by $$ S_F[\bar{\psi}, \psi, Q] = \sum_{\mathbf{n},\mathbf{m} \in \Lambda} \bar{\psi}(\mathbf{n}) M(\mathbf{n}, \mathbf{m}) \psi(\mathbf{m}). $$
Where, in index form: $$ M(\mathbf{n}, \mathbf{m}){\alpha \beta} = [m+\sqrt{g}Q(\mathbf{n})] \delta{\mathbf{n}, \mathbf{m}} \delta_{\alpha \beta} - (2a)^{-1}\sum_{\mu = 1} ^{2} (\gamma_{\mu}){\alpha \beta} [\delta{\mathbf{n}+\hat{\mu}, \mathbf{m}} - \delta_{\mathbf{n}-\hat{\mu}, \mathbf{m}}]. $$
Now, in order to avoid the doubling problem, we add to the last expression the Wilson term in the real space representation (omitting spin indices):
$$
M(\mathbf{n}, \mathbf{m}) - \frac{a}{2} \sum_{\mu = 1} ^{2} \frac{\delta_{\mathbf{n}+\hat{\mu}, \mathbf{m}} - 2\delta_{\mathbf{n}, \mathbf{m}}+ \delta_{\mathbf{n}-\hat{\mu}, \mathbf{m}}}{a^2} \longrightarrow M(\mathbf{n}, \mathbf{m}),
$$
Which implies:
$$
M(\mathbf{n}, \mathbf{m}){\alpha \beta} = [m+\sqrt{g}Q(\mathbf{n})] \delta{\mathbf{n}, \mathbf{m}} \delta_{\alpha \beta} - (2a)^{-1} \sum_{\mu = 1} ^{2} { -2\delta_{\mathbf{n}, \mathbf{m}} + \left(I+\gamma_{\mu}\right)\delta_{\mathbf{n}+\hat{\mu}, \mathbf{m}} + \left(I-\gamma_{\mu}\right)\delta_{\mathbf{n}-\hat{\mu}, \mathbf{m}} }
$$
Then, defining
For a
With the above labeling it is possible to associate a single
Let's start with finding the explicit form of
$$
\delta_{\mathbf{n}-\hat{2}, \mathbf{m}} = \delta_{[n_1,n_2-1],[m_1, m_2]} = \delta_{N , N' }; \hspace{4mm} \begin{cases}
r_2 = n_2 - 1 \mod N_1\
N = r_2 N_2 + n_1 \
N'= m_2 N_2 + m_1
\end{cases}
$$
The Modulo operations are introduced to take into account periodic boundary conditions. To make the p.b.c in the
Making the definition
-
\Gamma_{-1}\left[ \begin{matrix} 0 & -1 & 0 & 0\ 1 & 0 & 0 & 0\ 0 & 0 & 0 & -1\ 0 & 0 & 1 & 0\ \end{matrix} \right]
\Gamma_{2}\left[ \begin{matrix} 0 & 0 & 1 & 0\ 0 & 0 & 0 & 1\ 1 & 0 & 0 & 0\ 0 & 1 & 0 & 0\ \end{matrix} \right]
\Gamma_{-2}\left[ \begin{matrix} 0 & 0 & 1 & 0\ 0 & 0 & 0 & 1\ 1 & 0 & 0 & 0\ 0 & 1 & 0 & 0\ \end{matrix} \right] \ = \left[ \begin{matrix} Q'(0,0) & \Gamma_{-1}-\Gamma_{1} & -(\Gamma_{2}+\Gamma_{-2}) & 0\ \Gamma_{1}-\Gamma_{-1} & Q'(1,0) & 0 & -(\Gamma_{2}+\Gamma_{-2})\ -(\Gamma_{2}+\Gamma_{-2}) & 0 & Q'(0,1) & \Gamma_{-1}-\Gamma_{1})\ 0 & -(\Gamma_{2}+\Gamma_{-2}) & \Gamma_{1}-\Gamma_{-1} & Q'(1,1)\ \end{matrix} \right] $$
Now, taking into account the definitions of
And then: $$ \Gamma_{1} - \Gamma_{-1} = \left[\begin{matrix} 0 & 1\ 1 & 0\ \end{matrix} \right] $$
And therefore: $$ M(\mathbf{n}, \mathbf{m}) = \left[ \begin{matrix} \left[\begin{matrix}Q'(0,0) & 0\0 & Q'(0,0)\\end{matrix}\right] & \left[\begin{matrix} 0 & -1\-1 & 0\\end{matrix}\right] & \left[\begin{matrix} -1 & 0\0 & -1\\end{matrix}\right] & \left[\begin{matrix} 0 & 0\0 & 0\\end{matrix}\right]\ \left[\begin{matrix} 0 & 1\1 & 0\\end{matrix}\right] & \left[\begin{matrix}Q'(1,0) & 0\0 & Q'(1,0)\\end{matrix}\right] & \left[\begin{matrix} 0 & 0\0 & 0\\end{matrix}\right] & \left[\begin{matrix} -1 & 0\0 & -1\\end{matrix}\right]\ \left[\begin{matrix} -1 & 0\0 & -1\\end{matrix}\right] & \left[\begin{matrix} 0 & 0\0 & 0\\end{matrix}\right] & \left[\begin{matrix}Q'(0,1) & 0\0 & Q'(0,1)\\end{matrix}\right] &\left[\begin{matrix} 0 & -1\-1 & 0\\end{matrix}\right]\ \left[\begin{matrix} 0 & 0\0 & 0\\end{matrix}\right] & \left[\begin{matrix} -1 & 0\0 & -1\\end{matrix}\right] & \left[\begin{matrix} 0 & 1\1 & 0\\end{matrix}\right] & \left[\begin{matrix}Q'(1,1) & 0\0 & Q'(1,1)\\end{matrix}\right]\ \end{matrix} \right] \ =\left[\begin{matrix} Q'(0,0) & 0 & 0 & -1 & -1 & 0 & 0 & 0\ 0 & Q'(0,0) & -1 & 0 & 0 & -1 & 0 & 0\ 0 & 1 & Q'(1,0) & 0 & 0 & 0 & -1 & 0\ 1 & 0 & 0 & Q'(1,0) & 0 & 0 & 0 & -1\ -1 & 0 & 0 & 0 & Q'(0,1) & 0 & 0 & -1\ 0 & -1 & 0 & 0 & 0 & Q'(0,1) & -1 & 0\ 0 & 0 & -1 & 0 & 0 & 1 & Q'(1,1) & 0\ 0 & 0 & 0 & -1 & 1 & 0 & 0 & Q'(1,1)\ \end{matrix}\right] $$ The explicit forms of matrices M (n, m) corresponding to larger lattices are too extensive to be shown here. Because of this we will present examples of these bigger matrices using pictorial representations.
The observables are computed through Monte Carlo simulations using the Importance Sampling technique. Due to the extremely high dimensionality of the problem and the large amount of Monte Carlo steps required to obtain good results, a really high quality pseudo random number generator is needed. We've chosen Ranlux - Martin Lüscher - CERN , which exhibits despicable correlation and astronomical periods in the sequences it produces. This algorithm generates uniform random numbers between
void gaussian_random(double* gauss_rands,
const int M_dim,
const int seed)
{
int i{0};
double uniform_rands[2];
double u1, u2;
double gamma0, gamma1;
rlxd_init(1,seed);
while(i < M_dim/2)
{
ranlxd(uniform_rands, 2);
u1 = uniform_rands[0];
u2 = uniform_rands[1];
u1 = 2*u1-1;
u2 = 2*u2-1;
gamma0 = pow(u1, 2) + pow(u2, 2);
if (gamma0 < 1 && gamma0 != 0)
{
gamma1 = sqrt(-2*log(gamma0)/gamma0);
u1 = gamma1*u1;
u2 = gamma1*u2;
gauss_rands[2*i] = u1;
gauss_rands[2*i+1] = u2;
i += 1;
}
}
}
It is clear that this piece of code is crying for being optimized, specially because it is used inside the main Monte Carlo for-loop.
To calculate the determinants of equation umfpack library of suitesparse was used. This suite is specialized on sparse-matrix linear algebra.
Luckily there are a couple of cases with can solve exactly by analytical means. This is extremely useful because it allow us to test our numerical results. For the free case
$$ \langle \phi \rangle_{\phi} = \frac{\int^{+\infty}{-\infty} \dots \int^{+\infty}{-\infty}\prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2} \hspace{1.5mm} \left(\frac{1}{N_p} \sum_{\textbf m \in \Lambda}\phi(\textbf m)\right)}{\int^{+\infty}{-\infty} \dots \int^{+\infty}{-\infty} \prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2}} =0, $$
$$ \langle \phi \rangle_{\phi} = \frac{\int^{+\infty}{-\infty} \dots \int^{+\infty}{-\infty}\prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2} \hspace{1.5mm} \left(\frac{1}{N_p} \sum_{\textbf m \in \Lambda}\phi^2(\textbf m)\right)}{\int^{+\infty}{-\infty} \dots \int^{+\infty}{-\infty} \prod_{\textbf n \in \Lambda} \frac{d\phi(\textbf n)}{\sqrt{2\pi}} e^{-\phi^2(\textbf n)/2}} = \frac{\frac{2^{\frac{3}{2}}\Gamma(\frac{3}{2})}{\sqrt{2\pi}}}{\frac{2^{\frac{1}{2}}\Gamma(\frac{1}{2})}{\sqrt{2\pi}}} \ =\frac{(2\hspace{1mm}2^{1/2})(\frac{1}{2}\Gamma(\frac{1}{2}))}{(\hspace{1mm}2^{1/2})(\Gamma(\frac{1}{2}))} = \frac{\sqrt{2\pi}}{\sqrt{2\pi}} = 1, $$
the first result is simply because the integration argument is a odd function. To use our program to calculate the mean scalar field expected value for the free case , we use a input file of the form:
10 # Number of columns of the lattice
10 # Number of rows of the lattice
0 # Fermion mass
0 # g coupling constant
mean_scalar # Observable name
2000 # Number of Monte Carlo Steps
0 # print matrix (Not recommended)And we obtain:
Result: -4.18852e-05 + 2.83215e-37i Processes: 4 M-constr-time: 0.000921292 Comp-time: 1.92721 Comm-time: 0.0169848 Total-time: 1.94511Note
The square scalar case will need a input as the following:
10 # Number of columns of the lattice
10 # Number of rows of the lattice
0 # Fermion mass
0 # g coupling constant
mean_sq_scalar # Observable name
2000 # Number of Monte Carlo Steps
0 # print matrix (Not recommended)And it produces:
Result: 1.00077 + 2.36228e-32i Processes: 4 M-constr-time: 0.000912439 Comp-time: 1.94237 Comm-time: 0.000912096 Total-time: 1.94419With
The structure tree of the software is the following:
.
├── benchmarking
│ ├── plotter.sh
│ ├── runner.sh
│ └── TMP
│ ├── sdata.png
│ ├── str_avg_std_time.dat
│ ├── str_speed_up.dat
│ ├── wk_avg_std_time_mean_scalar.dat
│ ├── wk_avg_std_time_mean_sq_scalar.dat
│ ├── wk_efficiency_mean_scalar.dat
│ ├── wk_efficiency_mean_sq_scalar.dat
│ ├── wk_errors_mean_scalar.dat
│ └── wk_errors_mean_sq_scalar.dat
├── images
│ ├── 10x10.png
│ ├── 2x2.png
│ ├── 3x3.png
│ ├── 4x4.png
│ ├── 5x5.png
│ ├── 5x8.png
│ ├── 6x4.png
│ ├── error.png
│ └── space-time.png
├── include
│ ├── conditionalOStreamClass.h
│ ├── determinantWrapper.h
│ ├── gaussianRandomGenerator.h
│ ├── importanceSamplingIntegrator.h
│ ├── mpiCommManager.h
│ ├── observables.h
│ ├── ranlxd.h
│ ├── sparseDiracMatrixClass.h
│ ├── suiteSparseDet.h
│ └── utilities.h
├── input.inp
├── libref
│ ├── RANLUX_Guide.pdf
│ ├── UMFPACK_QuickStart.pdf
│ └── UMFPACK_UserGuide.pdf
├── main.cc
├── Makefile
├── python
│ ├── matrix.py
│ └── polinom_det.py
├── README.md
├── src
│ ├── determinantWrapper.cc
│ ├── gaussianRandomGenerator.cc
│ ├── importanceSamplingIntegrator.cc
│ ├── ranlxd.c
│ ├── sparseDiracMatrixClass.cc
│ ├── suiteSparseDet.c
│ └── utilities.cc
└── tests
8 directories, 49 files
To download the required dependencies run (only for debian systems):
$ sudo apt-get install libsuitesparse-dev
$ sudo apt-get install libsuitesparse-docTo compile it, just run:
$ makeTo run it:
$ mpirun -np [nprocs] gn [input file]The different parameters of the simulation are specified in the input file, for example:
10 # Number of columns of the lattice
10 # Number of rows of the lattice
0 # Fermion mass
0 # g coupling constant
mean_sq_scalar # Observable name
2000 # Number of Monte Carlo Steps
0 # print matrix (Not recommended)