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-2022 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 .

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 MODULE m_results_respfn
24 
25  use defs_basis
26  use m_errors
27  use m_abicore
28  use m_dtset
29 
30  implicit none
31 
32  private
33 
34 ! public procedures.
35  public :: init_results_respfn
36  public :: destroy_results_respfn

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

SOURCE

148 subroutine destroy_results_respfn(results_respfn)
149 
150 !Arguments ------------------------------------
151 !arrays
152  type(results_respfn_type),intent(inout) :: results_respfn
153 
154 !************************************************************************
155 
156  ABI_SFREE(results_respfn%gam_eig2nkq)
157 
158 end subroutine destroy_results_respfn
159 
160 !----------------------------------------------------------------------
161 
162 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

SOURCE

 86 subroutine init_results_respfn(dtsets,ndtset_alloc,results_respfn)
 87 
 88 !Arguments ------------------------------------
 89 !scalars
 90  integer,intent(in) :: ndtset_alloc
 91 !arrays
 92  type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc)
 93  type(results_respfn_type),intent(inout) :: results_respfn
 94 !Local variables-------------------------------
 95 !scalars
 96  integer :: getgam_eig2nkq,idtset,idtset_gam
 97  character(len=500) :: message
 98 
 99 !************************************************************************
100 
101  results_respfn%gam_jdtset=0
102 
103 !Possibly find the value of gam_jdtset
104  do idtset=1,ndtset_alloc
105    getgam_eig2nkq=dtsets(idtset)%getgam_eig2nkq
106    if(getgam_eig2nkq/=0)then
107      do idtset_gam=1,ndtset_alloc
108        if(dtsets(idtset_gam)%jdtset==getgam_eig2nkq)exit
109      enddo
110      if(dtsets(idtset_gam)%optdriver/=RUNL_RESPFN)then
111        write(message, '(a,i5,a,i5,2a,i5,3a)' )&
112          'For jdtset=',dtsets(idtset)%jdtset,', getgam_eig2nkq=',getgam_eig2nkq,ch10,&
113          'However this dataset with idtset=',idtset_gam, ' is not a phonon calculation;',ch10,&
114          'Action : correct the value of getgam_eig2nkq for that dataset.'
115        ABI_ERROR(message)
116      endif
117      if(results_respfn%gam_jdtset==0)then
118        results_respfn%gam_jdtset=-getgam_eig2nkq ! Store a negative value, indicating that it is expected.
119      else if(results_respfn%gam_jdtset/=-getgam_eig2nkq) then
120        write(message, '(a,i5,2a,i5,2a)' )&
121          'results_respfn%gam_jdtset=',results_respfn%gam_jdtset,ch10,&
122          'dtsets(idtset)%getgam_eig2nkq=',getgam_eig2nkq,ch10,&
123          'So, it seems that two gamma q point calculations should be stored, while this is not yet allowed.'
124        ABI_BUG(message)
125      endif
126    endif
127  enddo
128 
129 end subroutine init_results_respfn

m_results_respfn/results_respfn_type [ Types ]

[ Top ] [ 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

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