1
2! Copyright (C) 2002-2005 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
3! This file is distributed under the terms of the GNU General Public License.
4! See the file COPYING for license details.
5
6subroutine init2
7use modmain
8use modrdm
9use modphonon
10use modvars
11use modmpi
12implicit none
13! local variables
14logical lsym(48)
15integer isym,iv(3)
16real(8) boxl(3,0:3)
17real(8) ts0,ts1
18
19call timesec(ts0)
20
21!---------------------!
22!     q-point set     !
23!---------------------!
24! check if the system is an isolated molecule
25if (molecule) ngridq(:)=1
26! store the point group symmetries for reducing the q-point set
27if (reduceq.eq.0) then
28  nsymqpt=1
29  symqpt(:,:,1)=symlat(:,:,1)
30else
31  lsym(:)=.false.
32  do isym=1,nsymcrys
33    lsym(lsplsymc(isym))=.true.
34  end do
35  nsymqpt=0
36  do isym=1,nsymlat
37    if (lsym(isym)) then
38      nsymqpt=nsymqpt+1
39      symqpt(:,:,nsymqpt)=symlat(:,:,isym)
40    end if
41  end do
42end if
43if (any(task.eq.[105,180,185,320,330,331])) then
44! equal k- and q-point grids for nesting function, BSE and linear-reposnse TDDFT
45  ngridq(:)=ngridk(:)
46else if ((xctype(1).lt.0).or.(any(task.eq.[5,300,600,620,630,800,801]))) then
47! allow the q-point grid to be smaller than the k-point grid for OEP,
48! Hartree-Fock, RDMFT and GW
49  if (any(ngridq(:).le.0)) ngridq(:)=ngridk(:)
50else
51  ngridq(:)=abs(ngridq(:))
52end if
53! check that the q-point and k-point grids are commensurate for some tasks
54if ((xctype(1).lt.0).or.(any(task.eq.[5,205,240,241,300,600,620,630, &
55 800,801]))) then
56  iv(:)=mod(ngridk(:),ngridq(:))
57  if ((iv(1).ne.0).or.(iv(2).ne.0).or.(iv(3).ne.0)) then
58    write(*,*)
59    write(*,'("Error(init2): k-point grid incommensurate with q-point grid")')
60    write(*,'(" ngridk : ",3I6)') ngridk
61    write(*,'(" ngridq : ",3I6)') ngridq
62    write(*,*)
63    stop
64  end if
65end if
66! allocate the q-point arrays
67if (allocated(ivqiq)) deallocate(ivqiq)
68allocate(ivqiq(0:ngridq(1)-1,0:ngridq(2)-1,0:ngridq(3)-1))
69if (allocated(ivqiqnr)) deallocate(ivqiqnr)
70allocate(ivqiqnr(0:ngridq(1)-1,0:ngridq(2)-1,0:ngridq(3)-1))
71nqptnr=ngridq(1)*ngridq(2)*ngridq(3)
72if (allocated(ivq)) deallocate(ivq)
73allocate(ivq(3,nqptnr))
74if (allocated(vql)) deallocate(vql)
75allocate(vql(3,nqptnr))
76if (allocated(vqc)) deallocate(vqc)
77allocate(vqc(3,nqptnr))
78if (allocated(wqpt)) deallocate(wqpt)
79allocate(wqpt(nqptnr))
80! set up the q-point box (offset should always be zero)
81boxl(:,:)=0.d0
82boxl(1,1)=1.d0; boxl(2,2)=1.d0; boxl(3,3)=1.d0
83! generate the q-point set
84! (note that the vectors vql and vqc are in the first Brillouin zone)
85call genppts(.true.,nsymqpt,symqpt,ngridq,nqptnr,epslat,bvec,boxl,nqpt,ivqiq, &
86 ivqiqnr,ivq,vql,vqc,wqpt,wqptnr)
87! write the q-points to QPOINTS.OUT
88if (mp_mpi) call writeqpts
89! write to VARIABLES.OUT
90call writevars('nsymqpt',iv=nsymqpt)
91call writevars('symqpt',nv=9*nsymqpt,iva=symqpt)
92call writevars('ngridq',nv=3,iva=ngridq)
93call writevars('nqpt',iv=nqpt)
94call writevars('ivqiq',nv=nqptnr,iva=ivqiq)
95call writevars('ivq',nv=3*nqptnr,iva=ivq)
96call writevars('vql',nv=3*nqptnr,rva=vql)
97call writevars('wqpt',nv=nqpt,rva=wqpt)
98
99!--------------------------------------------------------!
100!     OEP, Hartree-Fock, RDMFT, BSE and GW variables     !
101!--------------------------------------------------------!
102if ((xctype(1).lt.0).or.(any(task.eq.[5,180,185,188,205,300,320,330,331,600, &
103 620,630,800,801]))) then
104! determine the regularised Coulomb Green's function for small q
105  call gengclq
106! output the Coulomb Green's function to GCLQ.OUT
107  if (mp_mpi) call writegclq
108! initialise OEP variables
109  if (xctype(1).lt.0) call initoep
110end if
111if (task.eq.300) then
112  if (allocated(vclmat)) deallocate(vclmat)
113  allocate(vclmat(nstsv,nstsv,nkpt))
114  if (allocated(dkdc)) deallocate(dkdc)
115  allocate(dkdc(nstsv,nstsv,nkpt))
116end if
117
118!-------------------------!
119!     phonon variables    !
120!-------------------------!
121if (allocated(wphq)) deallocate(wphq)
122if (task.eq.220) then
123  allocate(wphq(nbph,npp1d))
124end if
125if (any(task.eq.[240,241,250,270,271,280,285])) then
126  allocate(wphq(nbph,nqpt))
127end if
128
129call timesec(ts1)
130timeinit=timeinit+ts1-ts0
131
132end subroutine
133
134