File: vecplot.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (140 lines) | stat: -rw-r--r-- 4,120 bytes parent folder | download
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
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140

! Copyright (C) 2002-2006 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

!BOP
! !ROUTINE: vecplot
! !INTERFACE:
subroutine vecplot
! !DESCRIPTION:
!   Outputs a 2D or 3D vector field for plotting. The vector field can be the
!   magnetisation vector field, ${\bf m}$; the exchange-correlation magnetic
!   field, ${\bf B}_{\rm xc}$; or the electric field
!   ${\bf E}\equiv-\nabla V_{\rm C}$. The magnetisation is obtained from the
!   spin density matrix, $\rho_{\alpha\beta}$, by solving
!   $$ \rho_{\alpha\beta}({\bf r})=\frac{1}{2}\left(n({\bf r})
!    \delta_{\alpha\beta}+\sigma\cdot {\bf m}({\bf r})\right), $$
!   where $n\equiv\tr\rho_{\alpha\beta}$ is the total density. In the case of 2D
!   plots, the magnetisation vectors are still 3D, but are in the coordinate
!   system of the plane.
!
! !REVISION HISTORY:
!   Created August 2004 (JKD)
!   Included electric field plots, August 2006 (JKD)
!EOP
!BOC
use modmain
implicit none
! allocatable arrays
real(8), allocatable :: rvfmt(:,:,:),rvfir(:,:)
! initialise universal variables
call init0
if ((task.eq.72).or.(task.eq.73).or.(task.eq.82).or.(task.eq.83)) then
  if (.not.spinpol) then
    write(*,*)
    write(*,'("Error(vecplot): spin-unpolarised magnetisation/field is zero")')
    write(*,*)
    stop
  end if
end if
! read magnetisation and exchange-correlation magnetic field from file
call readstate
allocate(rvfmt(npmtmax,natmtot,3),rvfir(ngtot,3))
select case(task)
case(72,73)
! magnetisation
  if (ncmag) then
! non-collinear
    rvfmt(:,:,:)=magmt(:,:,:)
    rvfir(:,:)=magir(:,:)
  else
! collinear
    rvfmt(:,:,1:2)=0.d0
    rvfir(:,1:2)=0.d0
    rvfmt(:,:,3)=magmt(:,:,1)
    rvfir(:,3)=magir(:,1)
  end if
case(82,83)
! exchange-correlation magnetic field
  if (ncmag) then
! non-collinear
    rvfmt(:,:,:)=bxcmt(:,:,:)
    rvfir(:,:)=bxcir(:,:)
  else
! collinear
    rvfmt(:,:,1:2)=0.d0
    rvfir(:,1:2)=0.d0
    rvfmt(:,:,3)=bxcmt(:,:,1)
    rvfir(:,3)=bxcir(:,1)
  end if
case(142,143)
! electric field
  call gradrf(vclmt,vclir,rvfmt,rvfir)
! use the negative of the gradient
  rvfmt(:,:,:)=-rvfmt(:,:,:)
  rvfir(:,:)=-rvfir(:,:)
case(152,153)
  if (.not.ncmag) then
    write(*,*)
    write(*,'("Error(vecplot): collinear m(r) x B_xc(r) is zero")')
    write(*,*)
    stop
  end if
  call rvfcross(magmt,magir,bxcmt,bxcir,rvfmt,rvfir)
end select
select case(task)
case(72,82,142,152)
! project the magnetisation/field onto the plotting plane
  call proj2d(rvfmt,rvfir)
  if (task.eq.72) then
    open(50,file='MAG2D.OUT',form='FORMATTED')
  else if (task.eq.82) then
    open(50,file='BXC2D.OUT',form='FORMATTED')
  else if (task.eq.142) then
    open(50,file='EF2D.OUT',form='FORMATTED')
  else
    open(50,file='MCBXC2D.OUT',form='FORMATTED')
  end if
  call plot2d(50,3,rvfmt,rvfir)
  close(50)
  write(*,*)
  write(*,'("Info(vecplot):")')
  if (task.eq.72) then
    write(*,'(" 2D magnetisation density written to MAG2D.OUT")')
  else if (task.eq.82) then
    write(*,'(" 2D exchange-correlation field written to BXC2D.OUT")')
  else if (task.eq.142) then
    write(*,'(" 2D electric field written to EF2D.OUT")')
  else
    write(*,'(" 2D m(r) x B_xc(r) written to MCBXC2D.OUT")')
  end if
case(73,83,143,153)
  if (task.eq.73) then
    open(50,file='MAG3D.OUT',form='FORMATTED')
  else if (task.eq.83) then
    open(50,file='BXC3D.OUT',form='FORMATTED')
  else if (task.eq.143) then
    open(50,file='EF3D.OUT',form='FORMATTED')
  else
    open(50,file='MCBXC3D.OUT',form='FORMATTED')
  end if
  call plot3d(50,3,rvfmt,rvfir)
  close(50)
  write(*,*)
  write(*,'("Info(vecplot):")')
  if (task.eq.73) then
    write(*,'(" 3D magnetisation density written to MAG3D.OUT")')
  else if (task.eq.83) then
    write(*,'(" 3D exchange-correlation field written to BXC3D.OUT")')
  else if (task.eq.143) then
    write(*,'(" 3D electric field written to EF3D.OUT")')
  else
    write(*,'(" 3D m(r) x B_xc(r) written to MCBXC3D.OUT")')
  end if
end select
deallocate(rvfmt,rvfir)
return
end subroutine
!EOC