eos_quark_cfl6.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2020, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 /** \file eos_quark_cfl6.h
24  \brief File defining \ref o2scl::eos_quark_cfl6
25 */
26 #ifndef CFL6_EOS_H
27 #define CFL6_EOS_H
28 
29 #include <iostream>
30 #include <o2scl/test_mgr.h>
31 #include <o2scl/eos_quark_cfl.h>
32 
33 #ifndef DOXYGEN_NO_O2NS
34 namespace o2scl {
35 #endif
36 
37  /** \brief An EOS like \ref eos_quark_cfl but
38  with a color-superconducting 't Hooft interaction
39 
40  Beginning with the Lagrangian:
41  \f[
42  {\cal L} = {\cal L}_{Dirac} + {\cal L}_{NJL} +
43  {\cal L}_{'t Hooft} + {\cal L}_{SC} + {\cal L}_{SC6}
44  \f]
45  \f[
46  {\cal L}_{Dirac} = {\bar q} \left( i \partial -m -
47  \mu \gamma^0 \right) q
48  \f]
49  \f[
50  {\cal L}_{NJL} = G_S \sum_{a=0}^8
51  \left[ \left( {\bar q} \lambda^a q \right)^2
52  - \left( {\bar q} \lambda^a \gamma^5 q \right)^2 \right]
53  \f]
54  \f[
55  {\cal L}_{'t Hooft} = G_D \left[
56  \mathrm{det}_f {\bar q} \left(1-\gamma^5 \right) q
57  +\mathrm{det}_f {\bar q} \left(1+\gamma^5 \right) q
58  \right]
59  \f]
60  \f[
61  {\cal L}_{SC} = G_{DIQ}
62  \left( {\bar q}_{i \alpha} i \gamma^5
63  \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma}
64  q^C_{j \beta} \right)
65  \left( {\bar q}_{\ell \delta} i \gamma^5
66  \epsilon^{\ell m k}
67  \epsilon^{\delta \varepsilon \gamma}
68  q^C_{m \varepsilon} \right)
69  \f]
70  \f[
71  {\cal L}_{SC6} = K_D
72  \left( {\bar q}_{i \alpha} i \gamma^5
73  \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma}
74  q^C_{j \beta} \right)
75  \left( {\bar q}_{\ell \delta} i \gamma^5
76  \epsilon^{\ell m n}
77  \epsilon^{\delta \varepsilon \eta}
78  q^C_{m \varepsilon} \right)
79  \left( {\bar q}_{k \gamma} q_{n \eta} \right)
80  \f]
81 
82  We can simplify the relevant terms in \f${\cal L}_{NJL}\f$:
83  \f[
84  {\cal L}_{NJL} = G_S \left[
85  \left({\bar u} u\right)^2+
86  \left({\bar d} d\right)^2+
87  \left({\bar s} s\right)^2
88  \right]
89  \f]
90  and in \f${\cal L}_{'t Hooft}\f$:
91  \f[
92  {\cal L}_{NJL} = G_D \left(
93  {\bar u} u {\bar d} d {\bar s} s
94  \right)
95  \f]
96 
97  Using the definition:
98  \f[
99  \Delta^{k \gamma} = \left< {\bar q} i \gamma^5
100  \epsilon \epsilon q^C_{} \right>
101  \f]
102  and the ansatzes:
103  \f[
104  ({\bar q}_1 q_2) ({\bar q}_3 q_4) \rightarrow
105  {\bar q}_1 q_2 \left< {\bar q}_3 q_4 \right>
106  +{\bar q}_3 q_4 \left< {\bar q}_1 q_2 \right>
107  -\left< {\bar q}_1 q_2 \right> \left< {\bar q}_3 q_4 \right>
108  \f]
109  \f[
110  ({\bar q}_1 q_2) ({\bar q}_3 q_4) ({\bar q}_5 q_6) \rightarrow
111  {\bar q}_1 q_2 \left< {\bar q}_3 q_4 \right>
112  \left< {\bar q}_5 q_6 \right>
113  +{\bar q}_3 q_4 \left< {\bar q}_1 q_2 \right>
114  \left< {\bar q}_5 q_6 \right>
115  +{\bar q}_5 q_6 \left< {\bar q}_1 q_2 \right>
116  \left< {\bar q}_3 q_4 \right>
117  -2\left< {\bar q}_1 q_2 \right> \left< {\bar q}_3 q_4 \right>
118  \left< {\bar q}_5 q_6 \right>
119  \f]
120  for the mean field approximation, we can rewrite the Lagrangian
121  \f[
122  {\cal L}_{NJL} = 2 G_S \left[
123  \left( {\bar u} u \right) \left< {\bar u} u \right>
124  +\left( {\bar d} d \right) \left< {\bar d} d \right>
125  +\left( {\bar s} s \right) \left< {\bar s} s \right>
126  - \left< {\bar u} u \right>^2
127  - \left< {\bar d} d \right>^2
128  - \left< {\bar s} s \right>^2
129  \right]
130  \f]
131  \f[
132  {\cal L}_{'t Hooft} = - 2 G_D \left[
133  \left( {\bar u} u \right) \left< {\bar u} u \right>
134  \left< {\bar s} s \right>
135  + \left( {\bar d} d \right) \left< {\bar u} u \right>
136  \left< {\bar s} s \right>
137  + \left( {\bar s} s \right) \left< {\bar u} u \right>
138  \left< {\bar d} d \right>
139  - 2 \left< {\bar u} u \right>\left< {\bar d} d \right>
140  \left< {\bar s} s \right>
141  \right]
142  \f]
143  \f[
144  {\cal L}_{SC} = G_{DIQ} \left[
145  \Delta^{k \gamma}
146  \left( {\bar q}_{\ell \delta} i \gamma^5
147  \epsilon^{\ell m k}
148  \epsilon^{\delta \varepsilon \gamma}
149  q^C_{m \varepsilon} \right)
150  + \left( {\bar q}_{i \alpha} i \gamma^5
151  \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma}
152  q^C_{j \beta} \right)
153  \Delta^{k \gamma \dagger}
154  - \Delta^{k \gamma}
155  \Delta^{k \gamma \dagger}
156  \right]
157  \f]
158  \f[
159  {\cal L}_{SC6} = K_D \left[
160  \left( {\bar q}_{m \varepsilon} q_{n \eta} \right)
161  \Delta^{k \gamma} \Delta^{m \varepsilon \dagger}
162  + \left( {\bar q}_{i \alpha} i \gamma^5
163  \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma}
164  q^C_{j \beta} \right)
165  \Delta^{m \varepsilon \dagger}
166  \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
167  \right]
168  \f]
169  \f[
170  + K_D \left[\Delta^{k \gamma}
171  \left( {\bar q}_{\ell \delta} i \gamma^5
172  \epsilon^{\ell m n}
173  \epsilon^{\delta \varepsilon \eta}
174  q^C_{m \varepsilon} \right)
175  \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
176  -2
177  \Delta^{k \gamma} \Delta^{m \varepsilon \dagger}
178  \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
179  \right]
180  \f]
181 
182  If we make the definition \f$ {\tilde \Delta} =
183  2 G_{DIQ} \Delta \f$
184 
185  \hline
186  <b>References:</b>
187 
188  Created for \ref Steiner05.
189  */
190  class eos_quark_cfl6 : public eos_quark_cfl {
191  public:
192 
193  eos_quark_cfl6();
194 
195  virtual ~eos_quark_cfl6();
196 
197  /** \brief Calculate the EOS
198  \nothing
199 
200  Calculate the EOS from the quark condensates. Return the mass
201  gap equations in \c qq1, \c qq2, \c qq3, and the normal gap
202  equations in \c gap1, \c gap2, and \c gap3.
203 
204  Using \c fromqq=true as in eos_quark_njl and nambujl_temp_eos
205  does not work here and will return an error.
206 
207  If all of the gaps are less than gap_limit, then the
208  nambujl_temp_eos::calc_temp_p() is used, and \c gap1, \c gap2,
209  and \c gap3 are set to equal \c u.del, \c d.del, and \c s.del,
210  respectively.
211 
212  */
213  virtual int calc_eq_temp_p(quark &u, quark &d, quark &s,
214  double &qq1, double &qq2, double &qq3,
215  double &gap1, double &gap2, double &gap3,
216  double mu3, double mu8,
217  double &n3, double &n8, thermo &qb,
218  double temper);
219 
220  /// The momentum integrands
221  virtual int integrands(double p, double res[]);
222 
223  /// Check the derivatives specified by eigenvalues()
224  virtual int test_derivatives(double lmom, double mu3, double mu8,
225  test_mgr &t);
226 
227  /** \brief Calculate the energy eigenvalues and their derivatives
228  \nothing
229 
230  Given the momentum \c mom, and the chemical potentials
231  associated with the third and eighth gluons (\c mu3 and \c mu8),
232  this computes the eigenvalues of the inverse propagator and
233  the assocated derivatives.
234 
235  Note that this is not the same as eos_quark_cfl::eigenvalues()
236  which returns \c dedmu rather \c dedqqu.
237  */
238  virtual int eigenvalues6(double lmom, double mu3, double mu8,
239  double egv[36], double dedmuu[36],
240  double dedmud[36], double dedmus[36],
241  double dedmu[36], double dedmd[36],
242  double dedms[36], double dedu[36],
243  double dedd[36], double deds[36],
244  double dedmu3[36], double dedmu8[36]);
245 
246  /** \brief Construct the matrices, but don't solve the eigenvalue
247  problem
248  \nothing
249 
250  This is used by check_derivatives() to make sure that the derivative
251  entries are right.
252  */
253  virtual int make_matrices(double lmom, double mu3, double mu8,
254  double egv[36], double dedmuu[36],
255  double dedmud[36], double dedmus[36],
256  double dedmu[36], double dedmd[36],
257  double dedms[36], double dedu[36],
258  double dedd[36], double deds[36],
259  double dedmu3[36], double dedmu8[36]);
260 
261  /// The color superconducting 't Hooft coupling (default 0)
262  double KD;
263 
264  /// Return string denoting type ("eos_quark_cfl6")
265  virtual const char *type() { return "eos_quark_cfl6"; };
266 
267  /** \brief The absolute value below which the CSC 't Hooft coupling
268  is ignored(default \f$ 10^{-6} \f$)
269  */
270  double kdlimit;
271 
273  ubvector_complex;
275  ubmatrix_complex;
276 
277  protected:
278 
279 #ifndef DOXYGEN_INTERNAL
280 
281  /// Set the quark effective masses from the gaps and the condensates
282  int set_masses();
283 
284  /// The size of the matrix to be diagonalized
285  static const int mat_size=36;
286 
287  /// Storage for the inverse propagator
288  gsl_matrix_complex *iprop6;
289 
290  /// The eigenvectors
291  gsl_matrix_complex *eivec6;
292 
293  /// The derivative wrt the ds gap
295 
296  /// The derivative wrt the us gap
298 
299  /// The derivative wrt the ud gap
301 
302  /// The derivative wrt the up quark condensate
304 
305  /// The derivative wrt the down quark condensate
307 
308  /// The derivative wrt the strange quark condensate
310 
311  /// Storage for the eigenvalues
312  gsl_vector *eval6;
313 
314  /// GSL workspace for the eigenvalue computation
315  gsl_eigen_hermv_workspace *w6;
316 
317  private:
318 
320  eos_quark_cfl6& operator=(const eos_quark_cfl6&);
321 
322 #endif
323 
324  };
325 
326 #ifndef DOXYGEN_NO_O2NS
327 }
328 #endif
329 
330 #endif
boost::numeric::ublas::matrix
o2scl::eos_quark_cfl6::eigenvalues6
virtual int eigenvalues6(double lmom, double mu3, double mu8, double egv[36], double dedmuu[36], double dedmud[36], double dedmus[36], double dedmu[36], double dedmd[36], double dedms[36], double dedu[36], double dedd[36], double deds[36], double dedmu3[36], double dedmu8[36])
Calculate the energy eigenvalues and their derivatives.
boost::numeric::ublas::vector
o2scl::eos_quark_cfl6::integrands
virtual int integrands(double p, double res[])
The momentum integrands.
o2scl::eos_quark_cfl6::eivec6
gsl_matrix_complex * eivec6
The eigenvectors.
Definition: eos_quark_cfl6.h:291
o2scl::eos_quark_cfl6::kdlimit
double kdlimit
The absolute value below which the CSC 't Hooft coupling is ignored(default )
Definition: eos_quark_cfl6.h:265
o2scl::eos_quark_cfl6
An EOS like eos_quark_cfl but with a color-superconducting 't Hooft interaction.
Definition: eos_quark_cfl6.h:190
o2scl::eos_quark_cfl6::mat_size
static const int mat_size
The size of the matrix to be diagonalized.
Definition: eos_quark_cfl6.h:285
o2scl::eos_quark_cfl6::calc_eq_temp_p
virtual int calc_eq_temp_p(quark &u, quark &d, quark &s, double &qq1, double &qq2, double &qq3, double &gap1, double &gap2, double &gap3, double mu3, double mu8, double &n3, double &n8, thermo &qb, double temper)
Calculate the EOS.
o2scl::eos_quark_cfl6::dipdqqu
ubmatrix_complex dipdqqu
The derivative wrt the up quark condensate.
Definition: eos_quark_cfl6.h:303
thermo_tl< double >
o2scl::quark
o2scl::eos_quark_cfl6::type
virtual const char * type()
Return string denoting type ("eos_quark_cfl6")
Definition: eos_quark_cfl6.h:265
o2scl::eos_quark_cfl6::dipdgapu
ubmatrix_complex dipdgapu
The derivative wrt the ds gap.
Definition: eos_quark_cfl6.h:294
o2scl::eos_quark_cfl6::dipdgaps
ubmatrix_complex dipdgaps
The derivative wrt the ud gap.
Definition: eos_quark_cfl6.h:300
o2scl::eos_quark_cfl6::dipdgapd
ubmatrix_complex dipdgapd
The derivative wrt the us gap.
Definition: eos_quark_cfl6.h:297
o2scl::eos_quark_cfl6::eval6
gsl_vector * eval6
Storage for the eigenvalues.
Definition: eos_quark_cfl6.h:312
o2scl::eos_quark_cfl6::make_matrices
virtual int make_matrices(double lmom, double mu3, double mu8, double egv[36], double dedmuu[36], double dedmud[36], double dedmus[36], double dedmu[36], double dedmd[36], double dedms[36], double dedu[36], double dedd[36], double deds[36], double dedmu3[36], double dedmu8[36])
Construct the matrices, but don't solve the eigenvalue problem.
o2scl::test_mgr
o2scl::eos_quark_cfl::temper
double temper
Temperature.
Definition: eos_quark_cfl.h:440
o2scl::eos_quark_cfl6::iprop6
gsl_matrix_complex * iprop6
Storage for the inverse propagator.
Definition: eos_quark_cfl6.h:288
o2scl::eos_quark_cfl6::set_masses
int set_masses()
Set the quark effective masses from the gaps and the condensates.
o2scl::eos_quark_cfl6::dipdqqs
ubmatrix_complex dipdqqs
The derivative wrt the strange quark condensate.
Definition: eos_quark_cfl6.h:309
o2scl::eos_quark_cfl6::test_derivatives
virtual int test_derivatives(double lmom, double mu3, double mu8, test_mgr &t)
Check the derivatives specified by eigenvalues()
o2scl::eos_quark_cfl6::dipdqqd
ubmatrix_complex dipdqqd
The derivative wrt the down quark condensate.
Definition: eos_quark_cfl6.h:306
o2scl::eos_quark_cfl
Nambu Jona-Lasinio model with a schematic CFL di-quark interaction at finite temperature.
Definition: eos_quark_cfl.h:208
o2scl::eos_quark_cfl6::w6
gsl_eigen_hermv_workspace * w6
GSL workspace for the eigenvalue computation.
Definition: eos_quark_cfl6.h:315
o2scl::eos_quark_cfl6::KD
double KD
The color superconducting 't Hooft coupling (default 0)
Definition: eos_quark_cfl6.h:262

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).