1! 2! Copyright (C) 2001-2020 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8! 9!----------------------------------------------------------------------- 10SUBROUTINE hp_summary 11 !----------------------------------------------------------------------- 12 ! 13 ! This routine writes on output a summary of the variables which do not 14 ! depend on the q point. 15 ! 16 USE constants, ONLY : rytoev 17 USE ions_base, ONLY : nat, ityp, atm, tau, ntyp => nsp, amass 18 USE io_global, ONLY : stdout 19 USE cell_base, ONLY : at, bg, ibrav, alat, omega, celldm 20 USE gvecs, ONLY : dual 21 USE gvecw, ONLY : ecutwfc 22 USE funct, ONLY : write_dft_name 23 USE control_lr, ONLY : ethr_nscf 24 USE ldaU, ONLY : is_hubbard, Hubbard_U, lda_plus_u_kind, & 25 Hubbard_V, num_uc 26 USE ldaU_hp, ONLY : conv_thr_chi, skip_atom, skip_type, & 27 todo_atom, find_atpert, nath_pert 28 29 IMPLICIT NONE 30 ! 31 INTEGER :: i, ipol, apol, na, nb, nt, isymq, isym, ik, nsymtot, dimn 32 ! generic counter 33 ! counter on polarizations 34 ! counter on polarizations 35 ! counter on atoms 36 ! counter on atomic types 37 ! counter on symmetries 38 ! counter on symmetries 39 ! counter on k points 40 ! 41 WRITE( stdout, * ) 42 ! 43 WRITE( stdout, 100) ibrav, alat, omega, nat, ntyp, & 44 ecutwfc, ecutwfc * dual, ethr_nscf, conv_thr_chi 45100 FORMAT (/5x, 'bravais-lattice index = ',i12,/,5x, & 46 & 'lattice parameter (alat) = ',f12.4,' (a.u.)',/,5x, & 47 & 'unit-cell volume = ',f12.4,' (a.u.)^3',/,5x, & 48 & 'number of atoms/cell = ',i12,/,5x, & 49 & 'number of atomic types = ',i12,/,5x, & 50 & 'kinetic-energy cut-off = ',f12.2,' (Ry)',/,5x, & 51 & 'charge density cut-off = ',f12.2,' (Ry)',/,5x, & 52 & 'conv. thresh. for NSCF = ',3x,1pe9.1,/,5x, & 53 & 'conv. thresh. for chi = ',3x,1pe9.1) 54 ! 55 ! Info about Hubbard parameters 56 ! 57 WRITE (stdout,'(5x,a)') 'Input Hubbard parameters (in eV):' 58 IF (lda_plus_u_kind.EQ.0) THEN 59 DO nt = 1, ntyp 60 IF (is_hubbard(nt)) THEN 61 WRITE (stdout,'(7x,a,i2,a,21x,a,1x,1pe12.5)') & 62 & 'U (',nt,')','= ', Hubbard_U(nt)*rytoev 63 ENDIF 64 ENDDO 65 ELSEIF (lda_plus_u_kind.EQ.2) THEN 66 ! Number of atoms in the 3x3x3 supercell 67 dimn = num_uc * nat 68 DO na = 1, nat 69 DO nb = 1, dimn 70 IF ( ABS(Hubbard_V(na,nb,1)).GE.1.d-15 .OR. nb==na) THEN 71 WRITE (stdout,'(7x,a,i3,a,i4,a,2x,a,1x,f7.4)') & 72 & 'V (',na,',',nb,')','= ', Hubbard_V(na,nb,1)*rytoev 73 ENDIF 74 ENDDO 75 ENDDO 76 ENDIF 77 ! 78 ! Description of the unit cell 79 ! 80 WRITE( stdout, '(/,2(3x,3(2x,"celldm(",i1,") =",f9.5),/))') (i, celldm(i), i=1,6) 81 ! 82 ! Description of the direct and reciprocal lattice vectors 83 ! 84 WRITE( stdout, '(5x, "crystal axes: (cart. coord. in units of alat)",/, & 85 & 3(15x,"a(",i1,") = (",3f8.4," ) ",/ ) )') & 86 & (apol, (at(ipol,apol),ipol=1,3), apol=1,3) 87 WRITE( stdout, '(5x, "reciprocal axes: (cart. coord. in units 2 pi/alat)",/, & 88 & 3(15x,"b(",i1,") = (",3f8.4," ) ",/ ) )') & 89 & (apol, (bg(ipol,apol),ipol=1,3), apol=1,3) 90 ! 91 ! Description of the atoms inside the unit cell 92 ! 93 WRITE( stdout, '(5x,"Atoms inside the unit cell (Cartesian axes):")') 94 WRITE( stdout, '(5x,"site n. atom mass ", " positions (alat units)")') 95 WRITE( stdout, '(6x,i3,3x,a6,3x,f8.4," tau(",i3, ") = (",3f9.5," )")') & 96 & (na, atm(ityp(na)), amass(ityp(na)), na, (tau(ipol,na), ipol=1,3), na=1,nat) 97 ! 98 ! Print info about atoms which will not be perterbed (if requested from the input file) 99 ! 100 IF ( find_atpert.NE.1 .AND. ANY(skip_atom(:)) ) THEN 101 WRITE( stdout, '(/5x,"WARNING: Skipping perturbation of the following atom(s):")') 102 DO na = 1, nat 103 nt = ityp(na) 104 IF (skip_atom(na)) THEN 105 WRITE( stdout, '(7x,i2,3x,a6,3x,f8.4," tau(",i2, ") = (",3f9.5," )")') & 106 & na, atm(nt), amass(nt), na, (tau(ipol,na), ipol=1,3) 107 ENDIF 108 ENDDO 109 ENDIF 110 ! 111 IF (ANY(skip_type(:))) WRITE( stdout, '(/5x,"WARNING: Skipping perturbation of atoms with the following type:")') 112 DO nt = 1, ntyp 113 IF (skip_type(nt)) WRITE( stdout, '(8x,"Atomic type : ", a4)') atm(nt) 114 ENDDO 115 ! 116 IF ( nath_pert > 1 ) THEN 117 WRITE( stdout, '(/5x,"List of",i3,1x,"atoms which will be " & 118 & "perturbed (one at a time):",/)') nath_pert 119 ELSE 120 WRITE( stdout, '(/5x,"Atom which will be perturbed:",/)') 121 ENDIF 122 ! 123 DO na = 1, nat 124 IF (todo_atom(na)) THEN 125 WRITE( stdout, '(7x,i2,3x,a6,3x,f8.4," tau(",i2, ") = (",3f9.5," )")') & 126 & na, atm(ityp(na)), amass(ityp(na)), na, (tau(ipol,na), ipol=1,3) 127 ENDIF 128 ENDDO 129 ! 130 RETURN 131 ! 132END SUBROUTINE hp_summary 133