TABLE OF CONTENTS


ABINIT/m_ddb_hdr [ Modules ]

[ Top ] [ Modules ]

NAME

  m_ddb_hdr

FUNCTION

  This module contains the declaration of data types and methods
  to handle the header of the DDB files.

COPYRIGHT

 Copyright (C) 2011-2018 ABINIT group (MJV, XG, MT, MM, MVeithen, MG, PB, JCC, GA)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

SOURCE

 20 #if defined HAVE_CONFIG_H
 21 #include "config.h"
 22 #endif
 23 
 24 #include "abi_common.h"
 25 
 26 MODULE m_ddb_hdr
 27 
 28  use defs_basis
 29  use defs_datatypes
 30  use defs_abitypes
 31  use m_errors
 32  use m_abicore
 33  use m_xmpi
 34 
 35  use m_copy,      only : alloc_copy
 36  use m_pawtab,    only : pawtab_type, pawtab_nullify, pawtab_free !, pawtab_copy
 37  use m_psps,      only : psps_copy, psps_free
 38  use m_io_tools,  only : open_file
 39  use m_copy,      only : alloc_copy
 40  use m_fstrings,  only : sjoin
 41 
 42  implicit none
 43 
 44  private
 45 
 46  public :: ddb_getdims      ! Open a DDB file and read basic dimensions and variables.
 47  public :: ioddb8_in        ! Temporary
 48  public :: psddb8           ! Temporary
 49 
 50  type,public :: ddb_hdr_type
 51 
 52    integer :: ddb_version   ! Version of the DDB file
 53 
 54    integer :: matom
 55    integer :: mband
 56    integer :: mkpt
 57    integer :: msym
 58    integer :: mtypat
 59 
 60    integer :: intxc
 61    integer :: iscf
 62    integer :: ixc
 63    integer :: natom
 64    integer :: nkpt
 65    integer :: nspden
 66    integer :: nspinor
 67    integer :: nsppol
 68    integer :: nsym
 69    integer :: ntypat
 70    integer :: occopt
 71    integer :: usepaw
 72 
 73    integer :: nblok         ! Number of blocks in the ddb
 74    integer :: mblktyp       ! Max block type
 75                             ! 0 = Total energy
 76                             ! 1 = 2nd derivatives (non-stat.)
 77                             ! 2 = 2nd derivatives (stationary)
 78                             ! 3 = 3rd derivatives
 79                             ! 4 = 1st derivatives
 80                             ! 5 = 2nd eigenvalue derivatives
 81                             !
 82    integer :: fullinit      ! Whether the full info on the pseudo is present
 83                             ! TODO rename this variable
 84 
 85    real(dp) :: dilatmx
 86    real(dp) :: ecut
 87    real(dp) :: ecutsm
 88    real(dp) :: kptnrm
 89    real(dp) :: pawecutdg
 90    real(dp) :: dfpt_sciss
 91    real(dp) :: tolwfr
 92    real(dp) :: tphysel
 93    real(dp) :: tsmear
 94 
 95    character(len=fnlen) :: dscrpt
 96 
 97    integer :: ngfft(18)
 98    real(dp) :: acell(3)
 99    real(dp) :: rprim(3,3)
100 
101    integer,allocatable :: nband(:)
102    ! nband(mkpt*nsppol)
103 
104    integer,allocatable :: symafm(:)
105    ! symafm(msym)
106 
107    integer,allocatable :: symrel(:,:,:)
108    ! symrel(3,3,msym)
109 
110    integer,allocatable :: typat(:)
111    ! typat(matom)
112 
113    real(dp),allocatable :: amu(:)
114    ! amu(mtypat)
115 
116    real(dp),allocatable :: kpt(:,:)
117    ! kpt(3,mkpt)
118 
119    real(dp),allocatable :: occ(:)
120    ! occ(mband*mkpt*nsppol)
121 
122    real(dp),allocatable :: spinat(:,:)
123    ! spinat(3,matom)
124 
125    real(dp),allocatable :: tnons(:,:)
126    ! tnons(3,msym)
127 
128    real(dp),allocatable :: wtk(:)
129    ! wtk(mkpt)
130 
131    real(dp),allocatable :: xred(:,:)
132    ! xred(3,matom)
133 
134    real(dp),allocatable :: zion(:)
135    ! zion(mtypat)
136 
137    real(dp),allocatable :: znucl(:)
138    ! znucl(mtypat)
139 
140    type(pawtab_type),allocatable :: pawtab(:)
141    ! pawtab(psps%ntypat*psps%usepaw)
142 
143    type(pseudopotential_type) :: psps
144 
145 
146  end type ddb_hdr_type
147 
148  public :: ddb_hdr_init            ! Construct object
149  public :: ddb_hdr_malloc          ! Allocate dynamic memory.
150  public :: ddb_hdr_free            ! Free dynamic memory.
151  public :: ddb_hdr_open_write      ! Open the DDB file and write the header.
152  public :: ddb_hdr_open_read       ! Open the DDB file and read the header.
153  public :: ddb_hdr_compare         ! Compare two DDB headers.
154 
155 CONTAINS  !===========================================================

