Equation of State Sub-Library: Version 0.910
cfl6_eos.h
00001 /*
00002   -------------------------------------------------------------------
00003   
00004   Copyright (C) 2006-2012, 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 CFL6_EOS_H
00024 #define CFL6_EOS_H
00025 
00026 #include <iostream>
00027 #include <o2scl/test_mgr.h>
00028 #include <o2scl/cfl_njl_eos.h>
00029 
00030 #ifndef DOXYGENP
00031 namespace o2scl {
00032 #endif
00033 
00034   /** \brief An EOS like \ref cfl_njl_eos but 
00035       with a color-superconducting 't Hooft interaction
00036 
00037       Beginning with the Lagrangian:
00038       \f[
00039       {\cal L} = {\cal L}_{Dirac} + {\cal L}_{NJL} + 
00040       {\cal L}_{'t Hooft} + {\cal L}_{SC} + {\cal L}_{SC6} 
00041       \f]
00042       \f[
00043       {\cal L}_{Dirac} = {\bar q} \left( i \partial -m - 
00044       \mu \gamma^0 \right) q
00045       \f]
00046       \f[
00047       {\cal L}_{NJL} = G_S \sum_{a=0}^8 
00048       \left[ \left( {\bar q} \lambda^a q \right)^2
00049       - \left( {\bar q} \lambda^a \gamma^5 q \right)^2 \right]
00050       \f]
00051       \f[
00052       {\cal L}_{'t Hooft} = G_D \left[
00053       \mathrm{det}_f {\bar q} \left(1-\gamma^5 \right) q
00054       +\mathrm{det}_f {\bar q} \left(1+\gamma^5 \right) q
00055       \right]
00056       \f]
00057       \f[
00058       {\cal L}_{SC} = G_{DIQ}
00059       \left( {\bar q}_{i \alpha} i \gamma^5 
00060       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00061       q^C_{j \beta} \right)
00062       \left( {\bar q}_{\ell \delta} i \gamma^5 
00063       \epsilon^{\ell m k} 
00064       \epsilon^{\delta \varepsilon \gamma}
00065       q^C_{m \varepsilon} \right)
00066       \f]
00067       \f[
00068       {\cal L}_{SC6} = K_D
00069       \left( {\bar q}_{i \alpha} i \gamma^5 
00070       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00071       q^C_{j \beta} \right)
00072       \left( {\bar q}_{\ell \delta} i \gamma^5 
00073       \epsilon^{\ell m n} 
00074       \epsilon^{\delta \varepsilon \eta}
00075       q^C_{m \varepsilon} \right)
00076       \left( {\bar q}_{k \gamma} q_{n \eta} \right)
00077       \f]
00078 
00079       We can simplify the relevant terms in \f${\cal L}_{NJL}\f$:
00080       \f[
00081       {\cal L}_{NJL} = G_S \left[ 
00082       \left({\bar u} u\right)^2+
00083       \left({\bar d} d\right)^2+
00084       \left({\bar s} s\right)^2
00085       \right]
00086       \f]
00087       and in \f${\cal L}_{'t Hooft}\f$:
00088       \f[
00089       {\cal L}_{NJL} = G_D \left( 
00090       {\bar u} u {\bar d} d {\bar s} s
00091       \right)
00092       \f]
00093 
00094       Using the definition:
00095       \f[
00096       \Delta^{k \gamma} =  \left< {\bar q} i \gamma^5 
00097       \epsilon \epsilon q^C_{} \right>
00098       \f]
00099       and the ansatzes:
00100       \f[
00101       ({\bar q}_1 q_2) ({\bar q}_3 q_4) \rightarrow
00102       {\bar q}_1 q_2 \left< {\bar q}_3 q_4 \right>
00103       +{\bar q}_3 q_4 \left< {\bar q}_1 q_2 \right>
00104       -\left< {\bar q}_1 q_2 \right> \left< {\bar q}_3 q_4 \right>
00105       \f]
00106       \f[
00107       ({\bar q}_1 q_2) ({\bar q}_3 q_4) ({\bar q}_5 q_6) \rightarrow
00108       {\bar q}_1 q_2 \left< {\bar q}_3 q_4 \right> 
00109       \left< {\bar q}_5 q_6 \right>
00110       +{\bar q}_3 q_4 \left< {\bar q}_1 q_2 \right>
00111       \left< {\bar q}_5 q_6 \right>
00112       +{\bar q}_5 q_6 \left< {\bar q}_1 q_2 \right>
00113       \left< {\bar q}_3 q_4 \right>
00114       -2\left< {\bar q}_1 q_2 \right> \left< {\bar q}_3 q_4 \right>
00115       \left< {\bar q}_5 q_6 \right>
00116       \f]
00117       for the mean field approximation, we can rewrite the Lagrangian
00118       \f[
00119       {\cal L}_{NJL} = 2 G_S \left[ 
00120       \left( {\bar u} u \right) \left< {\bar u} u \right>
00121       +\left( {\bar d} d \right) \left< {\bar d} d \right>
00122       +\left( {\bar s} s \right) \left< {\bar s} s \right>
00123       - \left< {\bar u} u \right>^2
00124       - \left< {\bar d} d \right>^2
00125       - \left< {\bar s} s \right>^2
00126       \right]
00127       \f]
00128       \f[
00129       {\cal L}_{'t Hooft} = - 2 G_D \left[
00130       \left( {\bar u} u \right) \left< {\bar u} u \right>
00131       \left< {\bar s} s \right>
00132       + \left( {\bar d} d \right) \left< {\bar u} u \right>
00133       \left< {\bar s} s \right>
00134       + \left( {\bar s} s \right) \left< {\bar u} u \right>
00135       \left< {\bar d} d \right>
00136       - 2 \left< {\bar u} u \right>\left< {\bar d} d \right>
00137       \left< {\bar s} s \right>
00138       \right]
00139       \f]
00140       \f[
00141       {\cal L}_{SC} = G_{DIQ} \left[
00142       \Delta^{k \gamma}
00143       \left( {\bar q}_{\ell \delta} i \gamma^5 
00144       \epsilon^{\ell m k} 
00145       \epsilon^{\delta \varepsilon \gamma}
00146       q^C_{m \varepsilon} \right)
00147       + \left( {\bar q}_{i \alpha} i \gamma^5 
00148       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00149       q^C_{j \beta} \right)
00150       \Delta^{k \gamma \dagger}
00151       - \Delta^{k \gamma}
00152       \Delta^{k \gamma \dagger}
00153       \right]
00154       \f]
00155       \f[
00156       {\cal L}_{SC6} = K_D \left[
00157       \left( {\bar q}_{m \varepsilon} q_{n \eta} \right)
00158       \Delta^{k \gamma} \Delta^{m \varepsilon \dagger}
00159       + \left( {\bar q}_{i \alpha} i \gamma^5 
00160       \varepsilon^{i j k} \varepsilon^{\alpha \beta \gamma} 
00161       q^C_{j \beta} \right)
00162       \Delta^{m \varepsilon \dagger} 
00163       \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
00164       \right]
00165       \f]
00166       \f[
00167       + K_D \left[\Delta^{k \gamma}
00168       \left( {\bar q}_{\ell \delta} i \gamma^5 
00169       \epsilon^{\ell m n} 
00170       \epsilon^{\delta \varepsilon \eta}
00171       q^C_{m \varepsilon} \right)
00172       \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
00173       -2 
00174       \Delta^{k \gamma} \Delta^{m \varepsilon \dagger}
00175       \left< {\bar q}_{m \varepsilon} q_{n \eta} \right>
00176       \right]
00177       \f]
00178       
00179       If we make the definition \f$ {\tilde \Delta} =
00180       2 G_{DIQ} \Delta \f$
00181 
00182       \hline
00183       <b>References:</b>
00184 
00185       Created for \ref Steiner05.
00186   */
00187   class cfl6_eos : public cfl_njl_eos {
00188   public:
00189   
00190     cfl6_eos();
00191 
00192     virtual ~cfl6_eos();
00193   
00194     /** \brief Calculate the EOS
00195         \nothing
00196       
00197         Calculate the EOS from the quark condensates. Return the mass
00198         gap equations in \c qq1, \c qq2, \c qq3, and the normal gap
00199         equations in \c gap1, \c gap2, and \c gap3.
00200         
00201         Using \c fromqq=true as in nambujl_eos and nambujl_temp_eos
00202         does not work here and will return an error.
00203         
00204         If all of the gaps are less than gap_limit, then the
00205         nambujl_temp_eos::calc_temp_p() is used, and \c gap1, \c gap2,
00206         and \c gap3 are set to equal \c u.del, \c d.del, and \c s.del,
00207         respectively.
00208       
00209     */
00210     virtual int calc_eq_temp_p(quark &u, quark &d, quark &s, 
00211                             double &qq1, double &qq2, double &qq3, 
00212                             double &gap1, double &gap2, double &gap3, 
00213                             double mu3, double mu8,
00214                             double &n3, double &n8, thermo &qb, 
00215                             double temper);
00216 
00217     /// The momentum integrands
00218     virtual int integrands(double p, double res[]);
00219     
00220     /// Check the derivatives specified by eigenvalues()
00221     virtual int test_derivatives(double lmom, double mu3, double mu8,
00222                                  test_mgr &t);
00223 
00224     /** \brief Calculate the energy eigenvalues and their derivatives
00225         \nothing
00226 
00227         Given the momentum \c mom, and the chemical potentials
00228         associated with the third and eighth gluons (\c mu3 and \c mu8),
00229         this computes the eigenvalues of the inverse propagator and
00230         the assocated derivatives.
00231         
00232         Note that this is not the same as cfl_njl_eos::eigenvalues()
00233         which returns \c dedmu rather \c dedqqu.
00234     */
00235     virtual int eigenvalues6(double lmom, double mu3, double mu8, 
00236                              double egv[36], double dedmuu[36], 
00237                              double dedmud[36], double dedmus[36], 
00238                              double dedmu[36], double dedmd[36], 
00239                              double dedms[36], double dedu[36], 
00240                              double dedd[36], double deds[36], 
00241                              double dedmu3[36], double dedmu8[36]);
00242 
00243     /** \brief Construct the matrices, but don't solve the eigenvalue 
00244         problem
00245         \nothing
00246         
00247         This is used by check_derivatives() to make sure that the derivative
00248         entries are right.
00249      */
00250     virtual int make_matrices(double lmom, double mu3, double mu8, 
00251                               double egv[36], double dedmuu[36], 
00252                               double dedmud[36], double dedmus[36], 
00253                               double dedmu[36], double dedmd[36], 
00254                               double dedms[36], double dedu[36], 
00255                               double dedd[36], double deds[36], 
00256                               double dedmu3[36], double dedmu8[36]);
00257     
00258     /// The color superconducting 't Hooft coupling (default 0)
00259     double KD;
00260 
00261     /// Return string denoting type ("cfl6_eos")
00262     virtual const char *type() { return "cfl6_eos"; };
00263 
00264     /** \brief The absolute value below which the CSC 't Hooft coupling 
00265         is ignored(default \f$ 10^{-6} \f$)
00266     */
00267     double kdlimit;
00268 
00269   protected:
00270 
00271 #ifndef DOXYGEN_INTERNAL
00272 
00273     /// Set the quark effective masses from the gaps and the condensates
00274     int set_masses();
00275     
00276     /// The size of the matrix to be diagonalized
00277     static const int mat_size=36;
00278     
00279     /// Storage for the inverse propagator
00280     omatrix_cx iprop6;
00281 
00282     /// The eigenvectors
00283     omatrix_cx eivec6;
00284 
00285     /// The derivative wrt the ds gap
00286     omatrix_cx dipdgapu;
00287 
00288     /// The derivative wrt the us gap
00289     omatrix_cx dipdgapd;
00290 
00291     /// The derivative wrt the ud gap
00292     omatrix_cx dipdgaps;
00293 
00294     /// The derivative wrt the up quark condensate
00295     omatrix_cx dipdqqu;
00296 
00297     /// The derivative wrt the down quark condensate
00298     omatrix_cx dipdqqd;
00299 
00300     /// The derivative wrt the strange quark condensate
00301     omatrix_cx dipdqqs;
00302 
00303     /// Storage for the eigenvalues
00304     ovector eval6;
00305 
00306     /// GSL workspace for the eigenvalue computation
00307     gsl_eigen_hermv_workspace *w6;
00308 
00309   private:
00310 
00311     cfl6_eos(const cfl6_eos &);
00312     cfl6_eos& operator=(const cfl6_eos&);
00313 
00314 #endif
00315 
00316   };
00317 
00318 #ifndef DOXYGENP
00319 }
00320 #endif
00321 
00322 #endif
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Friends

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

Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads.