TABLE OF CONTENTS
- ABINIT/m_CtqmcoffdiagInterface
- ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_finalize
- ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_init
- ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_run
- ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_setOpts
- ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_setSweeps
- m_CtqmcoffdiagInterface/CtqmcoffdiagInterface
ABINIT/m_CtqmcoffdiagInterface [ Modules ]
NAME
m_CtqmcoffdiagInterface
FUNCTION
Manage a ctqmc simulation. friendly interface for the user
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder, B. Amadon, J. Denier) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
SOURCE
23 #include "defs.h" 24 MODULE m_CtqmcoffdiagInterface 25 USE m_Ctqmcoffdiag 26 use defs_basis 27 28 IMPLICIT NONE
ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_finalize [ Functions ]
NAME
CtqmcoffdiagInterface_finalize
FUNCTION
Destroy simulation
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
op=ctqmcinterface
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
428 SUBROUTINE CtqmcoffdiagInterface_finalize(op) 429 430 !Arguments ------------------------------------ 431 TYPE(CtqmcoffdiagInterface), INTENT(INOUT) :: op 432 433 !IF ( op%Hybrid%init .EQV. .TRUE. ) THEN 434 ! CALL Ctqmcoffdiag_printAll(op%Hybrid) 435 !write(6,*) "before ctqmc_destroy in Ctqmcoffdiaginterface_finalize" 436 CALL Ctqmcoffdiag_destroy(op%Hybrid) 437 !END IF 438 439 END SUBROUTINE CtqmcoffdiagInterface_finalize
ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_init [ Functions ]
NAME
CtqmcoffdiagInterface_init
FUNCTION
Initialize with permanent parameters
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
op=ctqmcinterface iseed=seed for rng sweeps=number of sweeps for the run thermalization=number of sweeps to thermalize measurements=how often we measure (modulo) flavors=number of orbitals (spin degenerated) samples=imaginary time slices beta=inverse temperature U=interaction parameter ostream=where to write output MPI_COMM=mpi communicator for the run nspinor=number of spinor
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
99 SUBROUTINE CtqmcoffdiagInterface_init(op,iseed,sweeps,thermalization,& 100 &measurements,flavors,samples,beta,U,ostream,MPI_COMM,opt_nondiag,nspinor) 101 102 !Arguments ------------------------------------ 103 TYPE(CtqmcoffdiagInterface), INTENT(INOUT) :: op 104 INTEGER, OPTIONAL, INTENT(IN) :: MPI_COMM 105 INTEGER, INTENT(IN) :: iseed 106 DOUBLE PRECISION, INTENT(IN) :: sweeps 107 INTEGER, INTENT(IN) :: thermalization 108 INTEGER, INTENT(IN) :: measurements 109 INTEGER, INTENT(IN) :: flavors 110 INTEGER, INTENT(IN) :: samples 111 !INTEGER, INTENT(IN) :: Wmax 112 INTEGER, INTENT(IN) :: ostream 113 INTEGER, INTENT(IN) :: opt_nondiag 114 DOUBLE PRECISION, INTENT(IN) :: beta 115 DOUBLE PRECISION, INTENT(IN) :: u 116 INTEGER, INTENT(IN) :: nspinor 117 !DOUBLE PRECISION, INTENT(IN) :: mu 118 !Local arguements ----------------------------- 119 INTEGER :: ifstream 120 DOUBLE PRECISION, DIMENSION(1:11) :: buffer 121 122 ifstream = 42 123 124 buffer(1)=DBLE(iseed) 125 buffer(2)=sweeps 126 buffer(3)=DBLE(thermalization) 127 buffer(4)=DBLE(measurements) 128 buffer(5)=DBLE(flavors) 129 buffer(6)=DBLE(samples) 130 buffer(7)=beta 131 buffer(8)=U 132 buffer(9)=GREENHYB_TAU 133 buffer(10)=DBLE(opt_nondiag) 134 buffer(11)=nspinor 135 !buffer(9)=0.d0!mu 136 !buffer(9)=DBLE(Wmax) 137 138 IF ( PRESENT( MPI_COMM ) ) THEN 139 CALL Ctqmcoffdiag_init(op%Hybrid, ostream, ifstream, .FALSE., MPI_COMM,buffer) 140 ELSE 141 CALL Ctqmcoffdiag_init(op%Hybrid, ostream, ifstream, .FALSE.,iBuffer=buffer) 142 END IF 143 op%opt_fk = 0 144 op%opt_order = 0 145 op%opt_histo = 0 146 op%opt_movie = 0 147 op%opt_analysis = 0 148 op%opt_check = 0 149 op%opt_noise = 0 150 op%opt_spectra = 0 151 END SUBROUTINE CtqmcoffdiagInterface_init
ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_run [ Functions ]
NAME
CtqmcoffdiagInterface_run
FUNCTION
run a ctqmc simu and get results
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
op=ctqmcinterface G0omega=Gw0 (according to opt_Fk) matU=interaction matrice opt_sym=weight factors to symmetrise G opt_levels=energy for each level (with respect to fermi level)
OUTPUT
Gtau=G(tau) Gw=fourier transform of Gtau D=full double occupancy E=Interaction energy Noise=Noise on E
SIDE EFFECTS
NOTES
SOURCE
260 SUBROUTINE CtqmcoffdiagInterface_run(op,G0omega, Gtau, Gw, D,E,Noise,matU,Docc,opt_sym,opt_levels,hybri_limit,Magmom_orb,& 261 &Magmom_spin,Magmom_tot,Iatom,fname) 262 263 !Arguments ------------------------------------ 264 TYPE(CtqmcoffdiagInterface), INTENT(INOUT) :: op 265 COMPLEX(KIND=8) , DIMENSION(:,:,:), INTENT(IN ) :: G0omega 266 DOUBLE PRECISION, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: Gtau 267 COMPLEX(KIND=8) , DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: Gw 268 DOUBLE PRECISION, OPTIONAL , INTENT(OUT) :: D 269 DOUBLE PRECISION, OPTIONAL , INTENT(OUT) :: E 270 DOUBLE PRECISION, OPTIONAL , INTENT(OUT) :: Noise 271 DOUBLE PRECISION, DIMENSION(:,:),OPTIONAL, INTENT(IN ) :: matU 272 DOUBLE PRECISION, DIMENSION(:,:),OPTIONAL, INTENT(OUT ) :: Docc 273 DOUBLE PRECISION, DIMENSION(:,:),OPTIONAL, INTENT(IN ) :: opt_sym 274 DOUBLE PRECISION, DIMENSION(:),OPTIONAL, INTENT(IN ) :: opt_levels 275 COMPLEX(KIND=8) , DIMENSION(:,:),OPTIONAL, INTENT(IN ) :: hybri_limit 276 DOUBLE PRECISION, DIMENSION(:,:),OPTIONAL, INTENT(IN ) :: Magmom_orb 277 DOUBLE PRECISION, DIMENSION(:,:),OPTIONAL, INTENT(IN ) :: Magmom_spin 278 DOUBLE PRECISION, DIMENSION(:,:),OPTIONAL, INTENT(IN ) :: Magmom_tot 279 INTEGER, INTENT(IN ) :: Iatom 280 character(len=fnlen), INTENT(INOUT) :: fname 281 !local variables-------------------------------- 282 ! INTEGER :: iflavor1,iflavor2 283 284 ! do iflavor1=1,10 285 ! do iflavor2=1,10 286 ! if(iflavor1==iflavor2) THEN 287 ! write(6,*) iflavor1, iflavor2, Magmom(iflavor1,iflavor2) 288 ! end if 289 ! end do 290 ! end do 291 292 CALL Ctqmcoffdiag_reset(op%Hybrid) 293 294 ! ifstream = 42 295 ! 296 ! OPEN(UNIT=ifstream, FILE="Gw.dat") 297 ! CALL Ctqmcoffdiag_setG0w(Hybrid, ifstream) 298 ! CLOSE(ifstream) 299 ! 300 301 IF ( PRESENT(opt_levels)) & 302 CALL Ctqmcoffdiag_setMu(op%Hybrid, opt_levels) 303 304 IF ( PRESENT(hybri_limit)) & 305 CALL Ctqmcoffdiag_sethybri_limit(op%Hybrid, hybri_limit) 306 307 !call xmpi_barrier(op%Hybrid%MY_COMM) 308 CALL Ctqmcoffdiag_setG0wTab(op%Hybrid, G0omega,op%opt_fk) 309 !call xmpi_barrier(op%Hybrid%MY_COMM) 310 311 IF ( PRESENT(matU) ) & 312 CALL Ctqmcoffdiag_setU(op%Hybrid, matU) 313 !call xmpi_barrier(op%Hybrid%MY_COMM) 314 315 IF ( PRESENT(Magmom_orb) ) & 316 CALL Ctqmcoffdiag_setMagmom(op%Hybrid, Magmom_orb, Magmom_spin, Magmom_tot) 317 !call xmpi_barrier(op%Hybrid%MY_COMM) 318 319 320 CALL Ctqmcoffdiag_run(op%Hybrid,opt_order=op%opt_order, & 321 opt_histo=op%opt_histo, & 322 opt_movie=op%opt_movie, & 323 opt_analysis=op%opt_analysis, & 324 opt_check=op%opt_check, & 325 opt_noise=op%opt_noise, & 326 opt_spectra=op%opt_spectra, & 327 opt_gMove=op%opt_gMove) 328 !call xmpi_barrier(op%Hybrid%MY_COMM) 329 330 ! write(6,*) "op%Hybrid%stats",op%Hybrid%stats 331 ! write(6,*) "opt_gMove",op%opt_gMove 332 333 CALL Ctqmcoffdiag_getResult(op%Hybrid,Iatom,fname) 334 335 IF ( PRESENT(opt_sym) ) THEN 336 CALL Ctqmcoffdiag_symmetrizeGreen(op%Hybrid,opt_sym) 337 END IF 338 339 ! write(6,*) "op%Hybrid%stats",op%Hybrid%stats 340 341 IF ( PRESENT(Gtau) .AND. PRESENT(Gw) ) THEN 342 CALL Ctqmcoffdiag_getGreen(op%Hybrid, Gtau=Gtau, Gw=Gw) 343 ! write(6,*) "size",size(Gw,dim=1),size(gw,dim=2) 344 ! IF ( hybrid%rank .EQ. 0 ) write(389,*) Gw(:,hybrid%flavors+1) 345 ! call flush(389) 346 ELSE IF ( PRESENT(Gtau) .AND. .NOT. PRESENT(Gw) ) THEN 347 CALL Ctqmcoffdiag_getGreen(op%Hybrid, Gtau=Gtau) 348 ELSE IF ( .NOT. PRESENT(Gtau) .AND. PRESENT(Gw) ) THEN 349 CALL Ctqmcoffdiag_getGreen(op%Hybrid, Gw=Gw) 350 END IF 351 352 !write(6,*) "op%Hybrid%stats",op%Hybrid%stats 353 354 IF ( PRESENT(D) ) & 355 CALL Ctqmcoffdiag_getD(op%Hybrid, D) 356 Docc=op%Hybrid%measDE 357 358 !write(6,*) "op%Hybrid%stats",op%Hybrid%stats 359 360 IF ( PRESENT(E) .AND. PRESENT(Noise) ) & 361 CALL Ctqmcoffdiag_getE(op%Hybrid, E, Noise) 362 363 CALL Ctqmcoffdiag_printAll(op%Hybrid) 364 !CALL Ctqmcoffdiag_printQMC(op%Hybrid) 365 366 END SUBROUTINE CtqmcoffdiagInterface_run
ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_setOpts [ Functions ]
NAME
CtqmcoffdiagInterface_setOpts
FUNCTION
Set and save options for many runs
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
op=ctqmcinterface opt_Fk=0 if we give us Gw0 and 1 if 1/Gw0+iwn opt_order=maximal perturbation order to scope(>0) opt_movie=print a latex file (0 or 1) opt_analysis=measure correlations (0 or 1) opt_check=check fast calculation :0 nothing 1 Impurity 2 Bath 3 Both opt_noise=calculate noise ofr green functions(0 ro 1) opt_spectra=fourier transform of time evolution of number of electrons (0 or 1) opt_gMove=number of global moves (>0)
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
190 SUBROUTINE CtqmcoffdiagInterface_setOpts(op,opt_Fk,opt_order,opt_histo,opt_movie,& 191 & opt_analysis,opt_check, opt_noise, opt_spectra, opt_gMove) 192 193 !Arguments ------------------------------------ 194 TYPE(CtqmcoffdiagInterface), INTENT(INOUT) :: op 195 INTEGER , OPTIONAL , INTENT(IN ) :: opt_Fk 196 INTEGER , OPTIONAL , INTENT(IN ) :: opt_order 197 INTEGER , OPTIONAL , INTENT(IN ) :: opt_histo 198 INTEGER , OPTIONAL , INTENT(IN ) :: opt_movie 199 INTEGER , OPTIONAL , INTENT(IN ) :: opt_analysis 200 INTEGER , OPTIONAL , INTENT(IN ) :: opt_check 201 INTEGER , OPTIONAL , INTENT(IN ) :: opt_noise 202 INTEGER , OPTIONAL , INTENT(IN ) :: opt_spectra 203 INTEGER , OPTIONAL , INTENT(IN ) :: opt_gMove 204 205 IF ( PRESENT(opt_Fk) ) & 206 op%opt_Fk = opt_fk 207 IF ( PRESENT(opt_order) ) & 208 op%opt_order = opt_order 209 IF ( PRESENT(opt_histo) ) & 210 op%opt_histo = opt_histo 211 IF ( PRESENT(opt_analysis) ) & 212 op%opt_analysis = opt_analysis 213 IF ( PRESENT(opt_check) ) & 214 op%opt_check = opt_check 215 IF ( PRESENT(opt_movie) ) & 216 op%opt_movie = opt_movie 217 IF ( PRESENT(opt_noise) ) & 218 op%opt_noise = opt_noise 219 IF ( PRESENT(opt_spectra) ) & 220 op%opt_spectra = opt_spectra 221 IF ( PRESENT(opt_gMove) ) & 222 op%opt_gMove = opt_gMove 223 224 END SUBROUTINE CtqmcoffdiagInterface_setOpts
ABINIT/m_CtqmcoffdiagInterface/CtqmcoffdiagInterface_setSweeps [ Functions ]
NAME
CtqmcoffdiagInterface_setSweeps
FUNCTION
change sweeps on the fly
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
op=ctqmcinterface sweeps=new number of sweeps
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
394 SUBROUTINE CtqmcoffdiagInterface_setSweeps(op, sweeps) 395 396 !Arguments ------------------------------------ 397 TYPE(CtqmcoffdiagInterface), INTENT(INOUT) :: op 398 DOUBLE PRECISION, INTENT(IN) :: sweeps 399 400 CALL Ctqmcoffdiag_setSweeps(op%Hybrid,sweeps) 401 END SUBROUTINE CtqmcoffdiagInterface_setSweeps
m_CtqmcoffdiagInterface/CtqmcoffdiagInterface [ Types ]
[ Top ] [ m_CtqmcoffdiagInterface ] [ Types ]
NAME
CtqmcoffdiagInterface
FUNCTION
This structured datatype contains the necessary data
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
47 TYPE CtqmcoffdiagInterface 48 TYPE(Ctqmcoffdiag) :: Hybrid 49 INTEGER :: opt_fk = 0 50 INTEGER :: opt_order = 0 51 INTEGER :: opt_histo = 0 52 INTEGER :: opt_movie = 0 53 INTEGER :: opt_analysis = 0 54 INTEGER :: opt_check = 0 55 INTEGER :: opt_spectra = 0 56 INTEGER :: opt_noise = 0 57 INTEGER :: opt_gMove = 0 58 END TYPE CtqmcoffdiagInterface