109 lines
2.8 KiB
C++
109 lines
2.8 KiB
C++
/**********************************************************
|
|
|
|
This software is part of J.-S. Caux's ABACUS library.
|
|
|
|
Copyright (c) J.-S. Caux.
|
|
|
|
-----------------------------------------------------------
|
|
|
|
File: LiebLin_Utils.cc
|
|
|
|
Purpose: Utilities for Lieb-Liniger gas
|
|
|
|
***********************************************************/
|
|
|
|
#include "ABACUS.h"
|
|
|
|
using namespace std;
|
|
|
|
namespace ABACUS {
|
|
|
|
DP LiebLin_dE0_dc (DP c_int, DP L, int N)
|
|
{
|
|
LiebLin_Bethe_State groundstate (c_int, L, N);
|
|
groundstate.Compute_All(true);
|
|
|
|
DP dc_int = 1.0e-5 * c_int;
|
|
LiebLin_Bethe_State groundstate_2 (c_int + dc_int, L, N);
|
|
groundstate_2.Compute_All(true);
|
|
|
|
return((groundstate_2.E - groundstate.E)/dc_int);
|
|
}
|
|
|
|
DP LiebLin_vs (DP c_int, DP L, int N)
|
|
{
|
|
LiebLin_Bethe_State gstate (c_int, L, N);
|
|
gstate.Compute_All(true);
|
|
|
|
DP Egs = gstate.E;
|
|
|
|
gstate.Ix2[N-1] += 2;
|
|
|
|
gstate.Compute_All(false);
|
|
|
|
return((gstate.E - Egs) * L/twoPI);
|
|
|
|
}
|
|
|
|
DP LiebLin_Dressed_Charge_N (DP c_int, DP L, int N)
|
|
{
|
|
LiebLin_Bethe_State gstate(c_int, L, N);
|
|
gstate.Compute_All(true);
|
|
return(twoPI/(c_int * L * (gstate.lambdaoc[N-1] - gstate.lambdaoc[N-2]))); // lambda in units of c within code
|
|
}
|
|
|
|
|
|
// For computation of threshold exponents:
|
|
|
|
int Momentum_Right_Excitations (LiebLin_Bethe_State& ScanState)
|
|
{
|
|
// Calculates the momentum of the excitations on the right of Fermi
|
|
// surface, discarding the rightmost excitation.
|
|
|
|
int nr = 0;
|
|
for (int i = ScanState.N - 2; i >= ScanState.N/2; --i)
|
|
if (ScanState.Ix2[i] >= 0) nr += (ScanState.Ix2[i] + ScanState.N - 1 - 2*i)/2;
|
|
|
|
return(nr);
|
|
}
|
|
|
|
int Momentum_Left_Excitations (LiebLin_Bethe_State& ScanState)
|
|
{
|
|
// Calculates the momentum of the excitations on the left of Fermi
|
|
// surface, discarding the rightmost excitation.
|
|
|
|
int nl = 0;
|
|
for (int i = 0; i < ScanState.N/2; ++i)
|
|
if (ScanState.Ix2[i] < 0) nl -= (ScanState.Ix2[i] + ScanState.N - 1 - 2*i)/2;
|
|
|
|
return(nl);
|
|
}
|
|
|
|
|
|
DP ln_Overlap_with_BEC (LiebLin_Bethe_State& lambda) // Note: factorial approximated with Lanczos
|
|
{
|
|
int N = lambda.N;
|
|
int L = lambda.L;
|
|
|
|
DP c_int = lambda.c_int;
|
|
|
|
// The overlap identically vanishes if the state is not parity invariant:
|
|
if (!lambda.Check_Symmetry()) return(-300.0);
|
|
|
|
SQMat_DP Gaudin( N);
|
|
|
|
SQMat_DP Gaudin_Quench( N/2);
|
|
|
|
lambda.Build_Reduced_Gaudin_Matrix ( Gaudin);
|
|
lambda.Build_Reduced_BEC_Quench_Gaudin_Matrix (Gaudin_Quench);
|
|
|
|
DP ln_prefactor = N/2.0 * log(2./(c_int*L) ) + 0.5 * real(ln_Gamma(complex<double>(N+1.0)));
|
|
|
|
for(int i =N/2; i<N; i ++)
|
|
ln_prefactor -= log(fabs(lambda.lambdaoc[i ])) + 0.5*log( 1. + 4. * lambda.lambdaoc[i ]*lambda.lambdaoc[i ]) ;
|
|
|
|
return (ln_prefactor + real(lndet_LU_dstry(Gaudin_Quench)) - 0.5 * real(lndet_LU_dstry(Gaudin)));
|
|
}
|
|
|
|
}
|