[SIESTA-L] Energy correction to charged defect calculations
Dear Users, Suppose I ran a SIESTA calculation with net charge q on a periodic system. The magnitude of the first order energy correction based on the multipole expansion is: delta_E=alpha*q^2/(2*epsilon*L), where alpha is the Madelung constant, epsilon is the dielectric constant and L is the supercell length. My question is: in order to obtain the corrected energy, should delta_E be added or subtract from the total energy in SIESTA? Thanks, FG -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Polarization Grid for Born Charge Calculation of a large slab
Dear Users, I have a large slab with dimensions Lx=Ly=23 Ang and Lz=88 Ang. I relaxed the slab using the Gamma point due to its large size and would like to compute the infrared spectra for the adsorbed species. I would like to seek advice from the developers and experienced users on how to choose an appropriate polarization grid for the Born Charge calculation of such a large slab (which will be used to compute infrared intensities via Vibra). Thank you. Regards, FG -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Constrained vibrational frequency calculation
Dear Users, I have a molecular species adsorbed on the surface of a big slab with 500+ atoms. Is it possible for me to compute the vibrational frequencies of only the surface layer atoms of the slab and the atoms of the adsorbed molecule? Thank you. Regards FG -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Band lines
Dear Users, Kindly forgive me if this question is too basic, but in the band lines below (in units of pi/a) for crystalline silicon, why is the first Gamma point [2 2 2] instead of the traditional [0 0 0]? Translational invariance perhaps? 1 2.000 2.000 2.000 \Gamma 15 2.000 0.000 0.000 X 25 0.000 0.000 0.000 \Gamma 20 1.000 1.000 1.000 L 20 2.000 0.000 0.000 X 15 2.000 1.000 0.000 W 20 1.000 1.000 1.000 L -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Libfdf error in SIESTA 5.0 compilation
Dear SIESTA developers, I am having trouble compiling SIESTA-5.0, specifically regarding libfdf. I get the following error message even though libfdf was present in the External folder: -- Searching for libfdf -- | Siesta_find_package[libfdf] METHODS | ALLOWED = source | cmake;pkgconf;source;fetch -- | source in folder: /tmp/test/SIESTA_5.0_LATEST/siesta-rel-5.0/External/libfdf -- | source in folder: /tmp/test/SIESTA_5.0_LATEST/siesta-rel-5.0/External/libfdf - not found/compatible -- Searching for libfdf - not found CMake Error at Config/cmake/SiestaFindPackage.cmake:383 (message): Required package libfdf cannot be found Call Stack (most recent call first): Config/cmake/Modules/FindCustomlibfdf.cmake:3 (Siesta_find_package) CMakeLists.txt:215 (find_package) I have outlined my installation steps below so that the developers can point out any potential flaws in my approach. Any assistance in fixing this problem would be greatly appreciated. Thank you. STEP 1: Installation of required libraries The latest versions of the following libraries were installed: libfdf, libgridxc, libpsml, libxc and xmlf90 STEP 2: SIESTA-5.0 set-up (a) The latest version of SIESTA-5.0 was downloaded (this is not the beta-5.0 version; it is the most up-to-date version). (b) After expanding the *.tar.gz file, I copied the libraries in STEP 1 into the folder labeled "External" in the main SIESTA directory. STEP 3: Compiling SIESTA with CMake #Define a few important directories MKLPATH=/opt/intel/oneapi/mkl/2023.2.0/lib/intel64 SDIR=/tmp/test/SIESTA_5.0_LATEST/siesta-rel-5.0 LIBXC_ROOT=${SDIR}/External/libxc GRIDXC_ROOT=${SDIR}/External/libgridxc XMLF90_ROOT=${SDIR}/External/xmlf90 PSML_ROOT=${SDIR}/External/libpsml LIBFDF_ROOT=${SDIR}/External/libfdf # Now compile following instructions in the installation guide FC=ftn cmake -S. -B_build \ -DCMAKE_INSTALL_PREFIX=${SDIR}/build \ -DWITH_LIBXC=ON-DLIBXC_ROOT=${LIBXC_ROOT} \ -DXMLF90_FIND_METHOD=source-DXMLF90_SOURCE_DIR=${XMLF90_ROOT} \ -DLIBPSML_FIND_METHOD=source -DLIBPSML_SOURCE_DIR=${PSML_ROOT} \ -DLIBGRIDXC_FIND_METHOD=source -DLIBGRIDXC_SOURCE_DIR=${GRIDXC_ROOT}\ -DLIBFDF_FIND_METHOD=source-DLIBFDF_SOURCE_DIR=${LIBFDF_ROOT} \ -DSIESTA_TOOLCHAIN=${SDIR}/Config/cmake/toolchains/intel.cmake \ -DWITH_NCDF=ON \ -DNetCDF_ROOT=/opt/cray/pe/netcdf/4.9.0.1/intel/19.0/ \ -DNetCDF_INCLUDE_DIR=/opt/cray/pe/netcdf/4.9.0.1/intel/19.0/include \ -DWITH_MPI=ON \ -DWITH_LUA=OFF \ -DWITH_FLOOK=OFF \ -DWITH_DFTD3=OFF \ -DBLAS_LIBRARY="-L${MKLPATH} ${MKLPATH}/libmkl_blas95_lp64.a -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -lmkl_blacs_intelmpi_lp64 -liomp5 -lpthread -lm -ldl" \ -DLAPACK_LIBRARY="-L${MKLPATH} ${MKLPATH}/libmkl_lapack95_lp64.a -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -lmkl_blacs_intelmpi_lp64 -liomp5 -lpthread -lm -ldl" \ -DSCALAPACK_LIBRARY="-L${MKLPATH} -lmkl_scalapack_lp64 -lmkl_intel_lp64 -lmkl_core -lmkl_intel_thread -lmkl_blacs_intelmpi_lp64 -lpthread -lm -ldl" \ -DBLAS_LIBRARY_DIR=${MKLPATH} \ -DLAPACK_LIBRARY_DIR=${MKLPATH} \ -DSCALAPACK_LIBRARY_DIR=${MKLPATH} \ STEP 4: Output of SIESTA compilation (up to error generation point) -- The Fortran compiler identification is IntelLLVM 2023.2.0 -- The C compiler identification is IntelLLVM 2023.2.0 -- Cray Programming Environment 2.7.21 Fortran -- Detecting Fortran compiler ABI info -- Detecting Fortran compiler ABI info - done -- Check for working Fortran compiler: /opt/cray/pe/craype/2.7.21/bin/ftn - skipped -- Cray Programming Environment 2.7.21 C -- Detecting C compiler ABI info -- Detecting C compiler ABI info - done -- Check for working C compiler: /opt/cray/pe/craype/2.7.21/bin/cc - skipped -- Detecting C compile features -- Detecting C compile features - done -- Using IntelLLVM compiler -- Using toolchain: /tmp/test/SIESTA_5.0_LATEST/siesta-rel-5.0/Config/cmake/toolchains/intel.cmake -- No build type selected. SIESTA will default to 'Release'. To override pass -DCMAKE_BUILD_TYPE= in order to configure SIESTA. Available options are: * -DCMAKE_BUILD_TYPE=Release - For an optimized build with no assertions or debug info. * -DCMAKE_BUILD_TYPE=Debug - For an unoptimized build with assertions and debug info. * -DCMAKE_BUILD_TYPE=Check - For an unoptimized build with assertions and debug info + code checks. * -DCMAKE_BUILD_TYPE=RelWithDebInfo - For an optimized build with no assertions but with debug info. * -DCMAKE_BUILD_TYPE=MinSizeRel - For a build optimized for size instead of speed. -- Flags for C-compiler (build type: Release): -O2 -ip -xHost -- Flags for Fortran-compiler (build type: Release): -O2 -ip -xHost -fp-model=strict -prec-div -prec-sqrt -- Found PkgConfig: /usr/bin/pkg-config (found version "0.29.2") -- Parsing BLAS options -- Locating BLAS library -- Searching in: /opt/intel/oneapi/mkl/2023.2.0/lib/intel64 -- Found CustomBlas: -L/opt/intel/oneapi/mkl/2023.2.0/lib/intel64 /opt/intel/oneapi/mkl/2023.2.0/lib/intel64/libmkl_blas95_lp64.a
[SIESTA-L] Ghost atom in vacancy
Dear Users, In a mono-vacancy formation energy calculation of a crystal, is it "acceptable" to place ghost basis orbitals centered at the vacancy site? I say it's acceptable because the ghost basis orbitals will enhance the description of the electron density in the vicinity of the vacancy (and effectively corrects BSSE). Excluding the ghost orbitals can lead to a poor description of the charge distribution around the vacancy and yield inaccurate vacancy formation energies. In spite of my reasoning, I would like to humbly get the views of experienced theorists on this matter. Thank you and best regards, FG -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
Re: [SIESTA-L] Passing constraint harmonic energy through constr.f
Dear Prof. Garcia Many thanks for your response. I looked up the MM.Potentials block but I couldn't figure out how to apply the harmonic potential to selected atoms. Any assistance would be appreciated. Thanks! Francisco On Mon, Nov 13, 2023 at 1:51 AM Alberto Garcia wrote: > Hi Francisco, > > You could indeed modify the constr.f file to include what you need, but > first check the section of the manual entitled "Auxiliary Force Field". In > there you can find a discussion of the MM.Potentials block, in which you > can define the parameters of, among others, a harmonic potential. > > Best regards, > >Alberto > > > - El 11 de Noviembre de 2023, a las 09:13, Francisco Garcia < > garcia.ff@gmail.com> escribió: > > Dear Users, > > The constraint subroutine takes the form: subroutine constr( cell, na, > isa, amass, xa, stress, fa, ntcon ). It doesn't contain information about > energy. > > Now I am using a harmonic potential with a stiff spring to constrain a > pair of atoms separated by a distance of r0. I should be able to pass the > constraint harmonic energy V=1/2 * k * (rij - r0)^2 to the subroutine to be > added to the total energy i.e., E + 1/2 * k * (rij - r0)^2 will be the > total energy. > > Which routines should I modify to allow the total energy to account for > the constraint energy? > > > > > -- > SIESTA is supported by the Spanish Research Agency (AEI) and by the > European H2020 MaX Centre of Excellence ( > https://urldefense.com/v3/__http://www.max-centre.eu/__;!!D9dNQwwGXtA!V_yhBOHCFFxdY5WMcg2Vv8Ymb5CUuYr4dpyQyrnB8zQRuVGasPatUDFCuulHsnFCbVfyyrstpp82eApPAg3lvMc$ > ) > > > -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Bond distance constraint
Dear Users, I want the distance between two specific atoms i and j to remain constant at all times during geometry relaxation : |r_i - r_j| = c. I also want the bond to be able to translate and rotate freely. Any assistance on how this may be achieved would be appreciated. Thanks! -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] How to define customized parameter is input FDF file
Dear Users, I am doing constrained relaxation using harmonic springs. I have been able to successfully modify constr.f. However, I currently have to set the value of each constraint distance D in constr.f hand, recompile siesta and run the relaxation. I was wondering if I can define a block in the input FDF file and simply set the value D. That way I only have to compile siesta once and use different values of D in the FDF input. How may I do that? Which routines do I have to modify and what special declarations should I include in constr.f? Thank you. Regards, Francisco -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Passing constraint harmonic energy through constr.f
Dear Users, The constraint subroutine takes the form: subroutine constr( cell, na, isa, amass, xa, stress, fa, ntcon ). It doesn't contain information about energy. Now I am using a harmonic potential with a stiff spring to constrain a pair of atoms separated by a distance of r0. I should be able to pass the constraint harmonic energy V=1/2 * k * (rij - r0)^2 to the subroutine to be added to the total energy i.e., E + 1/2 * k * (rij - r0)^2 will be the total energy. Which routines should I modify to allow the total energy to account for the constraint energy? -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Vacancy formation energy with Ghost orbitals
Dear Users, I computed the vacancy formation energy, E_vac_f, of a crystal as follows: E_vac_f = E(crystal with one vacant site) - ((n-1)/n)*E(pristine crystal). In computing E(crystal with one vacant site), I kept a floating/ghost orbital at the vacant site during the geometry relaxation to account for BSSE (i.e., the net number of basis functions in the defected and pristine crystals are the same). I would like to know from the experts if the procedure is correct. Thank you. Francisco -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Unfolding phonon dispersion
Dear Users I computed the phonon eigenvalues of a 108-atom FCC supercell (3x3x3 supercell built from a 4-atom FCC conventional unit cell). I was forced to do so because of a special magnetic configuration required to obtain the ground state. Is it possible to unfold the supercell phonon dispersion to obtain the primitive cell phonon dispersion in SIESTA? I know that electronic bands can be unfolded (SIESTA has a utility for this), but I couldn't find anything for phonons. Any advice or assistance would be highly appreciated. Thank you. Francisco. -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Gaussian basis vs numerical basis
Dear users, I would like to know why calculations with a Gaussian basis are in general more expensive than a numerical basis like the type used in SIESTA. Both Gaussian and numerical basis sets are atom-centered so one would expect their performance not to be too far apart. Is the SIESTA basis set faster because the radial part of the basis is stored on a log grid? Or is it because the confining potential shortens the radial extent of the basis set, thereby making it more compact? Any inputs from the experts would be appreciated. Thank you. -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Siesta-master Pseudopotential issues
Dear Users/Developers I noticed that siesta-master (the latest version) has a problem reading some pseudopotentials whereas siesta-v 4.1.5 has no such issue. For example: version 4.1.5 can read Ga psf with 3d10 4s2 4p1 or 4s2 4p1 4f0 as the reference state. However, the latest siesta-master fails to read the 3d10 4s2 4p1 psf but reads 4s2 4p1 4f0 successfully. The error message in siesta-master is: Total charge in occupied basis states different from valence charge. Is there any particular reason why one version of SIESTA reads both psfs and the other can't? Thanks! PS: I was forced to resort to siesta-master because I get the following error is siesta v 4.1.5 (the error occurs for a few systems; setting the optimization flag to -O0 and recompiling doesn't alleviate the problem is v 4.1.5): Image PCRoutineLineSource libpthread-2.31.s 14C32F8A88C0 Unknown Unknown Unknown siesta 00BD6729 gpfa2f329 fft1d.F siesta 00BD5F63 gpfa 194 fft1d.F siesta 004560F6 fft 193 fft.F siesta 00570A11 poison103 poison.F siesta 004623C1 dhscf_init406 dhscf.F siesta 0075E338 setup_h0 141 setup_H0.F siesta 006077FF siesta_forces 219 siesta_forces.F siesta 00BAF618 siesta 74 siesta.F siesta 0041C9BD Unknown Unknown Unknown libc-2.31.so 14C32A3712BD __libc_start_main Unknown Unknown siesta 0041C8EA Unknown Unknown Unknown forrtl: severe (174): SIGSEGV, segmentation fault occurred -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
Re: [SIESTA-L] Need help with Phonon Dispersion Band Lines in Vibra
Dear Prof. Postnikov, Many thanks and appreciation for your response. I believe I found a solution to my problem but I want to run it by you. First, an FCC cell with 2 unique atoms is equivalent to a tetragonal cell (this is the smallest unit cell to model antiferromagnetism). Using the website https://urldefense.com/v3/__https://www.materialscloud.org/work/tools/seekpath__;!!D9dNQwwGXtA!WbTOEoZkczNXgxFdjZz1Ho66DZKcaHdWandSB0IytM7jNKQONZgfMHtxCird0d0mH_PsbycR0mostpvReHj7Iw$ , the high symmetry points in the Brillouin zone are as follows (each set of points is scaled by the corresponding pi/a): Standard FCC primitive cell: Gamma (0,0,0), X(0,2,0), K(1.5,1.5,0), W(1,2,0), L(1,1,1) 2-atom tetragonal cell: Gamma(0,0,0), X(0,1,0), M(1,1,0), R(0,1,0.707107), A(1,1,0.707107), Z(0,0,0.707107). With this information, I believe the two Vibra inputs below, one for the primitive FCC cell and the other for 2-atom tetragonal cell, are formally equivalent (the last two k-points in each case, i.e. L and M, is what I'm a bit unsure about). Thank you very much for your kindness & happy holidays. (A) Primitive FCC cell: NumberOfAtoms 1 #Lattice parameters LatticeConstant 3.47 Ang %block LatticeVectors 0.50 0.50 0.00 0.50 0.00 0.50 0.00 0.50 0.50 %endblock LatticeVectors #Atomic positions AtomicCoordinatesFormat Fractional %block AtomicCoordinatesAndAtomicSpecies 0.00 0.00 0.00 1 54.938 %endblock AtomicCoordinatesAndAtomicSpecies #High symmetry Brillouin zones points scaled by pi/a: Gamma (0,0,0), X(0,2,0), K(1.5,1.5,0), W(1,2,0), L(1,1,1) BandLinesScale pi/a %block BandLines 1 0.000 0.000 0.000 \Gamma 30 0.000 2.000 0.000 X 30 2.000 2.000 2.000 \Gamma 30 1.000 1.000 1.000 L %endblock BandLines (B) 2-atom tetragonal cell to model antiferromagnetism (this is double the volume of the FCC primitive cell) NumberOfAtoms 2 #Lattice parameters LatticeConstant 2.453660531 Ang #[this is the FCC lattice constant divided by sqrt(2)] %block LatticeVectors 1.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 1.414213562 %endblock LatticeVectors #Atomic positions AtomicCoordinatesFormat Fractional %block AtomicCoordinatesAndAtomicSpecies 0.00 0.00 0.00 1 54.938 0.50 0.50 0.50 1 54.938 %endblock AtomicCoordinatesAndAtomicSpecies #High symmetry Brillouin zones points scaled by pi/a: Gamma(0,0,0), X(0,1,0), M(1,1,0), R(0,1,0.707107), A(1,1,0.707107), Z(0,0,0.707107) BandLinesScale pi/a %block BandLines 1 0.000 0.000 0.000 \Gamma 30 0.000 1.000 0.000 X 30 1.000 1.000 1.000 \Gamma 30 2.000 2.000 2.000 M %endblock BandLines On Thu, Dec 29, 2022 at 3:34 PM Andrei Postnikov < andrei.postni...@univ-lorraine.fr> wrote: > Dear Francisco, > it is difficult to give a useful advice on the basis of very limited > information you provide, > but my impression is that your problems are not obviously related with > Vibra. > Some questions: > 1. What (magnetic) structure are you modelling? How comes you have four > atoms per AFM unit cell? > Can there be two? > 2. Is electronic structure (and band dispersions) correct, prior to any > phonons? > 3. What means "incorrect phonon dispersion"? Do you have problems with > crystallography / > choosing the q-path, or is your calculation basically wrong? > 4. With 4 atoms as you use it so far, the Gamma phonon calculation would > yield > 9 modes, which would map genuine zone-center and zone-boundary modes. > Do they come out reasonably? > > To your problem: > "Basically I want to alter the band lines in input 2 so that they are > equivalent to the band lines in input 1" - > you have > BandLinesScale pi/a > in both inputs, the same lattice parameter, and the same definition of > path. > So if everything is correctly read, you must get the same Cartesian q-path > in both cases. > Either this is not so and there is something wrong with the input, > or the paths are identical but your problem is elsewhere. > > Best regards > > Andrei > > > > > > > > > - Le 29 Déc 22, à 0:40, garcia ff 000 a > écrit : > > Dear Users, > > I have appended 2 Vibra inputs below for computing the phonon dispersion > for FCC Mn. > > Input 1 works fine as it gives the expected band shapes for the dispersion > (but the frequencies are off). The main issue with input 1 is that it is > not suitable for antiferromagnetic calculations since there is only one Mn > atom in the primitive cell. > > This led me to consider input 2, which has 4 atoms in the unit cell and > can be used for antiferromagnetic calculations. The issue with input 2 is > that the bandlines yield an incorrect phonon dispersion. This is what I > need your help on. Basically I want to alter the band lines in input 2 so > that they are equivalent to the band lines in input 1. > > Any assistance with this, especially from the Vibra authors, would be >
Re: [SIESTA-L] Need help with Phonon Dispersion Band Lines in Vibra
Dear Prof. Postnikov, Special thanks again for your kind and prompt response. (i) The q-paths in input 2 must be rotated by 45 degrees to maintain consistency with the paths in input 1. I would like to kindly know the axis about which the rotation should be performed. (ii) Regarding your suggestion of using the lattice vectors [1/2 -1/2 0], [1/2 1/2 0], [0 0 1], I can't quite figure out why the q-paths will be the same as the q-paths for input 1. The primitive FCC lattice vectors in input 1 are [1/2 1/2 0], [1/2 0 1/2], and [0 1/2 1/2] (with an angle of 60 degrees between each pair of lattice vectors). The angle between the lattice vectors you proposed is 90 degrees for each pair. I'm having a difficulty grasping how this new set of lattice vectors yields the same q-paths as input 1. Any further explanation would be highly appreciated. Thank you very much Professor. On Fri, Dec 30, 2022 at 5:12 PM Andrei Postnikov < andrei.postni...@univ-lorraine.fr> wrote: > Dear Francisco, > so your AFM structure is of CuAu type. This is fine but of course > this is not the only AFM structure possible (and I don't know whether it is > realistic at all, but this of course depends on your objectives). > Now, if you want the q-path to be the same in your two settings, > you should consider that the second one is rotated by 45 degrees. > That is, if you choose the Gamma->X direction || [010] in the first > setting > it must be || [110] in the second setting, with the lattice vectors you > use. > Otherwise define the lattice vectors as [1/2 -1/2 0], [1/2 1/2 0], [0 0 1] > with the same lattice parameter as in the first setting > and enjoy the same coordinates of q points (cartesian, in terms of pi/a) > in both settings. > > Best regards > > Andrei > > > - Le 31 Déc 22, à 0:24, garcia ff 000 a > écrit : > > Dear Prof. Postnikov, > > Many thanks and appreciation for your response. I believe I found a > solution to my problem but I want to run it by you. > > First, an FCC cell with 2 unique atoms is equivalent to a tetragonal cell > (this is the smallest unit cell to model antiferromagnetism). > > Using the website > https://urldefense.com/v3/__https://www.materialscloud.org/work/tools/seekpath__;!!D9dNQwwGXtA!RtHnhiKf-1Tr0ZJZ17JGCt-WslBldkrwcANZ1KZuejXMskHtY6-_llZSytEsPLSiAktiz7DhWg1EBpjE-NmrWg$ > , the > high symmetry points in the Brillouin zone are as follows (each set of > points is scaled by the corresponding pi/a): > > Standard FCC primitive cell: Gamma (0,0,0), X(0,2,0), K(1.5,1.5,0), > W(1,2,0), L(1,1,1) > > 2-atom tetragonal cell: Gamma(0,0,0), X(0,1,0), M(1,1,0), R(0,1,0.707107), > A(1,1,0.707107), Z(0,0,0.707107). > > With this information, I believe the two Vibra inputs below, one for the > primitive FCC cell and the other for 2-atom tetragonal cell, are formally > equivalent (the last two k-points in each case, i.e. L and M, is what I'm a > bit unsure about). > > > Thank you very much for your kindness & happy holidays. > > > > (A) Primitive FCC cell: > > NumberOfAtoms 1 > > #Lattice parameters > LatticeConstant 3.47 Ang > %block LatticeVectors > 0.50 0.50 0.00 > 0.50 0.00 0.50 > 0.00 0.50 0.50 > %endblock LatticeVectors > > #Atomic positions > AtomicCoordinatesFormat Fractional > %block AtomicCoordinatesAndAtomicSpecies > 0.00 0.00 0.00 1 54.938 > %endblock AtomicCoordinatesAndAtomicSpecies > > #High symmetry Brillouin zones points scaled by pi/a: Gamma (0,0,0), > X(0,2,0), K(1.5,1.5,0), W(1,2,0), L(1,1,1) > > BandLinesScale pi/a > %block BandLines > 1 0.000 0.000 0.000 \Gamma > 30 0.000 2.000 0.000 X > 30 2.000 2.000 2.000 \Gamma > 30 1.000 1.000 1.000 L > %endblock BandLines > > > > (B) 2-atom tetragonal cell to model antiferromagnetism (this is double the > volume of the FCC primitive cell) > > NumberOfAtoms 2 > > #Lattice parameters > LatticeConstant 2.453660531 Ang #[this is the FCC lattice constant > divided by sqrt(2)] > %block LatticeVectors > 1.00 0.00 0.00 > 0.00 1.00 0.00 > 0.00 0.00 1.414213562 > %endblock LatticeVectors > > #Atomic positions > AtomicCoordinatesFormat Fractional > %block AtomicCoordinatesAndAtomicSpecies > 0.00 0.00 0.00 1 54.938 > 0.50 0.50 0.50 1 54.938 > %endblock AtomicCoordinatesAndAtomicSpecies > > #High symmetry Brillouin zones points scaled by pi/a: Gamma(0,0,0), > X(0,1,0), M(1,1,0), R(0,1,0.707107), A(1,1,0.707107), Z(0,0,0.707107) > > BandLinesScale pi/a > %block BandLines > 1 0.000 0.000 0.000 \Gamma > 30 0.000 1.000 0.000 X > 30 1.000 1.000 1.000 \Gamma > 30 2.000 2.000 2.000 M > %endblock BandLines > > > > On Thu, Dec 29, 2022 at 3:34 PM Andrei Postnikov < > andrei.postni...@univ-lorraine.fr> wrote: > >> Dear Francisco, >> it is difficult to give a useful advice on the basis of very limited >> information you provide, >> but my impression is that your
[SIESTA-L] Need help with Phonon Dispersion Band Lines in Vibra
Dear Users, I have appended 2 Vibra inputs below for computing the phonon dispersion for FCC Mn. Input 1 works fine as it gives the expected band shapes for the dispersion (but the frequencies are off). The main issue with input 1 is that it is not suitable for antiferromagnetic calculations since there is only one Mn atom in the primitive cell. This led me to consider input 2, which has 4 atoms in the unit cell and can be used for antiferromagnetic calculations. The issue with input 2 is that the bandlines yield an incorrect phonon dispersion. This is what I need your help on. Basically I want to alter the band lines in input 2 so that they are equivalent to the band lines in input 1. Any assistance with this, especially from the Vibra authors, would be greatly appreciated. Thank you very much for your kind assistance and God Bless! Francisco #INPUT 1 (1 atom in the FCC primitive cell; 125 atoms in Supercell) SystemName fccMn_1 SystemLabelfccMn_1 NumberOfAtoms 1 LatticeConstant3.47 Ang %block LatticeVectors 0.50 0.50 0.00 0.50 0.00 0.50 0.00 0.50 0.50 %endblock LatticeVectors AtomicCoordinatesFormat Fractional %block AtomicCoordinatesAndAtomicSpecies 0.00 0.00 0.00 1 54.938 %endblock AtomicCoordinatesAndAtomicSpecies SuperCell_1 2 SuperCell_2 2 SuperCell_3 2 AtomicDispl 0.04 Bohr BandLinesScale pi/a %block BandLines 1 0.000 0.000 0.000 \Gamma 30 2.000 0.000 0.000 X 30 2.000 2.000 2.000 \Gamma 30 1.000 1.000 1.000 L %endblock BandLines Eigenvectors True #INPUT 2 (4 atoms in the FCC conventional cell; 108 atoms in Supercell) SystemName fccMn_4 SystemLabelfccMn_4 NumberOfAtoms 4 LatticeConstant3.47 Ang %block LatticeVectors 1.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 1.00 %endblock LatticeVectors AtomicCoordinatesFormat Fractional %block AtomicCoordinatesAndAtomicSpecies 0.00 0.00 0.00 1 54.938 0.50 0.50 0.00 1 54.938 0.50 0.00 0.50 1 54.938 0.00 0.50 0.50 1 54.938 %endblock AtomicCoordinatesAndAtomicSpecies SuperCell_1 1 SuperCell_2 1 SuperCell_3 1 AtomicDispl 0.04 Bohr BandLinesScale pi/a # The band lines below are incorrect. %block BandLines 1 0.000 0.000 0.000 \Gamma 30 2.000 0.000 0.000 X 30 2.000 2.000 2.000 \Gamma 30 1.000 1.000 1.000 L %endblock BandLines Eigenvectors True -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Advice on GGA+U convergence issues
Dear Users, Can anyone give me some advice on how best to accelerate the convergence of GGA+U calculations for f-lectron systems. Starting the standard GGA calculations, which converged with no issues, I increased U in small steps: U=0.1, U=0,2, etcbut even U=0.1 couldn't converge. Any generic advice on how to alleviate this problem will be appreciated. Thanks! -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] No q-points in force constant matrix calculation
Dear Users, I am a tad confused about the calculation of the Force Constant (FC) matrix calculation in SIESTA. I would like to know why the FC matrix isn't computed for different q-points. Instead, a single FC matrix is computed. Wouldn't that affect the accuracy of the phonon spectra? Thank you. -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Compiling SIESTA on HPE Cray EX Supercomputer
Dear Users, Can someone kindly share his or her arch.make for compiling SIESTA on an HPE Cray EX Supercomputer? Compiling SIESTA on CRAY is turning out to be more problematic than I anticipated. Thank you. -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Wave vectors for Phonon Dispersion in SIESTA
Dear Users, I would appreciate it if you can suggest resources (textbooks, papers, review articles, software, etc) on how to generate the Brillouin zone wave vectors for computing the phonon dispersion curves for all the standard crystals. Inputs in SIESTA format for unit cells and supercells of different single-component or multicomponent crystals would be highly appreciated. Thank you. -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
Re: [SIESTA-L] Total energy
Thank you very much for the reference! I'll use E=1/2(Etot +FreeEng) for my E-V equation of state fit. I noticed that for semiconductors, Etot = FreeEng; it is in metallic systems that I see differences between Etot and FreeEng (with FreeEng being lower). In this regard, using the average of the two energies for T=0 makes sense. Thanks again! On Wed, Jun 22, 2022 at 3:01 PM RCP wrote: > Hello, > May I contribute my 0.01 ?. According to Eq.(21) in, > >Kresse & Furthmüller, Comp. Mat. Sci. 6 (1996) 15-50 > > a good estimate for the (T=0K) DFT energy is 1/2(Etot +FreeEng), and > this is the recipe I've always used. > > Regards, > Roberto > > On 22/06/2022 13:11, Francisco Garcia wrote: > > Thanks for your response Prof. Artacho! > > > > So if I have a situation like below for a large metallic system in which > FreeEng > > is lower than Etot, should I use Etot? Specifically, if I want to plot > the E-V > > curve from single point runs to obtain the equation of state I should > use Etot > > instead of FreeEng? Or if I want to compute the cohesive energy, Etot > should be > > used instead of FreeEng? > > > > siesta: Program's energy decomposition (eV): > > siesta: Ebs   =    -859.157108 > > siesta: Eions  =    9568.777238 > > siesta: Ena   =    277.651800 > > siesta: Ekin   =    4247.415878 > > siesta: Enl   =   -1698.956164 > > siesta: Eso   =     0.00 > > siesta: Edftu  =     0.00 > > siesta: DEna   =    162.636813 > > siesta: DUscf  =     26.981770 > > siesta: DUext  =     0.00 > > siesta: Exc   =   -6715.533109 > > siesta: eta*DQ  =     0.00 > > siesta: Emadel  =     0.00 > > siesta: Emeta  =     0.00 > > siesta: Emolmec =     0.00 > > siesta: Ekinion =     0.00 > > siesta: Eharris =    -13268.580292 > > siesta: Etot   =     -13268.580250 > > siesta: FreeEng =   -13269.911499 > > > > On Tue, Jun 21, 2022 at 2:07 PM Emilio Artacho > <mailto:e.arta...@nanogune.eu>> wrote: > > > > > >> On 20 Jun 2022, at 17:16, Francisco Garcia >> <mailto:garcia.ff@gmail.com>> wrote: > >> > >> Dear users, > >> > >> In metallic systems with a fairly sizable electronic smearing > temperature > >> T, is it accurate to claim that > >> > >> (i) in a single point calculation, the free energy is the > representative > >> energy of the system (due to the addition of -TS to the total > energy U) > > > > No, the -TS term for the phononic entropy is missing. > > > > It is the free energy defining the finite temperature equilibrium > for the > > purely electronic > > problem for fixed external potential (fixed nuclei), as in > Mermin’s finite > > -T DFT. > > > > > >> (ii) in a variable cell optimization, the enthalpy is the > representative > >> energy of the system (due to the addition of the PV term to the > energy). > > > > Yes, it is the free energy of the system (minimum defines > equilibrium) for T=0. > > > > Emilio > > > >> > >> Thanks! > >>  > >> > >> -- > >> SIESTA is supported by the Spanish Research Agency (AEI) and by the > >> European H2020 MaX Centre of Excellence (http://www.max-centre.eu/) > > > > -- > > Emilio Artacho > > > > Theory Group, Nanogune, 20018 San Sebastian, Spain, and > > Theory of Condensed Matter, Department of Physics, > > Cavendish Laboratory, University of Cambridge, Cambridge CB3 0HE, UK > > > > > > > > > > > > -- > > SIESTA is supported by the Spanish Research Agency (AEI) and by the > European > > H2020 MaX Centre of Excellence (http://www.max-centre.eu/) > > > > > > > > -- > Dr. Roberto C. Pasianot Phone: 54 11 4839 6709 > Gcia. Materiales, CAC-CNEA FAX : 54 11 6772 7362 > Avda. Gral. Paz 1499pasia...@cnea.gov.ar > 1650 San Martin, Buenos Aires > ARGENTINA > -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
Re: [SIESTA-L] Total energy
Thanks for your response Prof. Artacho! So if I have a situation like below for a large metallic system in which FreeEng is lower than Etot, should I use Etot? Specifically, if I want to plot the E-V curve from single point runs to obtain the equation of state I should use Etot instead of FreeEng? Or if I want to compute the cohesive energy, Etot should be used instead of FreeEng? siesta: Program's energy decomposition (eV): siesta: Ebs = -859.157108 siesta: Eions = 9568.777238 siesta: Ena = 277.651800 siesta: Ekin= 4247.415878 siesta: Enl = -1698.956164 siesta: Eso = 0.00 siesta: Edftu = 0.00 siesta: DEna= 162.636813 siesta: DUscf =26.981770 siesta: DUext = 0.00 siesta: Exc = -6715.533109 siesta: eta*DQ = 0.00 siesta: Emadel = 0.00 siesta: Emeta = 0.00 siesta: Emolmec = 0.00 siesta: Ekinion = 0.00 siesta: Eharris = -13268.580292 siesta: Etot=-13268.580250 siesta: FreeEng =-13269.911499 On Tue, Jun 21, 2022 at 2:07 PM Emilio Artacho wrote: > > On 20 Jun 2022, at 17:16, Francisco Garcia > wrote: > > Dear users, > > In metallic systems with a fairly sizable electronic smearing temperature > T, is it accurate to claim that > > (i) in a single point calculation, the free energy is the representative > energy of the system (due to the addition of -TS to the total energy U) > > > No, the -TS term for the phononic entropy is missing. > > It is the free energy defining the finite temperature equilibrium for the > purely electronic > problem for fixed external potential (fixed nuclei), as in Mermin’s finite > -T DFT. > > > (ii) in a variable cell optimization, the enthalpy is the representative > energy of the system (due to the addition of the PV term to the energy). > > > Yes, it is the free energy of the system (minimum defines equilibrium) for > T=0. > > Emilio > > > Thanks! > > > -- > SIESTA is supported by the Spanish Research Agency (AEI) and by the > European H2020 MaX Centre of Excellence (http://www.max-centre.eu/) > > > -- > Emilio Artacho > > Theory Group, Nanogune, 20018 San Sebastian, Spain, and > Theory of Condensed Matter, Department of Physics, > Cavendish Laboratory, University of Cambridge, Cambridge CB3 0HE, UK > > > > > > -- > SIESTA is supported by the Spanish Research Agency (AEI) and by the > European H2020 MaX Centre of Excellence (http://www.max-centre.eu/) > -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Total energy
Dear users, In metallic systems with a fairly sizable electronic smearing temperature T, is it accurate to claim that (i) in a single point calculation, the free energy is the representative energy of the system (due to the addition of -TS to the total energy U) (ii) in a variable cell optimization, the enthalpy is the representative energy of the system (due to the addition of the PV term to the energy). Thanks! -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Variable cell optimization spin-orbit inclusion in version 4.1.5
Dear users/developers, Is it possible to perform variable cell optimizations with spin-orbit coupling (SOC) in version 4.1.5? Does the warning message below imply that SOC shouldn't be used in version 4.1.5 (especially for strongly localized d & f electron systems): *"This spin-orbit implementation uses a local approximation. From 4.2 and onwards the full non-local approximation is implemented. You are strongly advised to use >=4.2 versions!"* I would like to ask the developers: (i) when will a version with fully functioning SOC will be released? (ii) Is the current implementation of SOC suitable for variable cell calculations in systems with strong SOC effects? Thank you very much in advance! -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] A question about basis set configurations
Dear Users, I have a question on basis sets. In the four cases below, is case 1 equivalent to case 2? Is case 3 equivalent to case 4? If not, then why? Thank you. Case 1 %block PAO.Basis O 2 n=202 0.0 0.0 1.0 1.0 n=213 P 1 0.0 0.0 0.0 1.0 1.0 1.0 %endblock PAO.Basis Case 2 %block PAO.Basis O 2 n=202 0.0 0.0 1.0 1.0 n=213 0.0 0.0 0.0 1.0 1.0 1.0 n=32 1 0.0 1.0 %endblock PAO.Basis Case 3 %block PAO.Basis Zn 2 n=4 0 2 P 1 0.0 0.0 1.0 1.0 n=3 2 1 0.0 1.0 %endblock PAO.Basis Case 4 %block PAO.Basis Zn 2 n=4 0 2 0.0 0.0 1.0 1.0 n=4 1 1 0.0 1.0 n=3 2 1 0.0 1.0 %endblock PAO.Basis -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
Re: [SIESTA-L] Forces on crystalline atoms are not zero
Thanks for your responses. I figured out my errors. The problem has been resolved. On Mon, Mar 28, 2022 at 3:50 PM Andrei Postnikov < andrei.postni...@univ-lorraine.fr> wrote: > Dear Francisco: > as a general statement, your expectation is justified. However, unexpected > things happen. > Having posted an (input) / output file would help to pinpoint a problem. > Without any additional information, I'd guess that your lattice vectors > are not exactly fcc, > or the atoms are not exactly at symmetric positions. (E.g., you think that > you scale > the atom coordinates with lattice parameter, but in fact you don't). > > Best regards > > Andrei Postnikov > > > - Francisco Garcia a écrit : > > > Dear Users > > I performed a series of single point energy calculations on an FCC crystal > by varying the lattice constant. I'm trying to generate data points to > optimize the lattice constant > > I was expecting the force on each atom to zero by virtue of crystalline > symmetry (no internal degree of freedom). But to my surprise, the forces > are not zero. Increasing the Mesh cutoff from 300 Ry to higher values > didn't help either. How is it possible for the atomic forces in a crystal > to be non-zero? > > > -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Forces on crystalline atoms are not zero
Dear Users I performed a series of single point energy calculations on an FCC crystal by varying the lattice constant. I'm trying to generate data points to optimize the lattice constant I was expecting the force on each atom to zero by virtue of crystalline symmetry (no internal degree of freedom). But to my surprise, the forces are not zero. Increasing the Mesh cutoff from 300 Ry to higher values didn't help either. How is it possible for the atomic forces in a crystal to be non-zero? -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Using prior calculation for spin-orbit calculation
Dear Users, Is it allowed for me to add spin-orbit coupling (SOC) to a prior regular calculation? Or do I have to start the SOC calculation from scratch? Thanks! -- SIESTA is supported by the Spanish Research Agency (AEI) and by the European H2020 MaX Centre of Excellence (http://www.max-centre.eu/)
[SIESTA-L] Blyp psf for Silicon
Dear users, Does anyone know where I can obtain a BLYP pseudopotential for Si? Thank you.
[SIESTA-L] Large forces in siesta molecular dynamics (data provided in email)
Hello users: I am running a molecular dynamics simulation with SIESTA and the job keep running into situations where the forces grow large all of a sudden and then it will go back to normal in the next time step. For example: take a look below at the maximum atomic force printed for the MD steps [this is not the entire trajaectory; I am just showing a tiny segment to illustrate this problem]. Max8.169457constrained Max7.725134constrained Max7.231393constrained Max6.698717constrained Max6.434356constrained Max 180.612893constrained <- large force of 180.61 eV/Ang here Max6.351384constrained Max6.500315constrained Max6.616196constrained Max6.705574constrained Max7.375029constrained Max7.971425constrained Max8.440149constrained Max8.768627constrained Now take a look at the SCF loop for the MD step where the large force was generated: iscf Eharris(eV)E_KS(eV) FreeEng(eV) dDmax Ef(eV) dHmax(eV) scf:1 -217188.442539 -217204.666859 -217206.302454 0.064243 -1.492048 5.165406 scf:2 -146020.001861 -193579.811302 -193580.989056 1.710060 -0.338743357.676349 scf:3 -217800.695121 -217208.716553 -217210.351765 1.701213 -1.495349 6.087548 scf:4 -217213.535597 -217211.812703 -217213.449302 0.032944 -1.494384 2.570971 scf:5 -217212.451595 -217212.518659 -217214.156524 0.027437 -1.493541 0.659157 scf:6 -217212.666366 -217212.617779 -217214.255836 0.005069 -1.493661 0.297730 scf:7 -217212.671543 -217212.651068 -217214.289179 0.002601 -1.493708 0.171617 scf:8 -217212.670437 -217212.664758 -217214.302857 0.002015 -1.493723 0.141605 scf:9 -217212.676500 -217212.671768 -217214.309790 0.001025 -1.493792 0.096767 scf: 10 -217212.679444 -217212.676034 -217214.314015 0.000846 -1.493810 0.083722 scf: 11 -217212.682014 -217212.679321 -217214.317264 0.000677 -1.493834 0.059076 scf: 12 -217212.684219 -217212.682043 -217214.319956 0.000602 -1.493864 0.037595 scf: 13 -217212.684718 -217212.683469 -217214.321359 0.000377 -1.493896 0.032728 scf: 14 -217212.685645 -217212.684620 -217214.322474 0.000429 -1.493929 0.029768 scf: 15 -217212.686201 -217212.685444 -217214.323249 0.000401 -1.493955 0.027512 scf: 16 -217212.687024 -217212.686259 -217214.323994 0.000466 -1.493982 0.024955 scf: 17 -217212.687903 -217212.687112 -217214.324756 0.000566 -1.494013 0.021878 scf: 18 -217212.688413 -217212.687788 -217214.325350 0.000508 -1.494041 0.018827 scf: 19 -217212.689186 -217212.688519 -217214.325984 0.000660 -1.494070 0.014994 scf: 20 -217212.689633 -217212.689113 -217214.326468 0.000768 -1.494101 0.010350 scf: 21 -217212.689427 -217212.689275 -217214.326592 0.000290 -1.494111 0.008499 scf: 22 -217212.689461 -217212.689369 -217214.326663 0.000185 -1.494118 0.006871 scf: 23 -217212.689521 -217212.689448 -217214.326723 0.000161 -1.494124 0.006698 scf: 24 -217212.689528 -217212.689492 -217214.326755 0.75 -1.494126 0.006319 scf: 25 -217212.689561 -217212.689528 -217214.326778 0.63 -1.494128 0.005641 scf: 26 -217212.689564 -217212.689544 -217214.326788 0.47 -1.494127 0.004544 scf: 27 -217212.689564 -217212.689554 -217214.326799 0.53 -1.494124 0.003604 scf: 28 -217212.689567 -217212.689561 -217214.326808 0.53 -1.494119 0.002869 scf: 29 -217212.689568 -217212.689565 -217214.326817 0.59 -1.494113 0.002762 scf: 30 -217212.689569 -217212.689567 -217214.326822 0.37 -1.494108 0.002618 scf: 31 -217212.689566 -217212.689567 -217214.326826 0.28 -1.494103 0.002409 scf: 32 -217212.689566 -217212.689567 -217214.326831 0.32 -1.494099 0.002048 scf: 33 -217212.689569 -217212.689568 -217214.326835 0.38 -1.494095 0.001620 scf: 34 -217212.689576 -217212.689572 -217214.326839 0.37 -1.494093 0.001183 scf: 35 -217212.689582 -217212.689577 -217214.326841 0.25 -1.494093 0.000872 SCF Convergence by DM+H criterion max |DM_out - DM_in| : 0.250455 max |H_out - H_in| (eV) : 0.0008717786 SCF cycle converged after 35 iterations Using DM_out to compute the final energy and forces No. of atoms with KB's overlaping orbs in proc 0. Max # of overlaps: 222 160 siesta: E_KS(eV) = -217212.6890 siesta: Atomic forces (eV/Ang): Tot -87.051127 -112.933365 -273.703901 Max 180.612893 Res 19.109629sqrt( Sum f_i^2 / 3N ) Max 180.612893constrained As far as I can tell, there is nothing wrong in the SCF convergence
[SIESTA-L] 75,000 atom system with orderN
Hello there, I am interested in optimizing a 75,000 disordered carbon system using the orderN approach. I would like to know if the calculations commence with full diagonalization, followed by orderN (not preferred) or the entire calculation commences with order N and stays orderN (preferred). I am asking because full diagonalization will be nearly impossible (even for a single geometry step). Thanks!
[SIESTA-L] Large spikes in pressure (and atomic forces during MD)
Dear all, I always encounter large spikes in the atomic forces in SIESTA molecular dynamics (and subsequently large cell pressure), especially when the MD simulations are quite long. I am curious as to why and how that happens. Once the forces spike up, the system is thrown out of equilibrium; sometimes the system recovers after a few steps and sometimes the system blows up. This happens in nearly all MD runs which run for several picoseconds. Can someone give me a few insights? Is it related my compilation optimization flags? Thank you.