TABLE OF CONTENTS


ABINIT/m_results_results [ Modules ]

[ Top ] [ Modules ]

NAME

  m_results_results

FUNCTION

  This module provides the definition of the results_respfn_type used
  to store results from RESPFN calculations.

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (MT,GG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

26 #if defined HAVE_CONFIG_H
27 #include "config.h"
28 #endif
29 
30 #include "abi_common.h"
31 
32 MODULE m_results_respfn
33 
34  use defs_basis
35  use defs_abitypes
36  use m_errors
37  use m_abicore
38 
39  implicit none
40 
41  private
42 
43 ! public procedures.
44  public :: init_results_respfn
45  public :: destroy_results_respfn

defs_datatypes/results_respfn_type [ Types ]

[ Top ] [ defs_datatypes ] [ Types ]

NAME

 results_respfn_type

FUNCTION

 This structured datatype contains a subset of the results of a RESPFN
 calculation, needed to compute the Diagonal Debye-Waller calculation.

SOURCE

58  type, public :: results_respfn_type
59 
60 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
61 ! declared in another part of ABINIT, that might need to take into account your modification.
62 
63 ! Integer scalar
64   integer :: gam_jdtset
65   ! jdtset if the results from a q=gamma wavevector calculation, with dataset jdtset have been stored
66   ! -jdtset if the results from a q=gamma wavevector calculation, with dataset jdtset should be stored
67   ! 0 if no q=gamma wavevector calculation should be stored
68 
69   real(dp), allocatable :: gam_eig2nkq(:,:,:,:,:,:,:)  ! one half second derivatives of the electronic eigenvalues
70   ! eig2nkq(2,mband*nsppol,nkpt,3,natom,3,natom)
71 
72  end type results_respfn_type

m_results_out/destroy_results_respfn [ Functions ]

[ Top ] [ m_results_out ] [ Functions ]

NAME

  destroy_results_respfn

FUNCTION

  Clean and destroy results_respfn datastructure

INPUTS

SIDE EFFECTS

  results_respfn(:)=<type(results_respfn_type)>=results_respfn datastructure

PARENTS

      driver

CHILDREN

SOURCE

176 subroutine destroy_results_respfn(results_respfn)
177 
178 
179 !This section has been created automatically by the script Abilint (TD).
180 !Do not modify the following lines by hand.
181 #undef ABI_FUNC
182 #define ABI_FUNC 'destroy_results_respfn'
183 !End of the abilint section
184 
185  implicit none
186 
187 !Arguments ------------------------------------
188 !arrays
189  type(results_respfn_type),intent(inout) :: results_respfn
190 !Local variables-------------------------------
191 !scalars
192 
193 !************************************************************************
194 
195  if(results_respfn%gam_jdtset>0)then
196 !DEBUG
197 !  write(std_out,*)' results_respfn%gam_jdtset=',results_respfn%gam_jdtset
198 !  write(std_out,*)' allocated(results_respfn%gam_eig2nkq)=',allocated(results_respfn%gam_eig2nkq)
199 !  call flush(6)
200 !ENDDEBUG
201    if(allocated(results_respfn%gam_eig2nkq))then
202      ABI_DEALLOCATE(results_respfn%gam_eig2nkq)
203    endif
204  endif
205 
206 end subroutine destroy_results_respfn
207 
208 !----------------------------------------------------------------------
209 
210 END MODULE m_results_respfn

m_results_out/init_results_respfn [ Functions ]

[ Top ] [ m_results_out ] [ Functions ]

NAME

  init_results_respfn

FUNCTION

  Init all scalars and nullify pointers in the structure.

INPUTS

  ndtset= number of datasets

SIDE EFFECTS

  results_respfn(:)=<type(results_respfn_type)>=results_respfn datastructure

PARENTS

      driver

CHILDREN

SOURCE

100 subroutine init_results_respfn(dtsets,ndtset_alloc,results_respfn)
101 
102 
103 !This section has been created automatically by the script Abilint (TD).
104 !Do not modify the following lines by hand.
105 #undef ABI_FUNC
106 #define ABI_FUNC 'init_results_respfn'
107 !End of the abilint section
108 
109  implicit none
110 
111 !Arguments ------------------------------------
112 !scalars
113  integer,intent(in) :: ndtset_alloc
114 !arrays
115  type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc)
116  type(results_respfn_type),intent(inout) :: results_respfn
117 !Local variables-------------------------------
118 !scalars
119  integer :: getgam_eig2nkq,idtset,idtset_gam
120  character(len=500) :: message
121 
122 !************************************************************************
123 
124  results_respfn%gam_jdtset=0
125 
126 !Possibly find the value of gam_jdtset
127  do idtset=1,ndtset_alloc
128    getgam_eig2nkq=dtsets(idtset)%getgam_eig2nkq
129    if(getgam_eig2nkq/=0)then
130      do idtset_gam=1,ndtset_alloc
131        if(dtsets(idtset_gam)%jdtset==getgam_eig2nkq)exit
132      enddo
133      if(dtsets(idtset_gam)%optdriver/=RUNL_RESPFN)then
134        write(message, '(a,i5,a,i5,2a,i5,3a)' )&
135 &        'For jdtset=',dtsets(idtset)%jdtset,', getgam_eig2nkq=',getgam_eig2nkq,ch10,&
136 &        'However this dataset with idtset=',idtset_gam, ' is not a phonon calculation;',ch10,&
137 &        'Action : correct the value of getgam_eig2nkq for that dataset.'
138        MSG_ERROR(message)
139      endif
140      if(results_respfn%gam_jdtset==0)then
141        results_respfn%gam_jdtset=-getgam_eig2nkq ! Store a negative value, indicating that it is expected.
142      else if(results_respfn%gam_jdtset/=-getgam_eig2nkq) then
143        write(message, '(a,i5,2a,i5,2a)' )&
144 &        'results_respfn%gam_jdtset=',results_respfn%gam_jdtset,ch10,&
145 &        'dtsets(idtset)%getgam_eig2nkq=',getgam_eig2nkq,ch10,&
146 &        'So, it seems that two gamma q point calculations should be stored, while this is not yet allowed.'
147        MSG_BUG(message)
148      endif
149    endif
150  enddo
151 
152 end subroutine init_results_respfn