A convenient implementation of the overlap between arbitrary Hartree-Fock-Bogoliubov vacua for projection
aa r X i v : . [ nu c l - t h ] A ug A convenient implementation of the overlap between arbitraryHartree-Fock-Bogoliubov vacua for projection
Zao-Chun Gao a, ∗ , Qing-Li Hu a , Y. S. Chen a a China Institute of Atomic Energy, P.O. Box 275 (10), Beijing 102413, P.R. China
Abstract
Overlap between Hartree-Fock-Bogoliubov(HFB) vacua is very important in the beyond mean-field calculations.However, in the HFB transformation, the U , V matrices are sometimes singular due to the exact emptiness ( v i = u i =
0) of some single-particle orbits. This singularity may cause some problem in evaluating theoverlap between HFB vacua through Pfa ffi an. We found that this problem can be well avoided by setting those zerooccupation numbers u i , v i to some tiny values denoted by ε ( > + ε = ε = − when using the double precision data type). This treatment does not change the HFB vacuum state because u i , v i = ε are numerically zero relative to 1. Therefore, for arbitrary HFB transformation, we say that the U , V matrices canalways be nonsingular. From this standpoint, we present a new convenient Pfa ffi an formula for the overlap betweenarbitrary HFB vacua, which is especially suitable for symmetry restoration. Testing calculations have been performedfor this new formula. It turns out that our method is reliable and accurate in evaluating the overlap between arbitraryHFB vacua. Keywords:
Hartree-Fock-Bogoliubov method, beyond mean-field method, Pfa ffi an
1. Introduction
The Hartree-Fock-Bogoliubov (HFB) approximationhas been a great success in understanding interactingmany-body quantum systems in all fields of physics.However, the beyond mean-field e ff ects (e.g., the nu-clear vibration and rotation) are missing in the HFB cal-culations. Methods that go beyond mean-field, such asthe Generator Coordinate Method(GCM) and the pro-jection method, are expected to take those missing ef-fects into consideration and present better description ofthe many-body quantum system. In the beyond mean-field calculations, operator matrix elements and over-laps between multi-quasiparticle HFB states are basicblocks. These matrix elements and overlaps can be eval-uated using the generalized Wick’s theorem (GWT)[1,2], or equivalently using Pfa ffi an [3, 4, 5, 6, 7], or us-ing the compact formula in Ref.[8]. However, in thee ffi cient calculations (e.g., see [5]), all of the matrix el-ements and overlaps require the value of the overlap be-tween HFB vacua. ∗ Corresponding author
Email address: [email protected] (Zao-Chun Gao)
Thus, the reliable and accurate evaluation of the over-lap between HFB vacua is very important for the stabil-ity and the e ffi ciency of the beyond mean-field calcu-lations. Especially in cases near to the Egido pole [9],the overlap between HFB vacua is very tiny, and a smallerror could lead to a large uncertainty of the matrix ele-ments. In the past, numerical calculations of the overlapwere performed with the Onishi formula [10]. Unfortu-nately, the Onishi formula leaves the sign of the overlapundefined due to the square root of a determinant. Sev-eral e ff orts have been made to overcome this sign prob-lem [11, 12, 13, 14, 15, 16]. In 2009, Robledo proposeda di ff erent overlap formula with the Pfa ffi an rather thanthe determinant [17]. This formula completely solvesthe sign problem but requires the inversion of the ma-trix U in the Bogoliubov transformation. To avoid thesingularity of U , the formula for the limit when severalorbits are fully occupied is given in Ref. [18]. Simul-taneously, the limit when some orbits are exact emptywas also considered to reduce the computational cost.Meanwhile, various Pfa ffi an formulae for the overlapbetween HFB vacua have been proposed by several au-thors [3, 6, 7]. In Ref. [7], the overlap formula does notrequire the inversion of U , but the empty orbits in theFock space should be omitted. Preprint submitted to Physics Letter B November 8, 2018 n practical calculation, one should first identify thesingularity of the matrices U and V in the Bogoli-ubov transformation. This can be easily tested with theBloch-Messiah theorem (see details in Ref. [19]). Thematrices U and V can be decomposed as U = D ¯ UC and V = D ∗ ¯ VC . Here, D and C are unitary matrices. ¯ U and¯ V refer to the BCS-transformation and are constructedfrom the occupation numbers u i , v i with 0 ≤ u i , v i ≤ u i + v i = u i = , v i = u i = , v i =
0) levels have been care-fully treated in Refs. [6, 7, 18] to avoid the collapse ofthe overlap computation.However, we note that in most realistic cases the v i ’scan be extremely close to 0 or 1 but not exact 0 or 1.Strictly speaking, these levels with such extreme empti-ness or occupation should be considered but may leadto exotic values ( extremely huge or extremely tiny) ofthe Pfa ffi an in the proposed formulae. What is worse,the Pfa ffi an values are easily out of the scope of thedouble precision data type and cause the computationcollapsed.Careful treatment must be made to avoid such dataoverflow. In this paper, we implement an accurateand reliable calculation for the overlap between arbi-trary HFB vacua in a unified way. For the cases of( u i = , v i =
0) and ( u i = , v i = u i = , v i = ε ) and ( u i = ε, v i = ε > ε should be numerical zero relative to 1 in the practi-cal calculation. In other words, ε should numericallysatisfy 1 + ε =
1. Under this condition, ε may be cho-sen as large as possible so that the calculated Pfa ffi anvalues are not necessarily too huge or too tiny. For in-stance, one can choose ε = − when using double pre-cision. Because v i ( u i ) = ε is actually zero relative to u i ( v i ) = u i , v i values are allowed to be extremely close to 0 or 1. Ide-ally, U , V are nonsingular in our assumption, and we canderive a new formula for the overlap between the HFBvacua based on the work of Bertsch and Robledo [7].This formula is especially convenient for the symmetryrestoration. Numerical calculations have been carriedout for heavy nuclear system to test the precision of thenew formula by comparing with the Onishi formula.In section 2, the formalism of the new overlap for-mula is given. Section 3 provides an example of numer-ical calculation. A summary is given in section 4.
2. The overlap between the HFB vacua
We denote ˆ c † i and ˆ c i as the creation and annihilationoperators defined in an M -dimensional Fock-space. TheHartree-Fock-Bogoliubov(HFB) transformation is ˆ β ˆ β † ! = U † V † V T U T ! ˆ c ˆ c † ! . (1)Here, we assume U and V are nonsingular matrices, andtheir shapes are M × M . The HFB vacuum (unnormal-ized) can be written as | φ i = ˆ β ˆ β ... ˆ β M |−i , (2)where |−i is the true vacuum. By definition, one hasˆ β i | φ i = ≤ i ≤ M . (3)The second HFB vacuum | φ ′ i is defined in the sameway, but the prime, ‘ ′ ’, is attached to the correspond-ing symbols to show di ff erence.The overlap between | φ i and | φ ′ i is given by h φ | φ ′ i = h−| ˆ β † M ˆ β † M − ... ˆ β † ˆ β ′ ˆ β ′ ... ˆ β ′ M |−i = s M h−| ˆ β † ˆ β † ... ˆ β † M ˆ β ′ ˆ β ′ ... ˆ β ′ M |−i , (4)where, s M = ( − [ M ( M − / . If M is even, s M = ( − M / .Following the technique of Bertsch and Robledo [7],one can obtain h φ | φ ′ i = s M pf V T U V T V ′∗ − V ′† V U ′† V ′∗ ! . (5)The shape of the matrix in Eq. (5) is 2 M × M , and noempty levels are omitted. For the norm overlap h φ | φ i , itis real and positive. From Eq.(5) and the Bloch-Messiahtheorem, one can get h φ | φ i = s M pf V T U V T V ∗ − V † V U † V ∗ ! = M / Y i = v i . (6)Denoting Q M / i = v i by N , the normalized quasi-particle vacuum, | ψ i , can be written as | ψ i = | φ i N . (7)Then, one finds that h ψ | ψ ′ i = s M NN ′ pf V T U V T V ′∗ − V ′† V U ′† V ′∗ ! . (8)In the symmetry restoration, the general rotationaloperator, involving the spin and particle number projec-tion, may be written asˆ R ( Ξ ) = ˆ R ( Ω ) e − i ˆ N φ n e − i ˆ Z φ p , (9)2here ˆ R ( Ω ) is the rotation operator, and Ω refers to thethree Euler angles α, β, γ . e − i ˆ N φ n and e − i ˆ Z φ p are ‘gauge’rotational operators induced by the neutron and protonnumber projection. ˆ N and ˆ Z are neutron and protonnumber operators, respectively. φ n and φ p are ”gauge”angles for neutron and proton, respectively. Ξ refers to( Ω , φ n , φ p ). The matrix element h ψ | ˆ R ( Ξ ) | ψ ′ i needs to becalculated. Let’s define the general rotation transforma-tion for symmetry restoration,ˆ R ( Ξ ) ˆ c ˆ c † ! ˆ R † ( Ξ ) = D † ( Ξ ) 00 D T ( Ξ ) ! ˆ c ˆ c † ! , (10)where D i j ( Ξ ) = h i | ˆ R ( Ξ ) | j i , and | i ( j ) i = ˆ c † i ( j ) |−i . The D ( Ξ ) matrix has the dimension M × M . One can getˆ R ( Ξ ) ˆ β ′ ˆ β ′† ! ˆ R † ( Ξ ) = D ( Ξ ) ˆ c ˆ c † ! , (11)where D ( Ξ ) = [ D ( Ξ ) U ′ ] † [ D ∗ ( Ξ ) V ′ ] † [ D ∗ ( Ξ ) V ′ ] T [ D ( Ξ ) U ′ ] T ! . (12)By comparing Eq.(11) with Eq.(1), one can obtain therotated overlap by replacing U ′ and V ′ in Eq.(8) with D ( Ξ ) U ′ and D ∗ ( Ξ ) V ′ , respectively. Thus N pf ( Ξ ) = h ψ | ˆ R ( Ξ ) | ψ ′ i = s M NN ′ pf[ M ( Ξ )] , (13)where M ( Ξ ) = V T U V T D ( Ξ ) V ′∗ − V ′† D T ( Ξ ) V U ′† V ′∗ ! . (14)This formula is essentially the same as the one proposedby Bertsch and Robledo [7], but we will transform itinto a new form. Supposing that there is a Ξ satisfying N pf ( Ξ ) ,
0, we have N pf ( Ξ ) N pf ( Ξ ) = pf[ M ( Ξ )]pf[ M ( Ξ )] = pf h P M ( Ξ ) P T i pf (cid:2) P M ( Ξ ) P T (cid:3) = pf[ W ( Ξ )]pf[ W ( Ξ )] , (15)where W ( Ξ ) = [ U ′ V ′− ] † − D T ( Ξ ) D ( Ξ ) UV − ! , (16)and P is P = V ′† ) − ( V T ) − ! . (17)Therefore, one can get N pf ( Ξ ) = C pf[ W ( Ξ )] , (18) where, the coe ffi cient C is actually independent of Ξ ,and can be written as C = N pf ( Ξ )pf[ W ( Ξ )] = s M NN ′ det P = s M ∆ NN ′ . (19)Here, ∆ is a phase determined by ∆ = det D ∗ det D ′ det C det C ′∗ . (20)In Eq.(18), we have used the Bloch-Messiah theoremand the following equationdet P = det[( V ′† ) − ]det[( V T ) − ] . (21)Eq.(18) looks more convenient to be implementedand may save some computing time in contrast toEq.(13), where extra evaluation of V T D ( Ξ ) V ′∗ is re-quired for each mesh point in the integral of projection.For comparison, let us present a brief introductionof the overlap of the Onishi formula [10]. The uni-tary transformation of the quasi-particles under rotationˆ R ( Ξ ) can be written asˆ R ( Ξ ) ˆ β ′ ˆ β ′† ! ˆ R † ( Ξ ) = X ( Ξ ) Y ( Ξ ) Y ∗ ( Ξ ) X ∗ ( Ξ ) ! ˆ β ˆ β † ! , (22)where X ( Ξ ) = U ′† D † ( Ξ ) U + V ′† D T ( Ξ ) V , Y ( Ξ ) = U ′† D † ( Ξ ) V ∗ + V ′† D T ( Ξ ) U ∗ . (23)The Onishi formula is then expressed as (see Ref.[20]), N Onishi ( Ξ ) = h ψ | ˆ R ( Ξ ) | ψ ′ i = ( ± ) p det[ X ( Ξ )] e − i ( M n φ n + M p φ p ) / , (24)where M n and M p are the numbers of neutron and pro-ton orbits in the Fock space, respectively. The valueof det[ X ( Ξ )] is a complex number, and the sign of thesquare root is left undefined. Extra e ff orts must be madeto determine the sign before the application of the On-ishi formula. For instance, in the Projected Shell Model[21] without particle number projection, the overlap be-tween the BCS vacua is real and positive, thus there isno sign ambiguity and the Onishi formula works.
3. Numerical test of the overlap formulae
Although the sign problem is solved in Eq.(13) andEq.(18), one can imagine that N , N ′ are extremely tinynumbers by definition. Thus pf[ M ( Ξ )] is also very tiny,but pf[ W ( Ξ )] should be huge. Numerical accuracy ofEqs. (13) and (18) needs to be carefully tested. It isbelieved that the Onishi formula is accurate except forits undetermined sign. So, it is helpful to compare the3umerical values of the overlaps using Eq.(13), Eq.(18)and Eq.(24).To demonstrate the accuracy and the reliability ofthe Eqs. (13) and (18), numerical calculations are per-formed for the typical example of the deformed heavynucleus Th. For projection, we should take | ψ i = | ψ ′ i , and then ∆ = U , V matrices are obtained from the Nils-son + BCS method. The single particle levels are gen-erated from the Nilsson Hamiltonian with the standardparameters [22]. The single-particle model space con-tains 5 neutron major shells with N = N = M n = M p = ǫ = .
2. Here, we only con-sider the axial symmetry for simplicity.In the no pairing case, the BCS vacuum becomes apure slater determinant, which is a challenge for Eq.(18)because all v i ’s above the Fermi surface are zero. Con-sequently, N = W ( Ξ ) is meaningless due to thesingularity of V . Here, we use the double precision datatype and set v i = ε = − for those v i = h φ n | φ n i = (10 − ) − = − , h φ p | φ p i = (10 − ) − = − , where, | φ n i and | φ p i are BCS vacua for neutrons andprotons, respectively, and | φ i = | φ n i| φ p i . The tiny num-bers 10 − and 10 − are too far out of the scope ofthe double precision data ( ∼ ± ). To avoid the dataoverflow, we multiply the tiny variable by 10 severaltimes until the scaled absolute value falls into the in-terval [10 − , ]. In other words, we use a number y and an integer number k to express a tiny number x through x = y × (10 − ) k . If x is a huge number, then k is negative.However, for the Onishi formula of Eq.(24), we donot need to change v i = v i = ε . The overlaps for theneutron part, calculated with Eq.(18) and Eq.(24), arecompared in Fig.1. The curves of Eq.(18) are continu-ous, but the sign uncertainty of Eq.(24) causes the dis-continuity. However, if one copies the sign of Eq.(18) toEq.(24), one can compare numerical di ff erence betweenEq.(18) and Eq.(24) using the following quantity, R , R = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) N Onishi ( φ n ) N pf ( φ n ) − (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) . (25)In all calculations, we found that R < − with dou-ble precision. This confirms that a small change of v i I m [ ( n ) ] n (in unit of )-0.4-0.20.00.20.4 R e [ ( n ) ] Onishi Formula Present Formula
Figure 1: (Color online) Overlaps of the ground state neutron slaterdeterminant for
Th as functions of φ n with Euler angles α = γ = ◦ , β = ◦ , calculated with present formula [Eq.(18)] and the Onishiformula [Eq.(24) with ‘ + ’ sign]. Re[ N ( φ n )] and Im[ N ( φ n )] are the realand imaginary parts of the overlap. from zero to ε almost does not a ff ect numerical accu-racy. However, it is crucial to keep Eq.(18) valid. Yetnotice that N pf ( Ξ ) in Eq.(18) is obtained from a productof tiny and giant numbers. The same calculations havealso been done with Eq.(13), and we also get R < − .Thus we have presented an alternative way of usingEq.(13), where we set v i = ε for those empty orbitsrather than omitting them[7].Once the overlap is available, it is straightforward toperform the symmetry restoration. The deformed BCSvacuum of Th has been projected onto good particlenumber and spin. Therefore, one can test how precisethe numerical calculations with Eq.(18) satisfy X N , Z , I h ψ | ˆ P N ˆ P Z ˆ P I | ψ i = , (26)where ˆ P N , ˆ P Z , and ˆ P IMK are neutron-number, proton-number, and spin projection operators, respectively. Forthe above vacuum state without pairing (i.e. the groundstate slater determinant), the particle numbers of bothneutrons and protons are good. Indeed, our particlenumber projection (using 16 mesh points in the integral)shows that h ψ n | ˆ P N | ψ n i = N = N , − . Calculations forthe protons also have the same accuracy. This againshows the reliability of Eq. (18). Angular momentumprojection is also performed on the same state in addi-tion to the particle number projection. The amplitudeof h ψ | ˆ P N ˆ P Z ˆ P I | ψ i with ( N = , Z =
70) is plotted asa function of spin I in Fig.2. In the integral of the spinprojection, 100 mesh points are taken, and the range ofspin is 0 ≤ I ≤
70, and we indeed reproduced Eq.(26)with numerical error around 10 − .4
20 40 60 8010 -12 -10 -8 -6 -4 -2 A m p li t ude Spin I Figure 2: The amplitude of projection, h ψ | ˆ P N ˆ P Z ˆ P I | ψ i , as a functionof spin I at N =
96 and Z =
70 using Eq.(18). | ψ i is the axiallydeformed BCS vacuum but without pairing. We also have tested Eq.(18) in the projection of thetriaxially deformed vacuum with normal pairing, whichseems more convenient to use Eq.(18). With the presentmethod, similar accuracy has also been achieved.
4. Summary
Following the strategy of Bertsch and Robledo [7],we have proposed a new formula of the overlap betweenHFB vacua by using the Pfa ffi an identity and assumingthat the inverse of the V matrix exists. This formulais especially convenient and e ffi cient in the symmetryrestoration, and has the same high accuracy as the On-ishi formula as well as the correct sign. The reliabilityof the present formula has been tested by carrying outthe calculations of the overlap and the quantum num-ber projection for the heavy nucleus Th. In the test-ing calculations, one has to be faced with two numericalproblems: (1) The extreme (huge or tiny) quantities arecertainly encountered, and we have properly treated thissituation to avoid data overflow (see the text). (2) Forthose empty orbits with v i =
0, which make Eq.(18) in-valid, one can change v i to a small quantity ε ( >
0) toavoid the singularity of V matrix. It turns out that suchtreatments work very well. Testing calculations haveconfirmed that the present formula is even applicable tothe pure slater determinant without losing the numeri-cal accuracy. Thus it is promising that Eq.(18) may beapplicable in evaluating the overlap between arbitraryHFB vacua. Acknowledgements
Z. G. thanks Prof. Y. Sun andDr. F. Q. Chen for the stimulating and fruitful discus-sions. The authors acknowledge support from the Na- tional Natural Science Foundation of China under Con-tract Nos. 11175258, 11021504 and 11275068.
References [1] R. Balian and E. Brezin, Nuovo Cimento B 64 (1969) 37.[2] K. Hara, S. Iwasaki, Nucl. Phys. A 332 (1979) 61.[3] M. Oi, T. Mizusaki, Phys. Lett. B 707 (2012) 305.[4] T. Mizusaki, M. Oi, Phys. Lett. B 715 (2012) 219.[5] T. Mizusaki, M. Oi, Fang-Qi Chen, Yang Sun, Phys. Lett. B 725(2013) 175.[6] B. Avez, M. Bender, Phys. Rev. C. 85 (2012) 034325.[7] G.F. Bertsch and L.M. Robledo, Phys. Rev. Lett. 108 (2012)042505.[8] S. Perez-Martin and L.M. Robledo, Phys. Rev. C 76 (2007)064314.[9] M. Anguiano, J.L. Egido, L.M. Robledo, Nucl. Phys. 696 (2001)467.[10] N. Onishi and S. Yoshida, Nucl. Phys. 80 (1966) 367.[11] K. Neergård and E. W¨ust, Nucl. Phys. A 402 (1983) 311.[12] Q. Haider and D. Gogny, J. Phys. G: Nucl. Part. Phys. 18 (1992)993.[13] F. D¨onau, Phys. Rev. C 58 (1998) 872.[14] M. Bender and Paul-Henri Heenen, Phys. Rev. C 78 (2008)024309.[15] K. Hara, A. Hayashi and P. Ring, Nucl. Phys. A 385 (1982) 14.[16] M. Oi and N. Tajima, Phys. Lett. B 606 (2005) 43.[17] L.M. Robledo, Phys. Rev. C 79 (2009) 021302(R).[18] L.M. Robledo, Phys. Rev. C 84 (2011) 014307.[19] P. Ring and P. Schuck, The Nuclear Many-Body Problem,Springer-Verlag, 1980.[20] K.W. Schmid, Prog. Part. Nucl. Phys. 52 (2004) 565.[21] K. Hara and Y. Sun, Int. J. Mod. Phys. E 4 (1995) 637.[22] T. Bengtsson and I. Ragnarsson, Nucl. Phys. A 436 (1985) 14.[1] R. Balian and E. Brezin, Nuovo Cimento B 64 (1969) 37.[2] K. Hara, S. Iwasaki, Nucl. Phys. A 332 (1979) 61.[3] M. Oi, T. Mizusaki, Phys. Lett. B 707 (2012) 305.[4] T. Mizusaki, M. Oi, Phys. Lett. B 715 (2012) 219.[5] T. Mizusaki, M. Oi, Fang-Qi Chen, Yang Sun, Phys. Lett. B 725(2013) 175.[6] B. Avez, M. Bender, Phys. Rev. C. 85 (2012) 034325.[7] G.F. Bertsch and L.M. Robledo, Phys. Rev. Lett. 108 (2012)042505.[8] S. Perez-Martin and L.M. Robledo, Phys. Rev. C 76 (2007)064314.[9] M. Anguiano, J.L. Egido, L.M. Robledo, Nucl. Phys. 696 (2001)467.[10] N. Onishi and S. Yoshida, Nucl. Phys. 80 (1966) 367.[11] K. Neergård and E. W¨ust, Nucl. Phys. A 402 (1983) 311.[12] Q. Haider and D. Gogny, J. Phys. G: Nucl. Part. Phys. 18 (1992)993.[13] F. D¨onau, Phys. Rev. C 58 (1998) 872.[14] M. Bender and Paul-Henri Heenen, Phys. Rev. C 78 (2008)024309.[15] K. Hara, A. Hayashi and P. Ring, Nucl. Phys. A 385 (1982) 14.[16] M. Oi and N. Tajima, Phys. Lett. B 606 (2005) 43.[17] L.M. Robledo, Phys. Rev. C 79 (2009) 021302(R).[18] L.M. Robledo, Phys. Rev. C 84 (2011) 014307.[19] P. Ring and P. Schuck, The Nuclear Many-Body Problem,Springer-Verlag, 1980.[20] K.W. Schmid, Prog. Part. Nucl. Phys. 52 (2004) 565.[21] K. Hara and Y. Sun, Int. J. Mod. Phys. E 4 (1995) 637.[22] T. Bengtsson and I. Ragnarsson, Nucl. Phys. A 436 (1985) 14.