forked from jautschbach/mcd-molcas
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathdiagonalize-magdip-gs.f90
103 lines (76 loc) · 3 KB
/
diagonalize-magdip-gs.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
subroutine diagonalize_magdip_gs (idir)
! -------------------------------------------------
! diagonalize the Zeeman Hamiltonian component idir
! within the degenerate GS
! -------------------------------------------------
use definitions
use constants_parameters
use namelist_module
use shared_variables
implicit none
integer(KINT), intent(in) :: idir
integer(KINT) :: jdir
complex(KREAL), dimension(:,:), allocatable :: eigv, tmpmat1, tmpmat2
allocate(eigv(degen,degen))
! the next call uses the complex lapack routine zheevd.
! upon return, the input matrix is replaced with the
! eigenvectors.
! we then transform the blocks
eigv = 0
eigv(1:degen,1:degen) = magdip(1:degen,1:degen,idir)
call diagonalize_matrix (degen, eigv)
allocate(tmpmat1(degen, degen))
tmpmat1 = transpose(conjg(eigv))
! transform the relevant blocks of the magnetic dipole matrix
do jdir = 1,3
allocate(tmpmat2(degen, nstates))
tmpmat2 = matmul(tmpmat1, magdip(1:degen,:,jdir))
magdip(1:degen,:, jdir) = tmpmat2
deallocate (tmpmat2)
allocate (tmpmat2(nstates, degen))
tmpmat2 = matmul(magdip(:, 1:degen, jdir), eigv)
magdip(:, 1:degen, jdir) = tmpmat2(:, 1:degen)
deallocate (tmpmat2)
end do
call print_rec_matrix(out, degen, real(magdip(1:degen,1:degen,idir)),&
& 'Transformed Zeeman Hamiltonian REAL part')
call print_rec_matrix(out, degen, aimag(magdip(1:degen,1:degen,idir)),&
& 'Transformed Zeeman Hamiltonian IMAG part')
! -------------------------------------------------------
! transform the dipole and quad matrix elements between the GS
! components and the ESs accordingly, store in *
! -------------------------------------------------------
do jdir = 1,3
allocate(tmpmat2(degen, nstates))
tmpmat2 = matmul(tmpmat1, eldip(1:degen,:,jdir))
eldip(1:degen,:, jdir) = tmpmat2
deallocate (tmpmat2)
allocate (tmpmat2(nstates, degen))
tmpmat2 = matmul(eldip(:, 1:degen, jdir), eigv)
eldip(:, 1:degen, jdir) = tmpmat2(:, 1:degen)
deallocate (tmpmat2)
end do ! jdir
if (havequad) then
do jdir = 1,6
allocate(tmpmat2(degen, nstates))
tmpmat2 = matmul(tmpmat1, elquad(1:degen,:,jdir))
elquad(1:degen,:, jdir) = tmpmat2
deallocate (tmpmat2)
allocate (tmpmat2(nstates, degen))
tmpmat2 = matmul(elquad(:, 1:degen, jdir), eigv)
elquad(:, 1:degen, jdir) = tmpmat2(:, 1:degen)
deallocate (tmpmat2)
end do ! jdir
end if ! havequad
deallocate(tmpmat1)
deallocate(eigv)
!!$ if (dbg>1) then
!!$ allocate (tmpmat1(nstates, nstates))
!!$ do jdir = 1,3
!!$ tmpmat1 = transpose(conjg(eldip(:,:,jdir)))
!!$ tmpmat1 = tmpmat1 - eldip(:,:,jdir)
!!$ write (out,*) jdir, pack(tmpmat1(:,:), abs(tmpmat1(:,:))>tiny)
!!$ end do ! jdir
!!$ deallocate (tmpmat1)
!!$ end if ! dbg
end subroutine diagonalize_magdip_gs