TABLE OF CONTENTS
- ABINIT/m_ddb_hdr
- m_ddb_hdr/chki8
- m_ddb_hdr/chkr8
- m_ddb_hdr/ddb_chkname
- m_ddb_hdr/ddb_getdims
- m_ddb_hdr/ddb_hdr_bcast
- m_ddb_hdr/ddb_hdr_bcast_dim
- m_ddb_hdr/ddb_hdr_close
- m_ddb_hdr/ddb_hdr_compare
- m_ddb_hdr/ddb_hdr_copy_missing_variables
- m_ddb_hdr/ddb_hdr_copy_psps_info
- m_ddb_hdr/ddb_hdr_free
- m_ddb_hdr/ddb_hdr_get_block_dims
- m_ddb_hdr/ddb_hdr_get_iomode
- m_ddb_hdr/ddb_hdr_init
- m_ddb_hdr/ddb_hdr_malloc
- m_ddb_hdr/ddb_hdr_open_read
- m_ddb_hdr/ddb_hdr_open_read_nc
- m_ddb_hdr/ddb_hdr_open_read_txt
- m_ddb_hdr/ddb_hdr_open_write
- m_ddb_hdr/ddb_hdr_open_write_nc
- m_ddb_hdr/ddb_hdr_open_write_txt
- m_ddb_hdr/ddb_hdr_print
- m_ddb_hdr/ddb_hdr_set_typ
- m_ddb_hdr/ddb_io_out
- m_ddb_hdr/inprep8
- m_ddb_hdr/ioddb8_in
- m_ddb_hdr/is_type_d0E
- m_ddb_hdr/is_type_d1E
- m_ddb_hdr/is_type_d2E
- m_ddb_hdr/is_type_d2eig
- m_ddb_hdr/is_type_d3E
- m_ddb_hdr/psddb8
ABINIT/m_ddb_hdr [ 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-2024 ABINIT group (GA, MJV, XG, MT, MM, MVeithen, MG, PB, JCC, MMignolet) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 MODULE m_ddb_hdr 24 25 use defs_basis 26 use m_errors 27 use m_abicore 28 use m_xmpi 29 use m_dtset 30 use m_crystal 31 use m_nctk 32 #ifdef HAVE_NETCDF 33 use netcdf 34 #endif 35 36 use defs_datatypes, only : pseudopotential_type 37 use m_copy, only : alloc_copy 38 use m_pawtab, only : pawtab_type, pawtab_nullify, pawtab_free, pawtab_bcast !, pawtab_copy 39 use m_psps, only : psps_copy, psps_free, psps_ncwrite, psps_ncread 40 use m_io_tools, only : open_file, get_unit, file_exists 41 use m_copy, only : alloc_copy 42 use m_fstrings, only : sjoin, endswith 43 use m_geometry, only : mkrdim 44 45 46 implicit none 47 48 private 49 50 public :: is_type_d0E ! Is this block type of a d0E kind? 51 public :: is_type_d1E ! Is this block type of a d1E kind? 52 public :: is_type_d2E ! Is this block type of a d2E kind? 53 public :: is_type_d3E ! Is this block type of a d3E kind? 54 public :: is_type_d2eig ! Is this block type of a d2eig kind? 55 56 ! GA: These 3 functions should no longer be public 57 public :: ddb_getdims ! Open a DDB file and read basic dimensions and variables. 58 public :: ioddb8_in ! Temporary 59 public :: psddb8 ! Temporary 60 61 62 ! Description of perturbations and block types 63 ! -------------------------------------------- 64 ! 65 ! GA: Variable names should not end with _d0E, _d1E, d2E, ... 66 ! Such names may affect the variable type, with certain compilers. 67 ! This is why I appended _xx, lack of a better suffix. 68 ! 69 ! 70 ! These parameters should be consistent with input variable rfmeth 71 ! and input parameter rftyp (e.g. ddb_get_block, ddb_read_block_txt). 72 ! Changing these values requires to change the documentation as well. 73 ! 74 integer,public,parameter :: BLKTYP_d0E_xx=0 ! Total energy 75 integer,public,parameter :: BLKTYP_d1E_xx=4 ! First-order derivatives of total energy 76 integer,public,parameter :: BLKTYP_d2E_ns=1 ! Second-order derivatives of total energy, non-stationary 77 integer,public,parameter :: BLKTYP_d2E_st=2 ! Second-order derivatives of total energy, stationary 78 integer,public,parameter :: BLKTYP_d2E_mbc=85 ! Molecular Berry curvature (MBC) 79 integer,public,parameter :: BLKTYP_d3E_xx=3 ! Third-order derivatives of total energy 80 integer,public,parameter :: BLKTYP_d3E_lw=33 ! Long-wave third-order derivatives of total energy 81 integer,public,parameter :: BLKTYP_d2eig_re=5 ! Second-order derivatives of eigenvalues 82 integer,public,parameter :: BLKTYP_d2eig_im=6 ! Static broadening of eigenvalues 83 84 integer,public,parameter :: descrlen=500 85 integer,public,parameter :: ntypes=9 ! The number of different block types 86 character(len=descrlen),public,parameter :: DESCR_d0E_xx = ' Total energy - ' 87 character(len=descrlen),public,parameter :: DESCR_d1E_xx = ' 1st derivatives - ' 88 character(len=descrlen),public,parameter :: DESCR_d2E_ns = ' 2nd derivatives (non-stat.) - ' 89 character(len=descrlen),public,parameter :: DESCR_d2E_st = ' 2nd derivatives (stationary) - ' 90 character(len=descrlen),public,parameter :: DESCR_d2E_mbc = ' 2nd derivatives (MBC) - ' 91 character(len=descrlen),public,parameter :: DESCR_d3E_xx = ' 3rd derivatives - ' 92 character(len=descrlen),public,parameter :: DESCR_d3E_lw = ' 3rd derivatives (long wave) - ' 93 character(len=descrlen),public,parameter :: DESCR_d2eig_re = ' 2nd eigenvalue derivatives - ' 94 character(len=descrlen),public,parameter :: DESCR_d2eig_im = ' 2nd eigenvalue derivatives (imaginary part) - ' 95 96 ! These describe the perturbation indices above natom 97 ! GA: It would be nice if different ddb instance could order perturbations 98 ! differently. This could potentially reduce the ddb size. 99 character(len=descrlen),public,parameter :: DESCR_ipert_0 = 'displacement of atom ' 100 character(len=descrlen),public,parameter :: DESCR_ipert_1 = 'derivative wrt k' 101 character(len=descrlen),public,parameter :: DESCR_ipert_2 = 'electric field' 102 character(len=descrlen),public,parameter :: DESCR_ipert_3 = 'strain' ! diagonal stress components 103 character(len=descrlen),public,parameter :: DESCR_ipert_4 = 'strain' ! non-diagonal stress 104 character(len=descrlen),public,parameter :: DESCR_ipert_5 = 'magnetic field' 105 character(len=descrlen),public,parameter :: DESCR_ipert_10 = '2nd derivative wrt to k' 106 character(len=descrlen),public,parameter :: DESCR_ipert_11 = '2nd derivative wrt to k and electric field' 107 108 integer,public,parameter :: DDB_VERSION=20230401 ! TODO: check if we should update this with new G matrix stuff 109 ! DDB Version number for text format. 110 111 integer,public,parameter :: DDB_VERSION_NC=20230219 ! TODO: check if we should update this with new G matrix stuff 112 ! DDB NetCDF version number. 113 114 type,public :: ddb_hdr_type 115 116 logical :: has_open_file_txt=.false. ! Has an open file in text format 117 logical :: has_open_file_nc=.false. ! Has an open file in netcdf format 118 119 integer :: iomode=IO_MODE_FORTRAN 120 121 integer :: unddb ! Unit for open file in text format 122 integer :: ncid ! Unit for open file in netcdf format 123 124 ! The variables with_psps / with_dfpt_vars 125 ! Must be identical in text format, but may differ in netcdf. 126 integer :: with_psps=0 ! Whether info on pseudopotentials is present. 127 integer :: with_dfpt_vars=0 ! Whether this ddb comes from a dfpt calculation 128 ! 0 -> comes from a ground state calculation. 129 ! 1 -> comes from a dfpt calculation 130 ! and k-point may or may not be reduced 131 ! by time reversal symmetry. 132 133 integer :: ddb_version ! Version of the DDB file 134 135 integer :: matom 136 integer :: mband 137 integer :: mkpt 138 integer :: msym ! GA: Do we need both msym and nsym ? 139 ! ddb_hdr reading is first called with msym=192 140 ! for dimensions only, then it is called with 141 ! actual number of sym. 142 integer :: mtypat 143 integer :: mpert 144 integer :: msize 145 146 integer :: intxc 147 integer :: iscf 148 integer :: ixc 149 integer :: natom 150 integer :: nkpt 151 integer :: nspden 152 integer :: nspinor 153 integer :: nsppol 154 integer :: nsym 155 integer :: ntypat 156 integer :: occopt 157 integer :: usepaw 158 159 integer :: nblok ! Number of blocks in the ddb 160 161 ! GA: FIXME mblktyp doesnt make sense. 162 ! Instead, we should have logical variables has_d2eig, has_d3E, etc. 163 ! MMig; I agree mblktyp does not make much sense and its not even always 164 ! treated properly throughout the code 165 ! The logical variables seems like a good idea to me 166 integer :: mblktyp ! Max block type 167 ! 0 = Total energy 168 ! 1 = 2nd derivatives (non-stat.) 169 ! 2 = 2nd derivatives (stationary) 170 ! 3 = 3rd derivatives 171 ! 4 = 1st derivatives 172 ! 5 = 2nd eigenvalue derivatives 173 real(dp) :: dilatmx 174 real(dp) :: ecut 175 real(dp) :: ecutsm 176 real(dp) :: kptnrm 177 real(dp) :: pawecutdg 178 real(dp) :: dfpt_sciss 179 real(dp) :: tolwfr 180 real(dp) :: tphysel 181 real(dp) :: tsmear 182 183 character(len=fnlen) :: dscrpt 184 185 integer :: ngfft(18) 186 real(dp) :: acell(3) 187 real(dp) :: rprim(3,3) 188 189 integer,allocatable :: nband(:) 190 ! nband(mkpt*nsppol) 191 192 integer,allocatable :: symafm(:) 193 ! symafm(msym) 194 195 integer,allocatable :: symrel(:,:,:) 196 ! symrel(3,3,msym) 197 198 integer,allocatable :: typat(:) 199 ! typat(matom) 200 201 integer,allocatable :: typ(:) 202 ! typ(nblok), type of the blocks: d2E_st, d3E_lw,... 203 204 real(dp),allocatable :: amu(:) 205 ! amu(mtypat) 206 207 real(dp),allocatable :: kpt(:,:) 208 ! kpt(3,mkpt) 209 210 real(dp),allocatable :: occ(:) 211 ! occ(mband*mkpt*nsppol) 212 213 real(dp),allocatable :: spinat(:,:) 214 ! spinat(3,matom) 215 216 real(dp),allocatable :: tnons(:,:) 217 ! tnons(3,msym) 218 219 real(dp),allocatable :: wtk(:) 220 ! wtk(mkpt) 221 222 real(dp),allocatable :: xred(:,:) 223 ! xred(3,matom) 224 225 real(dp),allocatable :: zion(:) 226 ! zion(mtypat) 227 228 real(dp),allocatable :: znucl(:) 229 ! znucl(mtypat) 230 231 type(pawtab_type),allocatable :: pawtab(:) 232 ! pawtab(psps%ntypat*psps%usepaw) 233 234 type(pseudopotential_type) :: psps 235 236 type(crystal_t) :: crystal 237 238 contains 239 240 procedure :: init => ddb_hdr_init 241 ! Initialize object 242 243 procedure :: malloc => ddb_hdr_malloc 244 ! Allocate dynamic memory. 245 246 procedure :: free => ddb_hdr_free 247 ! Free dynamic memory. 248 249 procedure :: bcast_dim => ddb_hdr_bcast_dim 250 ! Master broadcasts header dimensions. 251 252 procedure :: bcast => ddb_hdr_bcast 253 ! Master broadcasts header data. 254 255 procedure :: print => ddb_hdr_print 256 ! Print out the content of the header. 257 258 procedure :: compare => ddb_hdr_compare 259 ! Compare two DDB headers. 260 261 procedure :: copy_missing_variables => ddb_hdr_copy_missing_variables 262 ! Copy some missing variables from one header to another 263 264 procedure :: copy_psps_from => ddb_hdr_copy_psps_from 265 ! Copy the information on pseudopotentials from one header to another 266 267 procedure :: set_typ => ddb_hdr_set_typ 268 ! Set the type of each blok in the ddb 269 270 procedure :: get_block_dims => ddb_hdr_get_block_dims 271 ! Compute mpert and msize 272 273 procedure :: open_write => ddb_hdr_open_write 274 ! Open the DDB file and write the header. 275 276 procedure :: open_write_txt => ddb_hdr_open_write_txt 277 ! Open the DDB text file and write the header. 278 279 procedure :: open_write_nc => ddb_hdr_open_write_nc 280 ! Open the DDB NetCDF file and write the header. 281 282 procedure :: open_read => ddb_hdr_open_read 283 ! Open the DDB file and read the header. 284 285 procedure :: open_read_txt => ddb_hdr_open_read_txt 286 ! Open the DDB text file and read the header. 287 288 procedure :: open_read_nc => ddb_hdr_open_read_nc 289 ! Open the DDB NetCDF file and read the header. 290 291 procedure :: close => ddb_hdr_close 292 ! Close any open file 293 294 procedure :: get_iomode => ddb_hdr_get_iomode 295 ! Decide on iomode based on file name and current iomode parameter 296 297 end type ddb_hdr_type 298 299 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)
SOURCE
4707 subroutine chki8(inti,intt,name) 4708 4709 !Arguments ------------------------------- 4710 !scalars 4711 integer,intent(in) :: inti,intt 4712 character(len=6),intent(in) :: name 4713 4714 !Local variables------------------------------- 4715 !scalars 4716 character(len=500) :: message 4717 4718 ! ********************************************************************* 4719 4720 if(inti/=intt) then 4721 write(message, '(a,a,a,a,a,i10,a,a,a,i10,a,a,a)' )& 4722 'Comparing integers for variable',name,'.',ch10,& 4723 'Value from input DDB is',inti,' and',ch10,& 4724 'from transfer DDB is',intt,'.',ch10,& 4725 'Action: check your DDBs.' 4726 ABI_ERROR(message) 4727 end if 4728 4729 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)
SOURCE
4659 subroutine chkr8(reali,realt,name,tol) 4660 4661 !Arguments ------------------------------- 4662 !scalars 4663 real(dp),intent(in) :: reali,realt,tol 4664 character(len=6),intent(in) :: name 4665 4666 !Local variables------------------------------- 4667 !scalars 4668 character(len=500) :: message 4669 4670 ! ********************************************************************* 4671 4672 if(abs(reali-realt)>tol) then 4673 write(message, '(a,a,a,a,a,es16.6,a,a,a,es16.6,a,a,a)' )& 4674 'Comparing reals for variable',name,'.',ch10,& 4675 'Value from input DDB is',reali,' and',ch10,& 4676 'from transfer DDB is',realt,'.',ch10,& 4677 'Action: check your DDBs.' 4678 ABI_ERROR(message) 4679 end if 4680 4681 end subroutine chkr8
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
SOURCE
4553 subroutine ddb_chkname(nmfond,nmxpct,nmxpct2) 4554 4555 !Arguments ------------------------------- 4556 !scalars 4557 character(len=*),intent(in) :: nmfond,nmxpct 4558 character(len=*),intent(in),optional :: nmxpct2 4559 4560 !Local variables------------------------------- 4561 !scalars 4562 logical :: found 4563 character(len=500) :: nmfond_,nmxpct_,nmxpct2_ 4564 character(len=500) :: message 4565 4566 ! ********************************************************************* 4567 4568 nmxpct_ = trim(adjustl(nmxpct)) 4569 nmfond_ = trim(adjustl(nmfond)) 4570 4571 found = (nmxpct_ == nmfond_) 4572 4573 if (present(nmxpct2) .and. .not. found) then 4574 nmxpct2_ = trim(adjustl(nmxpct2)) 4575 found = (nmxpct2_==nmfond_) 4576 end if 4577 4578 if (.not. found) then 4579 write(message, '(a,a,a,a,a,a,a,a,a,a,a)' )& 4580 & 'Reading DDB, expected name was "',trim(nmxpct_),'"',ch10,& 4581 & ' and name found is "',trim(nmfond_),'"',ch10,& 4582 & 'Likely your DDB is incorrect.',ch10,& 4583 & 'Action: correct your DDB, or contact the ABINIT group.' 4584 ABI_ERROR(message) 4585 end if 4586 4587 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, and close the file.
INPUTS
character(len=*) filename: name of input or output file unddb=unit number for input or output
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.
SOURCE
3968 subroutine ddb_getdims(filename,comm,dimekb,lmnmax,mband,mblktyp,msym,natom,nblok,nkpt,nsppol,ntypat,usepaw) 3969 3970 !Arguments ------------------------------- 3971 !scalars 3972 character(len=*),intent(in) :: filename 3973 integer,intent(in) :: comm 3974 integer,intent(out) :: msym,dimekb,lmnmax,mband,mblktyp,natom,nblok,nkpt,ntypat,nsppol,usepaw 3975 3976 !Local variables------------------------------- 3977 !scalars 3978 integer,parameter :: master=0 3979 integer :: ierr,unddb 3980 3981 ! ********************************************************************* 3982 3983 ! Master node reads dims from file and then broadcast. 3984 if (xmpi_comm_rank(comm) == master) then 3985 unddb = get_unit() 3986 call inprep8(filename,unddb,dimekb,lmnmax,mband,mblktyp,msym,natom,nblok,nkpt,nsppol,ntypat,usepaw) 3987 end if 3988 3989 ! couldnt we use ddb_hdr_bcast_dim here? 3990 if (xmpi_comm_size(comm) > 1) then 3991 call xmpi_bcast(dimekb, master, comm, ierr) 3992 call xmpi_bcast(lmnmax, master, comm, ierr) 3993 call xmpi_bcast(mband, master, comm, ierr) 3994 call xmpi_bcast(mblktyp, master, comm, ierr) 3995 call xmpi_bcast(msym, master, comm, ierr) 3996 call xmpi_bcast(natom, master, comm, ierr) 3997 call xmpi_bcast(nblok, master, comm, ierr) 3998 call xmpi_bcast(nkpt, master, comm, ierr) 3999 call xmpi_bcast(ntypat, master, comm, ierr) 4000 call xmpi_bcast(nsppol, master, comm, ierr) 4001 call xmpi_bcast(usepaw, master, comm, ierr) 4002 end if 4003 4004 ! Maximum number of perturbations and size of matrix. 4005 !mpert=natom+6 4006 !msize=3*mpert*3*mpert; if (mblktyp==3) msize=msize*3*mpert 4007 4008 end subroutine ddb_getdims
m_ddb_hdr/ddb_hdr_bcast [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_bcast
FUNCTION
Master broadcasts the data and others allocate their arrays
INPUTS
comm=MPI communicator.
SOURCE
2693 subroutine ddb_hdr_bcast(ddb_hdr, comm) 2694 2695 !Arguments ------------------------------------ 2696 class(ddb_hdr_type),intent(inout) :: ddb_hdr 2697 integer, intent(in) :: comm 2698 2699 !Local variables ------------------------- 2700 integer, parameter :: master=0 2701 integer :: ierr 2702 integer :: ii,nn 2703 2704 ! ************************************************************************ 2705 2706 if (xmpi_comm_size(comm) == 1) return 2707 2708 DBG_ENTER("COLL") 2709 2710 call ddb_hdr%bcast_dim(comm) 2711 2712 ! Allocate arrays on the other nodes. 2713 if (xmpi_comm_rank(comm) /= master) then 2714 2715 call ddb_hdr%free() 2716 call ddb_hdr%malloc() 2717 2718 ABI_MALLOC(ddb_hdr%psps%indlmn,(6,ddb_hdr%psps%lmnmax,ddb_hdr%mtypat)) 2719 ABI_MALLOC(ddb_hdr%psps%pspso,(ddb_hdr%mtypat)) 2720 ABI_MALLOC(ddb_hdr%psps%ekb,(ddb_hdr%psps%dimekb,ddb_hdr%mtypat * (1 - ddb_hdr%psps%usepaw))) 2721 2722 ABI_MALLOC(ddb_hdr%pawtab,(ddb_hdr%psps%ntypat*ddb_hdr%psps%usepaw)) 2723 call pawtab_nullify(ddb_hdr%pawtab) 2724 2725 ddb_hdr%dscrpt = '' 2726 2727 end if 2728 2729 ! Floats 2730 call xmpi_bcast(ddb_hdr%dilatmx, master, comm, ierr) 2731 call xmpi_bcast(ddb_hdr%ecut, master, comm, ierr) 2732 call xmpi_bcast(ddb_hdr%kptnrm, master, comm, ierr) 2733 call xmpi_bcast(ddb_hdr%pawecutdg, master, comm, ierr) 2734 call xmpi_bcast(ddb_hdr%dfpt_sciss, master, comm, ierr) 2735 call xmpi_bcast(ddb_hdr%tolwfr, master, comm, ierr) 2736 call xmpi_bcast(ddb_hdr%tphysel, master, comm, ierr) 2737 call xmpi_bcast(ddb_hdr%tsmear, master, comm, ierr) 2738 2739 ! Arrays 2740 call xmpi_bcast(ddb_hdr%ngfft, master, comm, ierr) 2741 call xmpi_bcast(ddb_hdr%acell, master, comm, ierr) 2742 call xmpi_bcast(ddb_hdr%rprim, master, comm, ierr) 2743 2744 call xmpi_bcast(ddb_hdr%nband, master, comm, ierr) 2745 call xmpi_bcast(ddb_hdr%symafm, master, comm, ierr) 2746 call xmpi_bcast(ddb_hdr%symrel, master, comm, ierr) 2747 call xmpi_bcast(ddb_hdr%typat, master, comm, ierr) 2748 2749 call xmpi_bcast(ddb_hdr%amu, master, comm, ierr) 2750 call xmpi_bcast(ddb_hdr%kpt, master, comm, ierr) 2751 call xmpi_bcast(ddb_hdr%occ, master, comm, ierr) 2752 call xmpi_bcast(ddb_hdr%spinat, master, comm, ierr) 2753 call xmpi_bcast(ddb_hdr%tnons, master, comm, ierr) 2754 call xmpi_bcast(ddb_hdr%wtk, master, comm, ierr) 2755 call xmpi_bcast(ddb_hdr%xred, master, comm, ierr) 2756 call xmpi_bcast(ddb_hdr%zion, master, comm, ierr) 2757 call xmpi_bcast(ddb_hdr%znucl, master, comm, ierr) 2758 2759 call xmpi_bcast(ddb_hdr%psps%indlmn, master, comm, ierr) 2760 call xmpi_bcast(ddb_hdr%psps%pspso, master, comm, ierr) 2761 2762 if (ddb_hdr%psps%dimekb > 0 .and. ddb_hdr%psps%usepaw==0) then 2763 call xmpi_bcast(ddb_hdr%psps%ekb, master, comm, ierr) 2764 end if 2765 2766 nn = ddb_hdr%psps%ntypat * ddb_hdr%psps%usepaw 2767 if (nn.gt.0) then 2768 do ii=1,nn 2769 call pawtab_bcast(ddb_hdr%pawtab(ii), comm) 2770 end do 2771 end if 2772 2773 ! Strings 2774 !call xmpi_bcast(ddb_hdr%dscrpt, master, comm, ierr) 2775 2776 call ddb_hdr%crystal%bcast(comm) 2777 2778 end subroutine ddb_hdr_bcast
m_ddb_hdr/ddb_hdr_bcast_dim [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_bcast_dim
FUNCTION
Master broadcasts the header dimensions.
INPUTS
comm=MPI communicator.
SOURCE
2624 subroutine ddb_hdr_bcast_dim(ddb_hdr, comm) 2625 2626 !Arguments ------------------------------------ 2627 class(ddb_hdr_type),intent(inout) :: ddb_hdr 2628 integer, intent(in) :: comm 2629 2630 !Local variables ------------------------- 2631 integer, parameter :: master=0 2632 integer :: ierr 2633 2634 ! ************************************************************************ 2635 2636 if (xmpi_comm_size(comm) == 1) return 2637 2638 DBG_ENTER("COLL") 2639 2640 ! Integers 2641 call xmpi_bcast(ddb_hdr%ddb_version, master, comm, ierr) 2642 2643 call xmpi_bcast(ddb_hdr%matom, master, comm, ierr) 2644 call xmpi_bcast(ddb_hdr%mband, master, comm, ierr) 2645 call xmpi_bcast(ddb_hdr%mkpt, master, comm, ierr) 2646 call xmpi_bcast(ddb_hdr%msym, master, comm, ierr) 2647 call xmpi_bcast(ddb_hdr%mtypat, master, comm, ierr) 2648 2649 call xmpi_bcast(ddb_hdr%natom, master, comm, ierr) 2650 call xmpi_bcast(ddb_hdr%nkpt, master, comm, ierr) 2651 call xmpi_bcast(ddb_hdr%nsym, master, comm, ierr) 2652 call xmpi_bcast(ddb_hdr%ntypat, master, comm, ierr) 2653 2654 call xmpi_bcast(ddb_hdr%nsppol, master, comm, ierr) 2655 call xmpi_bcast(ddb_hdr%nspinor, master, comm, ierr) 2656 call xmpi_bcast(ddb_hdr%nspden, master, comm, ierr) 2657 2658 call xmpi_bcast(ddb_hdr%nblok, master, comm, ierr) 2659 call xmpi_bcast(ddb_hdr%mblktyp, master, comm, ierr) 2660 2661 call xmpi_bcast(ddb_hdr%mpert, master, comm, ierr) 2662 call xmpi_bcast(ddb_hdr%msize, master, comm, ierr) 2663 2664 call xmpi_bcast(ddb_hdr%with_psps, master, comm, ierr) 2665 call xmpi_bcast(ddb_hdr%with_dfpt_vars, master, comm, ierr) 2666 2667 call xmpi_bcast(ddb_hdr%usepaw, master, comm, ierr) 2668 call xmpi_bcast(ddb_hdr%psps%dimekb, master, comm, ierr) 2669 call xmpi_bcast(ddb_hdr%psps%ntypat, master, comm, ierr) 2670 call xmpi_bcast(ddb_hdr%psps%lmnmax, master, comm, ierr) 2671 call xmpi_bcast(ddb_hdr%psps%usepaw, master, comm, ierr) 2672 call xmpi_bcast(ddb_hdr%psps%useylm, master, comm, ierr) 2673 2674 DBG_EXIT("COLL") 2675 2676 end subroutine ddb_hdr_bcast_dim
m_ddb_hdr/ddb_hdr_close [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_close
FUNCTION
Close the open DDB file in either txt or nc format using the unit or ncid stored in the ddb_hdr object.
INPUTS
comm=MPI communicator.
SOURCE
2164 subroutine ddb_hdr_close(ddb_hdr, comm) 2165 !Arguments ------------------------------------ 2166 class(ddb_hdr_type),intent(inout) :: ddb_hdr 2167 integer,intent(in),optional :: comm 2168 2169 !Local variables ------------------------- 2170 !scalars 2171 integer,parameter :: master=0 2172 integer :: ncerr, ierr 2173 2174 ! ************************************************************************ 2175 2176 if (present(comm)) then 2177 if (xmpi_comm_rank(comm) /= master) return 2178 end if 2179 2180 if (ddb_hdr%has_open_file_txt) then 2181 close(ddb_hdr%unddb, iostat=ierr) 2182 ddb_hdr%has_open_file_txt = .false. 2183 end if 2184 if (ddb_hdr%has_open_file_nc) then 2185 #ifdef HAVE_NETCDF 2186 ncerr = nf90_close(ddb_hdr%ncid) 2187 ddb_hdr%has_open_file_nc = .false. 2188 #endif 2189 end if 2190 2191 end subroutine ddb_hdr_close
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
ddb_hdr2: ddb header to compare with.
OUTPUT
SOURCE
2211 subroutine ddb_hdr_compare(ddb_hdr1, ddb_hdr2) 2212 2213 !Arguments ------------------------------------ 2214 class(ddb_hdr_type),intent(inout) :: ddb_hdr1 2215 type(ddb_hdr_type),intent(inout) :: ddb_hdr2 2216 2217 !Local variables ------------------------- 2218 integer :: bantot,ii,ij,isym 2219 integer :: fullinit,fullinit2 2220 !integer :: itypat 2221 !real(dp) :: ekbcm8,ekbcmp 2222 real(dp),parameter :: tol=1.0d-6 ! Limited by the precision of text DDB 2223 character(len=500) :: msg 2224 2225 ! ************************************************************************ 2226 2227 fullinit = ddb_hdr1%with_psps * ddb_hdr1%with_dfpt_vars 2228 fullinit2 = ddb_hdr2%with_psps * ddb_hdr2%with_dfpt_vars 2229 2230 ! natom 2231 call chki8(ddb_hdr1%natom,ddb_hdr2%natom,' natom') 2232 2233 ! nkpt 2234 if(fullinit/=0)then 2235 if(ddb_hdr1%nkpt/=2*ddb_hdr2%nkpt .and. 2*ddb_hdr1%nkpt/=ddb_hdr2%nkpt)then 2236 ! GA: I've disabled this check, but I might want to put it back... 2237 !call chki8(nkpt,nkpt8,' nkpt') 2238 else 2239 write(std_out,*)' compar8 : assume that one of the DDB to be',' merged use Time-Reversal to' 2240 write(std_out,*)' decrease the number of k-points' 2241 end if 2242 end if 2243 2244 ! nband 2245 ! There can be the case of perturbation at Gamma, that 2246 ! only need half of the number of k points. 2247 if(fullinit==0 .or. ddb_hdr2%nkpt==2*ddb_hdr1%nkpt)then 2248 bantot=0 2249 do ii=1,ddb_hdr2%nkpt 2250 bantot=bantot+ddb_hdr1%nband(ii) 2251 end do 2252 else 2253 bantot=0 2254 do ii=1,ddb_hdr1%nkpt 2255 !if(ddb_hdr1%nkpt==ddb_hdr2%nkpt)then 2256 ! call chki8(ddb_hdr1%nband(ii),ddb_hdr2%nband(ii),' nband') 2257 !end if 2258 bantot=bantot+ddb_hdr1%nband(ii) 2259 end do 2260 end if 2261 2262 ! nsppol 2263 call chki8(ddb_hdr1%nsppol,ddb_hdr2%nsppol,'nsppol') 2264 2265 ! nsym 2266 !if(ddb_hdr1%nsym/=1 .and. ddb_hdr2%nsym/=1)then 2267 ! call chki8(ddb_hdr1%nsym,ddb_hdr2%nsym,' nsym') 2268 !end if 2269 2270 ! ntypat 2271 call chki8(ddb_hdr1%ntypat,ddb_hdr2%ntypat,'ntypat') 2272 2273 ! acell 2274 do ii=1,3 2275 call chkr8(ddb_hdr1%acell(ii),ddb_hdr2%acell(ii),' acell',tol) 2276 end do 2277 2278 ! amu 2279 do ii=1,ddb_hdr1%ntypat 2280 call chkr8(ddb_hdr1%amu(ii),ddb_hdr2%amu(ii),' amu',tol) 2281 end do 2282 2283 ! ecut 2284 call chkr8(ddb_hdr1%ecut,ddb_hdr2%ecut,' ecut',tol) 2285 2286 ! pawecutdg (PAW only) 2287 if (ddb_hdr1%usepaw==1) then 2288 call chkr8(ddb_hdr1%pawecutdg,ddb_hdr2%pawecutdg,' ecut',tol) 2289 end if 2290 2291 ! iscf 2292 if(fullinit/=0)then 2293 call chki8(ddb_hdr1%iscf,ddb_hdr2%iscf,' iscf') 2294 end if 2295 2296 ! ixc 2297 call chki8(ddb_hdr1%ixc,ddb_hdr2%ixc,' ixc') 2298 2299 ! kpt and 14. kptnrm 2300 if(ddb_hdr2%nkpt == ddb_hdr1%nkpt .and. fullinit/=0)then 2301 do ij=1,ddb_hdr1%nkpt 2302 do ii=1,3 2303 call chkr8(ddb_hdr1%kpt(ii,ij)/ddb_hdr1%kptnrm,ddb_hdr2%kpt(ii,ij)/ddb_hdr2%kptnrm,' kpt',tol) 2304 end do 2305 end do 2306 end if 2307 2308 ! ngfft 2309 ! Allow for different ngfft 2310 ! E.g. DFPT calculation might use a different ngfft from GS WFK. 2311 do ii=1,3 2312 if (ddb_hdr1%ngfft(ii) == ddb_hdr2%ngfft(ii)) cycle 2313 write(msg,'(3a,i10,3a,i10,a)') & 2314 'Comparing integers for variable ngfft.',ch10,& 2315 'Value from input DDB is',ddb_hdr1%ngfft(ii),' and',ch10,& 2316 'from transfer DDB is',ddb_hdr2%ngfft(ii),'.' 2317 ABI_WARNING(msg) 2318 end do 2319 2320 ! occ 2321 if(fullinit/=0 .and. ddb_hdr1%nkpt==ddb_hdr2%nkpt)then 2322 do ii=1,bantot 2323 call chkr8(ddb_hdr1%occ(ii),ddb_hdr2%occ(ii),' occ',tol) 2324 end do 2325 end if 2326 2327 ! rprim 2328 do ii=1,3 2329 do ij=1,3 2330 call chkr8(ddb_hdr1%rprim(ii,ij),ddb_hdr2%rprim(ii,ij),' rprim',tol) 2331 end do 2332 end do 2333 2334 ! dfpt_sciss 2335 if(fullinit/=0)then 2336 call chkr8(ddb_hdr1%dfpt_sciss,ddb_hdr2%dfpt_sciss,' dfpt_sciss',tol) 2337 end if 2338 2339 ! symrel and tnons 2340 if(ddb_hdr1%nsym==ddb_hdr2%nsym)then 2341 do isym=1,ddb_hdr1%nsym 2342 do ii=1,3 2343 call chkr8(ddb_hdr1%tnons(ii,isym),ddb_hdr2%tnons(ii,isym),' tnons',tol) 2344 do ij=1,3 2345 call chki8(ddb_hdr1%symrel(ii,ij,isym),ddb_hdr2%symrel(ii,ij,isym),'symrel') 2346 end do 2347 end do 2348 end do 2349 end if 2350 2351 ! typat 2352 do ii=1,ddb_hdr1%ntypat 2353 call chki8(ddb_hdr1%typat(ii),ddb_hdr2%typat(ii),' typat') 2354 end do 2355 2356 ! wtk 2357 if(ddb_hdr1%nkpt==ddb_hdr2%nkpt .and. fullinit/=0)then 2358 do ii=1,ddb_hdr1%nkpt 2359 call chkr8(ddb_hdr1%wtk(ii),ddb_hdr2%wtk(ii),' wtk',tol) 2360 end do 2361 end if 2362 2363 ! xred 2364 do ij=1,ddb_hdr1%natom 2365 do ii=1,3 2366 call chkr8(ddb_hdr1%xred(ii,ij),ddb_hdr2%xred(ii,ij),' xred',tol) 2367 end do 2368 end do 2369 ! zion 2370 do ii=1,ddb_hdr1%ntypat 2371 call chkr8(ddb_hdr1%zion(ii),ddb_hdr2%zion(ii),' zion',tol) 2372 end do 2373 2374 ! Now compare the NC pseudopotential information 2375 ! GA: I disabled this comparison because the ekb 2376 ! are not properly written to the DDB in text format. 2377 ! See comment in psddb8 2378 !if (ddb_hdr1%usepaw==0) then 2379 ! if(ddb_hdr1%psps%dimekb/=0 .and. fullinit/=0 .and. fullinit2/=0 )then 2380 ! do ii=1,dimekb 2381 ! do itypat=1,ntypat 2382 ! ekbcmp=ekb(ii,itypat) 2383 ! ekbcm8=ekb8(ii,itypat) 2384 ! !call chkr8(ekbcmp,ekbcm8,' ekb',tol) 2385 ! end do 2386 ! end do 2387 !end if 2388 2389 !Now compare several PAW dataset information 2390 ! if (ddb_hdr1%usepaw==1) then 2391 ! if (fullinit/=0 .and. fullinit2/=0) then 2392 ! do itypat=1,ddb_hdr1%ntypat 2393 ! call chki8(ddb_hdr1%pawtab(itypat)%basis_size,ddb_hdr2%pawtab(itypat)%basis_size,'bas_sz') 2394 ! call chki8(ddb_hdr1%pawtab(itypat)%lmn_size,ddb_hdr2%pawtab(itypat)%lmn_size,'lmn_sz') 2395 ! call chki8(ddb_hdr1%pawtab(itypat)%lmn2_size,ddb_hdr2%pawtab(itypat)%lmn2_size,'lmn2sz') 2396 ! call chkr8(ddb_hdr1%pawtab(itypat)%rpaw,ddb_hdr2%pawtab(itypat)%rpaw,' rpaw',tol3) 2397 ! call chkr8(ddb_hdr1%pawtab(itypat)%rshp,ddb_hdr2%pawtab(itypat)%rshp,'rshape',tol3) 2398 ! call chki8(ddb_hdr1%pawtab(itypat)%shape_type,ddb_hdr2%pawtab(itypat)%shape_type,'shp_tp') 2399 ! if (ddb_hdr1%pawtab(itypat)%lmn2_size>0) then 2400 ! do ii=1,ddb_hdr1%pawtab(itypat)%lmn2_size 2401 ! call chkr8(ddb_hdr1%pawtab(itypat)%dij0(ii),ddb_hdr2%pawtab(itypat)%dij0(ii),' dij0',tol) 2402 ! end do 2403 ! end if 2404 ! end do 2405 ! end if 2406 ! end if 2407 2408 2409 2410 ! Should also compare indlmn and pspso ... but suppose that 2411 ! the checking of ekb is enough for the psps. 2412 ! Should also compare many other variables ... 2413 2414 2415 end subroutine ddb_hdr_compare
m_ddb_hdr/ddb_hdr_copy_missing_variables [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_copy_missing_variables
FUNCTION
Copy missing variables from ddb_hdr2 to ddb_hdr1
INPUTS
OUTPUT
SOURCE
2433 subroutine ddb_hdr_copy_missing_variables(ddb_hdr1, ddb_hdr2) 2434 2435 !Arguments ------------------------------------ 2436 class(ddb_hdr_type),intent(inout) :: ddb_hdr1 2437 type(ddb_hdr_type),intent(inout) :: ddb_hdr2 2438 2439 !Local variables ------------------------- 2440 integer :: fullinit, fullinit2 2441 integer :: bantot,ii,ij,itypat 2442 2443 ! ************************************************************************ 2444 2445 fullinit = ddb_hdr1%with_psps * ddb_hdr1%with_dfpt_vars 2446 fullinit2 = ddb_hdr2%with_psps * ddb_hdr2%with_dfpt_vars 2447 2448 ! nkpt 2449 if(fullinit==0)then 2450 ddb_hdr1%nkpt = ddb_hdr2%nkpt 2451 end if 2452 2453 ! occopt 2454 ! Take here the most favorable case. 2455 if(ddb_hdr2%occopt==0) ddb_hdr1%occopt=0 2456 2457 ! nband 2458 ! There can be the case of perturbation at Gamma, that 2459 ! only need half of the number of k points. 2460 if(fullinit==0 .or. ddb_hdr2%nkpt==2*ddb_hdr1%nkpt)then 2461 bantot=0 2462 do ii=1,ddb_hdr2%nkpt 2463 ddb_hdr1%nband(ii) = ddb_hdr2%nband(ii) 2464 bantot=bantot+ddb_hdr1%nband(ii) 2465 end do 2466 else 2467 bantot=0 2468 do ii=1,ddb_hdr1%nkpt 2469 bantot=bantot+ddb_hdr1%nband(ii) 2470 end do 2471 end if 2472 2473 ! iscf 2474 if(fullinit==0)then 2475 ddb_hdr1%iscf = ddb_hdr2%iscf 2476 end if 2477 2478 ! kpt and kptnrm 2479 if(ddb_hdr2%nkpt == 2*ddb_hdr1%nkpt .or. fullinit==0) then 2480 2481 ! Copy the largest number of k points in the right place 2482 do ij=1,ddb_hdr2%nkpt 2483 do ii=1,3 2484 ddb_hdr1%kpt(ii,ij)=ddb_hdr2%kpt(ii,ij) 2485 end do 2486 end do 2487 2488 ddb_hdr1%kptnrm=ddb_hdr2%kptnrm 2489 2490 end if 2491 2492 ! occ 2493 if (fullinit==0 .or. ddb_hdr2%nkpt==2*ddb_hdr1%nkpt) then 2494 do ii=1,bantot 2495 ddb_hdr1%occ(ii)=ddb_hdr2%occ(ii) 2496 end do 2497 end if 2498 2499 ! dfpt_sciss 2500 if(fullinit==0)then 2501 ddb_hdr1%dfpt_sciss = ddb_hdr2%dfpt_sciss 2502 end if 2503 2504 ! symrel 2505 if(ddb_hdr1%nsym/=ddb_hdr2%nsym .and. ddb_hdr2%nsym/=1)then 2506 ! GA: It is suspicious to assume the array are allocated with the same size 2507 ! but this is only used in merge_ddb 2508 ddb_hdr1%nsym = ddb_hdr2%nsym 2509 ddb_hdr1%symrel(:,:,1:ddb_hdr2%nsym)=ddb_hdr2%symrel(:,:,1:ddb_hdr2%nsym) 2510 ddb_hdr1%tnons(:,1:ddb_hdr2%nsym)=ddb_hdr2%tnons(:,1:ddb_hdr2%nsym) 2511 end if 2512 2513 ! tolwfr 2514 !Take the less converged value... 2515 !ddb_hdr1%tolwfr=max(ddb_hdr1%tolwfr,ddb_hdr2%tolwfr) 2516 2517 ! wtk 2518 if(ddb_hdr2%nkpt==2*ddb_hdr1%nkpt .or. fullinit==0)then 2519 do ii=1,ddb_hdr2%nkpt 2520 ddb_hdr1%wtk(ii)=ddb_hdr2%wtk(ii) 2521 end do 2522 end if 2523 2524 ! Copy nkpt in the case of the use of the time-reversal symmetry 2525 if(2*ddb_hdr1%nkpt==ddb_hdr2%nkpt)then 2526 ddb_hdr1%nkpt=ddb_hdr2%nkpt 2527 end if 2528 2529 ! PAW dataset information 2530 if (ddb_hdr1%usepaw==1) then 2531 if (fullinit==0 .and. fullinit2/=0) then 2532 do itypat=1,ddb_hdr1%ntypat 2533 ddb_hdr1%pawtab(itypat)%basis_size = ddb_hdr2%pawtab(itypat)%basis_size 2534 ddb_hdr1%pawtab(itypat)%lmn_size = ddb_hdr2%pawtab(itypat)%lmn_size 2535 ddb_hdr1%pawtab(itypat)%rpaw = ddb_hdr2%pawtab(itypat)%rpaw 2536 ddb_hdr1%pawtab(itypat)%rshp = ddb_hdr2%pawtab(itypat)%rshp 2537 ddb_hdr1%pawtab(itypat)%shape_type = ddb_hdr2%pawtab(itypat)%shape_type 2538 if (ddb_hdr2%pawtab(itypat)%lmn2_size>0) then 2539 if (ddb_hdr1%pawtab(itypat)%lmn2_size==0) then 2540 ABI_MALLOC(ddb_hdr1%pawtab(itypat)%dij0,(ddb_hdr2%pawtab(itypat)%lmn2_size)) 2541 end if 2542 do ii=1,ddb_hdr2%pawtab(itypat)%lmn2_size 2543 ddb_hdr1%pawtab(itypat)%dij0(ii) = ddb_hdr2%pawtab(itypat)%dij0(ii) 2544 end do 2545 end if 2546 ddb_hdr1%pawtab(itypat)%lmn2_size = ddb_hdr2%pawtab(itypat)%lmn2_size 2547 end do 2548 end if 2549 end if 2550 end subroutine ddb_hdr_copy_missing_variables
m_ddb_hdr/ddb_hdr_copy_psps_info [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_compare
FUNCTION
Copy the information on pseudopotentials from one header to another.
INPUTS
ddb_hdr1: the ddb header into which info is copied ddb_hdr2: the ddb header from which info is copied
OUTPUT
SOURCE
2570 subroutine ddb_hdr_copy_psps_from(ddb_hdr1, ddb_hdr2) 2571 2572 !Arguments ------------------------------------ 2573 class(ddb_hdr_type),intent(inout) :: ddb_hdr1 2574 type(ddb_hdr_type),intent(inout) :: ddb_hdr2 2575 2576 !Local variables ------------------------- 2577 integer :: itypat 2578 2579 ! ************************************************************************ 2580 2581 if (ddb_hdr2%with_psps > 0) then 2582 ddb_hdr1%with_psps = ddb_hdr2%with_psps 2583 2584 call psps_free(ddb_hdr1%psps) 2585 call psps_copy(ddb_hdr2%psps, ddb_hdr1%psps) 2586 2587 end if 2588 2589 ! TODO: Need a function pawtab_copy 2590 if (ddb_hdr1%psps%usepaw > 0) then 2591 do itypat=1,ddb_hdr1%mtypat 2592 ddb_hdr1%pawtab(itypat)%basis_size = ddb_hdr2%pawtab(itypat)%basis_size 2593 ddb_hdr1%pawtab(itypat)%lmn_size = ddb_hdr2%pawtab(itypat)%lmn_size 2594 ddb_hdr1%pawtab(itypat)%lmn2_size = ddb_hdr2%pawtab(itypat)%lmn2_size 2595 ddb_hdr1%pawtab(itypat)%shape_type = ddb_hdr2%pawtab(itypat)%shape_type 2596 ddb_hdr1%pawtab(itypat)%rpaw = ddb_hdr2%pawtab(itypat)%rpaw 2597 ddb_hdr1%pawtab(itypat)%rshp = ddb_hdr2%pawtab(itypat)%rshp 2598 if (allocated(ddb_hdr2%pawtab(itypat)%dij0)) then 2599 ABI_SFREE(ddb_hdr1%pawtab(itypat)%dij0) 2600 ABI_MALLOC(ddb_hdr1%pawtab(itypat)%dij0, (ddb_hdr1%pawtab(itypat)%lmn2_size)) 2601 ddb_hdr1%pawtab(itypat)%dij0(:) = ddb_hdr2%pawtab(itypat)%dij0(:) 2602 end if 2603 end do 2604 end if 2605 2606 2607 end subroutine ddb_hdr_copy_psps_from
m_ddb_hdr/ddb_hdr_free [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_free
FUNCTION
Free memory.
SOURCE
634 subroutine ddb_hdr_free(ddb_hdr) 635 636 !Arguments ------------------------------------ 637 class(ddb_hdr_type),intent(inout) :: ddb_hdr 638 639 ! ************************************************************************ 640 641 ! integer 642 ABI_SFREE(ddb_hdr%nband) 643 ABI_SFREE(ddb_hdr%symafm) 644 ABI_SFREE(ddb_hdr%symrel) 645 ABI_SFREE(ddb_hdr%typat) 646 647 ! real 648 ABI_SFREE(ddb_hdr%amu) 649 ABI_SFREE(ddb_hdr%kpt) 650 ABI_SFREE(ddb_hdr%occ) 651 ABI_SFREE(ddb_hdr%spinat) 652 ABI_SFREE(ddb_hdr%tnons) 653 ABI_SFREE(ddb_hdr%wtk) 654 ABI_SFREE(ddb_hdr%xred) 655 ABI_SFREE(ddb_hdr%zion) 656 ABI_SFREE(ddb_hdr%znucl) 657 ABI_SFREE(ddb_hdr%typ) 658 659 ! types 660 call ddb_hdr%crystal%free() 661 call psps_free(ddb_hdr%psps) 662 663 if (allocated(ddb_hdr%pawtab)) then 664 call pawtab_free(ddb_hdr%pawtab) 665 ABI_FREE(ddb_hdr%pawtab) 666 end if 667 668 end subroutine ddb_hdr_free
m_ddb_hdr/ddb_hdr_get_block_dims [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_get_block_dims
FUNCTION
Get the maximum number of perturbations from natom, and get the maximum matrix size of one block, knowing the block types.
SOURCE
602 subroutine ddb_hdr_get_block_dims(ddb_hdr) 603 604 !Arguments ------------------------------------ 605 class(ddb_hdr_type),intent(inout) :: ddb_hdr 606 607 ! ************************************************************************ 608 609 ! Compute mpert 610 !ddb_hdr%mpert = ddb_hdr%natom+MPERT_MAX 611 ! GA: mpert is stored in netcdf format but not in text format. 612 613 ! Compute msize 614 if (is_type_d3E(ddb_hdr%mblktyp)) then 615 ddb_hdr%msize=3*ddb_hdr%mpert*3*ddb_hdr%mpert*3*ddb_hdr%mpert 616 else 617 ddb_hdr%msize=3*ddb_hdr%mpert*3*ddb_hdr%mpert 618 end if 619 620 end subroutine ddb_hdr_get_block_dims
m_ddb_hdr/ddb_hdr_get_iomode [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_get_iomode
FUNCTION
Decide whether to use text or netcdf IO
INPUTS
filename=name of the file to examine. io=specify whether we want to read (io=1) or write (io=2)
OUTPUT
iomode=The io mode. filenameout=The filename, possibly appended by '.nc'
SOURCE
2800 subroutine ddb_hdr_get_iomode(ddb_hdr, filenamein, io, iomode, filenameout) 2801 !Arguments ------------------------------------ 2802 class(ddb_hdr_type),intent(inout) :: ddb_hdr 2803 character(len=fnlen),intent(in) :: filenamein 2804 integer,intent(in) :: io 2805 integer,intent(out) :: iomode 2806 character(len=fnlen),intent(out) :: filenameout 2807 2808 !Local variables ------------------------------ 2809 character(len=fnlen) :: filenamenc 2810 2811 ! ************************************************************************ 2812 2813 ! Simple case: file name was given with netcdf extension 2814 #ifdef HAVE_NETCDF 2815 if (endswith(filenamein, ".nc")) then 2816 filenameout = filenamein 2817 iomode = IO_MODE_ETSF 2818 return 2819 else 2820 ! In most case, this function is called without filename extension. 2821 filenamenc = nctk_ncify(filenamein) 2822 end if 2823 #else 2824 if (endswith(filenamein, ".nc")) then 2825 ABI_ERROR("NETCDF support required to read or write DDB.nc file.") 2826 else 2827 filenameout = filenamein 2828 iomode = IO_MODE_FORTRAN 2829 return 2830 end if 2831 #endif 2832 2833 ! Reading 2834 if (io==1) then 2835 2836 ! First check for the existence of netcdf file 2837 if (file_exists(filenamenc)) then 2838 filenameout = filenamenc 2839 iomode = IO_MODE_ETSF 2840 return 2841 2842 ! Then check for the existence of text file 2843 else if (file_exists(filenamein)) then 2844 filenameout = filenamein 2845 iomode = IO_MODE_FORTRAN 2846 return 2847 2848 ! File is not available for reading. 2849 else 2850 2851 ! This doesnt always have to be an error 2852 !ABI_ERROR(sjoin("Cannot find DDB file:", filenamein)) 2853 ABI_WARNING(sjoin("Cannot find DDB file:", filenamein)) 2854 2855 ! Rely on default value of iomode. 2856 iomode = ddb_hdr%iomode 2857 if (iomode==IO_MODE_ETSF) then 2858 filenameout = filenamenc 2859 else 2860 filenameout = filenamein 2861 end if 2862 return 2863 end if 2864 2865 ! Writing 2866 else if (io==2) then 2867 2868 ! Rely on ddb_hdr%iomode, which should have been initialized from dtset 2869 if (ddb_hdr%iomode==IO_MODE_ETSF) then 2870 filenameout = filenamenc 2871 iomode = IO_MODE_ETSF 2872 return 2873 else if ((ddb_hdr%iomode==IO_MODE_FORTRAN).or.(ddb_hdr%iomode==IO_MODE_FORTRAN_MASTER).or.(ddb_hdr%iomode==IO_MODE_MPI)) then 2874 filenameout = filenamein 2875 iomode = IO_MODE_FORTRAN 2876 return 2877 else 2878 ABI_ERROR("Unexpected value for iomode.") 2879 end if 2880 else 2881 ABI_ERROR("Unexpected value for io.") 2882 end if 2883 2884 end subroutine ddb_hdr_get_iomode
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
ddb_hdr=the new ddb_hdr object dtset=dtset object of the current calculation psps=the pseudopotential objects pawtab=pawtab object dscrpt=string description of the ddb nblok=number of blocks xred=reduced coordinate positions of the atoms occ=occupation number of the bands ngfft=fft grid mband=maximum number of bands mkpt=maximum number of kpoints
OUTPUT
SOURCE
328 subroutine ddb_hdr_init(ddb_hdr, dtset, psps, pawtab, dscrpt, & 329 nblok, msize, mpert, xred, occ, ngfft,& 330 mband,nkpt,kpt) 331 332 !Arguments ------------------------------------ 333 class(ddb_hdr_type),intent(inout) :: ddb_hdr 334 type(dataset_type),intent(in) :: dtset 335 type(pseudopotential_type),intent(inout) :: psps 336 type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw) 337 character(len=*),intent(in) :: dscrpt ! TODO: Make this one optional. 338 integer,intent(in) :: nblok 339 !integer,intent(in),optional :: mpert ! GA: TODO 340 integer,intent(in),optional :: mband,nkpt,msize,mpert 341 integer,intent(in),optional :: ngfft(18) 342 real(dp),intent(in),optional :: xred(3,dtset%natom) 343 real(dp),intent(in),optional :: occ(dtset%mband*dtset%nkpt*dtset%nsppol) 344 real(dp),intent(in),optional :: kpt(:,:) 345 346 !Local variables ------------------------- 347 integer :: ii, nn, ikpt 348 349 ! ************************************************************************ 350 351 ddb_hdr%nblok = nblok 352 ddb_hdr%ddb_version = DDB_VERSION 353 ddb_hdr%dscrpt = trim(dscrpt) 354 355 call psps_copy(psps, ddb_hdr%psps) 356 357 ! Copy scalars from dtset 358 ddb_hdr%matom = dtset%natom 359 ddb_hdr%natom = dtset%natom 360 ddb_hdr%mband = dtset%mband 361 ddb_hdr%mkpt = dtset%nkpt 362 ddb_hdr%nkpt = dtset%nkpt 363 ddb_hdr%msym = dtset%nsym 364 ddb_hdr%nsym = dtset%nsym 365 ddb_hdr%mtypat = dtset%ntypat 366 ddb_hdr%ntypat = dtset%ntypat 367 ddb_hdr%nspden = dtset%nspden 368 ddb_hdr%nspinor = dtset%nspinor 369 ddb_hdr%nsppol = dtset%nsppol 370 371 if (present(mband)) then 372 ddb_hdr%mband = mband 373 end if 374 375 if (present(nkpt)) then 376 ddb_hdr%mkpt = nkpt 377 ddb_hdr%nkpt = nkpt 378 end if 379 380 if (present(mpert)) then 381 ddb_hdr%mpert = mpert 382 end if 383 384 if (present(msize)) then 385 ddb_hdr%msize = msize 386 else 387 if (present(mpert)) then 388 call ddb_hdr%get_block_dims() 389 end if 390 end if 391 392 ddb_hdr%occopt = dtset%occopt 393 ddb_hdr%usepaw = dtset%usepaw 394 395 ddb_hdr%intxc = dtset%intxc 396 ddb_hdr%ixc = dtset%ixc 397 ddb_hdr%iscf = dtset%iscf 398 399 ddb_hdr%dilatmx = dtset%dilatmx 400 ddb_hdr%ecut = dtset%ecut 401 ddb_hdr%ecutsm = dtset%ecutsm 402 ddb_hdr%pawecutdg = dtset%pawecutdg 403 ddb_hdr%kptnrm = dtset%kptnrm 404 ddb_hdr%dfpt_sciss = dtset%dfpt_sciss 405 ddb_hdr%tolwfr = one 406 ddb_hdr%tphysel = dtset%tphysel 407 ddb_hdr%tsmear = dtset%tsmear 408 409 ddb_hdr%with_psps = 1 410 ddb_hdr%with_dfpt_vars = 1 411 412 ddb_hdr%iomode = dtset%iomode 413 ddb_hdr%has_open_file_txt = .false. 414 ddb_hdr%has_open_file_nc = .false. 415 416 call ddb_hdr%malloc() 417 418 if (present(ngfft)) then 419 ddb_hdr%ngfft = ngfft 420 else 421 ddb_hdr%ngfft = dtset%ngfft 422 end if 423 424 ! Copy arrays from dtset 425 ddb_hdr%acell(:) = dtset%acell_orig(1:3,1) 426 ddb_hdr%rprim(:,:) = dtset%rprim_orig(1:3,1:3,1) 427 ddb_hdr%amu(:) = dtset%amu_orig(:,1) 428 ddb_hdr%nband(:) = dtset%nband(1:ddb_hdr%mkpt*ddb_hdr%nsppol) 429 ddb_hdr%symafm(:) = dtset%symafm(1:ddb_hdr%msym) 430 ddb_hdr%symrel(:,:,:) = dtset%symrel(1:3,1:3,1:ddb_hdr%msym) 431 ddb_hdr%typat(:) = dtset%typat(1:ddb_hdr%matom) 432 ddb_hdr%wtk(:) = dtset%wtk(1:ddb_hdr%mkpt) 433 ddb_hdr%spinat(:,:) = dtset%spinat(1:3,1:ddb_hdr%matom) 434 ddb_hdr%tnons(:,:) = dtset%tnons(1:3,1:ddb_hdr%msym) 435 ddb_hdr%zion(:) = dtset%ziontypat(1:ddb_hdr%mtypat) 436 ddb_hdr%znucl(:) = dtset%znucl(1:ddb_hdr%mtypat) 437 438 439 if (present(kpt)) then 440 do ikpt=1,ddb_hdr%nkpt 441 do ii = 1,3 442 ddb_hdr%kpt(ii,ikpt) = kpt(ii,ikpt) 443 end do 444 end do 445 else 446 ddb_hdr%kpt(:,:) = dtset%kpt(1:3,1:ddb_hdr%mkpt) 447 end if 448 449 if (present(mband)) then 450 do ii=1,ddb_hdr%nkpt 451 ddb_hdr%nband(ii) = mband 452 end do 453 end if 454 455 if (present(xred)) then 456 ddb_hdr%xred(:,:) = xred(1:3,1:ddb_hdr%matom) 457 else 458 ddb_hdr%xred(:,:) = dtset%xred_orig(1:3,1:ddb_hdr%matom,1) 459 end if 460 if (present(occ)) then 461 ddb_hdr%occ(:) = occ(1:ddb_hdr%mband*ddb_hdr%mkpt*ddb_hdr%nsppol) 462 else 463 ddb_hdr%occ(:) = dtset%occ_orig(1:ddb_hdr%mband*ddb_hdr%mkpt*ddb_hdr%nsppol,1) 464 end if 465 466 ! GA: I had way too many problems implementing pawtab_copy. 467 ! The script check-libpaw would report all sorts of errors. 468 ! Therefore, I do a cheap copy here, copying only the relevant info. 469 !call pawtab_copy(pawtab, ddb_hdr%pawtab) 470 nn=size(pawtab) 471 ABI_MALLOC(ddb_hdr%pawtab,(ddb_hdr%psps%ntypat*ddb_hdr%psps%usepaw)) 472 if (nn.gt.0) then 473 474 call pawtab_nullify(ddb_hdr%pawtab) 475 476 do ii=1,nn 477 ddb_hdr%pawtab(ii)%basis_size = pawtab(ii)%basis_size 478 ddb_hdr%pawtab(ii)%lmn_size = pawtab(ii)%lmn_size 479 ddb_hdr%pawtab(ii)%lmn2_size = pawtab(ii)%lmn2_size 480 ddb_hdr%pawtab(ii)%rpaw = pawtab(ii)%rpaw 481 ddb_hdr%pawtab(ii)%rshp = pawtab(ii)%rshp 482 ddb_hdr%pawtab(ii)%shape_type = pawtab(ii)%shape_type 483 if (allocated(pawtab(ii)%dij0)) then 484 call alloc_copy(pawtab(ii)%dij0, ddb_hdr%pawtab(ii)%dij0) 485 end if 486 end do 487 end if 488 489 call crystal_init(dtset%amu_orig(:,1), ddb_hdr%crystal, & 490 & dtset%spgroup, dtset%natom, dtset%npsp, psps%ntypat, & 491 & dtset%nsym, dtset%rprimd_orig(:,:,1), dtset%typat, & 492 & ddb_hdr%xred, dtset%ziontypat, dtset%znucl, 1, & 493 & dtset%nspden==2.and.dtset%nsppol==1, .false., '', & 494 & dtset%symrel, dtset%tnons, dtset%symafm) 495 496 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.
SOURCE
510 subroutine ddb_hdr_malloc(ddb_hdr) 511 512 !Arguments ------------------------------------ 513 class(ddb_hdr_type),intent(inout) :: ddb_hdr 514 515 ! ************************************************************************ 516 517 ! integer 518 ABI_MALLOC(ddb_hdr%nband,(ddb_hdr%mkpt*ddb_hdr%nsppol)) 519 ABI_MALLOC(ddb_hdr%symafm,(ddb_hdr%msym)) 520 ABI_MALLOC(ddb_hdr%symrel,(3,3,ddb_hdr%msym)) 521 ABI_MALLOC(ddb_hdr%typat,(ddb_hdr%matom)) 522 523 ! real 524 ABI_MALLOC(ddb_hdr%amu,(ddb_hdr%mtypat)) 525 ABI_MALLOC(ddb_hdr%kpt,(3, ddb_hdr%mkpt)) 526 ABI_MALLOC(ddb_hdr%occ,(ddb_hdr%mband*ddb_hdr%mkpt*ddb_hdr%nsppol)) 527 ABI_MALLOC(ddb_hdr%spinat,(3, ddb_hdr%matom)) 528 ABI_MALLOC(ddb_hdr%tnons,(3, ddb_hdr%msym)) 529 ABI_MALLOC(ddb_hdr%wtk,(ddb_hdr%mkpt)) 530 ABI_MALLOC(ddb_hdr%xred,(3, ddb_hdr%matom)) 531 ABI_MALLOC(ddb_hdr%zion,(ddb_hdr%mtypat)) 532 ABI_MALLOC(ddb_hdr%znucl,(ddb_hdr%mtypat)) 533 534 ! types 535 !ABI_MALLOC(ddb_hdr%pawtab,(ddb_hdr%psps%ntypat*ddb_hdr%psps%usepaw)) 536 !call pawtab_nullify(ddb_hdr%pawtab) 537 538 ! Do not allocate the pseudo 539 540 ddb_hdr%ngfft = zero 541 ddb_hdr%nband = zero 542 ddb_hdr%symafm = one 543 ddb_hdr%symrel = zero 544 ddb_hdr%typat = zero 545 ddb_hdr%amu = zero 546 ddb_hdr%kpt = zero 547 ddb_hdr%occ = zero 548 ddb_hdr%spinat = zero 549 ddb_hdr%tnons = zero 550 ddb_hdr%wtk = zero 551 ddb_hdr%xred = zero 552 ddb_hdr%zion = zero 553 ddb_hdr%znucl = zero 554 555 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
filename=name of the file being written (abo_DS*_DDB) comm=MPI communicator matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw=dimensions to be set. dimonly 1-> only read dimensions and close the file 0-> leave the file open
SOURCE
1516 subroutine ddb_hdr_open_read(ddb_hdr, filename, comm, & 1517 & matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw,dimonly) 1518 1519 !Arguments ------------------------------------ 1520 class(ddb_hdr_type),intent(inout) :: ddb_hdr 1521 character(len=*),intent(in) :: filename 1522 integer,intent(in) :: comm 1523 integer,intent(in),optional :: matom,mtypat,mband,mkpt,msym 1524 integer,intent(in),optional :: dimekb,lmnmax,usepaw 1525 integer,intent(in),optional :: dimonly 1526 1527 !Local variables------------------------------- 1528 !scalars 1529 integer,parameter :: master=0 1530 integer :: iomode 1531 character(len=fnlen) :: filename_ 1532 1533 1534 call ddb_hdr%get_iomode(filename, 1, iomode, filename_) 1535 1536 if (xmpi_comm_rank(comm) == master) then 1537 1538 ! GA: FIXME 1539 ! Here we keep the old output behavior just to make test pass. 1540 ! However, test should allow differences in filenames. 1541 !call wrtout(std_out, sjoin(" Opening DDB file:", filename_), 'COLL') 1542 call wrtout(std_out, sjoin(" Opening DDB file:", filename), 'COLL') 1543 end if 1544 1545 if (iomode==IO_MODE_ETSF) then 1546 1547 call ddb_hdr%open_read_nc(filename_, comm, & 1548 & matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw,dimonly) 1549 1550 else if (iomode==IO_MODE_FORTRAN) then 1551 1552 call ddb_hdr%open_read_txt(filename, comm, & 1553 & matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw,dimonly) 1554 1555 end if 1556 1557 1558 end subroutine ddb_hdr_open_read
m_ddb_hdr/ddb_hdr_open_read_nc [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_open_read
FUNCTION
Open the DDB.nc file and read the header.
INPUTS
filename=name of the file being written (abo_DS*_DDB.nc) comm=MPI communicator matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw=dimensions to be set. dimonly 1-> only read dimensions and close the file 0-> leave the file open
SOURCE
1833 subroutine ddb_hdr_open_read_nc(ddb_hdr, filename, comm, & 1834 & matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw,dimonly) 1835 1836 !Arguments ------------------------------------ 1837 class(ddb_hdr_type),intent(inout) :: ddb_hdr 1838 character(len=*),intent(in) :: filename 1839 integer,intent(in),optional :: comm 1840 integer,intent(in),optional :: matom,mtypat,mband,mkpt,msym 1841 integer,intent(in),optional :: dimekb,lmnmax,usepaw 1842 integer,intent(in),optional :: dimonly 1843 1844 !Local variables------------------------------- 1845 !scalars 1846 integer,parameter :: master=0 1847 integer :: ncid 1848 integer :: ncid_crystal, ncid_psps, ncid_pawtab 1849 integer :: ii,isppol,isym,ikpt,iband,itypat,iblok 1850 integer :: comm_ 1851 integer :: natom_file, ntypat_file, mband_file, nkpt_file, nsym_file 1852 integer :: mband_, nkpt_ 1853 integer :: ncerr 1854 !arrays 1855 integer :: ngfft(3) 1856 integer,allocatable :: nband(:,:) 1857 integer,allocatable :: pawtab_basis_size(:) 1858 integer,allocatable :: pawtab_lmn_size(:) 1859 integer,allocatable :: pawtab_lmn2_size(:) 1860 integer,allocatable :: pawtab_shape_type(:) 1861 real(dp),allocatable :: pawtab_rpaw(:) 1862 real(dp),allocatable :: pawtab_rshp(:) 1863 real(dp),allocatable :: occ(:,:,:) 1864 1865 ! ************************************************************************ 1866 1867 if (present(comm)) then 1868 comm_ = comm 1869 else 1870 comm_ = xmpi_comm_self 1871 end if 1872 1873 if (xmpi_comm_rank(comm_) == master) then 1874 NCF_CHECK(nctk_open_read(ncid, filename, xmpi_comm_self)) 1875 ddb_hdr%ncid = ncid 1876 end if 1877 1878 if (xmpi_comm_size(comm_) > 1) then 1879 call xmpi_bcast(ddb_hdr%ncid, master, comm_, ncerr) 1880 end if 1881 1882 ddb_hdr%has_open_file_nc = .True. 1883 1884 ! Flush allocated arrays 1885 call ddb_hdr%free() 1886 1887 ! =============== 1888 ! Read dimensions 1889 ! =============== 1890 if (xmpi_comm_rank(comm_) == master) then 1891 1892 ! No check is performed on the ddb version number, yet. 1893 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'ddb_version'), ddb_hdr%ddb_version)) 1894 1895 ! Scalar integers 1896 NCF_CHECK(nctk_get_dim(ncid, "number_of_atoms", natom_file)) 1897 NCF_CHECK(nctk_get_dim(ncid, "number_of_kpoints", nkpt_file)) 1898 ddb_hdr%natom = natom_file 1899 ddb_hdr%nkpt = nkpt_file 1900 1901 NCF_CHECK(nctk_get_dim(ncid, "number_of_spins", ddb_hdr%nsppol)) 1902 NCF_CHECK(nctk_get_dim(ncid, "number_of_spinor_components", ddb_hdr%nspinor)) 1903 NCF_CHECK(nctk_get_dim(ncid, "number_of_spin_densities", ddb_hdr%nspden)) 1904 1905 NCF_CHECK(nctk_get_dim(ncid, "number_of_symmetries", nsym_file)) 1906 NCF_CHECK(nctk_get_dim(ncid, "number_of_atom_types", ntypat_file)) 1907 ddb_hdr%nsym = nsym_file 1908 ddb_hdr%ntypat = ntypat_file 1909 1910 NCF_CHECK(nctk_get_dim(ncid, "maximum_number_of_bands", mband_file)) 1911 ddb_hdr%mband = mband_file 1912 1913 NCF_CHECK(nctk_get_dim(ncid, "number_of_blocks", ddb_hdr%nblok)) 1914 NCF_CHECK(nctk_get_dim(ncid, "number_of_perturbations", ddb_hdr%mpert)) 1915 1916 ddb_hdr%msym = ddb_hdr%nsym 1917 ddb_hdr%mkpt = ddb_hdr%nkpt 1918 ddb_hdr%matom = ddb_hdr%natom 1919 ddb_hdr%mtypat = ddb_hdr%ntypat 1920 1921 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "with_psps"), ddb_hdr%with_psps)) 1922 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "with_dfpt_vars"), ddb_hdr%with_dfpt_vars)) 1923 1924 ! Read pseudopotentials info 1925 ncid_psps = nctk_idgroup(ncid, 'Pseudopotentials') 1926 1927 call psps_ncread(ddb_hdr%psps, ncid_psps) 1928 ddb_hdr%usepaw = ddb_hdr%psps%usepaw 1929 1930 ABI_MALLOC(ddb_hdr%pawtab, (ddb_hdr%ntypat*ddb_hdr%usepaw)) 1931 if (ddb_hdr%usepaw > 0) then 1932 1933 ! Read PAW info 1934 ! ------------------------------------- 1935 ! GA: TODO There should be a function pawtab_ncread 1936 ! Here, we only read the info needed to print the ddb in text format. 1937 1938 ncid_pawtab = nctk_idgroup(ncid_psps, 'PAW') 1939 call pawtab_nullify(ddb_hdr%pawtab) 1940 1941 ABI_MALLOC(pawtab_basis_size, (ddb_hdr%ntypat)) 1942 ABI_MALLOC(pawtab_lmn_size, (ddb_hdr%ntypat)) 1943 ABI_MALLOC(pawtab_lmn2_size, (ddb_hdr%ntypat)) 1944 ABI_MALLOC(pawtab_shape_type, (ddb_hdr%ntypat)) 1945 1946 NCF_CHECK(nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'basis_size'), pawtab_basis_size)) 1947 NCF_CHECK(nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'lmn_size'), pawtab_lmn_size)) 1948 NCF_CHECK(nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'lmn2_size'), pawtab_lmn2_size)) 1949 NCF_CHECK(nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'shape_type'), pawtab_shape_type)) 1950 1951 ABI_MALLOC(pawtab_rpaw, (ddb_hdr%ntypat)) 1952 ncerr = nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'r_PAW'), pawtab_rpaw) 1953 NCF_CHECK(ncerr) 1954 ABI_MALLOC(pawtab_rshp, (ddb_hdr%ntypat)) 1955 ncerr = nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'r_shape'), pawtab_rshp) 1956 NCF_CHECK(ncerr) 1957 1958 do itypat=1,ddb_hdr%ntypat 1959 1960 ddb_hdr%pawtab(itypat)%basis_size = pawtab_basis_size(itypat) 1961 ddb_hdr%pawtab(itypat)%lmn_size = pawtab_lmn_size(itypat) 1962 ddb_hdr%pawtab(itypat)%lmn2_size = pawtab_lmn2_size(itypat) 1963 ddb_hdr%pawtab(itypat)%shape_type = pawtab_shape_type(itypat) 1964 1965 ddb_hdr%pawtab(itypat)%rpaw = pawtab_rpaw(itypat) 1966 ddb_hdr%pawtab(itypat)%rshp = pawtab_rshp(itypat) 1967 1968 ABI_MALLOC(ddb_hdr%pawtab(itypat)%dij0, (ddb_hdr%pawtab(itypat)%lmn2_size)) 1969 ncerr = nf90_get_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'Dij0'), & 1970 ddb_hdr%pawtab(itypat)%dij0, & 1971 start=[itypat,1], & 1972 count=[1,ddb_hdr%pawtab(itypat)%lmn2_size]) 1973 NCF_CHECK(ncerr) 1974 1975 end do 1976 1977 ABI_FREE(pawtab_rpaw) 1978 ABI_FREE(pawtab_rshp) 1979 ABI_FREE(pawtab_basis_size) 1980 ABI_FREE(pawtab_lmn_size) 1981 ABI_FREE(pawtab_lmn2_size) 1982 ABI_FREE(pawtab_shape_type) 1983 1984 ! ------------------------------------- 1985 1986 end if 1987 1988 ! Read block types and compute maximum block type 1989 ABI_SFREE(ddb_hdr%typ) 1990 ABI_MALLOC(ddb_hdr%typ,(ddb_hdr%nblok)) 1991 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'block_types'), ddb_hdr%typ)) 1992 1993 ! TODO: mblktyp should be removed 1994 ddb_hdr%mblktyp=0 1995 do iblok=1,ddb_hdr%nblok 1996 ddb_hdr%mblktyp = max(ddb_hdr%mblktyp, ddb_hdr%typ(iblok)) 1997 end do 1998 1999 ! Overwrite dimensions if requested 2000 ! GA: TODO Try to remove these optional arguments 2001 ! Should at least give a warning when inconsistencies are found 2002 if (present(matom)) ddb_hdr%matom = matom 2003 if (present(mtypat)) ddb_hdr%mtypat = mtypat 2004 if (present(mband)) ddb_hdr%mband = mband 2005 if (present(mkpt)) ddb_hdr%mkpt = mkpt 2006 if (present(msym)) ddb_hdr%msym = msym 2007 if (present(dimekb)) ddb_hdr%psps%dimekb = dimekb 2008 if (present(lmnmax)) ddb_hdr%psps%lmnmax = lmnmax 2009 if (present(usepaw)) ddb_hdr%usepaw = usepaw 2010 2011 ! Compute the block dimensions (mpert and msize) 2012 call ddb_hdr%get_block_dims() 2013 2014 end if 2015 2016 ddb_hdr%dscrpt = '' 2017 2018 if (present(dimonly)) then 2019 if (dimonly==1) then 2020 call ddb_hdr%bcast_dim(comm_) 2021 call ddb_hdr%close() 2022 return 2023 end if 2024 end if 2025 2026 ! =============== 2027 ! Allocate arrays 2028 ! =============== 2029 if (xmpi_comm_rank(comm_) == master) then 2030 2031 call ddb_hdr%malloc() 2032 2033 !! GA : Not sure we still need this. 2034 !ddb_hdr%symafm(:)=1 2035 !if(ddb_hdr%mtypat>=1)then 2036 ! ddb_hdr%znucl(:)=zero 2037 !end if 2038 !if(ddb_hdr%matom>=1)then 2039 ! ddb_hdr%spinat(:,:)=zero 2040 !end if 2041 2042 end if 2043 2044 ! ============== 2045 ! Read variables 2046 ! ============== 2047 if (xmpi_comm_rank(comm_) == master) then 2048 2049 ! ------- 2050 ! Scalars 2051 ! ------- 2052 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'intxc'), ddb_hdr%intxc)) 2053 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'iscf'), ddb_hdr%iscf)) 2054 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'ixc'), ddb_hdr%ixc)) 2055 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'occopt'), ddb_hdr%occopt)) 2056 2057 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'dilatmx'), ddb_hdr%dilatmx)) 2058 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'ecut'), ddb_hdr%ecut)) 2059 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'pawecutdg'), ddb_hdr%pawecutdg)) 2060 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'ecutsm'), ddb_hdr%ecutsm)) 2061 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'kptnrm'), ddb_hdr%kptnrm)) 2062 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'dfpt_sciss'), ddb_hdr%dfpt_sciss)) 2063 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'tolwfr'), ddb_hdr%tolwfr)) 2064 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'tphysel'), ddb_hdr%tphysel)) 2065 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'tsmear'), ddb_hdr%tsmear)) 2066 2067 ! ------ 2068 ! Arrays 2069 ! ------ 2070 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'acell'), ddb_hdr%acell)) 2071 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'ngfft'), ngfft)) 2072 ddb_hdr%ngfft = zero 2073 do ii=1,3 2074 ddb_hdr%ngfft(ii) = ngfft(ii) 2075 end do 2076 2077 ! Basis 2078 ! Note that internal value of mband might differ from file when merging ddbs. 2079 ABI_MALLOC(occ, (mband_file,nkpt_file,ddb_hdr%nsppol)) 2080 ABI_MALLOC(nband, (nkpt_file,ddb_hdr%nsppol)) 2081 2082 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'spinat'), ddb_hdr%spinat)) 2083 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'occupations'), occ)) 2084 ncerr = nf90_get_var(ncid, nctk_idname(ncid, 'reduced_coordinates_of_kpoints'), ddb_hdr%kpt, count=[3,nkpt_file]) 2085 NCF_CHECK(ncerr) 2086 ncerr = nf90_get_var(ncid, nctk_idname(ncid, 'kpoints_weights'), ddb_hdr%wtk, count=[nkpt_file]) 2087 NCF_CHECK(ncerr) 2088 NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, 'number_of_bands'), nband)) 2089 2090 2091 ! Flatten some arrays for internal storage 2092 ddb_hdr%occ = zero 2093 mband_ = min(mband_file, ddb_hdr%mband) 2094 nkpt_ = min(nkpt_file, ddb_hdr%nkpt) 2095 ii = 0 2096 do isppol=1,ddb_hdr%nsppol 2097 do ikpt=1,nkpt_ 2098 do iband=1,mband_ 2099 ii = ii + 1 2100 ddb_hdr%occ(ii) = occ(iband,ikpt,isppol) 2101 end do 2102 end do 2103 end do 2104 ddb_hdr%nband = zero 2105 ii = 0 2106 do isppol=1,ddb_hdr%nsppol 2107 do ikpt=1,nkpt_ 2108 ii = ii + 1 2109 ddb_hdr%nband(ii) = nband(ikpt,isppol) 2110 end do 2111 end do 2112 ABI_FREE(occ) 2113 ABI_FREE(nband) 2114 2115 ! ---------------------- 2116 ! Initialize the crystal 2117 ! ---------------------- 2118 ncid_crystal = nctk_idgroup(ncid, 'Crystal') 2119 call ddb_hdr%crystal%ncread(ncid_crystal) 2120 2121 ! Copy symmetries 2122 do isym=1,ddb_hdr%nsym 2123 ddb_hdr%symafm(isym) = ddb_hdr%crystal%symafm(isym) 2124 ddb_hdr%tnons(:,isym) = ddb_hdr%crystal%tnons(:,isym) 2125 ddb_hdr%symrel(:,:,isym) = ddb_hdr%crystal%symrel(:,:,isym) 2126 end do 2127 2128 ! Copy geometry 2129 do ii=1,3 2130 ddb_hdr%rprim(:,ii) = ddb_hdr%crystal%rprimd(:,ii) / ddb_hdr%acell(ii) 2131 end do 2132 2133 ddb_hdr%xred = ddb_hdr%crystal%xred 2134 ddb_hdr%typat = ddb_hdr%crystal%typat 2135 ddb_hdr%znucl = ddb_hdr%crystal%znucl 2136 ddb_hdr%zion = ddb_hdr%crystal%zion 2137 ddb_hdr%amu = ddb_hdr%crystal%amu 2138 2139 end if 2140 2141 ! ====================== 2142 ! Master broadcasts data 2143 ! ====================== 2144 call ddb_hdr%bcast(comm_) 2145 2146 end subroutine ddb_hdr_open_read_nc
m_ddb_hdr/ddb_hdr_open_read_txt [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_open_read_txt
FUNCTION
Open the DDB text file and read the header. If dimonly, only read dimensions and close the file, otherwise leave the file open.
INPUTS
filename=name of the file being written (abo_DS*_DDB) comm=MPI communicator matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw=dimensions to be set. dimonly 1-> only read dimensions and close the file 0-> leave the file open
SOURCE
1582 subroutine ddb_hdr_open_read_txt(ddb_hdr, filename, comm, & 1583 & matom,mtypat,mband,mkpt,msym,dimekb,lmnmax,usepaw,dimonly) 1584 1585 !Arguments ------------------------------------ 1586 class(ddb_hdr_type),intent(inout) :: ddb_hdr 1587 character(len=*),intent(in) :: filename 1588 integer,intent(in),optional :: comm 1589 integer,intent(in),optional :: matom,mtypat,mband,mkpt,msym 1590 integer,intent(in),optional :: dimekb,lmnmax,usepaw 1591 integer,intent(in),optional :: dimonly 1592 1593 !Local variables ------------------------- 1594 !scalars 1595 integer,parameter :: master=0 1596 integer :: unddb 1597 integer :: spgroup, timrev, npsp 1598 integer :: choice 1599 integer :: ii, mblktyp,nblok 1600 integer :: matom_l,mtypat_l,mband_l,mkpt_l,msym_l,nsppol_l 1601 integer :: dimekb_l,lmnmax_l,usepaw_l 1602 integer :: comm_l 1603 !arrays 1604 real(dp):: rprimd(3,3) 1605 1606 ! ************************************************************************ 1607 1608 if (present(comm)) then 1609 comm_l = comm 1610 else 1611 comm_l = xmpi_comm_self 1612 end if 1613 1614 ! ===================================================== 1615 ! Open the DDB once and read the dimensions from header 1616 ! ===================================================== 1617 1618 call ddb_getdims(filename,comm_l,dimekb_l,lmnmax_l,mband_l,mblktyp, & 1619 & msym_l,matom_l,nblok,mkpt_l,nsppol_l,mtypat_l,usepaw_l) 1620 1621 1622 ! When merging ddbs, we want to open a header with fixed dimensions. 1623 ! GA: Might not be necessary if we handle the merging better. 1624 ! 1625 if (present(matom)) matom_l = matom 1626 if (present(mtypat)) mtypat_l = mtypat 1627 if (present(mband)) mband_l = mband 1628 if (present(mkpt)) mkpt_l = mkpt 1629 if (present(msym)) msym_l = msym 1630 if (present(dimekb)) dimekb_l = dimekb 1631 if (present(lmnmax)) lmnmax_l = lmnmax 1632 if (present(usepaw)) usepaw_l = usepaw 1633 1634 ddb_hdr%ddb_version = DDB_VERSION 1635 ddb_hdr%usepaw = usepaw_l 1636 ddb_hdr%mband = mband_l 1637 ddb_hdr%matom = matom_l 1638 ddb_hdr%msym = msym_l 1639 ddb_hdr%mtypat = mtypat_l 1640 ddb_hdr%mkpt = mkpt_l 1641 ddb_hdr%nsppol = nsppol_l 1642 1643 ddb_hdr%nblok = nblok 1644 ddb_hdr%mblktyp = mblktyp 1645 1646 ddb_hdr%psps%dimekb = dimekb_l 1647 ddb_hdr%psps%ntypat = mtypat_l 1648 ddb_hdr%psps%lmnmax = lmnmax_l 1649 ddb_hdr%psps%usepaw = usepaw_l 1650 ddb_hdr%psps%useylm = usepaw_l ! yep. 1651 1652 ! GA: Eventually, we might want to initialize other entries of psps here. 1653 ! such as alchemical psps info. 1654 1655 ddb_hdr%natom = ddb_hdr%matom 1656 ddb_hdr%nkpt = ddb_hdr%mkpt 1657 ddb_hdr%nsym = ddb_hdr%msym 1658 ddb_hdr%ntypat = ddb_hdr%mtypat 1659 1660 npsp = ddb_hdr%mtypat 1661 1662 ! Set maximal value for mpert 1663 ddb_hdr%mpert = ddb_hdr%natom+MPERT_MAX 1664 1665 ! Compute the block dimensions 1666 call ddb_hdr%get_block_dims() 1667 1668 unddb = get_unit() 1669 ddb_hdr%unddb = unddb 1670 1671 if (present(dimonly)) then 1672 if (dimonly==1) return 1673 end if 1674 1675 1676 ! ========================================= 1677 ! Open the DDB a second time to read header 1678 ! ========================================= 1679 1680 ddb_hdr%has_open_file_txt = .true. 1681 1682 if (xmpi_comm_rank(comm) == master) then 1683 1684 ! Allocate the memory 1685 call ddb_hdr%malloc() 1686 1687 ABI_MALLOC(ddb_hdr%psps%indlmn,(6,ddb_hdr%psps%lmnmax,ddb_hdr%mtypat)) 1688 ABI_MALLOC(ddb_hdr%psps%pspso,(ddb_hdr%mtypat)) 1689 ABI_MALLOC(ddb_hdr%psps%ekb,(ddb_hdr%psps%dimekb,ddb_hdr%mtypat)) 1690 1691 ABI_MALLOC(ddb_hdr%pawtab,(ddb_hdr%psps%ntypat*ddb_hdr%psps%usepaw)) 1692 call pawtab_nullify(ddb_hdr%pawtab) 1693 1694 ! This is needed to read the DDBs in the old format 1695 ! GA : Not sure if we really need this. 1696 ddb_hdr%symafm(:)=1 1697 if(ddb_hdr%mtypat>=1)then 1698 ddb_hdr%psps%pspso(:)=0 1699 ddb_hdr%psps%ekb(:,:)=zero 1700 ddb_hdr%znucl(:)=zero 1701 end if 1702 if(ddb_hdr%matom>=1)then 1703 ddb_hdr%spinat(:,:)=zero 1704 end if 1705 1706 ! Note: the maximum parameters (matom, mkpt, etc.) are inputs to ioddb8_in 1707 ! wile the actual parameters (natom, nkpt, etc.) are outputs 1708 call ioddb8_in(filename,ddb_hdr%matom,ddb_hdr%mband,& 1709 & ddb_hdr%mkpt,ddb_hdr%msym,ddb_hdr%mtypat,unddb,& 1710 & ddb_hdr%acell,ddb_hdr%amu,ddb_hdr%ddb_version,ddb_hdr%dilatmx,ddb_hdr%ecut,ddb_hdr%ecutsm,& 1711 & ddb_hdr%intxc,ddb_hdr%iscf,ddb_hdr%ixc,ddb_hdr%kpt,ddb_hdr%kptnrm,& 1712 & ddb_hdr%natom,ddb_hdr%nband,ddb_hdr%ngfft,ddb_hdr%nkpt,ddb_hdr%nspden,& 1713 & ddb_hdr%nspinor,ddb_hdr%nsppol,ddb_hdr%nsym,ddb_hdr%ntypat,& 1714 & ddb_hdr%occ,ddb_hdr%occopt,ddb_hdr%pawecutdg,ddb_hdr%rprim,& 1715 & ddb_hdr%dfpt_sciss,ddb_hdr%spinat,ddb_hdr%symafm,ddb_hdr%symrel,& 1716 & ddb_hdr%tnons,ddb_hdr%tolwfr,ddb_hdr%tphysel,ddb_hdr%tsmear,& 1717 & ddb_hdr%typat,ddb_hdr%usepaw,ddb_hdr%wtk,ddb_hdr%xred,ddb_hdr%zion,& 1718 & ddb_hdr%znucl) 1719 1720 1721 ! Read the psp information of the input DDB 1722 choice=1 ! Read 1723 call psddb8(choice,ddb_hdr%psps%dimekb,ddb_hdr%psps%ekb,ddb_hdr%with_psps,& 1724 & ddb_hdr%psps%indlmn,ddb_hdr%psps%lmnmax,& 1725 & ddb_hdr%nblok,ddb_hdr%ntypat,unddb,ddb_hdr%pawtab,ddb_hdr%psps%pspso,& 1726 & ddb_hdr%psps%usepaw,ddb_hdr%psps%useylm) 1727 1728 1729 ! This is a limitation of the text file format. 1730 ! We might be able to get rid of it. 1731 ddb_hdr%with_dfpt_vars = ddb_hdr%with_psps 1732 1733 ddb_hdr%dscrpt = '' 1734 1735 ! Complete the initialization of the psps object 1736 ! so that we may copy it. 1737 ddb_hdr%psps%npsp = npsp 1738 ddb_hdr%psps%ntypat = ddb_hdr%ntypat 1739 ddb_hdr%psps%ntypalch = 0 ! Alchemical potential not supported 1740 ddb_hdr%psps%mtypalch = 0 ! Alchemical potential not supported 1741 ddb_hdr%psps%n1xccc = 0 1742 ddb_hdr%psps%usewvl = 0 1743 ddb_hdr%psps%mqgrid_vl = 0 1744 ddb_hdr%psps%mqgrid_ff = 0 1745 ddb_hdr%psps%mpspso = 1 ! FIXME 1746 do ii=1,npsp 1747 ddb_hdr%psps%mpspso = max(ddb_hdr%psps%mpspso, ddb_hdr%psps%pspso(ii)) 1748 end do 1749 ddb_hdr%psps%mpssoang = 1 1750 ddb_hdr%psps%mpsang = 1 1751 ddb_hdr%psps%mproj = 0 1752 ddb_hdr%psps%lnmax = 0 1753 ddb_hdr%psps%ntyppure = 0 1754 ddb_hdr%psps%npspalch = 0 1755 ddb_hdr%psps%nc_xccc_gspace = 0 1756 ddb_hdr%psps%optnlxccc = 0 1757 ddb_hdr%psps%positron = 0 1758 ddb_hdr%psps%vlspl_recipSpace = .false. 1759 1760 ABI_MALLOC(ddb_hdr%psps%ziontypat,(ddb_hdr%psps%ntypat)) 1761 ABI_MALLOC(ddb_hdr%psps%znucltypat,(ddb_hdr%psps%ntypat)) 1762 ABI_MALLOC(ddb_hdr%psps%filpsp,(ddb_hdr%psps%npsp)) 1763 ABI_MALLOC(ddb_hdr%psps%md5_pseudos,(ddb_hdr%psps%npsp)) 1764 ABI_MALLOC(ddb_hdr%psps%title,(ddb_hdr%psps%npsp)) 1765 ABI_MALLOC(ddb_hdr%psps%pspxc,(ddb_hdr%psps%ntypat)) 1766 1767 ABI_MALLOC(ddb_hdr%psps%pspcod,(ddb_hdr%psps%npsp)) 1768 ABI_MALLOC(ddb_hdr%psps%algalch,(ddb_hdr%psps%ntypalch)) 1769 ABI_MALLOC(ddb_hdr%psps%pspdat,(ddb_hdr%psps%ntypat)) 1770 ABI_MALLOC(ddb_hdr%psps%zionpsp,(ddb_hdr%psps%npsp)) 1771 ABI_MALLOC(ddb_hdr%psps%znuclpsp,(ddb_hdr%psps%npsp)) 1772 1773 ABI_MALLOC(ddb_hdr%psps%nctab,(ddb_hdr%psps%ntypat)) 1774 1775 ddb_hdr%psps%ziontypat(:) = ddb_hdr%zion(:) 1776 ddb_hdr%psps%znucltypat(:) = ddb_hdr%znucl(:) 1777 do ii=1,npsp 1778 ddb_hdr%psps%filpsp(ii) = '' 1779 ddb_hdr%psps%md5_pseudos(ii) = '' 1780 ddb_hdr%psps%title(ii) = '' 1781 !ddb_hdr%psps%pspxc(ii) = ddb_hdr%ixc 1782 end do 1783 1784 ddb_hdr%psps%pspxc = zero ! Unknown 1785 1786 ! ------------------------------- 1787 ! Initialize the crystal 1788 ! ------------------------------- 1789 1790 call mkrdim(ddb_hdr%acell,ddb_hdr%rprim,rprimd) 1791 1792 ! GA: space group and time reversal are not written in the text file 1793 ! but they are written in the netcdf version. 1794 ! It doesnt seem to make any difference to anaddb. 1795 spgroup = 1 1796 timrev = 2 1797 1798 call crystal_init(ddb_hdr%amu, ddb_hdr%crystal, & 1799 & spgroup, ddb_hdr%natom, npsp, ddb_hdr%ntypat, & 1800 & ddb_hdr%nsym, rprimd, ddb_hdr%typat, & 1801 & ddb_hdr%xred, ddb_hdr%zion, ddb_hdr%znucl, timrev, & 1802 & ddb_hdr%nspden==2.and.ddb_hdr%nsppol==1, .false., '', & 1803 & ddb_hdr%symrel, ddb_hdr%tnons, ddb_hdr%symafm) 1804 1805 1806 end if 1807 1808 ! Master broadcasts data 1809 call ddb_hdr%bcast(comm_l) 1810 1811 end subroutine ddb_hdr_open_read_txt
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
filename=name of the file being written
with_psps
1-> include information on pseudopoentials
0-> do not include information on pseudopoentials
comm=MPI communicator.
OUTPUT
SOURCE
1457 subroutine ddb_hdr_open_write(ddb_hdr, filename, with_psps, comm) 1458 1459 !Arguments ------------------------------------ 1460 class(ddb_hdr_type),intent(inout) :: ddb_hdr 1461 character(len=*),intent(in) :: filename 1462 integer,intent(in),optional :: with_psps 1463 integer,intent(in),optional :: comm 1464 1465 !Local variables------------------------------- 1466 !scalars 1467 integer,parameter :: master=0 1468 integer :: comm_ 1469 integer :: iomode 1470 character(len=fnlen) :: filename_ 1471 1472 ! ************************************************************************ 1473 1474 if (present(comm)) then 1475 comm_ = comm 1476 else 1477 comm_ = xmpi_comm_self 1478 end if 1479 1480 if (xmpi_comm_rank(comm_) /= master) return 1481 1482 call ddb_hdr%get_iomode(filename, 2, iomode, filename_) 1483 1484 if (iomode==IO_MODE_ETSF) then 1485 1486 call ddb_hdr%open_write_nc(nctk_ncify(filename), with_psps=with_psps) 1487 1488 else if (iomode==IO_MODE_FORTRAN) then 1489 1490 call ddb_hdr%open_write_txt(filename, with_psps) 1491 1492 end if 1493 1494 end subroutine ddb_hdr_open_write
m_ddb_hdr/ddb_hdr_open_write_nc [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_open_write_nc
FUNCTION
Open the DDB NetCDF file and write the header. Create the groups and define all variables.
INPUTS
filename=name of the file being written
with_psps
1-> include information on pseudopoentials
0-> do not include information on pseudopoentials
OUTPUT
ncid=netcdf identifier.
SOURCE
746 subroutine ddb_hdr_open_write_nc(ddb_hdr, filename, with_psps, with_dfpt_vars) 747 748 !Arguments ------------------------------------ 749 class(ddb_hdr_type),intent(inout) :: ddb_hdr 750 character(len=*),intent(in) :: filename 751 integer,intent(in),optional :: with_psps 752 integer,intent(in),optional :: with_dfpt_vars 753 754 !Local variables ------------------------- 755 !scalars 756 integer :: ncid 757 integer :: ncid_crystal, ncid_psps, ncid_pawtab, ncid_d0E, ncid_d1E, ncid_d2E, ncid_d3E, ncid_d2eig 758 integer :: ncerr,nkpt,bandtot,occopt,blktyp 759 integer :: ii,iat,itypat,isppol,ikpt,iband,iblok,giblok 760 integer :: nblok_d0E,nblok_d1E,nblok_d2E,nblok_d3E,nblok_d2eig 761 integer :: iblok_d0E,iblok_d1E,iblok_d2E,iblok_d3E,iblok_d2eig 762 integer :: basis_size, max_basis_size, lmn_size, lmn2_size, max_lmn2_size, shape_type 763 integer :: natom, mpert 764 !arrays 765 integer,allocatable :: available_block_types(:) 766 integer,allocatable :: blktyp_d0E(:),blktyp_d1E(:),blktyp_d2E(:),blktyp_d3E(:),blktyp_d2eig(:) 767 integer,allocatable :: nband(:,:) 768 character(len=descrlen) :: descr 769 real(dp),allocatable :: occ(:,:,:) 770 integer,allocatable :: pawtab_basis_size(:) 771 integer,allocatable :: pawtab_lmn_size(:) 772 integer,allocatable :: pawtab_lmn2_size(:) 773 integer,allocatable :: pawtab_shape_type(:) 774 real(dp),allocatable :: pawtab_rshp(:) 775 real(dp),allocatable :: pawtab_rpaw(:) 776 777 ! ************************************************************************ 778 779 if (present(with_psps)) then 780 ddb_hdr%with_psps = with_psps 781 end if 782 783 if (present(with_dfpt_vars)) then 784 ddb_hdr%with_dfpt_vars = with_dfpt_vars 785 else 786 ddb_hdr%with_dfpt_vars = ddb_hdr%with_psps 787 end if 788 789 ! Initialize NetCDF file. 790 NCF_CHECK(nctk_open_create(ncid, filename, xmpi_comm_self)) 791 ddb_hdr%ncid = ncid 792 ddb_hdr%has_open_file_nc = .True. 793 794 ! Write version number 795 NCF_CHECK(nctk_defnwrite_ivars(ncid, ['ddb_version'], [DDB_VERSION_NC])) 796 797 ! ---------- 798 ! Dimensions 799 ! ---------- 800 801 ! Dimensions specified by ETSF 802 ncerr = nctk_def_dims(ncid, [& 803 nctkdim_t('number_of_cartesian_directions', 3), & 804 nctkdim_t('number_of_atoms', ddb_hdr%natom), & 805 nctkdim_t('number_of_kpoints', ddb_hdr%nkpt), & 806 nctkdim_t('number_of_spins', ddb_hdr%nsppol), & 807 nctkdim_t('number_of_spinor_components', ddb_hdr%nspinor), & 808 nctkdim_t('number_of_spin_densities', ddb_hdr%nspden), & 809 nctkdim_t('one_dim', 1), & 810 nctkdim_t('two_dim', 2), & 811 nctkdim_t('three_dim', 3), & 812 nctkdim_t('cplex',2)], defmode=.True.) 813 NCF_CHECK(ncerr) 814 815 ! Dimensions relevant to the DDB 816 ncerr = nctk_def_dims(ncid, [& 817 nctkdim_t('number_of_symmetries', ddb_hdr%nsym), & 818 nctkdim_t('number_of_atom_types', ddb_hdr%ntypat), & 819 nctkdim_t('maximum_number_of_bands', ddb_hdr%mband) & 820 ], defmode=.True.) 821 NCF_CHECK(ncerr) 822 823 ! Compute total number of bands 824 occopt = ddb_hdr%occopt 825 nkpt = ddb_hdr%nkpt 826 if(ddb_hdr%occopt==2)then 827 bandtot=0 828 do ikpt=1,nkpt 829 bandtot=bandtot+ddb_hdr%nband(ikpt)*ddb_hdr%nsppol 830 end do 831 else 832 bandtot=nkpt*ddb_hdr%nband(1)*ddb_hdr%nsppol 833 end if 834 ncerr = nctk_def_dims(ncid, [& 835 nctkdim_t('total_number_of_states', bandtot) & 836 ], defmode=.True.) 837 NCF_CHECK(ncerr) 838 839 ! Construct local array of occupations and nband 840 ABI_MALLOC(occ, (ddb_hdr%mband,ddb_hdr%nkpt,ddb_hdr%nsppol)) 841 occ = zero 842 ii = 0 843 do isppol=1,ddb_hdr%nsppol 844 do ikpt=1,ddb_hdr%nkpt 845 do iband=1,ddb_hdr%mband 846 ii = ii + 1 847 occ(iband,ikpt,isppol) = ddb_hdr%occ(ii) 848 end do 849 end do 850 end do 851 ABI_MALLOC(nband, (ddb_hdr%nkpt,ddb_hdr%nsppol)) 852 nband = zero 853 ii = 0 854 do isppol=1,ddb_hdr%nsppol 855 do ikpt=1,ddb_hdr%nkpt 856 ii = ii + 1 857 nband(ikpt,isppol) = ddb_hdr%nband(ii) 858 end do 859 end do 860 861 862 ncerr = nctk_defnwrite_ivars(ncid, ["with_psps"], [ddb_hdr%with_psps]) 863 NCF_CHECK(ncerr) 864 ncerr = nctk_defnwrite_ivars(ncid, ["with_dfpt_vars"], [ddb_hdr%with_dfpt_vars]) 865 !ncerr = nctk_defnwrite_ivars(ncid, & 866 !& [ character(len=nctk_slen) :: & 867 !& "with_psps", "with_dfpt_vars"], & 868 !& [ddb_hdr%with_psps, ddb_hdr%with_dfpt_vars]) 869 NCF_CHECK(ncerr) 870 871 ! -------------------------- 872 ! Write objects in subgroups 873 ! -------------------------- 874 875 ! Crystal info 876 NCF_CHECK(nf90_def_grp(ncid, 'Crystal', ncid_crystal)) 877 NCF_CHECK(ddb_hdr%crystal%ncwrite(ncid_crystal)) 878 879 ! Pseudopotentials info 880 NCF_CHECK(nf90_def_grp(ncid, 'Pseudopotentials', ncid_psps)) 881 call psps_ncwrite(ddb_hdr%psps, ncid_psps) ! We still need to write usepaw, for example. 882 883 ! PAW info 884 ! Note that PAW shares the namespace with Pseudopotentials. 885 if (ddb_hdr%psps%usepaw > 0) then 886 887 ! ------------------------------------- 888 ! GA: TODO This should be a function pawtab_ncwrite 889 NCF_CHECK(nf90_def_grp(ncid_psps, 'PAW', ncid_pawtab)) 890 891 ! get dimensions 892 ABI_MALLOC(pawtab_basis_size, (ddb_hdr%ntypat)) 893 ABI_MALLOC(pawtab_lmn_size, (ddb_hdr%ntypat)) 894 ABI_MALLOC(pawtab_lmn2_size, (ddb_hdr%ntypat)) 895 ABI_MALLOC(pawtab_shape_type, (ddb_hdr%ntypat)) 896 897 pawtab_basis_size = zero 898 pawtab_lmn_size = zero 899 pawtab_lmn2_size = zero 900 pawtab_shape_type = zero 901 902 max_basis_size = 0 ; max_lmn2_size = 0 903 do itypat=1,ddb_hdr%ntypat 904 905 basis_size = ddb_hdr%pawtab(itypat)%basis_size 906 lmn_size = ddb_hdr%pawtab(itypat)%lmn_size 907 lmn2_size = ddb_hdr%pawtab(itypat)%lmn2_size 908 shape_type = ddb_hdr%pawtab(itypat)%shape_type 909 910 pawtab_basis_size(itypat) = basis_size 911 pawtab_lmn_size(itypat) = lmn_size 912 pawtab_lmn2_size(itypat) = lmn2_size 913 pawtab_shape_type(itypat) = shape_type 914 915 if (lmn2_size > max_lmn2_size) then 916 max_lmn2_size = lmn2_size 917 end if 918 if (basis_size > max_basis_size) then 919 max_basis_size = basis_size 920 end if 921 end do 922 923 ncerr = nctk_def_dims(ncid_pawtab, [& 924 nctkdim_t('maximum_basis_size', max_basis_size), & 925 nctkdim_t('maximum_lmn2_size', max_lmn2_size)], & 926 defmode=.True.) 927 NCF_CHECK(ncerr) 928 929 ncerr = nctk_def_arrays(ncid_pawtab, [& 930 nctkarr_t('basis_size', "i", 'number_of_atom_types'), & 931 nctkarr_t('lmn_size', "i", 'number_of_atom_types'), & 932 nctkarr_t('lmn2_size', "i", 'number_of_atom_types'), & 933 nctkarr_t('shape_type', "i", 'number_of_atom_types') & 934 &]) 935 NCF_CHECK(ncerr) 936 937 ncerr = nctk_def_arrays(ncid_pawtab, [& 938 nctkarr_t('Dij0', "dp", 'number_of_atom_types, maximum_lmn2_size'), & 939 nctkarr_t('r_PAW', "dp", 'number_of_atom_types'), & 940 nctkarr_t('r_shape', "dp", 'number_of_atom_types') & 941 &]) 942 NCF_CHECK(ncerr) 943 944 NCF_CHECK(nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'basis_size'), pawtab_basis_size)) 945 NCF_CHECK(nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'lmn_size'), pawtab_lmn_size)) 946 NCF_CHECK(nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'lmn2_size'), pawtab_lmn2_size)) 947 NCF_CHECK(nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'shape_type'), pawtab_shape_type)) 948 949 ABI_FREE(pawtab_basis_size) 950 ABI_FREE(pawtab_lmn_size) 951 ABI_FREE(pawtab_lmn2_size) 952 ABI_FREE(pawtab_shape_type) 953 954 ABI_MALLOC(pawtab_rpaw, (ddb_hdr%ntypat)) 955 ABI_MALLOC(pawtab_rshp, (ddb_hdr%ntypat)) 956 957 do itypat=1,ddb_hdr%ntypat 958 959 pawtab_rpaw(itypat) = ddb_hdr%pawtab(itypat)%rpaw 960 pawtab_rshp(itypat) = ddb_hdr%pawtab(itypat)%rshp 961 962 ncerr = nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'Dij0'), & 963 ddb_hdr%pawtab(itypat)%dij0, start=[itypat,1], & 964 count=[1,ddb_hdr%pawtab(itypat)%lmn2_size]) 965 NCF_CHECK(ncerr) 966 967 end do 968 969 ncerr = nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'r_PAW'), pawtab_rpaw) 970 NCF_CHECK(ncerr) 971 972 ncerr = nf90_put_var(ncid_pawtab, nctk_idname(ncid_pawtab, 'r_shape'), pawtab_rshp) 973 NCF_CHECK(ncerr) 974 975 ABI_FREE(pawtab_rpaw) 976 ABI_FREE(pawtab_rshp) 977 978 ! ------------------------------------- 979 980 end if 981 982 983 ! ---------------- 984 ! Scalar variables 985 ! ---------------- 986 987 ! Scalar integer variables 988 ncerr = nctk_defnwrite_ivars(ncid, & 989 [character(len=nctk_slen) :: 'intxc', 'iscf', 'ixc', 'occopt'], & 990 [ddb_hdr%intxc, ddb_hdr%iscf, ddb_hdr%ixc, ddb_hdr%occopt]) 991 NCF_CHECK(ncerr) 992 993 ! Scalar float variables 994 ncerr = nctk_defnwrite_dpvars(ncid, & 995 [character(len=nctk_slen) :: 'dilatmx', 'ecut', 'pawecutdg', 'ecutsm', 'kptnrm'], & 996 [ddb_hdr%dilatmx, ddb_hdr%ecut, ddb_hdr%pawecutdg, ddb_hdr%ecutsm, ddb_hdr%kptnrm]) 997 NCF_CHECK(ncerr) 998 999 ncerr = nctk_defnwrite_dpvars(ncid, & 1000 [character(len=nctk_slen) :: 'dfpt_sciss', 'tolwfr', 'tphysel', 'tsmear'], & 1001 [ddb_hdr%dfpt_sciss, ddb_hdr%tolwfr, ddb_hdr%tphysel, ddb_hdr%tsmear]) 1002 NCF_CHECK(ncerr) 1003 1004 ! ------ 1005 ! Arrays 1006 ! ------ 1007 1008 ! FFT grid 1009 ncerr = nctk_def_arrays(ncid, [nctkarr_t('ngfft', "i", 'number_of_cartesian_directions')]) 1010 NCF_CHECK(ncerr) 1011 1012 ! Geometry 1013 ncerr = nctk_def_arrays(ncid, [& 1014 nctkarr_t('acell', "dp", 'number_of_cartesian_directions') & 1015 &]) 1016 NCF_CHECK(ncerr) 1017 1018 ! Basis 1019 ! note: In case occopt!=2, only the first entry of number_of_bands matters 1020 ncerr = nctk_def_arrays(ncid, [& 1021 nctkarr_t('reduced_coordinates_of_kpoints', "dp", 'number_of_cartesian_directions, number_of_kpoints'), & 1022 nctkarr_t('number_of_bands', "i", 'number_of_kpoints, number_of_spins'), & 1023 nctkarr_t('spinat', "dp", 'number_of_cartesian_directions, number_of_atoms'), & 1024 nctkarr_t('occupations', "dp", 'maximum_number_of_bands, number_of_kpoints, number_of_spins'), & 1025 nctkarr_t('kpoints_weights', "dp", 'number_of_kpoints') & 1026 &]) 1027 NCF_CHECK(ncerr) 1028 1029 ! ---------------- 1030 ! Set array values 1031 ! ---------------- 1032 1033 ! FFT grid 1034 NCF_CHECK(nf90_put_var(ncid, vid('ngfft'), ddb_hdr%ngfft(1:3))) 1035 1036 ! Geometry 1037 NCF_CHECK(nf90_put_var(ncid, vid('acell'), ddb_hdr%acell(:))) 1038 1039 ! Basis 1040 NCF_CHECK(nf90_put_var(ncid, vid('spinat'), ddb_hdr%spinat(:,:))) 1041 NCF_CHECK(nf90_put_var(ncid, vid('occupations'), occ)) 1042 NCF_CHECK(nf90_put_var(ncid, vid('reduced_coordinates_of_kpoints'), ddb_hdr%kpt(:,:))) 1043 NCF_CHECK(nf90_put_var(ncid, vid('kpoints_weights'), ddb_hdr%wtk(:))) 1044 NCF_CHECK(nf90_put_var(ncid, vid('number_of_bands'), nband)) 1045 1046 ABI_FREE(occ) 1047 ABI_FREE(nband) 1048 1049 ! ============================================================= 1050 ! Write info on blocks, create groups, and define all variables 1051 ! ============================================================= 1052 1053 ! -------------------------------------------- 1054 ! Description of perturbations and block types 1055 ! -------------------------------------------- 1056 mpert = ddb_hdr%mpert 1057 natom = ddb_hdr%natom 1058 1059 ncerr = nctk_def_dims(ncid, [& 1060 nctkdim_t('number_of_perturbations', ddb_hdr%mpert), & 1061 nctkdim_t("number_of_available_block_types", ntypes), & 1062 nctkdim_t("description_length", descrlen) & 1063 &], defmode=.True.) 1064 NCF_CHECK(ncerr) 1065 1066 ncerr = nctk_def_arrays(ncid, [& 1067 nctkarr_t('available_block_types', "i", 'number_of_available_block_types'), & 1068 nctkarr_t('description_of_block_types', "char", 'description_length, number_of_available_block_types'), & 1069 nctkarr_t('description_of_perturbations', "char", 'description_length, number_of_perturbations') & 1070 &]) 1071 NCF_CHECK(ncerr) 1072 1073 1074 do iat=1,natom 1075 1076 write(descr, '(2a,i7)') trim(DESCR_ipert_0), ' ', iat 1077 1078 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1079 trim(descr) , start=[1, iat]) 1080 NCF_CHECK(ncerr) 1081 end do 1082 1083 1084 if (mpert >= natom+1) then 1085 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1086 trim(DESCR_ipert_1), start=[1, natom+1]) 1087 NCF_CHECK(ncerr) 1088 end if 1089 1090 if (mpert >= natom+2) then 1091 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1092 trim(DESCR_ipert_2), start=[1, natom+2]) 1093 NCF_CHECK(ncerr) 1094 end if 1095 1096 if (mpert >= natom+3) then 1097 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1098 trim(DESCR_ipert_3), start=[1, natom+3]) 1099 NCF_CHECK(ncerr) 1100 end if 1101 1102 if (mpert >= natom+4) then 1103 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1104 trim(DESCR_ipert_4), start=[1, natom+4]) 1105 NCF_CHECK(ncerr) 1106 end if 1107 1108 if (mpert >= natom+5) then 1109 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1110 trim(DESCR_ipert_5), start=[1, natom+5]) 1111 NCF_CHECK(ncerr) 1112 end if 1113 1114 ! GA: The numbering natom+5, natom+10 , natom+11 comes from dfpt_looppert 1115 ! However, we dont want to increase mpert needlessly. 1116 ! Hence, we use natom+5, natom+6, natom+7. 1117 if (mpert >= natom+6) then 1118 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1119 trim(DESCR_ipert_10), start=[1, natom+6]) 1120 NCF_CHECK(ncerr) 1121 end if 1122 1123 if (mpert >= natom+7) then 1124 ncerr = nf90_put_var(ncid, vid("description_of_perturbations"), & 1125 trim(DESCR_ipert_11), start=[1, natom+7]) 1126 NCF_CHECK(ncerr) 1127 end if 1128 1129 1130 ABI_MALLOC(available_block_types, (ntypes)) 1131 ! The order is unimportant, as long as it matches the description order 1132 1133 available_block_types(1) = BLKTYP_d0E_xx 1134 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1135 trim(DESCR_d0E_xx), start=[1, 1]) 1136 NCF_CHECK(ncerr) 1137 1138 available_block_types(2) = BLKTYP_d1E_xx 1139 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1140 trim(DESCR_d1E_xx), start=[1, 2]) 1141 NCF_CHECK(ncerr) 1142 1143 available_block_types(3) = BLKTYP_d2E_ns 1144 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1145 trim(DESCR_d2E_ns), start=[1, 3]) 1146 NCF_CHECK(ncerr) 1147 1148 available_block_types(4) = BLKTYP_d2E_st 1149 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1150 trim(DESCR_d2E_st), start=[1, 4]) 1151 NCF_CHECK(ncerr) 1152 1153 available_block_types(5) = BLKTYP_d2E_mbc 1154 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1155 trim(DESCR_d2E_mbc), start=[1, 5]) 1156 NCF_CHECK(ncerr) 1157 1158 available_block_types(6) = BLKTYP_d3E_xx 1159 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1160 trim(DESCR_d3E_xx), start=[1, 6]) 1161 NCF_CHECK(ncerr) 1162 1163 available_block_types(7) = BLKTYP_d3E_lw 1164 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1165 trim(DESCR_d3E_lw), start=[1, 7]) 1166 NCF_CHECK(ncerr) 1167 1168 available_block_types(8) = BLKTYP_d2eig_re 1169 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1170 trim(DESCR_d2eig_re), start=[1, 8]) 1171 NCF_CHECK(ncerr) 1172 1173 available_block_types(9) = BLKTYP_d2eig_im 1174 ncerr = nf90_put_var(ncid, vid("description_of_block_types"),& 1175 trim(DESCR_d2eig_im), start=[1, 9]) 1176 NCF_CHECK(ncerr) 1177 1178 ncerr = nf90_put_var(ncid, vid("available_block_types"), available_block_types) 1179 NCF_CHECK(ncerr) 1180 ABI_FREE(available_block_types) 1181 1182 ncerr = nctk_def_dims(ncid, [nctkdim_t('number_of_blocks', ddb_hdr%nblok)]) 1183 1184 ! Info on blocks 1185 ncerr = nctk_def_arrays(ncid, [& 1186 nctkarr_t('block_types', "i", 'number_of_blocks') & 1187 &]) 1188 NCF_CHECK(ncerr) 1189 1190 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, 'block_types'), ddb_hdr%typ)) 1191 1192 ! TODO Write a chart of block types 1193 1194 ! Count each types of blocks 1195 nblok_d0E = 0 1196 nblok_d1E = 0 1197 nblok_d2E = 0 1198 nblok_d3E = 0 1199 nblok_d2eig = 0 1200 1201 do iblok=1,ddb_hdr%nblok 1202 blktyp = ddb_hdr%typ(iblok) 1203 if (is_type_d0E(blktyp)) nblok_d0E = nblok_d0E + 1 1204 if (is_type_d1E(blktyp)) nblok_d1E = nblok_d1E + 1 1205 if (is_type_d2E(blktyp)) nblok_d2E = nblok_d2E + 1 1206 if (is_type_d3E(blktyp)) nblok_d3E = nblok_d3E + 1 1207 if (is_type_d2eig(blktyp)) nblok_d2eig = nblok_d2eig + 1 1208 end do 1209 1210 ! ------------------------------------------ 1211 ! Build local arrays of global block indices 1212 ! to keep track of the initial order. 1213 ! ------------------------------------------ 1214 1215 ABI_MALLOC(blktyp_d0E, (nblok_d0E)) 1216 ABI_MALLOC(blktyp_d1E, (nblok_d1E)) 1217 ABI_MALLOC(blktyp_d2E, (nblok_d2E)) 1218 ABI_MALLOC(blktyp_d3E, (nblok_d3E)) 1219 ABI_MALLOC(blktyp_d2eig, (nblok_d2eig)) 1220 1221 iblok_d0E = 0 1222 iblok_d1E = 0 1223 iblok_d2E = 0 1224 iblok_d3E = 0 1225 iblok_d2eig = 0 1226 1227 do giblok=1,ddb_hdr%nblok 1228 1229 blktyp = ddb_hdr%typ(giblok) 1230 1231 if (is_type_d0E(blktyp)) then 1232 iblok_d0E = iblok_d0E + 1 1233 blktyp_d0E(iblok_d0E) = blktyp 1234 1235 else if (is_type_d1E(blktyp)) then 1236 iblok_d1E = iblok_d1E + 1 1237 blktyp_d1E(iblok_d1E) = blktyp 1238 1239 else if (is_type_d2E(blktyp)) then 1240 iblok_d2E = iblok_d2E + 1 1241 blktyp_d2E(iblok_d2E) = blktyp 1242 1243 else if (is_type_d3E(blktyp)) then 1244 iblok_d3E = iblok_d3E + 1 1245 blktyp_d3E(iblok_d3E) = blktyp 1246 1247 else if (is_type_d2eig(blktyp)) then 1248 iblok_d2eig = iblok_d2eig + 1 1249 blktyp_d2eig(iblok_d2eig) = blktyp 1250 end if 1251 1252 end do 1253 1254 ! ------------------------ 1255 ! Zeroth-order derivatives 1256 ! ------------------------ 1257 1258 NCF_CHECK(nf90_def_grp(ncid, 'd0E', ncid_d0E)) 1259 1260 ncerr = nctk_def_dims(ncid_d0E, [nctkdim_t('number_of_d0E_blocks', nblok_d0E)], defmode=.True.) 1261 NCF_CHECK(ncerr) 1262 1263 ! Info on blocks and matrix values 1264 ncerr = nctk_def_arrays(ncid_d0E, [& 1265 nctkarr_t('d0E_block_types', "i", 'number_of_d0E_blocks'), & 1266 nctkarr_t('matrix_values', "dp", 'number_of_d0E_blocks'), & 1267 nctkarr_t('matrix_mask', "i", 'number_of_d0E_blocks') & 1268 &]) 1269 NCF_CHECK(ncerr) 1270 1271 NCF_CHECK(nf90_put_var(ncid_d0E, nctk_idname(ncid_d0E, 'd0E_block_types'), blktyp_d0E)) 1272 1273 ABI_FREE(blktyp_d0E) 1274 1275 ! ----------------------- 1276 ! First-order derivatives 1277 ! ----------------------- 1278 1279 NCF_CHECK(nf90_def_grp(ncid, 'd1E', ncid_d1E)) 1280 1281 ncerr = nctk_def_dims(ncid_d1E, [nctkdim_t('number_of_d1E_blocks', nblok_d1E)], defmode=.True.) 1282 NCF_CHECK(ncerr) 1283 1284 ! Info on blocks and matrix values 1285 ncerr = nctk_def_arrays(ncid_d1E, [& 1286 nctkarr_t('d1E_block_types', "i", 'number_of_d1E_blocks'),& 1287 nctkarr_t('matrix_values', "dp",& 1288 sjoin('cplex,',& 1289 'number_of_cartesian_directions, number_of_perturbations, ',& 1290 'number_of_d1E_blocks')),& 1291 nctkarr_t('matrix_mask', "i",& 1292 sjoin('number_of_cartesian_directions, number_of_perturbations, ',& 1293 'number_of_d1E_blocks'))& 1294 &]) 1295 NCF_CHECK(ncerr) 1296 1297 NCF_CHECK(nf90_put_var(ncid_d1E, nctk_idname(ncid_d1E, 'd1E_block_types'), blktyp_d1E)) 1298 1299 ABI_FREE(blktyp_d1E) 1300 1301 ! ------------------------ 1302 ! Second-order derivatives 1303 ! ------------------------ 1304 1305 NCF_CHECK(nf90_def_grp(ncid, 'd2E', ncid_d2E)) 1306 1307 ! Qpoints 1308 ncerr = nctk_def_dims(ncid_d2E, [& 1309 nctkdim_t('number_of_d2E_blocks', nblok_d2E) & 1310 ], defmode=.True.) 1311 NCF_CHECK(ncerr) 1312 1313 ncerr = nctk_def_arrays(ncid_d2E, [& 1314 nctkarr_t('reduced_coordinates_of_qpoints', "dp",& 1315 'number_of_cartesian_directions, number_of_d2E_blocks') & 1316 &]) 1317 NCF_CHECK(ncerr) 1318 1319 ncerr = nctk_def_arrays(ncid_d2E, [& 1320 nctkarr_t('qpoints_normalization', "dp", 'number_of_d2E_blocks') & 1321 &]) 1322 NCF_CHECK(ncerr) 1323 1324 ! Info on blocks and matrix values 1325 ncerr = nctk_def_arrays(ncid_d2E, [& 1326 nctkarr_t('matrix_values', "dp",& 1327 sjoin('cplex, ',& 1328 'number_of_cartesian_directions, number_of_perturbations, ',& 1329 'number_of_cartesian_directions, number_of_perturbations, ',& 1330 'number_of_d2E_blocks')), & 1331 nctkarr_t('matrix_mask', "i",& 1332 sjoin('number_of_cartesian_directions, number_of_perturbations, ',& 1333 'number_of_cartesian_directions, number_of_perturbations, ',& 1334 'number_of_d2E_blocks')), & 1335 nctkarr_t('d2E_block_types', "i", 'number_of_d2E_blocks') & 1336 &]) 1337 NCF_CHECK(ncerr) 1338 1339 ABI_FREE(blktyp_d2E) 1340 1341 ! ----------------------- 1342 ! Third-order derivatives 1343 ! ----------------------- 1344 1345 NCF_CHECK(nf90_def_grp(ncid, 'd3E', ncid_d3E)) 1346 1347 ! Qpoints 1348 ncerr = nctk_def_dims(ncid_d3E, [& 1349 nctkdim_t('number_of_d3E_blocks', nblok_d3E) & 1350 ], defmode=.True.) 1351 NCF_CHECK(ncerr) 1352 1353 ncerr = nctk_def_arrays(ncid_d3E, [& 1354 nctkarr_t('reduced_coordinates_of_qpoints', "dp",& 1355 'three_dim, number_of_cartesian_directions, number_of_d3E_blocks') & 1356 &]) 1357 NCF_CHECK(ncerr) 1358 1359 ncerr = nctk_def_arrays(ncid_d3E, [& 1360 nctkarr_t('qpoints_normalization', "dp",& 1361 'three_dim, number_of_d3E_blocks') & 1362 &]) 1363 NCF_CHECK(ncerr) 1364 1365 ! Info on blocks and matrix values 1366 ncerr = nctk_def_arrays(ncid_d3E, [& 1367 nctkarr_t('matrix_values', "dp",& 1368 sjoin('cplex, ',& 1369 'number_of_cartesian_directions, number_of_perturbations, ',& 1370 'number_of_cartesian_directions, number_of_perturbations, ',& 1371 'number_of_cartesian_directions, number_of_perturbations, ',& 1372 'number_of_d3E_blocks')), & 1373 nctkarr_t('matrix_mask', "i",& 1374 sjoin('number_of_cartesian_directions, number_of_perturbations, ',& 1375 'number_of_cartesian_directions, number_of_perturbations, ',& 1376 'number_of_cartesian_directions, number_of_perturbations, ',& 1377 'number_of_d3E_blocks')), & 1378 nctkarr_t('d3E_block_types', "i", 'number_of_d3E_blocks') & 1379 &]) 1380 NCF_CHECK(ncerr) 1381 1382 NCF_CHECK(nf90_put_var(ncid_d3E, nctk_idname(ncid_d3E, 'd3E_block_types'), blktyp_d3E)) 1383 1384 ABI_FREE(blktyp_d3E) 1385 1386 ! --------------------------------------- 1387 ! Second-order derivatives of eigenvalues 1388 ! --------------------------------------- 1389 1390 NCF_CHECK(nf90_def_grp(ncid, 'd2eig', ncid_d2eig)) 1391 1392 ! Qpoints 1393 ncerr = nctk_def_dims(ncid_d2eig, [nctkdim_t('number_of_d2eig_blocks', nblok_d2eig)], defmode=.True.) 1394 NCF_CHECK(ncerr) 1395 1396 ncerr = nctk_def_arrays(ncid_d2eig, [& 1397 nctkarr_t('reduced_coordinates_of_qpoints', "dp",& 1398 'number_of_cartesian_directions, number_of_d2eig_blocks') & 1399 &]) 1400 NCF_CHECK(ncerr) 1401 1402 ncerr = nctk_def_arrays(ncid_d2eig, [& 1403 nctkarr_t('qpoints_normalization', "dp", 'number_of_d2eig_blocks') & 1404 &]) 1405 NCF_CHECK(ncerr) 1406 1407 ! Info on blocks and matrix values 1408 ncerr = nctk_def_arrays(ncid_d2eig, [& 1409 nctkarr_t('matrix_values', "dp",& 1410 sjoin('cplex, ',& 1411 'number_of_cartesian_directions, number_of_perturbations, ',& 1412 'number_of_cartesian_directions, number_of_perturbations, ',& 1413 'maximum_number_of_bands, number_of_kpoints, number_of_spins, ',& 1414 'number_of_d2eig_blocks')), & 1415 nctkarr_t('matrix_mask', "i",& 1416 sjoin('number_of_cartesian_directions, number_of_perturbations, ',& 1417 'number_of_cartesian_directions, number_of_perturbations, ',& 1418 'number_of_d2eig_blocks')), & 1419 nctkarr_t('d2eig_block_types', "i", 'number_of_d2eig_blocks') & 1420 &]) 1421 NCF_CHECK(ncerr) 1422 1423 NCF_CHECK(nf90_put_var(ncid_d2eig, nctk_idname(ncid_d2eig, 'd2eig_block_types'), blktyp_d2eig)) 1424 1425 ABI_FREE(blktyp_d2eig) 1426 1427 1428 contains 1429 integer function vid(vname) 1430 character(len=*),intent(in) :: vname 1431 vid = nctk_idname(ncid, vname) 1432 end function vid 1433 1434 end subroutine ddb_hdr_open_write_nc
m_ddb_hdr/ddb_hdr_open_write_txt [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_open_write_txt
FUNCTION
Open the DDB text file and write the header.
INPUTS
filename=name of the file being written (abo_DS*_DDB)
with_psps
1-> include information on pseudopoentials
0-> do not include information on pseudopoentials
OUTPUT
unddb=unit of the open ddb file in text format.
SOURCE
691 subroutine ddb_hdr_open_write_txt(ddb_hdr, filename, with_psps) 692 693 !Arguments ------------------------------------ 694 class(ddb_hdr_type),intent(inout) :: ddb_hdr 695 character(len=*),intent(in) :: filename 696 integer,intent(in),optional :: with_psps 697 698 !Local variables ------------------------- 699 character(len=500) :: message 700 integer :: unddb 701 integer :: ierr 702 703 ! ************************************************************************ 704 705 if (present(with_psps)) then 706 ddb_hdr%with_psps = with_psps 707 ddb_hdr%with_dfpt_vars = with_psps 708 end if 709 710 unddb = get_unit() 711 ddb_hdr%unddb = unddb 712 ddb_hdr%has_open_file_txt = .true. 713 714 ! Open the output derivative database. 715 ierr = open_file(filename,message,unit=unddb,status='unknown',form='formatted') 716 if (ierr /= 0) then 717 ABI_ERROR(message) 718 end if 719 720 call ddb_hdr%print(unddb) 721 722 end subroutine ddb_hdr_open_write_txt
m_ddb_hdr/ddb_hdr_print [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_print
FUNCTION
Print out the content of the header.
INPUTS
unddb=unit to print out the content.
OUTPUT
SOURCE
4606 subroutine ddb_hdr_print(ddb_hdr, unddb) 4607 4608 !Arguments ------------------------------------ 4609 class(ddb_hdr_type),intent(inout) :: ddb_hdr 4610 integer,intent(in) :: unddb 4611 4612 !Local variables ------------------------- 4613 integer,parameter :: choice=2 4614 4615 ! ************************************************************************ 4616 4617 call ddb_io_out(unddb,ddb_hdr%dscrpt,ddb_hdr%matom,ddb_hdr%mband,& 4618 & ddb_hdr%mkpt,ddb_hdr%msym,ddb_hdr%mtypat,& 4619 & ddb_hdr%acell,ddb_hdr%amu,ddb_hdr%dilatmx,ddb_hdr%ecut,ddb_hdr%ecutsm,& 4620 & ddb_hdr%intxc,ddb_hdr%iscf,ddb_hdr%ixc,ddb_hdr%kpt,ddb_hdr%kptnrm,& 4621 & ddb_hdr%natom,ddb_hdr%nband,ddb_hdr%ngfft,ddb_hdr%nkpt,ddb_hdr%nspden,& 4622 & ddb_hdr%nspinor,ddb_hdr%nsppol,ddb_hdr%nsym,ddb_hdr%ntypat,ddb_hdr%occ,& 4623 & ddb_hdr%occopt,ddb_hdr%pawecutdg,ddb_hdr%rprim,ddb_hdr%dfpt_sciss,& 4624 & ddb_hdr%spinat,ddb_hdr%symafm,ddb_hdr%symrel,ddb_hdr%tnons,ddb_hdr%tolwfr,& 4625 & ddb_hdr%tphysel,ddb_hdr%tsmear,ddb_hdr%typat,ddb_hdr%usepaw,ddb_hdr%wtk,& 4626 & ddb_hdr%xred,ddb_hdr%zion,ddb_hdr%znucl) 4627 4628 call psddb8(choice,ddb_hdr%psps%dimekb,ddb_hdr%psps%ekb,ddb_hdr%with_psps,& 4629 & ddb_hdr%psps%indlmn,ddb_hdr%psps%lmnmax,ddb_hdr%nblok,ddb_hdr%ntypat,unddb,& 4630 & ddb_hdr%pawtab,ddb_hdr%psps%pspso,ddb_hdr%psps%usepaw,ddb_hdr%psps%useylm) 4631 4632 end subroutine ddb_hdr_print
m_ddb_hdr/ddb_hdr_set_typ [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
ddb_hdr_set_typ
FUNCTION
Set the type of each block.
INPUTS
nblok=number of blocks typ=type of each block
SOURCE
573 subroutine ddb_hdr_set_typ(ddb_hdr, nblok, typ) 574 575 !Arguments ------------------------------------ 576 class(ddb_hdr_type),intent(inout) :: ddb_hdr 577 integer,intent(in) :: nblok 578 integer,intent(in) :: typ(nblok) 579 580 ! ************************************************************************ 581 582 ddb_hdr%nblok = nblok 583 ABI_SFREE(ddb_hdr%typ) 584 ABI_MALLOC(ddb_hdr%typ,(nblok)) 585 ddb_hdr%typ(:) = typ(:) 586 587 end subroutine ddb_hdr_set_typ
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
unddb=unit number for output 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) filename: 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 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
OUTPUT
Only writing
SOURCE
4797 subroutine ddb_io_out (unddb,dscrpt,matom,mband,& 4798 & mkpt,msym,mtypat,& 4799 & acell,amu,dilatmx,ecut,ecutsm,intxc,iscf,ixc,kpt,kptnrm,& 4800 & natom,nband,ngfft,nkpt,nspden,nspinor,nsppol,nsym,ntypat,occ,occopt,& 4801 & pawecutdg,rprim,dfpt_sciss,spinat,symafm,symrel,tnons,tolwfr,tphysel,tsmear,& 4802 & typat,usepaw,wtk,xred,zion,znucl) 4803 4804 !Arguments ------------------------------- 4805 !scalars 4806 integer,intent(in) :: unddb,matom,mband,mkpt,msym,mtypat 4807 integer,intent(in) :: intxc,iscf,ixc,natom,nkpt,nspden,nspinor,nsppol,nsym 4808 integer,intent(in) :: ntypat,occopt,usepaw 4809 real(dp),intent(in) :: dilatmx,ecut,ecutsm,kptnrm,pawecutdg,dfpt_sciss,tolwfr,tphysel 4810 real(dp),intent(in) :: tsmear 4811 character(len=fnlen),intent(in) :: dscrpt 4812 !arrays 4813 integer,intent(in) :: nband(mkpt*nsppol),ngfft(18),symafm(msym),symrel(3,3,msym) 4814 integer,intent(in) :: typat(matom) 4815 real(dp),intent(in) :: acell(3),amu(mtypat),kpt(3,mkpt),occ(mband*mkpt*nsppol) 4816 real(dp),intent(in) :: rprim(3,3),spinat(3,matom),tnons(3,msym),wtk(mkpt) 4817 real(dp),intent(in) :: xred(3,matom),zion(mtypat),znucl(mtypat) 4818 4819 !Local variables ------------------------- 4820 !Set routine version number here: 4821 !scalars 4822 integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527 4823 integer :: bantot,ii,ij,ikpt,iline,im 4824 !arrays 4825 character(len=9) :: name(9) 4826 4827 ! ********************************************************************* 4828 4829 DBG_ENTER("COLL") 4830 4831 4832 !Write the heading 4833 write(unddb, '(/,a,/,a,i10,/,/,a,a,/)' ) & 4834 ' **** DERIVATIVE DATABASE **** ',& 4835 '+DDB, Version number',DDB_VERSION,' ',trim(dscrpt) 4836 4837 !Write the descriptive data 4838 !1. usepaw 4839 write(unddb, '(1x,a9,i10)' )' usepaw',usepaw 4840 !2. natom 4841 write(unddb, '(1x,a9,i10)' )' natom',natom 4842 !3. nkpt 4843 write(unddb, '(1x,a9,i10)' )' nkpt',nkpt 4844 !4. nsppol 4845 write(unddb, '(1x,a9,i10)' )' nsppol',nsppol 4846 !5. nsym 4847 write(unddb, '(1x,a9,i10)' )' nsym',nsym 4848 !6. ntypat 4849 write(unddb, '(1x,a9,i10)' )' ntypat',ntypat 4850 !7. occopt 4851 write(unddb, '(1x,a9,i10)' )' occopt',occopt 4852 !8. nband 4853 if(occopt==2)then 4854 im=12 4855 name(1)=' nband' 4856 do iline=1,(nkpt+11)/12 4857 if(iline==(nkpt+11)/12)im=nkpt-12*(iline-1) 4858 write(unddb, '(1x,a9,5x,12i5)' )name(1),(nband((iline-1)*12+ii),ii=1,im) 4859 name(1)=' ' 4860 end do 4861 bantot=0 4862 do ikpt=1,nkpt 4863 bantot=bantot+nband(ikpt) 4864 end do 4865 else 4866 write(unddb, '(1x,a9,i10)' )' nband',nband(1) 4867 bantot=nkpt*nband(1) 4868 end if 4869 4870 !9. acell 4871 write(unddb, '(1x,a9,3d22.14)' )' acell',acell 4872 !10. amu 4873 im=3 4874 name(1)=' amu' 4875 do iline=1,(ntypat+2)/3 4876 if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1) 4877 write (unddb, '(1x,a9,3d22.14)' )name(1),(amu((iline-1)*3+ii),ii=1,im) 4878 name(1)=' ' 4879 end do 4880 !11. dilatmx 4881 write(unddb, '(1x,a9,d22.14)' )' dilatmx',dilatmx 4882 !12. ecut 4883 write(unddb, '(1x,a9,d22.14)' )' ecut',ecut 4884 !12b. pawecutdg (PAW) 4885 if (usepaw==1) then 4886 write(unddb, '(1x,a9,d22.14)' )'pawecutdg',pawecutdg 4887 end if 4888 !13. ecutsm 4889 write(unddb, '(1x,a9,d22.14)' )' ecutsm',ecutsm 4890 !14. intxc 4891 write(unddb, '(1x,a9,i10)' )' intxc',intxc 4892 !15. iscf 4893 write(unddb, '(1x,a9,i10)' )' iscf',iscf 4894 !16. ixc 4895 write(unddb, '(1x,a9,i10)' )' ixc',ixc 4896 !17. kpt 4897 name(1)=' kpt' 4898 do iline=1,nkpt 4899 write (unddb, '(1x,a9,3d22.14)' )name(1),(kpt(ii,iline),ii=1,3) 4900 name(1)=' ' 4901 end do 4902 !18. kptnrm 4903 write(unddb, '(1x,a9,d22.14)' )' kptnrm',kptnrm 4904 !19. ngfft 4905 write(unddb, '(1x,a9,5x,3i5)' )' ngfft',ngfft(1:3) 4906 !20. nspden 4907 write(unddb, '(1x,a9,i10)' )' nspden',nspden 4908 !21. nspinor 4909 write(unddb, '(1x,a9,i10)' )' nspinor',nspinor 4910 !22. occ 4911 if(occopt==2)then 4912 im=3 4913 name(1)=' occ' 4914 do iline=1,(bantot+2)/3 4915 if(iline==(bantot+2)/3)im=bantot-3*(iline-1) 4916 write(unddb, '(1x,a9,3d22.14)' )name(1),(occ((iline-1)*3+ii),ii=1,im) 4917 name(1)=' ' 4918 end do 4919 else 4920 im=3 4921 name(1)=' occ' 4922 do iline=1,(nband(1)+2)/3 4923 if(iline==(nband(1)+2)/3)im=nband(1)-3*(iline-1) 4924 write(unddb, '(1x,a9,3d22.14)' )name(1),(occ((iline-1)*3+ii),ii=1,im) 4925 name(1)=' ' 4926 end do 4927 end if 4928 !23. rprim 4929 name(1)=' rprim' 4930 do iline=1,3 4931 write(unddb, '(1x,a9,3d22.14)' )name(1),(rprim(ii,iline),ii=1,3) 4932 name(1)=' ' 4933 end do 4934 !24. dfpt_sciss 4935 write(unddb, '(1x,a11,d22.14)' )' dfpt_sciss',dfpt_sciss 4936 !25. spinat 4937 name(1)=' spinat' 4938 do iline=1,natom 4939 write(unddb, '(1x,a9,3d22.14)' )name(1),(spinat(ii,iline),ii=1,3) 4940 name(1)=' ' 4941 end do 4942 !26. symafm 4943 im=12 4944 name(1)=' symafm' 4945 do iline=1,(nsym+11)/12 4946 if(iline==(nsym+11)/12)im=nsym-12*(iline-1) 4947 write(unddb, '(1x,a9,5x,12i5)' )name(1),(symafm((iline-1)*12+ii),ii=1,im) 4948 name(1)=' ' 4949 end do 4950 !27. symrel 4951 name(1)=' symrel' 4952 do iline=1,nsym 4953 write(unddb, '(1x,a9,5x,9i5)' )name(1),((symrel(ii,ij,iline),ii=1,3),ij=1,3) 4954 name(1)=' ' 4955 end do 4956 !28. tnons 4957 name(1)=' tnons' 4958 do iline=1,nsym 4959 write(unddb, '(1x,a9,3d22.14)' )name(1),(tnons(ii,iline),ii=1,3) 4960 name(1)=' ' 4961 end do 4962 !29. tolwfr 4963 write(unddb, '(1x,a9,d22.14)' )' tolwfr',tolwfr 4964 !30. tphysel 4965 write(unddb, '(1x,a9,d22.14)' )' tphysel',tphysel 4966 !31. tsmear 4967 write(unddb, '(1x,a9,d22.14)' )' tsmear',tsmear 4968 !32. typat 4969 im=12 4970 name(1)=' typat' 4971 do iline=1,(natom+11)/12 4972 if(iline==(natom+11)/12)im=natom-12*(iline-1) 4973 write(unddb, '(1x,a9,5x,12i5)' )name(1),(typat((iline-1)*12+ii),ii=1,im) 4974 name(1)=' ' 4975 end do 4976 !33. wtk 4977 name(1)=' wtk' 4978 im=3 4979 do iline=1,(nkpt+2)/3 4980 if(iline==(nkpt+2)/3)im=nkpt-3*(iline-1) 4981 write(unddb, '(1x,a9,3d22.14)' )name(1),(wtk((iline-1)*3+ii),ii=1,im) 4982 name(1)=' ' 4983 end do 4984 !34. xred 4985 name(1)=' xred' 4986 do iline=1,natom 4987 write(unddb, '(1x,a9,3d22.14)' )name(1),(xred(ii,iline),ii=1,3) 4988 name(1)=' ' 4989 end do 4990 !35. znucl 4991 name(1)=' znucl' 4992 im=3 4993 do iline=1,(ntypat+2)/3 4994 if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1) 4995 write(unddb, '(1x,a9,3d22.14)' )name(1),(znucl((iline-1)*3+ii),ii=1,im) 4996 name(1)=' ' 4997 end do 4998 !36. zion 4999 name(1)=' zion' 5000 im=3 5001 do iline=1,(ntypat+2)/3 5002 if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1) 5003 write(unddb, '(1x,a9,3d22.14)' )name(1),(zion((iline-1)*3+ii),ii=1,im) 5004 name(1)=' ' 5005 end do 5006 5007 DBG_EXIT("COLL") 5008 5009 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=*) filename: name of input or output file unddb=unit number for input or output
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
SOURCE
4042 subroutine inprep8 (filename,unddb,dimekb,lmnmax,mband,mblktyp,msym,natom,nblok,nkpt,& 4043 & nsppol,ntypat,usepaw) 4044 4045 !Arguments ------------------------------- 4046 !scalars 4047 integer,intent(in) :: unddb 4048 integer, intent(out) :: msym 4049 integer,intent(out) :: dimekb,lmnmax,mband,mblktyp,natom,nblok,nkpt,ntypat,nsppol,usepaw 4050 character(len=*),intent(in) :: filename 4051 4052 !Local variables ------------------------- 4053 !scalars 4054 !Set routine version number here: 4055 integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527 4056 integer,parameter :: cvrsio9=20230401,cvrsio8=20100401,cvrsio8_old=20010929,cvrsio8_old_old=19990527 4057 integer :: bantot,basis_size0,blktyp,ddbvrs,iband,iblok,iekb,ii,ikpt,iline,im,ios,iproj 4058 integer :: itypat,itypat0,jekb,lmn_size0,mproj,mpsang,nekb,ndig,nelmts 4059 integer :: occopt,pspso0,nsym 4060 logical :: ddbvrs_is_current_or_old,testn,testv 4061 character(len=12) :: string 4062 character(len=32) :: blkname 4063 character(len=500) :: message 4064 character(len=6) :: name_old, ddbvrs6 4065 character(len=80) :: rdstring 4066 character(len=8) :: ddbvrs8 4067 character(len=3) :: prefix 4068 !arrays 4069 integer,allocatable :: nband(:) 4070 character(len=12) :: name(9) 4071 4072 ! ********************************************************************* 4073 4074 !Open the input derivative database. 4075 if (open_file(filename,message,unit=unddb,form="formatted",status="old",action="read") /= 0) then 4076 ABI_ERROR(message) 4077 end if 4078 4079 !Check the compatibility of the input DDB with the DDB code 4080 read (unddb,*) 4081 read (unddb,*) 4082 read (unddb, '(20x,i10)' )ddbvrs 4083 4084 if (all(ddbvrs/= [cvrsio9, vrsio8, vrsio8_old, vrsio8_old_old]) )then 4085 write(message, '(a,i10,2a,4(a,i10))' )& 4086 & 'The input DDB version number=',ddbvrs,' does not agree',ch10,& 4087 & 'with the allowed code DDB version numbers,',cvrsio9,', ',vrsio8,', ',vrsio8_old,' and ',vrsio8_old_old 4088 ABI_ERROR(message) 4089 end if 4090 4091 !Convert older version to 8 digit format 4092 if (ddbvrs /= cvrsio9) then 4093 ndig= int(log10(real(ddbvrs))) + 1 4094 write(ddbvrs6,'(i0)') ddbvrs 4095 if (ddbvrs==vrsio8 .or.ddbvrs==vrsio8_old) then 4096 if (ndig==6) then 4097 write(prefix,'(i2)') 20 4098 else if (ndig==5) then 4099 write(prefix,'(i3)') 200 4100 end if 4101 else if (ddbvrs==vrsio8_old_old) then 4102 if (ndig==6) then 4103 write(prefix,'(i2)') 19 4104 else if (ndig==5) then 4105 write(prefix,'(i3)') 199 4106 end if 4107 end if 4108 ddbvrs8= trim(prefix) // trim(ddbvrs6) 4109 read(ddbvrs8,'(i8)') ddbvrs 4110 end if 4111 4112 !Read the 4 n-integers, also testing the names of data, 4113 !and checking that their value is acceptable. 4114 !This is important to insure that any array has a sufficient dimension. 4115 read (unddb,*) 4116 read (unddb,*) 4117 read (unddb,*) 4118 testn=.true. 4119 testv=.true. 4120 ! ddbvrs_is_current_or_old=(ddbvrs==vrsio8.or.ddbvrs==vrsio8_old) 4121 ddbvrs_is_current_or_old=(ddbvrs>=cvrsio8_old) 4122 4123 !1. usepaw 4124 if(ddbvrs>=cvrsio8)then 4125 read (unddb, '(1x,a9,i10)' )name(1),usepaw 4126 else 4127 usepaw=0;name(1)=' usepaw' 4128 end if 4129 if(name(1)/=' usepaw')testn=.false. 4130 !2. natom 4131 if(ddbvrs_is_current_or_old)then 4132 read (unddb, '(1x,a9,i10)' )name(2),natom 4133 else 4134 read (unddb, '(1x,a6,i10)' )name_old,natom ; name(2)=' '//name_old 4135 end if 4136 if(name(2)/=' natom')testn=.false. 4137 if(natom<=0)testv=.false. 4138 !3. nkpt 4139 if(ddbvrs_is_current_or_old)then 4140 read (unddb, '(1x,a9,i10)' )name(3),nkpt 4141 else 4142 read (unddb, '(1x,a6,i10)' )name_old,nkpt ; name(3)=' '//name_old 4143 end if 4144 if(name(3)/=' nkpt')testn=.false. 4145 if(nkpt <=0)testv=.false. 4146 !4. nsppol 4147 if(ddbvrs_is_current_or_old)then 4148 read (unddb, '(1x,a9,i10)' )name(4),nsppol 4149 else 4150 read (unddb, '(1x,a6,i10)' )name_old,nsppol ; name(4)=' '//name_old 4151 end if 4152 if(name(4)/=' nsppol')testn=.false. 4153 if(nsppol<=0.or.nsppol>2)testv=.false. 4154 !5. nsym 4155 if(ddbvrs_is_current_or_old)then 4156 read (unddb, '(1x,a9,i10)' )name(5),nsym 4157 else 4158 read (unddb, '(1x,a6,i10)' )name_old,nsym ; name(5)=' '//name_old 4159 end if 4160 if(name(5)/=' nsym')testn=.false. 4161 !MG FIXME Why this and why do we need msym? 4162 msym = 192 4163 if (nsym > msym) msym=nsym 4164 !if(nsym <=0.or.nsym >msym )testv=.false. 4165 !6. ntypat 4166 if(ddbvrs_is_current_or_old)then 4167 read (unddb, '(1x,a9,i10)' )name(6),ntypat 4168 else 4169 read (unddb, '(1x,a6,i10)' )name_old,ntypat ; name(6)=' '//name_old 4170 end if 4171 if(name(6)/=' ntypat' .and. name(6)/=' ntype')testn=.false. 4172 if(ntypat<=0)testv=.false. 4173 !7. occopt 4174 !Before reading nband, the last parameters that define 4175 !the dimension of some array, need to know what is their 4176 !representation, given by occopt 4177 if(ddbvrs_is_current_or_old)then 4178 read (unddb, '(1x,a9,i10)' )name(7),occopt 4179 else 4180 read (unddb, '(1x,a6,i10)' )name_old,occopt ; name(7)=' '//name_old 4181 end if 4182 if(name(7)/=' occopt')testn=.false. 4183 if(occopt<0.or.occopt>8)testv=.false. 4184 4185 !Message if the names or values are not right 4186 if (.not.testn.or..not.testv) then 4187 write(message, '(a,a)' )' inprep8 : An error has been found in the',' positive n-integers contained in the DDB : ' 4188 call wrtout(std_out,message,'COLL') 4189 write(message, '(a)' ) ' Expected Found ' 4190 call wrtout(std_out,message,'COLL') 4191 write(message, '(a,i9,a,a,a,i10)' )' natom , larger than',0,' ',trim(name(2)),' =',natom 4192 call wrtout(std_out,message,'COLL') 4193 write(message, '(a,i9,a,a,a,i10)' )' nkpt , larger than',0,' ',trim(name(3)),' =',nkpt 4194 call wrtout(std_out,message,'COLL') 4195 write(message, '(a,i1,a,a,a,i10)' )' nsppol, either 1 or ',2,' ',trim(name(4)),' =',nsppol 4196 call wrtout(std_out,message,'COLL') 4197 ! write(message, '(a,i10,a,a,a,i10)' )& ' nsym , lower than',msym,' ',trim(name(5)),' =',nsym 4198 call wrtout(std_out,message,'COLL') 4199 write(message, '(a,i9,a,a,a,i10)' )' ntypat , larger than',0,' ',trim(name(6)),' =',ntypat 4200 call wrtout(std_out,message,'COLL') 4201 write(message, '(a,a,a,i10)' )' occopt, equal to 0,1 or 2 ',trim(name(7)),' =',occopt 4202 call wrtout(std_out,message,'COLL') 4203 4204 ABI_ERROR('See the error message above.') 4205 end if 4206 4207 !One more set of parameters define the dimensions of the 4208 !array : nband. Morever, it depends on occopt and !nkpt, and has to be 4209 !tested after the test on nkpt is performed. 4210 4211 !8. nband 4212 ABI_MALLOC(nband,(nkpt)) 4213 if(occopt==2)then 4214 im=12 4215 do iline=1,(nkpt+11)/12 4216 if(iline==(nkpt+11)/12)im=nkpt-12*(iline-1) 4217 if(ddbvrs_is_current_or_old)then 4218 read (unddb, '(1x,a9,5x,12i5)' )name(1),(nband((iline-1)*12+ii),ii=1,im) 4219 else 4220 read (unddb, '(1x,a6,5x,12i5)' )name_old,& 4221 & (nband((iline-1)*12+ii),ii=1,im) ; name(1)=' '//name_old 4222 end if 4223 if (iline==1) then 4224 call ddb_chkname(name(1),' nband') 4225 else 4226 call ddb_chkname(name(1),' ') 4227 end if 4228 end do 4229 else 4230 if(ddbvrs_is_current_or_old)then 4231 read (unddb, '(1x,a9,i10)' )name(1),nband(1) 4232 else 4233 read (unddb, '(1x,a6,i10)' )name_old,nband(1) ; name(1)=' '//name_old 4234 end if 4235 call ddb_chkname(name(1),' nband') 4236 if(nkpt>1)then 4237 do ikpt=2,nkpt 4238 nband(ikpt)=nband(1) 4239 end do 4240 end if 4241 end if 4242 4243 !Check all nband values, and sum them 4244 bantot=0 4245 do ikpt=1,nkpt 4246 if(nband(ikpt)<0)then 4247 write(message, '(a,i0,a,i0,3a)' )& 4248 & 'For ikpt = ',ikpt,' nband = ',nband(ikpt),' is negative.',ch10,& 4249 & 'Action: correct your DDB.' 4250 ABI_ERROR(message) 4251 end if 4252 bantot=bantot+nband(ikpt) 4253 end do 4254 4255 mband=maxval(nband(:)) 4256 4257 !Skip the rest of variables 4258 !9. acell 4259 read (unddb,*) 4260 !10. amu 4261 do iline=1,(ntypat+2)/3 4262 read (unddb,*) 4263 end do 4264 !11. dilatmx 4265 if(ddbvrs_is_current_or_old) read (unddb,*) 4266 !12. ecut 4267 read (unddb,*) 4268 !12b. pawecutdg (PAW only) 4269 if(ddbvrs>=cvrsio8.and.usepaw==1) read (unddb,*) 4270 !13. ecutsm 4271 if(ddbvrs_is_current_or_old) read (unddb,*) 4272 !14. intxc 4273 if(ddbvrs_is_current_or_old) read (unddb,*) 4274 !15. iscf 4275 read (unddb,*) 4276 !16. ixc 4277 read (unddb,*) 4278 !17. kpt 4279 do iline=1,nkpt 4280 read (unddb,*) 4281 end do 4282 !18. kptnrm 4283 read (unddb,*) 4284 !19. ngfft 4285 read (unddb,*) 4286 !20. nspden 4287 if(ddbvrs_is_current_or_old) read (unddb,*) 4288 !21. nspinor 4289 if(ddbvrs_is_current_or_old) read (unddb,*) 4290 !22. occ 4291 if(occopt==2)then 4292 do iline=1,(bantot+2)/3 4293 read (unddb,*) 4294 end do 4295 else 4296 !write(message,*)' inprep8 : nband(1)=',nband(1) 4297 !call wrtout(std_out,message,'COLL') 4298 do iline=1,(nband(1)+2)/3 4299 read (unddb,'(a80)')rdstring 4300 !write(message,*)trim(rdstring) 4301 !call wrtout(std_out,message,'COLL') ! GA: why are we printing this? 4302 end do 4303 end if 4304 !23. rprim 4305 do iline=1,3 4306 read (unddb,*) 4307 end do 4308 !24. dfpt_sciss 4309 read (unddb,*) 4310 !25. spinat 4311 if(ddbvrs_is_current_or_old)then 4312 do iline=1,natom 4313 read (unddb,*) 4314 end do 4315 end if 4316 !26. symafm 4317 if(ddbvrs_is_current_or_old)then 4318 do iline=1,(nsym+11)/12 4319 read (unddb,*) 4320 end do 4321 end if 4322 !27. symrel 4323 do iline=1,nsym 4324 read (unddb,*) 4325 end do 4326 !28old. xred 4327 if(.not.ddbvrs_is_current_or_old)then 4328 do iline=1,natom 4329 read (unddb,*) 4330 end do 4331 end if 4332 !28. tnons 4333 do iline=1,nsym 4334 read (unddb,*) 4335 end do 4336 !29. tolwfr 4337 if(ddbvrs_is_current_or_old) read (unddb,*) 4338 !30. tphysel 4339 if(ddbvrs_is_current_or_old) read (unddb,*) 4340 !31. tsmear 4341 if(ddbvrs_is_current_or_old) read (unddb,*) 4342 !32. type 4343 do iline=1,(natom+11)/12 4344 read (unddb,*) 4345 end do 4346 !33old. tolwfr 4347 if(.not.ddbvrs_is_current_or_old) read (unddb,*) 4348 !33. wtk 4349 do iline=1,(nkpt+2)/3 4350 read (unddb,*) 4351 end do 4352 !34. xred 4353 if(ddbvrs_is_current_or_old)then 4354 do iline=1,natom 4355 read (unddb,*) 4356 end do 4357 end if 4358 !35. znucl 4359 if(ddbvrs_is_current_or_old)then 4360 do iline=1,(ntypat+2)/3 4361 read (unddb,*) 4362 end do 4363 end if 4364 !36. zion 4365 do iline=1,(ntypat+2)/3 4366 read (unddb,*) 4367 end do 4368 4369 read (unddb,*) 4370 4371 !Now, take care of the pseudopotentials 4372 read(unddb, '(a12)' )string 4373 4374 if(string==' Descriptio')then 4375 4376 read (unddb,*) 4377 if (ddbvrs==cvrsio8_old.or.ddbvrs==cvrsio8_old_old) then 4378 read (unddb, '(10x,i3,14x,i3,11x,i3)', iostat=ios )dimekb,lmnmax,usepaw 4379 if(ios/=0)then 4380 backspace(unddb) 4381 read (unddb, '(10x,i3,14x,i3)')dimekb,lmnmax 4382 usepaw=0 4383 end if 4384 else if (ddbvrs>=cvrsio8) then 4385 read (unddb, '(10x,i3)') usepaw 4386 if (usepaw==0) then 4387 read (unddb, '(10x,i3,14x,i3)' ) dimekb,lmnmax 4388 else 4389 dimekb=0;lmnmax=0 4390 end if 4391 end if 4392 if (usepaw==0) then 4393 do itypat=1,ntypat 4394 read(unddb, '(13x,i4,9x,i3,8x,i4)' )itypat0,pspso0,nekb 4395 read(unddb,*) 4396 do iekb=1,nekb 4397 do jekb=1,nekb,4 4398 read(unddb,*) 4399 end do 4400 end do 4401 end do 4402 else 4403 do itypat=1,ntypat 4404 read(unddb, '(12x,i4,12x,i3,12x,i5)' )itypat0,basis_size0,lmn_size0 4405 lmnmax=max(lmnmax,lmn_size0) 4406 read(unddb,*) 4407 read(unddb,*) 4408 read(unddb,'(24x,i3)') nekb 4409 read(unddb,*) 4410 do iekb=1,nekb,4 4411 read(unddb,*) 4412 end do 4413 end do 4414 end if 4415 4416 else if(string==' Description')then 4417 if (usepaw==1) then 4418 ABI_BUG('old DDB pspformat not compatible with PAW 1') 4419 end if 4420 4421 read (unddb, '(10x,i3,10x,i3)' )mproj,mpsang 4422 dimekb=mproj*mpsang 4423 usepaw=0 4424 do itypat=1,ntypat 4425 read (unddb,*) 4426 ! For f-electrons, one more line has been written 4427 do iproj=1,mproj*max(1,(mpsang+2)/3) 4428 read (unddb,*) 4429 end do 4430 end do 4431 4432 else if(string==' No informat')then 4433 4434 dimekb=0 4435 lmnmax=0 4436 !usepaw=0 ! GA: usepaw is also declared earlier in the header 4437 ! and it is that earlier value that usepaw will 4438 ! be compared to in ioddb8_in, so there is no reason 4439 ! to override the value here. 4440 4441 else 4442 write(message, '(a,a,a,a)' )& 4443 & 'Error when reading the psp information',ch10,& 4444 & 'String=',trim(string) 4445 ABI_BUG(message) 4446 end if 4447 4448 !Now, the number of blocks 4449 read(unddb,*) 4450 read(unddb,*) 4451 read(unddb, '(24x,i4)' )nblok 4452 4453 !Now, the type of each blok, in turn 4454 ! GA: Certain types of block are not expected to mix (3 and 5) 4455 mblktyp=1 4456 if(nblok>=1)then 4457 do iblok=1,nblok 4458 4459 read(unddb,*) 4460 read(unddb, '(a32,12x,i12)' )blkname,nelmts 4461 if(blkname==' 2nd derivatives (non-stat.) - ' .or. blkname==' 2rd derivatives (non-stat.) - ')then 4462 blktyp=BLKTYP_d2E_ns 4463 else if(blkname==' 2nd derivatives (stationary) - ' .or. blkname==' 2rd derivatives (stationary) - ')then 4464 blktyp=BLKTYP_d2E_st 4465 else if(blkname==' 2nd derivatives (MBC) - ') then 4466 ! wouldn't it be better to use conditions like: 4467 ! else if(blkname==DESCR_d2E_mbc) then 4468 blktyp=BLKTYP_d2E_mbc 4469 else if(blkname==' 3rd derivatives - ')then 4470 blktyp=BLKTYP_d3E_xx 4471 else if(blkname==' Total energy - ')then 4472 blktyp=BLKTYP_d0E_xx 4473 else if(blkname==' 1st derivatives - ')then 4474 blktyp=BLKTYP_d1E_xx 4475 else if(blkname==' 2nd eigenvalue derivatives - ' .or. blkname==' 2rd eigenvalue derivatives - ')then 4476 blktyp=BLKTYP_d2eig_re 4477 else if(blkname==' 3rd derivatives (long wave) - ')then 4478 blktyp=BLKTYP_d3E_lw 4479 else 4480 write(message, '(a,a,a,a,a,a,a,a,a)' )& 4481 & 'The following string appears in the DDB in place of',' the block type description :',ch10,blkname,ch10,& 4482 & 'Action: check your DDB.',ch10,& 4483 & 'Note: If you did use an abinit version prior to 6.12 to generate your DDB',& 4484 & 'pay attention to the change:: 2rd derivatives ==> 2nd derivatives' 4485 ABI_ERROR(message) 4486 end if 4487 4488 if (is_type_d2E(blktyp)) then 4489 ! Read the phonon wavevector 4490 read(unddb,*) 4491 else if (is_type_d3E(blktyp)) then 4492 ! Read the perturbation wavevectors 4493 read(unddb,*) 4494 read(unddb,*) 4495 read(unddb,*) 4496 mblktyp=blktyp 4497 else if (is_type_d2eig(blktyp)) then 4498 read(unddb,*) 4499 mblktyp=BLKTYP_d2eig_re 4500 end if 4501 4502 ! Read every element 4503 if (is_type_d2eig(blktyp)) then 4504 do ikpt=1,nkpt 4505 read(unddb,*) 4506 do iband=1,nband(ikpt) 4507 read(unddb,*) 4508 do ii=1,nelmts 4509 read(unddb,*) 4510 end do 4511 end do 4512 end do 4513 else 4514 do ii=1,nelmts 4515 read(unddb,*) 4516 end do 4517 end if 4518 4519 end do 4520 end if 4521 4522 ABI_FREE(nband) 4523 4524 !Close the DDB 4525 close(unddb) 4526 4527 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
OUTPUT
acell(3)=length scales of primitive translations (bohr) amu(mtypat)=mass of the atoms (atomic mass unit) ddb_version=version of the ddb file 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
SOURCE
3368 subroutine ioddb8_in(filename,matom,mband,mkpt,msym,mtypat,unddb,& 3369 & acell,amu,ddb_version,dilatmx,ecut,ecutsm,intxc,iscf,ixc,kpt,kptnrm,& 3370 & natom,nband,ngfft,nkpt,nspden,nspinor,nsppol,nsym,ntypat,occ,occopt,& 3371 & pawecutdg,rprim,dfpt_sciss,spinat,symafm,symrel,tnons,tolwfr,tphysel,tsmear,& 3372 & typat,usepaw,wtk,xred,zion,znucl) 3373 3374 !Arguments ------------------------------- 3375 !scalars 3376 integer,intent(in) :: matom,mband,mkpt,msym,mtypat,unddb,usepaw 3377 integer,intent(out) :: ddb_version,intxc,iscf,ixc,natom,nkpt,nspden,nspinor,nsppol,nsym,ntypat,occopt 3378 real(dp),intent(out) :: dilatmx,ecut,ecutsm,pawecutdg,kptnrm,dfpt_sciss,tolwfr,tphysel,tsmear 3379 character(len=*),intent(in) :: filename 3380 !arrays 3381 integer,intent(out) :: nband(mkpt),ngfft(18),symafm(msym),symrel(3,3,msym),typat(matom) 3382 real(dp),intent(out) :: acell(3),amu(mtypat),kpt(3,mkpt),occ(mband*mkpt) 3383 real(dp),intent(out) :: rprim(3,3),spinat(3,matom),tnons(3,msym),wtk(mkpt) 3384 real(dp),intent(out) :: xred(3,matom),zion(mtypat),znucl(mtypat) 3385 3386 !Local variables ------------------------- 3387 !Set routine version number here: 3388 !scalars 3389 integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527 ! should I modify this ? I guess not 3390 integer,parameter :: cvrsio9=20230401,cvrsio8=20100401,cvrsio8_old=20010929,cvrsio8_old_old=19990527 ! should I modify this ? I guess not 3391 integer :: bantot,ddbvrs,iband,ii,ij,ikpt,iline,im,ndig,usepaw0 3392 logical :: ddbvrs_is_current_or_old,testn,testv 3393 character(len=500) :: message 3394 character(len=6) :: name_old, ddbvrs6 3395 character(len=8) :: ddbvrs8 3396 character(len=3) :: prefix 3397 !arrays 3398 character(len=12) :: name(9) 3399 3400 ! ********************************************************************* 3401 3402 !Open the input derivative database. 3403 !write(message,'(a,a)')' About to open file ',TRIM(filename) 3404 !call wrtout(std_out,message,'COLL') 3405 if (open_file(filename,message,unit=unddb,form="formatted",status="old",action="read") /= 0) then 3406 ABI_ERROR(message) 3407 end if 3408 3409 !Check the compatibility of the input DDB with the DDB code 3410 read (unddb,*) 3411 read (unddb,*) 3412 read (unddb, '(20x,i10)' )ddbvrs 3413 3414 !write(std_out,'(a,i10)')' ddbvrs=',ddbvrs 3415 if(ddbvrs/=cvrsio9 .and. ddbvrs/=vrsio8 .and. ddbvrs/=vrsio8_old .and. ddbvrs/=vrsio8_old_old)then 3416 write(message, '(a,i10,2a,4(a,i10),a)' )& 3417 'The input DDB version number=',ddbvrs,' does not agree',ch10,& 3418 'with the allowed code DDB version numbers,',cvrsio9,', ',vrsio8,', ',vrsio8_old,' and ',vrsio8_old_old,' .' 3419 ABI_BUG(message) 3420 end if 3421 3422 !Convert older version to 8 digit format 3423 if (ddbvrs /= cvrsio9) then 3424 ndig= int(log10(real(ddbvrs))) + 1 3425 write(ddbvrs6,'(i0)') ddbvrs 3426 if (ddbvrs==vrsio8 .or.ddbvrs==vrsio8_old) then 3427 if (ndig==6) then 3428 write(prefix,'(i2)') 20 3429 else if (ndig==5) then 3430 write(prefix,'(i3)') 200 3431 end if 3432 else if (ddbvrs==vrsio8_old_old) then 3433 if (ndig==6) then 3434 write(prefix,'(i2)') 19 3435 else if (ndig==5) then 3436 write(prefix,'(i3)') 199 3437 end if 3438 end if 3439 ddbvrs8= trim(prefix) // trim(ddbvrs6) 3440 read(ddbvrs8,'(i8)') ddbvrs 3441 end if 3442 ddb_version=ddbvrs 3443 3444 !Read the 4 n-integers, also testing the names of data, and checking that their value is acceptable. 3445 !This is important to insure that any array has a sufficient dimension. 3446 read (unddb,*) 3447 read (unddb,*) 3448 read (unddb,*) 3449 testn=.true.; testv=.true. 3450 ! ddbvrs_is_current_or_old=(ddbvrs==vrsio8.or.ddbvrs==vrsio8_old) 3451 ddbvrs_is_current_or_old=(ddbvrs>=cvrsio8_old) 3452 3453 !1. usepaw 3454 if(ddbvrs>=cvrsio8)then 3455 read (unddb, '(1x,a9,i10)' )name(1),usepaw0 3456 else 3457 usepaw0=0;name(1)=' usepaw' 3458 end if 3459 if(name(1)/=' usepaw')testn=.false. 3460 if(usepaw0/=usepaw)testv=.false. 3461 !2. natom 3462 if(ddbvrs_is_current_or_old)then 3463 read (unddb, '(1x,a9,i10)' )name(2),natom 3464 else 3465 read (unddb, '(1x,a6,i10)' )name_old,natom ; name(2)=' '//name_old 3466 end if 3467 if(name(2)/=' natom')testn=.false. 3468 if(natom<=0.or.natom>matom)testv=.false. 3469 !3. nkpt 3470 if(ddbvrs_is_current_or_old)then 3471 read (unddb, '(1x,a9,i10)' )name(3),nkpt 3472 else 3473 read (unddb, '(1x,a6,i10)' )name_old,nkpt ; name(3)=' '//name_old 3474 end if 3475 if(name(3)/=' nkpt')testn=.false. 3476 if(nkpt <=0.or.nkpt >mkpt )testv=.false. 3477 !4. nsppol 3478 if(ddbvrs_is_current_or_old)then 3479 read (unddb, '(1x,a9,i10)' )name(4),nsppol 3480 else 3481 read (unddb, '(1x,a6,i10)' )name_old,nsppol ; name(4)=' '//name_old 3482 end if 3483 if(name(4)/=' nsppol')testn=.false. 3484 if(nsppol<=0.or.nsppol>2)testv=.false. 3485 !5. nsym 3486 if(ddbvrs_is_current_or_old)then 3487 read (unddb, '(1x,a9,i10)' )name(5),nsym 3488 else 3489 read (unddb, '(1x,a6,i10)' )name_old,nsym ; name(5)=' '//name_old 3490 end if 3491 if(name(5)/=' nsym')testn=.false. 3492 if(nsym <=0.or.nsym >msym )testv=.false. 3493 !6. ntypat 3494 if(ddbvrs_is_current_or_old)then 3495 read (unddb, '(1x,a9,i10)' )name(6),ntypat 3496 else 3497 read (unddb, '(1x,a6,i10)' )name_old,ntypat ; name(6)=' '//name_old 3498 end if 3499 if(name(6)/=' ntypat' .and. name(6)/=' ntype')testn=.false. 3500 if(ntypat<=0.or.ntypat>mtypat)testv=.false. 3501 !7. occopt 3502 !Before reading nband, the last parameters that define 3503 !the dimension of some array, need to know what is their 3504 !representation, given by occopt 3505 if(ddbvrs_is_current_or_old)then 3506 read (unddb, '(1x,a9,i10)' )name(7),occopt 3507 else 3508 read (unddb, '(1x,a6,i10)' )name_old,occopt ; name(7)=' '//name_old 3509 end if 3510 if(name(7)/=' occopt')testn=.false. 3511 if(occopt<0.or.occopt>8)testv=.false. 3512 !Message if the names or values are not right 3513 if (.not.testn.or..not.testv) then 3514 write(message, '(a,a,a)' )' ioddb8_in : An error has been found in one',ch10,& 3515 ' of the positive n-integers contained in the DDB : ' 3516 call wrtout(std_out,message,'COLL') 3517 write(message, '(a)' )& 3518 ' Expected Found ' 3519 call wrtout(std_out,message,'COLL') 3520 write(message, '(a,i10,a,a,a,i10)' )& 3521 ' usepaw equal to ',usepaw,' ',trim(name(1)),' =',usepaw0 3522 call wrtout(std_out,message,'COLL') 3523 write(message, '(a,i10,a,a,a,i10)' )& 3524 ' natom , lower than',matom+1,' ',trim(name(2)),' =',natom 3525 call wrtout(std_out,message,'COLL') 3526 write(message, '(a,i10,a,a,a,i10)' )& 3527 ' nkpt , lower than',mkpt+1 ,' ',trim(name(3)),' =',nkpt 3528 call wrtout(std_out,message,'COLL') 3529 write(message, '(a,i10,a,a,a,i10)' )& 3530 ' nsppol, lower than',3 ,' ',trim(name(4)),' =',nsppol 3531 call wrtout(std_out,message,'COLL') 3532 write(message, '(a,i10,a,a,a,i10)' )& 3533 ' nsym , lower than',msym+1 ,' ',trim(name(5)),' =',nsym 3534 call wrtout(std_out,message,'COLL') 3535 write(message, '(a,i10,a,a,a,i10)' )& 3536 ' ntypat, lower than',mtypat+1,' ',trim(name(6)),' =',ntypat 3537 call wrtout(std_out,message,'COLL') 3538 write(message, '(a,a,a,i10)' )& 3539 ' occopt, between 0 and 7 ',trim(name(7)),' =',occopt 3540 call wrtout(std_out,message,'COLL') 3541 3542 ABI_ERROR('See the error message above.') 3543 end if 3544 3545 !One more set of parameters define the dimensions of the 3546 !array : nband. Morever, it depends on occopt and nkpt, and has to be 3547 !tested after the test on nkpt is performed. 3548 !8. nband 3549 if(occopt==2)then 3550 im=12 3551 do iline=1,(nkpt+11)/12 3552 name(:) = "" ! reset all line strings 3553 if(iline==(nkpt+11)/12)im=nkpt-12*(iline-1) 3554 if(ddbvrs_is_current_or_old)then 3555 read (unddb, '(1x,a9,5x,12i5)' )name(1),(nband((iline-1)*12+ii),ii=1,im) 3556 else 3557 read (unddb, '(1x,a6,5x,12i5)' )name_old,(nband((iline-1)*12+ii),ii=1,im) ; name(1)=' '//name_old 3558 end if 3559 if (iline==1) then 3560 call ddb_chkname(name(1),' nband') 3561 else 3562 call ddb_chkname(name(1),' ') 3563 end if 3564 end do 3565 else 3566 name(:) = "" ! reset all line strings 3567 if(ddbvrs_is_current_or_old)then 3568 read (unddb, '(1x,a9,i10)' )name(1),nband(1) 3569 else 3570 read (unddb, '(1x,a6,i10)' )name_old,nband(1) ; name(1)=' '//name_old 3571 end if 3572 call ddb_chkname(name(1),' nband') 3573 if(nkpt>1)then 3574 do ikpt=2,nkpt 3575 nband(ikpt)=nband(1) 3576 end do 3577 end if 3578 end if 3579 3580 !check all nband values, and sum them 3581 bantot=0 3582 do ikpt=1,nkpt 3583 if(nband(ikpt)<0)then 3584 write(message, '(a,i4,a,i4,3a)' )& 3585 & 'For ikpt = ',ikpt,' nband = ',nband(ikpt),' is negative.',ch10,& 3586 & 'Action: correct your DDB.' 3587 ABI_ERROR(message) 3588 else if(nband(ikpt)>mband)then 3589 write(message, '(a,i4,a,i4,a,a,i4,3a)' )& 3590 & 'For ikpt = ',ikpt,', nband = ',nband(ikpt),ch10,& 3591 & 'is larger than mband = ',mband,'.',ch10,& 3592 & 'Action: recompile the calling code with a larger mband.' 3593 ABI_ERROR(message) 3594 end if 3595 bantot=bantot+nband(ikpt) 3596 end do 3597 3598 !Read the rest of variables, with check of the names 3599 !9. acell 3600 if(ddbvrs_is_current_or_old)then 3601 read (unddb, '(1x,a9,3d22.14)' )name(1),acell 3602 else 3603 read (unddb, '(1x,a6,3d22.14)' )name_old,acell ; name(1)=' '//name_old 3604 end if 3605 call ddb_chkname(name(1),' acell') 3606 !9. amu 3607 im=3 3608 do iline=1,(ntypat+2)/3 3609 if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1) 3610 if(ddbvrs_is_current_or_old)then 3611 read (unddb, '(1x,a9,3d22.14)' )name(1),(amu((iline-1)*3+ii),ii=1,im) 3612 else 3613 read (unddb, '(1x,a6,3d22.14)' )name_old,(amu((iline-1)*3+ii),ii=1,im) ; name(1)=' '//name_old 3614 end if 3615 if (iline==1) then 3616 call ddb_chkname(name(1),' amu') 3617 else 3618 call ddb_chkname(name(1),' ') 3619 end if 3620 end do 3621 !11. dilatmx 3622 if(ddbvrs_is_current_or_old)then 3623 read (unddb, '(1x,a9,d22.14)' )name(1),dilatmx 3624 call ddb_chkname(name(1),' dilatmx') 3625 else 3626 dilatmx=one 3627 end if 3628 !12. ecut 3629 if(ddbvrs_is_current_or_old)then 3630 read (unddb, '(1x,a9,d22.14)' )name(1),ecut 3631 else 3632 read (unddb, '(1x,a6,d22.14)' )name_old,ecut ; name(1)=' '//name_old 3633 end if 3634 call ddb_chkname(name(1),' ecut') 3635 !12b. pawecutdg (PAW only) 3636 if(ddbvrs>=cvrsio8.and.usepaw==1) then 3637 read (unddb, '(1x,a9,d22.14)' )name(1),pawecutdg 3638 else 3639 pawecutdg=ecut;name(1)='pawecutdg' 3640 end if 3641 call ddb_chkname(name(1),'pawecutdg') 3642 !13. ecutsm 3643 if(ddbvrs_is_current_or_old)then 3644 read (unddb, '(1x,a9,d22.14)' )name(1),ecutsm 3645 call ddb_chkname(name(1),' ecutsm') 3646 else 3647 ecutsm=zero 3648 end if 3649 !14. intxc 3650 if(ddbvrs_is_current_or_old)then 3651 read (unddb, '(1x,a9,i10)' )name(1),intxc 3652 call ddb_chkname(name(1),' intxc') 3653 else 3654 intxc=1 3655 end if 3656 !15. iscf 3657 if(ddbvrs_is_current_or_old)then 3658 read (unddb, '(1x,a9,i10)' )name(1),iscf 3659 else 3660 read (unddb, '(1x,a6,i10)' )name_old,iscf ; name(1)=' '//name_old 3661 end if 3662 call ddb_chkname(name(1),' iscf') 3663 !16. ixc 3664 if(ddbvrs_is_current_or_old)then 3665 read (unddb, '(1x,a9,i10)' )name(1),ixc 3666 else 3667 read (unddb, '(1x,a6,i10)' )name_old,ixc ; name(1)=' '//name_old 3668 end if 3669 call ddb_chkname(name(1),' ixc') 3670 !17. kpt 3671 do iline=1,nkpt 3672 if(ddbvrs_is_current_or_old)then 3673 read (unddb, '(1x,a9,3d22.14)' )name(1),(kpt(ii,iline),ii=1,3) 3674 else 3675 read (unddb, '(1x,a6,3d22.14)' )name_old,(kpt(ii,iline),ii=1,3) ; name(1)=' '//name_old 3676 end if 3677 if (iline==1) then 3678 call ddb_chkname(name(1),' kpt') 3679 else 3680 call ddb_chkname(name(1),' ') 3681 end if 3682 end do 3683 !18. kptnrm 3684 if(ddbvrs_is_current_or_old)then 3685 read (unddb, '(1x,a9,d22.14)' )name(1),kptnrm 3686 else 3687 read (unddb, '(1x,a6,d22.14)' )name_old,kptnrm ; name(1)=' '//name_old 3688 end if 3689 call ddb_chkname(name(1),' kptnrm') 3690 !19. ngfft 3691 if(ddbvrs_is_current_or_old)then 3692 read (unddb, '(1x,a9,5x,3i5)' )name(1),ngfft(1:3) 3693 else 3694 read (unddb, '(1x,a6,5x,3i5)' )name_old,ngfft(1:3) ; name(1)=' '//name_old 3695 end if 3696 !For the time being, do not check the validity of the name, 3697 !in order to accept both ng and ngfft 3698 !20. nspden 3699 if(ddbvrs_is_current_or_old)then 3700 read (unddb, '(1x,a9,i10)' )name(1),nspden 3701 call ddb_chkname(name(1),' nspden') 3702 else 3703 nspden=0 3704 end if 3705 !21. nspinor 3706 if(ddbvrs_is_current_or_old)then 3707 read (unddb, '(1x,a9,i10)' )name(1),nspinor 3708 call ddb_chkname(name(1),' nspinor') 3709 else 3710 nspinor=0 3711 end if 3712 !22. occ 3713 if(occopt==2)then 3714 im=3 3715 do iline=1,(bantot+2)/3 3716 if(iline==(bantot+2)/3)im=bantot-3*(iline-1) 3717 if(ddbvrs_is_current_or_old)then 3718 read (unddb, '(1x,a9,3d22.14)' )name(1),(occ((iline-1)*3+ii),ii=1,im) 3719 else 3720 read (unddb, '(1x,a6,3d22.14)' )name_old,(occ((iline-1)*3+ii),ii=1,im) ; name(1)=' '//name_old 3721 end if 3722 if (iline==1) then 3723 call ddb_chkname(name(1),' occ') 3724 else 3725 call ddb_chkname(name(1),' ') 3726 end if 3727 end do 3728 else 3729 im=3 3730 do iline=1,(nband(1)+2)/3 3731 if(iline==(nband(1)+2)/3)im=nband(1)-3*(iline-1) 3732 if(ddbvrs_is_current_or_old)then 3733 read (unddb, '(1x,a9,3d22.14)' )name(1),(occ((iline-1)*3+ii),ii=1,im) 3734 else 3735 read (unddb, '(1x,a6,3d22.14)' )name_old,(occ((iline-1)*3+ii),ii=1,im) ; name(1)=' '//name_old 3736 end if 3737 if (iline==1) then 3738 call ddb_chkname(name(1),' occ') 3739 else 3740 call ddb_chkname(name(1),' ') 3741 end if 3742 end do 3743 if(nkpt>1)then 3744 do ikpt=2,nkpt 3745 do iband=1,nband(1) 3746 occ(iband+nband(1)*(ikpt-1))=occ(iband) 3747 end do 3748 end do 3749 end if 3750 end if 3751 !23. rprim 3752 do iline=1,3 3753 if(ddbvrs_is_current_or_old)then 3754 read (unddb, '(1x,a9,3d22.14)' )name(1),(rprim(ii,iline),ii=1,3) 3755 else 3756 read (unddb, '(1x,a6,3d22.14)' )name_old,(rprim(ii,iline),ii=1,3) ; name(1)=' '//name_old 3757 end if 3758 if (iline==1) then 3759 call ddb_chkname(name(1),' rprim') 3760 else 3761 call ddb_chkname(name(1),' ') 3762 end if 3763 end do 3764 !24. dfpt_sciss 3765 if(ddbvrs_is_current_or_old)then 3766 read (unddb, '(1x,a11,d22.14)' )name(1),dfpt_sciss 3767 call ddb_chkname(name(1),'dfpt_sciss', 'sciss') 3768 else 3769 read (unddb, '(1x,a6,d22.14)' )name_old,dfpt_sciss ; name(1)=name_old 3770 call ddb_chkname(name(1),'sciss') 3771 end if 3772 !25. spinat 3773 if(ddbvrs_is_current_or_old)then 3774 do iline=1,natom 3775 read (unddb, '(1x,a9,3d22.14)' )name(1),(spinat(ii,iline),ii=1,3) 3776 if (iline==1) then 3777 call ddb_chkname(name(1),' spinat') 3778 else 3779 call ddb_chkname(name(1),' ') 3780 end if 3781 end do 3782 else 3783 ! spinat is set to zero by default in mrgddb.f 3784 ! spinat(:,1:natom)=zero 3785 end if 3786 !26. symafm 3787 if(ddbvrs_is_current_or_old)then 3788 im=12 3789 do iline=1,(nsym+11)/12 3790 if(iline==(nsym+11)/12)im=nsym-12*(iline-1) 3791 read (unddb, '(1x,a9,5x,12i5)' )name(1),(symafm((iline-1)*12+ii),ii=1,im) 3792 if (iline==1) then 3793 call ddb_chkname(name(1),' symafm') 3794 else 3795 call ddb_chkname(name(1),' ') 3796 end if 3797 end do 3798 else 3799 ! symafm is set to 1 by default in mrgddb.f 3800 ! symafm(1:nsym)=1 3801 end if 3802 !27. symrel 3803 do iline=1,nsym 3804 if(ddbvrs_is_current_or_old)then 3805 read (unddb, '(1x,a9,5x,9i5)' )name(1),((symrel(ii,ij,iline),ii=1,3),ij=1,3) 3806 else 3807 read (unddb, '(1x,a6,5x,9i5)' )name_old,& 3808 & ((symrel(ii,ij,iline),ii=1,3),ij=1,3) ; name(1)=' '//name_old 3809 end if 3810 if (iline==1) then 3811 call ddb_chkname(name(1),' symrel') 3812 else 3813 call ddb_chkname(name(1),' ') 3814 end if 3815 end do 3816 !28old. xred 3817 if(.not.ddbvrs_is_current_or_old)then 3818 do iline=1,natom 3819 read (unddb, '(1x,a6,3d22.14)' )name(1),(xred(ii,iline),ii=1,3) 3820 end do 3821 ! No check of name, to allow the old tn 3822 end if 3823 !28. tnons 3824 do iline=1,nsym 3825 if(ddbvrs_is_current_or_old)then 3826 read (unddb, '(1x,a9,3d22.14)' )name(1),(tnons(ii,iline),ii=1,3) 3827 else 3828 read (unddb, '(1x,a6,3d22.14)' )name_old,(tnons(ii,iline),ii=1,3) ; name(1)=' '//name_old 3829 end if 3830 if (iline==1) then 3831 call ddb_chkname(name(1),' tnons') 3832 else 3833 call ddb_chkname(name(1),' ') 3834 end if 3835 end do 3836 !29. tolwfr 3837 if(ddbvrs_is_current_or_old)then 3838 read (unddb, '(1x,a9,d22.14)' )name(1),tolwfr 3839 end if 3840 !Do not check the name, in order to allow both tolwfr and wftol 3841 !30. tphysel 3842 if(ddbvrs_is_current_or_old)then 3843 read (unddb, '(1x,a9,d22.14)' )name(1),tphysel 3844 call ddb_chkname(name(1),' tphysel') 3845 else 3846 tphysel=zero 3847 end if 3848 !31. tsmear 3849 if(ddbvrs_is_current_or_old)then 3850 read (unddb, '(1x,a9,d22.14)' )name(1),tsmear 3851 call ddb_chkname(name(1),' tsmear') 3852 else 3853 tsmear=zero 3854 end if 3855 !32. typat 3856 im=12 3857 do iline=1,(natom+11)/12 3858 if(iline==(natom+11)/12)im=natom-12*(iline-1) 3859 if(ddbvrs_is_current_or_old)then 3860 read (unddb, '(1x,a9,5x,12i5)' )name(1),(typat((iline-1)*12+ii),ii=1,im) 3861 else 3862 read (unddb, '(1x,a6,5x,12i5)' )name_old,(typat((iline-1)*12+ii),ii=1,im) ; name(1)=' '//name_old 3863 end if 3864 if (iline==1) then 3865 ! Both type and typat are allowed => no check 3866 ! call ddb_chkname(name(1),' typat') 3867 else 3868 call ddb_chkname(name(1),' ') 3869 end if 3870 end do 3871 !33old. tolwfr 3872 if(.not.ddbvrs_is_current_or_old)then 3873 read (unddb, '(1x,a6,d22.14)' )name(1),tolwfr 3874 end if 3875 !Do not check the name, in order to allow both tolwfr and wftol 3876 !33. wtk 3877 im=3 3878 do iline=1,(nkpt+2)/3 3879 if(iline==(nkpt+2)/3)im=nkpt-3*(iline-1) 3880 if(ddbvrs_is_current_or_old)then 3881 read (unddb, '(1x,a9,3d22.14)' )name(1),(wtk((iline-1)*3+ii),ii=1,im) 3882 else 3883 read (unddb, '(1x,a6,3d22.14)' )name_old,(wtk((iline-1)*3+ii),ii=1,im) ; name(1)=' '//name_old 3884 end if 3885 if (iline==1) then 3886 call ddb_chkname(name(1),' wtk') 3887 else 3888 call ddb_chkname(name(1),' ') 3889 end if 3890 end do 3891 !34. xred 3892 if(ddbvrs_is_current_or_old)then 3893 do iline=1,natom 3894 read (unddb, '(1x,a9,3d22.14)' )name(1),(xred(ii,iline),ii=1,3) 3895 if (iline==1) then 3896 call ddb_chkname(name(1),' xred') 3897 else 3898 call ddb_chkname(name(1),' ') 3899 end if 3900 end do 3901 end if 3902 !35. znucl 3903 if(ddbvrs_is_current_or_old)then 3904 im=3 3905 do iline=1,(ntypat+2)/3 3906 if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1) 3907 read (unddb, '(1x,a9,3d22.14)' )name(1),(znucl((iline-1)*3+ii),ii=1,im) 3908 if (iline==1) then 3909 call ddb_chkname(name(1),' znucl') 3910 else 3911 call ddb_chkname(name(1),' ') 3912 end if 3913 end do 3914 else 3915 ! znucl is set to zero by default in mrgddb.f 3916 ! znucl(:)=zero 3917 end if 3918 !36. zion 3919 im=3 3920 do iline=1,(ntypat+2)/3 3921 if(iline==(ntypat+2)/3)im=ntypat-3*(iline-1) 3922 if(ddbvrs_is_current_or_old)then 3923 read (unddb, '(1x,a9,3d22.14)' )name(1),(zion((iline-1)*3+ii),ii=1,im) 3924 else 3925 read (unddb, '(1x,a6,3d22.14)' )name_old,(zion((iline-1)*3+ii),ii=1,im) ; name(1)=' '//name_old 3926 end if 3927 if (iline==1) then 3928 ! Do not check the names, to allow both zion and znucl - the latter for 990527 format 3929 ! call ddb_chkname(name(1),' zion') 3930 else 3931 call ddb_chkname(name(1),' ') 3932 end if 3933 end do 3934 3935 end subroutine ioddb8_in
m_ddb_hdr/is_type_d0E [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
is_type_d0E
FUNCTION
Is this block type of a d0E kind?
SOURCE
5023 logical function is_type_d0E(blktyp) result(answer) 5024 5025 !Arguments ------------------------------------ 5026 integer, intent(in) :: blktyp 5027 5028 if (blktyp==BLKTYP_d0E_xx) then 5029 answer = .True. 5030 else 5031 answer = .False. 5032 end if 5033 5034 end function is_type_d0E
m_ddb_hdr/is_type_d1E [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
is_type_d1E
FUNCTION
Is this block type of a d1E kind?
SOURCE
5048 logical function is_type_d1E(blktyp) result(answer) 5049 5050 !Arguments ------------------------------------ 5051 integer, intent(in) :: blktyp 5052 5053 if (blktyp==BLKTYP_d1E_xx) then 5054 answer = .True. 5055 else 5056 answer = .False. 5057 end if 5058 5059 end function is_type_d1E
m_ddb_hdr/is_type_d2E [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
is_type_d2E
FUNCTION
Is this block type of a d2E kind?
SOURCE
5073 logical function is_type_d2E(blktyp) result(answer) 5074 5075 !Arguments ------------------------------------ 5076 integer, intent(in) :: blktyp 5077 5078 if ((blktyp==BLKTYP_d2E_ns).or.(blktyp==BLKTYP_d2E_st).or.(blktyp==BLKTYP_d2E_mbc)) then 5079 answer = .True. 5080 else 5081 answer = .False. 5082 end if 5083 5084 end function is_type_d2E
m_ddb_hdr/is_type_d2eig [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
is_type_d2eig
FUNCTION
Is this block type of a d2eig kind?
SOURCE
5123 logical function is_type_d2eig(blktyp) result(answer) 5124 5125 !Arguments ------------------------------------ 5126 integer, intent(in) :: blktyp 5127 5128 if (blktyp==BLKTYP_d2eig_re.or.blktyp==BLKTYP_d2eig_im) then 5129 answer = .True. 5130 else 5131 answer = .False. 5132 end if 5133 5134 end function is_type_d2eig
m_ddb_hdr/is_type_d3E [ Functions ]
[ Top ] [ m_ddb_hdr ] [ Functions ]
NAME
is_type_d3E
FUNCTION
Is this block type of a d3E kind?
SOURCE
5098 logical function is_type_d3E(blktyp) result(answer) 5099 5100 !Arguments ------------------------------------ 5101 integer, intent(in) :: blktyp 5102 5103 if ((blktyp==BLKTYP_d3E_xx).or.(blktyp==BLKTYP_d3E_lw)) then 5104 answer = .True. 5105 else 5106 answer = .False. 5107 end if 5108 5109 end function is_type_d3E
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
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
with_psps=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
SOURCE
2931 subroutine psddb8 (choice,dimekb,ekb,with_psps,indlmn,lmnmax,& 2932 & nblok,ntypat,nunit,pawtab,pspso,usepaw,useylm) 2933 2934 !Arguments ------------------------------- 2935 !scalars 2936 integer,intent(in) :: choice,dimekb,lmnmax,ntypat,nunit,usepaw,useylm 2937 integer,intent(inout) :: with_psps,nblok 2938 !arrays 2939 integer,intent(in) :: pspso(ntypat) 2940 integer,intent(inout) :: indlmn(6,lmnmax,ntypat) 2941 real(dp),intent(inout) :: ekb(dimekb,ntypat) 2942 type(pawtab_type),intent(inout) :: pawtab(ntypat*usepaw) 2943 2944 !Local variables ------------------------- 2945 !Set the version number 2946 !scalars 2947 integer,parameter :: vrsio8=100401,vrsio8_old=010929,vrsio8_old_old=990527 2948 integer :: basis_size0,dimekb0,iekb,ii,ij,il,ilm,ilmn,iln,iln0,im,ios,iproj,iproj0,itypat,itypat0 2949 integer :: jekb,jlmn,jln,lmnmax0,lmn_size0,lmn2_size0,lpsang,nekb,nproj,npsang,pspso0,shape_type0 2950 integer :: usepaw0,vrspsp8 2951 real(dp) :: rpaw0,rshape0 2952 character(len=12) :: string 2953 character(len=500) :: message 2954 !arrays 2955 integer,allocatable :: i1(:),i2(:),nprj(:),orbitals(:) 2956 real(dp),allocatable :: dij0(:),ekb0(:,:) 2957 2958 ! ********************************************************************* 2959 2960 !Check the value of choice 2961 if (choice<=0.or.choice>=3) then 2962 write(message, '(a,a,a,i10,a)' )& 2963 'The permitted values for choice are 1 or 2.',ch10,& 2964 'The calling routine asks ',choice,'.' 2965 ABI_BUG(message) 2966 end if 2967 2968 !================================================================================== 2969 !First option: read psp characteristic from file ================================= 2970 !================================================================================== 2971 if (choice==1) then 2972 2973 read(nunit,*) 2974 read(nunit, '(a12)' )string 2975 with_psps=1 ; if (lmnmax>0) indlmn(:,:,:)=0 2976 2977 ! -------------------------------------------- 2978 ! ----- NEW FORMAT (NCPP+PAW) --------------- 2979 ! -------------------------------------------- 2980 if (string==' Descriptio')then 2981 ! ============================== 2982 ! ======== Common data ========= 2983 ! ============================== 2984 read(nunit, '(32x,i6)' )vrspsp8 2985 if (vrspsp8==vrsio8_old.or.vrspsp8==vrsio8_old_old) then 2986 usepaw0=0 2987 read(nunit, '(10x,i3,14x,i3,11x,i3)', iostat=ios)dimekb0,lmnmax0,usepaw0 2988 if(ios/=0)then 2989 backspace(nunit) 2990 read (nunit, '(10x,i3,14x,i3)' )dimekb0,lmnmax0 2991 usepaw0=0 2992 end if 2993 else if (vrspsp8==vrsio8) then 2994 read(nunit, '(10x,i3)') usepaw0 2995 if (usepaw/=usepaw0) then 2996 write(message, '(a,i1,a,i1,a)' )'usepaw is announced to be ',usepaw,' but read usepaw is ',usepaw0,' !' 2997 ABI_ERROR(message) 2998 end if 2999 if (usepaw==0) then 3000 read (nunit, '(10x,i3,14x,i3)' )dimekb0,lmnmax0 3001 end if 3002 end if 3003 3004 ! ============================== 3005 ! === Norm-conserving psps ===== 3006 ! ============================== 3007 if (usepaw==0) then 3008 ekb(:,:)=zero 3009 ABI_MALLOC(ekb0,(dimekb,dimekb)) 3010 do itypat=1,ntypat 3011 read(nunit, '(13x,i4,9x,i3,8x,i4)' )itypat0,pspso0,nekb 3012 ! Check the compatibility with the main code dimensioning 3013 if(nekb>dimekb)then 3014 write(message, '(a,i8,a,a,a,i3,a)' )& 3015 ' ',nekb,' components of ekb are announced',ch10,'but dimekb=',dimekb,'.' 3016 ABI_BUG(message) 3017 end if 3018 read(nunit,*) 3019 ilmn=0;iproj0=0 3020 do iekb=1,nekb 3021 read(nunit, '(3i6,3x,8d15.7)' ) iln,lpsang,iproj,(ekb0(ii,iekb),ii=1,min(nekb,4)) 3022 if(nekb>4)then 3023 do jekb=5,nekb,4 3024 read(nunit, '(21x,8d15.7)' )(ekb0(ii,iekb),ii=jekb,min(nekb,jekb+3)) 3025 end do 3026 end if 3027 if (lpsang==0.and.iproj>iproj0) iproj0=iproj 3028 if (useylm==1) then 3029 do im=-lpsang,lpsang 3030 ilmn=ilmn+1 3031 indlmn(1,ilmn,itypat)=lpsang 3032 indlmn(2,ilmn,itypat)=im 3033 indlmn(3,ilmn,itypat)=iproj 3034 indlmn(4,ilmn,itypat)=lpsang**2+lpsang+1+im 3035 indlmn(5,ilmn,itypat)=iln 3036 indlmn(6,ilmn,itypat)=1 3037 if (pspso0/=1.and.iln>(nekb-iproj0)/2) indlmn(6,ilmn,itypat)=2 3038 end do 3039 else 3040 ilmn=ilmn+1 3041 indlmn(1,ilmn,itypat)=lpsang 3042 indlmn(2,ilmn,itypat)=lpsang 3043 indlmn(3,ilmn,itypat)=iproj 3044 indlmn(4,ilmn,itypat)=lpsang**2+lpsang+1 3045 indlmn(5,ilmn,itypat)=iln 3046 indlmn(6,ilmn,itypat)=1 3047 if (pspso0/=1.and.iln>(nekb-iproj0)/2) indlmn(6,ilmn,itypat)=2 3048 end if 3049 ! For the time being, only diagonal ekb are treated in abinit v3 3050 ekb(iekb,itypat)=ekb0(iekb,iekb) 3051 ! For non-diagonal ekb, one could use: 3052 ! do jekb=iekb to nekb 3053 ! ekb(jekb+iekb*(iekb-1)/2,itypat)=ekb0(jekb,iekb) 3054 ! end do 3055 end do 3056 end do 3057 ABI_FREE(ekb0) 3058 3059 ! ============================== 3060 ! ============ PAW ============= 3061 ! ============================== 3062 else 3063 do itypat=1,ntypat 3064 read(nunit, '(12x,i4,12x,i3,12x,i5)' )itypat0,basis_size0,lmn_size0 3065 lmn2_size0=lmn_size0*(lmn_size0+1)/2 3066 ABI_MALLOC(orbitals,(basis_size0)) 3067 read(nunit, '(20x,50i2)' ) orbitals(1:basis_size0) 3068 read(nunit, '(11x,f6.3,13x,i2,11x,f6.3)' ) rpaw0,shape_type0,rshape0 3069 read(nunit,'(24x,i3)') nekb 3070 read(nunit,*) 3071 ABI_MALLOC(dij0,(nekb)) 3072 ABI_MALLOC(i1,(nekb)) 3073 ABI_MALLOC(i2,(nekb)) 3074 do ii=1,nekb,4 3075 read(nunit,'(3x,4(1x,i4,1x,i4,1x,d12.5))') (i1(ij),i2(ij),dij0(ij),ij=ii,min(ii+3,nekb)) 3076 end do 3077 if (lmn_size0>lmnmax) then 3078 write(message, '(a,i5,3a,i5,a)' )& 3079 'max. value of ',lmnmax,' for lmn_size is announced',ch10,'but ',lmn_size0,' is read.' 3080 ABI_BUG(message) 3081 end if 3082 if (allocated(pawtab(itypat)%dij0)) then 3083 if (lmn_size0>pawtab(itypat)%lmn_size) then 3084 write(message, '(a,i5,3a,i5,a)' )& 3085 'lmn_size=,',pawtab(itypat)%lmn_size,' is announced',ch10,'but ',lmn_size0,' is read.' 3086 ABI_BUG(message) 3087 end if 3088 end if 3089 ABI_MALLOC(nprj,(0:maxval(orbitals))) 3090 ilmn=0;nprj=0 3091 do iln=1,basis_size0 3092 il=orbitals(iln) 3093 nprj(il)=nprj(il)+1 3094 do ilm=1,2*il+1 3095 indlmn(1,ilmn+ilm,itypat)=il 3096 indlmn(2,ilmn+ilm,itypat)=ilm-(il+1) 3097 indlmn(3,ilmn+ilm,itypat)=nprj(il) 3098 indlmn(4,ilmn+ilm,itypat)=il*il+ilm 3099 indlmn(5,ilmn+ilm,itypat)=iln 3100 indlmn(6,ilmn+ilm,itypat)=1 3101 end do 3102 ilmn=ilmn+2*il+1 3103 end do 3104 pawtab(itypat)%basis_size=basis_size0 3105 pawtab(itypat)%lmn_size =lmn_size0 3106 pawtab(itypat)%lmn2_size =lmn2_size0 3107 pawtab(itypat)%shape_type=shape_type0 3108 pawtab(itypat)%rpaw =rpaw0 3109 pawtab(itypat)%rshp =rshape0 3110 if (.not.allocated(pawtab(itypat)%dij0)) then 3111 ABI_MALLOC(pawtab(itypat)%dij0,(lmn2_size0)) 3112 end if 3113 pawtab(itypat)%dij0(1:lmn2_size0)=zero 3114 do ii=1,nekb 3115 ij=i1(ii)+i2(ii)*(i2(ii)-1)/2 3116 pawtab(itypat)%dij0(ij)=dij0(ii) 3117 end do 3118 ABI_FREE(nprj) 3119 ABI_FREE(orbitals) 3120 ABI_FREE(dij0) 3121 ABI_FREE(i1) 3122 ABI_FREE(i2) 3123 end do 3124 3125 end if ! NCPP or PAW 3126 3127 ! -------------------------------------------- 3128 ! ----- OLD FORMAT (NCPP only) -------------- 3129 ! -------------------------------------------- 3130 else if (string==' Description')then 3131 if (usepaw==1) then 3132 ABI_BUG("old DDB pspformat not compatible with PAW") 3133 end if 3134 3135 read (nunit, '(10x,i3,10x,i3)' )nproj,npsang 3136 nekb=nproj*npsang 3137 ! Check the compatibility with the main code dimensioning 3138 if(nekb>dimekb)then 3139 write(message, '(a,i8,a,a,a,i3,a)' )& 3140 ' ',nekb,' components of ekb are announced',ch10,'but the maximum is dimekb=',dimekb,'.' 3141 ABI_BUG(message) 3142 end if 3143 if(useylm/=0)then 3144 ABI_BUG('useylm must be 0 !') 3145 end if 3146 ! Read the data 3147 ABI_MALLOC(ekb0,(dimekb,dimekb)) 3148 ekb0(:,:)=zero 3149 do itypat=1,ntypat 3150 read (nunit, '(13x,i4)' )ij 3151 do iproj=1,nproj 3152 read (nunit, '(6x,3d22.14)' )(ekb0(iproj+nproj*(ii-1),iproj+nproj*(ii-1)),ii=1,min(npsang,3)) 3153 if(npsang>3)read (nunit, '(6x,3d22.14)' )(ekb0(iproj+nproj*(ii-1),iproj+nproj*(ii-1)),ii=4,npsang) 3154 do ii=1,npsang 3155 iekb=iproj+nproj*(ii-1) 3156 indlmn(1,iekb,itypat)=ii-1 3157 indlmn(2,iekb,itypat)=ii-1 3158 indlmn(3,iekb,itypat)=iproj 3159 indlmn(4,iekb,itypat)=ii**2-ii+1 3160 indlmn(5,iekb,itypat)=iekb 3161 indlmn(6,iekb,itypat)=1 3162 ! For the time being, only diagonal ekb are treated in abinit v3 3163 ekb(iekb,itypat)=ekb0(iekb,iekb) 3164 end do 3165 end do 3166 end do 3167 ABI_FREE(ekb0) 3168 3169 ! -------------------------------------------- 3170 ! ----- OTHER CASES ------------------------- 3171 ! -------------------------------------------- 3172 else if(string==' No informat')then 3173 with_psps=0 3174 else 3175 ABI_BUG('Error when reading the psp information') 3176 end if 3177 3178 ! Now, the number of blocks 3179 read(nunit,*) 3180 read(nunit,*) 3181 read(nunit, '(24x,i4)' )nblok 3182 3183 ! ================================================================================== 3184 ! Second option: write psp characteristic to file ================================== 3185 ! ================================================================================== 3186 else if(choice==2)then 3187 3188 write(nunit, '(a)' )' ' 3189 if (with_psps==0)then 3190 ! This possibility is used when the DDB is initialized, 3191 ! and the ekb s are not available from the GS input file... 3192 write(nunit, '(a)' )' No information on the potentials yet ' 3193 else 3194 3195 ! ============================== 3196 ! === Norm-conserving psps ===== 3197 ! ============================== 3198 ! GA: Something is wrong here. 3199 ! I noticed this with test tutorespfn[tnlo_2] 3200 ! where the ekb are reported both in the main output 3201 ! and in the ddb, but some components written to ddb are wrong. 3202 if (usepaw==0) then 3203 write(nunit, '(a)' )' Description of the potentials (KB energies)' 3204 write(nunit, '(a,i6)' )' vrsio8 (for pseudopotentials)=',vrsio8 3205 write(nunit, '(a,i3)' ) ' usepaw =',usepaw 3206 write(nunit, '(a,i3,a,i3,a,i3)' )' dimekb =',dimekb,' lmnmax=',lmnmax 3207 ABI_MALLOC(ekb0,(dimekb,dimekb)) 3208 do itypat=1,ntypat 3209 ! Compute nekb 3210 nekb=0 3211 do jlmn=1,lmnmax 3212 jln=indlmn(5,jlmn,itypat) 3213 if(jln>nekb)then 3214 nekb=jln 3215 end if 3216 end do 3217 write(nunit, '(a,i4,a,i3,a,i4)' )' Atom type= ',itypat,' pspso=',pspso(itypat),' nekb=',nekb 3218 write(nunit, '(a)' ) ' iln lpsang iproj ekb(:)' 3219 iln0=0 3220 ekb0(:,:)=zero 3221 do ilmn=1,lmnmax 3222 iln =indlmn(5,ilmn,itypat) 3223 if (iln>iln0) then 3224 iln0=iln 3225 lpsang=indlmn(1,ilmn,itypat) 3226 iproj=indlmn(3,ilmn,itypat) 3227 ! For the time being, only diagonal ekb are treated in abinit v3 3228 ekb0(iln,iln)=ekb(iln,itypat) 3229 ! For non-diagonal ekb, one could use: 3230 ! do ii=iln to nekb 3231 ! ekb0(ii,iln)=ekb(ii+iln*(iln-1)/2,itypat) 3232 ! end do 3233 write(nunit, '(3i6,3x,4es15.7)' ) iln,lpsang,iproj,(ekb0(ii,iln),ii=1,min(nekb,4)) 3234 if(nekb>4)then 3235 do iekb=5,nekb,4 3236 write(nunit, '(21x,4es15.7)' )(ekb0(ii,iekb),ii=iekb,min(nekb,iekb+3)) 3237 end do 3238 end if 3239 end if 3240 end do 3241 end do 3242 ABI_FREE(ekb0) 3243 3244 ! ============================== 3245 ! ============ PAW ============= 3246 ! ============================== 3247 else 3248 write(nunit, '(a)' )' Description of the PAW dataset(s)' 3249 write(nunit, '(a,i6)' )' vrsio8 (for pseudopotentials)=',vrsio8 3250 write(nunit, '(a,i3)' ) ' usepaw =',usepaw 3251 do itypat=1,ntypat 3252 iln0=0 3253 ABI_MALLOC(orbitals,(pawtab(itypat)%basis_size)) 3254 do ilmn=1,pawtab(itypat)%lmn_size 3255 iln =indlmn(5,ilmn,itypat) 3256 if (iln>iln0) then 3257 iln0=iln;orbitals(iln)=indlmn(1,ilmn,itypat) 3258 end if 3259 end do 3260 write(nunit, '(a,i4,a,i3,a,i5)' ) & 3261 ' Atom type=',itypat,' basis_size=',pawtab(itypat)%basis_size,& 3262 ' lmn_size=',pawtab(itypat)%lmn_size 3263 write(nunit, '(a,50i2)' ) & 3264 ' Basis functions=',orbitals(1:pawtab(itypat)%basis_size) 3265 write(nunit, '(a,f6.3,a,i2,a,f6.3)' ) & 3266 ' r_PAW= ',pawtab(itypat)%rpaw, & 3267 ' shape_type= ',pawtab(itypat)%shape_type,& 3268 ' r_shape= ',pawtab(itypat)%rshp 3269 nekb=0 3270 ABI_MALLOC(dij0,(pawtab(itypat)%lmn2_size)) 3271 ABI_MALLOC(i1,(pawtab(itypat)%lmn2_size)) 3272 ABI_MALLOC(i2,(pawtab(itypat)%lmn2_size)) 3273 do jlmn=1,pawtab(itypat)%lmn_size 3274 ij=jlmn*(jlmn-1)/2 3275 do ilmn=1,jlmn 3276 if (abs(pawtab(itypat)%dij0(ij+ilmn))>tol16) then 3277 nekb=nekb+1;i1(nekb)=ilmn;i2(nekb)=jlmn 3278 dij0(nekb)=pawtab(itypat)%dij0(ij+ilmn) 3279 end if 3280 end do 3281 end do 3282 write(nunit,'(a,i3,a)') ' Dij0= (only the ',nekb,' values different from zero)' 3283 write(nunit,'(2a)') ' i j Dij0 i j Dij0 ',& 3284 ' i j Dij0 i j Dij0' 3285 do ii=1,nekb,4 3286 write(nunit,'(3x,4(1x,i4,1x,i4,1x,es12.5))') (i1(ij),i2(ij),dij0(ij),ij=ii,min(ii+3,nekb)) 3287 end do 3288 ABI_FREE(dij0) 3289 ABI_FREE(i1) 3290 ABI_FREE(i2) 3291 ABI_FREE(orbitals) 3292 end do 3293 3294 end if ! NCPP or PAW 3295 end if ! with_psps==0 3296 3297 ! Now, write the number of blocks 3298 write(nunit, '(a)' )' ' 3299 write(nunit, '(a)' )' **** Database of total energy derivatives ****' 3300 write(nunit, '(a,i4)' ) ' Number of data blocks= ',nblok 3301 end if 3302 3303 end subroutine psddb8