cfl6_eos.h

00001 /*
00002   -------------------------------------------------------------------
00003   
00004   Copyright (C) 2006, 2007, Andrew W. Steiner
00005   
00006   This file is part of O2scl.
00007   
00008   O2scl is free software; you can redistribute it and/or modify
00009   it under the terms of the GNU General Public License as published by
00010   the Free Software Foundation; either version 3 of the License, or
00011   (at your option) any later version.
00012   
00013   O2scl is distributed in the hope that it will be useful,
00014   but WITHOUT ANY WARRANTY; without even the implied warranty of
00015   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00016   GNU General Public License for more details.
00017   
00018   You should have received a copy of the GNU General Public License
00019   along with O2scl. If not, see <http://www.gnu.org/licenses/>.
00020 
00021   -------------------------------------------------------------------
00022 */
00023 #ifndef O2SCL_CFL6_EOS_H
00024 #define O2SCL_CFL6_EOS_H
00025 
00026 #include <iostream>
00027 #include <o2scl/cfl_njl_eos.h>
00028 #include <o2scl/columnify.h>
00029 #include <o2scl/omatrix_cx_tlate.h>
00030 
00031 #ifndef DOXYGENP
00032 namespace o2scl {
00033 #endif
00034 
00035   /** 
00036       \brief CFL NJL EOS with a color-superconducting 't Hooft interaction
00037 
00038       Beginning with the Lagrangian:
00039       \f[
00040       {\cal L} = {\cal L}_{Dirac} + {\cal L}_{NJL} + 
00041       {\cal L}_{'t Hooft} + {\cal L}_{SC} + {\cal L}_{SC6} 
00042       \f]
00043       \f[
00044       {\cal L}_{Dirac} = {\bar q} \left( i \partial -m - 
00045       \mu \gamma^0 \right) q
00046       \f]
00047       \f[
00048       {\cal L}_{NJL} = G_S \sum_{a=0}^8 
00049       \left[ \left( {\bar q} \lambda^a q \right)^2
00050       - \left( {\bar q} \lambda^a \gamma^5 q \right)^2 \right]
00051       \f]
00052       \f[
00053       {\cal L}_{'t Hooft} = G_D \left[
00054       \mathrm{det}_f {\bar q} \left(1-\gamma^5 \right) q
00055       +\mathrm{det}_f {\bar q} \left(1+\gamma^5 \right) q
00056       \right]
00057       \f]
00058       \f[
00059       {\cal L}_{SC} = G_{DIQ}
00060       \left( {\bar q}_{i \alpha} i \gamma^5 
00061       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00062       q^C_{j \beta} \right)
00063       \left( {\bar q}_{\ell \delta} i \gamma^5 
00064       \epsilon^{\ell m k} 
00065       \epsilon^{\delta \varepsilon \gamma}
00066       q^C_{m \varepsilon} \right)
00067       \f]
00068       \f[
00069       {\cal L}_{SC6} = K_D
00070       \left( {\bar q}_{i \alpha} i \gamma^5 
00071       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00072       q^C_{j \beta} \right)
00073       \left( {\bar q}_{\ell \delta} i \gamma^5 
00074       \epsilon^{\ell m n} 
00075       \epsilon^{\delta \varepsilon \eta}
00076       q^C_{m \varepsilon} \right)
00077       \left( {\bar q}_{k \gamma} q_{n \eta} \right)
00078       \f]
00079 
00080       We can simplify the relevant terms in \f${\cal L}_{NJL}\f$:
00081       \f[
00082       {\cal L}_{NJL} = G_S \left[ 
00083       \left({\bar u} u\right)^2+
00084       \left({\bar d} d\right)^2+
00085       \left({\bar s} s\right)^2
00086       \right]
00087       \f]
00088       and in \f${\cal L}_{'t Hooft}\f$:
00089       \f[
00090       {\cal L}_{NJL} = G_D \left( 
00091       {\bar u} u {\bar d} d {\bar s} s
00092       \right)
00093       \f]
00094 
00095       Using the definition:
00096       \f[
00097       \Delta^{k \gamma} =  \left< {\bar q} i \gamma^5 
00098       \epsilon \epsilon q^C_{} \right>
00099       \f]
00100       and the ansatzes:
00101       \f[
00102       ({\bar q}_1 q_2) ({\bar q}_3 q_4) \rightarrow
00103       {\bar q}_1 q_2 \left< {\bar q}_3 q_4 \right>
00104       +{\bar q}_3 q_4 \left< {\bar q}_1 q_2 \right>
00105       -\left< {\bar q}_1 q_2 \right> \left< {\bar q}_3 q_4 \right>
00106       \f]
00107       \f[
00108       ({\bar q}_1 q_2) ({\bar q}_3 q_4) ({\bar q}_5 q_6) \rightarrow
00109       {\bar q}_1 q_2 \left< {\bar q}_3 q_4 \right> 
00110       \left< {\bar q}_5 q_6 \right>
00111       +{\bar q}_3 q_4 \left< {\bar q}_1 q_2 \right>
00112       \left< {\bar q}_5 q_6 \right>
00113       +{\bar q}_5 q_6 \left< {\bar q}_1 q_2 \right>
00114       \left< {\bar q}_3 q_4 \right>
00115       -2\left< {\bar q}_1 q_2 \right> \left< {\bar q}_3 q_4 \right>
00116       \left< {\bar q}_5 q_6 \right>
00117       \f]
00118       for the mean field approximation, we can rewrite the Lagrangian
00119       \f[
00120       {\cal L}_{NJL} = 2 G_S \left[ 
00121       \left( {\bar u} u \right) \left< {\bar u} u \right>
00122       +\left( {\bar d} d \right) \left< {\bar d} d \right>
00123       +\left( {\bar s} s \right) \left< {\bar s} s \right>
00124       - \left< {\bar u} u \right>^2
00125       - \left< {\bar d} d \right>^2
00126       - \left< {\bar s} s \right>^2
00127       \right]
00128       \f]
00129       \f[
00130       {\cal L}_{'t Hooft} = - 2 G_D \left[
00131       \left( {\bar u} u \right) \left< {\bar u} u \right>
00132       \left< {\bar s} s \right>
00133       + \left( {\bar d} d \right) \left< {\bar u} u \right>
00134       \left< {\bar s} s \right>
00135       + \left( {\bar s} s \right) \left< {\bar u} u \right>
00136       \left< {\bar d} d \right>
00137       - 2 \left< {\bar u} u \right>\left< {\bar d} d \right>
00138       \left< {\bar s} s \right>
00139       \right]
00140       \f]
00141       \f[
00142       {\cal L}_{SC} = G_{DIQ} \left[
00143       \Delta^{k \gamma}
00144       \left( {\bar q}_{\ell \delta} i \gamma^5 
00145       \epsilon^{\ell m k} 
00146       \epsilon^{\delta \varepsilon \gamma}
00147       q^C_{m \varepsilon} \right)
00148       + \left( {\bar q}_{i \alpha} i \gamma^5 
00149       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00150       q^C_{j \beta} \right)
00151       \Delta^{k \gamma \dagger}
00152       - \Delta^{k \gamma}
00153       \Delta^{k \gamma \dagger}
00154       \right]
00155       \f]
00156       \f[
00157       {\cal L}_{SC6} = K_D \left[
00158       \left( {\bar q}_{m \varepsilon} q_{n \eta} \right)
00159       \Delta^{k \gamma} \Delta^{m \varepsilon \dagger}
00160       + \left( {\bar q}_{i \alpha} i \gamma^5 
00161       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00162       q^C_{j \beta} \right)
00163       \Delta^{m \varepsilon \dagger} 
00164       \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
00165       \right]
00166       \f]
00167       \f[
00168       + K_D \left[\Delta^{k \gamma}
00169       \left( {\bar q}_{\ell \delta} i \gamma^5 
00170       \epsilon^{\ell m n} 
00171       \epsilon^{\delta \varepsilon \eta}
00172       q^C_{m \varepsilon} \right)
00173       \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
00174       -2 
00175       \Delta^{k \gamma} \Delta^{m \varepsilon \dagger}
00176       \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
00177       \right]
00178       \f]
00179       
00180       If we make the definition \f$ {\tilde \Delta} =
00181       2 G_{DIQ} \Delta \f$
00182 
00183   */
00184   class cfl6_eos : public cfl_njl_eos {
00185   public:
00186   
00187     cfl6_eos();
00188 
00189     virtual ~cfl6_eos();
00190   
00191     /** 
00192         \brief Calculate the EOS
00193       
00194         Calculate the EOS from the quark condensates.  Return the mass
00195         gap equations in \c qq1, \c qq2, \c qq3, and the normal gap
00196         equations in \c gap1, \c gap2, and \c gap3.
00197         
00198         Using \c fromqq=true as in nambujl_eos and nambujl_temp_eos does
00199         not work here and will return an error. The quarks must be set
00200         through quark_eos::quark_set() before use.
00201         
00202         If all of the gaps are less than gap_limit, then the
00203         nambujl_temp_eos::calc_temp_p() is used, and \c gap1, \c gap2,
00204         and \c gap3 are set to equal \c u.del, \c d.del, and \c s.del,
00205         respectively.
00206       
00207     */
00208     virtual int calc_eq_temp_p(quark &u, quark &d, quark &s, 
00209                             double &qq1, double &qq2, double &qq3, 
00210                             double &gap1, double &gap2, double &gap3, 
00211                             double mu3, double mu8,
00212                             double &n3, double &n8, thermo &qb, 
00213                             const double temper);
00214 
00215     /** \brief Direct calculation of the thermodynamic potential
00216      */
00217     virtual double thd_potential(quark &u, quark &d, quark &s,
00218                                  double mu3, double mu8, 
00219                                  const double ltemper);
00220 
00221 
00222     /** 
00223         \brief Calculate the energy eigenvalues as a function of the momentum
00224 
00225         Given the momentum 'mom', and the chemical potentials
00226         associated with the third and eighth gluons ('mu3' and 'mu8'),
00227         the energy eigenvalues are computed in egv[0] ... egv[35].  No
00228         space is allocated for the array by the function.
00229 
00230         \todo This function may make some inappropriate assumptions 
00231         on the vector egv.
00232     */
00233     virtual int eigenvalues(double mom, ovector_view &egv, double mu3, 
00234                             double mu8);
00235 
00236     /// The color superconducting 't Hooft coupling (default 0)
00237     double KD;
00238 
00239     /// Return string denoting type ("cfl6_eos")
00240     virtual const char *type() { return "cfl6_eos"; };
00241 
00242     /** \brief The absolute value below which the CSC 't Hooft coupling 
00243         is ignored(default \f$ 10^{-6} \f$)
00244     */
00245     double kdlimit;
00246 
00247   protected:
00248 
00249 #ifndef DOXYGEN_INTERNAL
00250 
00251     /// To clear all of the matrix entries
00252     gsl_complex zero;
00253     
00254     /// The size of the matrix to be diagonalized
00255     static const int size=36;
00256     
00257     /// Storage for the inverse propagator
00258     omatrix_cx m6;
00259 
00260     /// Storage for the eigenvalues
00261     ovector eval6;
00262 
00263     /// GSL workspace for the eigenvalue computation
00264     gsl_eigen_herm_workspace *w6;
00265 
00266     /** \brief The function used to take derivatives of the thermodynamic
00267         potential (used by calc_eq_temp_p() )
00268     */
00269     virtual double tpot(double var, void *&pa);
00270     
00271   private:
00272 
00273     cfl6_eos(const cfl6_eos &);
00274     cfl6_eos& operator=(const cfl6_eos&);
00275 
00276 #endif
00277 
00278   };
00279 
00280 #ifndef DOXYGENP
00281 }
00282 #endif
00283 
00284 #endif

Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.