TABLE OF CONTENTS


ABINIT/invars1m [ Functions ]

[ Top ] [ Functions ]

NAME

 invars1m

FUNCTION

 Initialisation phase : prepare the main input subroutine call by
 reading all the NO MULTI variables, as well as the dimensions
 needed for allocating the input arrays in abinit.

COPYRIGHT

 Copyright (C) 1999-2018 ABINIT group (XG, MKV)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt.

INPUTS

  iout=unit number of output file
  lenstr=actual length of string
  msym=default maximal number of symmetries
  mxnatom=maximal value of input natom for all the datasets
  mxnimage=maximal value of input nimage for all the datasets
  ndtset= number of datasets to be read; if 0, no multi-dataset mode
  ndtset_alloc=number of datasets, corrected for allocation of at least
               one data set.
  npsp= number of pseudopotential files
  string*(*)=string of characters containing all input variables and data
  zionpsp(npsp)= valence charge over all psps

OUTPUT

  dmatpuflag=flag controlling the use of an initial density matrix in PAW+U (max. value over datasets)
  mband_upper_(0:ndtset_alloc)=list of mband_upper values
  mxga_n_rules=maximal value of input ga_n_rules for all the datasets
  mxgw_nqlwl=maximal value of input gw_nqlwl for all the datasets
  mxlpawu=maximal value of input lpawu for all the datasets
  mxmband_upper=maximal value of input nband for all the datasets
  mxnatpawu=maximal value of number of atoms on which +U is applied for all the datasets
  mxnatsph=maximal value of input natsph for all the datasets
  mxnatsph_extra=maximal value of input natsph_extra for all the datasets
  mxnatvshift=maximal value of input natvshift for all the datasets
  mxnconeq=maximal value of input nconeq for all the datasets
  mxnkptgw=maximal value of input nkptgw for all the datasets
  mxnkpthf=maximal value of input nkpthf for all the datasets
  mxnkpt=maximal value of input nkpt for all the datasets
  mxnnos=maximal value of input nnos for all the datasets
  mxnqptdm=maximal value of input nqptdm for all the datasets
  mxnspinor=maximal value of input nspinor for all the datasets
  mxnsppol=maximal value of input nsppol for all the datasets
  mxnsym=maximum number of symmetries
  mxntypat=maximum number of types of atoms
  mxnzchempot=maximal value of input nzchempot for all the datasets

SIDE EFFECTS

  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables,
   some of which are initialized here (see invars1.f for more details on the
   initialized records)

PARENTS

      m_ab7_invars_f90

CHILDREN

      indefo1,invars1

SOURCE

 67 #if defined HAVE_CONFIG_H
 68 #include "config.h"
 69 #endif
 70 
 71 #include "abi_common.h"
 72 
 73 subroutine invars1m(dmatpuflag,dtsets,iout,lenstr,mband_upper_,&
 74 & msym,mxga_n_rules,mxgw_nqlwl,mxlpawu,mxmband_upper,mxnatom,&
 75 & mxnatpawu,mxnatsph,mxnatsph_extra,mxnatvshift,mxnconeq,&
 76 & mxnimage,mxn_efmas_dirs,mxnkpt,mxnkptgw,mxnkpthf,mxnnos,mxnqptdm,mxnspinor, &
 77 & mxnsppol,mxnsym,mxntypat,mxnimfrqs,mxnfreqsp,mxnzchempot,&
 78 & mxn_projection_frequencies,ndtset,ndtset_alloc,string,npsp,zionpsp)
 79 
 80  use defs_basis
 81  use defs_abitypes
 82  use m_profiling_abi
 83  use m_xmpi
 84 
 85 !This section has been created automatically by the script Abilint (TD).
 86 !Do not modify the following lines by hand.
 87 #undef ABI_FUNC
 88 #define ABI_FUNC 'invars1m'
 89  use interfaces_57_iovars, except_this_one => invars1m
 90 !End of the abilint section
 91 
 92  implicit none
 93 
 94 !Arguments ------------------------------------
 95 !scalars
 96  integer,intent(in) :: iout,lenstr,msym,mxnatom,mxnimage,ndtset,ndtset_alloc,npsp
 97  integer,intent(out) :: dmatpuflag,mxga_n_rules,mxgw_nqlwl,mxlpawu,mxmband_upper,mxnatpawu
 98  integer,intent(out) :: mxnatsph, mxnatsph_extra
 99  integer,intent(out) :: mxnatvshift,mxnconeq,mxn_efmas_dirs,mxnkpt,mxnkptgw,mxnkpthf,mxnnos
