1 /* $RCSfile$ 2 * $Author: hansonr $ 3 * 4 * Copyright (C) 2003-2005 Miguel, Jmol Development, www.jmol.org 5 * 6 * Contact: jmol-developers@lists.sf.net 7 * 8 * Copyright (C) 2009 Piero Canepa, University of Kent, UK 9 * 10 * Contact: pc229@kent.ac.uk or pieremanuele.canepa@gmail.com 11 * 12 * This library is free software; you can redistribute it and/or 13 * modify it under the terms of the GNU Lesser General Public 14 * License as published by the Free Software Foundation; either 15 * version 2.1 of the License, or (at your option) any later version. 16 * 17 * This library is distributed in the hope that it will be useful, 18 * but WITHOUT ANY WARRANTY; without even the implied warranty of 19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 20 * Lesser General Public License for more details. 21 * 22 * You should have received a copy of the GNU Lesser General Public 23 * License along with this library; if not, write to the Free Software 24 * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. 25 * 26 */ 27 28 package org.jmol.adapter.readers.xtal; 29 30 import javajs.util.DF; 31 import javajs.util.Lst; 32 import javajs.util.PT; 33 34 import org.jmol.adapter.smarter.AtomSetCollectionReader; 35 36 /** 37 * http://cms.mpi.univie.ac.at/vasp/ 38 * 39 * @author Pieremanuele Canepa, Room 104, FM Group School of Physical Sciences, 40 * Ingram Building, University of Kent, Canterbury, Kent, CT2 7NH United 41 * Kingdom, pc229@kent.ac.uk 42 * 43 * @version 1.0 44 */ 45 46 public class VaspOutcarReader extends AtomSetCollectionReader { 47 48 private String[] atomNames; 49 private boolean haveIonNames; 50 private int ac = 0; 51 private boolean inputOnly; 52 private boolean mDsimulation = false; //this is for MD simulations 53 private int vaspVersion; 54 55 56 @Override initializeReader()57 protected void initializeReader() { 58 isPrimitive = true; 59 setSpaceGroupName("P1"); 60 setFractionalCoordinates(true); 61 inputOnly = checkFilterKey("INPUT"); 62 } 63 64 @Override checkLine()65 protected boolean checkLine() throws Exception { 66 67 //reads if output is from vasp5 68 if (vaspVersion == 0 && line.contains(" vasp.")) { 69 readVersion(); 70 if (vaspVersion > 0) 71 this.appendLoadNote("VASP version " + vaspVersion + " " + line); 72 } else if (line.toUpperCase().startsWith(" POTCAR:")) { 73 //reads the kind of atoms namely H, Ca etc 74 readElementNames(); 75 } else if (line.contains("ions per type")) { 76 readAtomCountAndSetNames(); 77 } else if (line.contains("molecular dynamics for ions")) { 78 mDsimulation = true; 79 } else if (line.contains("direct lattice vectors")) { 80 readUnitCellVectors(); 81 } else if (ac > 0 && line.contains("position of ions in fractional coordinates")) { 82 readInitialCoordinates(); 83 if (inputOnly) 84 continuing = false; 85 } else if (line.contains("POSITION")) { 86 readPOSITION(); 87 return true; 88 } else if (line.startsWith(" FREE ENERGIE") && !mDsimulation) { 89 readEnergy(); 90 } else if (line.contains("ENERGIE OF THE ELECTRON-ION-THERMOSTAT") 91 && mDsimulation) { 92 readMdyn(); 93 } else if (line 94 .startsWith(" Eigenvectors and eigenvalues of the dynamical matrix")) { 95 readFrequency(); 96 } 97 return true; 98 } 99 readVersion()100 private void readVersion() { 101 String[] tokens = PT.split(line, "."); 102 vaspVersion = PT.parseInt(tokens[1]); 103 } 104 105 @Override finalizeSubclassReader()106 protected void finalizeSubclassReader() throws Exception { 107 setSymmetry(); 108 } 109 110 111 112 private Lst<String> elementNames = new Lst<String>(); 113 readElementNames()114 private void readElementNames() throws Exception { 115 // was: TITEL = PAW_PBE Al 04Jan2001 116 117 // now: POTCAR: PAW_PBE O 08Apr2002 118 // or: POTCAR: Ce: PAW 5s5p 6s5d4f valence 119 // or: POTCAR: PAW_PBE H_h 07Sep2000 120 121 line = PT.rep(line, " _ ", "_"); // Version 4 OUTCAR_H6Al2_2.dat has PAW _ PBE 122 String[] tokens = getTokens(); 123 int pt = tokens[1].indexOf(":"); 124 String name = (pt >= 0 ? tokens[1].substring(0, pt) : tokens[2]); 125 elementNames.addLast(name); 126 } 127 128 /* 129 Dimension of arrays: 130 k-Points NKPTS = 10 number of bands NBANDS= 16 131 number of dos NEDOS = 301 number of ions NIONS = 8 132 non local maximal LDIM = 4 non local SUM 2l+1 LMDIM = 8 133 total plane-waves NPLWV = 74088 134 max r-space proj IRMAX = 5763 max aug-charges IRDMAX= 42942 135 dimension x,y,z NGX = 42 NGY = 42 NGZ = 42 136 dimension x,y,z NGXF= 84 NGYF= 84 NGZF= 84 137 support grid NGXF= 84 NGYF= 84 NGZF= 84 138 ions per type = 6 2*/ 139 readAtomCountAndSetNames()140 private void readAtomCountAndSetNames() throws Exception { 141 int[] numofElement = new int[100]; 142 // readLine(); 143 String[] tokens = PT.getTokens(line.substring(line.indexOf("=") + 1)); 144 ac = 0; 145 for (int i = 0; i < tokens.length; i++) 146 ac += (numofElement[i] = parseIntStr(tokens[i])); 147 //this is to reconstruct the atomMappedarray containing the atom 148 atomNames = new String[ac]; 149 int nElements = elementNames.size(); 150 for (int pt = 0, i = 0; i < nElements; i++) 151 for (int j = 0; j < numofElement[i] && pt < ac; j++) 152 atomNames[pt++] = elementNames.get(i); 153 } 154 155 /*direct lattice vectors reciprocal lattice vectors 156 1.850000000 1.850000000 0.000000000 0.270270270 0.270270270 -0.270270270 157 0.000000000 1.850000000 1.850000000 -0.270270270 0.270270270 0.270270270 158 1.850000000 0.000000000 1.850000000 0.270270270 -0.270270270 0.270270270*/ 159 readUnitCellVectors()160 private void readUnitCellVectors() throws Exception { 161 if (asc.ac > 0) { 162 setSymmetry(); 163 asc.newAtomSet(); 164 setAtomSetInfo(); 165 } 166 float[] f = new float[3]; 167 for (int i = 0; i < 3; i++) 168 addExplicitLatticeVector(i, fillFloatArray(fixMinus(rd()), 0, f), 0); 169 } 170 fixMinus(String line)171 private String fixMinus(String line) { 172 return PT.rep(line, "-", " -"); 173 } 174 175 setSymmetry()176 private void setSymmetry() throws Exception { 177 applySymmetryAndSetTrajectory(); 178 setSpaceGroupName("P1"); 179 setFractionalCoordinates(false); 180 } 181 182 /* 183 position of ions in fractional coordinates (direct lattice) 184 0.87800000 0.62200000 0.25000000 185 0.25000000 0.87800000 0.62200000 186 0.62200000 0.25000000 0.87800000 187 0.12200000 0.37800000 0.75000000 188 0.75000000 0.12200000 0.37800000 189 0.37800000 0.75000000 0.12200000 190 0.00000000 0.00000000 0.00000000 191 0.50000000 0.50000000 0.50000000 192 193 position of ions in cartesian coordinates 194 */ 195 ///This is the initial geometry not the geometry during the geometry dump readInitialCoordinates()196 private void readInitialCoordinates() throws Exception { 197 int counter = 0; 198 while (rd() != null && line.length() > 10) { 199 addAtomXYZSymName(PT.getTokens(fixMinus(line)), 200 0, null, atomNames[counter++]); 201 } 202 asc.setAtomSetName("Initial Coordinates"); 203 } 204 205 /* 206 POSITION TOTAL-FORCE (eV/Angst) 207 ----------------------------------------------------------------------------------- 208 1.14298 0.83102 6.90311 -0.003060 0.001766 0.000000 209 -1.29117 0.57434 6.90311 0.000000 -0.003533 0.000000 210 0.14819 -1.40536 6.90311 0.003060 0.001766 0.000000 211 -1.14298 -0.83102 4.93079 0.003060 -0.001766 0.000000 212 1.29117 -0.57434 4.93079 0.000000 0.003533 0.000000 213 -0.14819 1.40536 4.93079 -0.003060 -0.001766 0.000000 214 0.00000 0.00000 0.00000 0.000000 0.000000 0.000000 215 0.00000 0.00000 5.91695 0.000000 0.000000 0.000000 216 ----------------------------------------------------------------------------------- 217 */ readPOSITION()218 private void readPOSITION() throws Exception { 219 int counter = 0; 220 readLines(1); 221 while (rd() != null && line.indexOf("----------") < 0) 222 addAtomXYZSymName(getTokens(), 0, null, atomNames[counter++]); 223 } 224 225 /* FREE ENERGIE OF THE ION-ELECTRON SYSTEM (eV) 226 --------------------------------------------------- 227 free energy TOTEN = -20.028155 eV 228 229 energy without entropy= -20.028155 energy(sigma->0) = -20.028155 230 */ 231 232 private Double gibbsEnergy, gibbsEntropy; 233 readEnergy()234 private void readEnergy() throws Exception { 235 rd(); 236 String[] tokens = PT.getTokens(rd()); 237 gibbsEnergy = Double.valueOf(Double.parseDouble(tokens[4])); 238 rd(); 239 tokens = PT.getTokens(rd()); 240 /* please double-check: 241 242 entropy T*S EENTRO = -0.01255935 243 eigenvalues EBANDS = 27.21110509 244 atomic energy EATOM = 181.97672381 245 --------------------------------------------------- 246 free energy TOTEN = 35.37614365 eV 247 248 energy without entropy = 35.38870300 energy(sigma->0) = 35.38242333 249 250 * G = H - TS, so TS = H - G 251 * 252 * My reading of this is that TOTEN is G, 253 * "energy without entropy" is H, and so the "T*S" line is 254 * actually -T*S, not T*S. 255 * 256 * Can that be? 257 * 258 * Bob 259 * 260 * 261 */ 262 double enthalpy = Double.parseDouble(tokens[3]); 263 gibbsEntropy = Double.valueOf(enthalpy - gibbsEnergy.doubleValue()); 264 } 265 setAtomSetInfo()266 private void setAtomSetInfo() { 267 if (gibbsEnergy == null) 268 return; 269 asc.setAtomSetEnergy("" + gibbsEnergy, 270 gibbsEnergy.floatValue()); 271 asc.setCurrentModelInfo("Energy", gibbsEnergy); 272 asc.setCurrentModelInfo("Entropy", gibbsEntropy); 273 asc.setInfo("Energy", gibbsEnergy); 274 asc 275 .setInfo("Entropy", gibbsEntropy); 276 asc.setAtomSetName("G = " + gibbsEnergy + " eV, T*S = " 277 + gibbsEntropy + " eV"); 278 } 279 280 private Double electronEne, kinEne, totEne; 281 private float temp; 282 readMdyn()283 private void readMdyn() throws Exception { 284 String[] tokens = getTokens(); 285 rd(); 286 tokens = PT.getTokens(rd()); 287 electronEne = Double.valueOf(Double.parseDouble(tokens[4])); 288 tokens = PT.getTokens(rd()); 289 kinEne = Double.valueOf(Double.parseDouble(tokens[4])); 290 temp = parseFloatStr(tokens[6]); 291 readLines(3); 292 tokens = PT.getTokens(rd()); 293 totEne = Double.valueOf(Double.parseDouble(tokens[4])); 294 setAtomSetInfoMd(); 295 } 296 setAtomSetInfoMd()297 private void setAtomSetInfoMd() { 298 asc.setAtomSetName("Temp. = " + DF.formatDecimal((temp), 2) 299 + " K, Energy = " + totEne + " eV"); 300 asc.setCurrentModelInfo("Energy", totEne); 301 asc.setInfo("Energy", totEne); 302 asc.setCurrentModelInfo("EleEnergy", kinEne); 303 asc.setInfo("EleEnergy", 304 electronEne); 305 asc.setCurrentModelInfo("Kinetic", electronEne); 306 asc.setInfo("Kinetic", kinEne); 307 asc.setCurrentModelInfo("Temperature", 308 DF.formatDecimal((temp), 2)); 309 asc.setInfo("Temperature", 310 DF.formatDecimal((temp), 2)); 311 } 312 313 /* 314 Eigenvectors after division by SQRT(mass) ///This is lines is not there in VASP5 315 316 Eigenvectors and eigenvalues of the dynamical matrix 317 ---------------------------------------------------- 318 319 320 1 f = 61.880092 THz 388.804082 2PiTHz 2064.097613 cm-1 255.915580 meV 321 X Y Z dx dy dz 322 1.154810 0.811470 6.887910 0.142328 0.284744 0.218505 323 -1.280150 0.594360 6.887910 -0.325017 -0.005914 0.256184 324 0.125350 -1.405820 6.887910 0.147813 -0.308571 0.250139 325 -1.154810 -0.811470 4.919940 -0.142322 -0.284738 -0.218502 326 1.280150 -0.594360 4.919940 0.325007 0.005917 -0.256181 327 -0.125350 1.405820 4.919940 -0.147808 0.308554 -0.250128 328 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 329 0.000000 0.000000 5.903930 0.000000 0.000000 0.000000 330 331 2 f = 56.226671 THz 353.282596 2PiTHz 1875.519821 cm-1 232.534905 meV 332 X Y Z dx dy dz 333 1.154810 0.811470 6.887910 0.005057 -0.193034 -0.074407 334 -1.280150 0.594360 6.887910 -0.215700 -0.046660 0.255430 335 0.125350 -1.405820 6.887910 -0.155522 0.437794 -0.342377 336 -1.154810 -0.811470 4.919940 0.005053 -0.193045 -0.074415 337 1.280150 -0.594360 4.919940 -0.215675 -0.046662 0.255412 338 -0.125350 1.405820 4.919940 -0.155521 0.437777 -0.342368 339 0.000000 0.000000 0.000000 0.015694 0.002977 0.005811 340 0.000000 0.000000 5.903930 0.011571 -0.016988 0.006533 341 342 343 344 */ readFrequency()345 private void readFrequency() throws Exception { 346 347 int pt = asc.iSet; 348 asc.baseSymmetryAtomCount = ac; 349 350 if (vaspVersion >= 5) { 351 readLines(3); 352 } else { 353 discardLinesUntilContains("Eigenvectors after division by SQRT(mass)"); 354 readLines(5); 355 } 356 357 boolean[] ignore = new boolean[1]; 358 while (rd() != null 359 && (line.contains("f = ") || line.contains("f/i= "))) { 360 applySymmetryAndSetTrajectory(); 361 int iAtom0 = asc.ac; 362 cloneLastAtomSet(ac, null); 363 if (!ignore[0]) { 364 asc.iSet = ++pt; 365 asc.setAtomSetFrequency(vibrationNumber, null, 366 null, line.substring(line.indexOf("2PiTHz") + 6, line.indexOf("c") - 1) 367 .trim(), null); 368 } 369 rd(); 370 fillFrequencyData(iAtom0, ac, ac, ignore, true, 35, 12, 371 null, 0, null); 372 rd(); 373 } 374 } 375 } 376