Branch data Line data Source code
1 : :
2 : : #if defined HAVE_CONFIG_H
3 : : #include "config.h"
4 : : #endif
5 : : !!****m* ABINIT/m_MatrixHybComplex
6 : : !! NAME
7 : : !! m_MatrixHybComplex
8 : : !!
9 : : !! FUNCTION
10 : : !! Module to deals with a matrix (mainly used for M matrix in m_BathOperatoroffdiagComplex).
11 : : !! Perform varius operation on matrices.
12 : : !!
13 : : !! COPYRIGHT
14 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
15 : : !! This file is distributed under the terms of the
16 : : !! GNU General Public License, see ~abinit/COPYING
17 : : !! or http://www.gnu.org/copyleft/gpl.txt .
18 : : !!
19 : : !! NOTES
20 : : !!
21 : : !! SOURCE
22 : :
23 : : #include "defs.h"
24 : : MODULE m_MatrixHybComplex
25 : : USE m_Global
26 : : IMPLICIT NONE
27 : :
28 : : !!***
29 : :
30 : : PRIVATE
31 : :
32 : : !!****t* m_MatrixHybComplex/MatrixHybComplex
33 : : !! NAME
34 : : !! MatrixHybComplex
35 : : !!
36 : : !! FUNCTION
37 : : !! This structured datatype contains the necessary data
38 : : !!
39 : : !! COPYRIGHT
40 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
41 : : !! This file is distributed under the terms of the
42 : : !! GNU General Public License, see ~abinit/COPYING
43 : : !! or http://www.gnu.org/copyleft/gpl.txt .
44 : : !!
45 : : !! SOURCE
46 : :
47 : : TYPE, PUBLIC :: MatrixHybComplex
48 : : INTEGER, PRIVATE :: size
49 : : ! size of the matrix
50 : :
51 : : INTEGER :: tail
52 : : ! the size of the matrix that is actually used.
53 : :
54 : : INTEGER, PRIVATE :: iTech = -1
55 : : ! precise if time or frequency values will be used.
56 : :
57 : : INTEGER, PRIVATE :: Wmax
58 : : ! size if the frequency grid for mat_omega
59 : :
60 : : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:,:) :: mat
61 : : ! matrix of size "size"
62 : :
63 : : INTEGER , ALLOCATABLE, DIMENSION(:,:) :: mat_tau
64 : : ! matrix of size "size"
65 : :
66 : : COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:,:,:) :: mat_omega
67 : : ! array of size (Wmax, size, size)
68 : :
69 : : END TYPE MatrixHybComplex
70 : : !!***
71 : :
72 : : !PUBLIC INTERFACE ASSIGNMENT (=)
73 : : ! MODULE PROCEDURE MatrixHybComplex_assign
74 : : !END INTERFACE
75 : :
76 : : PUBLIC :: MatrixHybComplex_init
77 : : PUBLIC :: MatrixHybComplex_setSize
78 : : PRIVATE :: MatrixHybComplex_enlarge
79 : : PUBLIC :: MatrixHybComplex_clear
80 : : PUBLIC :: MatrixHybComplex_assign
81 : : PUBLIC :: MatrixHybComplex_inverse
82 : : PRIVATE :: MatrixHybComplex_LU
83 : : PUBLIC :: MatrixHybComplex_getDet
84 : : PUBLIC :: MatrixHybComplex_print
85 : : PUBLIC :: MatrixHybComplex_destroy
86 : :
87 : : CONTAINS
88 : : !!***
89 : :
90 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_init
91 : : !! NAME
92 : : !! MatrixHybComplex_init
93 : : !!
94 : : !! FUNCTION
95 : : !! initialize
96 : : !!
97 : : !! COPYRIGHT
98 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
99 : : !! This file is distributed under the terms of the
100 : : !! GNU General Public License, see ~abinit/COPYING
101 : : !! or http://www.gnu.org/copyleft/gpl.txt .
102 : : !!
103 : : !! INPUTS
104 : : !! this=matrix
105 : : !! iTech=DO NOT USE => BUG
106 : : !! size=memory size for initialization
107 : : !! Wmax=maximum frequency number
108 : : !!
109 : : !! OUTPUT
110 : : !!
111 : : !! SIDE EFFECTS
112 : : !!
113 : : !! NOTES
114 : : !!
115 : : !! SOURCE
116 : :
117 : 8 : SUBROUTINE MatrixHybComplex_init(this, iTech, size, Wmax)
118 : :
119 : : !Arguments ------------------------------------
120 : : TYPE(MatrixHybComplex) , INTENT(INOUT) :: this
121 : : INTEGER , INTENT(IN ) :: iTech
122 : : INTEGER, OPTIONAL, INTENT(IN ) :: size
123 : : INTEGER, OPTIONAL, INTENT(IN ) :: Wmax
124 : : !Local variables ------------------------------
125 : : INTEGER :: size_val
126 : :
127 : 8 : size_val = Global_SIZE
128 [ + - ]: 8 : IF ( PRESENT(size) ) size_val = size
129 : 8 : this%size = size_val
130 [ - + ]: 8 : FREEIF(this%mat)
131 [ - + + - : 32 : MALLOC(this%mat,(1:size_val,1:size_val))
- + + - +
- + - - +
- + - + ]
132 : 8 : this%tail = 0
133 [ + + + + ]: 15691208 : this%mat = 0.d0
134 : 8 : this%iTech = iTech
135 : 8 : SELECT CASE(this%iTech)
136 : : CASE (GREENHYB_TAU)
137 [ - + ]: 8 : FREEIF(this%mat_tau)
138 [ - + + - : 32 : MALLOC(this%mat_tau,(1:size_val,1:size_val))
- + + - +
- + - - +
- + - + ]
139 : : CASE (GREENHYB_OMEGA)
140 [ # # # # ]: 0 : IF ( PRESENT(Wmax) .AND. Wmax .GT. 0 ) THEN
141 [ # # ]: 0 : FREEIF(this%mat_omega)
142 [ # # # # : 0 : MALLOC(this%mat_omega,(1:Wmax,1:size_val,1:size_val))
# # # # #
# # # # #
# # # # #
# # # ]
143 : 0 : this%Wmax=Wmax
144 : : ELSE
145 : 0 : CALL ERROR("MatrixHybComplex_init : Missing argument Wmax for measurement in omega")
146 : : END IF
147 : : CASE DEFAULT
148 [ + - - ]: 8 : CALL WARNALL("MatrixHybComplex_init : Wrong input argument iTech")
149 : : END SELECT
150 : 8 : END SUBROUTINE MatrixHybComplex_init
151 : : !!***
152 : :
153 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_setSize
154 : : !! NAME
155 : : !! MatrixHybComplex_setSize
156 : : !!
157 : : !! FUNCTION
158 : : !! impose size of the this
159 : : !!
160 : : !! COPYRIGHT
161 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
162 : : !! This file is distributed under the terms of the
163 : : !! GNU General Public License, see ~abinit/COPYING
164 : : !! or http://www.gnu.org/copyleft/gpl.txt .
165 : : !!
166 : : !! INPUTS
167 : : !! this=matrix
168 : : !! new_tail=new size
169 : : !!
170 : : !! OUTPUT
171 : : !!
172 : : !! SIDE EFFECTS
173 : : !!
174 : : !! NOTES
175 : : !!
176 : : !! SOURCE
177 : :
178 : 20815397 : SUBROUTINE MatrixHybComplex_setSize(this,new_tail)
179 : :
180 : : !Arguments ------------------------------------
181 : : TYPE(MatrixHybComplex), INTENT(INOUT) :: this
182 : : INTEGER , INTENT(IN ) :: new_tail
183 : : !Local variables ------------------------------
184 : : INTEGER :: size
185 : :
186 [ - + ]: 20815397 : IF ( .NOT. ALLOCATED(this%mat) ) &
187 : 0 : CALL MatrixHybComplex_init(this,this%iTech,Wmax=this%Wmax)
188 : 20815397 : size = this%size
189 [ - + ]: 20815397 : IF( new_tail .GT. size ) THEN
190 : 0 : CALL MatrixHybComplex_enlarge(this, MAX(Global_SIZE,new_tail-size))
191 : : END IF
192 : 20815397 : this%tail = new_tail
193 : 20815397 : END SUBROUTINE MatrixHybComplex_setSize
194 : : !!***
195 : :
196 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_enlarge
197 : : !! NAME
198 : : !! MatrixHybComplex_enlarge
199 : : !!
200 : : !! FUNCTION
201 : : !! This subroutine enlarges memory space
202 : : !!
203 : : !! COPYRIGHT
204 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
205 : : !! This file is distributed under the terms of the
206 : : !! GNU General Public License, see ~abinit/COPYING
207 : : !! or http://www.gnu.org/copyleft/gpl.txt .
208 : : !!
209 : : !! INPUTS
210 : : !! this=matrix
211 : : !! size=new memory size
212 : : !!
213 : : !! OUTPUT
214 : : !!
215 : : !! SIDE EFFECTS
216 : : !!
217 : : !! NOTES
218 : : !!
219 : : !! SOURCE
220 : :
221 : 0 : SUBROUTINE MatrixHybComplex_enlarge(this, size)
222 : :
223 : : !Arguments ------------------------------------
224 : : TYPE(MatrixHybComplex) , INTENT(INOUT) :: this
225 : : INTEGER, OPTIONAL, INTENT(IN ) :: size
226 : : !Local variables ------------------------------
227 : : INTEGER :: width
228 : : INTEGER :: tail
229 : 0 : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:,:) :: this_temp
230 : 0 : INTEGER , ALLOCATABLE, DIMENSION(:,:) :: this_temp_tau
231 : 0 : COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:,:,:) :: this_temp_omega
232 : : INTEGER :: size_val
233 : :
234 [ # # ]: 0 : IF ( ALLOCATED(this%mat) ) THEN
235 : : FREEIF(this_temp)
236 : : FREEIF(this_temp_tau)
237 : 0 : width = this%size
238 : 0 : tail = this%tail
239 : 0 : size_val = width
240 [ # # ]: 0 : IF ( PRESENT(size) ) size_val = size
241 : :
242 : : ! change size of mat
243 [ # # # # : 0 : MALLOC(this_temp,(1:tail,1:tail))
# # # # #
# # # # #
# # ]
244 [ # # # # ]: 0 : this_temp(1:tail,1:tail) = this%mat(1:tail,1:tail)
245 : 0 : FREE(this%mat)
246 : 0 : this%size = width + size_val
247 [ # # # # : 0 : MALLOC(this%mat,(1:this%size,1:this%size))
# # # # #
# # # # #
# # ]
248 [ # # # # ]: 0 : this%mat(1:tail,1:tail) = this_temp(1:tail,1:tail)
249 [ # # ]: 0 : FREE(this_temp)
250 : 0 : SELECT CASE(this%iTech)
251 : : CASE (GREENHYB_TAU)
252 : :
253 : : ! change size of mat_tau
254 [ # # # # : 0 : MALLOC(this_temp_tau,(1:tail,1:tail))
# # # # #
# # # # #
# # ]
255 [ # # # # ]: 0 : this_temp_tau(1:tail,1:tail) = this%mat_tau(1:tail,1:tail)
256 [ # # ]: 0 : FREE(this%mat_tau)
257 [ # # # # : 0 : MALLOC(this%mat_tau,(1:this%size,1:this%size))
# # # # #
# # # # #
# # ]
258 [ # # # # ]: 0 : this%mat_tau(1:tail,1:tail) = this_temp_tau(1:tail,1:tail)
259 [ # # ]: 0 : FREE(this_temp_tau)
260 : : CASE (GREENHYB_OMEGA)
261 : :
262 : : ! change size of mat_omega
263 [ # # # # : 0 : MALLOC(this_temp_omega,(1:this%Wmax,1:tail,1:tail))
# # # # #
# # # # #
# # # # #
# ]
264 [ # # # # : 0 : this_temp_omega(1:this%Wmax,1:tail,1:tail) = this%mat_omega(1:this%Wmax,1:tail,1:tail)
# # ]
265 [ # # ]: 0 : FREE(this%mat_omega)
266 [ # # # # : 0 : MALLOC(this%mat_omega,(1:this%Wmax,1:this%size,1:this%size))
# # # # #
# # # # #
# # # # #
# ]
267 [ # # # # : 0 : this%mat_omega(1:this%Wmax,1:tail,1:tail) = this_temp_omega(1:this%Wmax,1:tail,1:tail)
# # ]
268 [ # # # # : 0 : FREE(this_temp_omega)
# ]
269 : : END SELECT
270 : : ELSE
271 : 0 : CALL MatrixHybComplex_init(this, this%iTech, size=Global_SIZE, Wmax=this%Wmax)
272 : : END IF
273 : 0 : END SUBROUTINE MatrixHybComplex_enlarge
274 : : !!***
275 : :
276 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_clear
277 : : !! NAME
278 : : !! MatrixHybComplex_clear
279 : : !!
280 : : !! FUNCTION
281 : : !! Clear this
282 : : !!
283 : : !! COPYRIGHT
284 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
285 : : !! This file is distributed under the terms of the
286 : : !! GNU General Public License, see ~abinit/COPYING
287 : : !! or http://www.gnu.org/copyleft/gpl.txt .
288 : : !!
289 : : !! INPUTS
290 : : !! this=matrix
291 : : !!
292 : : !! OUTPUT
293 : : !!
294 : : !! SIDE EFFECTS
295 : : !!
296 : : !! NOTES
297 : : !!
298 : : !! SOURCE
299 : :
300 : 4 : SUBROUTINE MatrixHybComplex_clear(this)
301 : :
302 : : !Arguments ------------------------------------
303 : : TYPE(MatrixHybComplex), INTENT(INOUT) :: this
304 : 4 : this%tail = 0
305 : 4 : END SUBROUTINE MatrixHybComplex_clear
306 : : !!***
307 : :
308 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_assign
309 : : !! NAME
310 : : !! MatrixHybComplex_assign
311 : : !!
312 : : !! FUNCTION
313 : : !! assign this=matrix2
314 : : !!
315 : : !! COPYRIGHT
316 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
317 : : !! This file is distributed under the terms of the
318 : : !! GNU General Public License, see ~abinit/COPYING
319 : : !! or http://www.gnu.org/copyleft/gpl.txt .
320 : : !!
321 : : !! INPUTS
322 : : !! this=Matrix 1
323 : : !! Matrix2=Matrix 2
324 : : !!
325 : : !! OUTPUT
326 : : !!
327 : : !! SIDE EFFECTS
328 : : !!
329 : : !! NOTES
330 : : !!
331 : : !! SOURCE
332 : :
333 : 0 : SUBROUTINE MatrixHybComplex_assign(this, matrix)
334 : :
335 : : !Arguments ------------------------------------
336 : : TYPE(MatrixHybComplex), INTENT(INOUT) :: this
337 : : TYPE(MatrixHybComplex), INTENT(IN ) :: matrix
338 : : !Local variables ------------------------------
339 : : INTEGER :: tail
340 : :
341 : 0 : tail = matrix%tail
342 : 0 : CALL MatrixHybComplex_setSize(this, tail)
343 [ # # # # ]: 0 : this%mat(1:tail,1:tail) = matrix%mat(1:tail,1:tail)
344 [ # # ]: 0 : IF ( this%iTech .NE. matrix%iTech ) &
345 : 0 : CALL ERROR("MatrixHybComplex_assign : not compatible matrices")
346 : 0 : SELECT CASE(this%iTech)
347 : : CASE (GREENHYB_TAU)
348 [ # # # # ]: 0 : this%mat_tau(1:tail,1:tail) = matrix%mat_tau(1:tail,1:tail)
349 : : CASE (GREENHYB_OMEGA)
350 [ # # ]: 0 : IF ( this%Wmax .NE. matrix%Wmax ) &
351 : 0 : CALL ERROR("MatrixHybComplex_assig : different numbers of omega")
352 : : this%mat_omega(1:this%Wmax,1:tail,1:tail) = &
353 [ # # # # : 0 : matrix%mat_omega(1:this%Wmax,1:tail,1:tail)
# # # #
# ]
354 : : END SELECT
355 : :
356 : 0 : END SUBROUTINE MatrixHybComplex_assign
357 : : !!***
358 : :
359 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_inverse
360 : : !! NAME
361 : : !! MatrixHybComplex_inverse
362 : : !!
363 : : !! FUNCTION
364 : : !! inverse the matrix and compute the determinant
365 : : !!
366 : : !! COPYRIGHT
367 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
368 : : !! This file is distributed under the terms of the
369 : : !! GNU General Public License, see ~abinit/COPYING
370 : : !! or http://www.gnu.org/copyleft/gpl.txt .
371 : : !!
372 : : !! INPUTS
373 : : !! this=this
374 : : !!
375 : : !! OUTPUT
376 : : !! determinant=determinant of the this before inversion
377 : : !!
378 : : !! SIDE EFFECTS
379 : : !!
380 : : !! NOTES
381 : : !!
382 : : !! SOURCE
383 : :
384 : 0 : SUBROUTINE MatrixHybComplex_inverse(this,determinant)
385 : :
386 : : !Arguments ------------------------------------
387 : : TYPE(MatrixHybComplex), INTENT(INOUT) :: this
388 : : COMPLEX(KIND=8), OPTIONAL, INTENT(OUT) :: determinant
389 : 0 : COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:,:) :: invMatrix
390 : : !Local variables ------------------------------
391 : : !DOUBLE PRECISION :: reste
392 : : INTEGER :: ligne
393 : : INTEGER :: ligne_virtuelle
394 : : INTEGER :: colonne
395 : : INTEGER :: sys
396 : : INTEGER :: sys_vir
397 : : INTEGER :: tail
398 : 0 : INTEGER, DIMENSION(:), ALLOCATABLE :: pivot
399 : : COMPLEX(KIND=8) :: det
400 : :
401 : 0 : tail = this%tail
402 [ # # ]: 0 : IF ( tail .EQ. 0 ) THEN
403 [ # # ]: 0 : IF ( PRESENT(determinant) ) determinant = cmplx(1.d0,0.d0,kind=8)
404 : : RETURN
405 : : END IF
406 : :
407 [ # # # # : 0 : MALLOC(invMatrix,(1:tail,1:tail))
# # # # #
# # # # #
# # ]
408 : :
409 : : FREEIF(pivot)
410 [ # # # # : 0 : MALLOC(pivot,(1:tail))
# # # # #
# # # ]
411 : :
412 : 0 : CALL MatrixHybComplex_LU(this, pivot=pivot, determinant=det)
413 : :
414 : : !det = 1.d0
415 [ # # ]: 0 : DO sys = 1, tail
416 : 0 : sys_vir = pivot(sys)
417 : : !! DO ligne=1,tail
418 : : !! ligne_virtuelle = pivot(ligne)
419 : : !! reste = 0.d0
420 : : !! DO colonne = 1, ligne-1
421 : : !! reste = reste + this%mat(ligne_virtuelle,colonne)*invMatrix%mat(colonne,sys_vir)
422 : : !! END DO
423 : : !! IF ( ligne .EQ. sys ) THEN
424 : : !! invMatrix%mat(ligne,sys_vir) = 1.d0 - reste
425 : : !! ELSE
426 : : !! invMatrix%mat(ligne,sys_vir) = 0.d0 - reste
427 : : !! END IF
428 : : !! END DO
429 : : !!
430 : : !! DO ligne = tail, 1, -1
431 : : !! ligne_virtuelle=pivot(ligne)
432 : : !! reste = 0.d0
433 : : !! DO colonne = ligne+1, tail
434 : : !! reste= reste + this%mat(ligne_virtuelle,colonne)*invMatrix%mat(colonne,sys_vir)
435 : : !! END DO
436 : : !! invMatrix%mat(ligne, sys_vir) = (invMatrix%mat(ligne, sys_vir)-reste)/this%mat(ligne_virtuelle,ligne)
437 : : !! END DO
438 : : !det = det*this%mat(sys_vir,sys)
439 [ # # ]: 0 : invMatrix(:,sys_vir) = cmplx(0.d0,0.d0,kind=8)
440 : 0 : invMatrix(sys,sys_vir) = cmplx(1.d0,0.d0,kind=8)
441 [ # # ]: 0 : DO colonne = 1, tail-1
442 [ # # ]: 0 : DO ligne=colonne+1,tail
443 : 0 : ligne_virtuelle = pivot(ligne)
444 : : invMatrix(ligne, sys_vir) = invMatrix(ligne,sys_vir) &
445 : : - this%mat(ligne_virtuelle,colonne) &
446 : 0 : * invMatrix(colonne,sys_vir)
447 : : END DO
448 : : END DO
449 : :
450 [ # # ]: 0 : DO colonne = tail, 1, -1
451 : 0 : invMatrix(colonne, sys_vir) = invMatrix(colonne, sys_vir) / this%mat(pivot(colonne),colonne)
452 [ # # ]: 0 : DO ligne = 1, colonne-1
453 : 0 : ligne_virtuelle = pivot(ligne)
454 : : invMatrix(ligne, sys_vir) = invMatrix(ligne,sys_vir) &
455 : : - this%mat(ligne_virtuelle,colonne) &
456 : 0 : * invMatrix(colonne, sys_vir)
457 : : END DO
458 : : END DO
459 : :
460 : : END DO
461 [ # # ]: 0 : FREE(pivot)
462 [ # # # # ]: 0 : this%mat(1:tail,1:tail) = invMatrix(1:tail,1:tail)
463 [ # # ]: 0 : IF ( PRESENT(determinant) ) THEN
464 : 0 : determinant = det
465 : : END IF
466 [ # # ]: 0 : FREE(invMatrix)
467 : 0 : END SUBROUTINE MatrixHybComplex_inverse
468 : : !!***
469 : :
470 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_LU
471 : : !! NAME
472 : : !! MatrixHybComplex_LU
473 : : !!
474 : : !! FUNCTION
475 : : !! LU decomposition
476 : : !!
477 : : !! COPYRIGHT
478 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
479 : : !! This file is distributed under the terms of the
480 : : !! GNU General Public License, see ~abinit/COPYING
481 : : !! or http://www.gnu.org/copyleft/gpl.txt .
482 : : !!
483 : : !! INPUTS
484 : : !! this=matrix
485 : : !!
486 : : !! OUTPUT
487 : : !! pivot=gauss pivot
488 : : !! determinant=determinant of the this
489 : : !!
490 : : !! SIDE EFFECTS
491 : : !!
492 : : !! NOTES
493 : : !!
494 : : !! SOURCE
495 : :
496 : 0 : SUBROUTINE MatrixHybComplex_LU(this,pivot,determinant)
497 : :
498 : : !Arguments ------------------------------------
499 : : TYPE(MatrixHybComplex), INTENT(INOUT) :: this
500 : : INTEGER, DIMENSION(:), ALLOCATABLE, OPTIONAL, INTENT(INOUT) :: pivot
501 : : COMPLEX(KIND=8), OPTIONAL, INTENT(OUT) :: determinant
502 : : !Local variables ------------------------------
503 : : INTEGER :: ligne
504 : : INTEGER :: colonne
505 : : INTEGER :: colonne_de_1_ligne
506 : : INTEGER :: max_ind_lig
507 : : INTEGER :: ligne_virtuelle
508 : : INTEGER :: tail
509 : 0 : INTEGER, DIMENSION(:), ALLOCATABLE :: pivot_tmp
510 : : DOUBLE PRECISION :: max_col
511 : : COMPLEX(KIND=8) :: inverse_pivot
512 : : COMPLEX(KIND=8) :: coef
513 : : COMPLEX(KIND=8) :: det
514 : 0 : COMPLEX(KIND=8), DIMENSION(:,:), ALLOCATABLE :: mat_tmp
515 : :
516 : 0 : tail = this%tail
517 : 0 : det = cmplx(1.d0,0.d0,kind=8)
518 : : FREEIF(pivot_tmp)
519 [ # # # # : 0 : MALLOC(pivot_tmp,(1:tail))
# # # # #
# # # ]
520 : :
521 [ # # ]: 0 : DO ligne=1, tail
522 : 0 : pivot_tmp(ligne)=ligne
523 : : END DO
524 : :
525 [ # # # # : 0 : MALLOC(mat_tmp,(1:tail,1:tail))
# # # # #
# # # # #
# # ]
526 [ # # # # ]: 0 : mat_tmp(1:tail,1:tail) = this%mat(1:tail,1:tail)
527 [ # # ]: 0 : DO colonne = 1, tail-1
528 : 0 : max_col = ABS(mat_tmp(pivot_tmp(colonne),colonne))
529 : 0 : max_ind_lig = colonne
530 : :
531 [ # # ]: 0 : DO ligne = colonne+1, tail
532 : 0 : ligne_virtuelle = pivot_tmp(ligne)
533 [ # # ]: 0 : IF ( ABS( mat_tmp(ligne_virtuelle,colonne)).GT.max_col) then
534 : 0 : max_col = ABS(mat_tmp(ligne_virtuelle,colonne))
535 : 0 : max_ind_lig = ligne
536 : : ENDIF
537 : : END DO
538 : :
539 : 0 : ligne = pivot_tmp(colonne)
540 : 0 : pivot_tmp(colonne) = pivot_tmp(max_ind_lig)
541 : 0 : pivot_tmp(max_ind_lig) = ligne
542 [ # # ]: 0 : IF ( pivot_tmp(colonne) .NE. pivot_tmp(max_ind_lig) ) det = det * (-1.d0)
543 : 0 : inverse_pivot=1.d0 / mat_tmp(pivot_tmp(colonne),colonne)
544 : 0 : det = det * mat_tmp(pivot_tmp(colonne),colonne)
545 [ # # ]: 0 : DO ligne = colonne+1, tail
546 : 0 : ligne_virtuelle = pivot_tmp(ligne)
547 : 0 : coef = mat_tmp(ligne_virtuelle,colonne)*inverse_pivot
548 : 0 : mat_tmp(ligne_virtuelle,colonne) = coef
549 [ # # ]: 0 : DO colonne_de_1_ligne = colonne+1, tail
550 : : mat_tmp(ligne_virtuelle,colonne_de_1_ligne)= mat_tmp(ligne_virtuelle,colonne_de_1_ligne)&
551 : 0 : -coef * mat_tmp(pivot_tmp(colonne) ,colonne_de_1_ligne)
552 : : END DO
553 : : END DO
554 : : END DO
555 : 0 : det = det * mat_tmp(pivot_tmp(tail),tail)
556 [ # # ]: 0 : IF ( PRESENT(determinant) ) &
557 : 0 : determinant = det
558 [ # # ]: 0 : IF ( PRESENT(pivot) ) THEN
559 [ # # # # ]: 0 : this%mat(1:tail,1:tail) = mat_tmp(1:tail,1:tail)
560 [ # # # # ]: 0 : IF ( ALLOCATED(pivot) .AND. SIZE(pivot) .NE. tail ) THEN
561 : 0 : FREE(pivot)
562 [ # # # # : 0 : MALLOC(pivot,(1:tail))
# # # # ]
563 [ # # ]: 0 : ELSE IF ( .NOT. ALLOCATED(pivot) ) THEN
564 [ # # # # : 0 : MALLOC(pivot,(1:tail))
# # # # #
# ]
565 : : END IF
566 [ # # ]: 0 : pivot(1:tail)=pivot_tmp(1:tail)
567 : : END IF
568 : 0 : FREE(mat_tmp)
569 : 0 : FREE(pivot_tmp)
570 : 0 : END SUBROUTINE MatrixHybComplex_LU
571 : : !!***
572 : :
573 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_getDet
574 : : !! NAME
575 : : !! MatrixHybComplex_getDet
576 : : !!
577 : : !! FUNCTION
578 : : !! Just get the determinant
579 : : !!
580 : : !! COPYRIGHT
581 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
582 : : !! This file is distributed under the terms of the
583 : : !! GNU General Public License, see ~abinit/COPYING
584 : : !! or http://www.gnu.org/copyleft/gpl.txt .
585 : : !!
586 : : !! INPUTS
587 : : !! this=matrix
588 : : !!
589 : : !! OUTPUT
590 : : !! det=determinant
591 : : !!
592 : : !! SIDE EFFECTS
593 : : !!
594 : : !! NOTES
595 : : !!
596 : : !! SOURCE
597 : :
598 : 0 : SUBROUTINE MatrixHybComplex_getDet(this,det)
599 : :
600 : : !Arguments ------------------------------------
601 : : TYPE(MatrixHybComplex) , INTENT(INOUT) :: this
602 : : COMPLEX(KIND=8), INTENT( OUT) :: det
603 : :
604 [ # # ]: 0 : IF ( this%tail .EQ. 0 ) THEN
605 : 0 : det = cmplx(1.d0,0.d0,kind=8)
606 : 0 : RETURN
607 : : END IF
608 : 0 : CALL MatrixHybComplex_LU(this, determinant=det)
609 : : END SUBROUTINE MatrixHybComplex_getDet
610 : : !!***
611 : :
612 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_print
613 : : !! NAME
614 : : !! MatrixHybComplex_print
615 : : !!
616 : : !! FUNCTION
617 : : !! print Matrix
618 : : !!
619 : : !! COPYRIGHT
620 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
621 : : !! This file is distributed under the terms of the
622 : : !! GNU General Public License, see ~abinit/COPYING
623 : : !! or http://www.gnu.org/copyleft/gpl.txt .
624 : : !!
625 : : !! INPUTS
626 : : !! this=Matrix
627 : : !! ostream=file stream
628 : : !! opt_print=0 mat
629 : : !! 1 mat_tau
630 : : !! 2 both
631 : : !!
632 : : !! OUTPUT
633 : : !!
634 : : !! SIDE EFFECTS
635 : : !!
636 : : !! NOTES
637 : : !!
638 : : !! SOURCE
639 : :
640 : 0 : SUBROUTINE MatrixHybComplex_print(this,ostream,opt_print)
641 : :
642 : : !Arguments ------------------------------------
643 : : TYPE(MatrixHybComplex), INTENT(IN) :: this
644 : : INTEGER, OPTIONAL, INTENT(IN) :: ostream
645 : : INTEGER, OPTIONAL, INTENT(IN) :: opt_print
646 : : !Local variables ------------------------------
647 : : INTEGER :: ostream_val
648 : : INTEGER :: opt_val
649 : : INTEGER :: it1
650 : : INTEGER :: it2
651 : : CHARACTER(LEN=4 ) :: size
652 : : CHARACTER(LEN=22) :: string
653 : :
654 : 0 : ostream_val = 6
655 : 0 : opt_val=0
656 [ # # ]: 0 : IF ( PRESENT(ostream) ) ostream_val = ostream
657 [ # # ]: 0 : IF ( PRESENT(opt_print) ) opt_val = opt_print
658 : 0 : WRITE(size,'(I4)') this%tail
659 [ # # ]: 0 : IF ( MOD(opt_val,2) .EQ. 0 ) THEN
660 [ # # # # ]: 0 : string ='(1x,A,1x,'//TRIM(ADJUSTL(size))//'(ES10.2,1x),A)'
661 : 0 : WRITE(ostream_val,'(A)') "["
662 [ # # ]: 0 : DO it1 = 1, this%tail
663 [ # # # # : 0 : WRITE(ostream_val,string) "[",(/ (this%mat(it1,it2),it2=1,this%tail) /)," ]"
# # ]
664 : : END DO
665 : 0 : WRITE(ostream_val,'(A)') "]"
666 : : END IF
667 [ # # ]: 0 : IF ( opt_val .GE.1 ) THEN
668 [ # # # # ]: 0 : string ='(1x,A,1x,'//TRIM(ADJUSTL(size))//'(I4,1x),A)'
669 : 0 : WRITE(ostream_val,'(A)') "["
670 [ # # ]: 0 : DO it1 = 1, this%tail
671 [ # # # # : 0 : WRITE(ostream_val,string) "[",(/ (this%mat_tau(it1,it2),it2=1,this%tail) /),"]"
# # ]
672 : : END DO
673 : 0 : WRITE(ostream_val,'(A)') "]"
674 : : END IF
675 : 0 : END SUBROUTINE MatrixHybComplex_print
676 : : !!***
677 : :
678 : : !!****f* ABINIT/m_MatrixHybComplex/MatrixHybComplex_destroy
679 : : !! NAME
680 : : !! MatrixHybComplex_destroy
681 : : !!
682 : : !! FUNCTION
683 : : !! Destroy
684 : : !!
685 : : !! COPYRIGHT
686 : : !! Copyright (C) 2013-2025 ABINIT group (J. Bieder)
687 : : !! This file is distributed under the terms of the
688 : : !! GNU General Public License, see ~abinit/COPYING
689 : : !! or http://www.gnu.org/copyleft/gpl.txt .
690 : : !!
691 : : !! INPUTS
692 : : !! this=Matrix
693 : : !!
694 : : !! OUTPUT
695 : : !!
696 : : !! SIDE EFFECTS
697 : : !!
698 : : !! NOTES
699 : : !!
700 : : !! SOURCE
701 : :
702 : 8 : SUBROUTINE MatrixHybComplex_destroy(this)
703 : :
704 : : !Arguments ------------------------------------
705 : : TYPE(MatrixHybComplex), INTENT(INOUT) :: this
706 : :
707 [ + - ]: 8 : FREEIF(this%mat)
708 : 16 : SELECT CASE(this%iTech)
709 : : CASE (GREENHYB_TAU)
710 [ + - ]: 8 : FREEIF(this%mat_tau)
711 : : CASE (GREENHYB_OMEGA)
712 [ + - - - : 8 : FREEIF(this%mat_omega)
- ]
713 : : END SELECT
714 : :
715 : 8 : this%tail = 0
716 : 8 : this%size = 0
717 : 8 : this%iTech = -1
718 : 8 : END SUBROUTINE MatrixHybComplex_destroy
719 : : !!***
720 : :
721 [ # # # # : 0 : END MODULE m_MatrixHybComplex
# # # # #
# # # # #
# # # # ]
722 : : !!***
|
