Hello all, I'm trying to use SIESTA to calculate the phonon spectrum of two different 1D chains using the method perscribed in the vibra routine. I cannot, however, get past the calculation of the force matrix. For each system that I study, I use the relaxed coordinates (per each specific exchange-correlation functional) as input, create the FC.fdf with fcbuild with supercell_3=1 (1 above, 1 below in the z-direction, have also tried supercell_3=2 with the same results) and then run this file in SIESTA like the tutorials say to do:
SystemLabel LDA_111_68_vibra_PD NumberOfSpecies 2 %block ChemicalSpeciesLabel 1 14 Si 2 1 H %endblock ChemicalSpeciesLabel PAO.BasisSize DZP MeshCutoff 200.0 Ry %include FC.fdf I would use this file for the relaxed geometry corresponding to the stationary point found with an LDA xc-functional + PP and the PAO.BasisSize DZP in the original calculation. For other systems like a BLYP xc-functional and user-defined basis I would use those parameters in the input file (but with roughly the same results, which follow). Defaults are used for everything else, including the OrderN method since the full unit cell for the FC calculation is about 200 atoms. The FC calculation seems to run smoothly until I get something like this in the standard output: ********************************************************************************** cgwf: iter = 984 grad = -0.000241 Eb(Ry) = -203.356707 cgwf: iter = 985 grad = -0.000241 Eb(Ry) = -203.356720 cgwf: iter = 986 grad = -0.000241 Eb(Ry) = -203.356734 cgwf: iter = 987 grad = -0.000243 Eb(Ry) = -203.356748 cgwf: iter = 988 grad = -0.000248 Eb(Ry) = -203.356762 cgwf: iter = 989 grad = -0.000257 Eb(Ry) = -203.356775 cgwf: iter = 990 grad = -0.000260 Eb(Ry) = -203.356788 cgwf: iter = 991 grad = -0.000263 Eb(Ry) = -203.356801 cgwf: iter = 992 grad = -0.000257 Eb(Ry) = -203.356814 cgwf: iter = 993 grad = -0.000256 Eb(Ry) = -203.356827 cgwf: iter = 994 grad = -0.000250 Eb(Ry) = -203.356840 cgwf: iter = 995 grad = -0.000253 Eb(Ry) = -203.356853 cgwf: iter = 996 grad = -0.000252 Eb(Ry) = -203.356866 cgwf: iter = 997 grad = -0.000252 Eb(Ry) = -203.356880 cgwf: iter = 998 grad = -0.000251 Eb(Ry) = -203.356893 cgwf: iter = 999 grad = -0.000244 Eb(Ry) = -203.356906 cgwf: iter = 1000 grad = -0.000238 Eb(Ry) = -203.356920 cgwf: Maximum number of CG iterations reached denmat: qtot (before DM normalization) = 765.6153 ordern: qtot (after DM normalization) = 546.0000 siesta: iscf = 2 Eharris(eV) = -12344.2602 E_KS(eV) = -12020.4605 dDmax = 0.2649 ordern: enum = 546.0000 cgwf: iter = 1 grad = -67.320865 Eb(Ry) = -216.535666 cgwf: iter = 2 grad = -149.630916 Eb(Ry) = -224.149328 cgwf: iter = 3 grad = -165.080016 Eb(Ry) = -229.107766 cgwf: iter = 4 grad = -122.780882 Eb(Ry) = -235.802701 cgwf: iter = 5 grad = -128.457462 Eb(Ry) = -240.351352 cgwf: iter = 6 grad = -107.815766 Eb(Ry) = -244.409131 cgwf: iter = 7 grad = -105.320678 Eb(Ry) = -248.300033 cgwf: iter = 8 grad = -92.022410 Eb(Ry) = -251.093176 cgwf: iter = 9 grad = -88.452386 Eb(Ry) = -254.170805 cgwf: iter = 10 grad = -202.420032 Eb(Ry) = -257.806315 cgwf: iter = 11 grad = -2466.490254 Eb(Ry) = -222.354649 cgwf: iter = 12 grad = -163.364902 Eb(Ry) = -255.045915 cgwf: iter = 13 grad = -1168.957839 Eb(Ry) = -232.144564 cgwf: iter = 14 grad = -161.959255 Eb(Ry) = -254.849595 cgwf: iter = 15 grad = -1102.828538 Eb(Ry) = -232.803424 cgwf: iter = 16 grad = -157.425714 Eb(Ry) = -254.979384 cgwf: iter = 17 grad = -194.373222 Eb(Ry) = -255.299976 cgwf: iter = 18 grad = -182605.315981 Eb(Ry) = -100.555913 cgwf: iter = 19 grad = -217036.659307 Eb(Ry) = -99.721836 cgwf: iter = 20 grad = -249.395730 Eb(Ry) = -249.824155 cgwf: iter = 21 grad = -244.560067 Eb(Ry) = -253.941014 cgwf: iter = 22 grad = -152.248315 Eb(Ry) = -257.760953 cgwf: iter = 23 grad = -232.419799 Eb(Ry) = -261.955244 cgwf: iter = 24 grad = -96401477.679874 Eb(Ry) = -911.330200 cgwf: iter = 25 grad = -107154809.385566 Eb(Ry) = -912.719590 cgwf: iter = 26 grad = -1438.764936 Eb(Ry) = -262.289834 cgwf: iter = 27 grad = -1408.809784 Eb(Ry) = -262.325692 cgwf: iter = 28 grad = -186592.532696 Eb(Ry) = -109.282399 cgwf: iter = 29 grad = -217965.896792 Eb(Ry) = -108.576682 cgwf: iter = 30 grad = -302.647048 Eb(Ry) = -254.260707 cgwf: iter = 31 grad = -194.085227 Eb(Ry) = -259.234987 cgwf: iter = 32 grad = -148740.557077 Eb(Ry) = -279.409429 cgwf: iter = 33 grad = -2536890.389917 Eb(Ry) = -319.726178 cgwf: iter = 34 grad = -5780230.992347 Eb(Ry) = -331.297794 cgwf: iter = 35 grad = -132426997.296064 Eb(Ry) = -642.016407 cgwf: iter = 36 grad = -354656197.135856 Eb(Ry) = -740.210175 cgwf: iter = 37 grad = ****************** Eb(Ry) = -3046.280084 cgwf: iter = 38 grad = ****************** Eb(Ry) = -3113.114594 cgwf: iter = 39 grad = -2197752344.447697 Eb(Ry) = 19491.918692 cgwf: iter = 40 grad = -402.146008 Eb(Ry) = -261.851463 cgwf: iter = 41 grad = -430786270.291137 Eb(Ry) = 4171.581613 cgwf: iter = 42 grad = -430958482.835069 Eb(Ry) = 4171.584019 cgwf: iter = 43 grad = -82331696.109239 Eb(Ry) = -24934.797112 cgwf: iter = 44 grad = -2565001.988477 Eb(Ry) = 2255.634962 cgwf: iter = 45 grad = -80853.063402 Eb(Ry) = -174.667315 cgwf: iter = 46 grad = ****************** Eb(Ry) = 1621007.748789 cgwf: CG tolerance reached denmat: qtot (before DM normalization) = ************ ordern: qtot (after DM normalization) = 546.0000 siesta: iscf = 3 Eharris(eV) = -7881.8204 E_KS(eV) = -12044.3075 dDmax = 21.6098 ordern: enum = 546.0000 cgwf: iter = 1 grad = ****************** Eb(Ry) = ************** cgwf: iter = 2 grad = ****************** Eb(Ry) = ************** cgwf: iter = 3 grad = ****************** Eb(Ry) = ************** cgwf: iter = 4 grad = ****************** Eb(Ry) = ************** cgwf: iter = 5 grad = ****************** Eb(Ry) = ************** cgwf: iter = 6 grad = ****************** Eb(Ry) = ************** cgwf: iter = 7 grad = ****************** Eb(Ry) = ************** cgwf: iter = 8 grad = NaN Eb(Ry) = ************** cgwf: iter = 9 grad = NaN Eb(Ry) = NaN cgwf: iter = 10 grad = NaN Eb(Ry) = NaN etc. ****************************************************************************** The end result is that I get a force constants matrix full of NaN's, no matter which 1-D chain I run. I'm not sure what I'm doing wrong... has anybody run into this problem before? Could you suggest additional simulation parameters that help to converge this type of calculation? Thank you, wise SIESTA gurus. Abraham Hmiel Research Assistant College of Nanoscale Science and Engineering at SUNY Albany