100  integer,intent(out) :: mxnqptdm,mxnspinor,mxnsppol,mxnsym,mxntypat
101  integer,intent(out) :: mxnimfrqs,mxnfreqsp,mxnzchempot,mxn_projection_frequencies
102  character(len=*),intent(inout) :: string
103 !arrays
104  integer,intent(out) :: mband_upper_(0:ndtset_alloc)
105  type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc)
106  real(dp),intent(in) :: zionpsp(npsp)
107 
108 !Local variables-------------------------------
109 !scalars
110  integer :: idtset,ii,jdtset,lpawu,mband_upper,iatom,nat,nsp
111 !arrays
112  integer,allocatable :: symafm_(:,:),symrel_(:,:,:,:)
113  integer,allocatable :: symafm(:),symrel(:,:,:)
114  real(dp),allocatable :: tnons_(:,:,:),tnons(:,:)
115 
116 !******************************************************************
117 
118 !Here, allocation of the arrays that depend on msym.
119  ABI_ALLOCATE(symrel_,(3,3,msym,0:ndtset_alloc))
120  ABI_ALLOCATE(symafm_,(msym,0:ndtset_alloc))
121  ABI_ALLOCATE(tnons_,(3,msym,0:ndtset_alloc))
122  ABI_ALLOCATE(symafm,(msym))
123  ABI_ALLOCATE(symrel,(3,3,msym))
124  ABI_ALLOCATE(tnons,(3,msym))
125 
126 !Set up default values (note that the default acell, amu
127 !mkmem, mkmem1,mkqmem, and nkpt must be overcome
128 
129  do idtset=0,ndtset_alloc
130    call indefo1(dtsets(idtset))
131  end do
132 
133 !natom and nimage are already initialized in invars0
134  dtsets(0)%natom=-1
135  dtsets(0)%nimage=1
136 
137 !Initialization for parallelization data has changed
138 !these lines aim to keep old original default values
139  dtsets(0)%npimage=1
140  dtsets(0)%npkpt=1
141  dtsets(0)%npspinor=1
142  dtsets(0)%npfft=1
143  dtsets(0)%npband=1
144  dtsets(0)%bandpp=1
145 
146  symafm_(:,0)=1
147  symrel_(:,:,:,0)=0
148  symrel_(1,1,:,0)=1 ; symrel_(2,2,:,0)=1 ; symrel_(3,3,:,0)=1
149  tnons_(:,:,0)=0.0_dp
150 
151 !Loop on datasets
152  do idtset=1,ndtset_alloc
153    !write(std_out,'(2a,i0)') ch10,' invars1m : enter jdtset= ',jdtset
154    jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0
155 
156 !  Input default values
157    dtsets(idtset)%bravais(:)=0
158    symafm(:)=symafm_(:,0)
159    symrel(:,:,:)=symrel_(:,:,:,0)
160    tnons(:,:)=tnons_(:,:,0)
161 
162    call invars1(dtsets(idtset)%bravais,dtsets(idtset),iout,jdtset,lenstr,&
163 &   mband_upper,msym,npsp,string,symafm,symrel,tnons,zionpsp)
164 
165    mband_upper_ (idtset)=mband_upper
166    symafm_(:,idtset)=symafm(:)
167    symrel_(:,:,:,idtset)=symrel(:,:,:)
168    tnons_(:,:,idtset)=tnons(:,:)
169  end do
170 
171  mxmband_upper =maxval(mband_upper_ (1:ndtset_alloc))
172 
173  dmatpuflag=0;mxnatpawu=0;mxlpawu=0
174  mxnatsph=dtsets(1)%natsph
175  mxnatsph_extra=dtsets(1)%natsph_extra
176  mxnatvshift=dtsets(1)%natvshift
177  mxnconeq=dtsets(1)%nconeq
178  mxn_efmas_dirs=0
179  mxga_n_rules = dtsets(1)%ga_n_rules
180  mxgw_nqlwl = dtsets(1)%gw_nqlwl
181  mxnimfrqs = 0
182  mxnfreqsp = 0
183  mxn_projection_frequencies=0
184  mxnkpt  =dtsets(1)%nkpt
185  mxnkptgw=dtsets(1)%nkptgw
186  mxnkpthf=dtsets(1)%nkpthf
187  mxnnos  =dtsets(1)%nnos
188  mxnqptdm=dtsets(1)%nqptdm
189  mxnspinor=dtsets(1)%nspinor
190  mxnsppol=dtsets(1)%nsppol
191  mxntypat=dtsets(1)%ntypat
192  mxnzchempot=dtsets(1)%nzchempot
193 
194 !Get MAX dimension over datasets
195  do ii=1,ndtset_alloc
196    mxnatsph=max(dtsets(ii)%natsph,mxnatsph)
197    mxnatsph_extra=max(dtsets(ii)%natsph_extra,mxnatsph_extra)
198    mxnconeq=max(dtsets(ii)%nconeq,mxnconeq)
199    mxn_efmas_dirs=max(dtsets(ii)%efmas_n_dirs,mxn_efmas_dirs)
200    mxga_n_rules = max(dtsets(ii)%ga_n_rules,mxga_n_rules)
201    mxgw_nqlwl = max(dtsets(ii)%gw_nqlwl,mxgw_nqlwl)
202    mxnimfrqs = max(dtsets(ii)%cd_customnimfrqs,mxnimfrqs)
203    mxnfreqsp = max(dtsets(ii)%gw_customnfreqsp,mxnfreqsp)
204    mxn_projection_frequencies = max(dtsets(ii)%gwls_n_proj_freq,mxn_projection_frequencies)
205    mxnkpt  =max(dtsets(ii)%nkpt,mxnkpt)
206    mxnkptgw=max(dtsets(ii)%nkptgw,mxnkptgw)
207    mxnkpthf=max(dtsets(ii)%nkpthf,mxnkpthf)
208    mxnnos  =max(dtsets(ii)%nnos,mxnnos)
209    mxnqptdm=max(dtsets(ii)%nqptdm,mxnqptdm)
210    mxnspinor=max(dtsets(ii)%nspinor,mxnspinor)
211    mxnsppol=max(dtsets(ii)%nsppol,mxnsppol)
212    mxntypat=max(dtsets(ii)%ntypat,mxntypat)
213    mxnzchempot=max(dtsets(ii)%nzchempot,mxnzchempot)
214    if (dtsets(ii)%usepawu>0) then
215      if (dtsets(ii)%usedmatpu/=0) dmatpuflag=1
216      lpawu=maxval(dtsets(ii)%lpawu(:))
217      mxlpawu=max(lpawu,mxlpawu)
218      !dtsets(ii)%natpawu=count(dtsets(ii)%lpawu(dtsets(ii)%typat((/(i1,i1=1,dtsets(ii)%natom)/)))/=-1)
219      ! Old fashon way that should do fine
220      dtsets(ii)%natpawu = 0
221      do iatom=1, dtsets(ii)%natom
222        if (dtsets(ii)%lpawu(dtsets(ii)%typat(iatom)) /= -1 ) dtsets(ii)%natpawu = dtsets(ii)%natpawu + 1
223      end do
224      mxnatpawu=max(dtsets(ii)%natpawu,mxnatpawu)
225      if (dtsets(ii)%macro_uj/=0) dtsets(ii)%natvshift=lpawu*2+1
226    end if
227    mxnatvshift=max(dtsets(ii)%natvshift,mxnatvshift)
228  end do
229 
230 !mxnsym=maxval(dtsets(1:ndtset_alloc)%nsym) ! This might not work properly with HP compiler
231  mxnsym=dtsets(1)%nsym
232  do idtset=1,ndtset_alloc
233    mxnsym=max(dtsets(idtset)%nsym,mxnsym)
234  end do
235 
236  do idtset=0,ndtset_alloc
237    ABI_ALLOCATE(dtsets(idtset)%atvshift,(mxnatvshift,mxnsppol,mxnatom))
238    ABI_ALLOCATE(dtsets(idtset)%bs_loband,(mxnsppol))
239    ABI_ALLOCATE(dtsets(idtset)%bdgw,(2,mxnkptgw,mxnsppol))
240    ABI_ALLOCATE(dtsets(idtset)%cd_imfrqs,(mxnimfrqs))
241    ABI_ALLOCATE(dtsets(idtset)%chempot,(3,mxnzchempot,mxntypat))
242    nsp=max(mxnsppol,mxnspinor);nat=mxnatpawu*dmatpuflag
243    ABI_ALLOCATE(dtsets(idtset)%dmatpawu,(2*mxlpawu+1,2*mxlpawu+1,nsp,nat,mxnimage))
244    ABI_ALLOCATE(dtsets(idtset)%efmas_bands,(2,mxnkpt))
245    ABI_ALLOCATE(dtsets(idtset)%efmas_dirs,(3,mxn_efmas_dirs))
246    ABI_ALLOCATE(dtsets(idtset)%gw_freqsp,(mxnfreqsp))
247    ABI_ALLOCATE(dtsets(idtset)%gwls_list_proj_freq,(mxn_projection_frequencies))
248    ABI_ALLOCATE(dtsets(idtset)%gw_qlwl,(3,mxgw_nqlwl))
249    ABI_ALLOCATE(dtsets(idtset)%kpt,(3,mxnkpt))
250    ABI_ALLOCATE(dtsets(idtset)%kptgw,(3,mxnkptgw))
251    ABI_ALLOCATE(dtsets(idtset)%kptns,(3,mxnkpt))
252    ABI_ALLOCATE(dtsets(idtset)%kptns_hf,(3,mxnkpthf))
253    ABI_ALLOCATE(dtsets(idtset)%iatsph,(mxnatsph))
254    ABI_ALLOCATE(dtsets(idtset)%istwfk,(mxnkpt))
255    ABI_ALLOCATE(dtsets(idtset)%nband,(mxnkpt*mxnsppol))
256    ABI_ALLOCATE(dtsets(idtset)%occ_orig,(mxmband_upper*mxnkpt*mxnsppol))
257    ABI_ALLOCATE(dtsets(idtset)%qmass,(mxnnos))
258    ABI_ALLOCATE(dtsets(idtset)%qptdm,(3,mxnqptdm))
259    ABI_ALLOCATE(dtsets(idtset)%symafm,(mxnsym))
260    ABI_ALLOCATE(dtsets(idtset)%symrel,(3,3,mxnsym))
261    ABI_ALLOCATE(dtsets(idtset)%tnons,(3,mxnsym))
262    ABI_ALLOCATE(dtsets(idtset)%wtatcon,(3,mxnatom,mxnconeq))
263    ABI_ALLOCATE(dtsets(idtset)%wtk,(mxnkpt))
264    ABI_ALLOCATE(dtsets(idtset)%xredsph_extra,(3,mxnatsph_extra))
265    dtsets(idtset)%symrel(:,:,:)=symrel_(:,:,1:mxnsym,idtset)
266    dtsets(idtset)%symafm(:)    =symafm_(1:mxnsym,idtset)
267    dtsets(idtset)%tnons (:,:)  =tnons_ (:,1:mxnsym,idtset)
268  end do
269 
270  ABI_DEALLOCATE(symafm_)
271  ABI_DEALLOCATE(symrel_)
272  ABI_DEALLOCATE(tnons_)
273  ABI_DEALLOCATE(symafm)
274  ABI_DEALLOCATE(symrel)
275  ABI_DEALLOCATE(tnons)
276 
277 end subroutine invars1m