Branch data Line data Source code
1 : :
2 : : #if defined HAVE_CONFIG_H
3 : : #include "config.h"
4 : : #endif
5 : : !!****m* ABINIT/m_CtqmcoffdiagComplex
6 : : !! NAME
7 : : !! m_CtqmcoffdiagComplex
8 : : !!
9 : : !! FUNCTION
10 : : !! Manage and drive all the CTQMC
11 : : !! Should not be used if you don't know what you do
12 : : !! Please use CtqmcoffdiagComplexComplexInterface
13 : : !!
14 : : !! COPYRIGHT
15 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder, B. Amadon, J. Denier)
16 : : !! This file is distributed under the terms of the
17 : : !! GNU General Public License, see ~abinit/COPYING
18 : : !! or http://www.gnu.org/copyleft/gpl.txt .
19 : : !!
20 : : !! NOTES
21 : : !!
22 : : !! SOURCE
23 : :
24 : : #include "defs.h"
25 : : #define CTQMC_SLICE1 100
26 : : ! Coupe Sweeps en 100
27 : : #define CTQMC_SLICE2 100
28 : : ! Coupe modNoise1 en 2000
29 : : #define CTQMC_SEGME 1
30 : : #define CTQMC_ANTIS -2
31 : : #define CTQMC_ADDED 3
32 : : #define CTQMC_REMOV 4
33 : : #define CTQMC_DETSI 5
34 : : MODULE m_CtqmcoffdiagComplex
35 : :
36 : : USE m_Global
37 : : USE m_GreenHyboffdiagComplex
38 : : USE m_BathOperatoroffdiagComplex
39 : : USE m_ImpurityOperator
40 : : USE m_Stat
41 : : USE m_FFTHyb
42 : : USE m_OurRng
43 : : use defs_basis
44 : : use m_abicore
45 : : use m_io_tools, only : open_file
46 : :
47 : : #ifdef HAVE_MPI2
48 : : USE mpi
49 : : #endif
50 : : IMPLICIT NONE
51 : :
52 : : public :: CtqmcoffdiagComplex_init
53 : : public :: CtqmcoffdiagComplex_setParameters
54 : : public :: CtqmcoffdiagComplex_setSweeps
55 : : public :: CtqmcoffdiagComplex_setSeed
56 : : public :: CtqmcoffdiagComplex_allocateAll
57 : : public :: CtqmcoffdiagComplex_allocateOpt
58 : : public :: CtqmcoffdiagComplex_setG0wTab
59 : : public :: CtqmcoffdiagComplex_setU
60 : : public :: CtqmcoffdiagComplex_clear
61 : : public :: CtqmcoffdiagComplex_reset
62 : : public :: CtqmcoffdiagComplex_setMu
63 : : public :: CtqmcoffdiagComplex_computeF
64 : : public :: CtqmcoffdiagComplex_run
65 : : public :: CtqmcoffdiagComplex_tryAddRemove
66 : : public :: CtqmcoffdiagComplex_trySwap
67 : : public :: CtqmcoffdiagComplex_measN
68 : : public :: CtqmcoffdiagComplex_measCorrelation
69 : : public :: CtqmcoffdiagComplex_measPerturbation
70 : : public :: CtqmcoffdiagComplex_getResult
71 : : public :: CtqmcoffdiagComplex_symmetrizeGreen
72 : : public :: CtqmcoffdiagComplex_getGreen
73 : : public :: CtqmcoffdiagComplex_getD
74 : : public :: CtqmcoffdiagComplex_getE
75 : : public :: CtqmcoffdiagComplex_printAll
76 : : public :: CtqmcoffdiagComplex_printQMC
77 : : public :: CtqmcoffdiagComplex_printGreen
78 : : public :: CtqmcoffdiagComplex_printD
79 : : public :: CtqmcoffdiagComplex_printE
80 : : public :: CtqmcoffdiagComplex_printPerturbation
81 : : public :: CtqmcoffdiagComplex_printCorrelation
82 : : public :: CtqmcoffdiagComplex_printSpectra
83 : : public :: CtqmcoffdiagComplex_destroy
84 : : public :: CtqmcoffdiagComplex_setMagmom
85 : :
86 : : !!***
87 : :
88 : : !!****t* m_CtqmcoffdiagComplex/CtqmcoffdiagComplex
89 : : !! NAME
90 : : !! CtqmcoffdiagComplex
91 : : !!
92 : : !! FUNCTION
93 : : !! This structured datatype contains the necessary data
94 : : !!
95 : : !! COPYRIGHT
96 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
97 : : !! This file is distributed under the terms of the
98 : : !! GNU General Public License, see ~abinit/COPYING
99 : : !! or http://www.gnu.org/copyleft/gpl.txt .
100 : : !!
101 : : !! SOURCE
102 : :
103 : : TYPE CtqmcoffdiagComplex
104 : :
105 : : LOGICAL :: init = .FALSE.
106 : : ! Flag: is MC initialized
107 : :
108 : : LOGICAL :: set = .FALSE.
109 : : ! Flag: ??
110 : :
111 : : LOGICAL :: setU = .FALSE.
112 : : ! Flag: is U Set ?
113 : :
114 : : LOGICAL :: inF = .FALSE.
115 : : ! Flag: is hybridization fct in input ?
116 : :
117 : : LOGICAL :: done = .FALSE.
118 : : ! Flag: is MC terminated ?
119 : :
120 : : LOGICAL :: para = .FALSE.
121 : : ! Flag: do we have parameters in input
122 : :
123 : : LOGICAL :: have_MPI = .FALSE.
124 : : ! Flag:
125 : :
126 : : INTEGER :: opt_movie = 0
127 : : !
128 : :
129 : : INTEGER :: opt_analysis = 0
130 : : ! correlations
131 : :
132 : : INTEGER :: opt_check = 0
133 : : ! various check 0
134 : : ! various check 1 impurity
135 : : ! various check 2 bath
136 : : ! various check 3 both
137 : :
138 : : INTEGER :: opt_order = 0
139 : : ! nb of segments max for analysis
140 : :
141 : : INTEGER :: opt_histo = 0
142 : : ! Enable histo calc.
143 : :
144 : : INTEGER :: opt_noise = 0
145 : : ! compute noise
146 : :
147 : : INTEGER :: opt_spectra = 0
148 : : ! markov chain FT (correlation time)
149 : :
150 : : INTEGER :: opt_levels = 0
151 : : ! do we have energy levels
152 : :
153 : : INTEGER :: opt_hybri_limit = 0
154 : : ! do we have limit of hybridization (yes=1)
155 : :
156 : : INTEGER :: opt_nondiag = 0
157 : : ! if opt_nondiag = 1 F is non diagonal.
158 : :
159 : : INTEGER :: prtopt = 1
160 : : ! printing
161 : :
162 : : INTEGER :: flavors
163 : : ! number of flavors
164 : :
165 : : INTEGER :: nspinor
166 : : ! number of spinor
167 : :
168 : : INTEGER :: measurements
169 : : ! The modulo used to measure the interaction energy and the number of electrons. Example : 2 means the measure is perform every two sweeps.
170 : :
171 : : INTEGER :: samples
172 : : ! nb of L points (dmftqmc_l)
173 : :
174 : : ! INTEGER :: endDensity
175 : : ! second dimension of density allocatable array
176 : :
177 : : INTEGER(8) :: seed
178 : : !
179 : :
180 : : INTEGER :: sweeps
181 : : !
182 : :
183 : : INTEGER :: thermalization
184 : : !
185 : :
186 : : INTEGER :: ostream
187 : : ! output file
188 : :
189 : : INTEGER :: istream
190 : : ! input file
191 : :
192 : : INTEGER :: modNoise1
193 : : ! measure the noise each modNoise1
194 : :
195 : : INTEGER :: modNoise2
196 : : ! measure the noise each modNoise2
197 : :
198 : : INTEGER :: activeFlavor
199 : : ! orbital on which one do sth now
200 : :
201 : : INTEGER, DIMENSION(1:2) :: modGlobalMove
202 : : ! 1: global move each modglobalmove(1)
203 : : ! 2: we have done modglobalmove(2) for two different orbitals.
204 : :
205 : : INTEGER :: Wmax
206 : : ! Max freq for FT
207 : :
208 : : DOUBLE PRECISION, DIMENSION(1:6) :: stats
209 : : ! to now how many negative determinant, antisegments,seeme.e.twfs...j
210 : :
211 : : DOUBLE PRECISION :: swap
212 : : ! nb of successfull GM
213 : :
214 : : DOUBLE PRECISION :: signvalue
215 : :
216 : : COMPLEX(KIND=8) :: phasevalue
217 : :
218 : : INTEGER :: MY_COMM
219 : : !
220 : :
221 : : INTEGER :: rank
222 : : !
223 : :
224 : : INTEGER :: size
225 : : ! size of MY_COMM
226 : :
227 : : DOUBLE PRECISION :: runTime ! time for the run routine
228 : : !
229 : :
230 : : DOUBLE PRECISION :: beta
231 : : !
232 : :
233 : : DOUBLE PRECISION :: U
234 : :
235 : : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:) :: mu
236 : : ! levels
237 : :
238 : : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:,:) :: hybri_limit
239 : : ! coeff A such that F=-A/(iwn)
240 : :
241 : : TYPE(GreenHyboffdiagComplex) :: Greens
242 : : ! Green's function
243 : :
244 : : ! DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: measN
245 : : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:,:) :: measN
246 : : ! measure of occupations (3or4,flavor)
247 : :
248 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) :: measDE
249 : : ! (flavor,flavor) double occupancies
250 : : ! (1,1): total energy of correlation.
251 : :
252 : : DOUBLE PRECISION :: a_Noise
253 : : ! Noise a exp (-bx) for the noise
254 : :
255 : : DOUBLE PRECISION :: b_Noise
256 : : ! Noise a exp (-bx) for the noise
257 : :
258 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: abNoiseG !(ab,tau,flavor)
259 : : ! Noise but for G
260 : :
261 : : TYPE(Vector) , DIMENSION(1:2) :: measNoise
262 : : TYPE(Vector), ALLOCATABLE, DIMENSION(:,:,:) :: measNoiseG !(tau,flavor,mod)
263 : : ! accumulate each value relataed to measurenoise 1 2
264 : :
265 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
266 : : ! INTEGER :: order
267 : : DOUBLE PRECISION :: inv_dt
268 : : ! 1/(beta/L)
269 : :
270 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) :: measPerturbation
271 : : ! opt_order,nflavor
272 : :
273 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: occup_histo_time
274 : : ! nflavor
275 : :
276 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: occupconfig
277 : : ! 2**nflavor
278 : :
279 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: suscep
280 : : ! samples
281 : :
282 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: chi
283 : : ! samples
284 : :
285 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: chicharge
286 : : ! samples
287 : :
288 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: ntot
289 : : ! occupation total, t2g, eg
290 : :
291 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) :: meas_fullemptylines
292 : : ! opt_order,nflavor
293 : :
294 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: measCorrelation
295 : : ! segment,antisegment,nflavor,nflavor
296 : :
297 : : !#endif
298 : : !#ifdef CTCtqmcoffdiagComplex_CHECK
299 : : DOUBLE PRECISION :: errorImpurity
300 : : ! check
301 : :
302 : : DOUBLE PRECISION :: errorBath
303 : : ! for check
304 : :
305 : : !#endif
306 : : TYPE(BathOperatoroffdiagComplex) :: Bath
307 : :
308 : :
309 : : TYPE(ImpurityOperator) :: Impurity
310 : :
311 : : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: density
312 : :
313 : : END TYPE CtqmcoffdiagComplex
314 : : !!***
315 : :
316 : : !INTERFACE CtqmcoffdiagComplex_setG0w
317 : : ! MODULE PROCEDURE CtqmcoffdiagComplex_setG0wFile, CtqmcoffdiagComplex_setG0wTab
318 : : !END INTERFACE
319 : :
320 : : CONTAINS
321 : : !!***
322 : :
323 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_init
324 : : !! NAME
325 : : !! CtqmcoffdiagComplex_init
326 : : !!
327 : : !! FUNCTION
328 : : !! Initialize the type CtqmcoffdiagComplex
329 : : !! Allocate all the non optional variables
330 : : !!
331 : : !! COPYRIGHT
332 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
333 : : !! This file is distributed under the terms of the
334 : : !! GNU General Public License, see ~abinit/COPYING
335 : : !! or http://www.gnu.org/copyleft/gpl.txt .
336 : : !!
337 : : !! INPUTS
338 : : !! op=ctqmc
339 : : !! ostream=where to write
340 : : !! istream=where to read the input parameters if so
341 : : !! bFile=logical argument True if input is read from istream
342 : : !! MY_COMM=mpi communicator for the CTQMC
343 : : !! iBuffer=input parameters if bFile is false
344 : : !!
345 : : !! OUTPUT
346 : : !!
347 : : !! SIDE EFFECTS
348 : : !!
349 : : !! NOTES
350 : : !!
351 : : !! SOURCE
352 : :
353 : 0 : SUBROUTINE CtqmcoffdiagComplex_init(op, ostream, istream, bFile, MY_COMM, iBuffer)
354 : :
355 : :
356 : : #ifdef HAVE_MPI1
357 : : include 'mpif.h'
358 : : #endif
359 : : !Arguments ------------------------------------
360 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
361 : : INTEGER , INTENT(IN ) :: ostream
362 : : INTEGER , INTENT(IN ) :: istream
363 : : LOGICAL , INTENT(IN ) :: bFile
364 : : DOUBLE PRECISION, DIMENSION(1:11), OPTIONAL, INTENT(IN) :: iBuffer
365 : : INTEGER , OPTIONAL, INTENT(IN ) :: MY_COMM
366 : : !Local variables ------------------------------
367 : : #ifdef HAVE_MPI
368 : : INTEGER :: ierr
369 : : #endif
370 : : !INTEGER :: iflavor
371 : : #ifdef __GFORTRAN__
372 : : ! INTEGER :: pid
373 : : ! CHARACTER(LEN=5) :: Cpid
374 : : !
375 : : #endif
376 : : DOUBLE PRECISION, DIMENSION(1:11) :: buffer
377 : :
378 : 0 : op%ostream = ostream
379 : 0 : op%istream = istream
380 : :
381 : : ! --- RENICE ---
382 : : !#ifdef __GFORTRAN__
383 : : ! pid = GetPid()
384 : : ! WRITE(Cpid,'(I5)') pid
385 : : ! CALL SYSTEM('renice +19 '//TRIM(ADJUSTL(Cpid))//' > /dev/null')
386 : : !#endif
387 : : !! --- RENICE ---
388 : :
389 [ # # ]: 0 : IF ( PRESENT(MY_COMM)) THEN
390 : : #ifdef HAVE_MPI
391 : 0 : op%have_MPI = .TRUE.
392 : 0 : op%MY_COMM = MY_COMM
393 : 0 : CALL MPI_Comm_rank(op%MY_COMM, op%rank, ierr)
394 : 0 : CALL MPI_Comm_size(op%MY_COMM, op%size, ierr)
395 : : #else
396 : : CALL WARN("MPI is not used ")
397 : : op%have_MPI = .FALSE.
398 : : op%MY_COMM = -1
399 : : op%rank = 0
400 : : op%size = 1
401 : : #endif
402 : : ELSE
403 : 0 : op%have_MPI = .FALSE.
404 : 0 : op%MY_COMM = -1
405 : 0 : op%rank = 0
406 : 0 : op%size = 1
407 : : END IF
408 : :
409 : : !IF ( op%rank .EQ. 0 ) THEN
410 : : ! WRITE(ostream,'(A20)') 'Job reniced with +19'
411 : : !CALL FLUSH(ostream)
412 : : !END IF
413 : :
414 [ # # ]: 0 : IF ( bFile .EQV. .TRUE. ) THEN
415 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
416 : :
417 : 0 : READ(istream,*) buffer(1) !iseed
418 : 0 : READ(istream,*) buffer(2) !op%sweeps
419 : 0 : READ(istream,*) buffer(3) !op%thermalization
420 : 0 : READ(istream,*) buffer(4) !op%measurements
421 : 0 : READ(istream,*) buffer(5) !op%flavors
422 : 0 : READ(istream,*) buffer(6) !op%samples
423 : 0 : READ(istream,*) buffer(7) !op%beta
424 : 0 : READ(istream,*) buffer(8) !U
425 : 0 : READ(istream,*) buffer(9) !iTech
426 : 0 : READ(istream,*) buffer(11)!op%nspinor
427 : : !READ(istream,*) buffer(9) !Wmax
428 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
429 : : !READ(istream,*) buffer(10) !order
430 : : !#endif
431 : : END IF
432 : :
433 : : #ifdef HAVE_MPI
434 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE. ) &
435 : : CALL MPI_Bcast(buffer, 11, MPI_DOUBLE_PRECISION, 0, &
436 : 0 : op%MY_COMM, ierr)
437 : : #endif
438 [ # # ]: 0 : ELSE IF ( PRESENT(iBuffer) ) THEN
439 : 0 : buffer(1:11) = iBuffer(1:11)
440 : : ELSE
441 : 0 : CALL ERROR("CtqmcoffdiagComplex_init : No input parameters ")
442 : : END IF
443 : :
444 : 0 : CALL CtqmcoffdiagComplex_setParameters(op, buffer)
445 : :
446 : 0 : CALL CtqmcoffdiagComplex_allocateAll(op)
447 : :
448 : : CALL GreenHyboffdiagComplex_init(op%Greens,op%samples, op%beta,INT(buffer(5)), &
449 : 0 : iTech=INT(buffer(9)),MY_COMM=op%MY_COMM)
450 : :
451 : :
452 : : ! op%seg_added = 0.d0
453 : : ! op%anti_added = 0.d0
454 : : ! op%seg_removed = 0.d0
455 : : ! op%anti_removed = 0.d0
456 : : ! op%seg_sign = 0.d0
457 : : ! op%anti_sign = 0.d0
458 [ # # ]: 0 : op%stats(:) = 0.d0
459 : : ! write(std_out,*) "op%stats",op%stats
460 : 0 : op%signvalue = 1.d0
461 : 0 : op%phasevalue = cmplx(1.d0,0.d0,kind=8)
462 : : ! op%signvaluecurrent = 0.d0
463 : : ! op%signvaluemeas = 0.d0
464 : 0 : op%swap = 0.d0
465 : 0 : op%runTime = 0.d0
466 : :
467 : 0 : CALL Vector_init(op%measNoise(1),op%sweeps/op%modNoise1)
468 : 0 : CALL Vector_init(op%measNoise(2),(op%sweeps/op%modNoise1+1)*CTQMC_SLICE2)
469 : : !CALL Vector_init(op%measNoise(3),101)
470 : : !CALL Vector_init(op%measNoise(4),101)
471 : :
472 [ # # # # ]: 0 : op%set = op%para .AND. op%inF
473 : 0 : op%done = .FALSE.
474 : 0 : op%init = .TRUE.
475 : :
476 : : !#ifdef CTCtqmcoffdiagComplex_CHECK
477 : 0 : op%errorImpurity = 0.d0
478 : 0 : op%errorBath = 0.d0
479 : : !#endif
480 : 0 : END SUBROUTINE CtqmcoffdiagComplex_init
481 : : !!***
482 : :
483 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setParameters
484 : : !! NAME
485 : : !! CtqmcoffdiagComplex_setParameters
486 : : !!
487 : : !! FUNCTION
488 : : !! set all parameters and operators
489 : : !!
490 : : !! COPYRIGHT
491 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
492 : : !! This file is distributed under the terms of the
493 : : !! GNU General Public License, see ~abinit/COPYING
494 : : !! or http://www.gnu.org/copyleft/gpl.txt .
495 : : !!
496 : : !! INPUTS
497 : : !! op=ctqmc
498 : : !! buffer=input parameters
499 : : !!
500 : : !! OUTPUT
501 : : !!
502 : : !! SIDE EFFECTS
503 : : !!
504 : : !! NOTES
505 : : !!
506 : : !! SOURCE
507 : :
508 : 0 : SUBROUTINE CtqmcoffdiagComplex_setParameters(op,buffer)
509 : :
510 : : !Arguments ------------------------------------
511 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
512 : : DOUBLE PRECISION, DIMENSION(1:11), INTENT(IN ) :: buffer
513 : :
514 : :
515 : 0 : op%thermalization = INT(buffer(3)) !op%thermalization
516 : 0 : CALL CtqmcoffdiagComplex_setSeed(op,INT(buffer(1)))
517 : 0 : CALL CtqmcoffdiagComplex_setSweeps(op,buffer(2))
518 : :
519 : 0 : op%measurements = INT(buffer(4)) !op%measurements
520 : 0 : op%flavors = INT(buffer(5))
521 : 0 : op%samples = INT(buffer(6)) !op%samples
522 : 0 : op%beta = buffer(7) !op%beta
523 : 0 : op%U = buffer(8) !U
524 : 0 : op%opt_nondiag = INT(buffer(10))
525 : 0 : op%nspinor = INT(buffer(11)) !op%nspinor
526 : : ! op%mu = buffer(9) !op%mu
527 : : !op%Wmax = INT(buffer(9)) !Freq
528 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
529 : : ! op%order = INT(buffer(10)) ! order
530 : 0 : op%inv_dt = op%samples / op%beta
531 : : !#endif
532 : :
533 : : !CALL ImpurityOperator_init(op%Impurity,op%flavors,op%beta, op%samples)
534 : 0 : CALL ImpurityOperator_init(op%Impurity,op%flavors,op%beta)
535 [ # # ]: 0 : IF ( op%U .GE. 0.d0 ) THEN
536 : 0 : CALL ImpurityOperator_computeU(op%Impurity,op%U,0.d0)
537 : 0 : op%setU = .TRUE.
538 : : END IF
539 : : ! op%mu = op%mu + op%Impurity%shift_mu
540 : : !sui!write(std_out,*) "op%opt_nondiag",op%opt_nondiag
541 : 0 : CALL BathOperatoroffdiagComplex_init(op%Bath, op%flavors, op%samples, op%beta, INT(buffer(9)), op%opt_nondiag)
542 : :
543 : 0 : op%para = .TRUE.
544 : :
545 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setParameters
546 : : !!***
547 : :
548 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setSweeps
549 : : !! NAME
550 : : !! CtqmcoffdiagComplex_setSweeps
551 : : !!
552 : : !! FUNCTION
553 : : !! set the number of sweeps
554 : : !!
555 : : !! COPYRIGHT
556 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
557 : : !! This file is distributed under the terms of the
558 : : !! GNU General Public License, see ~abinit/COPYING
559 : : !! or http://www.gnu.org/copyleft/gpl.txt .
560 : : !!
561 : : !! INPUTS
562 : : !! op=ctqmc
563 : : !! sweeps=asked sweeps
564 : : !!
565 : : !! OUTPUT
566 : : !!
567 : : !! SIDE EFFECTS
568 : : !!
569 : : !! NOTES
570 : : !!
571 : : !! SOURCE
572 : :
573 : 0 : SUBROUTINE CtqmcoffdiagComplex_setSweeps(op,sweeps)
574 : :
575 : : !Arguments ------------------------------------
576 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
577 : : DOUBLE PRECISION , INTENT(IN ) :: sweeps
578 : :
579 : 0 : op%sweeps = NINT(sweeps / DBLE(op%size))
580 : : ! !write(std_out,*) op%sweeps,NINT(sweeps / DBLE(op%size)),ANINT(sweeps/DBLE(op%size))
581 [ # # ]: 0 : IF ( DBLE(op%sweeps) .NE. ANINT(sweeps/DBLE(op%size)) ) &
582 : 0 : CALL ERROR("CtqmcoffdiagComplex_setSweeps : sweeps is negative or too big ")
583 [ # # ]: 0 : IF ( op%sweeps .LT. 2*CTQMC_SLICE1 ) THEN !202
584 : 0 : CALL WARNALL("CtqmcoffdiagComplex_setSweeps : # sweeps automtically changed ")
585 : 0 : op%sweeps = 2*CTQMC_SLICE1
586 : : ! ELSE IF ( op%sweeps .LT. op%thermalization ) THEN
587 : : ! CALL WARNALL("CtqmcoffdiagComplex_setSweeps : Thermalization > sweeps / cpu -> auto fix")
588 : : ! op%sweeps = op%thermalization
589 : : END IF
590 : 0 : IF ( DBLE(NINT(DBLE(op%sweeps)*DBLE(op%size)/DBLE(CTQMC_SLICE1))) .NE. &
591 [ # # ]: 0 : ANINT(DBLE(op%sweeps)*DBLE(op%size)/DBLE(CTQMC_SLICE1)) ) THEN
592 : 0 : op%modNoise1 = op%sweeps
593 : : ELSE
594 : 0 : op%modNoise1 = MIN(op%sweeps,INT(DBLE(op%sweeps)*DBLE(op%size) / DBLE(CTQMC_SLICE1))) !101
595 : : END IF
596 : 0 : op%modNoise2 = MAX(op%modNoise1 / CTQMC_SLICE2, 1) ! 100
597 : : ! op%modGlobalMove(1) = op%thermalization / 10 + 1
598 : : ! op%modGlobalMove(2) = 0
599 : :
600 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setSweeps
601 : : !!***
602 : :
603 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setSeed
604 : : !! NAME
605 : : !! CtqmcoffdiagComplex_setSeed
606 : : !!
607 : : !! FUNCTION
608 : : !! initialize random number generator
609 : : !!
610 : : !! COPYRIGHT
611 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
612 : : !! This file is distributed under the terms of the
613 : : !! GNU General Public License, see ~abinit/COPYING
614 : : !! or http://www.gnu.org/copyleft/gpl.txt .
615 : : !!
616 : : !! INPUTS
617 : : !! op=ctqmc
618 : : !! iseed=seed from imput
619 : : !!
620 : : !! OUTPUT
621 : : !!
622 : : !! SIDE EFFECTS
623 : : !!
624 : : !! NOTES
625 : : !!
626 : : !! SOURCE
627 : :
628 : 0 : SUBROUTINE CtqmcoffdiagComplex_setSeed(op,iseed)
629 : :
630 : : !Arguments ------------------------------------
631 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
632 : : INTEGER , INTENT(IN ) :: iseed
633 : : !Local variables ------------------------------
634 : : !INTEGER :: n
635 : : !INTEGER :: i
636 : : !INTEGER, DIMENSION(:), ALLOCATABLE :: seed
637 : :
638 : :
639 : : !CALL RANDOM_SEED(size = n)
640 : : !MALLOC(seed,(n))
641 : : !seed = iseed + (/ (i - 1, i = 1, n) /)
642 : :
643 : : !CALL RANDOM_SEED(PUT = seed+op%rank)
644 : :
645 : : !FREE(seed)
646 : :
647 : 0 : op%seed=INT(iseed+op%rank,8)
648 : :
649 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setSeed
650 : : !!***
651 : :
652 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_allocateAll
653 : : !! NAME
654 : : !! CtqmcoffdiagComplex_allocateAll
655 : : !!
656 : : !! FUNCTION
657 : : !! Allocate all non option variables
658 : : !!
659 : : !! COPYRIGHT
660 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
661 : : !! This file is distributed under the terms of the
662 : : !! GNU General Public License, see ~abinit/COPYING
663 : : !! or http://www.gnu.org/copyleft/gpl.txt .
664 : : !!
665 : : !! INPUTS
666 : : !! op=ctqmc
667 : : !!
668 : : !! OUTPUT
669 : : !!
670 : : !! SIDE EFFECTS
671 : : !!
672 : : !! NOTES
673 : : !!
674 : : !! SOURCE
675 : :
676 : 0 : SUBROUTINE CtqmcoffdiagComplex_allocateAll(op)
677 : :
678 : : !Arguments ------------------------------------
679 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
680 : : !Local variables ------------------------------
681 : : INTEGER :: flavors
682 : :
683 [ # # ]: 0 : IF ( .NOT. op%para ) &
684 : 0 : CALL ERROR("CtqmcoffdiagComplex_allocateAll : CtqmcoffdiagComplex_setParameters never called ")
685 : :
686 : 0 : flavors = op%flavors
687 : :
688 : :
689 : : ! number of electrons
690 [ # # ]: 0 : FREEIF(op%measN)
691 [ # # # # : 0 : MALLOC(op%measN,(1:4,1:flavors))
# # # # #
# # # #
# ]
692 [ # # # # ]: 0 : op%measN = cmplx(0.d0,0.d0,kind=8)
693 : :
694 : : ! double occupancies
695 [ # # ]: 0 : FREEIF(op%measDE)
696 [ # # # # : 0 : MALLOC(op%measDE,(1:flavors,1:flavors) )
# # # # #
# # # # #
# # # # ]
697 [ # # # # ]: 0 : op%measDE = real(cmplx(0.d0,0.d0,kind=8))
698 : :
699 [ # # ]: 0 : FREEIF(op%mu)
700 : : #ifdef FC_LLVM
701 : : ! LLVM 16 doesn't recognize this macro here
702 : : MALLOC(op%mu, (1:flavors) )
703 : : #else
704 [ # # # # : 0 : MALLOC(op%mu, (1:flavors))
# # # # #
# # # #
# ]
705 : : #endif
706 [ # # ]: 0 : op%mu = 0.d0
707 [ # # ]: 0 : FREEIF(op%hybri_limit)
708 : : #ifdef FC_LLVM
709 : : ! LLVM 16 doesn't recognize this macro here
710 : : MALLOC(op%hybri_limit, (flavors,flavors) )
711 : : #else
712 [ # # # # : 0 : MALLOC(op%hybri_limit, (flavors,flavors))
# # # # #
# # # # #
# # # # ]
713 : : #endif
714 [ # # # # ]: 0 : op%hybri_limit = czero
715 : 0 : END SUBROUTINE CtqmcoffdiagComplex_allocateAll
716 : : !!***
717 : :
718 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
719 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_allocateOpt
720 : : !! NAME
721 : : !! CtqmcoffdiagComplex_allocateOpt
722 : : !!
723 : : !! FUNCTION
724 : : !! allocate all option variables
725 : : !!
726 : : !! COPYRIGHT
727 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
728 : : !! This file is distributed under the terms of the
729 : : !! GNU General Public License, see ~abinit/COPYING
730 : : !! or http://www.gnu.org/copyleft/gpl.txt .
731 : : !!
732 : : !! INPUTS
733 : : !! op=ctqmc
734 : : !!
735 : : !! OUTPUT
736 : : !!
737 : : !! SIDE EFFECTS
738 : : !!
739 : : !! NOTES
740 : : !!
741 : : !! SOURCE
742 : :
743 : 0 : SUBROUTINE CtqmcoffdiagComplex_allocateOpt(op)
744 : :
745 : : !Arguments ------------------------------------
746 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
747 : : !Local variables ------------------------------
748 : : INTEGER :: i
749 : : INTEGER :: j
750 : : INTEGER :: k
751 : :
752 [ # # ]: 0 : IF ( .NOT. op%para ) &
753 : 0 : CALL ERROR("CtqmcoffdiagComplex_allocateOpt : CtqmcoffdiagComplex_setParameters never called ")
754 : :
755 [ # # ]: 0 : IF ( op%opt_analysis .EQ. 1 ) THEN
756 [ # # ]: 0 : FREEIF(op%measCorrelation)
757 [ # # # # : 0 : MALLOC(op%measCorrelation,(1:op%samples+1,1:3,1:op%flavors))
# # # # #
# # # # #
# # # # ]
758 [ # # # # : 0 : op%measCorrelation = 0.d0
# # ]
759 : : END IF
760 : :
761 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
762 [ # # ]: 0 : FREEIF(op%measPerturbation)
763 [ # # # # : 0 : MALLOC(op%measPerturbation,(1:op%opt_order,1:op%flavors))
# # # # #
# # # # #
# # # # ]
764 [ # # # # ]: 0 : op%measPerturbation = 0.d0
765 [ # # ]: 0 : FREEIF(op%meas_fullemptylines)
766 [ # # # # : 0 : MALLOC(op%meas_fullemptylines,(2,1:op%flavors))
# # # # #
# # # #
# ]
767 [ # # # # ]: 0 : op%meas_fullemptylines = 0.d0
768 : : END IF
769 : :
770 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
771 [ # # ]: 0 : FREEIF(op%occup_histo_time)
772 [ # # # # : 0 : MALLOC(op%occup_histo_time,(1:op%flavors+1))
# # # # #
# # # #
# ]
773 [ # # ]: 0 : op%occup_histo_time= 0.d0
774 [ # # ]: 0 : FREEIF(op%occupconfig)
775 [ # # # # : 0 : MALLOC(op%occupconfig,(1:2**op%flavors))
# # # # #
# # # # #
# # # # ]
776 [ # # ]: 0 : op%occupconfig= 0.d0
777 [ # # ]: 0 : FREEIF(op%suscep)
778 [ # # # # : 0 : MALLOC(op%suscep,(1:3,1:op%samples))
# # # # #
# # # #
# ]
779 [ # # # # ]: 0 : op%suscep= 0.d0
780 [ # # ]: 0 : FREEIF(op%chi)
781 [ # # # # : 0 : MALLOC(op%chi,(1:3,1:op%samples))
# # # # #
# # # #
# ]
782 [ # # # # ]: 0 : op%chi= 0.d0
783 [ # # ]: 0 : FREEIF(op%chicharge)
784 [ # # # # : 0 : MALLOC(op%chicharge,(1:3,1:op%samples))
# # # # #
# # # #
# ]
785 [ # # # # ]: 0 : op%chicharge= 0.d0
786 [ # # ]: 0 : FREEIF(op%ntot)
787 [ # # # # ]: 0 : MALLOC(op%ntot,(1:3))
788 [ # # ]: 0 : op%ntot= 0.d0
789 : : END IF
790 : :
791 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
792 [ # # ]: 0 : IF ( ALLOCATED(op%measNoiseG) ) THEN
793 [ # # ]: 0 : DO i=1,2
794 [ # # ]: 0 : DO j = 1, op%flavors
795 [ # # ]: 0 : DO k= 1, op%samples+1
796 : 0 : CALL Vector_destroy(op%measNoiseG(k,j,i))
797 : : END DO
798 : : END DO
799 : : END DO
800 [ # # # # : 0 : DT_FREE(op%measNoiseG)
# # # # ]
801 : : END IF
802 [ # # # # : 0 : DT_MALLOC(op%measNoiseG,(1:op%samples+1,1:op%flavors,1:2))
# # # # #
# # # # #
# # # # #
# # # # #
# # # # ]
803 : : !DO i=1,2
804 [ # # ]: 0 : DO j = 1, op%flavors
805 [ # # ]: 0 : DO k= 1, op%samples+1
806 : 0 : CALL Vector_init(op%measNoiseG(k,j,1),CTQMC_SLICE1)
807 : : END DO
808 : : END DO
809 [ # # ]: 0 : DO j = 1, op%flavors
810 [ # # ]: 0 : DO k= 1, op%samples+1
811 : 0 : CALL Vector_init(op%measNoiseG(k,j,2),CTQMC_SLICE1*CTQMC_SLICE2+1) ! +1 pour etre remplacer ceil
812 : : END DO
813 : : END DO
814 : : !END DO
815 [ # # ]: 0 : FREEIF(op%abNoiseG)
816 [ # # # # : 0 : MALLOC(op%aBNoiseG,(1:2,1:op%samples+1,op%flavors))
# # # # #
# # # # #
# # # # ]
817 [ # # # # : 0 : op%abNoiseG = 0.d0
# # ]
818 : : END IF
819 : :
820 [ # # ]: 0 : IF (op%opt_spectra .GE. 1 ) THEN
821 [ # # ]: 0 : FREEIF(op%density)
822 : : !MALLOC(op%density,(1:op%thermalization,1:op%flavors))
823 [ # # ]: 0 : i = CEILING(DBLE(op%thermalization+op%sweeps)/DBLE(op%measurements*op%opt_spectra))
824 [ # # # # : 0 : MALLOC(op%density,(1:op%flavors+1,1:i))
# # # # #
# # # # #
# # # # ]
825 [ # # # # ]: 0 : op%density = 0.d0
826 : : ! op%endDensity=i
827 : : !ENDPROBLEM
828 : :
829 : : END IF
830 : : !#endif
831 : 0 : END SUBROUTINE CtqmcoffdiagComplex_allocateOpt
832 : : !!***
833 : :
834 : : !SUBROUTINE CtqmcoffdiagComplex_setG0wFile(op,istream,opt_fk)
835 : : !include 'mpif.h'
836 : : !#endif
837 : : ! TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
838 : : ! INTEGER , INTENT(IN ) :: istream
839 : : ! INTEGER , INTENT(IN ) :: opt_fk
840 : : !! DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: F
841 : : ! COMPLEX*16 , DIMENSION(:,:), ALLOCATABLE :: Gomega
842 : : ! INTEGER :: flavors
843 : : ! INTEGER :: iflavor
844 : : ! INTEGER :: iomega
845 : : ! TYPE(GreenHyboffdiag) :: temp
846 : : !#ifdef HAVE_MPI
847 : : ! INTEGER :: ierr
848 : : !#endif
849 : : !
850 : : ! IF ( .NOT. op%para ) &
851 : : ! CALL ERROR("CtqmcoffdiagComplex_setG0wFile : CtqmcoffdiagComplex_setParameters never called ")
852 : : !
853 : : ! flavors = op%flavors
854 : : !
855 : : !! MALLOC(F,(1:op%samples+1,1:flavors))
856 : : ! MALLOC(Gomega,(1:op%Wmax,1:flavors))
857 : : ! IF ( op%have_MPI .EQV. .TRUE. ) THEN
858 : : ! CALL GreenHyboffdiag_init(temp, op%samples, op%beta, op%MY_COMM)
859 : : ! ELSE
860 : : ! CALL GreenHyboffdiag_init(temp, op%samples, op%beta)
861 : : ! END IF
862 : : !
863 : : ! DO iflavor = 1, flavors
864 : : ! IF ( op%rank .EQ. 0 ) THEN
865 : : ! DO iomega=1, op%Wmax
866 : : ! READ(istream,*) Gomega(iomega,iflavor)
867 : : ! END DO
868 : : ! END IF
869 : : ! END DO
870 : : !
871 : : !#ifdef HAVE_MPI
872 : : ! CALL MPI_Bcast(Gomega, op%Wmax*flavors, MPI_DOUBLE_COMPLEX, 0, &
873 : : ! op%MY_COMM, ierr)
874 : : !#endif
875 : : !
876 : : !! CALL GreenHyboffdiag_backFourier(temp,Gomega(:,1))
877 : : !! WRITE(10,*) Gomega(:,1)
878 : : !! CALL GreenHyboffdiag_print(temp)
879 : : !! CALL GreenHyboffdiag_forFourier (temp,Gomega(:,1))
880 : : !! WRITE(11,*) Gomega(:,1)
881 : : ! CALL CtqmcoffdiagComplex_setG0wTab(op,Gomega,opt_fk)
882 : : !! CALL CtqmcoffdiagComplex_computeF(op,Gomega, op%Wmax, F) ! mu is changed
883 : : ! FREE(Gomega)
884 : : !! CALL BathOperatoroffdiag_setF(op%Bath, F)
885 : : !! CALL BathOperatoroffdiag_printF(op%Bath)
886 : : !! FREE(F)
887 : : !
888 : : !! stop
889 : : !! op%inF = .TRUE.
890 : : !! op%set = .TRUE.
891 : : !
892 : : !END SUBROUTINE CtqmcoffdiagComplex_setG0wFile
893 : : !!***
894 : :
895 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setG0wTab
896 : : !! NAME
897 : : !! CtqmcoffdiagComplex_setG0wTab
898 : : !!
899 : : !! FUNCTION
900 : : !! Set Gow from input array
901 : : !!
902 : : !! COPYRIGHT
903 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
904 : : !! This file is distributed under the terms of the
905 : : !! GNU General Public License, see ~abinit/COPYING
906 : : !! or http://www.gnu.org/copyleft/gpl.txt .
907 : : !!
908 : : !! INPUTS
909 : : !! op=ctqmc
910 : : !! Gomega=G0w
911 : : !! opt_fk=F is already inversed with out iwn
912 : : !!
913 : : !! OUTPUT
914 : : !!
915 : : !! SIDE EFFECTS
916 : : !!
917 : : !! NOTES
918 : : !!
919 : : !! SOURCE
920 : :
921 [ # # ]: 0 : SUBROUTINE CtqmcoffdiagComplex_setG0wTab(op,Gomega,opt_fk,Iatom,fname)
922 : :
923 : : !Arguments ------------------------------------
924 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
925 : : COMPLEX(KIND=8), DIMENSION(:,:,:), INTENT(IN ) :: Gomega
926 : : INTEGER , INTENT(IN ) :: opt_fk
927 : : INTEGER, INTENT(IN) :: Iatom
928 : : CHARACTER(LEN=fnlen), INTENT(INOUT) :: fname
929 : : !Local variable -------------------------------
930 : : COMPLEX(KIND=8), DIMENSION(:,:,:), ALLOCATABLE :: F
931 : :
932 [ # # ]: 0 : IF ( .NOT. op%para ) &
933 : 0 : CALL ERROR("CtqmcoffdiagComplex_setG0wTab : CtqmcoffdiagComplex_setParameters never called ")
934 : :
935 [ # # # # : 0 : MALLOC(F,(1:op%samples+1,1:op%flavors,1:op%flavors))
# # # # #
# # # # #
# # # # #
# ]
936 : 0 : CALL CtqmcoffdiagComplex_computeF(op,Gomega, F, opt_fk,Iatom,fname) ! mu is changed
937 : : !write(6,*) "eee111"
938 : 0 : CALL BathOperatoroffdiagComplex_setF(op%Bath, F)
939 : : ! CALL BathOperatoroffdiag_printF(op%Bath,333)
940 : : !write(6,*) "eee"
941 [ # # ]: 0 : FREE(F)
942 : :
943 : 0 : op%inF = .TRUE.
944 : 0 : op%set = .TRUE.
945 : :
946 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setG0wTab
947 : : !!***
948 : :
949 : : !SUBROUTINE CtqmcoffdiagComplex_setFwK(op,Gomega)
950 : : ! COMPLEX*16 , DIMENSION(:,:), INTENT(IN ) :: Gomega
951 : : ! DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: F
952 : : !
953 : : ! IF ( .NOT. op%para ) &
954 : : ! CALL ERROR("CtqmcoffdiagComplex_setG0wTab : CtqmcoffdiagComplex_setParameters never called ")
955 : : !
956 : : ! MALLOC(F,(1:op%samples+1,1:op%flavors))
957 : : ! CALL CtqmcoffdiagComplex_computeFK(op,Gomega, op%Wmax, F) ! mu is changed
958 : : ! CALL BathOperatoroffdiag_setF(op%Bath, F)
959 : : ! CALL BathOperatoroffdiag_printF(op%Bath)
960 : : ! FREE(F)
961 : : !
962 : : ! op%inF = .TRUE.
963 : : ! op%set = .TRUE.
964 : : !
965 : : !END SUBROUTINE CtqmcoffdiagComplex_setFwK
966 : : !!***
967 : :
968 : : !SUBROUTINE CtqmcoffdiagComplex_setBand(op, mu, U)
969 : : ! DOUBLE PRECISION, INTENT(IN ) :: mu
970 : : ! DOUBLE PRECISION, INTENT(IN ) :: U
971 : : !
972 : : ! IF ( .NOT. op%para ) &
973 : : ! CALL ERROR("CtqmcoffdiagComplex_setBand : CtqmcoffdiagComplex_setParameters never called ")
974 : : !
975 : : ! CALL ImpurityOperator_setU(op%Impurity, U, 0.d0)
976 : : ! !op%mu = mu + op%Impurity%shift_mu
977 : : !END SUBROUTINE CtqmcoffdiagComplex_setBand
978 : : !!***
979 : :
980 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setU
981 : : !! NAME
982 : : !! CtqmcoffdiagComplex_setU
983 : : !!
984 : : !! FUNCTION
985 : : !! set the interaction matrix
986 : : !!
987 : : !! COPYRIGHT
988 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
989 : : !! This file is distributed under the terms of the
990 : : !! GNU General Public License, see ~abinit/COPYING
991 : : !! or http://www.gnu.org/copyleft/gpl.txt .
992 : : !!
993 : : !! INPUTS
994 : : !! op=ctqmc
995 : : !! matU=interaction matrix
996 : : !!
997 : : !! OUTPUT
998 : : !!
999 : : !! SIDE EFFECTS
1000 : : !!
1001 : : !! NOTES
1002 : : !!
1003 : : !! SOURCE
1004 : :
1005 [ # # ]: 0 : SUBROUTINE CtqmcoffdiagComplex_setU(op,matU)
1006 : :
1007 : : !Arguments ------------------------------------
1008 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) ::op
1009 : : !Local variables ------------------------------
1010 : : COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN) :: matU
1011 : :
1012 [ # # # # ]: 0 : IF ( SIZE(matU) .NE. op%flavors*op%flavors ) &
1013 : 0 : CALL ERROR("CtqmcoffdiagComplex_setU : Wrong interaction matrix (size) ")
1014 : :
1015 : 0 : CALL ImpurityOperator_setUmatComplex(op%Impurity, matU)
1016 : 0 : op%setU = .TRUE.
1017 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setU
1018 : : !!***
1019 : :
1020 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_clear
1021 : : !! NAME
1022 : : !! CtqmcoffdiagComplex_clear
1023 : : !!
1024 : : !! FUNCTION
1025 : : !! clear a ctqmc run
1026 : : !!
1027 : : !! COPYRIGHT
1028 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1029 : : !! This file is distributed under the terms of the
1030 : : !! GNU General Public License, see ~abinit/COPYING
1031 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1032 : : !!
1033 : : !! INPUTS
1034 : : !! op=ctqmc
1035 : : !!
1036 : : !! OUTPUT
1037 : : !!
1038 : : !! SIDE EFFECTS
1039 : : !!
1040 : : !! NOTES
1041 : : !!
1042 : : !! SOURCE
1043 : :
1044 : 0 : SUBROUTINE CtqmcoffdiagComplex_clear(op)
1045 : :
1046 : : !Arguments ------------------------------------
1047 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
1048 : : !Local variables ------------------------------
1049 : : INTEGER :: i
1050 : : INTEGER :: j
1051 : : INTEGER :: k
1052 : :
1053 [ # # ]: 0 : op%measN(1,:) = cmplx(0.d0,0.d0,kind=8)
1054 [ # # ]: 0 : op%measN(2,:) = cmplx(0.d0,0.d0,kind=8)
1055 : : !Do not set measN(3,:) to 0 to avoid erasing N between therm and ctqmc
1056 [ # # ]: 0 : op%measN(4,:) = cmplx(0.d0,0.d0,kind=8)
1057 [ # # # # ]: 0 : op%measDE = real(cmplx(0.d0,0.d0,kind=8))
1058 : : ! op%seg_added = 0.d0
1059 : : ! op%anti_added = 0.d0
1060 : : ! op%seg_removed = 0.d0
1061 : : ! op%anti_removed = 0.d0
1062 : : ! op%seg_sign = 0.d0
1063 : : ! op%anti_sign = 0.d0
1064 [ # # ]: 0 : op%stats(:) = 0.d0
1065 : : ! op%signvaluecurrent = 0.d0
1066 : : ! op%signvaluemeas = 0.d0
1067 : 0 : op%swap = 0.d0
1068 : 0 : op%runTime = 0.d0
1069 : 0 : op%modGlobalMove(2) = 0
1070 : 0 : CALL Vector_clear(op%measNoise(1))
1071 : 0 : CALL Vector_clear(op%measNoise(2))
1072 : : !#ifdef CTCtqmcoffdiagComplex_CHECK
1073 : 0 : op%errorImpurity = 0.d0
1074 : 0 : op%errorBath = 0.d0
1075 : : !#endif
1076 : 0 : CALL GreenHyboffdiagComplex_clear(op%Greens)
1077 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
1078 [ # # # # ]: 0 : IF ( op%opt_analysis .EQ. 1 .AND. ALLOCATED(op%measCorrelation) ) &
1079 [ # # # # : 0 : op%measCorrelation = 0.d0
# # ]
1080 [ # # # # ]: 0 : IF ( op%opt_order .GT. 0 .AND. ALLOCATED(op%measPerturbation) ) &
1081 [ # # # # ]: 0 : op%measPerturbation = 0.d0
1082 [ # # # # ]: 0 : IF ( op%opt_order .GT. 0 .AND. ALLOCATED(op%meas_fullemptylines) ) &
1083 [ # # # # ]: 0 : op%meas_fullemptylines = 0.d0
1084 [ # # # # ]: 0 : IF ( op%opt_noise .EQ. 1 .AND. ALLOCATED(op%measNoiseG) ) THEN
1085 [ # # ]: 0 : DO i=1,2
1086 [ # # ]: 0 : DO j = 1, op%flavors
1087 [ # # ]: 0 : DO k= 1, op%samples+1
1088 : 0 : CALL Vector_clear(op%measNoiseG(k,j,i))
1089 : : END DO
1090 : : END DO
1091 : : END DO
1092 : : !DO j = 1, op%flavors
1093 : : ! CALL GreenHyboffdiag_clear(op%Greens(j))
1094 : : !END DO
1095 : : END IF
1096 : : !#endif
1097 : 0 : END SUBROUTINE CtqmcoffdiagComplex_clear
1098 : : !!***
1099 : :
1100 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_reset
1101 : : !! NAME
1102 : : !! CtqmcoffdiagComplex_reset
1103 : : !!
1104 : : !! FUNCTION
1105 : : !! reset a ctqmc simulation
1106 : : !!
1107 : : !! COPYRIGHT
1108 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1109 : : !! This file is distributed under the terms of the
1110 : : !! GNU General Public License, see ~abinit/COPYING
1111 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1112 : : !!
1113 : : !! INPUTS
1114 : : !! op=ctqmc
1115 : : !!
1116 : : !! OUTPUT
1117 : : !!
1118 : : !! SIDE EFFECTS
1119 : : !!
1120 : : !! NOTES
1121 : : !!
1122 : : !! SOURCE
1123 : :
1124 : 0 : SUBROUTINE CtqmcoffdiagComplex_reset(op)
1125 : :
1126 : : !Arguments ------------------------------------
1127 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
1128 : : !Local variables ------------------------------
1129 : : !INTEGER :: iflavor
1130 : : DOUBLE PRECISION :: sweeps
1131 : :
1132 : 0 : CALL GreenHyboffdiagComplex_reset(op%Greens)
1133 : 0 : CALL CtqmcoffdiagComplex_clear(op)
1134 : 0 : CALL ImpurityOperator_reset(op%Impurity)
1135 : 0 : CALL BathOperatoroffdiagComplex_reset (op%Bath)
1136 [ # # ]: 0 : op%measN(3,:) = cmplx(0.d0,0.d0,kind=8)
1137 : : !complete restart -> measN=0
1138 : 0 : op%done = .FALSE.
1139 : 0 : op%set = .FALSE.
1140 : 0 : op%inF = .FALSE.
1141 : 0 : op%opt_movie = 0
1142 : 0 : op%opt_analysis = 0
1143 : 0 : op%opt_order = 0
1144 : 0 : op%opt_check = 0
1145 : 0 : op%opt_noise = 0
1146 : 0 : op%opt_spectra = 0
1147 : 0 : op%opt_levels = 0
1148 : 0 : sweeps = DBLE(op%sweeps)*DBLE(op%size)
1149 : 0 : CALL CtqmcoffdiagComplex_setSweeps(op, sweeps)
1150 : : !#ifdef HAVE_MPI
1151 : : ! CALL MPI_BARRIER(op%MY_COMM,iflavor)
1152 : : ! IF ( op%rank .EQ. 0 ) &
1153 : : !#endif
1154 : : ! WRITE(op%ostream,'(A9)') "QMC reset"
1155 : : ! CALL FLUSH(op%ostream)
1156 : 0 : END SUBROUTINE CtqmcoffdiagComplex_reset
1157 : : !!***
1158 : :
1159 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setMu
1160 : : !! NAME
1161 : : !! CtqmcoffdiagComplex_setMu
1162 : : !!
1163 : : !! FUNCTION
1164 : : !! impose energy levels
1165 : : !!
1166 : : !! COPYRIGHT
1167 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1168 : : !! This file is distributed under the terms of the
1169 : : !! GNU General Public License, see ~abinit/COPYING
1170 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1171 : : !!
1172 : : !! INPUTS
1173 : : !! op=ctqmc
1174 : : !! levels=energy levels vector
1175 : : !!
1176 : : !! OUTPUT
1177 : : !! argout(sizeout)=description
1178 : : !!
1179 : : !! SIDE EFFECTS
1180 : : !!
1181 : : !! NOTES
1182 : : !!
1183 : : !! SOURCE
1184 : :
1185 [ # # ]: 0 : SUBROUTINE CtqmcoffdiagComplex_setMu(op, levels)
1186 : :
1187 : : !Arguments ------------------------------------
1188 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
1189 : : COMPLEX(KIND=8), DIMENSION(:), INTENT(IN ) :: levels
1190 : :
1191 [ # # ]: 0 : IF ( op%flavors .NE. SIZE(levels,1) ) &
1192 : 0 : CALL WARNALL("CtqmcoffdiagComplex_setMu : Taking energy levels from weiss G(iw)")
1193 : :
1194 [ # # ]: 0 : op%mu(:)=-levels(:) ! levels = \epsilon_j - \mu
1195 : : !op%mu =\tilde{\mu} = \mu -\epsilon_j
1196 : 0 : op%opt_levels = 1
1197 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setMu
1198 : : !!***
1199 : :
1200 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_sethybri_limit
1201 : : !! NAME
1202 : : !! CtqmcoffdiagComplex_sethybri_limit
1203 : : !!
1204 : : !! FUNCTION
1205 : : !! use coefficient A such that F=-A/(iwn) given by DMFT code.
1206 : : !!
1207 : : !! COPYRIGHT
1208 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1209 : : !! This file is distributed under the terms of the
1210 : : !! GNU General Public License, see ~abinit/COPYING
1211 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1212 : : !!
1213 : : !! INPUTS
1214 : : !! hybri_limit(nflavor,nflavor)=contains the limit for each couple of flavors
1215 : : !!
1216 : : !! OUTPUT
1217 : : !! argout(sizeout)=description
1218 : : !!
1219 : : !! SIDE EFFECTS
1220 : : !! op(CtqmcoffdiagComplex_type) = is the ctqmc main variable
1221 : : !! op&limit is now filled
1222 : : !! NOTES
1223 : : !!
1224 : : !! SOURCE
1225 : :
1226 [ # # ]: 0 : SUBROUTINE CtqmcoffdiagComplex_sethybri_limit(op, hybri_limit,opthybri)
1227 : :
1228 : : !Arguments ------------------------------------
1229 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
1230 : : COMPLEX(KIND=8) , DIMENSION(:,:), INTENT(IN ) :: hybri_limit
1231 : : INTEGER, INTENT(IN) :: opthybri
1232 : : !----------------------------------------------
1233 [ # # ]: 0 : IF ( op%flavors .NE. SIZE(hybri_limit,1) ) &
1234 : 0 : CALL ERROR("Error in sethybri_limit")
1235 : :
1236 [ # # # # ]: 0 : op%hybri_limit(:,:)=hybri_limit(:,:)
1237 : 0 : op%opt_hybri_limit = opthybri
1238 : 0 : END SUBROUTINE CtqmcoffdiagComplex_sethybri_limit
1239 : : !!***
1240 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_computeF
1241 : : !! NAME
1242 : : !! CtqmcoffdiagComplex_computeF
1243 : : !!
1244 : : !! FUNCTION
1245 : : !! Compute the hybridization function
1246 : : !!
1247 : : !! COPYRIGHT
1248 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1249 : : !! This file is distributed under the terms of the
1250 : : !! GNU General Public License, see ~abinit/COPYING
1251 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1252 : : !!
1253 : : !! INPUTS
1254 : : !! op=ctqmc
1255 : : !! Gomega=G0 to compute F
1256 : : !! opt_fk=What is Gomega
1257 : : !!
1258 : : !! OUTPUT
1259 : : !! F=hybridization function
1260 : : !!
1261 : : !! SIDE EFFECTS
1262 : : !!
1263 : : !! NOTES
1264 : : !!
1265 : : !! SOURCE
1266 : :
1267 [ # # # # ]: 0 : SUBROUTINE CtqmcoffdiagComplex_computeF(op, Gomega, F, opt_fk,Iatom,fname)
1268 : :
1269 : : use m_hide_lapack, only : xginv
1270 : : !Arguments ------------------------------------
1271 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
1272 : : COMPLEX(KIND=8), DIMENSION(:,:,:), INTENT(IN ) :: Gomega
1273 : : !INTEGER , INTENT(IN ) :: Wmax
1274 : : COMPLEX(KIND=8), DIMENSION(:,:,:), INTENT(INOUT) :: F
1275 : : INTEGER , INTENT(IN ) :: opt_fk
1276 : : INTEGER, INTENT(IN) :: Iatom
1277 : : CHARACTER(LEN=fnlen), INTENT(INOUT) :: fname
1278 : : !Local variables ------------------------------
1279 : : INTEGER :: flavors
1280 : : INTEGER :: samples
1281 : : INTEGER :: iflavor
1282 : : INTEGER :: iflavor2
1283 : : INTEGER :: iomega
1284 : : INTEGER :: itau,ioerr,debug
1285 : : DOUBLE PRECISION :: pi_invBeta
1286 : : COMPLEX(KIND=8) :: K
1287 : : !DOUBLE PRECISION :: re
1288 : : !DOUBLE PRECISION :: im
1289 : : !DOUBLE PRECISION :: det
1290 : 0 : COMPLEX(KIND=8), DIMENSION(:,:,:), ALLOCATABLE :: F_omega
1291 : 0 : COMPLEX(KIND=8), DIMENSION(:,:), ALLOCATABLE :: F_omega_inv
1292 : 0 : COMPLEX(KIND=8), DIMENSION(:,:,:), ALLOCATABLE :: Gomega_tmp
1293 : 0 : TYPE(GreenHyboffdiagComplex) :: F_tmp
1294 : : CHARACTER(LEN=200) :: message
1295 : : character(len=2) :: atomnb
1296 : 0 : DOUBLE PRECISION, allocatable :: x_r(:,:,:), x_i(:,:,:)
1297 : : CHARACTER(LEN=256) :: dummy
1298 : : DOUBLE PRECISION :: tautemp
1299 : : !character(len=30) :: tmpfil
1300 : : !INTEGER :: unitnb
1301 : :
1302 : : ABI_UNUSED((/opt_fk/))
1303 : :
1304 : 0 : flavors = op%flavors
1305 : :
1306 : 0 : samples = op%samples
1307 : 0 : pi_invBeta = ACOS(-1.d0) / op%beta
1308 : 0 : op%Wmax=SIZE(Gomega,1)
1309 : :
1310 : 0 : debug = 0 !set to one to build G0
1311 : :
1312 : : !=================================
1313 : : ! --- Initialize F_tmp
1314 : : !=================================
1315 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE. ) THEN
1316 : 0 : CALL GreenHyboffdiagComplex_init(F_tmp,samples,op%beta,flavors,MY_COMM=op%MY_COMM)
1317 : : ELSE
1318 : 0 : CALL GreenHyboffdiagComplex_init(F_tmp,samples,op%beta,flavors)
1319 : : END IF
1320 : : ! K = op%mu
1321 : :
1322 : : !=================================
1323 : : ! --- Allocate F_omega
1324 : : !=================================
1325 [ # # # # : 0 : MALLOC(F_omega,(1:op%Wmax,1:flavors,1:flavors))
# # # # #
# # # # #
# # # # #
# ]
1326 [ # # # # : 0 : MALLOC(F_omega_inv,(1:flavors,1:flavors))
# # # # #
# # # # #
# # ]
1327 [ # # # # : 0 : MALLOC(Gomega_tmp,(1:op%Wmax,1:flavors,1:flavors))
# # # # #
# # # # #
# # ]
1328 [ # # # # : 0 : Gomega_tmp=Gomega
# # # # #
# # # # #
# # ]
1329 : :
1330 : :
1331 : : !This section of the code is not used / Only building G0 so just skip it
1332 : : if(debug == 1) then
1333 : : !=============================================================================================
1334 : : ! --- Compute Bath Green's function from Hybridization function in imaginary time
1335 : : !=============================================================================================
1336 : : ! In the following we want to get G_0(iw_n) from F(iw_n).
1337 : : ! Be careful : G_0 is called F and F is called G_0 (Why ????)
1338 : : !-----------------------------------------------------------
1339 : : ! G_0(iw_n)^-1 = iw_n + mu - F(iw_n)
1340 : : !-----------------------------------
1341 : : do iomega=1,op%Wmax
1342 : : do iflavor=1,flavors
1343 : : do iflavor2=1,flavors
1344 : : if (iflavor==iflavor2) then
1345 : : F_omega_inv(iflavor,iflavor2)= (cmplx(0.d0,(2.d0*DBLE(iomega)-1.d0) * pi_invBeta,kind=8) &
1346 : : & + op%mu(iflavor)- Gomega_tmp(iomega,iflavor,iflavor2))
1347 : : else
1348 : : F_omega_inv(iflavor,iflavor2)= (- Gomega_tmp(iomega,iflavor,iflavor2))
1349 : : endif
1350 : : enddo
1351 : : enddo
1352 : :
1353 : : !Inverse of G_0(iw_n)^-1 to get G_0(iw_n)
1354 : : !--------------------
1355 : : call xginv(F_omega_inv,flavors)
1356 : :
1357 : : do iflavor=1,flavors
1358 : : do iflavor2=1,flavors
1359 : : F_omega(iomega,iflavor,iflavor2) = F_omega_inv(iflavor,iflavor2)
1360 : : enddo
1361 : : enddo
1362 : : enddo !iomega
1363 : :
1364 : : CALL GreenHyboffdiagComplex_setOperW(F_tmp,F_omega)
1365 : :
1366 : : CALL GreenHyboffdiagComplex_backFourier(F_tmp,func="green")
1367 : :
1368 : : !Put the result in F(tau) aka G_0(tau).
1369 : : !--------------------------------------
1370 : : DO iflavor = 1, flavors
1371 : : DO iflavor2 = 1, flavors
1372 : : DO itau=1,samples+1
1373 : : ! This symetrization is general and valid even with SOC
1374 : : ! Without SOC, it is useless ?
1375 : : F(itau,iflavor,iflavor2) = ( F_tmp%oper(itau,iflavor,iflavor2) + &
1376 : : & dconjg(F_tmp%oper(itau,iflavor2,iflavor)) )/2.d0
1377 : : END DO
1378 : : END DO
1379 : : END DO
1380 : :
1381 : : !Print G_0(tau) in file
1382 : : !------------------------
1383 : : !Real part
1384 : : open (unit=4367,file='G0tau_fromF_Re',status='unknown',form='formatted')
1385 : : rewind(4367)
1386 : : IF ( op%rank .EQ. 0 ) THEN
1387 : : DO iflavor = 1, flavors
1388 : : DO iflavor2 = 1, flavors
1389 : : write(4367,*) "#",iflavor,iflavor2
1390 : : do itau=1,op%samples+1
1391 : : write(4367,*) (itau-1)*op%beta/(op%samples),real(F(itau,iflavor,iflavor2))
1392 : : enddo
1393 : : write(4367,*)
1394 : : END DO
1395 : : END DO
1396 : : ENDIF
1397 : : !call flush(4367) nag compiler problrem
1398 : : close(4367)
1399 : :
1400 : : !Imag part
1401 : : open (unit=4367,file='G0tau_fromF_Im',status='unknown',form='formatted')
1402 : : rewind(4367)
1403 : : IF ( op%rank .EQ. 0 ) THEN
1404 : : DO iflavor = 1, flavors
1405 : : DO iflavor2 = 1, flavors
1406 : : write(4367,*) "#",iflavor,iflavor2
1407 : : do itau=1,op%samples+1
1408 : : write(4367,*) (itau-1)*op%beta/(op%samples),aimag(F(itau,iflavor,iflavor2))
1409 : : enddo
1410 : : write(4367,*)
1411 : : END DO
1412 : : END DO
1413 : : ENDIF
1414 : : !call flush(4367) !nag compiler problem
1415 : : close(4367)
1416 : :
1417 : : endif ! end of building G0
1418 : :
1419 : 0 : call xmpi_barrier(op%MY_COMM)
1420 : :
1421 : : !=============================================================================================
1422 : : ! --- Restore Hybridization in F_omega
1423 : : !=============================================================================================
1424 : :
1425 : : ! Restore Hybridization in F_omega for the following operations
1426 [ # # # # : 0 : F_omega = Gomega_tmp
# # # # #
# # # # #
# # ]
1427 : : ! From now F_omega is actually F(iw_n) calculated in m_forctqmc.f90
1428 : : ! and corresponding to F(iw_n) = -[G_0(iw_n^-1 - iwn - levels)]
1429 : :
1430 : : !==================================================================
1431 : : ! --- Full double loop on flavors to compute F (remove levels)
1432 : : !==================================================================
1433 [ # # ]: 0 : DO iflavor = 1, flavors
1434 [ # # ]: 0 : DO iflavor2 = 1, flavors
1435 : :
1436 : : ! --- Compute or use the levels for the diagonal hybridization (else K=0)
1437 [ # # ]: 0 : IF(iflavor==iflavor2) THEN
1438 [ # # ]: 0 : IF ( op%opt_levels .EQ. 1 ) THEN
1439 : : K = op%mu(iflavor)
1440 : : ELSE
1441 : 0 : K = -(F_omega(op%Wmax,iflavor,iflavor))
1442 : 0 : op%mu(iflavor) = K
1443 : : END IF
1444 : : ELSE
1445 : : K=cmplx(0.d0,0.d0,kind=8)
1446 : : ENDIF
1447 : :
1448 : : ! --- compute residual K (?)
1449 : : ! K corresponds to the first moment of hybridization (d1 in setMuD1)
1450 : : ! ------------------------------------------------------------------
1451 : 0 : K = CMPLX(0,(2.d0*DBLE(op%Wmax)-1.d0)*pi_invBeta,8)*F_omega(op%Wmax,iflavor,iflavor2)
1452 : 0 : CALL GreenHyboffdiagComplex_setMuD1(op%Greens,iflavor,iflavor2,op%mu(iflavor),K)
1453 : : END DO
1454 : : END DO
1455 : :
1456 : : ! --- Creates F_tmp%oper_w
1457 : 0 : CALL GreenHyboffdiagComplex_setOperW(F_tmp,F_omega)
1458 : :
1459 : : ! For all iflavor and iflavor2, do the Fourier transformation to have F(tau)
1460 : : ! ---------------------------------------------
1461 [ # # ]: 0 : if (op%opt_hybri_limit .eq. 0) then
1462 : 0 : write(message,'(a,a)') " == Not using the asymptotic limit of hybridization to enforce F(iw_n) -> -C_ij/iw_n (dmft_hybri_limit = 0)", ch10
1463 : : else
1464 : : ! Take into account asymptotic limit of hybridization function such that F(iw_n) -> -C_ij/iw_n
1465 : : ! with C_ij calculated in m_forctqmc.f90
1466 : : ! --------------------------------------
1467 : 0 : write(message,'(a,a)') " == Use asymptotic limit of hybridization function such that F(iw_n) -> -C_ij/iw_n (dmft_hybri_limit = 1)", ch10
1468 : : endif
1469 : 0 : CALL wrtout(std_out,message,'COLL')
1470 : 0 : CALL GreenHyboffdiagComplex_backFourierComplex(F_tmp,hybri_limit=op%hybri_limit,opt_hybri_limit=op%opt_hybri_limit)
1471 : :
1472 : : ! --- Put the result in F
1473 [ # # ]: 0 : DO iflavor = 1, flavors
1474 [ # # ]: 0 : DO iflavor2 = 1, flavors
1475 [ # # ]: 0 : DO itau=1,samples+1
1476 : 0 : F(itau,iflavor,iflavor2) = -F_tmp%oper(samples+2-itau,iflavor,iflavor2)
1477 : : END DO
1478 : : END DO
1479 : : END DO
1480 : :
1481 [ # # ]: 0 : DO iflavor = 1, flavors
1482 [ # # ]: 0 : DO iflavor2 = 1, flavors
1483 [ # # ]: 0 : DO itau=1,samples+1
1484 : : ! This symetrization is general and valid even with SOC
1485 : : ! Without SOC, it leads to zero.
1486 : : F(itau,iflavor,iflavor2) = -( (F_tmp%oper(samples+2-itau,iflavor,iflavor2)) + &
1487 : 0 : & dconjg(F_tmp%oper(samples+2-itau,iflavor2,iflavor)) )/2.d0
1488 : : END DO
1489 : : END DO
1490 : : END DO
1491 : :
1492 : : !debug
1493 : : if (3 .eq. 4) then
1494 : : !== Read Hybridization function from file ==
1495 : : !reset previous hybridization
1496 : : write(message,'(a)') "Warning: About to read Hybridization from file!"
1497 : : call wrtout(std_out,message,'COLL')
1498 : : DO iflavor = 1, flavors
1499 : : DO iflavor2 = 1, flavors
1500 : : DO itau=1,samples+1
1501 : : F(itau,iflavor,iflavor2) = czero
1502 : : END DO
1503 : : END DO
1504 : : END DO
1505 : :
1506 : : if (op%rank .eq. 0 ) then
1507 : : ABI_MALLOC(x_r,(op%samples+1,flavors,flavors))
1508 : : ABI_MALLOC(x_i,(op%samples+1,flavors,flavors))
1509 : : open(unit=735,file='Hybridization.dat',status='old',form='formatted',action='read',iostat=ioerr)
1510 : :
1511 : : read(735, '(a)') dummy
1512 : :
1513 : : do itau=1, op%samples+1
1514 : : read(735, '(2x,393(e25.17e3,2x))') tautemp, &
1515 : : & ((x_r(itau,iflavor,iflavor2), x_i(itau,iflavor,iflavor2), iflavor=1, flavors), iflavor2=1, flavors)
1516 : : end do
1517 : :
1518 : : close(735)
1519 : :
1520 : : ! Now assign to F array (reverse the sign flip and time ordering)
1521 : : do itau=1, op%samples+1
1522 : : do iflavor=1, flavors
1523 : : do iflavor2=1, flavors
1524 : : F(op%samples+2-itau, iflavor, iflavor2) = &
1525 : : CMPLX(-x_r(itau,iflavor,iflavor2), -x_i(itau,iflavor,iflavor2), kind=8)
1526 : : enddo
1527 : : enddo
1528 : : enddo
1529 : : ABI_FREE(x_r)
1530 : : ABI_FREE(x_i)
1531 : : endif
1532 : : end if ! end reading F from file
1533 : :
1534 : : !== Write Hybridization function in file ==
1535 [ # # ]: 0 : If (Iatom .lt. 10) then
1536 : 0 : write(atomnb,'("0",i1)') Iatom
1537 : : else
1538 : 0 : write(atomnb,'(i2)') Iatom
1539 : : endif
1540 : :
1541 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
1542 : 0 : open (unit=735,file=trim(fname)//'_Hybridization_iatom_'//atomnb//'.dat',status='unknown',form='formatted')
1543 : 0 : write(735,'(6a)') "# Real and Imaginary part of the Hybridization function Delta(tau) in the CTQMC basis"
1544 [ # # ]: 0 : do itau=1,op%samples+1
1545 : 0 : write(735,'(2x,393(e25.17e3,2x))') DBLE(itau-1)*(op%beta/op%samples),&
1546 [ # # # # : 0 : &((-1*real(F(op%samples+2-itau,iflavor,iflavor2)),-1*AIMAG(F(op%samples+2-itau,iflavor,iflavor2)),iflavor=1,flavors),iflavor2=1,flavors)
# # # # ]
1547 : : enddo
1548 : 0 : close(735)
1549 : : ENDIF
1550 : :
1551 : :
1552 : : ! call xmpi_barrier(op%MY_COMM)
1553 : : !write(6,*) "QQQQ3"
1554 [ # # ]: 0 : FREE(Gomega_tmp)
1555 [ # # ]: 0 : FREE(F_omega)
1556 [ # # ]: 0 : FREE(F_omega_inv)
1557 : : !write(6,*) "QQQQ4"
1558 : 0 : CALL GreenHyboffdiagComplex_destroy(F_tmp)
1559 : : !write(6,*) "QQQQ2"
1560 : :
1561 : :
1562 [ # # # # : 0 : END SUBROUTINE CtqmcoffdiagComplex_computeF
# # # # #
# # # # #
# # # # #
# # # ]
1563 : : !!***
1564 : :
1565 : : !SUBROUTINE CtqmcoffdiagComplex_computeFK(op, Gomega, Wmax, F)
1566 : : ! COMPLEX*16 , DIMENSION(:,:), INTENT(IN ) :: Gomega
1567 : : ! INTEGER , INTENT(IN ) :: Wmax
1568 : : ! DOUBLE PRECISION, DIMENSION(:,:), INTENT(INOUT) :: F
1569 : : ! INTEGER :: flavors
1570 : : ! INTEGER :: samples
1571 : : ! INTEGER :: iflavor
1572 : : ! INTEGER :: iomega
1573 : : ! INTEGER :: itau
1574 : : ! DOUBLE PRECISION :: pi_invBeta
1575 : : ! DOUBLE PRECISION :: K
1576 : : ! COMPLEX*16 , DIMENSION(:,:), ALLOCATABLE :: F_omega
1577 : : ! TYPE(GreenHyboffdiag) :: F_tmp
1578 : : !
1579 : : ! flavors = op%flavors
1580 : : !
1581 : : ! samples = op%samples
1582 : : ! pi_invBeta = ACOS(-1.d0) / op%beta
1583 : : !
1584 : : ! IF ( op%have_MPI .EQV. .TRUE. ) THEN
1585 : : ! CALL GreenHyboffdiag_init(F_tmp,samples,op%beta,op%MY_COMM)
1586 : : ! ELSE
1587 : : ! CALL GreenHyboffdiag_init(F_tmp,samples,op%beta)
1588 : : ! END IF
1589 : : !! K = op%mu
1590 : : !
1591 : : ! MALLOC(F_omega,(1:Wmax,1:flavors))
1592 : : !
1593 : : ! DO iflavor = 1, flavors
1594 : : ! K = REAL(Gomega(Wmax, iflavor))
1595 : : ! WRITE(op%ostream,*) "CTQMC K, op%mu = ",K,op%mu
1596 : : ! WRITE(op%ostream,*) "CTQMC beta = ",op%beta
1597 : : ! op%mu(iflavor) = K
1598 : : ! DO iomega = 1, Wmax
1599 : : ! F_omega(iomega,iflavor) = Gomega(iomega,iflavor) &
1600 : : ! - CMPLX(K, 0.d0, 8)
1601 : : ! !if(iflavor==1.and.op%rank==0) then
1602 : : ! !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor))
1603 : : ! !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor))
1604 : : ! !end if
1605 : : ! END DO
1606 : : ! CALL GreenHyboffdiag_backFourier(F_tmp,F_omega(:,iflavor))
1607 : : ! F(1:samples+1,iflavor) = (/ (-F_tmp%oper(samples+1-itau),itau=0,samples) /)
1608 : : ! END DO
1609 : : ! FREE(F_omega)
1610 : : ! CALL GreenHyboffdiag_destroy(F_tmp)
1611 : : !END SUBROUTINE CtqmcoffdiagComplex_computeFK
1612 : : !!***
1613 : :
1614 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_run
1615 : : !! NAME
1616 : : !! CtqmcoffdiagComplex_run
1617 : : !!
1618 : : !! FUNCTION
1619 : : !! set all options and run a simulation
1620 : : !!
1621 : : !! COPYRIGHT
1622 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1623 : : !! This file is distributed under the terms of the
1624 : : !! GNU General Public License, see ~abinit/COPYING
1625 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1626 : : !!
1627 : : !! INPUTS
1628 : : !! op=ctqmc
1629 : : !! opt_order=maximal perturbation order to scope
1630 : : !! opt_movie=draw a movie of the simulation
1631 : : !! opt_analysis=compute correlation functions
1632 : : !! opt_check=check fast calculations
1633 : : !! opt_noise=compute noise for green function
1634 : : !! opt_spectra=fourier transform of the time evolution of the number of electrons
1635 : : !! opt_gMove=steps without global move
1636 : : !!
1637 : : !! OUTPUT
1638 : : !!
1639 : : !! SIDE EFFECTS
1640 : : !!
1641 : : !! NOTES
1642 : : !!
1643 : : !! SOURCE
1644 : :
1645 : 0 : SUBROUTINE CtqmcoffdiagComplex_run(op,opt_order,opt_histo,opt_movie,opt_analysis,opt_check,opt_noise,opt_spectra,opt_gMove)
1646 : :
1647 : :
1648 : : #ifdef HAVE_MPI1
1649 : : include 'mpif.h'
1650 : : #endif
1651 : : !Arguments ------------------------------------
1652 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
1653 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_order
1654 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_histo
1655 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_movie
1656 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_analysis
1657 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_check
1658 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_noise
1659 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_spectra
1660 : : INTEGER, OPTIONAL, INTENT(IN ) :: opt_gMove
1661 : : !Local variables ------------------------------
1662 : : #ifdef HAVE_MPI
1663 : : INTEGER :: ierr
1664 : : DOUBLE PRECISION :: rtime(1)
1665 : : #endif
1666 : : !#ifdef CTCtqmcoffdiagComplex_MOVIE
1667 : : INTEGER :: ilatex
1668 : : CHARACTER(LEN=4) :: Cchar
1669 : : !#endif
1670 : : DOUBLE PRECISION :: estimatedTime
1671 : : CHARACTER(LEN=100) :: message
1672 : :
1673 [ # # ]: 0 : IF ( .NOT. op%set ) &
1674 : 0 : CALL ERROR("CtqmcoffdiagComplex_run : QMC not set up ")
1675 [ # # ]: 0 : IF ( .NOT. op%setU ) &
1676 : 0 : CALL ERROR("CtqmcoffdiagComplex_run : QMC does not have a U matrix ")
1677 : :
1678 : :
1679 : : ! OPTIONS of the run
1680 [ # # ]: 0 : IF ( PRESENT( opt_check ) ) THEN
1681 : 0 : op%opt_check = opt_check
1682 : 0 : CALL ImpurityOperator_doCheck(op%Impurity,opt_check)
1683 : 0 : CALL BathOperatoroffdiagComplex_doCheck(op%Bath,opt_check)
1684 : : END IF
1685 [ # # ]: 0 : IF ( PRESENT( opt_movie ) ) &
1686 : 0 : op%opt_movie = opt_movie
1687 [ # # ]: 0 : IF ( PRESENT( opt_analysis ) ) &
1688 : 0 : op%opt_analysis = opt_analysis
1689 [ # # ]: 0 : IF ( PRESENT ( opt_order ) ) &
1690 : 0 : op%opt_order = opt_order
1691 [ # # ]: 0 : IF ( PRESENT ( opt_histo ) ) &
1692 : 0 : op%opt_histo = opt_histo
1693 [ # # ]: 0 : IF ( PRESENT ( opt_noise ) ) THEN
1694 : 0 : op%opt_noise = opt_noise
1695 : : END IF
1696 [ # # ]: 0 : IF ( PRESENT ( opt_spectra ) ) &
1697 : 0 : op%opt_spectra = opt_spectra
1698 : :
1699 : 0 : op%modGlobalMove(1) = max(op%sweeps,op%thermalization)+1 ! No Global Move
1700 : : !!sui!write(std_out,*) "op%sweeps",op%thermalization,op%sweeps,opt_gMove
1701 : 0 : op%modGlobalMove(2) = 0
1702 [ # # ]: 0 : IF ( PRESENT ( opt_gMove ) ) THEN
1703 [ # # # # ]: 0 : IF ( opt_gMove .LE. 0 .OR. opt_gMove .GT. op%sweeps ) THEN
1704 : : ! op%modGlobalMove(1) = op%sweeps+1
1705 : : op%modGlobalMove(1) = max(op%sweeps,op%thermalization)+1 ! No Global Move
1706 : : !write(std_out,*) "op%sweeps",op%sweeps, op%modGlobalMove(1)
1707 : : !CALL WARNALL("CtqmcoffdiagComplex_run : global moves option is <= 0 or > sweeps/cpu -> No global Moves")
1708 : 0 : write(message,'(a,a)') " == No global moves are used in CT-QMC. Make sure it is correct", ch10
1709 : 0 : CALL wrtout(std_out,message,'COLL')
1710 : : ELSE
1711 : 0 : op%modGlobalMove(1) = opt_gMove
1712 : : END IF
1713 : : END IF
1714 : : !sui!write(std_out,*) "op%sweeps",op%thermalization,op%sweeps
1715 : :
1716 : 0 : CALL CtqmcoffdiagComplex_allocateOpt(op)
1717 : :
1718 : : !#ifdef CTCtqmcoffdiagComplex_MOVIE
1719 : 0 : ilatex = 0
1720 [ # # ]: 0 : IF ( op%opt_movie .EQ. 1 ) THEN
1721 : 0 : Cchar ="0000"
1722 : 0 : WRITE(Cchar,'(I4)') op%rank
1723 : 0 : ilatex = 87+op%rank
1724 [ # # ]: 0 : OPEN(UNIT=ilatex, FILE="Movie_"//TRIM(ADJUSTL(Cchar))//".tex")
1725 : 0 : WRITE(ilatex,'(A)') "\documentclass{beamer}"
1726 : 0 : WRITE(ilatex,'(A)') "\usepackage{color}"
1727 : 0 : WRITE(ilatex,'(A)') "\setbeamersize{sidebar width left=0pt}"
1728 : 0 : WRITE(ilatex,'(A)') "\setbeamersize{sidebar width right=0pt}"
1729 : 0 : WRITE(ilatex,'(A)') "\setbeamersize{text width left=0pt}"
1730 : 0 : WRITE(ilatex,'(A)') "\setbeamersize{text width right=0pt}"
1731 : 0 : WRITE(ilatex,*)
1732 : 0 : WRITE(ilatex,'(A)') "\begin{document}"
1733 : 0 : WRITE(ilatex,*)
1734 : : END IF
1735 : : !#endif
1736 : :
1737 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
1738 : 0 : WRITE(op%ostream,'(A25)') "== Starting QMC Solver"
1739 : : END IF
1740 : :
1741 : : !=================================
1742 : : ! STARTING THERMALIZATION
1743 : : !=================================
1744 : : !write(std_out,*) "sweeps before thermalization",op%sweeps
1745 : : !write(std_out,*) "op%stats",op%stats
1746 : 0 : CALL CtqmcoffdiagComplex_loop(op,op%thermalization,ilatex)
1747 : : !=================================
1748 : : ! ENDING THERMALIZATION
1749 : : !=================================
1750 : :
1751 : 0 : estimatedTime = op%runTime
1752 : : #ifdef HAVE_MPI
1753 [ # # ]: 0 : CALL MPI_REDUCE([op%runTime], rtime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, 0, op%MY_COMM, ierr)
1754 : 0 : estimatedTime=rtime(1)
1755 : : #endif
1756 : :
1757 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
1758 [ # # ]: 0 : WRITE(op%ostream,'(A26,I6,A11)') "Thermalization done in ", CEILING(estimatedTime), " seconds"
1759 : 0 : WRITE(op%ostream,'(A25,I7,A15,I5,A5)') "The QMC should run in ", &
1760 [ # # ]: 0 : CEILING(estimatedTime*DBLE(op%sweeps)/DBLE(op%thermalization)),&
1761 : 0 : " seconds on ", op%size, " CPUs"
1762 : : END IF
1763 : :
1764 : : !=================================
1765 : : ! CLEANING CTQMC
1766 : : !=================================
1767 : 0 : CALL CtqmcoffdiagComplex_clear(op)
1768 : :
1769 : : !=================================
1770 : : ! STARTING CTQMC
1771 : : !=================================
1772 : : !write(std_out,*) "sweeps before loop",op%sweeps
1773 : : !write(std_out,*) "op%stats",op%stats
1774 : 0 : CALL CtqmcoffdiagComplex_loop(op,op%sweeps,ilatex)
1775 : : !=================================
1776 : : ! ENDING CTQMC
1777 : : !=================================
1778 : :
1779 [ # # ]: 0 : IF ( op%opt_movie .EQ. 1 ) THEN
1780 : 0 : WRITE(ilatex,*) ""
1781 : 0 : WRITE(ilatex,'(A14)') "\end{document}"
1782 : 0 : CLOSE(ilatex)
1783 : : END IF
1784 : :
1785 : 0 : op%done = .TRUE.
1786 : : !sui!write(std_out,*) "op%stats en of ctqmc_run",op%stats
1787 : :
1788 : 0 : END SUBROUTINE CtqmcoffdiagComplex_run
1789 : : !!***
1790 : :
1791 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_loop
1792 : : !! NAME
1793 : : !! CtqmcoffdiagComplex_loop
1794 : : !!
1795 : : !! FUNCTION
1796 : : !! Definition the main loop of the CT-QMC
1797 : : !!
1798 : : !! COPYRIGHT
1799 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
1800 : : !! This file is distributed under the terms of the
1801 : : !! GNU General Public License, see ~abinit/COPYING
1802 : : !! or http://www.gnu.org/copyleft/gpl.txt .
1803 : : !!
1804 : : !! INPUTS
1805 : : !! op=ctqmc
1806 : : !! itotal=number of sweeps to perform : thermalization or sweeps
1807 : : !! ilatex=unit of file to write movie if so
1808 : : !!
1809 : : !! OUTPUT
1810 : : !!
1811 : : !! SIDE EFFECTS
1812 : : !!
1813 : : !! NOTES
1814 : : !!
1815 : : !! SOURCE
1816 : :
1817 : 0 : SUBROUTINE CtqmcoffdiagComplex_loop(op,itotal,ilatex)
1818 : :
1819 : : !Arguments ------------------------------------
1820 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
1821 : : INTEGER , INTENT(IN ) :: itotal
1822 : : INTEGER , INTENT(IN ) :: ilatex
1823 : : !Local variables ------------------------------
1824 : : LOGICAL :: updated
1825 : : LOGICAL :: updated_seg
1826 : 0 : LOGICAL, DIMENSION(:), ALLOCATABLE :: updated_swap
1827 : :
1828 : : INTEGER :: flavors
1829 : : INTEGER :: measurements
1830 : : INTEGER :: modNoise1
1831 : : INTEGER :: modNoise2
1832 : : INTEGER :: modGlobalMove
1833 : : INTEGER :: sp1
1834 : : INTEGER :: itau
1835 : : INTEGER :: ind
1836 : : INTEGER :: endDensity
1837 : : INTEGER :: indDensity
1838 : : INTEGER :: swapUpdate1
1839 : : INTEGER :: swapUpdate2
1840 : : INTEGER :: old_percent
1841 : : INTEGER :: new_percent
1842 : : INTEGER :: ipercent !,ii
1843 : : INTEGER :: iflavor,ifl1,iflavor_d
1844 : : INTEGER :: isweep
1845 : :
1846 : : DOUBLE PRECISION :: cpu_time1
1847 : : DOUBLE PRECISION :: cpu_time2
1848 : : DOUBLE PRECISION :: NRJ_old1
1849 : : DOUBLE PRECISION :: NRJ_old2
1850 : : DOUBLE PRECISION :: NRJ_new
1851 : : DOUBLE PRECISION :: total
1852 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old1
1853 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old2
1854 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_new
1855 : :
1856 : 0 : CALL CPU_TIME(cpu_time1)
1857 : :
1858 : 0 : flavors = op%flavors
1859 : 0 : measurements = op%measurements
1860 : 0 : modNoise1 = op%modNoise1
1861 : 0 : modNoise2 = op%modNoise2
1862 : 0 : modGlobalMove = op%modGlobalMove(1)
1863 : 0 : sp1 = op%samples+1
1864 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
1865 [ # # ]: 0 : op%occup_histo_time= 0.d0
1866 [ # # ]: 0 : op%occupconfig= 0.d0
1867 [ # # # # ]: 0 : op%suscep= 0.d0
1868 [ # # # # ]: 0 : op%chi= 0.d0
1869 [ # # # # ]: 0 : op%chicharge= 0.d0
1870 [ # # ]: 0 : op%ntot= 0.d0
1871 : : END IF
1872 : :
1873 : 0 : old_percent = 0
1874 : :
1875 [ # # # # : 0 : MALLOC(updated_swap,(1:flavors))
# # # # #
# # # ]
1876 [ # # ]: 0 : updated_swap(:) = .FALSE.
1877 : :
1878 : 0 : NRJ_old1 = 0.d0
1879 : 0 : NRJ_old2 = 0.d0
1880 : 0 : NRJ_new = 0.d0
1881 : :
1882 [ # # ]: 0 : MALLOC(gtmp_new,(1,1))
1883 [ # # # # ]: 0 : gtmp_new = 0.d0
1884 [ # # ]: 0 : MALLOC(gtmp_old1,(1,1))
1885 [ # # # # ]: 0 : gtmp_old1 = 0.d0
1886 [ # # ]: 0 : MALLOC(gtmp_old2,(1,1))
1887 [ # # # # ]: 0 : gtmp_old2 = 0.d0
1888 : :
1889 : : !PROBLEM eos_gnu_13.2_mpich . %endDensity was introduced throughout
1890 : 0 : endDensity = SIZE(op%density,2)
1891 : : ! endDensity=op%endDensity
1892 : : !ENDPROBLEM
1893 : :
1894 : :
1895 [ # # ]: 0 : IF ( op%opt_noise .GT. 0 ) THEN
1896 [ # # ]: 0 : FREEIF(gtmp_new)
1897 [ # # # # : 0 : MALLOC(gtmp_new,(1:sp1,1:flavors))
# # # # #
# # # # #
# # # # ]
1898 [ # # ]: 0 : FREEIF(gtmp_old1)
1899 [ # # # # : 0 : MALLOC(gtmp_old1,(1:sp1,1:flavors))
# # # # #
# # # #
# ]
1900 [ # # ]: 0 : FREEIF(gtmp_old2)
1901 [ # # # # : 0 : MALLOC(gtmp_old2,(1:sp1,1:flavors))
# # # # #
# # # #
# ]
1902 : : END IF
1903 : :
1904 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
1905 : : WRITE(op%ostream, '(1x,103A)') &
1906 : 0 : "|----------------------------------------------------------------------------------------------------|"
1907 : 0 : WRITE(op%ostream,'(1x,A)', ADVANCE="NO") "|"
1908 : : END IF
1909 : :
1910 : 0 : total = DBLE(itotal)
1911 : : !write(std_out,*) "itotal",itotal
1912 : 0 : indDensity = 1
1913 : : !write(std_out,*) "op%stats",op%stats
1914 [ # # ]: 0 : DO isweep = 1, itotal
1915 : : !ii if(op%prtopt==1) write(std_out,*) "======== Isweep = ",isweep
1916 : : !updated_seg=.FALSE.
1917 [ # # ]: 0 : DO iflavor = 1, flavors
1918 : : ! if(isweep==itotal) write(std_out,*) " Iflavor = ",iflavor,op%Impurity%Particles(iflavor)%tail
1919 : : !ii if(op%prtopt==1) write(std_out,*) " ===Iflavor = ",iflavor
1920 : 0 : op%Impurity%activeFlavor=iflavor
1921 : 0 : op%Bath%activeFlavor=iflavor ; op%Bath%MAddFlag= .FALSE. ; op%Bath%MRemoveFlag = .FALSE.
1922 : :
1923 : : !write(std_out,*) "before tryaddremove"
1924 : :
1925 : : ! For iflavor, Try a move
1926 : : !==========================
1927 : 0 : CALL CtqmcoffdiagComplex_tryAddRemove(op,updated_seg)
1928 : : !sui!write(std_out,*) "after tryaddremove",updated_seg
1929 : :
1930 [ # # # # : 0 : updated = updated_seg .OR. updated_swap(iflavor).OR.(isweep==1)
# # ]
1931 : 0 : updated_swap(iflavor) = .FALSE.
1932 [ # # ]: 0 : if ( op%opt_nondiag >0 ) iflavor_d=0
1933 [ # # ]: 0 : if ( op%opt_nondiag==0 ) iflavor_d=iflavor
1934 : 0 : CALL GreenHyboffdiagComplex_measHybrid(op%Greens, op%Bath%M, op%Impurity%Particles, updated,op%signvalue,op%phasevalue,iflavor_d)
1935 : :
1936 : 0 : CALL CtqmcoffdiagComplex_measN (op, iflavor, updated)
1937 [ # # ]: 0 : IF ( op%opt_analysis .EQ. 1 ) &
1938 : 0 : CALL CtqmcoffdiagComplex_measCorrelation (op, iflavor)
1939 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) &
1940 : 0 : CALL CtqmcoffdiagComplex_measPerturbation(op, iflavor)
1941 : : END DO
1942 : : !CALL GreenHyboffdiag_measHybrid(op%Greens, op%Bath%M, op%Impurity%Particles, updated,op%signvalue,iflavor_d)
1943 : : !DO iflavor = 1,flavors
1944 : : ! CALL CtqmcoffdiagComplex_measN (op, iflavor, updated)
1945 : : !END DO
1946 : :
1947 [ # # ]: 0 : IF ( MOD(isweep,modGlobalMove) .EQ. 0 ) THEN
1948 : : ! !sui!write(std_out,*) "isweep,modGlobalMove,inside",isweep,modGlobalMove
1949 : 0 : CALL CtqmcoffdiagComplex_trySwap(op,swapUpdate1, swapUpdate2)
1950 : : ! !write(std_out,*) "no global move yet for non diag hybridization"
1951 [ # # # # ]: 0 : IF ( swapUpdate1 .NE. 0 .AND. swapUpdate2 .NE. 0 ) THEN
1952 : 0 : updated_swap(swapUpdate1) = .TRUE.
1953 : 0 : updated_swap(swapUpdate2) = .TRUE.
1954 : : END IF
1955 : : END IF
1956 : :
1957 [ # # ]: 0 : IF ( MOD(isweep,measurements) .EQ. 0 ) THEN ! default is always
1958 : 0 : CALL ImpurityOperator_measDE(op%Impurity,op%measDE)
1959 [ # # # # ]: 0 : IF ( op%opt_spectra .GE. 1 .AND. MOD(isweep,measurements*op%opt_spectra) .EQ. 0 ) THEN
1960 [ # # ]: 0 : op%density(1:flavors,indDensity) = real(op%measN(3,1:flavors))
1961 : 0 : indDensity = indDensity+1
1962 : : END IF
1963 : : END IF
1964 : :
1965 [ # # ]: 0 : IF ( MOD(isweep,measurements) .EQ. 0 ) THEN
1966 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
1967 : : CALL ImpurityOperator_occup_histo_time(op%Impurity,op%occup_histo_time,op%occupconfig,op%suscep,op%samples,op%chi,&
1968 : 0 : & op%chicharge,op%ntot,op%opt_histo,op%nspinor)
1969 : : END IF
1970 : : ENDIF
1971 : :
1972 [ # # ]: 0 : IF ( MOD(isweep, modNoise1) .EQ. 0 ) THEN
1973 : : !modNext = isweep + modNoise2
1974 : 0 : NRJ_new = op%measDE(1,1)
1975 : 0 : CALL Vector_pushBack(op%measNoise(1),NRJ_new - NRJ_old1)
1976 : 0 : NRJ_old1 = NRJ_new
1977 : :
1978 : : !! Try to limit accumulation error
1979 : 0 : CALL ImpurityOperator_cleanOverlaps(op%Impurity)
1980 : :
1981 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
1982 [ # # ]: 0 : DO ifl1 = 1, flavors
1983 [ # # ]: 0 : DO ind = 1, op%Greens%map(ifl1,ifl1)%tail
1984 : 0 : itau = op%Greens%map(ifl1,ifl1)%listINT(ind)
1985 : : gtmp_new(itau,ifl1) = op%Greens%oper(itau,ifl1,ifl1) &
1986 : 0 : +op%Greens%map(ifl1,ifl1)%listDBLE(ind)*DBLE(op%Greens%factor)
1987 : : END DO
1988 [ # # ]: 0 : DO itau = 1, sp1
1989 : 0 : CALL Vector_pushBack(op%measNoiseG(itau,ifl1,1), gtmp_new(itau,ifl1) - gtmp_old1(itau,ifl1))
1990 : 0 : gtmp_old1(itau,ifl1) = gtmp_new(itau,ifl1)
1991 : : END DO
1992 : : END DO
1993 : : END IF
1994 : : END IF
1995 : :
1996 [ # # ]: 0 : IF ( MOD(isweep,modNoise2) .EQ. 0 ) THEN
1997 : 0 : NRJ_new = op%measDE(1,1)
1998 : 0 : CALL Vector_pushBack(op%measNoise(2),NRJ_new - NRJ_old2)
1999 : 0 : NRJ_old2 = NRJ_new
2000 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
2001 [ # # ]: 0 : DO ifl1 = 1, flavors
2002 [ # # ]: 0 : DO ind = 1, op%Greens%map(ifl1,ifl1)%tail
2003 : 0 : itau = op%Greens%map(ifl1,ifl1)%listINT(ind)
2004 : : gtmp_new(itau,ifl1) = op%Greens%oper(itau,ifl1,ifl1) &
2005 : 0 : +op%Greens%map(ifl1,ifl1)%listDBLE(ind)*op%Greens%factor
2006 : : END DO
2007 [ # # ]: 0 : DO itau = 1, sp1
2008 : 0 : CALL Vector_pushBack(op%measNoiseG(itau,ifl1,2), gtmp_new(itau,ifl1) - gtmp_old2(itau,ifl1))
2009 : 0 : gtmp_old2(itau,ifl1) = gtmp_new(itau,ifl1)
2010 : : END DO
2011 : : END DO
2012 : : END IF
2013 : :
2014 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
2015 [ # # ]: 0 : new_percent = CEILING(DBLE(isweep)*100.d0/DBLE(itotal))
2016 [ # # ]: 0 : DO ipercent = old_percent+1, new_percent
2017 : 0 : WRITE(op%ostream,'(A)',ADVANCE="NO") "-"
2018 : : END DO
2019 : : old_percent = new_percent
2020 : : END IF
2021 : : END IF
2022 : :
2023 [ # # ]: 0 : IF ( op%opt_movie .EQ. 1 ) THEN
2024 : 0 : WRITE(ilatex,'(A11,I9)') "%iteration ", isweep
2025 : 0 : CALL ImpurityOperator_printLatex(op%Impurity,ilatex,isweep)
2026 : : END IF
2027 : :
2028 : : END DO
2029 : :
2030 [ # # ]: 0 : IF ( op%rank .EQ. 0 ) THEN
2031 [ # # ]: 0 : DO ipercent = old_percent+1, 100
2032 : 0 : WRITE(op%ostream,'(A)',ADVANCE="NO") "-"
2033 : : END DO
2034 : 0 : WRITE(op%ostream,'(A)') "|"
2035 : : END IF
2036 : :
2037 [ # # ]: 0 : FREE(gtmp_new)
2038 [ # # ]: 0 : FREE(gtmp_old1)
2039 [ # # ]: 0 : FREE(gtmp_old2)
2040 [ # # ]: 0 : FREE(updated_swap)
2041 : :
2042 [ # # # # ]: 0 : IF ( op%opt_spectra .GE. 1 .AND. itotal .EQ. op%sweeps ) THEN
2043 [ # # ]: 0 : IF ( endDensity .NE. indDensity-1 ) THEN
2044 [ # # ]: 0 : op%density(:,endDensity) = -1.d0
2045 : : END IF
2046 : : END IF
2047 : :
2048 : 0 : CALL CPU_TIME(cpu_time2)
2049 : :
2050 : 0 : op%runTime = (cpu_time2 - cpu_time1)*1.05d0 ! facteur arbitraire de correction
2051 : 0 : END SUBROUTINE CtqmcoffdiagComplex_loop
2052 : : !!***
2053 : :
2054 : :
2055 : :
2056 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_tryAddRemove
2057 : : !! NAME
2058 : : !! CtqmcoffdiagComplex_tryAddRemove
2059 : : !!
2060 : : !! FUNCTION
2061 : : !! Try to add or remove a segment and an anti-segment
2062 : : !!
2063 : : !! COPYRIGHT
2064 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
2065 : : !! This file is distributed under the terms of the
2066 : : !! GNU General Public License, see ~abinit/COPYING
2067 : : !! or http://www.gnu.org/copyleft/gpl.txt .
2068 : : !!
2069 : : !! INPUTS
2070 : : !! op=ctqmc
2071 : : !!
2072 : : !! OUTPUT
2073 : : !! updated=something changed
2074 : : !!
2075 : : !! SIDE EFFECTS
2076 : : !!
2077 : : !! NOTES
2078 : : !!
2079 : : !! SOURCE
2080 : :
2081 : 0 : SUBROUTINE CtqmcoffdiagComplex_tryAddRemove(op,updated)
2082 : :
2083 : : !Arguments ------------------------------------
2084 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
2085 : : ! TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: Bath
2086 : : ! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity
2087 : : LOGICAL , INTENT( OUT) :: updated
2088 : : !Local variables ------------------------------
2089 : : INTEGER :: position
2090 : : INTEGER , DIMENSION(1:2) :: nature ! -2 for antiseg and 1 for seg
2091 : : INTEGER :: i! -2 for antiseg and 1 for seg
2092 : : !INTEGER :: it,it1 !ii,
2093 : : DOUBLE PRECISION :: action
2094 : : DOUBLE PRECISION :: beta
2095 : : DOUBLE PRECISION :: time1
2096 : : DOUBLE PRECISION :: time2
2097 : : DOUBLE PRECISION :: time_avail
2098 : : COMPLEX(KIND=8) :: det_ratio
2099 : : DOUBLE PRECISION :: sign_det_ratio
2100 : : DOUBLE PRECISION :: overlap
2101 : : DOUBLE PRECISION :: length
2102 : : DOUBLE PRECISION :: signe
2103 : : DOUBLE PRECISION :: tail
2104 : : INTEGER :: tailint
2105 : : DOUBLE PRECISION :: signdet, signdetprev
2106 : : DOUBLE PRECISION, DIMENSION(1:2) :: CdagC_1
2107 : :
2108 [ # # ]: 0 : IF ( .NOT. op%set ) &
2109 : 0 : CALL ERROR("CtqmcoffdiagComplex_trySegment : QMC not set ")
2110 : :
2111 : : !write(std_out,*) " TryAddRemove start"
2112 : 0 : nature(1) = CTQMC_SEGME
2113 : 0 : nature(2) = CTQMC_ANTIS
2114 : 0 : beta = op%beta
2115 : :
2116 : 0 : updated = .FALSE.
2117 : : !tail is the number of segments for a given flavor
2118 : 0 : tailint = (op%Impurity%particles(op%Impurity%activeFlavor)%tail)
2119 : 0 : tail = DBLE(tailint)
2120 : :
2121 : :
2122 : : !=====================================
2123 : : ! First choose segment or antisegment
2124 : : !=====================================
2125 [ # # ]: 0 : DO i = 1, 2
2126 : 0 : signe = SIGN(1.d0,DBLE(nature(i)))
2127 : : ! ----- 1: segment signe= 1 ( CTQMC_SEGME = 1 )
2128 : : ! ----- 2: antisegment signe=-1 ( CTQMC_ANTIS = -2 )
2129 : :
2130 : 0 : tailint = (op%Impurity%particles(op%Impurity%activeFlavor)%tail)
2131 : :
2132 : : !CALL RANDOM_NUMBER(action)
2133 : 0 : CALL OurRng(op%seed,action)
2134 : :
2135 : : !==========================
2136 : : ! Add segment/antisegment
2137 : : !==========================
2138 [ # # ]: 0 : IF ( action .LT. .5d0 ) THEN ! Add a segment or antisegment
2139 : :
2140 : : ! Select time1 (>0) in [0,beta]
2141 : : !==============================
2142 : : !CALL RANDOM_NUMBER(time1)
2143 : 0 : CALL OurRng(op%seed,time1)
2144 : 0 : time1 = time1 * beta
2145 : :
2146 : : ! time_avail is the distance between between time1 and
2147 : : ! - the next start of a segment for a segment addition
2148 : : ! - the next end of a segment for an antisegment addition
2149 : : ! ImpurityOperator_getAvailableTime > 0 for a segment (signe>0) -> time_avail>0
2150 : : ! ImpurityOperator_getAvailableTime < 0 for an antisegment (signe<0) -> time_avail>0
2151 : : !====================================================================
2152 : 0 : time_avail = ImpurityOperator_getAvailableTime(op%Impurity,time1,position) * signe
2153 : :
2154 [ # # ]: 0 : IF ( time_avail .GT. 0.d0 ) THEN
2155 : :
2156 : : ! Time2 is the length of the proposed new (anti)segment
2157 : : !=======================================================
2158 : 0 : CALL OurRng(op%seed,time2)
2159 [ # # ]: 0 : IF ( time2 .EQ. 0.d0 ) CALL OurRng(op%seed,time2) ! Prevent null segment
2160 : :
2161 : : ! Now time2 is the time at the end of the proposed new (anti) segment
2162 : : ! time2 > time1
2163 : : !====================================================================
2164 : 0 : time2 = time1 + time2 * time_avail
2165 : :
2166 : : ! CdagC_1 gives the stard/end times for the proposed new segment/antisegment
2167 : : ! CdagC1(C_) can be above beta.
2168 : : ! For a segment CdagC_1(Cdag_) = time1 < CdagC_1(C_) = time2, l=time2-time1 > 0
2169 : : ! For a anti segment CdagC_1(Cdag_) = time2 > CdagC_1(C_) = time1, l=time1-time2 < 0
2170 : : ! time2 can be above beta and thus for a segment CdagC_1(C_ ) > beta
2171 : : ! time2 can be above beta and thus for an antisegment CdagC_1(Cdag_) > beta
2172 : : ! length > 0 for segment
2173 : : ! length < 0 for antisegment
2174 : : !====================================================================================
2175 : 0 : CdagC_1(Cdag_) = ((1.d0+signe)*time1+(1.d0-signe)*time2)*0.5d0
2176 : 0 : CdagC_1(C_ ) = ((1.d0+signe)*time2+(1.d0-signe)*time1)*0.5d0
2177 : 0 : length = CdagC_1(C_ ) - CdagC_1(Cdag_)
2178 : :
2179 : : ! ----- Computes the determinant ratio
2180 : 0 : det_ratio = BathOperatoroffdiagComplex_getDetAdd(op%Bath,CdagC_1,position,op%Impurity%particles)
2181 : :
2182 : : ! ----- Computes the overlap
2183 : 0 : overlap = ImpurityOperator_getNewOverlap(op%Impurity,CdagC_1)
2184 : 0 : signdetprev = ImpurityOperator_getsign(op%Impurity, time2, i, action, position)
2185 : : !write(std_out,*)"det_ratio",det_ratio,"signdetprev",signdetprev
2186 : :
2187 : : !Remove negative sign of det_ration
2188 : 0 : det_ratio=det_ratio*signdetprev
2189 : :
2190 [ # # ]: 0 : IF ( real(det_ratio) .LT. 0.d0 ) THEN
2191 : 0 : det_ratio = - det_ratio
2192 : 0 : sign_det_ratio=-1
2193 : 0 : op%stats(nature(i)+CTQMC_DETSI) = op%stats(nature(i)+CTQMC_DETSI) + 1.d0
2194 : : ELSE
2195 : : sign_det_ratio=1
2196 : : END IF
2197 : :
2198 : 0 : CALL OurRng(op%seed,time1)
2199 : : !write(std_out,*) " .................",(time1 * (tail + 1.d0 )),beta * time_avail * abs(det_ratio) * DEXP(real(op%mu(op%Impurity%activeFlavor))*length + overlap)
2200 : : !write(std_out,*) " .................",beta , time_avail , op%mu(op%Impurity%activeFlavor),op%Impurity%activeFlavor
2201 : :
2202 [ # # ]: 0 : IF ( (time1 * (tail + 1.d0 )) &
2203 : : .LT. (beta * time_avail * abs(det_ratio) * DEXP(real(op%mu(op%Impurity%activeFlavor))*length + overlap) ) ) THEN
2204 : 0 : CALL ImpurityOperator_add(op%Impurity,CdagC_1,position)
2205 : : ! write(*,*) "after "
2206 : : ! CALL ListCdagC_print(op%Impurity%particles(op%Impurity%activeFlavor),6)
2207 : 0 : CALL BathOperatoroffdiagComplex_setMAdd(op%bath,op%Impurity%particles)
2208 : 0 : op%stats(nature(i)+CTQMC_ADDED) = op%stats(nature(i)+CTQMC_ADDED) + 1.d0
2209 : 0 : updated = .TRUE. .OR. updated
2210 : 0 : tail = tail + 1.d0
2211 : 0 : tailint = tailint + 1
2212 : :
2213 [ # # ]: 0 : IF ( sign_det_ratio .LT. 0.d0 ) op%signvalue=-op%signvalue
2214 : 0 : op%phasevalue=op%phasevalue*det_ratio/abs(det_ratio)*sign_det_ratio
2215 : : END IF
2216 : :
2217 : : END IF ! time_avail > 0
2218 : :
2219 : : !========================================
2220 : : ! Remove segment/antisegment
2221 : : !========================================
2222 : : ELSE ! Remove a segment among the segment of the flavor activeflavor
2223 : : !ii if(op%prtopt==1) write(6,*) " =try: Segment removed of type",i
2224 [ # # ]: 0 : IF ( tail .GT. 0.d0 ) THEN
2225 : 0 : CALL OurRng(op%seed,time1)
2226 : 0 : position = INT(((time1 * tail) + 1.d0) * signe )
2227 : : !prt!if(op%prtopt==1) write(6,*) " position",position
2228 : 0 : time_avail = ImpurityOperator_getAvailedTime(op%Impurity,position)
2229 : 0 : det_ratio = BathOperatoroffdiagComplex_getDetRemove(op%Bath,position)
2230 : : !write(6,*) " det_ratio", det_ratio
2231 : 0 : CdagC_1 = ImpurityOperator_getSegment(op%Impurity,position)
2232 : : ! length = CdagC_length(CdagC_1)
2233 : 0 : length = CdagC_1(C_) - CdagC_1(Cdag_)
2234 : : !write(6,*) " length ", length
2235 : 0 : overlap = ImpurityOperator_getNewOverlap(op%Impurity,CdagC_1)
2236 : : !write(6,*) " overlap ", overlap
2237 : 0 : CALL OurRng(op%seed,time1)
2238 : : !write(6,*) " Random ",time1
2239 : 0 : signdetprev = ImpurityOperator_getsign(op%Impurity, time2, i, action, position)
2240 : 0 : det_ratio=det_ratio*signdetprev
2241 : 0 : signdet=1.d0
2242 : :
2243 [ # # ]: 0 : IF ( real(det_ratio) .LT. 0.d0 ) THEN
2244 : : !sui!write(6,*) " NEGATIVE DET",det_ratio,signdetprev
2245 : 0 : det_ratio = -det_ratio
2246 : 0 : sign_det_ratio=-1
2247 : 0 : op%stats(nature(i)+CTQMC_DETSI) = op%stats(nature(i)+CTQMC_DETSI) + 1.d0
2248 : 0 : signdet=-1.d0
2249 : : ELSE
2250 : : sign_det_ratio=1
2251 : : END IF
2252 : :
2253 [ # # ]: 0 : IF ( (time1 * beta * time_avail * DEXP(real(op%mu(op%Impurity%activeFlavor))*length+overlap)) &
2254 : : .LT. (tail * abs(det_ratio) ) ) THEN
2255 : 0 : CALL ImpurityOperator_remove(op%Impurity,position)
2256 : 0 : CALL BathOperatoroffdiagComplex_setMRemove(op%Bath,op%Impurity%particles)
2257 : : !op%seg_removed = op%seg_removed + 1.d0
2258 : 0 : op%stats(nature(i)+CTQMC_REMOV) = op%stats(nature(i)+CTQMC_REMOV) + 1.d0
2259 : 0 : updated = .TRUE. .OR. updated
2260 : 0 : tail = tail -1.d0
2261 : 0 : tailint = tailint -1
2262 : :
2263 [ # # ]: 0 : IF ( sign_det_ratio .LT. 0.d0 ) op%signvalue=-op%signvalue
2264 : 0 : op%phasevalue=op%phasevalue*det_ratio/abs(det_ratio)*sign_det_ratio
2265 : : END IF
2266 : : END IF !tail > O
2267 : : END IF ! action
2268 : : !========================================
2269 : : ! End Add/Remove Antisegment
2270 : : !========================================
2271 : : END DO
2272 : 0 : END SUBROUTINE CtqmcoffdiagComplex_tryAddRemove
2273 : : !!***
2274 : :
2275 : : !!SUBROUTINE CtqmcoffdiagComplex_trySegment(op,updated)
2276 : : !!! TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: Bath
2277 : : !!! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity
2278 : : !! LOGICAL , INTENT(INOUT) :: updated
2279 : : !! INTEGER :: position
2280 : : !! DOUBLE PRECISION :: action
2281 : : !! DOUBLE PRECISION :: beta
2282 : : !! DOUBLE PRECISION :: time1
2283 : : !! DOUBLE PRECISION :: time2
2284 : : !! DOUBLE PRECISION :: time_avail
2285 : : !! DOUBLE PRECISION :: det_ratio
2286 : : !! DOUBLE PRECISION :: overlap
2287 : : !! DOUBLE PRECISION :: length
2288 : : !! DOUBLE PRECISION :: tail
2289 : : !! DOUBLE PRECISION, DIMENSION(1:2) :: CdagC_1
2290 : : !!
2291 : : !! IF ( .NOT. op%set ) &
2292 : : !! CALL ERROR("CtqmcoffdiagComplex_trySegment : QMC not set ")
2293 : : !!
2294 : : !! beta = op%beta
2295 : : !! tail = DBLE(op%Impurity%particles(op%Impurity%activeFlavor)%tail)
2296 : : !!
2297 : : !! CALL RANDOM_NUMBER(action)
2298 : : !!
2299 : : !! updated = .FALSE.
2300 : : !!
2301 : : !! IF ( action .LT. .5d0 ) THEN ! Ajout de segment
2302 : : !! CALL RANDOM_NUMBER(time1)
2303 : : !! time1 = time1 * beta
2304 : : !! time_avail = ImpurityOperator_getAvailableTime(op%Impurity,time1,position)
2305 : : !! IF ( time_avail .GT. 0.d0 ) THEN
2306 : : !! CALL RANDOM_NUMBER(time2)
2307 : : !! time2 = time1 + time2 * time_avail
2308 : : !!! CALL CdagC_init(CdagC_1,time1,time2)
2309 : : !! CdagC_1(Cdag_) = time1
2310 : : !! CdagC_1(C_ ) = time2
2311 : : !!! length = CdagC_length(CdagC_1)
2312 : : !! length = time2 - time1
2313 : : !! det_ratio = BathOperatoroffdiag_getDetAdd(op%Bath,CdagC_1,position,op%Impurity%particles(op%Impurity%activeFlavor))
2314 : : !! overlap = ImpurityOperator_getNewOverlap(op%Impurity,CdagC_1)
2315 : : !! CALL RANDOM_NUMBER(time1)
2316 : : !! IF ( det_ratio .LT. 0.d0 ) THEN
2317 : : !! det_ratio = -det_ratio
2318 : : !! op%stats(CTQMC_SEGME+CTQMC_DETSI) = op%stats(CTQMC_SEGME+CTQMC_DETSI) + SIGN(1.d0,det_ratio)
2319 : : !! END IF
2320 : : !! IF ( (time1 * (tail + 1.d0 )) &
2321 : : !! .LT. (beta * time_avail * det_ratio * DEXP(op%mu(op%Impurity%activeFlavor)*length + overlap) ) ) THEN
2322 : : !! write(*,*) position
2323 : : !! CALL ImpurityOperator_add(op%Impurity,CdagC_1,position)
2324 : : !! CALL BathOperatoroffdiag_setMAdd(op%bath,op%Impurity%particles(op%Impurity%activeFlavor))
2325 : : !! op%stats(CTQMC_SEGME+CTQMC_ADDED) = op%stats(CTQMC_SEGME+CTQMC_ADDED) + 1.d0
2326 : : !! updated = .TRUE.
2327 : : !! END IF
2328 : : !! END IF
2329 : : !!
2330 : : !! ELSE ! Supprimer un segment
2331 : : !! IF ( tail .GT. 0.d0 ) THEN
2332 : : !! CALL RANDOM_NUMBER(time1)
2333 : : !! position = INT(time1 * tail) + 1
2334 : : !! time_avail = ImpurityOperator_getAvailedTime(op%Impurity,position)
2335 : : !! det_ratio = BathOperatoroffdiag_getDetRemove(op%Bath,position)
2336 : : !! CdagC_1 = ImpurityOperator_getSegment(op%Impurity,position)
2337 : : !!! length = CdagC_length(CdagC_1)
2338 : : !! length = CdagC_1(C_) - CdagC_1(Cdag_)
2339 : : !! overlap = ImpurityOperator_getNewOverlap(op%Impurity,CdagC_1)
2340 : : !! CALL RANDOM_NUMBER(time1)
2341 : : !! IF ( det_ratio .LT. 0.d0 ) THEN
2342 : : !! det_ratio = -det_ratio
2343 : : !!! op%seg_sign = op%seg_sign + 1.d0
2344 : : !! op%stats(CTQMC_SEGME+CTQMC_DETSI) = op%stats(CTQMC_SEGME+CTQMC_DETSI) + SIGN(1.d0,det_ratio)
2345 : : !! END IF
2346 : : !! IF ( (time1 * beta * time_avail * DEXP(op%mu(op%Impurity%activeFlavor)*length+overlap)) &
2347 : : !! .LT. (tail * det_ratio ) ) THEN
2348 : : !! write(*,*) position
2349 : : !! CALL ImpurityOperator_remove(op%Impurity,position)
2350 : : !! CALL BathOperatoroffdiag_setMRemove(op%Bath,op%Impurity%particles(op%Impurity%activeFlavor))
2351 : : !! !op%seg_removed = op%seg_removed + 1.d0
2352 : : !! op%stats(CTQMC_SEGME+CTQMC_REMOV) = op%stats(CTQMC_SEGME+CTQMC_REMOV) + 1.d0
2353 : : !! updated = .TRUE.
2354 : : !! END IF
2355 : : !! END IF
2356 : : !! END IF
2357 : : !!END SUBROUTINE CtqmcoffdiagComplex_trySegment
2358 : : !!***
2359 : : !!
2360 : : !!SUBROUTINE CtqmcoffdiagComplex_tryAntiSeg(op, updated)
2361 : : !!! TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: Bath
2362 : : !!! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity
2363 : : !! LOGICAL , INTENT(INOUT) :: updated
2364 : : !! INTEGER :: position
2365 : : !! DOUBLE PRECISION :: action
2366 : : !! DOUBLE PRECISION :: beta
2367 : : !! DOUBLE PRECISION :: time1
2368 : : !! DOUBLE PRECISION :: time2
2369 : : !! DOUBLE PRECISION :: time_avail
2370 : : !! DOUBLE PRECISION :: det_ratio
2371 : : !! DOUBLE PRECISION :: overlap
2372 : : !! DOUBLE PRECISION :: length
2373 : : !! DOUBLE PRECISION :: tail
2374 : : !! DOUBLE PRECISION, DIMENSION(1:2) :: CdagC_1
2375 : : !!
2376 : : !! IF ( .NOT. op%set ) &
2377 : : !! CALL ERROR("CtqmcoffdiagComplex_trySegment : QMC not set ")
2378 : : !!
2379 : : !! beta = op%beta
2380 : : !! tail = DBLE(op%Impurity%particles(op%Impurity%activeFlavor)%tail)
2381 : : !!
2382 : : !! CALL RANDOM_NUMBER(action)
2383 : : !!
2384 : : !! updated = .FALSE.
2385 : : !!
2386 : : !! IF ( action .LT. .5d0 ) THEN ! Ajout d'un antiseg
2387 : : !! CALL RANDOM_NUMBER(time1)
2388 : : !! time1 = time1 * beta
2389 : : !! time_avail = ImpurityOperator_getAvailableTime(op%Impurity,time1,position)
2390 : : !! IF ( time_avail .LT. 0.d0 ) THEN
2391 : : !! CALL RANDOM_NUMBER(time2)
2392 : : !! time2 = time1 - time2 * time_avail
2393 : : !!! CALL CdagC_init(CdagC_1,time2,time1)
2394 : : !! CdagC_1(Cdag_) = time2
2395 : : !! CdagC_1(C_ ) = time1
2396 : : !!! length = CdagC_length(CdagC_1) ! /!\ length is negative
2397 : : !! length = time1 - time2
2398 : : !! det_ratio = BathOperatoroffdiag_getDetAdd(op%Bath,CdagC_1,position,op%Impurity%particles(op%Impurity%activeFlavor))
2399 : : !! overlap = ImpurityOperator_getNewOverlap(op%Impurity,CdagC_1) ! OK
2400 : : !! CALL RANDOM_NUMBER(time1)
2401 : : !! IF ( det_ratio .LT. 0.d0 ) THEN
2402 : : !! det_ratio = -det_ratio
2403 : : !!! op%anti_sign = op%anti_sign + 1.d0
2404 : : !! op%stats(CTQMC_ANTIS+CTQMC_DETSI) = op%stats(CTQMC_ANTIS+CTQMC_DETSI) + SIGN(1.d0,det_ratio)
2405 : : !! END IF
2406 : : !! IF ( (time1 * (tail + 1.d0 )) &
2407 : : !! .LT. (beta * ABS(time_avail) * det_ratio * DEXP(op%mu(op%Impurity%activeFlavor)*length + overlap) ) ) THEN
2408 : : !! CALL ImpurityOperator_add(op%Impurity,CdagC_1,position)
2409 : : !! !write(*,*) position
2410 : : !! CALL BathOperatoroffdiag_setMAdd(op%bath,op%Impurity%particles(op%Impurity%activeFlavor))
2411 : : !! !op%anti_added = op%anti_added + 1.d0
2412 : : !! op%stats(CTQMC_ANTIS+CTQMC_ADDED) = op%stats(CTQMC_ANTIS+CTQMC_ADDED) + 1.d0
2413 : : !! updated = .TRUE.
2414 : : !! END IF
2415 : : !! END IF
2416 : : !! ELSE ! Supprimer un antiseg
2417 : : !! IF ( tail .GT. 0.d0 ) THEN
2418 : : !! CALL RANDOM_NUMBER(time1)
2419 : : !! position = -(INT(time1 * tail) + 1)
2420 : : !! time_avail = ImpurityOperator_getAvailedTime(op%Impurity,position)!OK
2421 : : !! det_ratio = BathOperatoroffdiag_getDetRemove(op%Bath,position)!OK
2422 : : !! CdagC_1 = ImpurityOperator_getSegment(op%Impurity,position)!OK
2423 : : !!! length = CdagC_length(CdagC_1) ! /!\ length is negative
2424 : : !! length = CdagC_1(C_) - CdagC_1(Cdag_)
2425 : : !! overlap = ImpurityOperator_getNewOverlap(op%Impurity,CdagC_1) !OK
2426 : : !! CALL RANDOM_NUMBER(time1)
2427 : : !! IF ( det_ratio .LT. 0.d0 ) THEN
2428 : : !! det_ratio = -det_ratio
2429 : : !!! op%anti_sign = op%anti_sign + 1.d0
2430 : : !! op%stats(CTQMC_ANTIS+CTQMC_DETSI) = op%stats(CTQMC_ANTIS+CTQMC_DETSI) + SIGN(1.d0,det_ratio)
2431 : : !! END IF
2432 : : !! IF ( (time1 * beta * time_avail * DEXP(op%mu(op%Impurity%activeFlavor)*length+overlap)) &
2433 : : !! .LT. (tail * det_ratio ) ) THEN
2434 : : !! CALL ImpurityOperator_remove(op%Impurity,position)
2435 : : !! !write(*,*) position
2436 : : !! CALL BathOperatoroffdiag_setMRemove(op%Bath,op%Impurity%particles(op%Impurity%activeFlavor))
2437 : : !! !op%anti_removed = op%anti_removed + 1.d0
2438 : : !! op%stats(CTQMC_ANTIS+CTQMC_REMOV) = op%stats(CTQMC_ANTIS+CTQMC_REMOV) + 1.d0
2439 : : !! updated = .TRUE.
2440 : : !! END IF
2441 : : !! END IF
2442 : : !! END IF
2443 : : !!END SUBROUTINE CtqmcoffdiagComplex_tryAntiSeg
2444 : : !!***
2445 : :
2446 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_trySwap
2447 : : !! NAME
2448 : : !! CtqmcoffdiagComplex_trySwap
2449 : : !!
2450 : : !! FUNCTION
2451 : : !! try a global move (swap to flavors)
2452 : : !!
2453 : : !! COPYRIGHT
2454 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
2455 : : !! This file is distributed under the terms of the
2456 : : !! GNU General Public License, see ~abinit/COPYING
2457 : : !! or http://www.gnu.org/copyleft/gpl.txt .
2458 : : !!
2459 : : !! INPUTS
2460 : : !! op=ctqmc
2461 : : !!
2462 : : !! OUTPUT
2463 : : !! flav_i=first flavor swaped
2464 : : !! flav_j=second flavor swaped
2465 : : !!
2466 : : !! SIDE EFFECTS
2467 : : !!
2468 : : !! NOTES
2469 : : !!
2470 : : !! SOURCE
2471 : :
2472 : 0 : SUBROUTINE CtqmcoffdiagComplex_trySwap(op,flav_i,flav_j)
2473 : :
2474 : : !Arguments ------------------------------------
2475 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
2476 : : ! TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: Bath
2477 : : ! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity
2478 : : INTEGER , INTENT( OUT) :: flav_i
2479 : : INTEGER , INTENT( OUT) :: flav_j
2480 : : !Local variables ------------------------------
2481 : : INTEGER :: flavor_i
2482 : : INTEGER :: flavor_j !,ii,it,it1 !,iflavor
2483 : : DOUBLE PRECISION :: rnd
2484 : : DOUBLE PRECISION :: lengthi
2485 : : DOUBLE PRECISION :: lengthj
2486 : : DOUBLE PRECISION :: overlapic1
2487 : : DOUBLE PRECISION :: overlapjc1
2488 : : DOUBLE PRECISION :: overlapic2
2489 : : DOUBLE PRECISION :: overlapjc2
2490 : : !DOUBLE PRECISION :: detic1
2491 : : !DOUBLE PRECISION :: detjc1
2492 : : !DOUBLE PRECISION :: detic2
2493 : : !DOUBLE PRECISION :: detjc2
2494 : : DOUBLE PRECISION :: det_ratio,detnew,detold
2495 : : DOUBLE PRECISION :: local_ratio
2496 : : ! TYPE(BathOperatoroffdiag) :: Bathnew
2497 : :
2498 : :
2499 : : !CALL RANDOM_NUMBER(rnd)
2500 : 0 : CALL OurRng(op%seed,rnd)
2501 : 0 : flavor_i = NINT(rnd*DBLE(op%flavors-1.d0))+1
2502 : : !CALL RANDOM_NUMBER(rnd)
2503 : 0 : CALL OurRng(op%seed,rnd)
2504 : 0 : flavor_j = NINT(rnd*DBLE(op%flavors-1.d0))+1
2505 : : !ii write(6,'(a,2i4)') "--------------- new swap --------------------------------",flavor_i,flavor_j
2506 : :
2507 : 0 : flav_i = 0
2508 : 0 : flav_j = 0
2509 : : !ii do iflavor=1,op%flavors
2510 : : !ii write(6,*) "BEFORE GMOVE For flavor", iflavor,"size is",op%Impurity%particles(iflavor)%tail," and Conf is :"
2511 : : !ii do ii=1, op%Impurity%Particles(iflavor)%tail
2512 : : !ii write(6,'(i4,100f12.3)') ii, op%Impurity%Particles(iflavor)%list(ii,1),&
2513 : : !ii & op%Impurity%Particles(iflavor)%list(ii,2)
2514 : : !ii enddo
2515 : : !ii enddo
2516 : : !ii write(6,*) " = M Matrix"
2517 : : !ii write(6,'(a,2x,100(i12))') "Flavor=",((iflavor,it=1,op%Impurity%particles(iflavor)%tail),iflavor=1,op%flavors)
2518 : : !ii write(6,'(i21,100i12)') ((it,it=1,op%Impurity%particles(iflavor)%tail),iflavor=1,op%flavors)
2519 : : !ii do it=1,op%Bath%sumtails
2520 : : !ii write(6,'(a,100f12.3)') " M before",(op%Bath%M%mat(it,it1),it1=1,op%Bath%sumtails)
2521 : : !ii enddo
2522 : :
2523 : : ! todoba this part
2524 [ # # ]: 0 : IF ( flavor_i .NE. flavor_j ) THEN
2525 : : !CALL BathOperatoroffdiag_init(Bathnew, op%flavors, op%samples, op%beta, 0, op%opt_nondiag)
2526 : : ! On tente d'intervertir i et j
2527 : : ! Configuration actuelle :
2528 : :
2529 : 0 : op%modGlobalMove(2) = op%modGlobalMove(2)+1
2530 : : ! ===========================================
2531 : : ! First use M matrix to compute determinant
2532 : : ! ===========================================
2533 : 0 : detold = BathOperatoroffdiagComplex_getDetF(op%Bath) ! use op%Bath%M
2534 : :
2535 : : ! ===========================================
2536 : : ! Second build M_update matrix to compute determinant after.
2537 : : ! ===========================================
2538 : : !CALL ListCdagCoffdiag_print(particle)
2539 : 0 : call BathOperatoroffdiagComplex_recomputeM(op%Bath,op%impurity%particles,flavor_i,flavor_j) ! compute op%Bath%M_update
2540 : 0 : detnew = BathOperatoroffdiagComplex_getDetF(op%Bath,option=1) ! use op%Bath%M_update
2541 : :
2542 : 0 : lengthi = real(ImpurityOperator_measN(op%Impurity,flavor_i))
2543 : 0 : lengthj = real(ImpurityOperator_measN(op%Impurity,flavor_j))
2544 : 0 : overlapic1 = ImpurityOperator_overlapFlavor(op%Impurity,flavor_i)
2545 : 0 : overlapjc1 = ImpurityOperator_overlapFlavor(op%Impurity,flavor_j)
2546 : : ! lengths unchanged
2547 : 0 : overlapic2 = ImpurityOperator_overlapSwap(op%Impurity,flavor_i,flavor_j)
2548 : 0 : overlapjc2 = ImpurityOperator_overlapSwap(op%Impurity,flavor_j,flavor_i)
2549 : :
2550 : : ! IF ( detic1*detjc1 .EQ. detic2*detjc2 ) THEN
2551 : : ! det_ratio = 1.d0
2552 : : ! ELSE IF ( detic1*detjc1 .EQ. 0.d0 ) THEN
2553 : : ! det_ratio = detic2*detjc2 ! evite de diviser par 0 si pas de segment
2554 : : ! ELSE
2555 : :
2556 : 0 : det_ratio = detnew/detold ! because the determinant is the determinant of F
2557 : : !ii write(6,*) "det_ratio, detold,detnew",det_ratio, detold,detnew, detold/detnew
2558 : :
2559 : : ! END IF
2560 : : local_ratio = DEXP(-overlapic2*overlapjc2+overlapic1*overlapjc1 &
2561 : 0 : +(lengthj-lengthi)*real((op%mu(flavor_i))-real(op%mu(flavor_j))))
2562 : : !ii write(6,*) "local_ratio",local_ratio
2563 : :
2564 : : ! Wloc = exp(muN-Uo)
2565 : : !CALL RANDOM_NUMBER(rnd)
2566 : 0 : CALL OurRng(op%seed,rnd)
2567 [ # # ]: 0 : IF ( rnd .LT. local_ratio*det_ratio ) THEN ! swap accepted
2568 : : !ii write(6,*) " = M Matrix before swap"
2569 : : !ii write(6,'(a,2x,100(i12))') "Flavor=",((iflavor,it=1,op%Impurity%particles(iflavor)%tail),iflavor=1,op%flavors)
2570 : : !ii write(6,'(i21,100i12)') ((it,it=1,op%Impurity%particles(iflavor)%tail),iflavor=1,op%flavors)
2571 : : !ii do it=1,op%Bath%sumtails
2572 : : !ii write(6,'(a,100f12.3)') " M after ",(op%Bath%M%mat(it,it1),it1=1,op%Bath%sumtails)
2573 : : !ii enddo
2574 : : !ii do it=1,op%Bath%sumtails
2575 : : !ii write(6,'(a,100f12.3)') " update M after ",(op%Bath%M_update%mat(it,it1),it1=1,op%Bath%sumtails)
2576 : : !ii enddo
2577 : : !ii write(6,*) "Gmove accepted",rnd,local_ratio*det_ratio
2578 : 0 : CALL ImpurityOperator_swap(op%Impurity, flavor_i,flavor_j)
2579 : 0 : CALL BathOperatoroffdiagComplex_swap (op%Bath , flavor_i,flavor_j) ! use op%Bath%M_update to built new op%Bath%M
2580 : :
2581 : 0 : op%swap = op%swap + 1.d0
2582 : 0 : flav_i = flavor_i
2583 : 0 : flav_j = flavor_j
2584 : : ELSE
2585 : : !ii write(6,*) "Gmove refused",rnd,local_ratio*det_ratio
2586 : : ! CALL WARN("Swap refused")
2587 : : ! WRITE(op%ostream,'(6E24.14)') local_ratio, det_ratio, detic1, detjc1, detic2, detjc2
2588 : : END IF
2589 : : ! CALL BathOperatoroffdiag_destroy(Bathnew)
2590 : : END IF
2591 : : !ii do iflavor=1,op%flavors
2592 : : !ii write(6,*) "AFTER GMOVE For flavor", iflavor,"size is",op%Impurity%particles(iflavor)%tail," and Conf is :"
2593 : : !ii do ii=1, op%Impurity%Particles(iflavor)%tail
2594 : : !ii write(6,'(15x,i4,100f12.3)') ii, op%Impurity%Particles(iflavor)%list(ii,1),&
2595 : : !ii & op%Impurity%Particles(iflavor)%list(ii,2)
2596 : : !ii enddo
2597 : : !ii enddo
2598 : : !ii write(6,*) " = M Matrix"
2599 : : !ii write(6,'(a,2x,100(i12))') "Flavor=",((iflavor,it=1,op%Impurity%particles(iflavor)%tail),iflavor=1,op%flavors)
2600 : : !ii write(6,'(i21,100i12)') ((it,it=1,op%Impurity%particles(iflavor)%tail),iflavor=1,op%flavors)
2601 : : !ii do it=1,op%Bath%sumtails
2602 : : !ii write(6,'(a,100f12.3)') " M after ",(op%Bath%M%mat(it,it1),it1=1,op%Bath%sumtails)
2603 : : !ii enddo
2604 : : !ii do it=1,op%Bath%sumtails
2605 : : !ii write(6,'(a,100f12.3)') " update M after ",(op%Bath%M_update%mat(it,it1),it1=1,op%Bath%sumtails)
2606 : : !ii enddo
2607 : :
2608 : 0 : END SUBROUTINE CtqmcoffdiagComplex_trySwap
2609 : : !!***
2610 : :
2611 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_measN
2612 : : !! NAME
2613 : : !! CtqmcoffdiagComplex_measN
2614 : : !!
2615 : : !! FUNCTION
2616 : : !! measures the number of electron
2617 : : !! by taking into account the value for the move before before this one
2618 : : !! with the correct weight.
2619 : : !!
2620 : : !! COPYRIGHT
2621 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
2622 : : !! This file is distributed under the terms of the
2623 : : !! GNU General Public License, see ~abinit/COPYING
2624 : : !! or http://www.gnu.org/copyleft/gpl.txt .
2625 : : !!
2626 : : !! INPUTS
2627 : : !! op=ctqmc
2628 : : !! iflavor=which flavor to measure
2629 : : !! updated=something has changed since last time
2630 : : !!
2631 : : !! OUTPUT
2632 : : !!
2633 : : !! SIDE EFFECTS
2634 : : !!
2635 : : !! NOTES
2636 : : !!
2637 : : !! SOURCE
2638 : :
2639 : 0 : SUBROUTINE CtqmcoffdiagComplex_measN(op, iflavor, updated)
2640 : :
2641 : : !Arguments ------------------------------------
2642 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
2643 : : !TYPE(ImpurityOperator), INTENT(IN ) :: impurity
2644 : : INTEGER , INTENT(IN ) :: iflavor
2645 : : LOGICAL , INTENT(IN ) :: updated
2646 : :
2647 : : ! IF ( .NOT. op%set ) &
2648 : : ! CALL ERROR("CtqmcoffdiagComplex_measN : QMC not set ")
2649 : :
2650 : :
2651 [ # # ]: 0 : IF ( updated .EQV. .TRUE. ) THEN
2652 : : ! --- accumulate occupations with values op%measN(3,iflavor) from the last measurements with the corresponding weight
2653 : : ! --- op*measN(4,iflavor)
2654 : 0 : op%measN(1,iflavor) = op%measN(1,iflavor) + op%measN(3,iflavor)*op%measN(4,iflavor)
2655 : : ! write(6,*) "Cllll42"
2656 : :
2657 : : ! --- Compute total number of new measurements
2658 : 0 : op%measN(2,iflavor) = op%measN(2,iflavor) + op%measN(4,iflavor)
2659 : :
2660 : : ! write(6,*) "Allll42"
2661 : : ! --- Compute the occupation for this configuration (will be put in
2662 : : ! --- op%measN(1,iflavor) at the next occurence of updated=.true.), with
2663 : : ! --- the corresponding weight op%measN(4,iflavor) (we do not now it yet)
2664 : 0 : op%measN(3,iflavor) = ImpurityOperator_measN(op%impurity)
2665 : :
2666 : : ! --- set weight: as update=true, it is a new measurement , so put it to one
2667 : 0 : op%measN(4,iflavor) = 1.d0
2668 : :
2669 : : ELSE
2670 : : ! --- increased the count so that at new move, we will be able to update measN(1) correctly.
2671 : 0 : op%measN(4,iflavor) = op%measN(4,iflavor) + 1.d0
2672 : : ! write(6,*) "Bllll42"
2673 : : END IF
2674 : 0 : END SUBROUTINE CtqmcoffdiagComplex_measN
2675 : : !!***
2676 : :
2677 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
2678 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_measCorrelation
2679 : : !! NAME
2680 : : !! CtqmcoffdiagComplex_measCorrelation
2681 : : !!
2682 : : !! FUNCTION
2683 : : !! measure all correlations in times for a flavor
2684 : : !!
2685 : : !! COPYRIGHT
2686 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
2687 : : !! This file is distributed under the terms of the
2688 : : !! GNU General Public License, see ~abinit/COPYING
2689 : : !! or http://www.gnu.org/copyleft/gpl.txt .
2690 : : !!
2691 : : !! INPUTS
2692 : : !! op=ctqmc
2693 : : !! iflavor=the flavor to measure
2694 : : !!
2695 : : !! OUTPUT
2696 : : !!
2697 : : !! SIDE EFFECTS
2698 : : !!
2699 : : !! NOTES
2700 : : !!
2701 : : !! SOURCE
2702 : :
2703 : 0 : SUBROUTINE CtqmcoffdiagComplex_measCorrelation(op, iflavor)
2704 : :
2705 : : !Arguments ------------------------------------
2706 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
2707 : : !TYPE(ImpurityOperator), INTENT(IN ) :: impurity
2708 : : INTEGER , INTENT(IN ) :: iflavor
2709 : : !Local variables ------------------------------
2710 : : INTEGER :: iCdag
2711 : : INTEGER :: iCdagBeta
2712 : : INTEGER :: iC
2713 : : INTEGER :: index
2714 : : INTEGER :: size
2715 : : DOUBLE PRECISION :: tC
2716 : : DOUBLE PRECISION :: tCdag
2717 : : !DOUBLE PRECISION :: time
2718 : : DOUBLE PRECISION :: inv_dt
2719 : : DOUBLE PRECISION :: beta
2720 : :
2721 [ # # ]: 0 : IF ( .NOT. op%set ) &
2722 : 0 : CALL ERROR("CtqmcoffdiagComplex_measCorrelation : QMC not set ")
2723 : : !write(6,*) "not available"
2724 : 0 : stop
2725 : :
2726 : : size = op%impurity%particles(op%impurity%activeFlavor)%tail
2727 : : beta = op%beta
2728 : :
2729 : : IF ( size .EQ. 0 ) RETURN
2730 : :
2731 : : inv_dt = op%inv_dt
2732 : :
2733 : : DO iCdag = 1, size ! first segments
2734 : : tCdag = op%impurity%particles(op%impurity%activeFlavor)%list(iCdag,Cdag_)
2735 : : tC = op%impurity%particles(op%impurity%activeFlavor)%list(iCdag,C_ )
2736 : : index = INT( ( (tC - tCdag) * inv_dt ) + .5d0 ) + 1
2737 : : op%measCorrelation(index,1,iflavor) = op%measCorrelation(index,1,iflavor) + 1.d0
2738 : : MODCYCLE(iCdag+1,size,iCdagBeta)
2739 : : index = INT( ( ( &
2740 : : op%impurity%particles(op%impurity%activeFlavor)%list(iCdagBeta,Cdag_) - tC &
2741 : : + AINT(DBLE(iCdag)/DBLE(size))*beta &
2742 : : ) * inv_dt ) + .5d0 ) + 1
2743 : : IF ( index .LT. 1 .OR. index .GT. op%samples+1 ) THEN
2744 : : CALL WARN("CtqmcoffdiagComplex_measCorrelation : bad index line 1095 ")
2745 : : ELSE
2746 : : op%measCorrelation(index,2,iflavor) = op%measCorrelation(index,2,iflavor) + 1.d0
2747 : : END IF
2748 : : ! DO iC = 1, size
2749 : : ! tC = impurity%particles(impurity%activeFlavor)%list(C_,iC)
2750 : : ! time = tC - tCdag
2751 : : ! IF ( time .LT. 0.d0 ) time = time + beta
2752 : : ! index = INT( ( time * inv_dt ) + .5d0 ) + 1
2753 : : ! op%measCorrelation(index,3,iflavor) = op%measCorrelation(index,3,iflavor) + 1.d0
2754 : : ! END DO
2755 : : DO iC = 1, size! op%Greens(iflavor)%index_old%tail
2756 : : !todoba op%measCorrelation(op%Greens(iflavor)%map%listINT(iC+(iCdag-1)*size),3,iflavor) = &
2757 : : !todoba op%measCorrelation(op%Greens(iflavor)%map%listINT(iC+(iCdag-1)*size),3,iflavor) + 1.d0
2758 : : END DO
2759 : : END DO
2760 : :
2761 : : END SUBROUTINE CtqmcoffdiagComplex_measCorrelation
2762 : : !!***
2763 : :
2764 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_measPerturbation
2765 : : !! NAME
2766 : : !! CtqmcoffdiagComplex_measPerturbation
2767 : : !!
2768 : : !! FUNCTION
2769 : : !! measure perturbation order
2770 : : !!
2771 : : !! COPYRIGHT
2772 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
2773 : : !! This file is distributed under the terms of the
2774 : : !! GNU General Public License, see ~abinit/COPYING
2775 : : !! or http://www.gnu.org/copyleft/gpl.txt .
2776 : : !!
2777 : : !! INPUTS
2778 : : !! op=ctqmc
2779 : : !! iflavor=the flavor to measure
2780 : : !!
2781 : : !! OUTPUT
2782 : : !!
2783 : : !! SIDE EFFECTS
2784 : : !!
2785 : : !! NOTES
2786 : : !!
2787 : : !! SOURCE
2788 : :
2789 : 0 : SUBROUTINE CtqmcoffdiagComplex_measPerturbation(op, iflavor)
2790 : :
2791 : : !Arguments ------------------------------------
2792 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
2793 : : !TYPE(ImpurityOperator), INTENT(IN ) :: impurity
2794 : : INTEGER , INTENT(IN ) :: iflavor
2795 : : !Local variables ------------------------------
2796 : : INTEGER :: index
2797 : :
2798 [ # # ]: 0 : IF ( .NOT. op%set ) &
2799 : 0 : CALL ERROR("CtqmcoffdiagComplex_measiPerturbation : QMC not set ")
2800 : :
2801 : 0 : index = op%impurity%particles(op%impurity%activeFlavor)%tail + 1
2802 [ # # ]: 0 : IF ( index .LE. op%opt_order ) &
2803 : 0 : op%measPerturbation(index,iflavor) = op%measPerturbation(index,iflavor) + 1.d0
2804 [ # # ]: 0 : IF ( index == 1 ) THEN
2805 [ # # ]: 0 : IF (op%impurity%particles(iflavor)%list(0,C_) < op%impurity%particles(iflavor)%list(0,Cdag_) ) THEN
2806 : 0 : op%meas_fullemptylines(1,iflavor) = op%meas_fullemptylines(1,iflavor) + 1.d0
2807 : : ELSE
2808 : 0 : op%meas_fullemptylines(2,iflavor) = op%meas_fullemptylines(2,iflavor) + 1.d0
2809 : : ENDIF
2810 : : ENDIF
2811 : :
2812 : 0 : END SUBROUTINE CtqmcoffdiagComplex_measPerturbation
2813 : : !!***
2814 : : !#endif
2815 : :
2816 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_getResult
2817 : : !! NAME
2818 : : !! CtqmcoffdiagComplex_getResult
2819 : : !!
2820 : : !! FUNCTION
2821 : : !! reduce everything to get the result of the simulation
2822 : : !!
2823 : : !! COPYRIGHT
2824 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder,F. Gendron)
2825 : : !! This file is distributed under the terms of the
2826 : : !! GNU General Public License, see ~abinit/COPYING
2827 : : !! or http://www.gnu.org/copyleft/gpl.txt .
2828 : : !!
2829 : : !! INPUTS
2830 : : !! op=ctqmc
2831 : : !!
2832 : : !! OUTPUT
2833 : : !!
2834 : : !! SIDE EFFECTS
2835 : : !!
2836 : : !! NOTES
2837 : : !!
2838 : : !! SOURCE
2839 : :
2840 : 0 : SUBROUTINE CtqmcoffdiagComplex_getResult(op,Iatom,fname)
2841 : :
2842 : :
2843 : : #ifdef HAVE_MPI1
2844 : : include 'mpif.h'
2845 : : #endif
2846 : : !Arguments ------------------------------------
2847 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
2848 : : INTEGER, INTENT(IN ) :: Iatom
2849 : : character(len=fnlen), INTENT(INOUT) :: fname
2850 : : !Local variables ------------------------------
2851 : : INTEGER :: iflavor
2852 : : ! INTEGER :: iflavor1
2853 : : ! INTEGER :: iflavor2
2854 : : INTEGER :: flavors
2855 : : INTEGER :: itau
2856 : : INTEGER :: endDensity
2857 : : character(len=2) :: atomnb
2858 : : DOUBLE PRECISION :: inv_flavors
2859 : : DOUBLE PRECISION :: a
2860 : : DOUBLE PRECISION :: b
2861 : : DOUBLE PRECISION :: r
2862 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: alpha
2863 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: beta
2864 : 0 : COMPLEX(kind=8), ALLOCATABLE, DIMENSION(:) :: measN_1
2865 : : DOUBLE PRECISION, DIMENSION(1:2) :: TabX
2866 : : DOUBLE PRECISION, DIMENSION(1:2) :: TabY
2867 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: freqs
2868 : 0 : INTEGER, ALLOCATABLE, DIMENSION(:) :: counts
2869 : 0 : INTEGER, ALLOCATABLE, DIMENSION(:) :: displs
2870 : 0 : INTEGER, ALLOCATABLE, DIMENSION(:) :: occtot
2871 : 0 : INTEGER, ALLOCATABLE, DIMENSION(:) :: spintot
2872 : 0 : INTEGER, ALLOCATABLE, DIMENSION(:,:) :: occ
2873 : : INTEGER :: sp1,spinmax,spinmin,dspin,nelec,spin
2874 : : INTEGER :: spAll
2875 : : INTEGER :: last
2876 : : INTEGER :: n1
2877 : : INTEGER :: n2,n3,quotient,remainder,signe
2878 : : INTEGER :: debut
2879 : : DOUBLE PRECISION :: signvaluemeassum
2880 : : COMPLEX(KIND=8) :: phasevaluemeassum
2881 : : ! INTEGER :: fin
2882 : : #ifdef HAVE_MPI
2883 : : INTEGER :: ierr
2884 : : DOUBLE PRECISION, DIMENSION(1) :: arr
2885 : : #endif
2886 : : INTEGER :: sizeoper,nbprocs,myrank
2887 : : DOUBLE PRECISION :: inv_size,sumh,sumtot
2888 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: buffer
2889 : 0 : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:,:,:) :: buffer2,buffer2s
2890 : 0 : DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: fullempty
2891 : 0 : TYPE(FFTHyb) :: FFTmrka
2892 : :
2893 : : #if defined HAVE_MPI && !defined HAVE_MPI2_INPLACE
2894 : : DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:) :: buffer1_out,freqs_buf
2895 : : DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:,:) :: buffer2_out
2896 : : #endif
2897 : :
2898 [ # # ]: 0 : IF ( .NOT. op%done ) &
2899 : 0 : CALL ERROR("CtqmcoffdiagComplex_getResult : Simulation not run ")
2900 : :
2901 : 0 : flavors = op%flavors
2902 : 0 : inv_flavors = 1.d0 / DBLE(flavors)
2903 : :
2904 : :
2905 : 0 : inv_size = 1.d0 / DBLE(op%size)
2906 : 0 : sp1 = 0
2907 : 0 : spAll = 0
2908 : :
2909 : : !#ifdef CTCtqmcoffdiagComplex_CHECK
2910 [ # # ]: 0 : IF ( op%opt_check .GT. 0 ) THEN
2911 : 0 : op%errorImpurity = ImpurityOperator_getError(op%Impurity) * inv_flavors
2912 : 0 : op%errorBath = BathOperatoroffdiagComplex_getError(op%Bath) * inv_flavors
2913 : : END IF
2914 : : !#endif
2915 : :
2916 [ # # ]: 0 : MALLOC(alpha,(1,1))
2917 [ # # ]: 0 : MALLOC(beta,(1,1))
2918 [ # # ]: 0 : MALLOC(buffer,(1,1))
2919 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1) THEN
2920 [ # # ]: 0 : FREEIF(alpha)
2921 [ # # # # : 0 : MALLOC(alpha,(1:op%samples+1,1:flavors))
# # # # #
# # # # #
# # # # ]
2922 [ # # ]: 0 : FREEIF(beta)
2923 [ # # # # : 0 : MALLOC(beta,(1:op%samples+1,1:flavors))
# # # # #
# # # #
# ]
2924 : : END IF
2925 : :
2926 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE.) THEN
2927 : 0 : sp1 = 0
2928 : 0 : spAll = sp1 + flavors + 6
2929 : :
2930 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
2931 : : IF ( op%opt_analysis .EQ. 1 ) &
2932 : : spAll = spAll + 3*sp1
2933 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) &
2934 : 0 : spAll = spAll + op%opt_order
2935 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) &
2936 : 0 : spAll = spAll + 2*(op%samples + 1)
2937 : : !#endif
2938 : :
2939 [ # # ]: 0 : FREEIF(buffer)
2940 [ # # # # : 0 : MALLOC(buffer,(1:spAll,1:MAX(2,flavors)))
# # # # #
# # # # #
# # # # ]
2941 : : END IF
2942 : :
2943 : : ! op%seg_added = op%seg_added * inv_flavors
2944 : : ! op%seg_removed = op%seg_removed * inv_flavors
2945 : : ! op%seg_sign = op%seg_sign * inv_flavors
2946 : : ! op%anti_added = op%anti_added * inv_flavors
2947 : : ! op%anti_removed = op%anti_removed * inv_flavors
2948 : : ! op%anti_sign = op%anti_sign * inv_flavors
2949 [ # # ]: 0 : op%stats(:) = op%stats(:) * inv_flavors
2950 : :
2951 [ # # ]: 0 : DO iflavor = 1, flavors
2952 : : ! Accumulate last values of N (see also ctqmc_measn)
2953 : 0 : op%measN(1,iflavor) = op%measN(1,iflavor) + op%measN(3,iflavor)*op%measN(4,iflavor)
2954 : 0 : op%measN(2,iflavor) = op%measN(2,iflavor) + op%measN(4,iflavor)
2955 : : ! Reduction
2956 : 0 : op%measN(1,iflavor) = op%measN(1,iflavor) / ( op%measN(2,iflavor) * op%beta )
2957 : : ! Correction
2958 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
2959 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) &
2960 : : op%measPerturbation(: ,iflavor) = op%measPerturbation(:,iflavor) &
2961 [ # # # # ]: 0 : / SUM(op%measPerturbation(:,iflavor))
2962 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) &
2963 : : op%meas_fullemptylines(: ,iflavor) = op%meas_fullemptylines(:,iflavor) &
2964 [ # # # # ]: 0 : / SUM(op%meas_fullemptylines(:,iflavor))
2965 : : !write(6,*) "sum fullempty",iflavor,op%meas_fullemptylines(:,iflavor)
2966 : :
2967 [ # # ]: 0 : IF ( op%opt_analysis .EQ. 1 ) THEN
2968 : : op%measCorrelation (:,1,iflavor) = op%measCorrelation (:,1,iflavor) &
2969 : : / SUM(op%measCorrelation (:,1,iflavor)) &
2970 [ # # # # ]: 0 : * op%inv_dt
2971 : : op%measCorrelation (:,2,iflavor) = op%measCorrelation (:,2,iflavor) &
2972 : : / SUM(op%measCorrelation (:,2,iflavor)) &
2973 [ # # # # ]: 0 : * op%inv_dt
2974 : : op%measCorrelation (:,3,iflavor) = op%measCorrelation (:,3,iflavor) &
2975 : : / SUM(op%measCorrelation (:,3,iflavor)) &
2976 [ # # # # ]: 0 : * op%inv_dt
2977 : : END IF
2978 : : !#endif
2979 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
2980 : 0 : TabX(1) = DBLE(op%modNoise2)
2981 : 0 : TabX(2) = DBLE(op%modNoise1)
2982 [ # # ]: 0 : DO itau = 1, op%samples+1
2983 : : op%measNoiseG(itau,iflavor,2)%vec = -op%measNoiseG(itau,iflavor,2)%vec*op%inv_dt &
2984 [ # # ]: 0 : /(op%beta*DBLE(op%modNoise2))
2985 : : op%measNoiseG(itau,iflavor,1)%vec = -op%measNoiseG(itau,iflavor,1)%vec*op%inv_dt &
2986 [ # # ]: 0 : /(op%beta*DBLE(op%modNoise1))
2987 : 0 : n2 = op%measNoiseG(itau,iflavor,2)%tail
2988 : 0 : TabY(1) = Stat_deviation(op%measNoiseG(itau,iflavor,2)%vec(1:n2))!*SQRT(n2/(n2-1))
2989 : 0 : n1 = op%measNoiseG(itau,iflavor,1)%tail
2990 : 0 : TabY(2) = Stat_deviation(op%measNoiseG(itau,iflavor,1)%vec(1:n1))!*SQRT(n1/(n1-1))
2991 [ # # ]: 0 : CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,alpha(itau,iflavor),beta(itau,iflavor),r)
2992 : : ! ecart type -> 60%
2993 : : ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma
2994 : : END DO
2995 : : END IF
2996 : :
2997 : : END DO
2998 : : !sui!write(6,*) "getresults"
2999 : 0 : CALL GreenHyboffdiagComplex_measHybrid(op%Greens, op%Bath%M, op%Impurity%Particles, .TRUE.,op%signvalue,op%phasevalue)
3000 : 0 : CALL GreenHyboffdiagComplex_getHybrid(op%Greens)
3001 : : !write(6,*) "op%measN",op%measN(1,:)
3002 [ # # # # : 0 : MALLOC(measN_1,(flavors))
# # # # #
# # # ]
3003 [ # # ]: 0 : do iflavor=1,flavors
3004 : 0 : measN_1(iflavor)=op%measN(1,iflavor)
3005 : : enddo
3006 : 0 : CALL GreenHyboffdiagComplex_setN(op%Greens, measN_1(:))
3007 [ # # ]: 0 : FREE(measN_1)
3008 : :
3009 : : ! todoab case _nd and _d are not completely described.
3010 : : FREEIF(buffer2)
3011 : : FREEIF(buffer2s)
3012 : 0 : sizeoper=size(op%Greens%oper,1)
3013 : : !write(6,*) "sss",size(op%Greens%oper,1),sizeoper
3014 : : !write(6,*) "sss",size(op%Greens%oper,2),flavors
3015 : : !write(6,*) "sss",size(op%Greens%oper,3),flavors
3016 [ # # # # : 0 : MALLOC(buffer2,(1:sizeoper,flavors,flavors))
# # # # #
# # # # #
# # # # #
# ]
3017 [ # # # # : 0 : MALLOC(buffer2s,(1:sizeoper,flavors,flavors))
# # # # #
# # # # #
# # ]
3018 [ # # # # : 0 : MALLOC(fullempty,(2,flavors))
# # # # #
# # # ]
3019 : : !sui!write(6,*) "greens1"
3020 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE. ) THEN
3021 : : !sui!write(6,*) "greens2"
3022 [ # # # # ]: 0 : fullempty=0.d0
3023 [ # # # # : 0 : buffer2 = op%Greens%oper
# # # # #
# # # # #
# # # # #
# # # ]
3024 : : !write(6,*) "buffer2",(op%Greens%oper(1,n1,n1),n1=1,flavors)
3025 [ # # # # : 0 : buffer2s= 0.d0
# # ]
3026 : : do iflavor=1,flavors
3027 : : do itau=1,sizeoper
3028 : : !sui!write(6,*) "greens",iflavor,itau,op%Greens%oper(itau,iflavor,iflavor)
3029 : : enddo
3030 : : enddo
3031 : : !write(6,*) "beforempi",op%Greens%oper(1,1,1) ,buffer2(1,1,1)
3032 : : #ifdef HAVE_MPI
3033 : 0 : CALL MPI_COMM_SIZE(op%MY_COMM,nbprocs,ierr)
3034 : 0 : CALL MPI_COMM_RANK(op%MY_COMM,myrank,ierr)
3035 : : #endif
3036 : : !write(6,*) "procs",nbprocs,myrank
3037 : : END IF
3038 : 0 : last = sp1
3039 : :
3040 [ # # # # ]: 0 : op%measDE(:,:) = op%measDE(:,:) * DBLE(op%measurements) /(DBLE(op%sweeps)*op%beta)
3041 : :
3042 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
3043 [ # # ]: 0 : op%occup_histo_time(:) = op%occup_histo_time(:) / INT(op%sweeps/op%measurements)
3044 [ # # ]: 0 : op%occupconfig(:) = op%occupconfig(:) / INT(op%sweeps/op%measurements)
3045 [ # # # # ]: 0 : op%suscep(:,:) = op%suscep(:,:) / INT(op%sweeps/op%measurements)
3046 [ # # # # ]: 0 : op%chi(:,:) = op%chi(:,:) / INT(op%sweeps/op%measurements)
3047 [ # # # # ]: 0 : op%chicharge(:,:) = op%chicharge(:,:) / INT(op%sweeps/op%measurements)
3048 [ # # ]: 0 : op%ntot(:) = op%ntot(:) / INT(op%sweeps/op%measurements)
3049 : : END IF
3050 : : ! HISTO before MPI_SUM
3051 : : ! write(6,*) "=== Histogram of occupations for complete simulation 3 ====",INT(op%sweeps/op%measurements)
3052 : : ! sumh=0
3053 : : ! do n1=1,op%flavors+1
3054 : : ! write(6,'(i4,f10.4)') n1-1, op%occup_histo_time(n1)
3055 : : ! sumh=sumh+op%occup_histo_time(n1)
3056 : : ! enddo
3057 : : ! write(6,*) "=================================",sumh
3058 : :
3059 : 0 : n1 = op%measNoise(1)%tail
3060 : 0 : n2 = op%measNoise(2)%tail
3061 : :
3062 : : ! On utilise freqs comme tableau de regroupement
3063 : : ! Gather de Noise1
3064 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE. ) THEN
3065 [ # # # # : 0 : MALLOC(counts,(1:op%size))
# # # # #
# # # ]
3066 [ # # # # : 0 : MALLOC(displs,(1:op%size))
# # # # ]
3067 : : FREEIF(freqs)
3068 [ # # # # : 0 : MALLOC(freqs,(1:op%size*n1))
# # # # #
# # # ]
3069 [ # # ]: 0 : freqs = 0.d0
3070 [ # # ]: 0 : counts(:) = n1
3071 [ # # # # : 0 : displs(:) = (/ ( iflavor*n1, iflavor=0, op%size-1 ) /)
# # ]
3072 : : #ifdef HAVE_MPI
3073 : : #if defined HAVE_MPI2_INPLACE
3074 [ # # ]: 0 : freqs(n1*op%rank+1:n1*(op%rank+1)) = op%measNoise(1)%vec(1:n1)
3075 : : CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, &
3076 : : freqs, counts, displs, &
3077 : 0 : MPI_DOUBLE_PRECISION, op%MY_COMM, ierr)
3078 : : #else
3079 : : MALLOC(freqs_buf,(n1))
3080 : : freqs_buf(1:n2)=op%measNoise(2)%vec(1:n1)
3081 : : CALL MPI_ALLGATHERV(freqs_buf, n1, MPI_DOUBLE_PRECISION, &
3082 : : freqs, counts, displs, &
3083 : : MPI_DOUBLE_PRECISION, op%MY_COMM, ierr)
3084 : : FREE(freqs_buf)
3085 : : #endif
3086 : : #endif
3087 : 0 : n1 = op%size*n1
3088 : 0 : CALL Vector_setSize(op%measNoise(1),n1)
3089 [ # # ]: 0 : op%measNoise(1)%vec(1:n1) = freqs(:)
3090 : : ! Gather de Noise2
3091 [ # # ]: 0 : FREE(freqs)
3092 [ # # # # : 0 : MALLOC(freqs,(1:op%size*n2))
# # # # #
# # # ]
3093 [ # # ]: 0 : freqs = 0.d0
3094 [ # # ]: 0 : counts(:) = n2
3095 [ # # # # ]: 0 : displs(:) = (/ ( iflavor*n2, iflavor=0, op%size-1 ) /)
3096 : : #ifdef HAVE_MPI
3097 : : #if defined HAVE_MPI2_INPLACE
3098 [ # # ]: 0 : freqs(n2*op%rank+1:n2*(op%rank+1)) = op%measNoise(2)%vec(1:n2)
3099 : : CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, &
3100 : : freqs, counts, displs, &
3101 : 0 : MPI_DOUBLE_PRECISION, op%MY_COMM, ierr)
3102 : : #else
3103 : : MALLOC(freqs_buf,(n2))
3104 : : freqs_buf(1:n2)=op%measNoise(2)%vec(1:n2)
3105 : : CALL MPI_ALLGATHERV(freqs_buf, n2, MPI_DOUBLE_PRECISION, &
3106 : : freqs, counts, displs, &
3107 : : MPI_DOUBLE_PRECISION, op%MY_COMM, ierr)
3108 : : FREE(freqs_buf)
3109 : : #endif
3110 : : #endif
3111 : 0 : n2 = op%size*n2
3112 : 0 : CALL Vector_setSize(op%measNoise(2),n2)
3113 [ # # ]: 0 : op%measNoise(2)%vec(1:n2) = freqs(:)
3114 [ # # ]: 0 : FREE(counts)
3115 [ # # ]: 0 : FREE(displs)
3116 [ # # ]: 0 : FREE(freqs)
3117 : : END IF
3118 : : !n1 = op%measNoise(1)%tail
3119 : : !n2 = op%measNoise(2)%tail
3120 : :
3121 : : ! Transformation des paquets pour que ca fit a CTQMC_SLICE(1|2)
3122 [ # # ]: 0 : IF ( n1 .GT. CTQMC_SLICE1 ) THEN
3123 : 0 : itau = n1/CTQMC_SLICE1
3124 [ # # # # : 0 : MALLOC(freqs,(1:n1/itau))
# # # # #
# # # #
# ]
3125 [ # # ]: 0 : DO debut=1, n1/itau
3126 [ # # ]: 0 : freqs(debut)=SUM(op%measNoise(1)%vec((debut-1)*itau+1:itau*debut))
3127 : : END DO
3128 [ # # ]: 0 : freqs(:) = freqs(:)/DBLE(itau)
3129 : 0 : op%modNoise1 = op%modNoise1*itau
3130 : 0 : n1 = n1/itau
3131 : 0 : CALL Vector_setSize(op%measNoise(1),n1)
3132 [ # # ]: 0 : op%measNoise(1)%vec(1:n1) = freqs(:)
3133 [ # # ]: 0 : FREE(freqs)
3134 : : END IF
3135 [ # # ]: 0 : IF ( n2 .GT. CTQMC_SLICE1*CTQMC_SLICE2 ) THEN
3136 : 0 : itau = n2/(CTQMC_SLICE1*CTQMC_SLICE2)
3137 [ # # # # : 0 : MALLOC(freqs,(1:n2/itau))
# # # # #
# # # #
# ]
3138 [ # # ]: 0 : DO debut=1, n2/itau
3139 [ # # ]: 0 : freqs(debut)=SUM(op%measNoise(2)%vec((debut-1)*itau+1:itau*debut))
3140 : : END DO
3141 [ # # ]: 0 : freqs(:) = freqs(:)/DBLE(itau)
3142 : 0 : op%modNoise2 = op%modNoise2*itau
3143 : 0 : n2 = n2/itau
3144 : 0 : CALL Vector_setSize(op%measNoise(2),n2)
3145 [ # # ]: 0 : op%measNoise(2)%vec(1:n2) = freqs(:)
3146 [ # # ]: 0 : FREE(freqs)
3147 : : END IF
3148 : : ! On peut s'amuser avec nos valeur d'energies
3149 : : !MALLOC(TabX,(1:20))
3150 : : !MALLOC(TabY,(1:20))
3151 : :
3152 : 0 : TabX(1) = DBLE(op%modNoise2)
3153 : 0 : TabX(2) = DBLE(op%modNoise1)
3154 : :
3155 : : ! Il faut calculer pour chaque modulo 10 ecarts type sur les donnes acquises
3156 [ # # ]: 0 : op%measNoise(1)%vec(1:n1) = op%measNoise(1)%vec(1:n1)/(op%beta*DBLE(op%modNoise1))*DBLE(op%measurements)
3157 [ # # ]: 0 : op%measNoise(2)%vec(1:n2) = op%measNoise(2)%vec(1:n2)/(op%beta*DBLE(op%modNoise2))*DBLE(op%measurements)
3158 : : ! CALL Vector_print(op%measNoise(1),op%rank+70)
3159 : : ! CALL Vector_print(op%measNoise(2),op%rank+50)
3160 : : ! DO iflavor=1,10
3161 : : ! debut = (iflavor-1)*n2/10+1
3162 : : ! fin = iflavor*n2/10
3163 : : ! TabY(iflavor) = Stat_deviation(op%measNoise(2)%vec(debut:fin))
3164 : : ! debut = (iflavor-1)*n1/10+1
3165 : : ! fin = iflavor*n1/10
3166 : : ! TabY(10+iflavor) = Stat_deviation(op%measNoise(1)%vec(debut:fin))
3167 : : ! END DO
3168 : : !! TabY(1:n) = (op%measNoise(2)%vec(1:n) &
3169 : : !! )
3170 : : !! !/(op%beta*DBLE(op%modNoise2))*DBLE(op%measurements) &
3171 : : !! !- op%measDE(1,1))
3172 : : !! TabY(op%measNoise(2)%tail+1:n+op%measNoise(2)%tail) = (op%measNoise(1)%vec(1:n) &
3173 : : !! )
3174 : : !! ! /(op%beta*DBLE(op%modNoise1))*DBLE(op%measurements) &
3175 : : !! ! - op%measDE(1,1))
3176 : : ! IF ( op%rank .EQ. 0 ) THEN
3177 : : ! DO iflavor=1,20
3178 : : ! write(45,*) TabX(iflavor), TabY(iflavor)
3179 : : ! END DO
3180 : : ! END IF
3181 : : !
3182 : :
3183 : :
3184 : 0 : TabY(1) = Stat_deviation(op%measNoise(2)%vec(1:n2))!*SQRT(n2/(n2-1))
3185 : : !! write(op%rank+10,*) TabX(2)
3186 : : !! write(op%rank+40,*) TabX(1)
3187 : : !! CALL Vector_print(op%measNoise(1),op%rank+10)
3188 : : !! CALL Vector_print(op%measNoise(2),op%rank+40)
3189 : : !! CLOSE(op%rank+10)
3190 : : !! CLOSE(op%rank+40)
3191 : 0 : TabY(2) = Stat_deviation(op%measNoise(1)%vec(1:n1))!*SQRT(n1/(n1-1))
3192 : : !! ! Ecart carre moyen ~ ecart type mais non biaise. Serait moins precis. Aucun
3193 : : ! impact sur la pente, juste sur l'ordonnee a l'origine.
3194 : :
3195 [ # # ]: 0 : CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,a,b,r)
3196 : : ! FREE(TabX)
3197 : : ! FREE(TabY)
3198 : : ! ecart type -> 60%
3199 : : ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma
3200 : :
3201 : : !op%measDE(1,1) = SUM(op%measNoise(1)%vec(1:op%measNoise(1)%tail))/(DBLE(op%measNoise(1)%tail*op%modNoise1)*op%beta)
3202 : : !op%measDE(2:flavors,1:flavors) = op%measDE(2:flavors,1:flavors) /(DBLE(op%sweeps)*op%beta)
3203 : 0 : CALL ImpurityOperator_getErrorOverlap(op%Impurity,op%measDE)
3204 : : ! Add the difference between true calculation and quick calculation of the
3205 : : ! last sweep overlap to measDE(2,2)
3206 : : !op%measDE = op%measDE * DBLE(op%measurements)
3207 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE. ) THEN
3208 [ # # ]: 0 : IF ( op%opt_analysis .EQ. 1 ) THEN
3209 [ # # # # ]: 0 : buffer(last+1:last+sp1,:) = op%measCorrelation(:,1,:)
3210 : 0 : last = last + sp1
3211 [ # # # # ]: 0 : buffer(last+1:last+sp1,:) = op%measCorrelation(:,2,:)
3212 : 0 : last = last + sp1
3213 [ # # # # ]: 0 : buffer(last+1:last+sp1,:) = op%measCorrelation(:,3,:)
3214 : : last = last + sp1
3215 : : END IF
3216 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
3217 [ # # # # ]: 0 : buffer(last+1:last+op%opt_order, :) = op%measPerturbation(:,:)
3218 : : last = last + op%opt_order
3219 : : END IF
3220 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
3221 [ # # # # ]: 0 : buffer(last+1:last+op%samples+1,:) = alpha(:,:)
3222 : 0 : last = last + op%samples + 1
3223 [ # # # # ]: 0 : buffer(last+1:last+op%samples+1,:) = beta(:,:)
3224 : 0 : last = last + op%samples + 1
3225 : : END IF
3226 : : ! op%measDE(2,2) = a*EXP(b*LOG(DBLE(op%sweeps*op%size)))
3227 [ # # # # ]: 0 : buffer(spall-(flavors+5):spAll-6,:) = op%measDE(:,:)
3228 : : ! buffer(spAll ,1) = op%seg_added
3229 : : ! buffer(spAll-1,1) = op%seg_removed
3230 : : ! buffer(spAll-2,1) = op%seg_sign
3231 : : ! buffer(spAll ,2) = op%anti_added
3232 : : ! buffer(spAll-1,2) = op%anti_removed
3233 : : ! buffer(spAll-2,2) = op%anti_sign
3234 : 0 : buffer(spAll ,1) = op%stats(1)
3235 : 0 : buffer(spAll-1,1) = op%stats(2)
3236 : 0 : buffer(spAll-2,1) = op%stats(3)
3237 : 0 : buffer(spAll ,2) = op%stats(4)
3238 : 0 : buffer(spAll-1,2) = op%stats(5)
3239 : 0 : buffer(spAll-2,2) = op%stats(6)
3240 : 0 : buffer(spAll-3,1) = op%swap
3241 : 0 : buffer(spAll-3,2) = DBLE(op%modGlobalMove(2))
3242 : 0 : buffer(spAll-4,1) = a
3243 : 0 : buffer(spAll-4,2) = b
3244 : : !#ifdef CTCtqmcoffdiagComplex_CHECK
3245 : 0 : buffer(spAll-5,1) = op%errorImpurity
3246 : 0 : buffer(spAll-5,2) = op%errorBath
3247 : 0 : signvaluemeassum = 0
3248 : 0 : phasevaluemeassum = 0
3249 : : !#endif
3250 : :
3251 : : #ifdef HAVE_MPI
3252 : : #if defined HAVE_MPI2_INPLACE
3253 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, buffer, spAll*flavors, &
3254 : 0 : MPI_DOUBLE_PRECISION, MPI_SUM, op%MY_COMM, ierr)
3255 : : #else
3256 : : MALLOC(buffer2_out,(spAll,flavors))
3257 : : CALL MPI_ALLREDUCE(buffer, buffer2_out, spAll*flavors, &
3258 : : MPI_DOUBLE_PRECISION, MPI_SUM, op%MY_COMM, ierr)
3259 : : buffer(1:spAll,1:flavors)=buffer2_out(1:spAll,1:flavors)
3260 : : FREE(buffer2_out)
3261 : : #endif
3262 : : CALL MPI_ALLREDUCE(buffer2, buffer2s, sizeoper*flavors*flavors, &
3263 : 0 : MPI_DOUBLE_COMPLEX, MPI_SUM, op%MY_COMM, ierr)
3264 [ # # ]: 0 : CALL MPI_ALLREDUCE([op%runtime], arr, 1, MPI_DOUBLE_PRECISION, MPI_MAX, op%MY_COMM, ierr)
3265 : 0 : op%runtime=arr(1)
3266 [ # # ]: 0 : CALL MPI_ALLREDUCE([op%Greens%signvaluemeas], arr, 1, MPI_DOUBLE_PRECISION, MPI_SUM, op%MY_COMM, ierr)
3267 : 0 : signvaluemeassum=arr(1)
3268 [ # # ]: 0 : CALL MPI_ALLREDUCE([op%Greens%phasevaluemeas], arr , 1, MPI_DOUBLE_COMPLEX, MPI_SUM, op%MY_COMM, ierr)
3269 : 0 : phasevaluemeassum=arr(1)
3270 : : #if defined HAVE_MPI2_INPLACE
3271 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
3272 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE,op%occup_histo_time, flavors+1, MPI_DOUBLE_PRECISION, MPI_SUM, &
3273 : 0 : op%MY_COMM, ierr)
3274 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, op%occupconfig, 2**flavors, MPI_DOUBLE_PRECISION, MPI_SUM, &
3275 [ # # # # ]: 0 : op%MY_COMM, ierr)
3276 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, op%suscep, 3*op%samples, MPI_DOUBLE_PRECISION, MPI_SUM, &
3277 : 0 : op%MY_COMM, ierr)
3278 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, op%chi, 3*op%samples, MPI_DOUBLE_PRECISION, MPI_SUM, &
3279 : 0 : op%MY_COMM, ierr)
3280 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, op%chicharge, 3*op%samples, MPI_DOUBLE_PRECISION, MPI_SUM, &
3281 : 0 : op%MY_COMM, ierr)
3282 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, op%ntot, 3, MPI_DOUBLE_PRECISION, MPI_SUM, &
3283 : 0 : op%MY_COMM, ierr)
3284 : : END IF
3285 : : CALL MPI_ALLREDUCE(MPI_IN_PLACE, sumh, 1, MPI_DOUBLE_PRECISION, MPI_SUM, &
3286 : 0 : op%MY_COMM, ierr)
3287 : : #else
3288 : : MALLOC(buffer1_out,(flavors+1))
3289 : : CALL MPI_ALLREDUCE(op%occup_histo_time, buffer1_out, flavors+1, MPI_DOUBLE_PRECISION, MPI_SUM, &
3290 : : op%MY_COMM, ierr)
3291 : : op%occup_histo_time(1:flavors+1) =buffer1_out(1:flavors+1)
3292 : : FREE(buffer1_out)
3293 : : #endif
3294 [ # # ]: 0 : CALL MPI_ALLREDUCE([sumh], arr, 1, MPI_DOUBLE_PRECISION, MPI_SUM, op%MY_COMM, ierr)
3295 : 0 : sumh=arr(1)
3296 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
3297 : : CALL MPI_ALLREDUCE(op%meas_fullemptylines, fullempty, 2*flavors, MPI_DOUBLE_PRECISION, MPI_SUM, &
3298 : 0 : op%MY_COMM, ierr)
3299 : : ENDIF
3300 : : #endif
3301 : :
3302 : :
3303 [ # # # # ]: 0 : buffer = buffer * inv_size
3304 [ # # # # ]: 0 : op%measDE(:,:) = buffer(spall-(flavors+5):spAll-6,:)
3305 : : ! op%seg_added = buffer(spAll ,1)
3306 : : ! op%seg_removed = buffer(spAll-1,1)
3307 : : ! op%seg_sign = buffer(spAll-2,1)
3308 : : ! op%anti_added = buffer(spAll ,2)
3309 : : ! op%anti_removed = buffer(spAll-1,2)
3310 : : ! op%anti_sign = buffer(spAll-2,2)
3311 : 0 : op%stats(1) = buffer(spAll ,1)
3312 : 0 : op%stats(2) = buffer(spAll-1,1)
3313 : 0 : op%stats(3) = buffer(spAll-2,1)
3314 : 0 : op%stats(4) = buffer(spAll ,2)
3315 : 0 : op%stats(5) = buffer(spAll-1,2)
3316 : 0 : op%stats(6) = buffer(spAll-2,2)
3317 : 0 : op%swap = buffer(spAll-3,1)
3318 : 0 : op%modGlobalMove(2) = NINT(buffer(spAll-3,2))
3319 : 0 : a = buffer(spAll-4,1)
3320 : 0 : b = buffer(spAll-4,2)
3321 : : !!#ifdef CTCtqmcoffdiagComplex_CHECK
3322 : 0 : op%errorImpurity= buffer(spAll-5,1)
3323 : 0 : op%errorBath = buffer(spAll-5,2)
3324 : : !#endif
3325 : :
3326 : : ! DO iflavor = 1, flavors
3327 : : ! op%Greens(iflavor)%oper = buffer(1:sp1 , iflavor)
3328 : : ! END DO
3329 [ # # # # : 0 : op%Greens%oper = buffer2s/float(nbprocs)
# # # # #
# # # # #
# # ]
3330 : : ! write(6,*) "buffer2s",(op%Greens%oper(1,n1,n1),n1=1,flavors)
3331 : 0 : op%Greens%signvaluemeas = signvaluemeassum/float(nbprocs)
3332 : 0 : op%Greens%phasevaluemeas = phasevaluemeassum/float(nbprocs)
3333 : : !sui!write(6,*) "nbprocs",nbprocs,op%Greens%signvaluemeas
3334 [ # # # # : 0 : op%Greens%oper = op%Greens%oper / op%Greens%phasevaluemeas
# # ]
3335 : : ! write(6,*) "buffer3s",(op%Greens%oper(1,n1,n1),n1=1,flavors)
3336 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
3337 [ # # # # : 0 : op%meas_fullemptylines= fullempty/float(nbprocs)
# # # # #
# # # ]
3338 : : ENDIF
3339 : 0 : do iflavor=1,flavors
3340 : : do itau=1,sizeoper
3341 : : !sui!write(6,*) "greens_av",iflavor,itau,op%Greens%oper(itau,iflavor,iflavor)
3342 : : enddo
3343 : : enddo
3344 : : !write(6,*) "aftermpi",op%Greens%oper(1,1,1) ,buffer2s(1,1,1)
3345 : 0 : last = sp1
3346 : : IF ( op%opt_analysis .EQ. 1 ) THEN
3347 : : op%measCorrelation(:,1,:) = buffer(last+1:last+sp1,:)
3348 : : last = last + sp1
3349 : : op%measCorrelation(:,2,:) = buffer(last+1:last+sp1,:)
3350 : : last = last + sp1
3351 : : op%measCorrelation(:,3,:) = buffer(last+1:last+sp1,:)
3352 : : last = last + sp1
3353 : : END IF
3354 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
3355 [ # # # # ]: 0 : op%measPerturbation(:,:) = buffer(last+1:last+op%opt_order, :)
3356 : : last = last + op%opt_order
3357 : : END IF
3358 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
3359 [ # # # # ]: 0 : alpha(:,:) = buffer(last+1:last+op%samples+1,:)
3360 : 0 : last = last + op%samples + 1
3361 [ # # # # ]: 0 : beta(:,:) = buffer(last+1:last+op%samples+1,:)
3362 : 0 : last = last + op%samples + 1
3363 : : END IF
3364 : : END IF
3365 [ # # ]: 0 : DO iflavor = 1, flavors
3366 : : ! complete DE matrix
3367 [ # # ]: 0 : op%measDE(iflavor, iflavor+1:flavors) = op%measDE(iflavor+1:flavors,iflavor)
3368 : : END DO
3369 [ # # ]: 0 : FREE(buffer)
3370 [ # # ]: 0 : FREE(buffer2)
3371 [ # # ]: 0 : FREE(buffer2s)
3372 [ # # ]: 0 : FREE(fullempty)
3373 : :
3374 [ # # ]: 0 : IF ( op%opt_spectra .GE. 1 ) THEN
3375 : : !PROBLEM eos_gnu_13.2_mpich . %endDensity was introduced throughout
3376 : 0 : endDensity = SIZE(op%density,2)
3377 : : ! endDensity=op%endDensity
3378 : : !ENDPROBLEM
3379 : :
3380 [ # # ]: 0 : IF ( op%density(1,endDensity) .EQ. -1.d0 ) &
3381 : 0 : endDensity = endDensity - 1
3382 : 0 : CALL FFTHyb_init(FFTmrka,endDensity,DBLE(op%thermalization)/DBLE(op%measurements*op%opt_spectra))
3383 : : ! Not very Beauty
3384 [ # # # # : 0 : MALLOC(freqs,(1:FFTmrka%size/2))
# # # # #
# # # #
# ]
3385 [ # # ]: 0 : DO iflavor = 1, flavors
3386 : : ! mean value is removed to supress the continue composent
3387 [ # # # # ]: 0 : CALL FFTHyb_setData(FFTmrka,op%density(iflavor,1:endDensity)/op%beta+real(op%Greens%oper(op%samples+1,iflavor,iflavor)))
3388 : 0 : CALL FFTHyb_run(FFTmrka,1)
3389 : 0 : CALL FFTHyb_getData(FFTmrka,endDensity,op%density(iflavor,:),freqs)
3390 : : END DO
3391 [ # # ]: 0 : op%density(flavors+1,:) = -1.d0
3392 [ # # ]: 0 : op%density(flavors+1,1:FFTmrka%size/2) = freqs
3393 : 0 : CALL FFTHyb_destroy(FFTmrka)
3394 [ # # ]: 0 : FREE(freqs)
3395 : : END IF
3396 : :
3397 : 0 : op%a_Noise = a
3398 : 0 : op%b_Noise = b
3399 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1 ) THEN
3400 [ # # # # ]: 0 : op%abNoiseG(1,:,:) = alpha
3401 [ # # # # ]: 0 : op%abNoiseG(2,:,:) = beta
3402 : : END IF
3403 [ # # ]: 0 : FREE(alpha)
3404 [ # # ]: 0 : FREE(beta)
3405 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
3406 : 0 : write(op%ostream,*) "=== Histogram of occupations for complete simulation ===="
3407 : : ! write(6,*) "sumh over procs", sumh
3408 : 0 : sumh=0
3409 [ # # ]: 0 : do n1=1,op%flavors+1
3410 : 0 : write(op%ostream,'(i4,f10.4)') n1-1, op%occup_histo_time(n1)/float(nbprocs)
3411 : 0 : sumh=sumh+op%occup_histo_time(n1)/float(nbprocs)
3412 : : enddo
3413 : 0 : write(op%ostream,'(a,f10.4)') " all" , sumh
3414 : 0 : write(op%ostream,*) "================================="
3415 : :
3416 [ # # # # : 0 : MALLOC(occ,(2**op%flavors,1:flavors))
# # # # #
# # # # #
# # # # #
# # # # #
# # # # ]
3417 : : #ifdef FC_LLVM
3418 : : MALLOC(occtot,(2**op%flavors) )
3419 : : MALLOC(spintot,(2**op%flavors) )
3420 : : #else
3421 [ # # # # : 0 : MALLOC(occtot,(2**op%flavors))
# # # # #
# # # # #
# # # # #
# # # ]
3422 [ # # # # : 0 : MALLOC(spintot,(2**op%flavors))
# # # # #
# # # # #
# # # # #
# # # ]
3423 : : #endif
3424 [ # # # # : 0 : do n1=1,2**op%flavors
# # ]
3425 : : ! Compute occupations of individual Orbitals
3426 : 0 : n3=n1-1
3427 : 0 : occtot(n1)=0
3428 : 0 : spintot(n1)=0
3429 : 0 : signe=1
3430 [ # # ]: 0 : do n2=1,op%flavors
3431 : 0 : remainder=modulo(n3,2)
3432 : 0 : quotient=(n3-remainder)/2
3433 : 0 : occ(n1,n2)=remainder
3434 : 0 : n3=quotient
3435 : 0 : occtot(n1)=occtot(n1)+occ(n1,n2)
3436 [ # # ]: 0 : if(n2>=7) signe =-1
3437 : : !if(n2>=7) signe =0
3438 : 0 : spintot(n1)=spintot(n1)+occ(n1,n2)*signe
3439 : : enddo
3440 : 0 : op%occupconfig(n1)=op%occupconfig(n1)/float(nbprocs)
3441 : : enddo
3442 : :
3443 : 0 : write(op%ostream,*) "=== Histogram of occupations of configurations for complete simulation ===="
3444 : 0 : sumh=0
3445 [ # # ]: 0 : if(op%flavors==14) then
3446 [ # # ]: 0 : do n1=1,2**op%flavors
3447 : 0 : write(op%ostream,'(i4,14i2,f20.2)') n1, (occ(n1,n2),n2=1,op%flavors),op%occupconfig(n1)
3448 : 0 : sumh=sumh+op%occupconfig(n1)
3449 : : enddo
3450 [ # # ]: 0 : else if(op%flavors==10) then
3451 [ # # ]: 0 : do n1=1,2**op%flavors
3452 : 0 : write(op%ostream,'(i4,10i2,f20.2)') n1, (occ(n1,n2),n2=1,op%flavors),op%occupconfig(n1)
3453 : 0 : sumh=sumh+op%occupconfig(n1)
3454 : : enddo
3455 : : end if
3456 : 0 : write(op%ostream,'(a,f10.4)') " all" , sumh
3457 : :
3458 : 0 : sumtot=0
3459 [ # # ]: 0 : do nelec=0,10
3460 : 0 : spinmin=modulo(nelec,2)
3461 : : if(nelec<=5) spinmax=nelec
3462 : : if(nelec>=6) spinmax=10-nelec
3463 : 0 : spinmax=6
3464 : 0 : dspin=1
3465 : 0 : write(op%ostream,*) "=== Histogram of occupations of configurations for total number of electrons",nelec
3466 [ # # ]: 0 : do spin=spinmin,spinmax,dspin
3467 : 0 : sumh=0
3468 [ # # # # : 0 : do n1=1,2**op%flavors
# # ]
3469 [ # # # # ]: 0 : if(occtot(n1)==nelec.and.abs(spintot(n1))==spin) then
3470 : 0 : sumh=sumh+op%occupconfig(n1)
3471 [ # # ]: 0 : if(op%flavors==14) then
3472 : 0 : write(op%ostream,'(i8,14i2,a,i2,i3,f10.4)') n1,(occ(n1,n2),n2=1,op%flavors)," ",occtot(n1),spintot(n1),&
3473 : 0 : & op%occupconfig(n1)
3474 [ # # ]: 0 : else if(op%flavors==10) then
3475 : 0 : write(op%ostream,'(i8,10i2,a,i2,i3,f10.4)') n1,(occ(n1,n2),n2=1,op%flavors)," ",occtot(n1),spintot(n1),&
3476 : 0 : op%occupconfig(n1)
3477 : : end if
3478 : : endif
3479 : : enddo
3480 : 0 : write(op%ostream,'(a,i4,a,i4,a,f10.4)') " === Sum of weights for",nelec," electrons and spin",spin," is ",sumh
3481 [ # # ]: 0 : sumtot=sumtot+sumh
3482 : : enddo
3483 : : enddo
3484 : 0 : write(op%ostream,'(a,f10.4)') "Full sum is",sumtot
3485 : :
3486 [ # # ]: 0 : FREE(occ)
3487 [ # # ]: 0 : FREE(occtot)
3488 [ # # ]: 0 : FREE(spintot)
3489 : : !END IF
3490 : :
3491 : : !==============================
3492 : : ! Write Susceptibilities
3493 : : !==============================
3494 [ # # ]: 0 : if(Iatom .lt. 10) then
3495 : 0 : write(atomnb, '("0",i1)') Iatom
3496 : : else
3497 : 0 : write(atomnb, '(i2)') Iatom
3498 : : end if
3499 : : ! Local Magnetic Susceptibility
3500 [ # # ]: 0 : if(op%opt_histo .gt. 1) then
3501 : : ! Scalar
3502 [ # # ]: 0 : if(op%nspinor .eq. 1) then
3503 : 0 : open (unit=735,file=trim(fname)//'_LocalSpinSuscept_atom_'//atomnb//'.dat',status='unknown',form='formatted')
3504 : 0 : write(735,*) '#Tau Total t2g eg'
3505 [ # # ]: 0 : do n1=1,op%samples
3506 [ # # ]: 0 : op%suscep(:,n1)=op%suscep(:,n1)/float(nbprocs)/float(op%samples)
3507 [ # # # # ]: 0 : write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (n1-1)*op%beta/op%samples,(op%suscep(n2,n1),n2=1,3)
3508 : : enddo
3509 : : !add tau=beta
3510 [ # # # # ]: 0 : write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (op%samples)*op%beta/op%samples,(op%suscep(n2,1),n2=1,3)
3511 : :
3512 : : else
3513 : : ! SOC
3514 : 0 : open (unit=735,file=trim(fname)//'_LocalMagnSuscept_atom_'//atomnb//'.dat',status='unknown',form='formatted')
3515 : 0 : write(735,*) '#Tau Total Orbital Spin'
3516 [ # # ]: 0 : do n1=1,op%samples
3517 [ # # ]: 0 : op%chi(:,n1) = op%chi(:,n1)/float(nbprocs)/float(op%samples)
3518 [ # # # # ]: 0 : write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (n1-1)*op%beta/op%samples,(op%chi(n2,n1),n2=1,3)
3519 : : end do
3520 : : !add tau=beta
3521 [ # # # # ]: 0 : write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (op%samples)*op%beta/op%samples,(op%chi(n2,1),n2=1,3)
3522 : : endif
3523 : 0 : close(unit=735)
3524 : : endif
3525 : :
3526 : : ! Local Charge Susceptiblity
3527 [ # # ]: 0 : if(op%opt_histo .gt. 2) then
3528 [ # # ]: 0 : op%ntot(:) = op%ntot(:)/float(nbprocs)/float(op%samples)
3529 : 0 : open (unit=735,file=trim(fname)//'_LocalChargeSuscept_atom_'//atomnb//'.dat',status='unknown',form='formatted')
3530 : 0 : write(735,*) '#Tau Total <ntot>'
3531 [ # # ]: 0 : do n1=1,op%samples
3532 : 0 : op%chicharge(1,n1)=(op%chicharge(1,n1)/float(nbprocs)/float(op%samples))-(op%ntot(1)*op%ntot(1))
3533 : : !op%chicharge(2,n1)=(op%chicharge(2,n1)/float(nbprocs)/float(op%samples))-(op%ntot(2)*op%ntot(2))
3534 : : !op%chicharge(3,n1)=(op%chicharge(3,n1)/float(nbprocs)/float(op%samples))-(op%ntot(3)*op%ntot(3))
3535 : 0 : write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (n1-1)*op%beta/op%samples,(op%chicharge(1,n1)),op%ntot(1)
3536 : : enddo
3537 : : !add tau=beta
3538 : 0 : write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (op%samples)*op%beta/op%samples,(op%chicharge(1,1)),op%ntot(1)
3539 : 0 : close(unit=735)
3540 : : endif
3541 : :
3542 : : END IF
3543 : :
3544 : :
3545 [ # # # # : 0 : END SUBROUTINE CtqmcoffdiagComplex_getResult
# # ]
3546 : : !!***
3547 : :
3548 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_symmetrizeGreen
3549 : : !! NAME
3550 : : !! CtqmcoffdiagComplex_symmetrizeGreen
3551 : : !!
3552 : : !! FUNCTION
3553 : : !! optionnaly symmetrize the green functions
3554 : : !!
3555 : : !! COPYRIGHT
3556 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
3557 : : !! This file is distributed under the terms of the
3558 : : !! GNU General Public License, see ~abinit/COPYING
3559 : : !! or http://www.gnu.org/copyleft/gpl.txt .
3560 : : !!
3561 : : !! INPUTS
3562 : : !! op=ctqmc
3563 : : !! syms=weight factors
3564 : : !!
3565 : : !! OUTPUT
3566 : : !! argout(sizeout)=description
3567 : : !!
3568 : : !! SIDE EFFECTS
3569 : : !!
3570 : : !! NOTES
3571 : : !!
3572 : : !! SOURCE
3573 : :
3574 : 0 : SUBROUTINE CtqmcoffdiagComplex_symmetrizeGreen(op, syms)
3575 : :
3576 : : !Arguments ------------------------------------
3577 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
3578 : : DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN ) :: syms
3579 : : !Local variables ------------------------------
3580 : : !INTEGER :: iflavor1
3581 : : !INTEGER :: iflavor2
3582 : : !INTEGER :: flavors
3583 : : !DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: green_tmp
3584 : : !DOUBLE PRECISION, ALLOCATABLE, DIMENSION(: ) :: n_tmp
3585 : :
3586 : : ABI_UNUSED((/syms(1,1), op%swap/))
3587 : :
3588 : : ! flavors = op%flavors
3589 : : ! IF ( SIZE(syms,1) .NE. flavors .OR. SIZE(syms,2) .NE. flavors ) THEN
3590 : : ! CALL WARNALL("CtqmcoffdiagComplex_symmetrizeGreen : wrong opt_sym -> not symmetrizing")
3591 : : ! RETURN
3592 : : ! END IF
3593 : : !
3594 : : ! MALLOC(green_tmp,(1:op%samples+1,flavors))
3595 : : ! green_tmp(:,:) = 0.d0
3596 : : ! MALLOC(n_tmp,(1:flavors))
3597 : : ! n_tmp(:) = 0.d0
3598 : : ! DO iflavor1=1, flavors
3599 : : ! DO iflavor2=1,flavors
3600 : : ! green_tmp(:,iflavor1) = green_tmp(:,iflavor1) &
3601 : : ! + syms(iflavor2,iflavor1) * op%Greens(iflavor2)%oper(:)
3602 : : ! n_tmp(iflavor1) = n_tmp(iflavor1) &
3603 : : ! + syms(iflavor2,iflavor1) * op%measN(1,iflavor2)
3604 : : ! END DO
3605 : : ! END DO
3606 : : ! DO iflavor1=1, flavors
3607 : : ! op%Greens(iflavor1)%oper(:) = green_tmp(:,iflavor1)
3608 : : ! op%measN(1,iflavor1) = n_tmp(iflavor1)
3609 : : ! END DO
3610 : : ! FREE(green_tmp)
3611 : : ! FREE(n_tmp)
3612 : 0 : END SUBROUTINE CtqmcoffdiagComplex_symmetrizeGreen
3613 : : !!***
3614 : :
3615 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_getGreen
3616 : : !! NAME
3617 : : !! CtqmcoffdiagComplex_getGreen
3618 : : !!
3619 : : !! FUNCTION
3620 : : !! Get the full green functions in time and/or frequency
3621 : : !!
3622 : : !! COPYRIGHT
3623 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
3624 : : !! This file is distributed under the terms of the
3625 : : !! GNU General Public License, see ~abinit/COPYING
3626 : : !! or http://www.gnu.org/copyleft/gpl.txt .
3627 : : !!
3628 : : !! INPUTS
3629 : : !! op=ctqmc
3630 : : !!
3631 : : !! OUTPUT
3632 : : !! Gtau=green function in time
3633 : : !! Gw=green function in frequency
3634 : : !!
3635 : : !! SIDE EFFECTS
3636 : : !!
3637 : : !! NOTES
3638 : : !!
3639 : : !! SOURCE
3640 : :
3641 [ # # # # : 0 : SUBROUTINE CtqmcoffdiagComplex_getGreen(op, Gtau, Gw)
# # # # #
# # # ]
3642 : :
3643 : : !Arguments ------------------------------------
3644 : : USE m_GreenHyboffdiagComplex
3645 : : TYPE(CtqmcoffdiagComplex) , INTENT(INOUT) :: op
3646 : : COMPLEX(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: Gtau
3647 : : COMPLEX(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: Gw
3648 : : !Local variables ------------------------------
3649 : : !INTEGER :: itime
3650 : : INTEGER :: iflavor1
3651 : : INTEGER :: iflavor1b !,iflavor,iflavorbis
3652 : : INTEGER :: iflavor2
3653 : : INTEGER :: iflavor3
3654 : : INTEGER :: flavors,tail
3655 : : INTEGER :: ifreq,itime
3656 : : COMPLEX(KIND=8) :: u1
3657 : : COMPLEX(KIND=8) :: u2
3658 : : COMPLEX(KIND=8) :: u3
3659 : : COMPLEX(KIND=8) :: Un
3660 : : COMPLEX(KIND=8) :: UUnn
3661 : : DOUBLE PRECISION :: iw !omega,
3662 : : CHARACTER(LEN=4) :: cflavors
3663 : : CHARACTER(LEN=50) :: string
3664 : 0 : TYPE(GreenHyboffdiagComplex) :: F_tmp
3665 : :
3666 : 0 : flavors = op%flavors
3667 [ # # ]: 0 : DO iflavor1 = 1, flavors
3668 : 0 : u1 = cmplx(0.d0,0.d0,kind=8)
3669 : 0 : u2 = cmplx(0.d0,0.d0,kind=8)
3670 : 0 : u3 = cmplx(0.d0,0.d0,kind=8)
3671 [ # # ]: 0 : DO iflavor2 = 1, flavors
3672 [ # # ]: 0 : IF ( iflavor2 .EQ. iflavor1 ) CYCLE
3673 : 0 : Un = op%Impurity%mat_U(iflavor2,iflavor1) * op%measN(1,iflavor2)
3674 : : ! Un = op%Impurity%mat_U(iflavor2,iflavor1) * (op%Greens%oper(1,iflavor2,iflavor2) + 1.d0)
3675 : : !write(6,*) "forsetmoments",iflavor1,iflavor2,(op%Greens%oper(1,iflavor2,iflavor2) + 1.d0), Un
3676 : 0 : u1 = u1 + Un
3677 : 0 : u2 = u2 + Un*op%Impurity%mat_U(iflavor2,iflavor1)
3678 [ # # ]: 0 : DO iflavor3 = 1, flavors
3679 [ # # ]: 0 : IF ( iflavor3 .EQ. iflavor2 .OR. iflavor3 .EQ. iflavor1 ) CYCLE
3680 : : UUnn = (op%Impurity%mat_U(iflavor2,iflavor1)*op%Impurity%mat_U(iflavor3,iflavor1)) * &
3681 : 0 : & op%measDE(iflavor2,iflavor3)
3682 : 0 : u2 = u2 + UUnn
3683 : : END DO
3684 : : END DO
3685 : : ! write(6,*) "u1,u2",u1,u2
3686 : :
3687 [ # # ]: 0 : DO iflavor1b = 1, flavors
3688 : 0 : u3 =-(op%Impurity%mat_U(iflavor1,iflavor1b))*op%Greens%oper(1,iflavor1,iflavor1b)
3689 : : ! u3=U_{1,1b}*G_{1,1b}
3690 : 0 : CALL GreenHyboffdiagComplex_setMoments(op%Greens,iflavor1,iflavor1b,u1,u2,u3)
3691 : : END DO ! iflavor1b
3692 : :
3693 : : END DO ! iflavor1
3694 : :
3695 [ # # ]: 0 : IF ( PRESENT( Gtau ) ) THEN
3696 [ # # ]: 0 : DO iflavor1 = 1, flavors
3697 [ # # ]: 0 : DO iflavor2 = 1, flavors
3698 [ # # ]: 0 : Gtau(1:op%samples,iflavor1,iflavor2) = op%Greens%oper(1:op%samples,iflavor1,iflavor2)
3699 : : END DO
3700 : : END DO ! iflavor1
3701 : : END IF
3702 : : ! !--------- Write Occupation matrix before Gtau
3703 : : ! write(ostream,'(17x,a)') "Occupation matrix"
3704 : : ! write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
3705 : : ! write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
3706 : : ! do iflavor=1, op%flavors
3707 : : ! write(ostream,'(7x,i10,a,30f10.4)') iflavor,"|",(-op%Greens%oper(op%samples,iflavor,iflavorbis),iflavorbis=1,op%flavors)
3708 : : ! enddo
3709 : : ! write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors),ch10
3710 : : ! !------------------------------------------------------------------------------------------
3711 : : ! !--------- Write Occupation matrix Gtau
3712 : : ! write(ostream,'(17x,a)') "Occupation matrix"
3713 : : ! write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
3714 : : ! write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
3715 : : ! do iflavor=1, op%flavors
3716 : : ! write(ostream,'(7x,i10,a,30f10.4)') iflavor,"|",(Gtau(op%samples,iflavor,iflavorbis),iflavorbis=1,op%flavors)
3717 : : ! enddo
3718 : : ! write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors),ch10
3719 : : ! !------------------------------------------------------------------------------------------
3720 : :
3721 : : !================================================
3722 : : if(3==4) then
3723 : : !================================================
3724 : : DO iflavor1 = 1, flavors
3725 : : DO iflavor1b = 1, flavors
3726 : : !call nfourier3(op%Greens%oper(1:op%samples,iflavor1,iflavor1b),Gw(1:op%samples,iflavor1,iflavor1b),iflavor1==iflavor1b,op%Greens%samples,op%Greens%samples-1,op%Greens%beta,1.d0,op%Greens%Mk(iflavor1,iflavor1b,1),op%Greens%Mk(iflavor1,iflavor1b,2),op%Greens%Mk(iflavor1,iflavor1b,3))
3727 : : END DO
3728 : : END DO ! iflavor1
3729 : : ! ============== write Gomega_nd.dat
3730 : : if(op%rank==0) then
3731 : : OPEN(UNIT=44, FILE="Gomega_nd_nfourier2.dat")
3732 : : WRITE(cflavors,'(I4)') 2*(flavors*flavors+1)
3733 : : string = '(1x,'//TRIM(ADJUSTL(cflavors))//'E15.5)'
3734 : : !write(6,*) " op%Greens%Wmax", op%Greens%Wmax
3735 : : do iflavor1=1, flavors
3736 : : do iflavor1b=1, flavors
3737 : : write(44,*) "#op%Greens%Mk(iflavor1,iflavor2,1",op%Greens%Mk(iflavor1,iflavor1b,:)
3738 : : DO ifreq = 1, op%samples
3739 : : ! !write(6,string) (DBLE(ifreq)*2-1)*3.1415/op%Greens%beta, &
3740 : : ! (/ ((real(Gw(ifreq,iflavor1,iflavor1b)),imag(Gw(ifreq,iflavor1,iflavor1b)), iflavor1=1, flavors),iflavor1b=1,flavors) /)
3741 : : ! WRITE(44,string) (DBLE(ifreq)*2.d0-1.d0)*3.1415926/op%Greens%beta, &
3742 : : iw=aimag(Gw(ifreq,op%flavors,op%flavors+1))
3743 : : WRITE(44,string) aimag(Gw(ifreq,op%flavors,op%flavors+1)),&
3744 : : real(Gw(ifreq,iflavor1,iflavor1b)),aimag(Gw(ifreq,iflavor1,iflavor1b)),&
3745 : : ( -op%Greens%Mk(iflavor1,iflavor1b,2) )/(iw*iw) , (op%Greens%Mk(iflavor1,iflavor1b,1))/iw!-op%Greens%Mk(iflavor1,iflavor1b,3)/(iw*iw))/iw
3746 : : ! WRITE(102,*) aimag(Gw(ifreq,op%flavors,op%flavors+1)), (op%Greens%Mk(iflavor1,iflavor1b,1))/iw,op%Greens%Mk(iflavor1,iflavor1b,1),iw
3747 : : END DO
3748 : : WRITE(44,*)
3749 : : END DO
3750 : : END DO
3751 : : close(44)
3752 : : endif
3753 : : !================================================
3754 : : endif
3755 : : !================================================
3756 : : !!write(6,*) "present gw", present(gw)
3757 [ # # ]: 0 : IF ( PRESENT( Gw ) ) THEN
3758 : : !!write(6,*) "size gw",SIZE(Gw,DIM=2) ,flavors+1
3759 [ # # ]: 0 : IF ( SIZE(Gw,DIM=3) .EQ. flavors+1 ) THEN
3760 : : ! CALL GreenHyboffdiag_forFourier(op%Greens, Gomega=Gw, omega=Gw(:,op%flavors,op%flavors+1))
3761 : 0 : CALL GreenHyboffdiagComplex_forFourierComplex(op%Greens, Gomega=Gw, omega=Gw(:,op%flavors,op%flavors+1))
3762 : : !write(6,*) "1"
3763 : : !IF ( op%rank .EQ. 0 ) write(20,*) Gw(:,iflavor1)
3764 [ # # ]: 0 : ELSE IF ( SIZE(Gw,DIM=3) .EQ. flavors ) THEN
3765 : 0 : CALL GreenHyboffdiagComplex_forFourierComplex(op%Greens,Gomega=Gw)
3766 : : !write(6,*) "2"
3767 : : ELSE
3768 : 0 : CALL WARNALL("CtqmcoffdiagComplex_getGreen : Gw is not valid ")
3769 : 0 : CALL GreenHyboffdiagComplex_forFourierComplex(op%Greens,Wmax=op%Wmax)
3770 : : !write(6,*) "3"
3771 : : END IF
3772 : : ELSE
3773 : 0 : CALL GreenHyboffdiagComplex_forFourierComplex(op%Greens,Wmax=op%Wmax)
3774 : : END IF
3775 : : ! ============== write Gomega_nd.dat
3776 : : !================================================
3777 : : ! if(3==4) then
3778 : : !================================================
3779 : : if(op%rank==0.and.3==4) then
3780 : : OPEN(UNIT=44, FILE="Gomega_nd.dat")
3781 : : WRITE(cflavors,'(I4)') 2*(flavors*flavors+1)
3782 : : string = '(1x,'//TRIM(ADJUSTL(cflavors))//'E15.5)'
3783 : : !write(6,*) " op%Greens%Wmax", op%Greens%Wmax
3784 : : do iflavor1=1, flavors
3785 : : do iflavor1b=1, flavors
3786 : : write(44,*) "#op%Greens%Mk(iflavor1,iflavor2,1",op%Greens%Mk(iflavor1,iflavor1b,:)
3787 : : DO ifreq = 1, SIZE(Gw,1)
3788 : : ! !write(6,string) (DBLE(ifreq)*2-1)*3.1415/op%Greens%beta, &
3789 : : ! (/ ((real(Gw(ifreq,iflavor1,iflavor1b)),imag(Gw(ifreq,iflavor1,iflavor1b)), iflavor1=1, flavors),iflavor1b=1,flavors) /)
3790 : : ! WRITE(44,string) (DBLE(ifreq)*2.d0-1.d0)*3.1415926/op%Greens%beta, &
3791 : : iw=aimag(Gw(ifreq,op%flavors,op%flavors+1))
3792 : : WRITE(44,string) aimag(Gw(ifreq,op%flavors,op%flavors+1)),&
3793 : : real(Gw(ifreq,iflavor1,iflavor1b)),aimag(Gw(ifreq,iflavor1,iflavor1b)),&
3794 : : ( -op%Greens%Mk(iflavor1,iflavor1b,2) )/(iw*iw) , &
3795 : : & (op%Greens%Mk(iflavor1,iflavor1b,1)-op%Greens%Mk(iflavor1,iflavor1b,3)/(iw*iw))/iw
3796 : : END DO
3797 : : WRITE(44,*)
3798 : : END DO
3799 : : END DO
3800 : : endif
3801 : : !================================================
3802 : : ! endif
3803 : : !================================================
3804 : :
3805 : :
3806 : : ! ==============================
3807 : : ! --- Initialize F_tmp
3808 : : !write(6,*) "10"
3809 : :
3810 [ # # ]: 0 : IF ( op%have_MPI .EQV. .TRUE. ) THEN
3811 : : !CALL GreenHyboffdiag_init(F_tmp,op%samples,op%beta,op%flavors,MY_COMM=op%MY_COMM)
3812 : 0 : CALL GreenHyboffdiagComplex_init(F_tmp,op%samples,op%beta,flavors)
3813 : : !write(6,*) "10a"
3814 : : ELSE
3815 : 0 : CALL GreenHyboffdiagComplex_init(F_tmp,op%samples,op%beta,flavors)
3816 : : !write(6,*) "10b"
3817 : : END IF
3818 : :
3819 : : !write(6,*) "11"
3820 : : ! CALL GreenHyboffdiag_setOperW(F_tmp,Gw)
3821 : :
3822 : 0 : tail = op%samples
3823 : 0 : F_tmp%Wmax=op%samples ! backFourier only works for linear freq: calculation of A and etc..
3824 [ # # # # : 0 : MALLOC(F_tmp%oper_w,(1:tail,op%flavors,op%flavors))
# # # # #
# # # # #
# # # # #
# # # ]
3825 [ # # # # : 0 : F_tmp%oper_w(1:tail,1:F_tmp%nflavors,1:F_tmp%nflavors) = Gw(1:tail,1:F_tmp%nflavors,1:F_tmp%nflavors)
# # ]
3826 : : !write(6,*) "example",F_tmp%oper_w(1,1,1)
3827 : : !write(6,*) "example",Gw(1,1,1)
3828 : 0 : F_tmp%setW = .TRUE.
3829 : : !write(6,*) size(F_tmp%oper_w,1)
3830 : : !write(6,*) size(F_tmp%oper_w,2)
3831 : : !write(6,*) size(F_tmp%oper_w,3)
3832 : : !write(6,*) size(Gw,1)
3833 : : !write(6,*) size(Gw,2)
3834 : : !write(6,*) size(Gw,3)
3835 : :
3836 : : !write(6,*) "eee", (2.d0*DBLE(ifreq)-1.d0) * 3.1415/op%beta,real(F_tmp%oper_w(1,1,1)),imag(F_tmp%oper_w(1,1,1))
3837 : : !================================================
3838 : : if(3==4) then
3839 : : !================================================
3840 : : OPEN(UNIT=3337, FILE="Gomega_nd2.dat")
3841 : : do iflavor1=1, flavors
3842 : : do iflavor1b=1, flavors
3843 : : do ifreq=1, tail
3844 : : ! write(3337,*) (2.d0*DBLE(ifreq)-1.d0) * 3.1415/op%beta,real(F_tmp%oper_w(ifreq,iflavor1,iflavor1b)),&
3845 : : ! & imag(F_tmp%oper_w(ifreq,iflavor1,iflavor1b))
3846 : : write(3337,*) aimag(Gw(ifreq,op%flavors,op%flavors+1)), real(F_tmp%oper_w(ifreq,iflavor1,iflavor1b)),&
3847 : : & aimag(F_tmp%oper_w(ifreq,iflavor1,iflavor1b))
3848 : :
3849 : : ! omega=(2.d0*DBLE(ifreq)-1.d0) * 3.1415/op%beta
3850 : : ! F_tmp%oper_w(ifreq,iflavor1,iflavor1b)=0.1**2/Gw(ifreq,op%flavors,op%flavors+1)
3851 : : enddo
3852 : : write(3337,*)
3853 : : enddo
3854 : : enddo
3855 : : close(3337)
3856 : : !================================================
3857 : : endif
3858 : : !================================================
3859 : :
3860 : : !write(6,*) "12",F_tmp%Wmax
3861 : :
3862 : : ! CALL GreenHyboffdiag_backFourier(F_tmp,func="green")
3863 : :
3864 : : !write(6,*) "13"
3865 : :
3866 : : !================================================
3867 : : if(3==4) then
3868 : : !================================================
3869 : : OPEN(UNIT=48, FILE="Gtau_nd_2.dat")
3870 : : ! --- Print full non diagonal Gtau in Gtau_nd.dat
3871 : : WRITE(cflavors,'(I4)') flavors*flavors+1
3872 : : string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)'
3873 : : DO itime = 1, op%samples+1
3874 : : WRITE(48,string) DBLE(itime-1)*op%beta/DBLE(op%samples), &
3875 : : & ((F_tmp%oper(itime,iflavor1,iflavor1b), iflavor1=1, flavors),iflavor1b=1,flavors)
3876 : : END DO
3877 : : !================================================
3878 : : endif
3879 : : !================================================
3880 : :
3881 : 0 : CALL GreenHyboffdiagComplex_destroy(F_tmp)
3882 : :
3883 : : !FREE(F_tmp%oper_w)
3884 : : ! ==============================
3885 [ # # # # : 0 : END SUBROUTINE CtqmcoffdiagComplex_getGreen
# # # # #
# # # # #
# # # # #
# # # ]
3886 : : !!***
3887 : :
3888 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_getD
3889 : : !! NAME
3890 : : !! CtqmcoffdiagComplex_getD
3891 : : !!
3892 : : !! FUNCTION
3893 : : !! get double occupation
3894 : : !!
3895 : : !! COPYRIGHT
3896 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
3897 : : !! This file is distributed under the terms of the
3898 : : !! GNU General Public License, see ~abinit/COPYING
3899 : : !! or http://www.gnu.org/copyleft/gpl.txt .
3900 : : !!
3901 : : !! INPUTS
3902 : : !! op=ctqmc
3903 : : !!
3904 : : !! OUTPUT
3905 : : !! D=full double occupation
3906 : : !!
3907 : : !! SIDE EFFECTS
3908 : : !!
3909 : : !! NOTES
3910 : : !!
3911 : : !! SOURCE
3912 : :
3913 : 0 : SUBROUTINE CtqmcoffdiagComplex_getD(op, D)
3914 : :
3915 : : !Arguments ------------------------------------
3916 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN ) :: op
3917 : : DOUBLE PRECISION, INTENT(OUT) :: D
3918 : : !Local variables ------------------------------
3919 : : INTEGER :: iflavor1
3920 : : INTEGER :: iflavor2
3921 : :
3922 : 0 : D = 0.d0
3923 : :
3924 [ # # ]: 0 : DO iflavor1 = 1, op%flavors
3925 [ # # ]: 0 : DO iflavor2 = iflavor1+1, op%flavors
3926 : 0 : D = D + op%measDE(iflavor2,iflavor1)
3927 : : END DO
3928 : : END DO
3929 : : !IF ( op%rank .EQ. 0 ) THEN
3930 : : ! DO iflavor1 = 1, op%flavors
3931 : : ! DO iflavor2 = iflavor1+1, op%flavors
3932 : : ! write(4533,*) op%measDE(iflavor2,iflavor1)k
3933 : : ! write(4534,*) op%Impurity%mat_U(iflavor2,iflavor1)k
3934 : : ! END DO
3935 : : ! END DO
3936 : :
3937 : : !ENDIF
3938 : :
3939 : 0 : END SUBROUTINE CtqmcoffdiagComplex_getD
3940 : : !!***
3941 : :
3942 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_getE
3943 : : !! NAME
3944 : : !! CtqmcoffdiagComplex_getE
3945 : : !!
3946 : : !! FUNCTION
3947 : : !! get interaction energy and noise on it
3948 : : !!
3949 : : !! COPYRIGHT
3950 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
3951 : : !! This file is distributed under the terms of the
3952 : : !! GNU General Public License, see ~abinit/COPYING
3953 : : !! or http://www.gnu.org/copyleft/gpl.txt .
3954 : : !!
3955 : : !! INPUTS
3956 : : !! op=ctqmc
3957 : : !!
3958 : : !! OUTPUT
3959 : : !! E=interaction energy
3960 : : !! noise=noise on this value
3961 : : !!
3962 : : !! SIDE EFFECTS
3963 : : !!
3964 : : !! NOTES
3965 : : !!
3966 : : !! SOURCE
3967 : :
3968 : 0 : SUBROUTINE CtqmcoffdiagComplex_getE(op,E,noise)
3969 : :
3970 : : !Arguments ------------------------------------
3971 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN ) :: op
3972 : : DOUBLE PRECISION, INTENT(OUT) :: E
3973 : : DOUBLE PRECISION, INTENT(OUT) :: Noise
3974 : :
3975 : 0 : E = op%measDE(1,1)
3976 : 0 : Noise = op%a_Noise*(DBLE(op%sweeps)*DBLE(op%size))**op%b_Noise
3977 : 0 : END SUBROUTINE CtqmcoffdiagComplex_getE
3978 : : !!***
3979 : :
3980 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printAll
3981 : : !! NAME
3982 : : !! CtqmcoffdiagComplex_printAll
3983 : : !!
3984 : : !! FUNCTION
3985 : : !! print different functions computed during the simulation
3986 : : !!
3987 : : !! COPYRIGHT
3988 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
3989 : : !! This file is distributed under the terms of the
3990 : : !! GNU General Public License, see ~abinit/COPYING
3991 : : !! or http://www.gnu.org/copyleft/gpl.txt .
3992 : : !!
3993 : : !! INPUTS
3994 : : !! op=ctqmc
3995 : : !!
3996 : : !! OUTPUT
3997 : : !!
3998 : : !! SIDE EFFECTS
3999 : : !!
4000 : : !! NOTES
4001 : : !!
4002 : : !! SOURCE
4003 : :
4004 : 0 : SUBROUTINE CtqmcoffdiagComplex_printAll(op)
4005 : :
4006 : : !Arguments ------------------------------------
4007 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
4008 : :
4009 [ # # ]: 0 : IF ( .NOT. op%done ) &
4010 : 0 : CALL WARNALL("CtqmcoffdiagComplex_printAll : Simulation not run ")
4011 : :
4012 : : !sui!write(6,*) "op%stats",op%stats
4013 : 0 : CALL CtqmcoffdiagComplex_printQMC(op)
4014 : :
4015 : 0 : CALL CtqmcoffdiagComplex_printGreen(op)
4016 : :
4017 : 0 : CALL CtqmcoffdiagComplex_printD(op)
4018 : :
4019 : : ! CALL CtqmcoffdiagComplex_printE(op)
4020 : :
4021 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
4022 : : ! CALL CtqmcoffdiagComplex_printPerturbation(op)
4023 : :
4024 : 0 : CALL CtqmcoffdiagComplex_printCorrelation(op)
4025 : : !#endif
4026 : :
4027 : 0 : CALL CtqmcoffdiagComplex_printSpectra(op)
4028 : :
4029 : 0 : END SUBROUTINE CtqmcoffdiagComplex_printAll
4030 : : !!***
4031 : :
4032 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printQMC
4033 : : !! NAME
4034 : : !! CtqmcoffdiagComplex_printQMC
4035 : : !!
4036 : : !! FUNCTION
4037 : : !! print ctqmc statistics
4038 : : !!
4039 : : !! COPYRIGHT
4040 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4041 : : !! This file is distributed under the terms of the
4042 : : !! GNU General Public License, see ~abinit/COPYING
4043 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4044 : : !!
4045 : : !! INPUTS
4046 : : !! op=ctqmc
4047 : : !!
4048 : : !! OUTPUT
4049 : : !!
4050 : : !! SIDE EFFECTS
4051 : : !!
4052 : : !! NOTES
4053 : : !!
4054 : : !! SOURCE
4055 : :
4056 : 0 : SUBROUTINE CtqmcoffdiagComplex_printQMC(op)
4057 : :
4058 : : !Arguments ------------------------------------
4059 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
4060 : : !Local variables ------------------------------
4061 : : INTEGER :: ostream
4062 : : INTEGER :: iflavor,iflavorbis,iorder
4063 : : DOUBLE PRECISION :: sweeps
4064 : : DOUBLE PRECISION :: invSweeps
4065 : : CHARACTER(LEN=2) :: a
4066 : : CHARACTER(LEN=15) :: string
4067 : :
4068 : : !IF ( op%rank .NE. 0) RETURN
4069 [ # # ]: 0 : IF ( op%rank .NE. MOD(op%size,op%size)) RETURN
4070 : :
4071 : 0 : ostream = op%ostream
4072 : 0 : sweeps = DBLE(op%sweeps)
4073 : 0 : invSweeps = 1.d0/sweeps
4074 : :
4075 : 0 : WRITE(ostream,'(a,F13.0,A11,F10.2,A12,I5,A5)') "The QMC Solver did ",sweeps*DBLE(op%size), " sweeps in ", op%runTime, &
4076 : 0 : " seconds on ", op%size, " CPUs"
4077 : 0 : WRITE(ostream,'(A28,F6.2)') "Segments added [%] : ", op%stats(4)*invSweeps*100.d0
4078 : 0 : WRITE(ostream,'(A28,F6.2)') "Segments removed [%] : ", op%stats(5)*invSweeps*100.d0
4079 : 0 : WRITE(ostream,'(A28,F6.2)') "Segments <0 sign [%] : ", op%stats(6)*invSweeps*100.d0
4080 : : !WRITE(ostream,'(A28,F12.2)') "Number of meas [%] : ", op%stats(6)
4081 : 0 : WRITE(ostream,'(A28,F6.2)') "Anti-segments added [%] : ", op%stats(1)*invSweeps*100.d0
4082 : 0 : WRITE(ostream,'(A28,F6.2)') "Anti-segments removed [%] : ", op%stats(2)*invSweeps*100.d0
4083 : 0 : WRITE(ostream,'(A28,F6.2)') "Anti-segments <0 sign [%] : ", op%stats(3)*invSweeps*100.d0
4084 : : !WRITE(ostream,'(A28,F12.2)') "Sum of sign [%] : ", op%stats(3)
4085 : 0 : WRITE(ostream,'(A28,F13.2)') "Signe value : ", op%Greens%signvaluemeas
4086 : 0 : WRITE(ostream,'(A28,F13.2,F13.2)') "Phase value : ", op%Greens%phasevaluemeas
4087 [ # # ]: 0 : IF ( op%modGlobalMove(1) .LT. op%sweeps + 1 ) THEN
4088 : 0 : WRITE(ostream,'(A28,F6.2)') "Global Move [%] : ", op%swap *invSweeps*100.d0*op%modGlobalMove(1)
4089 : 0 : WRITE(ostream,'(A28,F6.2)') "Global Move Reduced [%] : ", op%swap / DBLE(op%modGlobalMove(2))*100.d0
4090 : : END IF
4091 : : !#ifdef CTCtqmcoffdiagComplex_CHECK
4092 [ # # ]: 0 : IF ( op%opt_check .EQ. 1 .OR. op%opt_check .EQ. 3 ) &
4093 : 0 : WRITE(ostream,'(A28,E22.14)') "Impurity test [%] : ", op%errorImpurity*100.d0
4094 [ # # ]: 0 : IF ( op%opt_check .GE. 2 ) &
4095 : 0 : WRITE(ostream,'(A28,E22.14)') "Bath test [%] : ", op%errorBath *100.d0
4096 : : !#endif
4097 : 0 : WRITE(ostream,'(A28,ES22.14,A5,ES21.14)') "<Epot> [U] : ", op%measDE(1,1), " +/- ",&
4098 : : !#ifdef HAVE_MPI
4099 : 0 : op%a_Noise*(sweeps*DBLE(op%size))**op%b_Noise
4100 : : !#else
4101 : : ! op%a_Noise*(sweeps)**op%b_Noise
4102 : : !#endif
4103 : : !--------- Write double occupation between all pairs of orbitals --------------------------
4104 : 0 : write(ostream,'(17x,a)') "Double occupation between pairs of orbitals"
4105 [ # # # # ]: 0 : write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
4106 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
4107 [ # # ]: 0 : do iflavor=1, op%flavors
4108 [ # # # # ]: 0 : write(ostream,'(7x,i10,a,30f10.4)') iflavor,"|",(op%measDE(iflavor,iflavorbis),iflavorbis=1,op%flavors)
4109 : : enddo
4110 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors),ch10
4111 : : !------------------------------------------------------------------------------------------
4112 : :
4113 : : !--------- Write number of segments for each orbitals
4114 : : ! write(ostream,'(a)') "Number of segments for each orbitals"
4115 : : ! write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
4116 : : ! write(ostream,'(17x,30a)') ("----------",iflavorbis=1,op%flavors)
4117 : : ! do iflavor=1, op%flavors
4118 : : ! write(ostream,'(i17,a,30f10.4)') iflavor,"|",(op%Impurity%particles(IT)%tail)
4119 : : ! enddo
4120 : : ! write(ostream,'(17x,30a)') ("----------",iflavorbis=1,op%flavors)
4121 : : !------------------------------------------------------------------------------------------
4122 : : !--------- Write G(L)
4123 : 0 : write(ostream,'(17x,a)') "G(tau = L)"
4124 [ # # # # ]: 0 : write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
4125 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
4126 [ # # ]: 0 : do iflavor=1, op%flavors
4127 [ # # # # ]: 0 : write(ostream,'(7x,i10,a,30f10.4)') iflavor,"|",(real(op%Greens%oper(op%samples,iflavor,iflavorbis)),iflavorbis=1,op%flavors)
4128 : : enddo
4129 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors),ch10
4130 : : !------------------------------------------------------------------------------------------
4131 : : !--------- Write G(1)
4132 : 0 : write(ostream,'(17x,a)') "G(tau = 1)"
4133 [ # # # # ]: 0 : write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
4134 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
4135 [ # # ]: 0 : do iflavor=1, op%flavors
4136 [ # # # # ]: 0 : write(ostream,'(7x,i10,a,30f10.4)') iflavor,"|",(real(op%Greens%oper(1,iflavor,iflavorbis)),iflavorbis=1,op%flavors)
4137 : : enddo
4138 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors),ch10
4139 : : !------------------------------------------------------------------------------------------
4140 : :
4141 : 0 : WRITE(ostream,'(A28,F8.4,A3,F7.4)') "Noise [U] : ", op%a_Noise, " x^", op%b_Noise
4142 : 0 : WRITE(ostream,'(A28,E10.2)') "Niquist puls. [/beta] : ", ACOS(-1.d0)*op%inv_dt
4143 : 0 : WRITE(ostream,'(A28,E22.14)') "Max Acc. Epot Error [U] : ", op%measDE(2,2)/(op%beta*op%modNoise1*2.d0)*sweeps
4144 : :
4145 : : !WRITE(ostream,'(A28,F7.4,A3,F7.4,A4,E20.14)') "Noise [G(tau)] : ", op%a_Noise(2), "x^", op%b_Noise(2), " -> ", &
4146 : : !op%a_Noise(2)*(sweeps*DBLE(op%size))**op%b_Noise(2)
4147 : : !----- PERTURBATION ORDER------------------------------------------------------------------
4148 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
4149 : 0 : write(ostream,*)
4150 : 0 : WRITE(a,'(I2)') op%flavors
4151 [ # # # # ]: 0 : string = '(A28,'//TRIM(ADJUSTL(a))//'(1x,I3))'
4152 [ # # # # : 0 : WRITE(ostream,string) "Perturbation orders : ",(/ (MAXLOC(op%measPerturbation(:, iflavor))-1, iflavor=1, op%flavors) /)
# # # # #
# # # # #
# # # # ]
4153 : 0 : write(ostream,'(17x,a)') "order of Perturbation for flavors"
4154 [ # # # # ]: 0 : write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
4155 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
4156 [ # # # # : 0 : write(ostream,'(12x,a,30i10)') " max ",(/ (MAXLOC(op%measPerturbation(:, iflavor))-1, iflavor=1, op%flavors) /)
# # # # #
# # # # #
# # # # ]
4157 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
4158 [ # # ]: 0 : do iorder=0, op%opt_order-1
4159 [ # # # # ]: 0 : write(ostream,'(7x,i10,a,30f10.4)') iorder,"|",(op%measPerturbation(iorder+1,iflavor),iflavor=1,op%flavors)
4160 : : enddo
4161 : : END IF
4162 : : !------------------------------------------------------------------------------------------
4163 : : !----- PERTURBATION ORDER------------------------------------------------------------------
4164 [ # # ]: 0 : IF ( op%opt_order .GT. 0 ) THEN
4165 : 0 : write(ostream,*)
4166 : 0 : write(ostream,'(17x,a)') "Proportion of full and empty orbital for order 0"
4167 [ # # # # ]: 0 : write(ostream,'(17x,30i10)') (iflavorbis,iflavorbis=1,op%flavors)
4168 [ # # # # ]: 0 : write(ostream,'(18x,30a)') ("----------",iflavorbis=1,op%flavors)
4169 [ # # # # ]: 0 : write(ostream,'(2x,a,30f10.4)') " full orbital |",(op%meas_fullemptylines(1,iflavor),iflavor=1,op%flavors)
4170 [ # # # # ]: 0 : write(ostream,'(2x,a,30f10.4)') " empty orbital |",(op%meas_fullemptylines(2,iflavor),iflavor=1,op%flavors)
4171 : : END IF
4172 : : !------------------------------------------------------------------------------------------
4173 : : !CALL FLUSH(op%ostream)
4174 : : IF ( ABS(((op%stats(4) *invSweeps*100.d0) / (op%stats(5) *invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 &
4175 [ # # # # ]: 0 : .OR. ABS(((op%stats(1)*invSweeps*100.d0) / (op%stats(2)*invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 ) &
4176 : : THEN
4177 : 0 : CALL WARNALL("CtqmcoffdiagComplex_printQMC : bad statistic according to moves. Increase sweeps")
4178 : : END IF
4179 [ # # ]: 0 : IF ( ABS(op%b_Noise+0.5)/0.5d0 .GE. 0.05d0 ) &
4180 : 0 : CALL WARNALL("CtqmcoffdiagComplex_printQMC : bad statistic according to Noise. Increase sweeps")
4181 : : ! IF ( ISNAN(op%a_Noise) .OR. ISNAN(op%a_Noise) ) &
4182 : : ! CALL WARNALL("CtqmcoffdiagComplex_printQMC : NaN appeared. Increase sweeps ")
4183 : :
4184 : :
4185 : : END SUBROUTINE CtqmcoffdiagComplex_printQMC
4186 : : !!***
4187 : :
4188 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printGreen
4189 : : !! NAME
4190 : : !! CtqmcoffdiagComplex_printGreen
4191 : : !!
4192 : : !! FUNCTION
4193 : : !! print green functions
4194 : : !!
4195 : : !! COPYRIGHT
4196 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4197 : : !! This file is distributed under the terms of the
4198 : : !! GNU General Public License, see ~abinit/COPYING
4199 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4200 : : !!
4201 : : !! INPUTS
4202 : : !! op=ctqmc
4203 : : !! oFileIn=file stream
4204 : : !!
4205 : : !! OUTPUT
4206 : : !!
4207 : : !! SIDE EFFECTS
4208 : : !!
4209 : : !! NOTES
4210 : : !!
4211 : : !! SOURCE
4212 : :
4213 : 0 : SUBROUTINE CtqmcoffdiagComplex_printGreen(op, oFileIn)
4214 : :
4215 : : !Arguments ------------------------------------
4216 : : use m_io_tools, only : flush_unit
4217 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN) :: op
4218 : : INTEGER , OPTIONAL, INTENT(IN) :: oFileIn
4219 : : !Local variables ------------------------------
4220 : : INTEGER :: oFile
4221 : : INTEGER :: itime
4222 : : INTEGER :: sp1
4223 : : INTEGER :: iflavor,iflavorb
4224 : : INTEGER :: flavors !, iflavor2 !,iflavor1,
4225 : : CHARACTER(LEN=4) :: cflavors,hcflavors
4226 : : CHARACTER(LEN=50) :: string,hstring
4227 : : DOUBLE PRECISION :: dt
4228 : : DOUBLE PRECISION :: sweeps
4229 : :
4230 : : !IF ( op%rank .NE. MOD(1,op%size)) RETURN
4231 [ # # ]: 0 : IF ( op%rank .NE. MOD(op%size+1,op%size)) RETURN
4232 : :
4233 : 0 : oFile = 40
4234 [ # # ]: 0 : IF ( PRESENT(oFileIn) ) THEN
4235 : 0 : oFile = oFileIn
4236 : : ELSE
4237 : 0 : OPEN(UNIT=oFile, FILE="Gtau.dat")
4238 : : END IF
4239 : 0 : OPEN(UNIT=43, FILE="Gtau_nd.dat")
4240 : 0 : rewind(43)
4241 : 0 : OPEN(UNIT=143, FILE="Gtau_nd_Re_ctqmc.dat")
4242 : 0 : rewind(143)
4243 : 0 : OPEN(UNIT=243, FILE="Gtau_nd_Im_ctqmc.dat")
4244 : 0 : rewind(243)
4245 : 0 : sp1 = op%samples
4246 : 0 : dt = op%beta / DBLE(sp1)
4247 : 0 : sp1 = sp1 + 1
4248 : 0 : flavors = op%flavors
4249 : 0 : sweeps = DBLE(op%sweeps)*DBLE(op%size)
4250 : :
4251 [ # # ]: 0 : IF ( op%opt_noise .EQ. 1) THEN
4252 : 0 : WRITE(cflavors,'(I4)') (2*flavors+1)*2
4253 [ # # # # ]: 0 : string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)'
4254 [ # # ]: 0 : DO itime = 1, sp1
4255 : 0 : WRITE(oFile,string) DBLE(itime-1)*dt, &
4256 [ # # # # : 0 : (/ (op%Greens%oper(itime,iflavor,iflavor), iflavor=1, flavors) /), &
# # ]
4257 [ # # # # : 0 : (/ (op%abNoiseG(1,itime,iflavor)*(sweeps)**op%abNoiseG(2,itime,iflavor), iflavor=1, flavors) /)
# # ]
4258 : : END DO
4259 : : ELSE
4260 : 0 : WRITE(cflavors,'(I4)') (flavors+1)*2
4261 [ # # # # ]: 0 : string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)'
4262 [ # # ]: 0 : DO itime = 1, sp1
4263 : : ! WRITE(45,string) DBLE(itime-1)*dt, &
4264 : : ! (/ (op%Greens%oper(itime,iflavor,iflavor), iflavor=1, flavors) /)
4265 : 0 : WRITE(oFile,string) DBLE(itime-1)*dt, &
4266 [ # # # # : 0 : (/ (op%Greens%oper(itime,iflavor,iflavor), iflavor=1, flavors) /)
# # ]
4267 : : ! WRITE(46,*) DBLE(itime-1)*dt, &
4268 : : ! & (/ ((op%Greens%oper(itime,iflavor,iflavorb), iflavor=1, flavors),iflavorb=1,flavors) /)
4269 : : END DO
4270 : : ! DO itime = 1, sp1
4271 : : ! WRITE(47,*) DBLE(itime-1)*dt, &
4272 : : ! & (/ ((op%Greens%oper(itime,iflavor,iflavorb), iflavor=1, flavors),iflavorb=1,flavors) /)
4273 : : ! END DO
4274 : : ! --- Print full non diagonal Gtau in Gtau_nd.dat
4275 : 0 : WRITE(cflavors,'(I4)') (flavors*flavors+1)
4276 : 0 : WRITE(hcflavors,'(I4)') (flavors*flavors+1)
4277 : : ! write(47,*) "cflavors",cflavors
4278 [ # # # # ]: 0 : string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)'
4279 [ # # # # ]: 0 : hstring = '(1x,'//TRIM(ADJUSTL(hcflavors))//'ES22.14)'
4280 : : ! write(47,*) string
4281 [ # # ]: 0 : DO itime = 1, sp1
4282 : 0 : WRITE(43,string) DBLE(itime-1)*dt, &
4283 [ # # # # : 0 : & (/ ((op%Greens%oper(itime,iflavor,iflavorb), iflavorb=1, flavors),iflavor=1,flavors) /)
# # # # ]
4284 : 0 : WRITE(143,hstring) DBLE(itime-1)*dt, &
4285 [ # # # # : 0 : & (/ ((real(op%Greens%oper(itime,iflavor,iflavorb)), iflavorb=1, flavors),iflavor=1,flavors) /)
# # # # ]
4286 : 0 : WRITE(243,hstring) DBLE(itime-1)*dt, &
4287 [ # # # # : 0 : & (/ ((aimag(op%Greens%oper(itime,iflavor,iflavorb)), iflavorb=1, flavors),iflavor=1,flavors) /)
# # # # ]
4288 : : ! WRITE(44,*) DBLE(itime-1)*dt, &
4289 : : ! & (op%Greens%oper(itime,iflavor,iflavor), iflavor=1, flavors)
4290 : : ! WRITE(44,string) DBLE(itime-1)*dt, &
4291 : : ! & (op%Greens%oper(itime,iflavor,iflavor), iflavor=1, flavors)
4292 : : END DO
4293 : 0 : WRITE(43,*)
4294 : 0 : WRITE(143,*)
4295 : 0 : WRITE(243,*)
4296 : : END IF
4297 : : ! DO iflavor = 1, flavors
4298 : : ! DO iflavor2 = 1, flavors
4299 : : ! write(4436,*) "#",iflavor,iflavor2
4300 : : ! do itime=1,sp1
4301 : : ! write(4436,*) DBLE(itime-1)*dt,real(op%Greens%oper(itime,iflavor,iflavor2))
4302 : : ! enddo
4303 : : ! write(4436,*)
4304 : : ! END DO
4305 : : ! END DO
4306 : : ! close(4436)
4307 : :
4308 [ # # ]: 0 : IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
4309 : 0 : CLOSE(43)
4310 : 0 : CLOSE(143)
4311 : 0 : CLOSE(243)
4312 : : ! CLOSE(46)
4313 : : ! CLOSE(47)
4314 : : !call flush_unit(43)
4315 : :
4316 : : END SUBROUTINE CtqmcoffdiagComplex_printGreen
4317 : : !!***
4318 : :
4319 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printD
4320 : : !! NAME
4321 : : !! CtqmcoffdiagComplex_printD
4322 : : !!
4323 : : !! FUNCTION
4324 : : !! print individual double occupancy
4325 : : !!
4326 : : !! COPYRIGHT
4327 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4328 : : !! This file is distributed under the terms of the
4329 : : !! GNU General Public License, see ~abinit/COPYING
4330 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4331 : : !!
4332 : : !! INPUTS
4333 : : !! op=ctqmc
4334 : : !! oFileIn=file stream
4335 : : !!
4336 : : !! OUTPUT
4337 : : !!
4338 : : !! SIDE EFFECTS
4339 : : !!
4340 : : !! NOTES
4341 : : !!
4342 : : !! SOURCE
4343 : :
4344 : 0 : SUBROUTINE CtqmcoffdiagComplex_printD(op,oFileIn)
4345 : :
4346 : : !Arguments ------------------------------------
4347 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN) :: op
4348 : : INTEGER , OPTIONAL, INTENT(IN) :: oFileIn
4349 : : !Local variables ------------------------------
4350 : : INTEGER :: oFile
4351 : : INTEGER :: iflavor1
4352 : : INTEGER :: iflavor2
4353 : :
4354 : : !IF ( op%rank .NE. MOD(2,op%size)) RETURN
4355 [ # # ]: 0 : IF ( op%rank .NE. MOD(op%size+2,op%size)) RETURN
4356 : :
4357 : 0 : oFile = 41
4358 [ # # ]: 0 : IF ( PRESENT(oFileIn) ) THEN
4359 : 0 : oFile = oFileIn
4360 : : ELSE
4361 : 0 : OPEN(UNIT=oFile, FILE="D.dat")
4362 : : END IF
4363 : :
4364 [ # # ]: 0 : DO iflavor1 = 1, op%flavors
4365 [ # # ]: 0 : DO iflavor2 = iflavor1+1, op%flavors
4366 : 0 : WRITE(oFile,'(1x,A8,I4,A1,I4,A3,ES21.14)') "Orbitals", iflavor1, "-", iflavor2, " : ", op%measDE(iflavor2,iflavor1)
4367 : : END DO
4368 : : END DO
4369 : :
4370 [ # # ]: 0 : IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
4371 : :
4372 : : END SUBROUTINE CtqmcoffdiagComplex_printD
4373 : : !!***
4374 : :
4375 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printE
4376 : : !! NAME
4377 : : !! CtqmcoffdiagComplex_printE
4378 : : !!
4379 : : !! FUNCTION
4380 : : !! print energy and noise
4381 : : !!
4382 : : !! COPYRIGHT
4383 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4384 : : !! This file is distributed under the terms of the
4385 : : !! GNU General Public License, see ~abinit/COPYING
4386 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4387 : : !!
4388 : : !! INPUTS
4389 : : !! op=ctqmc
4390 : : !! oFileIn=file stream
4391 : : !!
4392 : : !! OUTPUT
4393 : : !!
4394 : : !! SIDE EFFECTS
4395 : : !!
4396 : : !! NOTES
4397 : : !!
4398 : : !! SOURCE
4399 : :
4400 : 0 : SUBROUTINE CtqmcoffdiagComplex_printE(op,oFileIn)
4401 : :
4402 : : !Arguments ------------------------------------
4403 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN) :: op
4404 : : INTEGER , OPTIONAL, INTENT(IN) :: oFileIn
4405 : : !Local variables ------------------------------
4406 : : INTEGER :: oFile
4407 : : DOUBLE PRECISION :: E
4408 : : DOUBLE PRECISION :: Noise
4409 : :
4410 : : !IF ( op%rank .NE. MOD(3,op%size)) RETURN
4411 [ # # ]: 0 : IF ( op%rank .NE. MOD(op%size+3,op%size)) RETURN
4412 : :
4413 : 0 : oFile = 42
4414 [ # # ]: 0 : IF ( PRESENT(oFileIn) ) THEN
4415 : 0 : oFile = oFileIn
4416 : : ELSE
4417 : 0 : OPEN(UNIT=oFile, FILE="BetaENoise.dat")
4418 : : END IF
4419 : :
4420 : 0 : CALL CtqmcoffdiagComplex_getE(op,E,Noise)
4421 : :
4422 : 0 : WRITE(oFile,'(1x,F3.2,A2,ES21.14,A2,ES21.14)') op%beta, " ", E, " ", Noise
4423 : :
4424 [ # # ]: 0 : IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
4425 : :
4426 : : END SUBROUTINE CtqmcoffdiagComplex_printE
4427 : : !!***
4428 : :
4429 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
4430 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printPerturbation
4431 : : !! NAME
4432 : : !! CtqmcoffdiagComplex_printPerturbation
4433 : : !!
4434 : : !! FUNCTION
4435 : : !! print perturbation order
4436 : : !!
4437 : : !! COPYRIGHT
4438 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4439 : : !! This file is distributed under the terms of the
4440 : : !! GNU General Public License, see ~abinit/COPYING
4441 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4442 : : !!
4443 : : !! INPUTS
4444 : : !! op=ctqmc
4445 : : !! oFileIn=file stream
4446 : : !!
4447 : : !! OUTPUT
4448 : : !! argout(sizeout)=description
4449 : : !!
4450 : : !! SIDE EFFECTS
4451 : : !!
4452 : : !! NOTES
4453 : : !!
4454 : : !! SOURCE
4455 : :
4456 : 0 : SUBROUTINE CtqmcoffdiagComplex_printPerturbation(op, oFileIn)
4457 : :
4458 : : !Arguments ------------------------------------
4459 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN) :: op
4460 : : INTEGER , OPTIONAL, INTENT(IN) :: oFileIn
4461 : : !Local variables-------------------------------
4462 : : INTEGER :: oFile
4463 : : INTEGER :: iorder
4464 : : INTEGER :: order
4465 : : INTEGER :: iflavor
4466 : : INTEGER :: flavors
4467 : : CHARACTER(LEN=2) :: a
4468 : : CHARACTER(LEN=50) :: string
4469 : :
4470 : : !IF ( op%rank .NE. MOD(4,op%size)) RETURN
4471 [ # # ]: 0 : IF ( op%rank .NE. MOD(op%size+4,op%size)) RETURN
4472 [ # # ]: 0 : IF ( op%opt_order .LE. 0 ) RETURN
4473 : :
4474 : 0 : oFile = 43
4475 [ # # ]: 0 : IF ( PRESENT(oFileIn) ) THEN
4476 : 0 : oFile = oFileIn
4477 : : ELSE
4478 : 0 : OPEN(UNIT=oFile, FILE="Perturbation.dat")
4479 : : END IF
4480 : :
4481 : 0 : order = op%opt_order
4482 : 0 : flavors = op%flavors
4483 : :
4484 : 0 : WRITE(a,'(I2)') flavors
4485 [ # # # # ]: 0 : string = '(I5,'//TRIM(ADJUSTL(a))//'F19.15)'
4486 [ # # ]: 0 : DO iorder = 1, order
4487 : 0 : WRITE(oFile,string) iorder-1, &
4488 [ # # # # : 0 : (/ (op%measPerturbation(iorder, iflavor), iflavor=1, flavors) /)
# # ]
4489 : : END DO
4490 : :
4491 [ # # ]: 0 : IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
4492 : : END SUBROUTINE CtqmcoffdiagComplex_printPerturbation
4493 : : !!***
4494 : :
4495 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printCorrelation
4496 : : !! NAME
4497 : : !! CtqmcoffdiagComplex_printCorrelation
4498 : : !!
4499 : : !! FUNCTION
4500 : : !! print correlation fonctions
4501 : : !!
4502 : : !! COPYRIGHT
4503 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4504 : : !! This file is distributed under the terms of the
4505 : : !! GNU General Public License, see ~abinit/COPYING
4506 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4507 : : !!
4508 : : !! INPUTS
4509 : : !! op=ctqmc
4510 : : !! oFileIn=file stream
4511 : : !!
4512 : : !! OUTPUT
4513 : : !!
4514 : : !! SIDE EFFECTS
4515 : : !!
4516 : : !! NOTES
4517 : : !!
4518 : : !! SOURCE
4519 : :
4520 : 0 : SUBROUTINE CtqmcoffdiagComplex_printCorrelation(op, oFileIn)
4521 : :
4522 : : !Arguments ------------------------------------
4523 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN) :: op
4524 : : INTEGER , OPTIONAL, INTENT(IN) :: oFileIn
4525 : : !Local variables ------------------------------
4526 : : INTEGER :: oFile
4527 : : INTEGER :: itime
4528 : : INTEGER :: sp1
4529 : : INTEGER :: iflavor
4530 : : INTEGER :: i
4531 : : INTEGER :: flavors
4532 : : CHARACTER(LEN=2) :: a
4533 : : CHARACTER(LEN=50) :: string
4534 : : DOUBLE PRECISION :: dt
4535 : :
4536 : : !IF ( op%rank .NE. MOD(5,op%size)) RETURN
4537 [ # # ]: 0 : IF ( op%rank .NE. MOD(op%size+5,op%size)) RETURN
4538 [ # # ]: 0 : IF ( op%opt_analysis .NE. 1 ) RETURN
4539 : :
4540 : 0 : oFile = 44
4541 [ # # ]: 0 : IF ( PRESENT(oFileIn) ) THEN
4542 : 0 : oFile = oFileIn
4543 : : ELSE
4544 : 0 : OPEN(UNIT=oFile, FILE="Correlation.dat")
4545 : : END IF
4546 : :
4547 : 0 : sp1 = op%samples
4548 : 0 : dt = op%beta / sp1
4549 : 0 : sp1 = sp1 + 1
4550 : 0 : flavors = op%flavors
4551 : :
4552 : 0 : i = 3*flavors + 1
4553 : 0 : WRITE(a,'(I2)') i
4554 : 0 : WRITE(oFile,*) "# time (/ (segement, antiseg, correl), i=1, flavor/)"
4555 [ # # # # ]: 0 : string = '(1x,'//TRIM(ADJUSTL(a))//'F19.15)'
4556 [ # # ]: 0 : DO itime = 1, sp1
4557 : 0 : WRITE(oFile,string) DBLE(itime-1)*dt, &
4558 : : (/ ( &
4559 : : (/ ( op%measCorrelation(itime, i, iflavor), i=1,3) /) &
4560 [ # # # # : 0 : , iflavor=1, flavors) /)
# # ]
4561 : : END DO
4562 : :
4563 [ # # ]: 0 : IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
4564 : :
4565 : : END SUBROUTINE CtqmcoffdiagComplex_printCorrelation
4566 : : !!***
4567 : : !#endif
4568 : :
4569 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_printSpectra
4570 : : !! NAME
4571 : : !! CtqmcoffdiagComplex_printSpectra
4572 : : !!
4573 : : !! FUNCTION
4574 : : !! print fourier transform of time evolution of number of electrons
4575 : : !!
4576 : : !! COPYRIGHT
4577 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4578 : : !! This file is distributed under the terms of the
4579 : : !! GNU General Public License, see ~abinit/COPYING
4580 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4581 : : !!
4582 : : !! INPUTS
4583 : : !! op=ctqmc
4584 : : !! oFileIn=file stream
4585 : : !!
4586 : : !! OUTPUT
4587 : : !!
4588 : : !! SIDE EFFECTS
4589 : : !!
4590 : : !! NOTES
4591 : : !!
4592 : : !! SOURCE
4593 : :
4594 : 0 : SUBROUTINE CtqmcoffdiagComplex_printSpectra(op, oFileIn)
4595 : :
4596 : : !Arguments ------------------------------------
4597 : : TYPE(CtqmcoffdiagComplex) , INTENT(IN) :: op
4598 : : INTEGER , OPTIONAL, INTENT(IN) :: oFileIn
4599 : : !Local variables ------------------------------
4600 : : INTEGER :: oFile
4601 : : INTEGER :: flavors
4602 : : INTEGER :: indDensity
4603 : : INTEGER :: endDensity
4604 : : CHARACTER(LEN=4) :: a
4605 : : CHARACTER(LEN=16) :: formatSpectra
4606 : :
4607 : : !IF ( op%rank .NE. MOD(6,op%size)) RETURN
4608 [ # # ]: 0 : IF ( op%opt_spectra .LT. 1 ) RETURN
4609 : :
4610 : 0 : oFile = 45+op%rank
4611 : 0 : a ="0000"
4612 : 0 : WRITE(a,'(I4)') op%rank
4613 [ # # ]: 0 : IF ( PRESENT(oFileIn) ) THEN
4614 : 0 : oFile = oFileIn
4615 : : ELSE
4616 [ # # ]: 0 : OPEN(UNIT=oFile, FILE="Markov_"//TRIM(ADJUSTL(a))//".dat")
4617 : : END IF
4618 : :
4619 : 0 : flavors = op%flavors
4620 : 0 : WRITE(a,'(I4)') flavors+1
4621 [ # # # # ]: 0 : formatSpectra ='(1x,'//TRIM(ADJUSTL(a))//'ES22.14)'
4622 : 0 : WRITE(oFile,*) "# freq[/hermalization] FFT"
4623 : :
4624 : : !PROBLEM eos_gnu_13.2_mpich . %endDensity was introduced throughout
4625 : 0 : endDensity = SIZE(op%density,2)
4626 : : ! endDensity=op%endDensity
4627 : : !ENDPROBLEM
4628 [ # # ]: 0 : DO WHILE ( op%density(flavors+1,endDensity) .EQ. -1 )
4629 : 0 : endDensity = endDensity -1
4630 : : END DO
4631 : :
4632 [ # # ]: 0 : DO indDensity = 1, endDensity
4633 : 0 : WRITE(oFile,formatSpectra) op%density(flavors+1,indDensity), op%density(1:flavors,indDensity)
4634 : : END DO
4635 : :
4636 [ # # ]: 0 : IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
4637 : :
4638 : : END SUBROUTINE CtqmcoffdiagComplex_printSpectra
4639 : : !!***
4640 : :
4641 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_destroy
4642 : : !! NAME
4643 : : !! CtqmcoffdiagComplex_destroy
4644 : : !!
4645 : : !! FUNCTION
4646 : : !! destroy and deallocate all variables
4647 : : !!
4648 : : !! COPYRIGHT
4649 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
4650 : : !! This file is distributed under the terms of the
4651 : : !! GNU General Public License, see ~abinit/COPYING
4652 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4653 : : !!
4654 : : !! INPUTS
4655 : : !! op=ctqmc
4656 : : !!
4657 : : !! OUTPUT
4658 : : !!
4659 : : !! SIDE EFFECTS
4660 : : !!
4661 : : !! NOTES
4662 : : !!
4663 : : !! SOURCE
4664 : :
4665 : 0 : SUBROUTINE CtqmcoffdiagComplex_destroy(op)
4666 : :
4667 : : !Arguments ------------------------------------
4668 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) :: op
4669 : : !Local variables ------------------------------
4670 : : !INTEGER :: iflavor
4671 : : INTEGER :: flavors
4672 : : INTEGER :: i
4673 : : INTEGER :: j
4674 : : INTEGER :: k
4675 : :
4676 : 0 : flavors = op%flavors
4677 : :
4678 : 0 : CALL ImpurityOperator_destroy(op%Impurity)
4679 : 0 : CALL BathOperatoroffdiagComplex_destroy(op%Bath)
4680 : 0 : CALL Vector_destroy(op%measNoise(1))
4681 : 0 : CALL Vector_destroy(op%measNoise(2))
4682 : :
4683 : : !sui!write(6,*) "before greenhyb_destroy in ctmqc_destroy"
4684 : 0 : CALL GreenHyboffdiagComplex_destroy(op%Greens)
4685 : : !#ifdef CTCtqmcoffdiagComplex_ANALYSIS
4686 [ # # ]: 0 : FREEIF(op%measCorrelation)
4687 [ # # ]: 0 : FREEIF(op%measPerturbation)
4688 [ # # ]: 0 : FREEIF(op%meas_fullemptylines)
4689 [ # # ]: 0 : FREEIF(op%measN)
4690 [ # # ]: 0 : IF ( op%opt_histo .GT. 0 ) THEN
4691 [ # # ]: 0 : FREEIF(op%occup_histo_time)
4692 [ # # ]: 0 : FREEIF(op%suscep)
4693 [ # # ]: 0 : FREEIF(op%occupconfig)
4694 [ # # ]: 0 : FREEIF(op%chi)
4695 [ # # ]: 0 : FREEIF(op%chicharge)
4696 [ # # ]: 0 : FREEIF(op%ntot)
4697 : : END IF
4698 [ # # ]: 0 : FREEIF(op%measDE)
4699 [ # # ]: 0 : FREEIF(op%mu)
4700 [ # # ]: 0 : FREEIF(op%hybri_limit)
4701 [ # # ]: 0 : FREEIF(op%abNoiseG)
4702 [ # # ]: 0 : IF ( ALLOCATED(op%measNoiseG) ) THEN
4703 [ # # ]: 0 : DO i=1,2
4704 [ # # ]: 0 : DO j = 1, op%flavors
4705 [ # # ]: 0 : DO k= 1, op%samples+1
4706 : 0 : CALL Vector_destroy(op%measNoiseG(k,j,i))
4707 : : END DO
4708 : : END DO
4709 : : END DO
4710 [ # # # # : 0 : DT_FREE(op%measNoiseG)
# # # # ]
4711 : : END IF
4712 [ # # ]: 0 : FREEIF(op%density)
4713 : : !#endif
4714 : 0 : op%ostream = 0
4715 : 0 : op%istream = 0
4716 : :
4717 : 0 : op%sweeps = 0
4718 : 0 : op%thermalization = 0
4719 : 0 : op%flavors = 0
4720 : : ! op%endDensity = 0
4721 : 0 : op%samples = 0
4722 : 0 : op%beta = 0.d0
4723 : : ! op%seg_added = 0.d0
4724 : : ! op%anti_added = 0.d0
4725 : : ! op%seg_removed = 0.d0
4726 : : ! op%anti_removed = 0.d0
4727 : : ! op%seg_sign = 0.d0
4728 : : ! op%anti_sign = 0.d0
4729 [ # # ]: 0 : op%stats = 0.d0
4730 : 0 : op%swap = 0.d0
4731 : :
4732 : :
4733 : 0 : op%set = .FALSE.
4734 : 0 : op%done = .FALSE.
4735 : 0 : op%init = .FALSE.
4736 : 0 : END SUBROUTINE CtqmcoffdiagComplex_destroy
4737 : : !!***
4738 : :
4739 : : !!****f* ABINIT/m_CtqmcoffdiagComplex/CtqmcoffdiagComplex_setMagmom
4740 : : !! NAME
4741 : : !! CtqmcoffdiagComplex_setMagmom
4742 : : !!
4743 : : !! FUNCTION
4744 : : !! set the Magnetic moment matrix for susceptibility
4745 : : !!
4746 : : !! COPYRIGHT
4747 : : !! Copyright (C) 2013-2025 ABINIT group (F. Gendron)
4748 : : !! This file is distributed under the terms of the
4749 : : !! GNU General Public License, see ~abinit/COPYING
4750 : : !! or http://www.gnu.org/copyleft/gpl.txt .
4751 : : !!
4752 : : !! INPUTS
4753 : : !!
4754 : : !! OUTPUT
4755 : : !!
4756 : : !! SIDE EFFECTS
4757 : : !!
4758 : : !! NOTES
4759 : : !!
4760 : : !! PARENTS
4761 : : !! Will be filled automatically by the parent script
4762 : : !!
4763 : : !! CHILDREN
4764 : : !! Will be filled automatically by the parent script
4765 : : !!
4766 : : !! SOURCE
4767 : :
4768 [ # # # # : 0 : SUBROUTINE CtqmcoffdiagComplex_setMagmom(op,Magmom_orb, Magmom_spin, Magmom_tot)
# # ]
4769 : :
4770 : : !Arguments ------------------------------------
4771 : : TYPE(CtqmcoffdiagComplex), INTENT(INOUT) ::op
4772 : : DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: Magmom_orb
4773 : : DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: Magmom_spin
4774 : : DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: Magmom_tot
4775 : : !Local variables ------------------------------
4776 : : ! INTEGER :: iflavor1,iflavor2
4777 : :
4778 : : ! do iflavor1=1,10
4779 : : ! do iflavor2=1,10
4780 : : ! if(iflavor1==iflavor2) THEN
4781 : : ! write(6,*) iflavor1, iflavor2, Magmom(iflavor1,iflavor2)
4782 : : ! end if
4783 : : ! end do
4784 : : ! end do
4785 : :
4786 [ # # # # ]: 0 : IF ( SIZE(Magmom_orb) .NE. op%flavors*op%flavors ) &
4787 : 0 : CALL ERROR("CtqmcoffdiagComplex_setU : Wrong Magnetic Moment matrix (size) ")
4788 : :
4789 : 0 : CALL ImpurityOperator_setMagmommat(op%Impurity, Magmom_orb, Magmom_spin, Magmom_tot)
4790 : :
4791 : 0 : END SUBROUTINE CtqmcoffdiagComplex_setMagmom
4792 : : !!***
4793 : :
4794 [ # # # # : 0 : END MODULE m_CtqmcoffdiagComplex
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # ]
4795 : : !!***
4796 : :
4797 : : #undef CTQMC_SLICE1
4798 : : #undef CTQMC_SLICE2
4799 : : #undef CTQMC_SEGME
4800 : : #undef CTQMC_ANTIS
4801 : : #undef CTQMC_ADDED
4802 : : #undef CTQMC_REMOV
4803 : : #undef CTQMC_DETSI
|