m_ddb_hdr/chki8 [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 chki8

FUNCTION

 This small subroutine check the identity of inti and intt,
 who are integers, and eventually send a message and stop
 if they are found unequal

INPUTS

 inti=first integer
 intt=second integer
 character(len=6) name=name of the variable in the calling routine,
                       to be echoed

OUTPUT

  (only checking)

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

2984 subroutine chki8(inti,intt,name)
2985 
2986 
2987 !This section has been created automatically by the script Abilint (TD).
2988 !Do not modify the following lines by hand.
2989 #undef ABI_FUNC
2990 #define ABI_FUNC 'chki8'
2991 !End of the abilint section
2992 
2993  implicit none
2994 
2995 !Arguments -------------------------------
2996 !scalars
2997  integer,intent(in) :: inti,intt
2998  character(len=6),intent(in) :: name
2999 
3000 !Local variables-------------------------------
3001 !scalars
3002  character(len=500) :: message
3003 
3004 ! *********************************************************************
3005 
3006  if(inti/=intt) then
3007    write(message, '(a,a,a,a,a,i10,a,a,a,i10,a,a,a)' )&
3008    'Comparing integers for variable',name,'.',ch10,&
3009    'Value from input DDB is',inti,' and',ch10,&
3010    'from transfer DDB is',intt,'.',ch10,&
3011    'Action: check your DDBs.'
3012    MSG_ERROR(message)
3013  end if
3014 
3015  end subroutine chki8

m_ddb_hdr/chkr8 [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 chkr8

FUNCTION

 This small subroutine check the identity of reali and realt,
 who are integers, and eventually send a message and stop
 if they are found unequal by more than tol

INPUTS

 reali=first real number
 intt=second  real number
 character(len=6) name=name of the variable in the calling routine,
                       to be echoed
 tol=tolerance

OUTPUT

  (only checking)

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

2922 subroutine chkr8(reali,realt,name,tol)
2923 
2924 
2925 !This section has been created automatically by the script Abilint (TD).
2926 !Do not modify the following lines by hand.
2927 #undef ABI_FUNC
2928 #define ABI_FUNC 'chkr8'
2929 !End of the abilint section
2930 
2931  implicit none
2932 
2933 !Arguments -------------------------------
2934 !scalars
2935  real(dp),intent(in) :: reali,realt,tol
2936  character(len=6),intent(in) :: name
2937 
2938 !Local variables-------------------------------
2939 !scalars
2940  character(len=500) :: message
2941 
2942 ! *********************************************************************
2943 
2944  if(abs(reali-realt)>tol) then
2945    write(message, '(a,a,a,a,a,es16.6,a,a,a,es16.6,a,a,a)' )&
2946    'Comparing reals for variable',name,'.',ch10,&
2947    'Value from input DDB is',reali,' and',ch10,&
2948    'from transfer DDB is',realt,'.',ch10,&
2949    'Action: check your DDBs.'
2950    MSG_ERROR(message)
2951  end if
2952 
2953  end subroutine chkr8

m_ddb_hdr/compare_ddb_variables [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 compare_ddb_variables

FUNCTION

 Compare the temporary DDB and input DDB preliminary information,
 as well as psp information.

COPYRIGHT

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

NOTES

 1. All the variables have their usual meaning.
 2. Note that fullinit==0  means that the input DDB has been
 initialized by a ground state input file. Some comparison are
 then not required.
 3. All variables with 8 appended are from the new input DDB

INPUTS

  acell, acell8 = lattice parameters
  amu, amu8 = atomic masses
  dimekb = dimension of KB projector set (only used for NCPP)
  ecut, ecut8 = cutoff energy
  ekb, ekb8 = KB energies for pseudopotentials
  fullinit, fullmrgddb_init = flags (see notes)
  iscf, iscf8 = SCF algorithm
  ixc, ixc8 = XC functional
  kpt, kpt8 = kpoint array
  kptnrm, kptnrm8 = normalization factor for kpt
  natom, natom8 = number of atoms
  nband, nband8 = number of bands at each kpt
  ngfft, ngfft8 = FFT grid sizes
  nkpt, nkpt8 = number of kpoints
  nsppol, nsppol8 = number of spin polarization (1 or 2)
  nsym, nsym8 = number of symmetry operations
  ntypat, ntypat8 = number of types of atoms
  occ, occ8 = occupation numbers
  occopt, occop8 = occupation style (metal, insulator, smearing...)
  pawecutdg,pawecutdg8= cutoff energy used for the fine "double grid" (PAW only)
  pawtab,pawtab8= PAW tabulated data (PAW dataset)
  rprim, rprim8 = primitive vectors of unit cell (cartesian coordinates)
  dfpt_sciss, dfpt_sciss8 = scissor correction (Ha)
  symrel, symrel8 = symmetry operations in reciprocal space
  tnons, tnons8 = translations associated to symrel
  tolwfr, tolwfr8 = tolerance on convergence of wavefunctions
  typat, typat8 = array of atom types
  usepaw = flag for utilization of PAW
  wtk, wtk8 = weights of kpoints
  xred, xred8 = reduced coordinates of atoms
  zion, zion8 = ionic charges of nuclei

 OUTPUT (corresponding values, checked and/or set)
  acell, amu, dimekb, ecut, ekb, fullinit, iscf, ixc, kpt, kptnrm,
  natom, nband, ngfft, nkpt, nsppol, nsym, ntypat, occ, occopt,
  rprim, dfpt_sciss, symrel, tnons, tolwfr, typat, usepaw, wtk, xred, zion

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

2574 subroutine compare_ddb_variables(&
2575 & matom, matom8, mtypat, mtypat8, mkpt, mkpt8,&
2576 & mband, mband8, msym, msym8,&
2577 & acell,acell8,amu,amu8,dimekb,ecut,ecut8,ekb,ekb8,&
2578 & fullinit,fullmrgddb_init,iscf,iscf8,ixc,ixc8,kpt,kpt8,&
2579 & kptnrm,kptnrm8,&
2580 & natom,natom8,nband,nband8,ngfft,ngfft8,nkpt,nkpt8,&
2581 & nsppol,nsppol8,nsym,nsym8,ntypat,ntypat8,occ,occ8,&
2582 & occopt,occop8,pawecutdg,pawecutdg8,pawtab,pawtab8,&
2583 & rprim,rprim8,dfpt_sciss,dfpt_sciss8,symrel,symrel8,&
2584 & tnons,tnons8,tolwfr,tolwfr8,typat,typat8,usepaw,wtk,wtk8,&
2585 & xred,xred8,zion,zion8)
2586 
2587 
2588 !This section has been created automatically by the script Abilint (TD).
2589 !Do not modify the following lines by hand.
2590 #undef ABI_FUNC
2591 #define ABI_FUNC 'compare_ddb_variables'
2592 !End of the abilint section
2593 
2594  implicit none
2595 
2596 !Arguments -------------------------------
2597 !scalars
2598  integer,intent(in) :: matom, matom8, mtypat, mtypat8, mkpt, mkpt8
2599  integer,intent(in) :: mband, mband8, msym, msym8
2600  integer,intent(in) :: dimekb,fullmrgddb_init,iscf8,ixc8,natom8,nkpt8,nsppol8
2601  integer,intent(in) :: nsym8,ntypat8,occop8,usepaw
2602  integer,intent(inout) :: fullinit,iscf,ixc,natom,nkpt,nsppol,nsym,ntypat
2603  integer,intent(inout) :: occopt
2604  real(dp),intent(in) :: ecut8,kptnrm8,pawecutdg8,dfpt_sciss8,tolwfr8
2605  real(dp),intent(inout) :: ecut,kptnrm,pawecutdg,dfpt_sciss,tolwfr
2606 !arrays
2607  integer,intent(in) :: nband8(mkpt8*nsppol8),ngfft8(18)
2608  integer,intent(inout) :: nband(mkpt*nsppol),ngfft(18)
2609  integer,intent(in) :: symrel8(3,3,msym8),typat8(matom8)
2610  integer,intent(inout) :: symrel(3,3,msym),typat(matom)
2611  real(dp),intent(in) :: acell8(3),amu8(mtypat8),ekb8(dimekb,mtypat8)
2612  real(dp),intent(inout) :: acell(3),amu(mtypat),ekb(dimekb,mtypat)
2613  real(dp),intent(in) :: kpt8(3,mkpt8),occ8(mband8*mkpt8*nsppol8)
2614  real(dp),intent(inout) :: kpt(3,mkpt),occ(mband*mkpt*nsppol)
2615  real(dp),intent(in) :: rprim8(3,3),tnons8(3,msym8)
2616  real(dp),intent(inout) :: rprim(3,3),tnons(3,msym)
2617  real(dp),intent(in) :: wtk8(mkpt8),xred8(3,matom8),zion8(mtypat8)
2618  real(dp),intent(inout) :: wtk(mkpt),xred(3,matom),zion(mtypat)
2619  type(pawtab_type),intent(in) :: pawtab8(ntypat8*usepaw)
2620  type(pawtab_type),intent(inout) :: pawtab(ntypat*usepaw)
2621 
2622 !Local variables -------------------------
2623 !scalars
2624  integer :: bantot,ii,ij,isym,itypat
2625  real(dp) :: ekbcm8,ekbcmp
2626  real(dp),parameter :: tol=2.0d-14
2627  character(len=500) :: msg
2628 
2629 ! *********************************************************************
2630 
2631 
2632 !Compare all the preliminary information
2633 !1. natom
2634  call chki8(natom,natom8,' natom')
2635 !2. nkpt
2636 !Compares the input and transfer values only if the input has not
2637 !been initialized by a ground state input file
2638 !There can also be the case of perturbation at Gamma, that
2639 !only need half of the number of k points.
2640  if(fullinit/=0)then
2641    if(nkpt/=2*nkpt8 .and. 2*nkpt/=nkpt8)then
2642 
2643      ! GKA: We don't always need this variable to be consistent
2644      !      For example, we might have reduced the number of k-points
2645      !      with TRS only for certain q-points.
2646      !call chki8(nkpt,nkpt8,'  nkpt')
2647    else
2648      write(std_out,*)' compar8 : assume that one of the DDB to be',&
2649 &     ' merged use Time-Reversal to'
2650      write(std_out,*)' decrease the number of k-points'
2651    end if
2652  else
2653 !  Otherwise, takes the meaningful value
2654    nkpt=nkpt8
2655  end if
2656 !3a. occopt
2657 !Because the program will stop if the bloks
2658 !do not compare well, take here the most favorable case.
2659  if(occop8==0)occopt=0
2660 !3b. nband
2661 !Compares the input and transfer values only if the input has not
2662 !been initialized by a ground state input file
2663 !There can also be the case of perturbation at Gamma, that
2664 !only need half of the number of k points.
2665  if(fullinit==0 .or. nkpt8==2*nkpt)then
2666    bantot=0
2667    do ii=1,nkpt8
2668      nband(ii)=nband8(ii)
2669      bantot=bantot+nband(ii)
2670    end do
2671  else
2672    bantot=0
2673    do ii=1,nkpt
2674      if(nkpt==nkpt8)then
2675        call chki8(nband(ii),nband8(ii),' nband')
2676      end if
2677      bantot=bantot+nband(ii)
2678    end do
2679  end if
2680 !9. nsppol
2681  call chki8(nsppol,nsppol8,'nsppol')
2682 !4. nsym
2683  if(nsym/=1 .and. nsym8/=1)then
2684    call chki8(nsym,nsym8,'  nsym')
2685  end if
2686 !5. ntypat
2687  call chki8(ntypat,ntypat8,'ntypat')
2688 !6. acell
2689  do ii=1,3
2690    call chkr8(acell(ii),acell8(ii),' acell',tol)
2691  end do
2692 !7. amu
2693  do ii=1,ntypat
2694    call chkr8(amu(ii),amu8(ii),'   amu',tol)
2695  end do
2696 !9. date
2697 !10. ecut
2698  call chkr8(ecut,ecut8,'  ecut',tol)
2699 !10b. pawecutdg (PAW only)
2700  if (usepaw==1) then
2701    call chkr8(pawecutdg,pawecutdg8,'  ecut',tol)
2702  end if
2703 !11. iscf
2704 !Compares the input and transfer values only if the input has not
2705 !been initialized by a ground state input file
2706  if(fullinit/=0)then
2707    call chki8(iscf,iscf8,'  iscf')
2708  else
2709 !  Otherwise, takes the meaningful value
2710    iscf=iscf8
2711  end if
2712 !12. ixc
2713  call chki8(ixc,ixc8,'   ixc')
2714 !13. kpt and 14. kptnrm
2715 !Compares the input and transfer values only if the input
2716 !has not been initialized by a ground state input file
2717 !and if the number of k points is identical
2718  if(nkpt8 == 2*nkpt .or. fullinit==0)then
2719 !  Copy the largest number of k points in the right place
2720    do ij=1,nkpt8
2721      do ii=1,3
2722        kpt(ii,ij)=kpt8(ii,ij)
2723      end do
2724    end do
2725    kptnrm=kptnrm8
2726  else if (nkpt==nkpt8)then
2727    do ij=1,nkpt
2728      do ii=1,3
2729 !      Compares the input and transfer values only if the input
2730 !      has not been initialized by a ground state input file
2731        call chkr8(kpt(ii,ij)/kptnrm,kpt8(ii,ij)/kptnrm8,'   kpt',tol)
2732      end do
2733    end do
2734  end if
2735 !16. ngfft
2736 !MT dec 2013: deactivate the stop on ngfft to allow for
2737 ! (nfft-converged) DFPT calculations with GS WFK obtained with a different ngfft
2738  do ii=1,3
2739    if (ngfft(ii) == ngfft8(ii)) cycle
2740    write(msg,'(3a,i10,3a,i10,a)') &
2741 &   'Comparing integers for variable ngfft.',ch10,&
2742 &   'Value from input DDB is',ngfft(ii),' and',ch10,&
2743 &   'from transfer DDB is',ngfft8(ii),'.'
2744    MSG_WARNING(msg)
2745  end do
2746 !17. occ
2747 !Compares the input and transfer values only if the input has not
2748 !been inititialized by a ground state input file
2749  do ii=1,bantot
2750    if (fullinit==0 .or. nkpt8==2*nkpt) then
2751      occ(ii)=occ8(ii)
2752    else if(nkpt==nkpt8)then
2753      call chkr8(occ(ii),occ8(ii),'   occ',tol)
2754    end if
2755  end do
2756 !18. rprim
2757  do ii=1,3
2758    do ij=1,3
2759      call chkr8(rprim(ii,ij),rprim8(ii,ij),' rprim',tol)
2760    end do
2761  end do
2762 !19. dfpt_sciss
2763 !Compares the input and transfer values only if the input has not
2764 !been inititialized by a ground state input file
2765  if(fullinit/=0)then
2766    call chkr8(dfpt_sciss,dfpt_sciss8,' dfpt_sciss',tol)
2767  else
2768 !  Otherwise, takes the meaningful value
2769    dfpt_sciss=dfpt_sciss8
2770  end if
2771 !20. symrel
2772 !If nsym == nsym8, compares the symmetry operations,
2773 !otherwise, one of nsym or nsym8 is 1, and thus take the
2774 !symrel corresponding to the largest set.
2775 !nsym will be changed later
2776  if(nsym==nsym8)then
2777    do isym=1,nsym
2778      do ii=1,3
2779        do ij=1,3
2780          call chki8(symrel(ii,ij,isym),symrel8(ii,ij,isym),'symrel')
2781        end do
2782      end do
2783    end do
2784  else if(nsym8/=1)then
2785    symrel(:,:,1:nsym8)=symrel8(:,:,1:nsym8)
2786  end if
2787 !21. tnons (see symrel)
2788  if(nsym==nsym8)then
2789    do isym=1,nsym
2790      do ii=1,3
2791        call chkr8(tnons(ii,isym),tnons8(ii,isym),' tnons',tol)
2792      end do
2793    end do
2794  else if(nsym8/=1)then
2795    tnons(:,1:nsym8)=tnons8(:,1:nsym8)
2796    nsym=nsym8
2797  end if
2798 !22. tolwfr
2799 !Take the less converged value...
2800  tolwfr=max(tolwfr,tolwfr8)
2801 !23. typat
2802  do ii=1,ntypat
2803    call chki8(typat(ii),typat8(ii),' typat')
2804  end do
2805 !24. wtk
2806 !Compares the input and transfer values only if the input has not
2807 !been initialized by a ground state input file and the
2808 !number of k-points is identical.
2809  if(nkpt8==2*nkpt .or. fullinit==0)then
2810    do ii=1,nkpt8
2811      wtk(ii)=wtk8(ii)
2812    end do
2813  else if(nkpt==nkpt8)then
2814    do ii=1,nkpt
2815      call chkr8(wtk(ii),wtk8(ii),'   wtk',tol)
2816    end do
2817  end if
2818 !25.xred
2819  do ij=1,natom
2820    do ii=1,3
2821      call chkr8(xred(ii,ij),xred8(ii,ij),'  xred',tol)
2822    end do
2823  end do
2824 !26. zion
2825  do ii=1,ntypat
2826    call chkr8(zion(ii),zion8(ii),'  zion',tol)
2827  end do
2828 
2829 !Finally, put the correct value of nkpt in the case
2830 !of the use of the time-reversal symmetry
2831  if(2*nkpt==nkpt8)then
2832    nkpt=nkpt8
2833  end if
2834 
2835 !Now compare the NC pseudopotential information
2836  if (usepaw==0) then
2837    if(dimekb/=0 .and. fullinit/=0 .and. fullmrgddb_init/=0 )then
2838      do ii=1,dimekb
2839        do itypat=1,ntypat
2840          ekbcmp=ekb(ii,itypat)
2841          ekbcm8=ekb8(ii,itypat)
2842          call chkr8(ekbcmp,ekbcm8,'   ekb',tol)
2843        end do
2844      end do
2845    else if(dimekb/=0 .and. fullmrgddb_init/=0)then
2846      do ii=1,dimekb
2847        do itypat=1,ntypat
2848          ekb(ii,itypat)=ekb8(ii,itypat)
2849        end do
2850      end do
2851    end if
2852  end if
2853 
2854 !Now compare several PAW dataset information
2855  if (usepaw==1) then
2856    if (fullinit/=0 .and. fullmrgddb_init/=0) then
2857      do itypat=1,ntypat
2858        call chki8(pawtab(itypat)%basis_size,pawtab8(itypat)%basis_size,'bas_sz')
2859        call chki8(pawtab(itypat)%lmn_size,pawtab8(itypat)%lmn_size,'lmn_sz')
2860        call chki8(pawtab(itypat)%lmn2_size,pawtab8(itypat)%lmn2_size,'lmn2sz')
2861        call chkr8(pawtab(itypat)%rpaw,pawtab8(itypat)%rpaw,'  rpaw',tol3)
2862        call chkr8(pawtab(itypat)%rshp,pawtab8(itypat)%rshp,'rshape',tol3)
2863        call chki8(pawtab(itypat)%shape_type,pawtab8(itypat)%shape_type,'shp_tp')
2864        if (pawtab(itypat)%lmn2_size>0) then
2865          do ii=1,pawtab(itypat)%lmn2_size
2866            call chkr8(pawtab(itypat)%dij0(ii),pawtab8(itypat)%dij0(ii),'  dij0',tol)
2867          end do
2868        end if
2869      end do
2870    else if (fullmrgddb_init/=0) then
2871      do itypat=1,ntypat
2872        pawtab(itypat)%basis_size =pawtab8(itypat)%basis_size
2873        pawtab(itypat)%lmn_size   =pawtab8(itypat)%lmn_size
2874        pawtab(itypat)%rpaw       =pawtab8(itypat)%rpaw
2875        pawtab(itypat)%rshp       =pawtab8(itypat)%rshp
2876        pawtab(itypat)%shape_type =pawtab8(itypat)%shape_type
2877        if (pawtab8(itypat)%lmn2_size>0) then
2878          if (pawtab(itypat)%lmn2_size==0)  then
2879            ABI_ALLOCATE(pawtab(itypat)%dij0,(pawtab8(itypat)%lmn2_size))
2880          end if
2881          do ii=1,pawtab8(itypat)%lmn2_size
2882            pawtab(itypat)%dij0(ii)=pawtab8(itypat)%dij0(ii)
2883          end do
2884        end if
2885        pawtab(itypat)%lmn2_size  =pawtab8(itypat)%lmn2_size
2886      end do
2887    end if
2888  end if
2889 
2890 end subroutine compare_ddb_variables

m_ddb_hdr/ddb_chkname [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_chkname

FUNCTION

 This small subroutine check the identity of its argument,
 who are a6 names, and eventually send a message and stop
 if they are found unequal

INPUTS

 nmfond= name which has to be checked
 nmxpct= name expected for nmfond
 nmxpct2= eventual second optional name (backward compatibility)

OUTPUT

TODO

 Describe the inputs

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

2457 subroutine ddb_chkname(nmfond,nmxpct,nmxpct2)
2458 
2459 
2460 !This section has been created automatically by the script Abilint (TD).
2461 !Do not modify the following lines by hand.
2462 #undef ABI_FUNC
2463 #define ABI_FUNC 'ddb_chkname'
2464 !End of the abilint section
2465 
2466  implicit none
2467 
2468 !Arguments -------------------------------
2469 !scalars
2470  character(len=*),intent(in) :: nmfond,nmxpct
2471  character(len=*),intent(in),optional :: nmxpct2
2472 
2473 !Local variables-------------------------------
2474 !scalars
2475  logical :: found
2476  character(len=500) :: nmfond_,nmxpct_,nmxpct2_
2477  character(len=500) :: message
2478 
2479 ! *********************************************************************
2480 
2481  nmxpct_ = trim(adjustl(nmxpct))
2482  nmfond_ = trim(adjustl(nmfond))
2483 
2484  found = (nmxpct_ == nmfond_)
2485 
2486  if (present(nmxpct2) .and. .not. found) then
2487    nmxpct2_ = trim(adjustl(nmxpct2))
2488    found = (nmxpct2_==nmfond_)
2489  end if
2490 
2491  if (.not. found) then
2492    write(message, '(a,a,a,a,a,a,a,a,a,a,a)' )&
2493 &   'Reading DDB, expected name was "',trim(nmxpct_),'"',ch10,&
2494 &   '             and name found is "',trim(nmfond_),'"',ch10,&
2495 &   'Likely your DDB is incorrect.',ch10,&
2496 &   'Action: correct your DDB, or contact the ABINIT group.'
2497    MSG_ERROR(message)
2498  end if
2499 
2500 end subroutine ddb_chkname

m_ddb_hdr/ddb_getdims [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_getdims

FUNCTION

 Open Derivative DataBase, then reads the variables that
 must be known in order to dimension the arrays before complete reading

INPUTS

 character(len=*) filnam: name of input or output file
 unddb=unit number for input or output
 vrsddb=6 digit integer giving date, in form yymmdd for month=mm(1-12),
  day=dd(1-31), and year=yy(90-99 for 1990 to 1999,00-89 for 2000 to 2089),
  of current DDB version.

OUTPUT

 dimekb=dimension of ekb (only used for norm-conserving psps)
 lmnmax=if useylm=1, max number of (l,m,n) comp. over all type of psps
       =if useylm=0, max number of (l,n)   comp. over all type of psps
 mband=maximum number of bands
 mblktyp=largest block type
 msym=maximum number of symmetries
 natom=number of atoms
 nblok=number of bloks in the DDB
 nkpt=number of k points
 ntypat=number of atom types
 usepaw= 0 for non paw calculation; =1 for paw calculation
 comm=MPI communicator.

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

1867 subroutine ddb_getdims(dimekb,filnam,lmnmax,mband,mblktyp,msym,natom,nblok,nkpt,ntypat,unddb,usepaw,vrsddb,comm)
1868 
1869 
1870 !This section has been created automatically by the script Abilint (TD).
1871 !Do not modify the following lines by hand.
1872 #undef ABI_FUNC
1873 #define ABI_FUNC 'ddb_getdims'
1874 !End of the abilint section
1875 
1876  implicit none
1877 
1878 !Arguments -------------------------------
1879 !scalars
1880  integer,intent(in) :: unddb,vrsddb,comm
1881  integer,intent(out) :: msym,dimekb,lmnmax,mband,mblktyp,natom,nblok,nkpt,ntypat,usepaw
1882  character(len=*),intent(in) :: filnam
1883 
1884 !Local variables-------------------------------
1885 !scalars
1886  integer,parameter :: master=0
1887  integer :: ierr
1888  !integer :: mpert,msize
1889 
1890 ! *********************************************************************
1891 
1892  ! Master node reads dims from file and then broadcast.
1893  if (xmpi_comm_rank(comm) == master) then
1894    call inprep8(dimekb,filnam,lmnmax,mband,mblktyp,msym,natom,nblok,nkpt,ntypat,unddb,usepaw,vrsddb)
1895  end if
1896 
1897  if (xmpi_comm_size(comm) > 1) then
1898    call xmpi_bcast(dimekb, master, comm, ierr)
1899    call xmpi_bcast(lmnmax, master, comm, ierr)
1900    call xmpi_bcast(mband, master, comm, ierr)
1901    call xmpi_bcast(mblktyp, master, comm, ierr)
1902    call xmpi_bcast(msym, master, comm, ierr)
1903    call xmpi_bcast(natom, master, comm, ierr)
1904    call xmpi_bcast(nblok, master, comm, ierr)
1905    call xmpi_bcast(nkpt, master, comm, ierr)
1906    call xmpi_bcast(ntypat, master, comm, ierr)
1907    call xmpi_bcast(usepaw, master, comm, ierr)
1908  end if
1909 
1910  ! Maximum number of perturbations and size of matrix.
1911  !mpert=natom+6
1912  !msize=3*mpert*3*mpert; if (mblktyp==3) msize=msize*3*mpert
1913 
1914 end subroutine ddb_getdims

m_ddb_hdr/ddb_hdr_compare [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_hdr_compare

FUNCTION

  Compare two DDB headers and raise error if they differ.
  Also, complete psps information if one has more info than the other.

INPUTS

OUTPUT

PARENTS

      mblktyp1,mblktyp5

CHILDREN

SOURCE

690 subroutine ddb_hdr_compare(ddb_hdr1, ddb_hdr2)
691 
692 
693 !This section has been created automatically by the script Abilint (TD).
694 !Do not modify the following lines by hand.
695 #undef ABI_FUNC
696 #define ABI_FUNC 'ddb_hdr_compare'
697 !End of the abilint section
698 
699  implicit none
700 
701 !Arguments ------------------------------------
702  type(ddb_hdr_type),intent(inout) :: ddb_hdr1, ddb_hdr2
703 
704 !Local variables -------------------------
705 
706 ! ************************************************************************
707 
708 !    Should also compare indlmn and pspso ... but suppose that
709 !    the checking of ekb is enough for the psps.
710 !    Should also compare many other variables ... this is still
711 !    to be done ...
712  call compare_ddb_variables(&
713 &     ddb_hdr1%matom, ddb_hdr2%matom,&
714 &     ddb_hdr1%mtypat, ddb_hdr2%mtypat,&
715 &     ddb_hdr1%mkpt, ddb_hdr2%mkpt,&
716 &     ddb_hdr1%mband, ddb_hdr2%mband,&
717 &     ddb_hdr1%msym, ddb_hdr2%msym,&
718 &     ddb_hdr1%acell, ddb_hdr2%acell,&
719 &     ddb_hdr1%amu, ddb_hdr2%amu,&
720 &     ddb_hdr1%psps%dimekb,&
721 &     ddb_hdr1%ecut, ddb_hdr2%ecut,&
722 &     ddb_hdr1%psps%ekb, ddb_hdr2%psps%ekb,&
723 &     ddb_hdr1%fullinit, ddb_hdr2%fullinit,&
724 &     ddb_hdr1%iscf, ddb_hdr2%iscf,&
725 &     ddb_hdr1%ixc, ddb_hdr2%ixc,&
726 &     ddb_hdr1%kpt, ddb_hdr2%kpt,&
727 &     ddb_hdr1%kptnrm, ddb_hdr2%kptnrm,&
728 &     ddb_hdr1%natom, ddb_hdr2%natom,&
729 &     ddb_hdr1%nband, ddb_hdr2%nband,&
730 &     ddb_hdr1%ngfft, ddb_hdr2%ngfft,&
731 &     ddb_hdr1%nkpt, ddb_hdr2%nkpt,&
732 &     ddb_hdr1%nsppol, ddb_hdr2%nsppol,&
733 &     ddb_hdr1%nsym, ddb_hdr2%nsym,&
734 &     ddb_hdr1%ntypat, ddb_hdr2%ntypat,&
735 &     ddb_hdr1%occ, ddb_hdr2%occ,&
736 &     ddb_hdr1%occopt, ddb_hdr2%occopt,&
737 &     ddb_hdr1%pawecutdg, ddb_hdr2%pawecutdg,&
738 &     ddb_hdr1%pawtab, ddb_hdr2%pawtab,&
739 &     ddb_hdr1%rprim, ddb_hdr2%rprim,&
740 &     ddb_hdr1%dfpt_sciss, ddb_hdr2%dfpt_sciss,&
741 &     ddb_hdr1%symrel, ddb_hdr2%symrel,&
742 &     ddb_hdr1%tnons, ddb_hdr2%tnons,&
743 &     ddb_hdr1%tolwfr, ddb_hdr2%tolwfr,&
744 &     ddb_hdr1%typat, ddb_hdr2%typat,&
745 &     ddb_hdr1%usepaw,&
746 &     ddb_hdr1%wtk, ddb_hdr2%wtk,&
747 &     ddb_hdr1%xred, ddb_hdr2%xred,&
748 &     ddb_hdr1%zion, ddb_hdr2%zion)
749 
750 end subroutine ddb_hdr_compare

m_ddb_hdr/ddb_hdr_free [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_hdr_free

FUNCTION

INPUTS

OUTPUT

PARENTS

      anaddb,dfpt_looppert,eig2tot,gstate,m_ddb,m_effective_potential_file
      m_gruneisen,mblktyp1,mblktyp5,mrgddb,nonlinear,respfn,thmeig

CHILDREN

SOURCE

386 subroutine ddb_hdr_free(ddb_hdr)
387 
388 
389 !This section has been created automatically by the script Abilint (TD).
390 !Do not modify the following lines by hand.
391 #undef ABI_FUNC
392 #define ABI_FUNC 'ddb_hdr_free'
393 !End of the abilint section
394 
395  implicit none
396 
397 !Arguments ------------------------------------
398  type(ddb_hdr_type),intent(inout) :: ddb_hdr
399 
400 ! ************************************************************************
401 
402  ! integer
403  if (allocated(ddb_hdr%nband)) then
404    ABI_FREE(ddb_hdr%nband)
405  end if
406  if (allocated(ddb_hdr%symafm)) then
407    ABI_FREE(ddb_hdr%symafm)
408  end if
409  if (allocated(ddb_hdr%symrel)) then
410    ABI_FREE(ddb_hdr%symrel)
411  end if
412  if (allocated(ddb_hdr%typat)) then
413    ABI_FREE(ddb_hdr%typat)
414  end if
415 
416  ! real
417  if (allocated(ddb_hdr%amu)) then
418    ABI_FREE(ddb_hdr%amu)
419  end if
420  if (allocated(ddb_hdr%kpt)) then
421    ABI_FREE(ddb_hdr%kpt)
422  end if
423  if (allocated(ddb_hdr%occ)) then
424    ABI_FREE(ddb_hdr%occ)
425  end if
426  if (allocated(ddb_hdr%spinat)) then
427    ABI_FREE(ddb_hdr%spinat)
428  end if
429  if (allocated(ddb_hdr%tnons)) then
430    ABI_FREE(ddb_hdr%tnons)
431  end if
432  if (allocated(ddb_hdr%wtk)) then
433    ABI_FREE(ddb_hdr%wtk)
434  end if
435  if (allocated(ddb_hdr%xred)) then
436    ABI_FREE(ddb_hdr%xred)
437  end if
438  if (allocated(ddb_hdr%zion)) then
439    ABI_FREE(ddb_hdr%zion)
440  end if
441  if (allocated(ddb_hdr%znucl)) then
442    ABI_FREE(ddb_hdr%znucl)
443  end if
444 
445  ! types
446  call psps_free(ddb_hdr%psps)
447 
448  if (allocated(ddb_hdr%pawtab)) then
449    call pawtab_free(ddb_hdr%pawtab)
450    ABI_DATATYPE_DEALLOCATE(ddb_hdr%pawtab)
451  end if
452 
453 end subroutine ddb_hdr_free

m_ddb_hdr/ddb_hdr_init [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_hdr_init

FUNCTION

   Initialize a ddb_hdr object from a dataset.

INPUTS

OUTPUT

PARENTS

      dfpt_looppert,eig2tot,gstate,nonlinear,respfn

CHILDREN

SOURCE

178 subroutine ddb_hdr_init(ddb_hdr, dtset, psps, pawtab, ddb_version, dscrpt, &
179                         nblok, xred, occ, ngfft)
180 
181 
182 !This section has been created automatically by the script Abilint (TD).
183 !Do not modify the following lines by hand.
184 #undef ABI_FUNC
185 #define ABI_FUNC 'ddb_hdr_init'
186 !End of the abilint section
187 
188  implicit none
189 
190 !Arguments ------------------------------------
191  type(ddb_hdr_type),intent(out) :: ddb_hdr
192  type(dataset_type),intent(in) :: dtset
193  type(pseudopotential_type),intent(in) :: psps
194  type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw)
195  character(len=*),intent(in) :: dscrpt
196  integer,intent(in) :: ddb_version
197  integer,intent(in) :: nblok
198  integer,intent(in),optional :: ngfft(18)
199  real(dp),intent(in),optional :: xred(3,dtset%natom)
200  real(dp),intent(in),optional :: occ(dtset%mband*dtset%nkpt*dtset%nsppol)
201 
202 !Local variables -------------------------
203  integer :: ii, nn
204 
205 
206 ! ************************************************************************
207 
208  ddb_hdr%nblok = nblok
209  ddb_hdr%ddb_version = ddb_version
210  ddb_hdr%dscrpt = trim(dscrpt)
211  if (present(ngfft)) then
212    ddb_hdr%ngfft = ngfft
213  else
214    ddb_hdr%ngfft = dtset%ngfft
215  end if
216 
217  call psps_copy(psps, ddb_hdr%psps)
218 
219  ! Copy scalars from dtset
220  ddb_hdr%matom = dtset%natom
221  ddb_hdr%natom = dtset%natom
222  ddb_hdr%mband = dtset%mband
223  !ddb_hdr%mkpt = dtset%maxnkpt
224  ddb_hdr%mkpt = dtset%nkpt
225  ddb_hdr%nkpt = dtset%nkpt
226  !ddb_hdr%msym = dtset%maxnsym
227  ddb_hdr%msym = dtset%nsym
228  ddb_hdr%nsym = dtset%nsym
229  ddb_hdr%mtypat = dtset%ntypat
230  ddb_hdr%ntypat = dtset%ntypat
231 
232  ddb_hdr%nspden = dtset%nspden
233  ddb_hdr%nspinor = dtset%nspinor
234  ddb_hdr%nsppol = dtset%nsppol
235 
236  ddb_hdr%occopt = dtset%occopt
237  ddb_hdr%usepaw = dtset%usepaw
238 
239  ddb_hdr%intxc = dtset%intxc
240  ddb_hdr%ixc = dtset%ixc
241  ddb_hdr%iscf = dtset%iscf
242 
243  ddb_hdr%dilatmx = dtset%dilatmx
244  ddb_hdr%ecut = dtset%ecut
245  ddb_hdr%ecutsm = dtset%ecutsm
246  ddb_hdr%pawecutdg = dtset%pawecutdg
247  ddb_hdr%kptnrm = dtset%kptnrm
248  ddb_hdr%dfpt_sciss = dtset%dfpt_sciss
249  ddb_hdr%tolwfr = 1.0_dp  ! dummy
250  ddb_hdr%tphysel = dtset%tphysel
251  ddb_hdr%tsmear = dtset%tsmear
252 
253  ddb_hdr%fullinit = 1
254 
255  call ddb_hdr_malloc(ddb_hdr)
256 
257  ! Copy arrays from dtset
258  ddb_hdr%acell(:) = dtset%acell_orig(1:3,1)
259  ddb_hdr%rprim(:,:) = dtset%rprim_orig(1:3,1:3,1)
260  ddb_hdr%amu(:) = dtset%amu_orig(:,1)
261  ddb_hdr%nband(:) = dtset%nband(1:ddb_hdr%mkpt*ddb_hdr%nsppol)
262  ddb_hdr%symafm(:) = dtset%symafm(1:ddb_hdr%msym)
263  ddb_hdr%symrel(:,:,:) = dtset%symrel(1:3,1:3,1:ddb_hdr%msym)
264  ddb_hdr%typat(:) = dtset%typat(1:ddb_hdr%matom)
265  ddb_hdr%kpt(:,:) = dtset%kpt(1:3,1:ddb_hdr%mkpt)
266  ddb_hdr%wtk(:) = dtset%wtk(1:ddb_hdr%mkpt)
267  ddb_hdr%spinat(:,:) = dtset%spinat(1:3,1:ddb_hdr%matom)
268  ddb_hdr%tnons(:,:) = dtset%tnons(1:3,1:ddb_hdr%msym)
269  ddb_hdr%zion(:) = dtset%ziontypat(1:ddb_hdr%mtypat)
270  ddb_hdr%znucl(:) = dtset%znucl(1:ddb_hdr%mtypat)
271 
272  if (present(xred)) then
273    ddb_hdr%xred(:,:) = xred(1:3,1:ddb_hdr%matom)
274  else
275    ddb_hdr%xred(:,:) = dtset%xred_orig(1:3,1:ddb_hdr%matom,1)
276  end if
277  if (present(occ)) then
278    ddb_hdr%occ(:) = occ(1:ddb_hdr%mband*ddb_hdr%mkpt*ddb_hdr%nsppol)
279  else
280    ddb_hdr%occ(:) = dtset%occ_orig(1:ddb_hdr%mband*ddb_hdr%mkpt*ddb_hdr%nsppol,1)
281  end if
282 
283 
284  ! GA: I had way too much problems implementing pawtab_copy.
285  !     The script check-libpaw would report all sorts of errors.
286  !     Therefore, I do a cheap copy here, copying only the relevant info.
287  !call pawtab_copy(pawtab, ddb_hdr%pawtab)
288  nn=size(pawtab)
289  if (nn.gt.0) then
290    do ii=1,nn
291     ddb_hdr%pawtab(ii)%basis_size = pawtab(ii)%basis_size
292     ddb_hdr%pawtab(ii)%lmn_size = pawtab(ii)%lmn_size
293     ddb_hdr%pawtab(ii)%lmn2_size = pawtab(ii)%lmn2_size
294     ddb_hdr%pawtab(ii)%rpaw = pawtab(ii)%rpaw
295     ddb_hdr%pawtab(ii)%rshp = pawtab(ii)%rshp
296     ddb_hdr%pawtab(ii)%shape_type = pawtab(ii)%shape_type
297     if (allocated(pawtab(ii)%dij0)) then
298       call alloc_copy(pawtab(ii)%dij0, ddb_hdr%pawtab(ii)%dij0)
299     end if
300    end do
301  end if
302 
303 end subroutine ddb_hdr_init

m_ddb_hdr/ddb_hdr_malloc [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_hdr_malloc

FUNCTION

  Allocate dynamic memory.

INPUTS

OUTPUT

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

326 subroutine ddb_hdr_malloc(ddb_hdr)
327 
328 
329 !This section has been created automatically by the script Abilint (TD).
330 !Do not modify the following lines by hand.
331 #undef ABI_FUNC
332 #define ABI_FUNC 'ddb_hdr_malloc'
333 !End of the abilint section
334 
335  implicit none
336 
337 !Arguments ------------------------------------
338  type(ddb_hdr_type),intent(inout) :: ddb_hdr
339 
340 ! ************************************************************************
341 
342  ! integer
343  ABI_MALLOC(ddb_hdr%nband,(ddb_hdr%mkpt*ddb_hdr%nsppol))
344  ABI_MALLOC(ddb_hdr%symafm,(ddb_hdr%msym))
345  ABI_MALLOC(ddb_hdr%symrel,(3,3,ddb_hdr%msym))
346  ABI_MALLOC(ddb_hdr%typat,(ddb_hdr%matom))
347 
348  ! real
349  ABI_MALLOC(ddb_hdr%amu,(ddb_hdr%mtypat))
350  ABI_MALLOC(ddb_hdr%kpt,(3, ddb_hdr%mkpt))
351  ABI_MALLOC(ddb_hdr%occ,(ddb_hdr%mband*ddb_hdr%mkpt*ddb_hdr%nsppol))
352  ABI_MALLOC(ddb_hdr%spinat,(3, ddb_hdr%matom))
353  ABI_MALLOC(ddb_hdr%tnons,(3, ddb_hdr%msym))
354  ABI_MALLOC(ddb_hdr%wtk,(ddb_hdr%mkpt))
355  ABI_MALLOC(ddb_hdr%xred,(3, ddb_hdr%matom))
356  ABI_MALLOC(ddb_hdr%zion,(ddb_hdr%mtypat))
357  ABI_MALLOC(ddb_hdr%znucl,(ddb_hdr%mtypat))
358 
359  ! types
360  ABI_DATATYPE_ALLOCATE(ddb_hdr%pawtab,(ddb_hdr%psps%ntypat*ddb_hdr%psps%usepaw))
361  call pawtab_nullify(ddb_hdr%pawtab)
362 
363 end subroutine ddb_hdr_malloc

m_ddb_hdr/ddb_hdr_open_read [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_hdr_open_read

FUNCTION

  Open the DDB file and read the header.

INPUTS

OUTPUT

PARENTS

      anaddb,m_ddb,m_effective_potential_file,m_gruneisen,mblktyp1,mblktyp5
      mrgddb,thmeig

CHILDREN

SOURCE

544 subroutine ddb_hdr_open_read(ddb_hdr, filnam, unddb, ddb_version, comm, &
545 &        matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw,dimonly)
546 
547 
548 !This section has been created automatically by the script Abilint (TD).
549 !Do not modify the following lines by hand.
550 #undef ABI_FUNC
551 #define ABI_FUNC 'ddb_hdr_open_read'
552 !End of the abilint section
553 
554  implicit none
555 
556 !Arguments ------------------------------------
557  type(ddb_hdr_type),intent(inout) :: ddb_hdr
558  character(len=*),intent(in) :: filnam
559  integer,intent(in) :: unddb
560  integer,intent(in) :: ddb_version
561  integer,intent(in),optional :: comm
562  integer,intent(in),optional :: matom,mtypat,mband,mkpt,msym
563  integer,intent(in),optional :: dimekb,lmnmax,usepaw
564  integer,intent(in),optional :: dimonly
565 
566 !Local variables -------------------------
567  integer :: choice
568  integer :: mblktyp,nblok
569  integer :: matom_l,mtypat_l,mband_l,mkpt_l,msym_l
570  integer :: dimekb_l,lmnmax_l,usepaw_l
571  integer :: comm_l
572 
573 ! ************************************************************************
574 
575  if (present(comm)) then
576    comm_l = comm
577  else
578    comm_l = xmpi_comm_self
579  end if
580 
581  ! Read the dimensions from the DDB
582  call ddb_getdims(dimekb_l,filnam,lmnmax_l,mband_l,mblktyp, &
583 &       msym_l,matom_l,nblok,mkpt_l,mtypat_l,unddb,usepaw_l, &
584 &       ddb_version,comm_l)
585 
586  close(unddb)
587 
588  if (present(matom)) matom_l = matom
589  if (present(mtypat)) mtypat_l = mtypat
590  if (present(mband)) mband_l = mband
591  if (present(mkpt)) mkpt_l = mkpt
592  if (present(msym)) msym_l = msym
593  if (present(dimekb)) dimekb_l = dimekb
594  if (present(lmnmax)) lmnmax_l = lmnmax
595  if (present(usepaw)) usepaw_l = usepaw
596 
597  ddb_hdr%ddb_version = ddb_version
598  ddb_hdr%usepaw = usepaw_l
599  ddb_hdr%mband = mband_l
600  ddb_hdr%matom = matom_l
601  ddb_hdr%msym = msym_l
602  ddb_hdr%mtypat = mtypat_l
603  ddb_hdr%mkpt = mkpt_l
604  ddb_hdr%nsppol = 1     ! GA: Is nsppol not read?? Have to fix this...
605 
606  ddb_hdr%nblok = nblok
607  ddb_hdr%mblktyp = mblktyp
608 
609  ddb_hdr%psps%dimekb = dimekb_l
610  ddb_hdr%psps%ntypat = mtypat_l
611  ddb_hdr%psps%lmnmax = lmnmax_l
612  ddb_hdr%psps%usepaw = usepaw_l
613  ddb_hdr%psps%useylm = usepaw_l  ! yep...
614 
615  ddb_hdr%natom = ddb_hdr%matom
616  ddb_hdr%nkpt = ddb_hdr%mkpt
617  ddb_hdr%nsym = ddb_hdr%msym
618  ddb_hdr%ntypat = ddb_hdr%mtypat
619 
620  if (present(dimonly)) then
621    if (dimonly==1) return
622  end if
623 
624  ! Allocate the memory
625  call ddb_hdr_malloc(ddb_hdr)
626 
627  ABI_ALLOCATE(ddb_hdr%psps%indlmn,(6,ddb_hdr%psps%lmnmax,ddb_hdr%mtypat))
628  ABI_ALLOCATE(ddb_hdr%psps%pspso,(ddb_hdr%mtypat))
629  ABI_ALLOCATE(ddb_hdr%psps%ekb,(ddb_hdr%psps%dimekb,ddb_hdr%mtypat))
630 
631  ! This is needed to read the DDBs in the old format
632  ! GA : Not sure why this is required
633  ddb_hdr%symafm(:)=1
634  if(ddb_hdr%mtypat>=1)then
635    ddb_hdr%psps%pspso(:)=0
636    ddb_hdr%znucl(:)=zero
637    ddb_hdr%psps%ekb(:,:)=zero
638  end if
639  if(ddb_hdr%matom>=1)then
640    ddb_hdr%spinat(:,:)=zero
641  end if
642 
643 
644  ! Note: the maximum parameters (matom, mkpt, etc.) are inputs to ioddb8_in
645  !       wile the actual parameters (natom, nkpt, etc.) are outputs
646  call ioddb8_in(filnam,ddb_hdr%matom,ddb_hdr%mband,&
647 &       ddb_hdr%mkpt,ddb_hdr%msym,ddb_hdr%mtypat,unddb,ddb_version,&
648 &       ddb_hdr%acell,ddb_hdr%amu,ddb_hdr%dilatmx,ddb_hdr%ecut,ddb_hdr%ecutsm,&
649 &       ddb_hdr%intxc,ddb_hdr%iscf,ddb_hdr%ixc,ddb_hdr%kpt,ddb_hdr%kptnrm,&
650 &       ddb_hdr%natom,ddb_hdr%nband,ddb_hdr%ngfft,ddb_hdr%nkpt,ddb_hdr%nspden,&
651 &       ddb_hdr%nspinor,ddb_hdr%nsppol,ddb_hdr%nsym,ddb_hdr%ntypat,&
652 &       ddb_hdr%occ,ddb_hdr%occopt,ddb_hdr%pawecutdg,ddb_hdr%rprim,&
653 &       ddb_hdr%dfpt_sciss,ddb_hdr%spinat,ddb_hdr%symafm,ddb_hdr%symrel,&
654 &       ddb_hdr%tnons,ddb_hdr%tolwfr,ddb_hdr%tphysel,ddb_hdr%tsmear,&
655 &       ddb_hdr%typat,ddb_hdr%usepaw,ddb_hdr%wtk,ddb_hdr%xred,ddb_hdr%zion,&
656 &       ddb_hdr%znucl)
657 
658 
659 !  Read the psp information of the input DDB
660    choice=1  ! Read
661    call psddb8(choice,ddb_hdr%psps%dimekb,ddb_hdr%psps%ekb,ddb_hdr%fullinit,&
662 &       ddb_hdr%psps%indlmn,ddb_hdr%psps%lmnmax,&
663 &       ddb_hdr%nblok,ddb_hdr%ntypat,unddb,ddb_hdr%pawtab,ddb_hdr%psps%pspso,&
664 &       ddb_hdr%psps%usepaw,ddb_hdr%psps%useylm,ddb_version)
665 
666 end subroutine ddb_hdr_open_read

m_ddb_hdr/ddb_hdr_open_write [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_hdr_open_write

FUNCTION

  Open the DDB file and write the header.

INPUTS

OUTPUT

PARENTS

      ddb_interpolate,dfpt_looppert,eig2tot,gstate,mblktyp1,mblktyp5
      nonlinear,respfn

CHILDREN

SOURCE

477 subroutine ddb_hdr_open_write(ddb_hdr, filnam, unddb, fullinit)
478 
479 
480 !This section has been created automatically by the script Abilint (TD).
481 !Do not modify the following lines by hand.
482 #undef ABI_FUNC
483 #define ABI_FUNC 'ddb_hdr_open_write'
484 !End of the abilint section
485 
486  implicit none
487 
488 !Arguments ------------------------------------
489  type(ddb_hdr_type),intent(inout) :: ddb_hdr
490  character(len=*),intent(in) :: filnam
491  integer,intent(in) :: unddb
492  integer,intent(in),optional :: fullinit
493 
494 !Local variables -------------------------
495  integer :: choice
496 
497 ! ************************************************************************
498 
499  if (present(fullinit)) ddb_hdr%fullinit = fullinit
500 
501  call ddb_io_out(ddb_hdr%dscrpt,filnam,ddb_hdr%matom,ddb_hdr%mband,&
502 &  ddb_hdr%mkpt,ddb_hdr%msym,ddb_hdr%mtypat,unddb,ddb_hdr%ddb_version,&
503 &  ddb_hdr%acell,ddb_hdr%amu,ddb_hdr%dilatmx,ddb_hdr%ecut,ddb_hdr%ecutsm,&
504 &  ddb_hdr%intxc,ddb_hdr%iscf,ddb_hdr%ixc,ddb_hdr%kpt,ddb_hdr%kptnrm,&
505 &  ddb_hdr%natom,ddb_hdr%nband,ddb_hdr%ngfft,ddb_hdr%nkpt,ddb_hdr%nspden,&
506 &  ddb_hdr%nspinor,ddb_hdr%nsppol,ddb_hdr%nsym,ddb_hdr%ntypat,ddb_hdr%occ,&
507 &  ddb_hdr%occopt,ddb_hdr%pawecutdg,ddb_hdr%rprim,ddb_hdr%dfpt_sciss,&
508 &  ddb_hdr%spinat,ddb_hdr%symafm,ddb_hdr%symrel,ddb_hdr%tnons,ddb_hdr%tolwfr,&
509 &  ddb_hdr%tphysel,ddb_hdr%tsmear,ddb_hdr%typat,ddb_hdr%usepaw,ddb_hdr%wtk,&
510 &  ddb_hdr%xred,ddb_hdr%zion,ddb_hdr%znucl)
511 
512 
513  choice=2  ! Write
514  call psddb8(choice,ddb_hdr%psps%dimekb,ddb_hdr%psps%ekb,ddb_hdr%fullinit,&
515 &  ddb_hdr%psps%indlmn,ddb_hdr%psps%lmnmax,ddb_hdr%nblok,ddb_hdr%ntypat,unddb,&
516 &  ddb_hdr%pawtab,ddb_hdr%psps%pspso,ddb_hdr%psps%usepaw,ddb_hdr%psps%useylm,&
517 &  ddb_hdr%ddb_version)
518 
519 
520 end subroutine ddb_hdr_open_write

m_ddb_hdr/ddb_io_out [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ddb_io_out

FUNCTION

 Open Derivative DataBase, then
 write Derivative DataBase preliminary information.
 Note: only one processor writes the DDB.

INPUTS

 acell(3)=length scales of primitive translations (bohr)
 amu(mtypat)=mass of the atoms (atomic mass unit)
 dilatmx=the maximal dilatation factor
 character(len=fnlen) dscrpt:string that describe the output database
 ecut=kinetic energy planewave cutoff (hartree)
 ecutsm=smearing energy for plane wave kinetic energy (Ha)
 character(len=fnlen) filnam: name of output file
 intxc=control xc quadrature
 iscf=parameter controlling scf or non-scf choice
 ixc=exchange-correlation choice parameter
 kpt(3,mkpt)=k point set (reduced coordinates)
 kptnrm=normalisation of k points
 matom=maximum number of atoms
 mband=maximum number of bands
 mkpt=maximum number of special points
 msym=maximum number of symetries
 mtypat=maximum number of atom types
 natom=number of atoms in the unit cell
 nband(mkpt)=number of bands at each k point, for each polarization
 ngfft(18)=contain all needed information about 3D FFT,
        see ~abinit/doc/variables/vargs.htm#ngfft
 nkpt=number of k points
 nspden=number of spin-density components
 nspinor=number of spinorial components of the wavefunctions
 nsppol=1 for unpolarized, 2 for spin-polarized
 nsym=number of symmetry elements in space group
 ntypat=number of atom types
 occ(mband*mkpt)=occupation number for each band and k
 occopt=option for occupancies
 pawecutdg=cut-off for fine "double grid" used in PAW calculations (unused for NCPP)
 rprim(3,3)=dimensionless primitive translations in real space
 dfpt_sciss=scissor shift (Ha)
 spinat(3,matom)=initial spin of each atom, in unit of hbar/2
 symafm(msym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,msym)=symmetry operations in real space
 tnons(3,msym)=nonsymmorphic translations for symmetry operations
 tolwfr=tolerance on largest wf residual
 tphysel="physical" electronic temperature with FD occupations
 tsmear=smearing width (or temperature) in Hartree
 typat(matom)=type of each atom
 unddb=unit number for output
 usepaw=flag for PAW
 vrsddb=6 digit integer giving date, in form yymmdd for month=mm(1-12),
  day=dd(1-31), and year=yy(90-99 for 1990 to 1999,00-89 for 2000 to 2089),
  of current DDB version.
 wtk(mkpt)=weight assigned to each k point
 xred(3,matom)=reduced atomic coordinates
 zion(mtypat)=valence charge of each type of atom
 znucl(mtypat)=atomic number of atom type

OUTPUT

  Only writing

PARENTS

      m_ddb_hdr

CHILDREN

      wrtout

SOURCE

3092 subroutine ddb_io_out (dscrpt,filnam,matom,mband,&
3093 &  mkpt,msym,mtypat,unddb,vrsddb,&
3094 &  acell,amu,dilatmx,ecut,ecutsm,intxc,iscf,ixc,kpt,kptnrm,&
3095 &  natom,nband,ngfft,nkpt,nspden,nspinor,nsppol,nsym,ntypat,occ,occopt,&
3096 &  pawecutdg,rprim,dfpt_sciss,spinat,symafm,symrel,tnons,tolwfr,tphysel,tsmear,&
3097 &  typat,usepaw,wtk,xred,zion,znucl)
3098 
3099 
3100 !This section has been created automatically by the script Abilint (TD).
3101 !Do not modify the following lines by hand.
3102 #undef ABI_FUNC
3103 #define ABI_FUNC 'ddb_io_out'
3104 !End of the abilint section
3105 
3106  implicit none
3107 
3108 !Arguments -------------------------------
3109 !scalars
3110  integer,intent(in) :: matom,mband,mkpt,msym,mtypat,unddb,vrsddb
3111  integer,intent(in) :: intxc,iscf,ixc,natom,nkpt,nspden,nspinor,nsppol,nsym
3112  integer,intent(in) :: ntypat,occopt,usepaw
3113  real(dp),intent(in) :: dilatmx,ecut,ecutsm,kptnrm,pawecutdg,dfpt_sciss,tolwfr,tphysel
3114  real(dp),intent(in) :: tsmear
3115  character(len=fnlen),intent(in) :: dscrpt,filnam
3116 !arrays
3117  integer,intent(in) :: nband(mkpt*nsppol),ngfft(18),symafm(msym),symrel(3,3,msym)
3118  integer,intent(in) :: typat(matom)
3119  real(dp),intent(in) :: acell(3),amu(mtypat),kpt(3,mkpt),occ(mband*mkpt*nsppol)
3120  real(dp),intent(in) :: rprim(3,3),spinat(3,matom),tnons(3,msym),wtk(mkpt)
3121  real(dp),intent(in) :: xred(3,matom),zion(mtypat),znucl(mtypat)
3122 
3123 !Local variables -------------------------
3124 !Set routine version number here:
3125 !scalars
3126  integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527
3127  integer :: bantot,ii,ij,ikpt,iline,im,ierr
3128  character(len=500) :: message
3129 !arrays
3130  character(len=9) :: name(9)
3131 
3132 ! *********************************************************************
3133 
3134  DBG_ENTER("COLL")
3135 
3136 
3137 !Check ioddb8 version number (vrsio8) against mkddb version number
3138 !(vrsddb)
3139  if (vrsio8/=vrsddb) then
3140    write(message, '(a,a,a,i10,a,a,i10,a)' )&
3141 &   ' ddb_io_out: WARNING -',ch10,&
3142 &   '  The input/output DDB version number=',vrsio8,ch10,&
3143 &   '  is not equal to the DDB version number=',vrsddb,'.'
3144    call wrtout(std_out,message,'COLL')
3145  end if
3146 
3147 !Open the output derivative database.
3148 !(version 2.1. : changed because of a bug in a Perl script
3149 !should set up a name checking procedure, with change of name
3150 !like for the output file)
3151  ierr = open_file(filnam,message,unit=unddb,status='unknown',form='formatted')
3152  if (ierr /= 0) then
3153    MSG_ERROR(message)
3154  end if
3155 
3156 !Write the heading
3157  write(unddb, '(/,a,/,a,i10,/,/,a,a,/)' ) &
3158 & ' **** DERIVATIVE DATABASE ****    ',&
3159 & '+DDB, Version number',vrsddb,' ',trim(dscrpt)
3160 
3161 !Write the descriptive data
3162 !1. usepaw
3163  write(unddb, '(1x,a9,i10)' )'   usepaw',usepaw
3164 !2. natom
3165  write(unddb, '(1x,a9,i10)' )'    natom',natom
3166 !3. nkpt
3167  write(unddb, '(1x,a9,i10)' )'     nkpt',nkpt
3168 !4. nsppol
3169  write(unddb, '(1x,a9,i10)' )'   nsppol',nsppol
3170 !5. nsym
3171  write(unddb, '(1x,a9,i10)' )'     nsym',nsym
3172 !6. ntypat
3173  write(unddb, '(1x,a9,i10)' )'   ntypat',ntypat
3174 !7. occopt
3175  write(unddb, '(1x,a9,i10)' )'   occopt',occopt
3176 !8. nband
3177  if(occopt==2)then
3178    im=12
3179    name(1)='    nband'
3180    do iline=1,(nkpt+11)/12
3181      if(iline==(nkpt+11)/12)im=nkpt-12*(iline-1)
3182      write(unddb, '(1x,a9,5x,12i5)' )name(1),&
3183 &     (nband((iline-1)*12+ii),ii=1,im)
3184      name(1)='         '
3185    end do
3186    bantot=0
3187    do ikpt=1,nkpt
3188      bantot=bantot+nband(ikpt)
3189    end do
3190  else
3191    write(unddb, '(1x,a9,i10)' )'    nband',nband(1)
3192    bantot=nkpt*nband(1)
3193  end if
3194 
3195 !9. acell
3196  write(unddb, '(1x,a9,3d22.14)' )'    acell',acell
3197 !10. amu
3198  im=3
3199  name(1)='      amu'
3200  do iline=1,(ntypat+2)/3
3201    if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1)
3202    write (unddb, '(1x,a9,3d22.14)' )name(1),&
3203 &   (amu((iline-1)*3+ii),ii=1,im)
3204    name(1)='         '
3205  end do
3206 !11. dilatmx
3207  write(unddb, '(1x,a9,d22.14)' )'  dilatmx',dilatmx
3208 !12. ecut
3209  write(unddb, '(1x,a9,d22.14)' )'     ecut',ecut
3210 !12b. pawecutdg (PAW)
3211  if (usepaw==1) then
3212    write(unddb, '(1x,a9,d22.14)' )'pawecutdg',pawecutdg
3213  end if
3214 !13. ecutsm
3215  write(unddb, '(1x,a9,d22.14)' )'   ecutsm',ecutsm
3216 !14. intxc
3217  write(unddb, '(1x,a9,i10)' )'    intxc',intxc
3218 !15. iscf
3219  write(unddb, '(1x,a9,i10)' )'     iscf',iscf
3220 !16. ixc
3221  write(unddb, '(1x,a9,i10)' )'      ixc',ixc
3222 !17. kpt
3223  name(1)='      kpt'
3224  do iline=1,nkpt
3225    write (unddb, '(1x,a9,3d22.14)' )name(1),&
3226 &   (kpt(ii,iline),ii=1,3)
3227    name(1)='      '
3228  end do
3229 !18. kptnrm
3230  write(unddb, '(1x,a9,d22.14)' )'   kptnrm',kptnrm
3231 !19. ngfft
3232  write(unddb, '(1x,a9,5x,3i5)' )'    ngfft',ngfft(1:3)
3233 !20. nspden
3234  write(unddb, '(1x,a9,i10)' )'   nspden',nspden
3235 !21. nspinor
3236  write(unddb, '(1x,a9,i10)' )'  nspinor',nspinor
3237 !22. occ
3238  if(occopt==2)then
3239    im=3
3240    name(1)='      occ'
3241    do iline=1,(bantot+2)/3
3242      if(iline==(bantot+2)/3)im=bantot-3*(iline-1)
3243      write(unddb, '(1x,a9,3d22.14)' )name(1),&
3244 &     (occ((iline-1)*3+ii),ii=1,im)
3245      name(1)='         '
3246    end do
3247  else
3248    im=3
3249    name(1)='      occ'
3250    do iline=1,(nband(1)+2)/3
3251      if(iline==(nband(1)+2)/3)im=nband(1)-3*(iline-1)
3252      write(unddb, '(1x,a9,3d22.14)' )name(1),&
3253 &     (occ((iline-1)*3+ii),ii=1,im)
3254      name(1)='         '
3255    end do
3256  end if
3257 !23. rprim
3258  name(1)='    rprim'
3259  do iline=1,3
3260    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3261 &   (rprim(ii,iline),ii=1,3)
3262    name(1)='      '
3263  end do
3264 !24. dfpt_sciss
3265  write(unddb, '(1x,a11,d22.14)' )' dfpt_sciss',dfpt_sciss
3266 !25. spinat
3267  name(1)='   spinat'
3268  do iline=1,natom
3269    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3270 &   (spinat(ii,iline),ii=1,3)
3271    name(1)='         '
3272  end do
3273 !26. symafm
3274  im=12
3275  name(1)='   symafm'
3276  do iline=1,(nsym+11)/12
3277    if(iline==(nsym+11)/12)im=nsym-12*(iline-1)
3278    write(unddb, '(1x,a9,5x,12i5)' )name(1),&
3279 &   (symafm((iline-1)*12+ii),ii=1,im)
3280    name(1)='         '
3281  end do
3282 !27. symrel
3283  name(1)='   symrel'
3284  do iline=1,nsym
3285    write(unddb, '(1x,a9,5x,9i5)' )name(1),&
3286 &   ((symrel(ii,ij,iline),ii=1,3),ij=1,3)
3287    name(1)='         '
3288  end do
3289 !28. tnons
3290  name(1)='    tnons'
3291  do iline=1,nsym
3292    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3293 &   (tnons(ii,iline),ii=1,3)
3294    name(1)='         '
3295  end do
3296 !29. tolwfr
3297  write(unddb, '(1x,a9,d22.14)' )'   tolwfr',tolwfr
3298 !30. tphysel
3299  write(unddb, '(1x,a9,d22.14)' )'  tphysel',tphysel
3300 !31. tsmear
3301  write(unddb, '(1x,a9,d22.14)' )'   tsmear',tsmear
3302 !32. typat
3303  im=12
3304  name(1)='    typat'
3305  do iline=1,(natom+11)/12
3306    if(iline==(natom+11)/12)im=natom-12*(iline-1)
3307    write(unddb, '(1x,a9,5x,12i5)' )name(1),&
3308 &   (typat((iline-1)*12+ii),ii=1,im)
3309    name(1)='         '
3310  end do
3311 !33. wtk
3312  name(1)='      wtk'
3313  im=3
3314  do iline=1,(nkpt+2)/3
3315    if(iline==(nkpt+2)/3)im=nkpt-3*(iline-1)
3316    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3317 &   (wtk((iline-1)*3+ii),ii=1,im)
3318    name(1)='         '
3319  end do
3320 !34. xred
3321  name(1)='     xred'
3322  do iline=1,natom
3323    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3324 &   (xred(ii,iline),ii=1,3)
3325    name(1)='         '
3326  end do
3327 !35. znucl
3328  name(1)='    znucl'
3329  im=3
3330  do iline=1,(ntypat+2)/3
3331    if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1)
3332    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3333 &   (znucl((iline-1)*3+ii),ii=1,im)
3334    name(1)='         '
3335  end do
3336 !36. zion
3337  name(1)='     zion'
3338  im=3
3339  do iline=1,(ntypat+2)/3
3340    if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1)
3341    write(unddb, '(1x,a9,3d22.14)' )name(1),&
3342 &   (zion((iline-1)*3+ii),ii=1,im)
3343    name(1)='         '
3344  end do
3345 
3346  DBG_EXIT("COLL")
3347 
3348 end subroutine ddb_io_out

m_ddb_hdr/inprep8 [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 inprep8

FUNCTION

 Open Derivative DataBase, then reads the variables that
 must be known in order to dimension the arrays before complete reading
 Note: only one processor read or write the DDB.

INPUTS

 character(len=*) filnam: name of input or output file
 unddb=unit number for input or output
 vrsddb=6 digit integer giving date, in form yymmdd for month=mm(1-12),
  day=dd(1-31), and year=yy(90-99 for 1990 to 1999,00-89 for 2000 to 2089), of current DDB version.

OUTPUT

 dimekb=dimension of ekb (only used for norm-conserving psps)
 lmnmax=if useylm=1, max number of (l,m,n) comp. over all type of psps
       =if useylm=0, max number of (l,n)   comp. over all type of psps
 mband=maximum number of bands
 mblktyp=largest block type
 msym=maximum number of symmetries
 natom=number of atoms
 nblok=number of bloks in the DDB
 nkpt=number of k points
 ntypat=number of atom types
 usepaw= 0 for non paw calculation; =1 for paw calculation

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

1955 subroutine inprep8 (dimekb,filnam,lmnmax,mband,mblktyp,msym,natom,nblok,nkpt,&
1956 & ntypat,unddb,usepaw,vrsddb)
1957 
1958 
1959 !This section has been created automatically by the script Abilint (TD).
1960 !Do not modify the following lines by hand.
1961 #undef ABI_FUNC
1962 #define ABI_FUNC 'inprep8'
1963 !End of the abilint section
1964 
1965  implicit none
1966 
1967 !Arguments -------------------------------
1968 !scalars
1969  integer,intent(in) :: unddb,vrsddb
1970  integer, intent(out) :: msym
1971  integer,intent(out) :: dimekb,lmnmax,mband,mblktyp,natom,nblok,nkpt,ntypat,usepaw
1972  character(len=*),intent(in) :: filnam
1973 
1974 !Local variables -------------------------
1975 !scalars
1976 !Set routine version number here:
1977  integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527
1978  integer :: bantot,basis_size0,blktyp,ddbvrs,iband,iblok,iekb,ii,ikpt,iline,im,ios,iproj
1979  integer :: itypat,itypat0,jekb,lmn_size0,mproj,mpsang,nekb,nelmts,nsppol
1980  integer :: occopt,pspso0,nsym
1981  logical :: ddbvrs_is_current_or_old,testn,testv
1982  character(len=12) :: string
1983  character(len=32) :: blkname
1984  character(len=500) :: message
1985  character(len=6) :: name_old
1986  character(len=80) :: rdstring
1987 !arrays
1988  integer,allocatable :: nband(:)
1989  character(len=12) :: name(9)
1990 
1991 ! *********************************************************************
1992 
1993 !Check inprep8 version number (vrsio8) against mkddb version number (vrsddb)
1994  if (vrsio8/=vrsddb) then
1995    write(message, '(a,i0,2a,i0)' )&
1996 &   'The input/output DDB version number= ',vrsio8,ch10,&
1997 &   'is not equal to the DDB version number= ',vrsddb
1998    MSG_BUG(message)
1999  end if
2000 
2001 !Open the input derivative database.
2002  call wrtout(std_out, sjoin(" Opening DDB file:", filnam), 'COLL')
2003  if (open_file(filnam,message,unit=unddb,form="formatted",status="old",action="read") /= 0) then
2004    MSG_ERROR(message)
2005  end if
2006 
2007 !Check the compatibility of the input DDB with the DDB code
2008  read (unddb,*)
2009  read (unddb,*)
2010  read (unddb, '(20x,i10)' )ddbvrs
2011 
2012  if (all(ddbvrs/= [vrsio8, vrsio8_old, vrsio8_old_old]) )then
2013    write(message, '(a,i10,2a,3(a,i10))' )&
2014 &   'The input DDB version number=',ddbvrs,' does not agree',ch10,&
2015 &   'with the allowed code DDB version numbers,',vrsio8,', ',vrsio8_old,' and ',vrsio8_old_old
2016    MSG_ERROR(message)
2017  end if
2018 
2019 !Read the 4 n-integers, also testing the names of data,
2020 !and checking that their value is acceptable.
2021 !This is important to insure that any array has a sufficient dimension.
2022  read (unddb,*)
2023  read (unddb,*)
2024  read (unddb,*)
2025  testn=.true.
2026  testv=.true.
2027  ddbvrs_is_current_or_old=(ddbvrs==vrsio8.or.ddbvrs==vrsio8_old)
2028 
2029 !1. usepaw
2030  if(ddbvrs==vrsio8)then
2031    read (unddb, '(1x,a9,i10)' )name(1),usepaw
2032  else
2033    usepaw=0;name(1)='   usepaw'
2034  end if
2035  if(name(1)/='   usepaw')testn=.false.
2036 !2. natom
2037  if(ddbvrs_is_current_or_old)then
2038    read (unddb, '(1x,a9,i10)' )name(2),natom
2039  else
2040    read (unddb, '(1x,a6,i10)' )name_old,natom ; name(2)='   '//name_old
2041  end if
2042  if(name(2)/='    natom')testn=.false.
2043  if(natom<=0)testv=.false.
2044 !3. nkpt
2045  if(ddbvrs_is_current_or_old)then
2046    read (unddb, '(1x,a9,i10)' )name(3),nkpt
2047  else
2048    read (unddb, '(1x,a6,i10)' )name_old,nkpt ; name(3)='   '//name_old
2049  end if
2050  if(name(3)/='     nkpt')testn=.false.
2051  if(nkpt <=0)testv=.false.
2052 !4. nsppol
2053  if(ddbvrs_is_current_or_old)then
2054    read (unddb, '(1x,a9,i10)' )name(4),nsppol
2055  else
2056    read (unddb, '(1x,a6,i10)' )name_old,nsppol ; name(4)='   '//name_old
2057  end if
2058  if(name(4)/='   nsppol')testn=.false.
2059  if(nsppol<=0.or.nsppol>2)testv=.false.
2060 !5. nsym
2061  if(ddbvrs_is_current_or_old)then
2062    read (unddb, '(1x,a9,i10)' )name(5),nsym
2063  else
2064    read (unddb, '(1x,a6,i10)' )name_old,nsym ; name(5)='   '//name_old
2065  end if
2066  if(name(5)/='     nsym')testn=.false.
2067 !MG FIXME Why this and why do we need msym?
2068  msym = 192
2069  if (nsym > msym) msym=nsym
2070 !if(nsym <=0.or.nsym >msym )testv=.false.
2071 !6. ntypat
2072  if(ddbvrs_is_current_or_old)then
2073    read (unddb, '(1x,a9,i10)' )name(6),ntypat
2074  else
2075    read (unddb, '(1x,a6,i10)' )name_old,ntypat ; name(6)='   '//name_old
2076  end if
2077  if(name(6)/='   ntypat' .and. name(6)/='    ntype')testn=.false.
2078  if(ntypat<=0)testv=.false.
2079 !7. occopt
2080 !Before reading nband, the last parameters that define
2081 !the dimension of some array, need to know what is their
2082 !representation, given by occopt
2083  if(ddbvrs_is_current_or_old)then
2084    read (unddb, '(1x,a9,i10)' )name(7),occopt
2085  else
2086    read (unddb, '(1x,a6,i10)' )name_old,occopt ; name(7)='   '//name_old
2087  end if
2088  if(name(7)/='   occopt')testn=.false.
2089  if(occopt<0.or.occopt>8)testv=.false.
2090 
2091 !Message if the names or values are not right
2092  if (.not.testn.or..not.testv) then
2093    write(message, '(a,a)' )' inprep8 : An error has been found in the',' positive n-integers contained in the DDB : '
2094    call wrtout(std_out,message,'COLL')
2095    write(message, '(a)' )   '     Expected                      Found     '
2096    call wrtout(std_out,message,'COLL')
2097    write(message, '(a,i9,a,a,a,i10)' )'    natom , larger than',0,'    ',trim(name(2)),' =',natom
2098    call wrtout(std_out,message,'COLL')
2099    write(message, '(a,i9,a,a,a,i10)' )'    nkpt  , larger than',0,'    ',trim(name(3)),' =',nkpt
2100    call wrtout(std_out,message,'COLL')
2101    write(message, '(a,i1,a,a,a,i10)' )'    nsppol, either    1 or     ',2,'    ',trim(name(4)),' =',nsppol
2102    call wrtout(std_out,message,'COLL')
2103 !  write(message, '(a,i10,a,a,a,i10)' )&   '    nsym  , lower than',msym,'    ',trim(name(5)),' =',nsym
2104    call wrtout(std_out,message,'COLL')
2105    write(message, '(a,i9,a,a,a,i10)' )'    ntypat , larger than',0,'   ',trim(name(6)),' =',ntypat
2106    call wrtout(std_out,message,'COLL')
2107    write(message, '(a,a,a,i10)' )'    occopt,     equal to 0,1 or 2   ',trim(name(7)),' =',occopt
2108    call wrtout(std_out,message,'COLL')
2109 
2110    MSG_ERROR('See the error message above.')
2111  end if
2112 
2113 !One more set of parameters define the dimensions of the
2114 !array : nband. Morever, it depends on occopt and !nkpt, and has to be
2115 !tested after the test on nkpt is performed.
2116 
2117 !8. nband
2118  ABI_MALLOC(nband,(nkpt))
2119  if(occopt==2)then
2120    im=12
2121    do iline=1,(nkpt+11)/12
2122      if(iline==(nkpt+11)/12)im=nkpt-12*(iline-1)
2123      if(ddbvrs_is_current_or_old)then
2124        read (unddb, '(1x,a9,5x,12i5)' )name(1),(nband((iline-1)*12+ii),ii=1,im)
2125      else
2126        read (unddb, '(1x,a6,5x,12i5)' )name_old,&
2127 &       (nband((iline-1)*12+ii),ii=1,im) ; name(1)='   '//name_old
2128      end if
2129      if (iline==1) then
2130        call ddb_chkname(name(1),'    nband')
2131      else
2132        call ddb_chkname(name(1),'         ')
2133      end if
2134    end do
2135  else
2136    if(ddbvrs_is_current_or_old)then
2137      read (unddb, '(1x,a9,i10)' )name(1),nband(1)
2138    else
2139      read (unddb, '(1x,a6,i10)' )name_old,nband(1) ; name(1)='   '//name_old
2140    end if
2141    call ddb_chkname(name(1),'    nband')
2142    if(nkpt>1)then
2143      do ikpt=2,nkpt
2144        nband(ikpt)=nband(1)
2145      end do
2146    end if
2147  end if
2148 
2149 !Check all nband values, and sum them
2150  bantot=0
2151  do ikpt=1,nkpt
2152    if(nband(ikpt)<0)then
2153      write(message, '(a,i0,a,i0,3a)' )&
2154 &     'For ikpt = ',ikpt,'  nband = ',nband(ikpt),' is negative.',ch10,&
2155 &     'Action: correct your DDB.'
2156      MSG_ERROR(message)
2157    end if
2158    bantot=bantot+nband(ikpt)
2159  end do
2160 
2161  mband=maxval(nband(:))
2162 
2163 !Skip the rest of variables
2164 !9. acell
2165  read (unddb,*)
2166 !10. amu
2167  do iline=1,(ntypat+2)/3
2168    read (unddb,*)
2169  end do
2170 !11. dilatmx
2171  if(ddbvrs_is_current_or_old) read (unddb,*)
2172 !12. ecut
2173  read (unddb,*)
2174 !12b. pawecutdg (PAW only)
2175  if(ddbvrs==vrsio8.and.usepaw==1) read (unddb,*)
2176 !13. ecutsm
2177  if(ddbvrs_is_current_or_old) read (unddb,*)
2178 !14. intxc
2179  if(ddbvrs_is_current_or_old) read (unddb,*)
2180 !15. iscf
2181  read (unddb,*)
2182 !16. ixc
2183  read (unddb,*)
2184 !17. kpt
2185  do iline=1,nkpt
2186    read (unddb,*)
2187  end do
2188 !18. kptnrm
2189  read (unddb,*)
2190 !19. ngfft
2191  read (unddb,*)
2192 !20. nspden
2193  if(ddbvrs_is_current_or_old) read (unddb,*)
2194 !21. nspinor
2195  if(ddbvrs_is_current_or_old) read (unddb,*)
2196 !22. occ
2197  if(occopt==2)then
2198    do iline=1,(bantot+2)/3
2199      read (unddb,*)
2200    end do
2201  else
2202    write(message,*)' inprep8 : nband(1)=',nband(1)
2203    call wrtout(std_out,message,'COLL')
2204    do iline=1,(nband(1)+2)/3
2205      read (unddb,'(a80)')rdstring
2206      write(message,*)trim(rdstring)
2207      call wrtout(std_out,message,'COLL')  ! GA: why are we printing this?
2208    end do
2209  end if
2210 !23. rprim
2211  do iline=1,3
2212    read (unddb,*)
2213  end do
2214 !24. dfpt_sciss
2215  read (unddb,*)
2216 !25. spinat
2217  if(ddbvrs_is_current_or_old)then
2218    do iline=1,natom
2219      read (unddb,*)
2220    end do
2221  end if
2222 !26. symafm
2223  if(ddbvrs_is_current_or_old)then
2224    do iline=1,(nsym+11)/12
2225      read (unddb,*)
2226    end do
2227  end if
2228 !27. symrel
2229  do iline=1,nsym
2230    read (unddb,*)
2231  end do
2232 !28old. xred
2233  if(.not.ddbvrs_is_current_or_old)then
2234    do iline=1,natom
2235      read (unddb,*)
2236    end do
2237  end if
2238 !28. tnons
2239  do iline=1,nsym
2240    read (unddb,*)
2241  end do
2242 !29. tolwfr
2243  if(ddbvrs_is_current_or_old) read (unddb,*)
2244 !30. tphysel
2245  if(ddbvrs_is_current_or_old) read (unddb,*)
2246 !31. tsmear
2247  if(ddbvrs_is_current_or_old) read (unddb,*)
2248 !32. type
2249  do iline=1,(natom+11)/12
2250    read (unddb,*)
2251  end do
2252 !33old. tolwfr
2253  if(.not.ddbvrs_is_current_or_old) read (unddb,*)
2254 !33. wtk
2255  do iline=1,(nkpt+2)/3
2256    read (unddb,*)
2257  end do
2258 !34. xred
2259  if(ddbvrs_is_current_or_old)then
2260    do iline=1,natom
2261      read (unddb,*)
2262    end do
2263  end if
2264 !35. znucl
2265  if(ddbvrs_is_current_or_old)then
2266    do iline=1,(ntypat+2)/3
2267      read (unddb,*)
2268    end do
2269  end if
2270 !36. zion
2271  do iline=1,(ntypat+2)/3
2272    read (unddb,*)
2273  end do
2274 
2275  read (unddb,*)
2276 
2277 !Now, take care of the pseudopotentials
2278  read(unddb, '(a12)' )string
2279 
2280  if(string=='  Descriptio')then
2281 
2282    read (unddb,*)
2283    if (ddbvrs==vrsio8_old.or.ddbvrs==vrsio8_old_old) then
2284      read (unddb, '(10x,i3,14x,i3,11x,i3)', iostat=ios )dimekb,lmnmax,usepaw
2285      if(ios/=0)then
2286        backspace(unddb)
2287        read (unddb, '(10x,i3,14x,i3)')dimekb,lmnmax
2288        usepaw=0
2289      end if
2290    else if (ddbvrs==vrsio8) then
2291      read (unddb, '(10x,i3)') usepaw
2292      if (usepaw==0) then
2293        read (unddb, '(10x,i3,14x,i3)' ) dimekb,lmnmax
2294      else
2295        dimekb=0;lmnmax=0
2296      end if
2297    end if
2298    if (usepaw==0) then
2299      do itypat=1,ntypat
2300        read(unddb, '(13x,i4,9x,i3,8x,i4)' )itypat0,pspso0,nekb
2301        read(unddb,*)
2302        do iekb=1,nekb
2303          do jekb=1,nekb,4
2304            read(unddb,*)
2305          end do
2306        end do
2307      end do
2308    else
2309      do itypat=1,ntypat
2310        read(unddb, '(12x,i4,12x,i3,12x,i5)' )itypat0,basis_size0,lmn_size0
2311        lmnmax=max(lmnmax,lmn_size0)
2312        read(unddb,*)
2313        read(unddb,*)
2314        read(unddb,'(24x,i3)') nekb
2315        read(unddb,*)
2316        do iekb=1,nekb,4
2317          read(unddb,*)
2318        end do
2319      end do
2320    end if
2321 
2322  else if(string==' Description')then
2323    if (usepaw==1) then
2324      MSG_BUG('old DDB pspformat not compatible with PAW 1')
2325    end if
2326 
2327    read (unddb, '(10x,i3,10x,i3)' )mproj,mpsang
2328    dimekb=mproj*mpsang
2329    usepaw=0
2330    do itypat=1,ntypat
2331      read (unddb,*)
2332 !    For f-electrons, one more line has been written
2333      do iproj=1,mproj*max(1,(mpsang+2)/3)
2334        read (unddb,*)
2335      end do
2336    end do
2337 
2338  else if(string==' No informat')then
2339 
2340    dimekb=0
2341    lmnmax=0
2342    !usepaw=0   ! GA: usepaw is also declared earlier in the header
2343                !     and it is that earlier value that usepaw will
2344                !     be compared in ioddb8_in, so there is no reason
2345                !     to override the value here.
2346 
2347  else
2348    write(message, '(a,a,a,a)' )&
2349 &   'Error when reading the psp information',ch10,&
2350 &   'String=',trim(string)
2351    MSG_BUG(message)
2352  end if
2353 
2354 !Now, the number of blocks
2355  read(unddb,*)
2356  read(unddb,*)
2357  read(unddb, '(24x,i4)' )nblok
2358 
2359 !Now, the type of each blok, in turn
2360  mblktyp=1
2361  if(nblok>=1)then
2362    do iblok=1,nblok
2363 
2364      read(unddb,*)
2365      read(unddb, '(a32,12x,i8)' )blkname,nelmts
2366      if(blkname==' 2nd derivatives (non-stat.)  - ' .or.  blkname==' 2rd derivatives (non-stat.)  - ')then
2367        blktyp=1
2368      else if(blkname==' 2nd derivatives (stationary) - ' .or. blkname==' 2rd derivatives (stationary) - ')then
2369        blktyp=2
2370      else if(blkname==' 3rd derivatives              - ')then
2371        blktyp=3
2372      else if(blkname==' Total energy                 - ')then
2373        blktyp=0
2374      else if(blkname==' 1st derivatives              - ')then
2375        blktyp=4
2376      else if(blkname==' 2nd eigenvalue derivatives   - ' .or. blkname==' 2rd eigenvalue derivatives   - ')then
2377        blktyp=5
2378      else
2379        write(message, '(a,a,a,a,a,a,a,a,a)' )&
2380 &       'The following string appears in the DDB in place of',' the block type description :',ch10,blkname,ch10,&
2381 &       'Action: check your DDB.',ch10,&
2382 &       'Note: If you did use an abinit version prior to 6.12 to generate your DDB',&
2383 &       'pay attention to the change:: 2rd derivatives ==> 2nd derivatives'
2384        MSG_ERROR(message)
2385      end if
2386 
2387      if(blktyp==1.or.blktyp==2)then
2388 !      Read the phonon wavevector
2389        read(unddb,*)
2390      else if(blktyp==3)then
2391 !      Read the perturbation wavevectors
2392        read(unddb,*)
2393        read(unddb,*)
2394        read(unddb,*)
2395        mblktyp=3
2396      else if(blktyp==5)then
2397        read(unddb,*)
2398        mblktyp=5
2399      end if
2400 
2401 !    Read every element
2402      if(blktyp==5)then
2403        do ikpt=1,nkpt
2404          read(unddb,*)
2405          do iband=1,nband(ikpt)
2406            read(unddb,*)
2407            do ii=1,nelmts
2408              read(unddb,*)
2409            end do
2410          end do
2411        end do
2412      else
2413        do ii=1,nelmts
2414          read(unddb,*)
2415        end do
2416      end if
2417 
2418    end do
2419  end if
2420 
2421  ABI_FREE(nband)
2422 
2423 !Close the DDB
2424  close(unddb)
2425 
2426 end subroutine inprep8

m_ddb_hdr/ioddb8_in [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 ioddb8_in

FUNCTION

 Open Derivative DataBase, and read preliminary information.
 Note: only one processor reads the DDB.

INPUTS

 character(len=*)= name of input file
 matom=maximum number of atoms
 mband=maximum number of bands
 mkpt=maximum number of special points
 msym=maximum number of symetries
 mtypat=maximum number of atom types
 unddb=unit number for input
 vrsddb=6 digit integer giving date, in form yymmdd for month=mm(1-12),
  day=dd(1-31), and year=yy(90-99 for 1990 to 1999,00-89 for 2000 to 2089),
  of current DDB version.

OUTPUT

 acell(3)=length scales of primitive translations (bohr)
 amu(mtypat)=mass of the atoms (atomic mass unit)
 dilatmx=the maximal dilatation factor
 ecut=kinetic energy planewave cutoff (hartree)
 ecutsm=smearing energy for plane wave kinetic energy (Ha)
 intxc=control xc quadrature
 iscf=parameter controlling scf or non-scf choice
 ixc=exchange-correlation choice parameter
 kpt(3,mkpt)=k point set (reduced coordinates)
 kptnrm=normalisation of k points
 natom=number of atoms in the unit cell
 nband(mkpt)=number of bands at each k point, for each polarization
 ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
 nkpt=number of k points
 nspden=number of spin-density components
 nspinor=number of spinorial components of the wavefunctions
 nsppol=1 for unpolarized, 2 for spin-polarized
 nsym=number of symmetry elements in space group
 ntypat=number of atom types
 occ(mband*mkpt)=occupation number for each band and k
 occopt=option for occupancies
 pawecutdg=cut-off for fine "double grid" used in PAW calculations (unused for NCPP)
 rprim(3,3)=dimensionless primitive translations in real space
 dfpt_sciss=scissor shift (Ha)
 spinat(3,matom)=initial spin of each atom, in unit of hbar/2
 symafm(msym)=(anti)ferromagnetic part of symmetry operations
 symrel(3,3,msym)=symmetry operations in real space
 tnons(3,msym)=nonsymmorphic translations for symmetry operations
 tolwfr=tolerance on largest wf residual
 tphysel="physical" electronic temperature with FD occupations
 tsmear=smearing width (or temperature) in Hartree
 typat(matom)=type of each atom
 usepaw=flag for PAW
 wtk(mkpt)=weight assigned to each k point
 xred(3,matom)=reduced atomic coordinates
 zion(mtypat)=valence charge of each type of atom
 znucl(mtypat)=atomic number of atom type

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

1268 subroutine ioddb8_in(filnam,matom,mband,mkpt,msym,mtypat,unddb,vrsddb,&
1269 &  acell,amu,dilatmx,ecut,ecutsm,intxc,iscf,ixc,kpt,kptnrm,&
1270 &  natom,nband,ngfft,nkpt,nspden,nspinor,nsppol,nsym,ntypat,occ,occopt,&
1271 &  pawecutdg,rprim,dfpt_sciss,spinat,symafm,symrel,tnons,tolwfr,tphysel,tsmear,&
1272 &  typat,usepaw,wtk,xred,zion,znucl)
1273 
1274 
1275 !This section has been created automatically by the script Abilint (TD).
1276 !Do not modify the following lines by hand.
1277 #undef ABI_FUNC
1278 #define ABI_FUNC 'ioddb8_in'
1279 !End of the abilint section
1280 
1281  implicit none
1282 
1283 !Arguments -------------------------------
1284 !scalars
1285  integer,intent(in) :: matom,mband,mkpt,msym,mtypat,unddb,vrsddb,usepaw
1286  integer,intent(out) :: intxc,iscf,ixc,natom,nkpt,nspden,nspinor,nsppol,nsym,ntypat,occopt
1287  real(dp),intent(out) :: dilatmx,ecut,ecutsm,pawecutdg,kptnrm,dfpt_sciss,tolwfr,tphysel,tsmear
1288  character(len=*),intent(in) :: filnam
1289 !arrays
1290  integer,intent(out) :: nband(mkpt),ngfft(18),symafm(msym),symrel(3,3,msym),typat(matom)
1291  real(dp),intent(out) :: acell(3),amu(mtypat),kpt(3,mkpt),occ(mband*mkpt)
1292  real(dp),intent(out) :: rprim(3,3),spinat(3,matom),tnons(3,msym),wtk(mkpt)
1293  real(dp),intent(out) :: xred(3,matom),zion(mtypat),znucl(mtypat)
1294 
1295 !Local variables -------------------------
1296 !Set routine version number here:
1297 !scalars
1298  integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527
1299  integer :: bantot,ddbvrs,iband,ii,ij,ikpt,iline,im,usepaw0
1300  logical :: ddbvrs_is_current_or_old,testn,testv
1301  character(len=500) :: message
1302  character(len=6) :: name_old
1303 !arrays
1304  character(len=12) :: name(9)
1305 
1306 ! *********************************************************************
1307 
1308 !Check ioddb8 version number (vrsio8) against mkddb version number (vrsddb)
1309  if (vrsio8/=vrsddb) then
1310    write(message, '(a,i10,a,a,i10,a)' )&
1311 &   'The input/output DDB version number=',vrsio8,ch10,&
1312 &   'is not equal to the DDB version number=',vrsddb,'.'
1313    MSG_WARNING(message)
1314  end if
1315 
1316 !Open the input derivative database.
1317  write(message,'(a,a)')' About to open file ',TRIM(filnam)
1318  call wrtout(std_out,message,'COLL')
1319  if (open_file(filnam,message,unit=unddb,form="formatted",status="old",action="read") /= 0) then
1320    MSG_ERROR(message)
1321  end if
1322 
1323 !Check the compatibility of the input DDB with the DDB code
1324  read (unddb,*)
1325  read (unddb,*)
1326  read (unddb, '(20x,i10)' )ddbvrs
1327 
1328  !write(std_out,'(a,i10)')' ddbvrs=',ddbvrs
1329  if(ddbvrs/=vrsio8 .and. ddbvrs/=vrsio8_old .and. ddbvrs/=vrsio8_old_old)then
1330    write(message, '(a,i10,2a,3(a,i10),a)' )&
1331 &   'The input DDB version number=',ddbvrs,' does not agree',ch10,&
1332 &   'with the allowed code DDB version numbers,',vrsio8,', ',vrsio8_old,' and ',vrsio8_old_old,' .'
1333    MSG_BUG(message)
1334  end if
1335 
1336 !Read the 4 n-integers, also testing the names of data, and checking that their value is acceptable.
1337 !This is important to insure that any array has a sufficient dimension.
1338  read (unddb,*)
1339  read (unddb,*)
1340  read (unddb,*)
1341  testn=.true.; testv=.true.
1342  ddbvrs_is_current_or_old=(ddbvrs==vrsio8.or.ddbvrs==vrsio8_old)
1343 
1344 !1. usepaw
1345  if(ddbvrs==vrsio8)then
1346    read (unddb, '(1x,a9,i10)' )name(1),usepaw0
1347  else
1348    usepaw0=0;name(1)='   usepaw'
1349  end if
1350  if(name(1)/='   usepaw')testn=.false.
1351  if(usepaw0/=usepaw)testv=.false.
1352 !2. natom
1353  if(ddbvrs_is_current_or_old)then
1354    read (unddb, '(1x,a9,i10)' )name(2),natom
1355  else
1356    read (unddb, '(1x,a6,i10)' )name_old,natom ; name(2)='   '//name_old
1357  end if
1358  if(name(2)/='    natom')testn=.false.
1359  if(natom<=0.or.natom>matom)testv=.false.
1360 !3. nkpt
1361  if(ddbvrs_is_current_or_old)then
1362    read (unddb, '(1x,a9,i10)' )name(3),nkpt
1363  else
1364    read (unddb, '(1x,a6,i10)' )name_old,nkpt ; name(3)='   '//name_old
1365  end if
1366  if(name(3)/='     nkpt')testn=.false.
1367  if(nkpt <=0.or.nkpt >mkpt )testv=.false.
1368 !4. nsppol
1369  if(ddbvrs_is_current_or_old)then
1370    read (unddb, '(1x,a9,i10)' )name(4),nsppol
1371  else
1372    read (unddb, '(1x,a6,i10)' )name_old,nsppol ; name(4)='   '//name_old
1373  end if
1374  if(name(4)/='   nsppol')testn=.false.
1375  if(nsppol<=0.or.nsppol>2)testv=.false.
1376 !5. nsym
1377  if(ddbvrs_is_current_or_old)then
1378    read (unddb, '(1x,a9,i10)' )name(5),nsym
1379  else
1380    read (unddb, '(1x,a6,i10)' )name_old,nsym ; name(5)='   '//name_old
1381  end if
1382  if(name(5)/='     nsym')testn=.false.
1383  if(nsym <=0.or.nsym >msym )testv=.false.
1384 !6. ntypat
1385  if(ddbvrs_is_current_or_old)then
1386    read (unddb, '(1x,a9,i10)' )name(6),ntypat
1387  else
1388    read (unddb, '(1x,a6,i10)' )name_old,ntypat ; name(6)='   '//name_old
1389  end if
1390  if(name(6)/='   ntypat' .and. name(6)/='    ntype')testn=.false.
1391  if(ntypat<=0.or.ntypat>mtypat)testv=.false.
1392 !7. occopt
1393 !Before reading nband, the last parameters that define
1394 !the dimension of some array, need to know what is their
1395 !representation, given by occopt
1396  if(ddbvrs_is_current_or_old)then
1397    read (unddb, '(1x,a9,i10)' )name(7),occopt
1398  else
1399    read (unddb, '(1x,a6,i10)' )name_old,occopt ; name(7)='   '//name_old
1400  end if
1401  if(name(7)/='   occopt')testn=.false.
1402  if(occopt<0.or.occopt>8)testv=.false.
1403 !Message if the names or values are not right
1404  if (.not.testn.or..not.testv) then
1405    write(message, '(a,a,a)' )' ioddb8_in : An error has been found in one',ch10,&
1406 &   ' of the positive n-integers contained in the DDB : '
1407    call wrtout(std_out,message,'COLL')
1408    write(message, '(a)' )&
1409 &   '               Expected                      Found     '
1410    call wrtout(std_out,message,'COLL')
1411    write(message, '(a,i10,a,a,a,i10)' )&
1412 &   '    usepaw equal to   ',usepaw,'    ',trim(name(1)),' =',usepaw0
1413    call wrtout(std_out,message,'COLL')
1414    write(message, '(a,i10,a,a,a,i10)' )&
1415 &   '    natom , lower than',matom+1,'    ',trim(name(2)),' =',natom
1416    call wrtout(std_out,message,'COLL')
1417    write(message, '(a,i10,a,a,a,i10)' )&
1418 &   '    nkpt  , lower than',mkpt+1 ,'    ',trim(name(3)),' =',nkpt
1419    call wrtout(std_out,message,'COLL')
1420    write(message, '(a,i10,a,a,a,i10)' )&
1421 &   '    nsppol, lower than',3      ,'    ',trim(name(4)),' =',nsppol
1422    call wrtout(std_out,message,'COLL')
1423    write(message, '(a,i10,a,a,a,i10)' )&
1424 &   '    nsym  , lower than',msym+1 ,'    ',trim(name(5)),' =',nsym
1425    call wrtout(std_out,message,'COLL')
1426    write(message, '(a,i10,a,a,a,i10)' )&
1427 &   '    ntypat, lower than',mtypat+1,'   ',trim(name(6)),' =',ntypat
1428    call wrtout(std_out,message,'COLL')
1429    write(message, '(a,a,a,i10)' )&
1430 &   '    occopt,  between 0 and 7        ',trim(name(7)),' =',occopt
1431    call wrtout(std_out,message,'COLL')
1432 
1433    MSG_ERROR('See the error message above.')
1434  end if
1435 
1436 !One more set of parameters define the dimensions of the
1437 !array : nband. Morever, it depends on occopt and nkpt, and has to be
1438 !tested after the test on nkpt is performed.
1439 !8. nband
1440  if(occopt==2)then
1441    im=12
1442    do iline=1,(nkpt+11)/12
1443      name(:) = "" ! reset all line strings
1444      if(iline==(nkpt+11)/12)im=nkpt-12*(iline-1)
1445      if(ddbvrs_is_current_or_old)then
1446        read (unddb, '(1x,a9,5x,12i5)' )name(1),(nband((iline-1)*12+ii),ii=1,im)
1447      else
1448        read (unddb, '(1x,a6,5x,12i5)' )name_old,(nband((iline-1)*12+ii),ii=1,im) ; name(1)='   '//name_old
1449      end if
1450      if (iline==1) then
1451        call ddb_chkname(name(1),'    nband')
1452      else
1453        call ddb_chkname(name(1),'         ')
1454      end if
1455    end do
1456  else
1457    name(:) = "" ! reset all line strings
1458    if(ddbvrs_is_current_or_old)then
1459      read (unddb, '(1x,a9,i10)' )name(1),nband(1)
1460    else
1461      read (unddb, '(1x,a6,i10)' )name_old,nband(1) ; name(1)='   '//name_old
1462    end if
1463    call ddb_chkname(name(1),'    nband')
1464    if(nkpt>1)then
1465      do ikpt=2,nkpt
1466        nband(ikpt)=nband(1)
1467      end do
1468    end if
1469  end if
1470 
1471 !check all nband values, and sum them
1472  bantot=0
1473  do ikpt=1,nkpt
1474    if(nband(ikpt)<0)then
1475      write(message, '(a,i4,a,i4,3a)' )&
1476 &     'For ikpt = ',ikpt,'  nband = ',nband(ikpt),' is negative.',ch10,&
1477 &     'Action: correct your DDB.'
1478      MSG_ERROR(message)
1479    else if(nband(ikpt)>mband)then
1480      write(message, '(a,i4,a,i4,a,a,i4,3a)' )&
1481 &     'For ikpt = ',ikpt,', nband = ',nband(ikpt),ch10,&
1482 &     'is larger than mband = ',mband,'.',ch10,&
1483 &     'Action: recompile the calling code with a larger mband.'
1484      MSG_ERROR(message)
1485    end if
1486    bantot=bantot+nband(ikpt)
1487  end do
1488 
1489 !Read the rest of variables, with check of the names
1490 !9. acell
1491  if(ddbvrs_is_current_or_old)then
1492    read (unddb, '(1x,a9,3d22.14)' )name(1),acell
1493  else
1494    read (unddb, '(1x,a6,3d22.14)' )name_old,acell ; name(1)='   '//name_old
1495  end if
1496  call ddb_chkname(name(1),'    acell')
1497 !9. amu
1498  im=3
1499  do iline=1,(ntypat+2)/3
1500    if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1)
1501    if(ddbvrs_is_current_or_old)then
1502      read (unddb, '(1x,a9,3d22.14)' )name(1),(amu((iline-1)*3+ii),ii=1,im)
1503    else
1504      read (unddb, '(1x,a6,3d22.14)' )name_old,(amu((iline-1)*3+ii),ii=1,im) ; name(1)='   '//name_old
1505    end if
1506    if (iline==1) then
1507      call ddb_chkname(name(1),'      amu')
1508    else
1509      call ddb_chkname(name(1),'         ')
1510    end if
1511  end do
1512 !11. dilatmx
1513  if(ddbvrs_is_current_or_old)then
1514    read (unddb, '(1x,a9,d22.14)' )name(1),dilatmx
1515    call ddb_chkname(name(1),'  dilatmx')
1516  else
1517    dilatmx=one
1518  end if
1519 !12. ecut
1520  if(ddbvrs_is_current_or_old)then
1521    read (unddb, '(1x,a9,d22.14)' )name(1),ecut
1522  else
1523    read (unddb, '(1x,a6,d22.14)' )name_old,ecut ; name(1)='   '//name_old
1524  end if
1525  call ddb_chkname(name(1),'     ecut')
1526 !12b. pawecutdg (PAW only)
1527  if(ddbvrs==vrsio8.and.usepaw==1) then
1528    read (unddb, '(1x,a9,d22.14)' )name(1),pawecutdg
1529  else
1530    pawecutdg=ecut;name(1)='pawecutdg'
1531  end if
1532  call ddb_chkname(name(1),'pawecutdg')
1533 !13. ecutsm
1534  if(ddbvrs_is_current_or_old)then
1535    read (unddb, '(1x,a9,d22.14)' )name(1),ecutsm
1536    call ddb_chkname(name(1),'   ecutsm')
1537  else
1538    ecutsm=zero
1539  end if
1540 !14. intxc
1541  if(ddbvrs_is_current_or_old)then
1542    read (unddb, '(1x,a9,i10)' )name(1),intxc
1543    call ddb_chkname(name(1),'    intxc')
1544  else
1545    intxc=1
1546  end if
1547 !15. iscf
1548  if(ddbvrs_is_current_or_old)then
1549    read (unddb, '(1x,a9,i10)' )name(1),iscf
1550  else
1551    read (unddb, '(1x,a6,i10)' )name_old,iscf ; name(1)='   '//name_old
1552  end if
1553  call ddb_chkname(name(1),'     iscf')
1554 !16. ixc
1555  if(ddbvrs_is_current_or_old)then
1556    read (unddb, '(1x,a9,i10)' )name(1),ixc
1557  else
1558    read (unddb, '(1x,a6,i10)' )name_old,ixc ; name(1)='   '//name_old
1559  end if
1560  call ddb_chkname(name(1),'      ixc')
1561 !17. kpt
1562  do iline=1,nkpt
1563    if(ddbvrs_is_current_or_old)then
1564      read (unddb, '(1x,a9,3d22.14)' )name(1),(kpt(ii,iline),ii=1,3)
1565    else
1566      read (unddb, '(1x,a6,3d22.14)' )name_old,(kpt(ii,iline),ii=1,3) ; name(1)='   '//name_old
1567    end if
1568    if (iline==1) then
1569      call ddb_chkname(name(1),'      kpt')
1570    else
1571      call ddb_chkname(name(1),'         ')
1572    end if
1573  end do
1574 !18. kptnrm
1575  if(ddbvrs_is_current_or_old)then
1576    read (unddb, '(1x,a9,d22.14)' )name(1),kptnrm
1577  else
1578    read (unddb, '(1x,a6,d22.14)' )name_old,kptnrm ; name(1)='   '//name_old
1579  end if
1580  call ddb_chkname(name(1),'   kptnrm')
1581 !19. ngfft
1582  if(ddbvrs_is_current_or_old)then
1583    read (unddb, '(1x,a9,5x,3i5)' )name(1),ngfft(1:3)
1584  else
1585    read (unddb, '(1x,a6,5x,3i5)' )name_old,ngfft(1:3) ; name(1)='   '//name_old
1586  end if
1587 !For the time being, do not check the validity of the name,
1588 !in order to accept both ng and ngfft
1589 !20. nspden
1590  if(ddbvrs_is_current_or_old)then
1591    read (unddb, '(1x,a9,i10)' )name(1),nspden
1592    call ddb_chkname(name(1),'   nspden')
1593  else
1594    nspden=0
1595  end if
1596 !21. nspinor
1597  if(ddbvrs_is_current_or_old)then
1598    read (unddb, '(1x,a9,i10)' )name(1),nspinor
1599    call ddb_chkname(name(1),'  nspinor')
1600  else
1601    nspinor=0
1602  end if
1603 !22. occ
1604  if(occopt==2)then
1605    im=3
1606    do iline=1,(bantot+2)/3
1607      if(iline==(bantot+2)/3)im=bantot-3*(iline-1)
1608      if(ddbvrs_is_current_or_old)then
1609        read (unddb, '(1x,a9,3d22.14)' )name(1),(occ((iline-1)*3+ii),ii=1,im)
1610      else
1611        read (unddb, '(1x,a6,3d22.14)' )name_old,(occ((iline-1)*3+ii),ii=1,im) ; name(1)='   '//name_old
1612      end if
1613      if (iline==1) then
1614        call ddb_chkname(name(1),'      occ')
1615      else
1616        call ddb_chkname(name(1),'         ')
1617      end if
1618    end do
1619  else
1620    im=3
1621    do iline=1,(nband(1)+2)/3
1622      if(iline==(nband(1)+2)/3)im=nband(1)-3*(iline-1)
1623      if(ddbvrs_is_current_or_old)then
1624        read (unddb, '(1x,a9,3d22.14)' )name(1),(occ((iline-1)*3+ii),ii=1,im)
1625      else
1626        read (unddb, '(1x,a6,3d22.14)' )name_old,(occ((iline-1)*3+ii),ii=1,im) ; name(1)='   '//name_old
1627      end if
1628      if (iline==1) then
1629        call ddb_chkname(name(1),'      occ')
1630      else
1631        call ddb_chkname(name(1),'         ')
1632      end if
1633    end do
1634    if(nkpt>1)then
1635      do ikpt=2,nkpt
1636        do iband=1,nband(1)
1637          occ(iband+nband(1)*(ikpt-1))=occ(iband)
1638        end do
1639      end do
1640    end if
1641  end if
1642 !23. rprim
1643  do iline=1,3
1644    if(ddbvrs_is_current_or_old)then
1645      read (unddb, '(1x,a9,3d22.14)' )name(1),(rprim(ii,iline),ii=1,3)
1646    else
1647      read (unddb, '(1x,a6,3d22.14)' )name_old,(rprim(ii,iline),ii=1,3) ; name(1)='   '//name_old
1648    end if
1649    if (iline==1) then
1650      call ddb_chkname(name(1),'    rprim')
1651    else
1652      call ddb_chkname(name(1),'         ')
1653    end if
1654  end do
1655 !24. dfpt_sciss
1656  if(ddbvrs_is_current_or_old)then
1657    read (unddb, '(1x,a11,d22.14)' )name(1),dfpt_sciss
1658    call ddb_chkname(name(1),'dfpt_sciss', 'sciss')
1659  else
1660    read (unddb, '(1x,a6,d22.14)' )name_old,dfpt_sciss ; name(1)=name_old
1661    call ddb_chkname(name(1),'sciss')
1662  end if
1663 !25. spinat
1664  if(ddbvrs_is_current_or_old)then
1665    do iline=1,natom
1666      read (unddb, '(1x,a9,3d22.14)' )name(1),(spinat(ii,iline),ii=1,3)
1667      if (iline==1) then
1668        call ddb_chkname(name(1),'   spinat')
1669      else
1670        call ddb_chkname(name(1),'         ')
1671      end if
1672    end do
1673  else
1674 !  spinat is set to zero by default in mrgddb.f
1675 !  spinat(:,1:natom)=zero
1676  end if
1677 !26. symafm
1678  if(ddbvrs_is_current_or_old)then
1679    im=12
1680    do iline=1,(nsym+11)/12
1681      if(iline==(nsym+11)/12)im=nsym-12*(iline-1)
1682      read (unddb, '(1x,a9,5x,12i5)' )name(1),(symafm((iline-1)*12+ii),ii=1,im)
1683      if (iline==1) then
1684        call ddb_chkname(name(1),'   symafm')
1685      else
1686        call ddb_chkname(name(1),'         ')
1687      end if
1688    end do
1689  else
1690 !  symafm is set to 1 by default in mrgddb.f
1691 !  symafm(1:nsym)=1
1692  end if
1693 !27. symrel
1694  do iline=1,nsym
1695    if(ddbvrs_is_current_or_old)then
1696      read (unddb, '(1x,a9,5x,9i5)' )name(1),((symrel(ii,ij,iline),ii=1,3),ij=1,3)
1697    else
1698      read (unddb, '(1x,a6,5x,9i5)' )name_old,&
1699 &     ((symrel(ii,ij,iline),ii=1,3),ij=1,3) ; name(1)='   '//name_old
1700    end if
1701    if (iline==1) then
1702      call ddb_chkname(name(1),'   symrel')
1703    else
1704      call ddb_chkname(name(1),'         ')
1705    end if
1706  end do
1707 !28old. xred
1708  if(.not.ddbvrs_is_current_or_old)then
1709    do iline=1,natom
1710      read (unddb, '(1x,a6,3d22.14)' )name(1),(xred(ii,iline),ii=1,3)
1711    end do
1712 !  No check of name, to allow the old tn
1713  end if
1714 !28. tnons
1715  do iline=1,nsym
1716    if(ddbvrs_is_current_or_old)then
1717      read (unddb, '(1x,a9,3d22.14)' )name(1),(tnons(ii,iline),ii=1,3)
1718    else
1719      read (unddb, '(1x,a6,3d22.14)' )name_old,(tnons(ii,iline),ii=1,3) ; name(1)='   '//name_old
1720    end if
1721    if (iline==1) then
1722      call ddb_chkname(name(1),'    tnons')
1723    else
1724      call ddb_chkname(name(1),'         ')
1725    end if
1726  end do
1727 !29. tolwfr
1728  if(ddbvrs_is_current_or_old)then
1729    read (unddb, '(1x,a9,d22.14)' )name(1),tolwfr
1730  end if
1731 !Do not check the name, in order to allow both tolwfr and wftol
1732 !30. tphysel
1733  if(ddbvrs_is_current_or_old)then
1734    read (unddb, '(1x,a9,d22.14)' )name(1),tphysel
1735    call ddb_chkname(name(1),'  tphysel')
1736  else
1737    tphysel=zero
1738  end if
1739 !31. tsmear
1740  if(ddbvrs_is_current_or_old)then
1741    read (unddb, '(1x,a9,d22.14)' )name(1),tsmear
1742    call ddb_chkname(name(1),'   tsmear')
1743  else
1744    tsmear=zero
1745  end if
1746 !32. typat
1747  im=12
1748  do iline=1,(natom+11)/12
1749    if(iline==(natom+11)/12)im=natom-12*(iline-1)
1750    if(ddbvrs_is_current_or_old)then
1751      read (unddb, '(1x,a9,5x,12i5)' )name(1),(typat((iline-1)*12+ii),ii=1,im)
1752    else
1753      read (unddb, '(1x,a6,5x,12i5)' )name_old,(typat((iline-1)*12+ii),ii=1,im) ; name(1)='   '//name_old
1754    end if
1755    if (iline==1) then
1756 !    Both type and typat are allowed => no check
1757 !    call ddb_chkname(name(1),'    typat')
1758    else
1759      call ddb_chkname(name(1),'         ')
1760    end if
1761  end do
1762 !33old. tolwfr
1763  if(.not.ddbvrs_is_current_or_old)then
1764    read (unddb, '(1x,a6,d22.14)' )name(1),tolwfr
1765  end if
1766 !Do not check the name, in order to allow both tolwfr and wftol
1767 !33. wtk
1768  im=3
1769  do iline=1,(nkpt+2)/3
1770    if(iline==(nkpt+2)/3)im=nkpt-3*(iline-1)
1771    if(ddbvrs_is_current_or_old)then
1772      read (unddb, '(1x,a9,3d22.14)' )name(1),(wtk((iline-1)*3+ii),ii=1,im)
1773    else
1774      read (unddb, '(1x,a6,3d22.14)' )name_old,(wtk((iline-1)*3+ii),ii=1,im) ; name(1)='   '//name_old
1775    end if
1776    if (iline==1) then
1777      call ddb_chkname(name(1),'      wtk')
1778    else
1779      call ddb_chkname(name(1),'         ')
1780    end if
1781  end do
1782 !34. xred
1783  if(ddbvrs_is_current_or_old)then
1784    do iline=1,natom
1785      read (unddb, '(1x,a9,3d22.14)' )name(1),(xred(ii,iline),ii=1,3)
1786      if (iline==1) then
1787        call ddb_chkname(name(1),'     xred')
1788      else
1789        call ddb_chkname(name(1),'         ')
1790      end if
1791    end do
1792  end if
1793 !35. znucl
1794  if(ddbvrs_is_current_or_old)then
1795    im=3
1796    do iline=1,(ntypat+2)/3
1797      if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1)
1798      read (unddb, '(1x,a9,3d22.14)' )name(1),(znucl((iline-1)*3+ii),ii=1,im)
1799      if (iline==1) then
1800        call ddb_chkname(name(1),'    znucl')
1801      else
1802        call ddb_chkname(name(1),'         ')
1803      end if
1804    end do
1805  else
1806 !  znucl is set to zero by default in mrgddb.f
1807 !  znucl(:)=zero
1808  end if
1809 !36. zion
1810  im=3
1811  do iline=1,(ntypat+2)/3
1812    if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1)
1813    if(ddbvrs_is_current_or_old)then
1814      read (unddb, '(1x,a9,3d22.14)' )name(1),(zion((iline-1)*3+ii),ii=1,im)
1815    else
1816      read (unddb, '(1x,a6,3d22.14)' )name_old,(zion((iline-1)*3+ii),ii=1,im) ; name(1)='   '//name_old
1817    end if
1818    if (iline==1) then
1819 !    Do not check the names, to allow both zion and znucl - the latter for 990527 format
1820 !    call ddb_chkname(name(1),'     zion')
1821    else
1822      call ddb_chkname(name(1),'         ')
1823    end if
1824  end do
1825 
1826 end subroutine ioddb8_in

m_ddb_hdr/psddb8 [ Functions ]

[ Top ] [ m_ddb_hdr ] [ Functions ]

NAME

 psddb8

FUNCTION

 Take care of the i/o of pseudopotentials for the
 Derivative DataBase, and also the number of data blocks.

INPUTS

  choice=(1 => read), (2=> write)
  dimekb=dimension of ekb (contains Kleimann-Bylander energies)
         used only for norm-conserving pseudopotentials
  lmnmax=if useylm=1, max number of (l,m,n) comp. over all type of psps
        =if useylm=0, max number of (l,n)   comp. over all type of psps
  nunit=unit number for the Derivative DataBase.
  ntypat=number of atom types
  pspso(ntypat)=For each type of psp, 1 if no spin-orbit component is taken
     into account, 2 if a spin-orbit component is used
  usepaw= 0 for non paw calculation; =1 for paw calculation
  useylm=governs the way the nonlocal operator is to be applied:
         1=using Ylm, 0=using Legendre polynomials
  vrsddb=Derivative Database version, for check of compatibility.

OUTPUT

  (see side effects)

SIDE EFFECTS

  indlmn(6,i,ntypat)= array giving l,m,n,lm,ln,s for i=ln  (if useylm=0)
                                                 or i=lmn  (if useylm=1)
  ekb(dimekb,ntypat)= (norm-conserving psps only) (Real) Kleinman-Bylander energies (hartree)
                      Presently the only characteristics of the psp
  fullinit=0 if the ekb are not available, at input as well as at output
  pawtab(ntypat*usepaw)= (PAW only) PAW datasets characteristics
                  Presently only pawtab%basis_size,pawtab%lmn_size,pawtab%shape_type
                  pawtab%rpaw,pawtab%rshp,pawtab%dij0  are used
  nblok=number of blocks

NOTES

 Only executed by one processor

PARENTS

      m_ddb_hdr

CHILDREN

SOURCE

 803 subroutine psddb8 (choice,dimekb,ekb,fullinit,indlmn,lmnmax,&
 804 &          nblok,ntypat,nunit,pawtab,pspso,usepaw,useylm,vrsddb)
 805 
 806 
 807 !This section has been created automatically by the script Abilint (TD).
 808 !Do not modify the following lines by hand.
 809 #undef ABI_FUNC
 810 #define ABI_FUNC 'psddb8'
 811 !End of the abilint section
 812 
 813  implicit none
 814 
 815 !Arguments -------------------------------
 816 !scalars
 817  integer,intent(in) :: choice,dimekb,lmnmax,ntypat,nunit,usepaw,useylm
 818  integer,intent(in) :: vrsddb
 819  integer,intent(inout) :: fullinit,nblok
 820 !arrays
 821  integer,intent(in) :: pspso(ntypat)
 822  integer,intent(inout) :: indlmn(6,lmnmax,ntypat)
 823  real(dp),intent(inout) :: ekb(dimekb,ntypat)
 824  type(pawtab_type),intent(inout) :: pawtab(ntypat*usepaw)
 825 
 826 !Local variables -------------------------
 827 !Set the version number
 828 !scalars
 829  integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527
 830  integer :: basis_size0,dimekb0,iekb,ii,ij,il,ilm,ilmn,iln,iln0,im,ios,iproj,iproj0,itypat,itypat0
 831  integer :: jekb,jlmn,jln,lmnmax0,lmn_size0,lmn2_size0,lpsang,nekb,nproj,npsang,pspso0,shape_type0
 832  integer :: usepaw0,vrspsp8
 833  real(dp) :: rpaw0,rshape0
 834  character(len=12) :: string
 835  character(len=500) :: message
 836 !arrays
 837  integer,allocatable :: i1(:),i2(:),nprj(:),orbitals(:)
 838  real(dp),allocatable :: dij0(:),ekb0(:,:)
 839 
 840 ! *********************************************************************
 841 
 842 !Check psddb8 version number (vrsio8) against DDB version number (vrsddb)
 843  if (vrsio8/=vrsddb) then
 844    write(message, '(a,i10,a,a,i10,a)' )&
 845 &   'the psddb8 DDB version number=',vrsio8,ch10,&
 846 &   'is not equal to the calling code DDB version number=',vrsddb,'.'
 847    MSG_WARNING(message)
 848  end if
 849 
 850 !Check the value of choice
 851  if (choice<=0.or.choice>=3) then
 852    write(message, '(a,a,a,i10,a)' )&
 853 &   'The permitted values for choice are 1 or 2.',ch10,&
 854 &   'The calling routine asks ',choice,'.'
 855    MSG_BUG(message)
 856  end if
 857 
 858 !==================================================================================
 859 !First option: read psp characteristic from file =================================
 860 !==================================================================================
 861  if (choice==1) then
 862 
 863    read(nunit,*)
 864    read(nunit, '(a12)' )string
 865    fullinit=1 ; if (lmnmax>0) indlmn(:,:,:)=0
 866 
 867 !  --------------------------------------------
 868 !  -----  NEW FORMAT (NCPP+PAW) ---------------
 869 !  --------------------------------------------
 870    if (string=='  Descriptio')then
 871 !    ==============================
 872 !    ======== Common data =========
 873 !    ==============================
 874      read(nunit, '(32x,i6)' )vrspsp8
 875      if (vrspsp8==vrsio8_old.or.vrspsp8==vrsio8_old_old) then
 876        usepaw0=0
 877        read(nunit, '(10x,i3,14x,i3,11x,i3)', iostat=ios)dimekb0,lmnmax0,usepaw0
 878        if(ios/=0)then
 879          backspace(nunit)
 880          read (nunit, '(10x,i3,14x,i3)' )dimekb0,lmnmax0
 881          usepaw0=0
 882        end if
 883      else if (vrspsp8==vrsio8) then
 884        read(nunit, '(10x,i3)') usepaw0
 885        if (usepaw/=usepaw0) then
 886          write(message, '(a,i1,a,i1,a)' )&
 887 &         'usepaw is announced to be ',usepaw,' but read usepaw is ',usepaw0,' !'
 888          MSG_ERROR(message)
 889        end if
 890        if (usepaw==0) then
 891          read (nunit, '(10x,i3,14x,i3)' )dimekb0,lmnmax0
 892        end if
 893      end if
 894 
 895 !    ==============================
 896 !    === Norm-conserving psps =====
 897 !    ==============================
 898      if (usepaw==0) then
 899        ekb(:,:)=zero
 900        ABI_MALLOC(ekb0,(dimekb,dimekb))
 901        do itypat=1,ntypat
 902          read(nunit, '(13x,i4,9x,i3,8x,i4)' )itypat0,pspso0,nekb
 903 !        Check the compatibility with the main code dimensioning
 904          if(nekb>dimekb)then
 905            write(message, '(a,i8,a,a,a,i3,a)' )&
 906 &           '  ',nekb,' components of ekb are announced',ch10,&
 907 &           'but dimekb=',dimekb,'.'
 908            MSG_BUG(message)
 909          end if
 910          read(nunit,*)
 911          ilmn=0;iproj0=0
 912          do iekb=1,nekb
 913            read(nunit, '(3i6,3x,8d15.7)' ) iln,lpsang,iproj,(ekb0(ii,iekb),ii=1,min(nekb,4))
 914            if(nekb>4)then
 915              do jekb=5,nekb,4
 916                read(nunit, '(21x,8d15.7)' )(ekb0(ii,iekb),ii=jekb,min(nekb,jekb+3))
 917              end do
 918            end if
 919            if (lpsang==0.and.iproj>iproj0) iproj0=iproj
 920            if (useylm==1) then
 921              do im=-lpsang,lpsang
 922                ilmn=ilmn+1
 923                indlmn(1,ilmn,itypat)=lpsang
 924                indlmn(2,ilmn,itypat)=im
 925                indlmn(3,ilmn,itypat)=iproj
 926                indlmn(4,ilmn,itypat)=lpsang**2+lpsang+1+im
 927                indlmn(5,ilmn,itypat)=iln
 928                indlmn(6,ilmn,itypat)=1
 929                if (pspso0/=1.and.iln>(nekb-iproj0)/2) indlmn(6,ilmn,itypat)=2
 930              end do
 931            else
 932              ilmn=ilmn+1
 933              indlmn(1,ilmn,itypat)=lpsang
 934              indlmn(2,ilmn,itypat)=lpsang
 935              indlmn(3,ilmn,itypat)=iproj
 936              indlmn(4,ilmn,itypat)=lpsang**2+lpsang+1
 937              indlmn(5,ilmn,itypat)=iln
 938              indlmn(6,ilmn,itypat)=1
 939              if (pspso0/=1.and.iln>(nekb-iproj0)/2) indlmn(6,ilmn,itypat)=2
 940            end if
 941 !          For the time being, only diagonal ekb are treated in abinit v3
 942            ekb(iekb,itypat)=ekb0(iekb,iekb)
 943 !          For non-diagonal ekb, one could use:
 944 !          do jekb=iekb to nekb
 945 !          ekb(jekb+iekb*(iekb-1)/2,itypat)=ekb0(jekb,iekb)
 946 !          end do
 947          end do
 948        end do
 949        ABI_FREE(ekb0)
 950 
 951 !      ==============================
 952 !      ============ PAW =============
 953 !      ==============================
 954      else
 955        do itypat=1,ntypat
 956          read(nunit, '(12x,i4,12x,i3,12x,i5)' )itypat0,basis_size0,lmn_size0
 957          lmn2_size0=lmn_size0*(lmn_size0+1)/2
 958          ABI_MALLOC(orbitals,(basis_size0))
 959          read(nunit, '(20x,50i2)' ) orbitals(1:basis_size0)
 960          read(nunit, '(11x,f6.3,13x,i2,11x,f6.3)' ) rpaw0,shape_type0,rshape0
 961          read(nunit,'(24x,i3)') nekb
 962          read(nunit,*)
 963          ABI_MALLOC(dij0,(nekb))
 964          ABI_MALLOC(i1,(nekb))
 965          ABI_MALLOC(i2,(nekb))
 966          do ii=1,nekb,4
 967            read(nunit,'(3x,4(1x,i4,1x,i4,1x,d12.5))') (i1(ij),i2(ij),dij0(ij),ij=ii,min(ii+3,nekb))
 968          end do
 969          if (lmn_size0>lmnmax) then
 970            write(message, '(a,i5,3a,i5,a)' )&
 971 &           'max. value of ',lmnmax,' for lmn_size is announced',ch10,&
 972 &           'but ',lmn_size0,' is read.'
 973            MSG_BUG(message)
 974          end if
 975          if (allocated(pawtab(itypat)%dij0)) then
 976            if (lmn_size0>pawtab(itypat)%lmn_size) then
 977              write(message, '(a,i5,3a,i5,a)' )&
 978 &             'lmn_size=,',pawtab(itypat)%lmn_size,' is announced',ch10,&
 979 &             'but ',lmn_size0,' is read.'
 980              MSG_BUG(message)
 981            end if
 982          end if
 983          ABI_MALLOC(nprj,(0:maxval(orbitals)))
 984          ilmn=0;nprj=0
 985          do iln=1,basis_size0
 986            il=orbitals(iln)
 987            nprj(il)=nprj(il)+1
 988            do ilm=1,2*il+1
 989              indlmn(1,ilmn+ilm,itypat)=il
 990              indlmn(2,ilmn+ilm,itypat)=ilm-(il+1)
 991              indlmn(3,ilmn+ilm,itypat)=nprj(il)
 992              indlmn(4,ilmn+ilm,itypat)=il*il+ilm
 993              indlmn(5,ilmn+ilm,itypat)=iln
 994              indlmn(6,ilmn+ilm,itypat)=1
 995            end do
 996            ilmn=ilmn+2*il+1
 997          end do
 998          pawtab(itypat)%basis_size=basis_size0
 999          pawtab(itypat)%lmn_size  =lmn_size0
1000          pawtab(itypat)%lmn2_size =lmn2_size0
1001          pawtab(itypat)%shape_type=shape_type0
1002          pawtab(itypat)%rpaw      =rpaw0
1003          pawtab(itypat)%rshp      =rshape0
1004          if (.not.allocated(pawtab(itypat)%dij0))  then
1005            ABI_MALLOC(pawtab(itypat)%dij0,(lmn2_size0))
1006          end if
1007          pawtab(itypat)%dij0(1:lmn2_size0)=zero
1008          do ii=1,nekb
1009            ij=i1(ii)+i2(ii)*(i2(ii)-1)/2
1010            pawtab(itypat)%dij0(ij)=dij0(ii)
1011          end do
1012          ABI_FREE(nprj)
1013          ABI_FREE(orbitals)
1014          ABI_FREE(dij0)
1015          ABI_FREE(i1)
1016          ABI_FREE(i2)
1017        end do
1018 
1019      end if ! NCPP or PAW
1020 
1021 !    --------------------------------------------
1022 !    -----  OLD FORMAT (NCPP only) --------------
1023 !    --------------------------------------------
1024    else if (string==' Description')then
1025      if (usepaw==1) then
1026        MSG_BUG("old DDB pspformat not compatible with PAW")
1027      end if
1028 
1029      read (nunit, '(10x,i3,10x,i3)' )nproj,npsang
1030      nekb=nproj*npsang
1031 !    Check the compatibility with the main code dimensioning
1032      if(nekb>dimekb)then
1033        write(message, '(a,i8,a,a,a,i3,a)' )&
1034 &       '  ',nekb,' components of ekb are announced',ch10,&
1035 &       'but the maximum is dimekb=',dimekb,'.'
1036        MSG_BUG(message)
1037      end if
1038      if(useylm/=0)then
1039        MSG_BUG('useylm must be 0 !')
1040      end if
1041 !    Read the data
1042      ABI_MALLOC(ekb0,(dimekb,dimekb))
1043      ekb0(:,:)=zero
1044      do itypat=1,ntypat
1045        read (nunit, '(13x,i4)' )ij
1046        do iproj=1,nproj
1047          read (nunit, '(6x,3d22.14)' )&
1048 &         (ekb0(iproj+nproj*(ii-1),iproj+nproj*(ii-1)),ii=1,min(npsang,3))
1049          if(npsang>3)read (nunit, '(6x,3d22.14)' )&
1050 &         (ekb0(iproj+nproj*(ii-1),iproj+nproj*(ii-1)),ii=4,npsang)
1051          do ii=1,npsang
1052            iekb=iproj+nproj*(ii-1)
1053            indlmn(1,iekb,itypat)=ii-1
1054            indlmn(2,iekb,itypat)=ii-1
1055            indlmn(3,iekb,itypat)=iproj
1056            indlmn(4,iekb,itypat)=ii**2-ii+1
1057            indlmn(5,iekb,itypat)=iekb
1058            indlmn(6,iekb,itypat)=1
1059 !          For the time being, only diagonal ekb are treated in abinit v3
1060            ekb(iekb,itypat)=ekb0(iekb,iekb)
1061          end do
1062        end do
1063      end do
1064      ABI_FREE(ekb0)
1065 
1066 !    --------------------------------------------
1067 !    -----  OTHER CASES -------------------------
1068 !    --------------------------------------------
1069    else if(string==' No informat')then
1070      fullinit=0
1071    else
1072      MSG_BUG('Error when reading the psp information')
1073    end if
1074 
1075 !  Now, the number of blocks
1076    read(nunit,*)
1077    read(nunit,*)
1078    read(nunit, '(24x,i4)' )nblok
1079 
1080 !  ==================================================================================
1081 !  Second option: read psp characteristic from file ================================
1082 !  ==================================================================================
1083  else if(choice==2)then
1084 
1085    write(nunit, '(a)' )' '
1086    if (fullinit==0)then
1087 !    This possibility is used when the DDB is initialized,
1088 !    and the ekb s are not available from the GS input file...
1089      write(nunit, '(a)' )' No information on the potentials yet '
1090    else
1091 
1092 !    ==============================
1093 !    === Norm-conserving psps =====
1094 !    ==============================
1095      if (usepaw==0) then
1096        write(nunit, '(a)' )'  Description of the potentials (KB energies)'
1097        write(nunit, '(a,i6)' )'  vrsio8 (for pseudopotentials)=',vrsio8
1098        write(nunit, '(a,i3)' ) '  usepaw =',usepaw
1099        write(nunit, '(a,i3,a,i3,a,i3)' )'  dimekb =',dimekb,'       lmnmax=',lmnmax
1100        ABI_MALLOC(ekb0,(dimekb,dimekb))
1101        do itypat=1,ntypat
1102 !        Compute nekb
1103          nekb=0
1104          do jlmn=1,lmnmax
1105            jln=indlmn(5,jlmn,itypat)
1106            if(jln>nekb)then
1107              nekb=jln
1108            end if
1109          end do
1110          write(nunit, '(a,i4,a,i3,a,i4)' ) &
1111 &         '  Atom type= ',itypat,'   pspso=',pspso(itypat),'   nekb=',nekb
1112          write(nunit, '(a)' ) '  iln lpsang iproj  ekb(:)'
1113          iln0=0
1114          ekb0(:,:)=zero
1115          do ilmn=1,lmnmax
1116            iln =indlmn(5,ilmn,itypat)
1117            if (iln>iln0) then
1118              iln0=iln
1119              lpsang=indlmn(1,ilmn,itypat)
1120              iproj=indlmn(3,ilmn,itypat)
1121 !            For the time being, only diagonal ekb are treated in abinit v3
1122              ekb0(iln,iln)=ekb(iln,itypat)
1123 !            For non-diagonal ekb, one could use:
1124 !            do ii=iln to nekb
1125 !            ekb0(ii,iln)=ekb(ii+iln*(iln-1)/2,itypat)
1126 !            end do
1127              write(nunit, '(3i6,3x,4es15.7)' ) iln,lpsang,iproj,(ekb0(ii,iln),ii=1,min(nekb,4))
1128              if(nekb>4)then
1129                do iekb=5,nekb,4
1130                  write(nunit, '(21x,4es15.7)' )(ekb0(ii,iekb),ii=iekb,min(nekb,iekb+3))
1131                end do
1132              end if
1133            end if
1134          end do
1135        end do
1136        ABI_FREE(ekb0)
1137 
1138 !      ==============================
1139 !      ============ PAW =============
1140 !      ==============================
1141      else
1142        write(nunit, '(a)' )'  Description of the PAW dataset(s)'
1143        write(nunit, '(a,i6)' )'  vrsio8 (for pseudopotentials)=',vrsio8
1144        write(nunit, '(a,i3)' ) '  usepaw =',usepaw
1145        do itypat=1,ntypat
1146          iln0=0
1147          ABI_MALLOC(orbitals,(pawtab(itypat)%basis_size))
1148          do ilmn=1,pawtab(itypat)%lmn_size
1149            iln =indlmn(5,ilmn,itypat)
1150            if (iln>iln0) then
1151              iln0=iln;orbitals(iln)=indlmn(1,ilmn,itypat)
1152            end if
1153          end do
1154          write(nunit, '(a,i4,a,i3,a,i5)' ) &
1155 &         '  Atom type=',itypat,' basis_size=',pawtab(itypat)%basis_size,&
1156 &         '   lmn_size=',pawtab(itypat)%lmn_size
1157          write(nunit, '(a,50i2)' ) &
1158 &         '    Basis functions=',orbitals(1:pawtab(itypat)%basis_size)
1159          write(nunit, '(a,f6.3,a,i2,a,f6.3)' ) &
1160 &         '    r_PAW= ',pawtab(itypat)%rpaw,' shape_type= ',pawtab(itypat)%shape_type,&
1161 &         '  r_shape= ',pawtab(itypat)%rshp
1162          nekb=0
1163          ABI_MALLOC(dij0,(pawtab(itypat)%lmn2_size))
1164          ABI_MALLOC(i1,(pawtab(itypat)%lmn2_size))
1165          ABI_MALLOC(i2,(pawtab(itypat)%lmn2_size))
1166          do jlmn=1,pawtab(itypat)%lmn_size
1167            ij=jlmn*(jlmn-1)/2
1168            do ilmn=1,jlmn
1169              if (abs(pawtab(itypat)%dij0(ij+ilmn))>tol16) then
1170                nekb=nekb+1;i1(nekb)=ilmn;i2(nekb)=jlmn
1171                dij0(nekb)=pawtab(itypat)%dij0(ij+ilmn)
1172              end if
1173            end do
1174          end do
1175          write(nunit,'(a,i3,a)') '    Dij0=     (only the ',nekb,' values different from zero)'
1176          write(nunit,'(2a)') '       i    j     Dij0        i    j     Dij0 ',&
1177 &         '       i    j     Dij0        i    j     Dij0'
1178          do ii=1,nekb,4
1179            write(nunit,'(3x,4(1x,i4,1x,i4,1x,es12.5))') (i1(ij),i2(ij),dij0(ij),ij=ii,min(ii+3,nekb))
1180          end do
1181          ABI_FREE(dij0)
1182          ABI_FREE(i1)
1183          ABI_FREE(i2)
1184          ABI_FREE(orbitals)
1185        end do
1186 
1187      end if ! NCPP or PAW
1188    end if ! fullinit==0
1189 
1190 !  Now, write the number of blocks
1191    write(nunit, '(a)' )' '
1192    write(nunit, '(a)' )' **** Database of total energy derivatives ****'
1193    write(nunit, '(a,i4)' ) ' Number of data blocks= ',nblok
1194  end if
1195 
1196 end subroutine psddb8