TABLE OF CONTENTS
- ABINIT/m_CtqmcInterface
- ABINIT/m_CtqmcInterface/CtqmcInterface_finalize
- ABINIT/m_CtqmcInterface/CtqmcInterface_init
- ABINIT/m_CtqmcInterface/CtqmcInterface_run
- ABINIT/m_CtqmcInterface/CtqmcInterface_setOpts
- ABINIT/m_CtqmcInterface/CtqmcInterface_setSweeps
- m_CtqmcInterface/CtqmcInterface
ABINIT/m_CtqmcInterface [ Modules ]
NAME
m_CtqmcInterface
FUNCTION
Manage a ctqmc simulation. friendly interface for the user
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 .
NOTES
SOURCE
23 #include "defs.h" 24 MODULE m_CtqmcInterface 25 USE m_Ctqmc 26 use defs_basis 27 28 IMPLICIT NONE
ABINIT/m_CtqmcInterface/CtqmcInterface_finalize [ Functions ]
NAME
CtqmcInterface_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
this=ctqmcinterface
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
408 SUBROUTINE CtqmcInterface_finalize(this) 409 410 !Arguments ------------------------------------ 411 TYPE(CtqmcInterface), INTENT(INOUT) :: this 412 413 !IF ( this%Hybrid%init .EQV. .TRUE. ) THEN 414 ! CALL Ctqmc_printAll(this%Hybrid) 415 CALL Ctqmc_destroy(this%Hybrid) 416 !END IF 417 418 END SUBROUTINE CtqmcInterface_finalize
ABINIT/m_CtqmcInterface/CtqmcInterface_init [ Functions ]
NAME
CtqmcInterface_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
this=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
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
106 SUBROUTINE CtqmcInterface_init(this,iseed,sweeps,thermalization,measurements,flavors,samples,beta,U,ostream,MPI_COMM,nspinor) 107 108 !Arguments ------------------------------------ 109 TYPE(CtqmcInterface), INTENT(INOUT) :: this 110 INTEGER, OPTIONAL, INTENT(IN) :: MPI_COMM 111 INTEGER, INTENT(IN) :: iseed 112 DOUBLE PRECISION, INTENT(IN) :: sweeps 113 INTEGER, INTENT(IN) :: thermalization 114 INTEGER, INTENT(IN) :: measurements 115 INTEGER, INTENT(IN) :: flavors 116 INTEGER, INTENT(IN) :: samples 117 !INTEGER, INTENT(IN) :: Wmax 118 INTEGER, INTENT(IN) :: ostream 119 INTEGER, INTENT(IN) :: nspinor 120 DOUBLE PRECISION, INTENT(IN) :: beta 121 DOUBLE PRECISION, INTENT(IN) :: u 122 !DOUBLE PRECISION, INTENT(IN) :: mu 123 !Local arguements ----------------------------- 124 INTEGER :: ifstream!,opt_nondiag 125 DOUBLE PRECISION, DIMENSION(1:10) :: buffer 126 ! opt_nondiag=0 127 128 ifstream = 42 129 130 buffer(1)=DBLE(iseed) 131 buffer(2)=sweeps 132 buffer(3)=DBLE(thermalization) 133 buffer(4)=DBLE(measurements) 134 buffer(5)=DBLE(flavors) 135 buffer(6)=DBLE(samples) 136 buffer(7)=beta 137 buffer(8)=U 138 buffer(9)=GREENHYB_TAU 139 buffer(10)=nspinor 140 ! buffer(10)=DBLE(opt_nondiag) 141 !buffer(9)=0.d0!mu 142 !buffer(9)=DBLE(Wmax) 143 144 IF ( PRESENT( MPI_COMM ) ) THEN 145 CALL Ctqmc_init(this%Hybrid, ostream, ifstream, .FALSE., MY_COMM=MPI_COMM,iBuffer=buffer) 146 ELSE 147 CALL Ctqmc_init(this%Hybrid, ostream, ifstream, .FALSE.,iBuffer=buffer) 148 END IF 149 this%opt_fk = 0 150 this%opt_order = 0 151 this%opt_histo = 0 152 this%opt_movie = 0 153 this%opt_analysis = 0 154 this%opt_check = 0 155 this%opt_noise = 0 156 this%opt_spectra = 0 157 END SUBROUTINE CtqmcInterface_init
ABINIT/m_CtqmcInterface/CtqmcInterface_run [ Functions ]
NAME
CtqmcInterface_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
this=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
266 SUBROUTINE CtqmcInterface_run(this,G0omega, Gtau, Gw, D,E,Noise,matU,opt_sym,opt_levels,Magmom_orb,Magmom_spin,Magmom_tot,Iatom, & 267 &fname) 268 269 !Arguments ------------------------------------ 270 TYPE(CtqmcInterface), INTENT(INOUT) :: this 271 COMPLEX(KIND=8) , DIMENSION(:,:) , INTENT(IN ) :: G0omega 272 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT( OUT) :: Gtau 273 COMPLEX(KIND=8) , DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gw 274 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT( OUT) :: D 275 DOUBLE PRECISION , OPTIONAL, INTENT( OUT) :: E 276 DOUBLE PRECISION , OPTIONAL, INTENT( OUT) :: Noise 277 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(IN ) :: matU 278 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(IN ) :: opt_sym 279 DOUBLE PRECISION, DIMENSION(: ), OPTIONAL, INTENT(IN ) :: opt_levels 280 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(IN ) :: Magmom_orb 281 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(IN ) :: Magmom_spin 282 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(IN ) :: Magmom_tot 283 INTEGER, INTENT(IN ) :: Iatom 284 character(len=fnlen), INTENT(INOUT) :: fname 285 CALL Ctqmc_reset(this%Hybrid) 286 287 ! ifstream = 42 288 ! 289 ! OPEN(UNIT=ifstream, FILE="Gw.dat") 290 ! CALL Ctqmc_setG0w(Hybrid, ifstream) 291 ! CLOSE(ifstream) 292 ! 293 294 IF ( PRESENT(opt_levels) ) & 295 CALL Ctqmc_setMu(this%Hybrid, opt_levels) 296 297 CALL Ctqmc_setG0wTab(this%Hybrid, G0omega,this%opt_fk) 298 299 IF ( PRESENT(matU) ) & 300 CALL Ctqmc_setU(this%Hybrid, matU) 301 302 IF ( PRESENT(Magmom_orb) ) & 303 CALL Ctqmc_setMagmom(this%Hybrid, Magmom_orb, Magmom_spin, Magmom_tot) 304 305 CALL Ctqmc_run(this%Hybrid,opt_order=this%opt_order, & 306 opt_histo=this%opt_histo, & 307 opt_movie=this%opt_movie, & 308 opt_analysis=this%opt_analysis, & 309 opt_check=this%opt_check, & 310 opt_noise=this%opt_noise, & 311 opt_spectra=this%opt_spectra, & 312 opt_gMove=this%opt_gMove) 313 314 CALL Ctqmc_getResult(this%Hybrid,Iatom,fname) 315 316 IF ( PRESENT(opt_sym) ) THEN 317 CALL Ctqmc_symmetrizeGreen(this%Hybrid,opt_sym) 318 END IF 319 320 IF ( PRESENT(Gtau) .AND. PRESENT(Gw) ) THEN 321 CALL Ctqmc_getGreen(this%Hybrid, Gtau=Gtau, Gw=Gw) 322 ! write(6,*) "size",size(Gw,dim=1),size(gw,dim=2) 323 ! IF ( hybrid%rank .EQ. 0 ) write(389,*) Gw(:,hybrid%flavors+1) 324 ! call flush(389) 325 ELSE IF ( PRESENT(Gtau) .AND. .NOT. PRESENT(Gw) ) THEN 326 CALL Ctqmc_getGreen(this%Hybrid, Gtau=Gtau) 327 ELSE IF ( .NOT. PRESENT(Gtau) .AND. PRESENT(Gw) ) THEN 328 CALL Ctqmc_getGreen(this%Hybrid, Gw=Gw) 329 END IF 330 331 IF ( PRESENT(D) ) & 332 CALL Ctqmc_getD(this%Hybrid, D) 333 334 IF ( PRESENT(E) .AND. PRESENT(Noise) ) THEN 335 CALL Ctqmc_getE(this%Hybrid, E=E, Noise=Noise) 336 ELSE IF ( PRESENT(E) .AND. .NOT. PRESENT(Noise) ) THEN 337 CALL Ctqmc_getE(this%Hybrid, E=E) 338 ELSE IF ( .NOT. PRESENT(E) .AND. PRESENT(Noise) ) THEN 339 CALL Ctqmc_getE(this%Hybrid, Noise=noise) 340 END IF 341 342 343 CALL Ctqmc_printAll(this%Hybrid) 344 !CALL Ctqmc_printQMC(this%Hybrid) 345 346 END SUBROUTINE CtqmcInterface_run
ABINIT/m_CtqmcInterface/CtqmcInterface_setOpts [ Functions ]
NAME
CtqmcInterface_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
this=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
196 SUBROUTINE CtqmcInterface_setOpts(this,opt_Fk,opt_order,opt_histo,opt_movie,& 197 & opt_analysis,opt_check, opt_noise, opt_spectra, opt_gMove) 198 199 !Arguments ------------------------------------ 200 TYPE(CtqmcInterface), INTENT(INOUT) :: this 201 INTEGER , OPTIONAL , INTENT(IN ) :: opt_Fk 202 INTEGER , OPTIONAL , INTENT(IN ) :: opt_order 203 INTEGER , OPTIONAL , INTENT(IN ) :: opt_histo 204 INTEGER , OPTIONAL , INTENT(IN ) :: opt_movie 205 INTEGER , OPTIONAL , INTENT(IN ) :: opt_analysis 206 INTEGER , OPTIONAL , INTENT(IN ) :: opt_check 207 INTEGER , OPTIONAL , INTENT(IN ) :: opt_noise 208 INTEGER , OPTIONAL , INTENT(IN ) :: opt_spectra 209 INTEGER , OPTIONAL , INTENT(IN ) :: opt_gMove 210 211 IF ( PRESENT(opt_Fk) ) & 212 this%opt_Fk = opt_fk 213 IF ( PRESENT(opt_order) ) & 214 this%opt_order = opt_order 215 IF ( PRESENT(opt_histo) ) & 216 this%opt_histo = opt_histo 217 IF ( PRESENT(opt_analysis) ) & 218 this%opt_analysis = opt_analysis 219 IF ( PRESENT(opt_check) ) & 220 this%opt_check = opt_check 221 IF ( PRESENT(opt_movie) ) & 222 this%opt_movie = opt_movie 223 IF ( PRESENT(opt_noise) ) & 224 this%opt_noise = opt_noise 225 IF ( PRESENT(opt_spectra) ) & 226 this%opt_spectra = opt_spectra 227 IF ( PRESENT(opt_gMove) ) & 228 this%opt_gMove = opt_gMove 229 230 END SUBROUTINE CtqmcInterface_setOpts
ABINIT/m_CtqmcInterface/CtqmcInterface_setSweeps [ Functions ]
NAME
CtqmcInterface_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
this=ctqmcinterface sweeps=new number of sweeps
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
374 SUBROUTINE CtqmcInterface_setSweeps(this, sweeps) 375 376 !Arguments ------------------------------------ 377 TYPE(CtqmcInterface), INTENT(INOUT) :: this 378 DOUBLE PRECISION, INTENT(IN) :: sweeps 379 380 CALL Ctqmc_setSweeps(this%Hybrid,sweeps) 381 END SUBROUTINE CtqmcInterface_setSweeps
m_CtqmcInterface/CtqmcInterface [ Types ]
[ Top ] [ m_CtqmcInterface ] [ Types ]
NAME
CtqmcInterface
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
49 TYPE, PUBLIC :: CtqmcInterface 50 TYPE(Ctqmc) :: Hybrid 51 INTEGER _PRIVATE :: opt_fk = 0 52 INTEGER _PRIVATE :: opt_order = 0 53 INTEGER _PRIVATE :: opt_histo = 0 54 INTEGER _PRIVATE :: opt_movie = 0 55 INTEGER _PRIVATE :: opt_analysis = 0 56 INTEGER _PRIVATE :: opt_check = 0 57 INTEGER _PRIVATE :: opt_spectra = 0 58 INTEGER _PRIVATE :: opt_noise = 0 59 INTEGER _PRIVATE :: opt_gMove = 0 60 END TYPE CtqmcInterface