From bfda49d0396d32ae900d7c7a7a0185269549e297 Mon Sep 17 00:00:00 2001 From: paulatto Date: Tue, 17 Nov 2009 10:51:35 +0000 Subject: [PATCH] An additional file is now generated by NEB calculations; it is called $prefix.crd and contains the path coordinates in pw.x input format. Documentation updated (actually written) accordingly, I hope in the right file. LP git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@6161 c92efa57-630b-4861-b058-cf58834340f0 --- Doc/user_guide.tex | 30 +++++++++++++- Modules/io_files.f90 | 2 + Modules/path_base.f90 | 3 +- Modules/path_io_routines.f90 | 76 +++++++++++++++++++++++++++++++++--- 4 files changed, 103 insertions(+), 8 deletions(-) diff --git a/Doc/user_guide.tex b/Doc/user_guide.tex index b33728315..692e44747 100644 --- a/Doc/user_guide.tex +++ b/Doc/user_guide.tex @@ -1818,6 +1818,32 @@ of the initial and final image of the elastic band have to be specified in the ATOMIC POSITIONS card. A detailed description of all input variables is contained in the file Doc/INPUT PW. See also Example 17. +An NEB calculation will produce a number of files in the current directory +(i.e. in the directory were the code is run) containing additional information +on the minimum-energy path. The files are organized as following, where \$prefix +is the prefix of the calculation (specified in the input file): +\begin{description} +\item[\$prefix.dat] +is a three-column file containig the position of each image on the reaction +coordinate (arb. units), its energy in eV relative to the energy of the first image +and the residual error for the image in eV/$a_0$. +\item[\$prefix.int] +contains an interpolation of the path energy profile that pass exactly through each +image; it is computed using both the image energies and their derivatives +\item[\$prefix.path] +information used by QE to restart a path calculation, its format depends on the input +details and is undocumented +\item[\$prefix.axsf] +atomic positions of all path images in the XCrysDen animation format, +{\texttt xcrysden -\--axsf \$prefix.axsf} to open it +\item[\$prefix.xyz] +atomic positions of all path images in the generic xyz format, used by +many quantum-chemistry softwares +\item[\$prefix.crd] +path information in the input format used by {\texttt pw.x}, suitable for a manual +restart of the calculation +\end{description} + "NEB calculation are a bit tricky in general and require extreme care to be setup correctly. NEB also takes easily hunders of iteration to converge, of course depending on the number of atoms and of images. Here is some @@ -1834,7 +1860,9 @@ Carefully choose the initial path. Remember that QE assumes continuity between the first and the last image at the initial condition. In other words, periodic images are NOT used; you may have to manually translate an atom by one or more unit cell base vectors in order to have a meaningful -initial path. +initial path. You can visualize NEB input files with XCrysDen as animations, +take some time to check if any atoms overlap or get very close in the initial +path (you will have to add intermediate images, in this case). \item Try to start the NEB process with most atomic positions fixed, in order to converge the more "problematic" ones, before leaving diff --git a/Modules/io_files.f90 b/Modules/io_files.f90 index 52b3afb98..56d00b3c3 100644 --- a/Modules/io_files.f90 +++ b/Modules/io_files.f90 @@ -60,6 +60,7 @@ MODULE io_files CHARACTER (LEN=256) :: & dat_file = 'os.dat', &! file containing the enegy profile int_file = 'os.int', &! file containing the interpolated energy profile + crd_file = 'os.crd', &! file containing path coordinates in pw.x input format path_file = 'os.path', &! file containing informations needed to restart a path simulation xyz_file = 'os.xyz', &! file containing coordinates of all images in xyz format axsf_file = 'os.axsf', &! file containing coordinates of all images in axsf format @@ -119,6 +120,7 @@ MODULE io_files INTEGER :: iunxyz = 24 ! unit for saving coordinates ( xyz format ) INTEGER :: iunaxsf = 25 ! unit for saving coordinates ( axsf format ) INTEGER :: iunbroy = 26 ! unit for saving broyden's history + INTEGER :: iuncrd = 27 ! unit for saving coordinates in pw.x input format ! ! ... meta-dynamics ! diff --git a/Modules/path_base.f90 b/Modules/path_base.f90 index ec7b6aa73..23989cf2a 100644 --- a/Modules/path_base.f90 +++ b/Modules/path_base.f90 @@ -64,7 +64,7 @@ MODULE path_base USE control_flags, ONLY : conv_elec, lcoarsegrained USE ions_base, ONLY : nat, amass, ityp, if_pos USE metadyn_vars, ONLY : ncolvar - USE io_files, ONLY : prefix, tmp_dir, path_file, dat_file, & + USE io_files, ONLY : prefix, tmp_dir, path_file, dat_file, crd_file, & int_file, xyz_file, axsf_file, broy_file USE path_variables, ONLY : climbing, pos, istep_path, nstep_path, & dim1, num_of_images, pes, grad_pes, mass, & @@ -86,6 +86,7 @@ MODULE path_base path_file = TRIM( prefix ) // ".path" dat_file = TRIM( prefix ) // ".dat" int_file = TRIM( prefix ) // ".int" + crd_file = TRIM( prefix ) // ".crd" xyz_file = TRIM( prefix ) // ".xyz" axsf_file = TRIM( prefix ) // ".axsf" ! diff --git a/Modules/path_io_routines.f90 b/Modules/path_io_routines.f90 index 2c3dd06e1..173349839 100644 --- a/Modules/path_io_routines.f90 +++ b/Modules/path_io_routines.f90 @@ -650,20 +650,20 @@ MODULE path_io_routines !----------------------------------------------------------------------- ! USE constants, ONLY : pi - USE input_parameters, ONLY : atom_label + USE input_parameters, ONLY : atom_label, atomic_positions USE control_flags, ONLY : lcoarsegrained - USE cell_base, ONLY : alat, at - USE ions_base, ONLY : ityp, nat + USE cell_base, ONLY : alat, at, bg + USE ions_base, ONLY : ityp, nat, if_pos USE path_formats, ONLY : dat_fmt, int_fmt, xyz_fmt, axsf_fmt USE path_variables, ONLY : pos, grad_pes, pes, num_of_images, & tangent, dim1, error - USE io_files, ONLY : iundat, iunint, iunxyz, iunaxsf, & - dat_file, int_file, xyz_file, axsf_file + USE io_files, ONLY : iundat, iunint, iunxyz, iuncrd, iunaxsf, & + dat_file, int_file, xyz_file, axsf_file, crd_file ! IMPLICIT NONE ! REAL(DP) :: r, delta, x - REAL(DP), ALLOCATABLE :: a(:), b(:), c(:), d(:), f(:), s(:) + REAL(DP), ALLOCATABLE :: a(:), b(:), c(:), d(:), f(:), s(:), tau_out(:,:,:) REAL(DP) :: ener, ener_0 INTEGER :: i, j, ia INTEGER, PARAMETER :: max_i = 250 @@ -771,6 +771,70 @@ MODULE path_io_routines ! CLOSE( UNIT = iunxyz ) ! + ! ... the *.xyz file is written here + ! + OPEN( UNIT = iuncrd, FILE = crd_file, STATUS = "UNKNOWN", & + ACTION = "WRITE" ) + ALLOCATE( tau_out(3,nat,num_of_images) ) + ! + DO i = 1, num_of_images + DO ia = 1,nat + tau_out(1,ia,i) = pos(3*ia-2,i) + tau_out(2,ia,i) = pos(3*ia-1,i) + tau_out(3,ia,i) = pos(3*ia-0,i) + ENDDO + ENDDO + ! + SELECT CASE( atomic_positions ) + ! + ! ... convert output atomic positions from internally used format + ! ... (bohr units, for path) to the same format used in input + ! + CASE( 'alat' ) + WRITE( iuncrd, '(/"ATOMIC_POSITIONS (alat)")' ) + tau_out(:,:,:) = tau_out(:,:,:) / alat + CASE( 'bohr' ) + WRITE( iuncrd, '(/"ATOMIC_POSITIONS (bohr)")' ) + CASE( 'crystal' ) + WRITE( iuncrd, '(/"ATOMIC_POSITIONS (crystal)")' ) + tau_out(:,:,:) = tau_out(:,:,:) / alat + DO i = 1, num_of_images + call cryst_to_cart( nat, tau_out(1,1,i), bg, -1 ) + ENDDO + CASE( 'angstrom' ) + WRITE( iuncrd, '(/"ATOMIC_POSITIONS (angstrom)")' ) + tau_out(:,:,:) = tau_out(:,:,:) * bohr_radius_angs + CASE DEFAULT + WRITE( iuncrd, '(/"ATOMIC_POSITIONS")' ) + END SELECT + DO i = 1, num_of_images + ! Ad the image label + IF ( i == 1) THEN + WRITE( UNIT = iuncrd, FMT='("first_image")' ) + ELSE IF ( i == num_of_images) THEN + WRITE( UNIT = iuncrd, FMT='("last_image")' ) + ELSE + WRITE( UNIT = iuncrd, FMT='("intermediate_image", i5)') i-1 + ENDIF + ! + DO ia = 1, nat + ! + IF ( i == 1 .and. ANY(if_pos(:,ia) /= 1) ) THEN + WRITE( UNIT = iuncrd, FMT = '(x,a4,3f18.10,3i2)' ) & + TRIM( atom_label( ityp( ia ) ) ), & + tau_out(1:3,ia,i), if_pos(1:3,ia) + ELSE + WRITE( UNIT = iuncrd, FMT = '(x,a4,3f18.10)' ) & + TRIM( atom_label( ityp( ia ) ) ), & + tau_out(1:3,ia,i) + ENDIF + ! + END DO + ! + END DO + ! + CLOSE( UNIT = iuncrd ) + ! ! ... the *.axsf file is written here ! OPEN( UNIT = iunaxsf, FILE = axsf_file, STATUS = "UNKNOWN", &