TABLE OF CONTENTS
ABINIT/fock2ACE [ Functions ]
NAME
fock2ACE
FUNCTION
Compute nonlocal contribution to the Fock part of the hamiltonian in the ACE formalism. optionally contribution to Fock forces
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (FJ,XG,MT) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
cg(2,mcg)=wavefunctions (may be read from disk file) cprj(natom,mcprj*usecprj)=<p_lmn|Cnk> coefficients for each WF |Cnk> and each NL proj |p_lmn> fock <type(fock_type)>= quantities to calculate Fock exact exchange istwfk(nkpt)=input option parameter that describes the storage of wfs kg(3,mpw*mkmem)=reduced coordinates (integers) of G vecs in basis kpt(3,nkpt)=k points in reduced coordinates mband=maximum number of bands mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol mcprj=size of projected wave-functions array (cprj) =nspinor*mband*mkmem*nsppol mgfft=maximum size of 1D FFTs mkmem=number of k points treated by this node. mpi_enreg=informations about MPI parallelization mpsang= mpw= maximum number of plane waves my_natom=number of atoms treated by current processor natom=number of atoms in cell. nband(nkpt)=number of bands at each k point nfft=number of FFT grid points ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft nkpt=number of k points in Brillouin zone nloalg(3)=governs the choice of the algorithm for non-local operator. npwarr(nkpt)=number of planewaves in basis and boundary at each k nspden=Number of spin Density components nspinor=number of spinorial components of the wavefunctions nsppol=1 for unpolarized, 2 for spin-polarized ntypat=number of types of atoms occ(mband*nkpt*nsppol)=occupation numbers for each band over all k points optfor=1 if computation of forces is required paw_ij(my_natom*usepaw) <type(paw_ij_type)>=paw arrays given on (i,j) channels pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data ph1d(2,3*(2*mgfft+1)*natom)=one-dimensional structure factor information psps <type(pseudopotential_type)>=variables related to pseudopotentials rprimd(3,3)=dimensional primitive translations in real space (bohr) typat(natom)=type of each atom usecprj=1 if cprj datastructure has been allocated use_gpu_cuda= 0 or 1 to know if we use cuda for nonlop call wtk(nkpt)=weight associated with each k point xred(3,natom)=reduced dimensionless atomic coordinates ylm(mpw*mkmem,mpsang*mpsang*useylm)= real spherical harmonics for each G and k point
OUTPUT
fock%fockACE(ikpt,isppol)%xi if optfor=1, fock%fock_common%forces
NOTES
PARENTS
scfcv
CHILDREN
bandfft_kpt_restoretabs,bandfft_kpt_savetabs,destroy_hamiltonian dotprod_g,fock_getghc,init_hamiltonian,load_k_hamiltonian load_spin_hamiltonian,mkffnl,mkkpg,pawcprj_alloc,pawcprj_free pawcprj_get,pawcprj_reorder,prep_bandfft_tabs,timab,xmpi_sum,zpotrf ztrtrs
SOURCE
77 #if defined HAVE_CONFIG_H 78 #include "config.h" 79 #endif 80 81 #include "abi_common.h" 82 83 subroutine fock2ACE(cg,cprj,fock,istwfk,kg,kpt,mband,mcg,mcprj,mgfft,mkmem,mpi_enreg,mpsang,& 84 & mpw,my_natom,natom,nband,nfft,ngfft,nkpt,nloalg,npwarr,nspden,nspinor,nsppol,& 85 & ntypat,occ,optfor,paw_ij,pawtab,ph1d,psps,rprimd,typat,usecprj,use_gpu_cuda,wtk,xred,ylm) 86 87 use defs_basis 88 use defs_datatypes 89 use defs_abitypes 90 use m_profiling_abi 91 use m_xmpi 92 use m_errors 93 use m_fock 94 use m_hamiltonian, only : init_hamiltonian,destroy_hamiltonian,load_spin_hamiltonian,& 95 & load_k_hamiltonian,gs_hamiltonian_type 96 use m_bandfft_kpt, only : bandfft_kpt,bandfft_kpt_type,& 97 & bandfft_kpt_savetabs,bandfft_kpt_restoretabs 98 use m_pawtab, only : pawtab_type 99 use m_paw_ij, only : paw_ij_type 100 use m_pawcprj, only : pawcprj_type,pawcprj_alloc,pawcprj_free,pawcprj_get,pawcprj_reorder 101 use m_cgtools 102 103 !This section has been created automatically by the script Abilint (TD). 104 !Do not modify the following lines by hand. 105 #undef ABI_FUNC 106 #define ABI_FUNC 'fock2ACE' 107 use interfaces_18_timing 108 use interfaces_32_util 109 use interfaces_66_nonlocal 110 use interfaces_66_wfs, except_this_one => fock2ACE 111 !End of the abilint section 112 113 implicit none 114 115 !Arguments ------------------------------------ 116 !scalars 117 integer,intent(in) :: mband,mcg,mcprj,mgfft,mkmem,mpsang,mpw,my_natom,natom,nfft,nkpt 118 integer,intent(in) :: nspden,nsppol,nspinor,ntypat,optfor 119 integer,intent(in) :: usecprj,use_gpu_cuda 120 type(MPI_type),intent(inout) :: mpi_enreg 121 type(pseudopotential_type),intent(in) :: psps 122 !arrays 123 integer,intent(in) :: istwfk(nkpt),kg(3,mpw*mkmem),nband(nkpt*nsppol) 124 integer,intent(in) :: ngfft(18),nloalg(3),npwarr(nkpt) 125 integer,intent(in) :: typat(natom) 126 real(dp),intent(in) :: cg(2,mcg) 127 real(dp),intent(in) :: kpt(3,nkpt) 128 real(dp),intent(in) :: occ(mband*nkpt*nsppol),ph1d(2,3*(2*mgfft+1)*natom) 129 real(dp),intent(in) :: rprimd(3,3),wtk(nkpt),xred(3,natom) 130 real(dp),intent(in) :: ylm(mpw*mkmem,mpsang*mpsang*psps%useylm) 131 type(pawcprj_type),intent(inout) :: cprj(natom,mcprj*usecprj) 132 type(paw_ij_type),intent(in) :: paw_ij(my_natom*psps%usepaw) 133 type(pawtab_type),intent(in) :: pawtab(ntypat*psps%usepaw) 134 type(fock_type),pointer, intent(inout) :: fock 135 !Local variables------------------------------- 136 !scalars 137 integer :: bandpp,bdtot_index,dimffnl,iband,iband_cprj,iband_last,ibg,icg,ider 138 integer :: idir,ierr,ikg,ikpt,ilm,ipw,isppol,istwf_k,kk,ll 139 integer :: mband_cprj,me_distrb,my_ikpt,my_nspinor,nband_k,nband_cprj_k,ndat,nkpg 140 integer :: npw_k,spaceComm 141 integer :: use_ACE_old 142 integer :: blocksize,iblock,jblock,iblocksize,jblocksize,nblockbd 143 !integer, save :: counter=0 144 type(gs_hamiltonian_type) :: gs_hamk 145 logical :: compute_gbound 146 character(len=500) :: msg 147 type(fock_common_type),pointer :: fockcommon 148 !arrays 149 integer,allocatable :: kg_k(:,:) 150 real(dp) :: kpoint(3),rmet(3,3),tsec(2) 151 real(dp),allocatable :: bb(:,:,:),cwavef(:,:),cwavefk(:,:),ffnl_sav(:,:,:,:) 152 real(dp),allocatable :: kpg_k(:,:),kpg_k_sav(:,:) 153 real(dp),allocatable :: mkl(:,:,:),occblock(:),ph3d(:,:,:),ph3d_sav(:,:,:) 154 real(dp),allocatable :: wi(:,:,:),weight(:),ylm_k(:,:),ylmgr_k(:,:,:) 155 real(dp),allocatable,target :: ffnl(:,:,:,:) 156 type(bandfft_kpt_type),pointer :: my_bandfft_kpt => null() 157 type(pawcprj_type),target,allocatable :: cwaveprj(:,:) 158 type(pawcprj_type),pointer :: cwaveprj_idat(:,:) 159 160 !************************************************************************* 161 162 call timab(920,1,tsec) 163 call timab(921,1,tsec) 164 165 !DEBUG 166 !if(counter>0)return 167 !counter=counter+1 168 !ENDDEBUG 169 170 !Init mpicomm and me 171 if(mpi_enreg%paral_kgb==1)then 172 spaceComm=mpi_enreg%comm_kpt 173 me_distrb=mpi_enreg%me_kpt 174 else 175 !* In case of HF calculation 176 if (mpi_enreg%paral_hf==1) then 177 spaceComm=mpi_enreg%comm_kpt 178 me_distrb=mpi_enreg%me_kpt 179 else 180 spaceComm=mpi_enreg%comm_cell 181 me_distrb=mpi_enreg%me_cell 182 end if 183 end if 184 185 !Some initializations 186 my_nspinor=max(1,nspinor/mpi_enreg%nproc_spinor) 187 compute_gbound=.true. 188 fockcommon => fock%fock_common 189 use_ACE_old=fockcommon%use_ACE 190 fockcommon%use_ACE=0 191 192 !Initialize Hamiltonian (k- and spin-independent terms) 193 194 call init_hamiltonian(gs_hamk,psps,pawtab,nspinor,nsppol,nspden,natom,& 195 & typat,xred,nfft,mgfft,ngfft,rprimd,nloalg,usecprj=usecprj,& 196 & comm_atom=mpi_enreg%comm_atom,mpi_atmtab=mpi_enreg%my_atmtab,mpi_spintab=mpi_enreg%my_isppoltab,& 197 & paw_ij=paw_ij,ph1d=ph1d,fock=fock,& 198 & use_gpu_cuda=use_gpu_cuda) 199 rmet = MATMUL(TRANSPOSE(rprimd),rprimd) 200 fockcommon%use_ACE=use_ACE_old 201 call timab(921,2,tsec) 202 203 !need to reorder cprj=<p_lmn|Cnk> (from unsorted to atom-sorted) 204 if (psps%usepaw==1) then 205 call pawcprj_reorder(cprj,gs_hamk%atindx) 206 end if 207 208 !LOOP OVER SPINS 209 bdtot_index=0;ibg=0;icg=0 210 211 do isppol=1,nsppol 212 fockcommon%isppol=isppol 213 ! Continue to initialize the Hamiltonian (PAW DIJ coefficients) 214 call load_spin_hamiltonian(gs_hamk,isppol,with_nonlocal=.true.) 215 216 ! Loop over k points 217 ikg=0 218 do ikpt=1,nkpt 219 fockcommon%ikpt=ikpt 220 nband_k=nband(ikpt+(isppol-1)*nkpt) 221 npw_k=npwarr(ikpt) 222 kpoint(:)=kpt(:,ikpt) 223 istwf_k=istwfk(ikpt) 224 if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,isppol,me_distrb)) then 225 bdtot_index=bdtot_index+nband_k 226 cycle 227 end if 228 229 call timab(922,1,tsec) 230 231 ! Parallelism over FFT and/or bands: define sizes and tabs 232 if (mpi_enreg%paral_kgb==1) then 233 my_ikpt=mpi_enreg%my_kpttab(ikpt) 234 nblockbd=nband_k/(mpi_enreg%nproc_band*mpi_enreg%bandpp) 235 bandpp=mpi_enreg%bandpp 236 my_bandfft_kpt => bandfft_kpt(my_ikpt) 237 else 238 my_ikpt=ikpt 239 bandpp=mpi_enreg%bandpp 240 nblockbd=nband_k/bandpp 241 end if 242 blocksize=nband_k/nblockbd 243 mband_cprj=mband/mpi_enreg%nproc_band 244 nband_cprj_k=nband_k/mpi_enreg%nproc_band 245 246 ABI_ALLOCATE(cwavef,(2,npw_k*my_nspinor*blocksize)) 247 if (psps%usepaw==1) then 248 ABI_DATATYPE_ALLOCATE(cwaveprj,(natom,my_nspinor*bandpp)) 249 call pawcprj_alloc(cwaveprj,0,gs_hamk%dimcprj) 250 else 251 ABI_DATATYPE_ALLOCATE(cwaveprj,(0,0)) 252 end if 253 254 ABI_ALLOCATE(kg_k,(3,mpw)) 255 !$OMP PARALLEL DO 256 do ipw=1,npw_k 257 kg_k(:,ipw)=kg(:,ipw+ikg) 258 end do 259 260 ABI_ALLOCATE(ylm_k,(npw_k,mpsang*mpsang*psps%useylm)) 261 ABI_ALLOCATE(ylmgr_k,(0,0,0)) 262 if (psps%useylm==1) then 263 !$OMP PARALLEL DO COLLAPSE(2) 264 do ilm=1,mpsang*mpsang 265 do ipw=1,npw_k 266 ylm_k(ipw,ilm)=ylm(ipw+ikg,ilm) 267 end do 268 end do 269 end if 270 271 ! Compute (k+G) vectors 272 nkpg=3*nloalg(3) 273 ABI_ALLOCATE(kpg_k,(npw_k,nkpg)) 274 if (nkpg>0) then 275 call mkkpg(kg_k,kpg_k,kpoint,nkpg,npw_k) 276 end if 277 278 279 ! Compute nonlocal form factors ffnl at all (k+G) 280 ider=0;idir=0;dimffnl=1 281 ABI_ALLOCATE(ffnl,(npw_k,dimffnl,psps%lmnmax,ntypat)) 282 call mkffnl(psps%dimekb,dimffnl,psps%ekb,ffnl,psps%ffspl,gs_hamk%gmet,gs_hamk%gprimd,& 283 & ider,idir,psps%indlmn,kg_k,kpg_k,kpoint,psps%lmnmax,psps%lnmax,psps%mpsang,psps%mqgrid_ff,& 284 & nkpg,npw_k,ntypat,psps%pspso,psps%qgrid_ff,rmet,psps%usepaw,psps%useylm,ylm_k,ylmgr_k) 285 286 ! Load k-dependent part in the Hamiltonian datastructure 287 ! - Compute 3D phase factors 288 ! - Prepare various tabs in case of band-FFT parallelism 289 ! - Load k-dependent quantities in the Hamiltonian 290 291 ABI_ALLOCATE(ph3d,(2,npw_k,gs_hamk%matblk)) 292 call load_k_hamiltonian(gs_hamk,kpt_k=kpoint,istwf_k=istwf_k,npw_k=npw_k,& 293 & kg_k=kg_k,kpg_k=kpg_k,ffnl_k=ffnl,ph3d_k=ph3d,compute_gbound=compute_gbound,compute_ph3d=.true.) 294 295 ! Load band-FFT tabs (transposed k-dependent arrays) 296 if (mpi_enreg%paral_kgb==1) then 297 call bandfft_kpt_savetabs(my_bandfft_kpt,ffnl=ffnl_sav,ph3d=ph3d_sav,kpg=kpg_k_sav) 298 call prep_bandfft_tabs(gs_hamk,ikpt,mkmem,mpi_enreg) 299 call load_k_hamiltonian(gs_hamk,npw_fft_k=my_bandfft_kpt%ndatarecv, & 300 & kg_k =my_bandfft_kpt%kg_k_gather, & 301 & kpg_k =my_bandfft_kpt%kpg_k_gather, & 302 ffnl_k =my_bandfft_kpt%ffnl_gather, & 303 ph3d_k =my_bandfft_kpt%ph3d_gather,compute_gbound=compute_gbound) 304 end if 305 306 call timab(922,2,tsec) 307 308 ! The following is now wrong. In sequential, nblockbd=nband_k/bandpp 309 ! blocksize= bandpp (JB 2016/04/16) 310 ! Note that in sequential mode iblock=iband, nblockbd=nband_k and blocksize=1 311 ! 312 ABI_ALLOCATE(occblock,(blocksize)) 313 ABI_ALLOCATE(weight,(blocksize)) 314 occblock=zero;weight=zero 315 316 if (fockcommon%optfor) then 317 fockcommon%forces_ikpt=zero 318 end if 319 320 ABI_ALLOCATE(wi,(2,npw_k*my_nspinor*blocksize,nblockbd)) 321 wi=zero 322 ABI_ALLOCATE(mkl,(2,nband_k,nband_k)) 323 mkl=zero 324 ! Calculate all the Wi for the current k-point 325 326 do iblock=1,nblockbd 327 328 iband=(iblock-1)*blocksize+1;iband_last=min(iband+blocksize-1,nband_k) 329 iband_cprj=(iblock-1)*bandpp+1 330 if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,iband,iband_last,isppol,me_distrb)) cycle 331 332 ! Select occupied bandsddk 333 occblock(:)=occ(1+(iblock-1)*blocksize+bdtot_index:iblock*blocksize+bdtot_index) 334 call timab(923,1,tsec) 335 weight(:)=wtk(ikpt)*occblock(:) 336 337 ! Load contribution from n,k 338 cwavef(:,1:npw_k*my_nspinor*blocksize)=& 339 & cg(:,1+(iblock-1)*npw_k*my_nspinor*blocksize+icg:iblock*npw_k*my_nspinor*blocksize+icg) 340 if (psps%usepaw==1) then 341 call pawcprj_get(gs_hamk%atindx1,cwaveprj,cprj,natom,iband_cprj,ibg,ikpt,0,isppol,& 342 & mband_cprj,mkmem,natom,bandpp,nband_cprj_k,my_nspinor,nsppol,0,& 343 & mpicomm=mpi_enreg%comm_kpt,proc_distrb=mpi_enreg%proc_distrb) 344 end if 345 346 call timab(926,2,tsec) 347 348 if (mpi_enreg%paral_kgb==1) then 349 msg='fock2ACE: Paral_kgb is not yet implemented for fock calculations' 350 MSG_BUG(msg) 351 end if 352 ndat=mpi_enreg%bandpp 353 if (gs_hamk%usepaw==0) cwaveprj_idat => cwaveprj 354 355 do iblocksize=1,blocksize 356 fockcommon%ieigen=(iblock-1)*blocksize+iblocksize 357 fockcommon%iband=(iblock-1)*blocksize+iblocksize 358 if (gs_hamk%usepaw==1) then 359 cwaveprj_idat => cwaveprj(:,(iblocksize-1)*my_nspinor+1:iblocksize*my_nspinor) 360 end if 361 call fock_getghc(cwavef(:,1+(iblocksize-1)*npw_k*my_nspinor:iblocksize*npw_k*my_nspinor),cwaveprj_idat,& 362 & wi(:,1+(iblocksize-1)*npw_k*my_nspinor:iblocksize*npw_k*my_nspinor,iblock),gs_hamk,mpi_enreg) 363 mkl(1,fockcommon%ieigen,fockcommon%ieigen)=fockcommon%eigen_ikpt(fockcommon%ieigen) 364 if (fockcommon%optfor) then 365 fockcommon%forces(:,:)=fockcommon%forces(:,:)+weight(iblocksize)*fockcommon%forces_ikpt(:,:,fockcommon%ieigen) 366 end if 367 end do 368 369 370 end do ! End of loop on block of bands 371 372 ! Calculate Mkl for the current k-point 373 ABI_ALLOCATE(cwavefk,(2,npw_k*my_nspinor)) 374 do iblock=1,nblockbd 375 cwavef(:,1:npw_k*my_nspinor*blocksize)=& 376 & cg(:,1+(iblock-1)*npw_k*my_nspinor*blocksize+icg:iblock*npw_k*my_nspinor*blocksize+icg) 377 do iblocksize=1,blocksize 378 kk=(iblock-1)*blocksize+iblocksize 379 cwavefk(:,:)=cwavef(:,1+(iblocksize-1)*npw_k*my_nspinor:iblocksize*npw_k*my_nspinor) 380 do jblock=1,iblock 381 do jblocksize=1,blocksize 382 ll=(jblock-1)*blocksize+jblocksize 383 if (ll<kk) then 384 call dotprod_g(mkl(1,kk,ll),mkl(2,kk,ll),gs_hamk%istwf_k,npw_k,2,wi(:,1+(jblocksize-1)*npw_k*my_nspinor:& 385 & jblocksize*npw_k*my_nspinor,jblock),cwavefk,mpi_enreg%me_g0,mpi_enreg%comm_fft) 386 end if 387 end do 388 end do 389 end do 390 end do ! End of loop on block of bands 391 392 ABI_DEALLOCATE(cwavefk) 393 mkl=-mkl 394 395 ! Cholesky factorisation of -mkl=Lx(trans(L)*. On output mkl=L 396 call zpotrf("L",nband_k,mkl,nband_k,ierr) 397 398 ! calculate trans(L-1) 399 ABI_ALLOCATE(bb,(2,nband_k,nband_k)) 400 bb=zero 401 do kk=1,nband_k 402 bb(1,kk,kk)=one 403 end do 404 call ztrtrs("L","T","N",nband_k,nband_k,mkl,nband_k,bb,nband_k,ierr) 405 fock%fockACE(ikpt,isppol)%xi=zero 406 407 ! Calculate ksi 408 do kk=1,nband_k 409 do jblock=1,nblockbd 410 do jblocksize=1,blocksize 411 ll=(jblock-1)*blocksize+jblocksize 412 fock%fockACE(ikpt,isppol)%xi(1,:,kk)=fock%fockACE(ikpt,isppol)%xi(1,:,kk)+bb(1,ll,kk)*wi(1,1+(jblocksize-1)*& 413 & npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock)-& 414 & bb(2,ll,kk)*wi(2,1+(jblocksize-1)*npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock) 415 fock%fockACE(ikpt,isppol)%xi(2,:,kk)=fock%fockACE(ikpt,isppol)%xi(2,:,kk)+bb(1,ll,kk)*wi(2,1+(jblocksize-1)*& 416 npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock)+& 417 & bb(2,ll,kk)*wi(1,1+(jblocksize-1)*npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock) 418 end do 419 end do 420 end do 421 422 ! DEBUG 423 ! fock%fockACE(ikpt,isppol)%xi=zero 424 ! ENDDEBUG 425 426 ABI_DEALLOCATE(wi) 427 ABI_DEALLOCATE(mkl) 428 429 ! Restore the bandfft tabs 430 if (mpi_enreg%paral_kgb==1) then 431 call bandfft_kpt_restoretabs(my_bandfft_kpt,ffnl=ffnl_sav,ph3d=ph3d_sav,kpg=kpg_k_sav) 432 end if 433 434 ! Increment indices 435 bdtot_index=bdtot_index+nband_k 436 if (mkmem/=0) then 437 ibg=ibg+my_nspinor*nband_cprj_k 438 icg=icg+npw_k*my_nspinor*nband_k 439 ikg=ikg+npw_k 440 end if 441 442 if (psps%usepaw==1) then 443 call pawcprj_free(cwaveprj) 444 end if 445 ABI_DATATYPE_DEALLOCATE(cwaveprj) 446 ABI_DEALLOCATE(cwavef) 447 ABI_DEALLOCATE(bb) 448 ABI_DEALLOCATE(occblock) 449 ABI_DEALLOCATE(weight) 450 ABI_DEALLOCATE(ffnl) 451 ABI_DEALLOCATE(kg_k) 452 ABI_DEALLOCATE(kpg_k) 453 ABI_DEALLOCATE(ylm_k) 454 ABI_DEALLOCATE(ylmgr_k) 455 ABI_DEALLOCATE(ph3d) 456 end do ! End k point loop 457 end do ! End loop over spins 458 459 !Parallel case: accumulate (n,k) contributions 460 if (xmpi_paral==1) then 461 ! Forces 462 if (optfor==1) then 463 call timab(65,2,tsec) 464 if (psps%usepaw==1) then 465 call xmpi_sum(fockcommon%forces,spaceComm,ierr) 466 end if 467 end if 468 end if 469 470 call timab(925,1,tsec) 471 472 !need to reorder cprj=<p_lmn|Cnk> (from atom-sorted to unsorted) 473 if (psps%usepaw==1) then 474 call pawcprj_reorder(cprj,gs_hamk%atindx1) 475 end if 476 !Deallocate temporary space 477 call destroy_hamiltonian(gs_hamk) 478 479 call timab(925,2,tsec) 480 call timab(920,2,tsec) 481 482 end subroutine fock2ACE