Coverage for sparc/sparc_parsers/static.py: 87%
120 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-18 16:19 +0000
1"""
2Created on Thu Oct 18 14:16:21 2018
4Ben Comer (Georgia Tech)
6This file has been heavily modified since SPARC 0.1
8TODO: more descriptions about this file io parser
9"""
10from warnings import warn
12import numpy as np
13from ase.units import Bohr, GPa, Hartree
15# Safe wrappers for both string and fd
16from ase.utils import reader, writer
18from ..api import SparcAPI
19from .utils import strip_comments
21# TODO: should allow user to select the api
22# defaultAPI = SparcAPI()
25@reader
26def _read_static(fileobj):
27 """
28 Read the .static file content
30 Each .static file should only host 1 or more images
31 (if socket mode is enabled), but the output may vary
32 a lot depending on the flags (e.g. PRINT_ATOMS, PRINT_FORCES etc)
33 """
34 contents = fileobj.read()
35 data, comments = strip_comments(contents)
36 # Most static files should containe the Atom positions lines
37 # this can happen for both a single- or multi-image file
38 step_bounds = [i for i, x in enumerate(data) if "Atom positions" in x] + [len(data)]
39 raw_static_steps = [
40 data[start:end] for start, end in zip(step_bounds[:-1], step_bounds[1:])
41 ]
42 # In some cases (e.g. PRINT_ATOMS=0), the static file may not contain Atom positions
43 # All existing lines will be regarded as one step
44 if len(raw_static_steps) == 0:
45 raw_static_steps = [data]
46 static_steps = [_read_static_step(step) for step in raw_static_steps]
47 return {"static": static_steps}
50def _read_static_block(raw_block):
51 """Parse ONE static data block, this will return dict with keys:
52 {"name": PARAM, "value": value}
54 Arguments:
55 raw_block: un-parsed block as list of strings
56 """
57 header, body = raw_block[0], raw_block[1:]
58 header_name, header_rest = header.split(":")
59 if len(header_rest.strip()) > 0:
60 body = [header_rest.strip()] + body
62 try:
63 value = np.genfromtxt(body, dtype=float)
64 if np.isnan(value).any():
65 warn(
66 (
67 f"Field contains data that are not parsable by numpy!\n"
68 f"Contents are: {body}"
69 )
70 )
71 except Exception:
72 value = body
74 # name = None
75 if "Total free energy" in header_name:
76 name = "free energy"
77 elif "Atomic forces" in header_name:
78 name = "forces"
79 elif "Net magnetization" in header_name:
80 name = "net_magnetization"
81 elif "Atomic magnetization" in header_name:
82 name = "atomic_magnetization"
83 if type(value) is list:
84 try:
85 value = np.array(value)
86 except:
87 print("Cannot convert atomic magnetization to array")
88 elif "Stress (GPa)" in header_name:
89 name = "stress"
90 elif "Stress equiv." in header_name:
91 name = "stress_equiv"
92 elif "Stress (Ha/Bohr)" in header_name:
93 name = "stress_1d"
94 elif "Stress (Ha/Bohr**2)" in header_name:
95 name = "stress_2d"
96 elif "Fractional coordinates" in header_name:
97 # Fractional coordinates of Si -- > name=coord_frac symbol="Si"
98 name = "coord_frac"
99 symbol = header_name.split("of")[1].strip()
100 clean_array = value.reshape((-1, 3))
101 value = {"value": clean_array, "symbol": symbol}
102 # Exclusive to the socket mode
103 elif "Lattice (Bohr)" in header_name:
104 name = "lattice"
105 value = value.reshape((3, 3))
106 else:
107 name = header_name.strip()
109 return {"name": name, "value": value}
112def _read_static_step(step):
113 """Parse all the lines in one step and compose a dict containing sanitized blocks
114 Args:
115 step (list): Lines of raw lines in one step
116 """
117 separator = "*" * 60 # Make the separator long enough
118 # Clean up boundary lines
119 data = [
120 line
121 for line in step
122 if ("Atom positions" not in line) and (separator not in line)
123 ]
124 block_bounds = [i for i, x in enumerate(data) if ":" in x] + [len(data)]
125 raw_blocks = [
126 data[start:end] for start, end in zip(block_bounds[:-1], block_bounds[1:])
127 ]
128 # import pdb; pdb.set_trace()
129 block_dict = {}
130 coord_dict = {}
131 block_contents = [_read_static_block(block) for block in raw_blocks]
132 for bc in block_contents:
133 name, raw_value = bc["name"], bc["value"]
134 # Coord frac needs to be collected in all positions
135 if name == "coord_frac":
136 value = None
137 symbol, coord_f = raw_value["symbol"], raw_value["value"]
138 pos_count = len(coord_f)
139 if coord_dict == {}:
140 coord_dict.update(
141 {
142 "symbols": [
143 symbol,
144 ]
145 * pos_count,
146 "coord_frac": coord_f,
147 }
148 )
149 else:
150 coord_dict["symbols"] += [
151 symbol,
152 ] * pos_count
153 # import pdb; pdb.set_trace()
154 coord_dict["coord_frac"] = np.vstack(
155 [coord_dict["coord_frac"], coord_f]
156 )
158 elif name == "free energy":
159 value = raw_value * Hartree
160 elif name == "forces":
161 value = raw_value * Hartree / Bohr
162 elif name == "atomic_magnetization":
163 value = raw_value
164 elif name == "net_magnetization":
165 value = raw_value
166 elif name == "stress":
167 # Stress is in eV/Ang^3, may need to convert to Virial later when cell is known
168 # For low-dimension stress info, use stress_equiv
169 stress_ev_a3 = raw_value * GPa
170 if stress_ev_a3.shape != (3, 3):
171 raise ValueError("Stress from static file is not a 3x3 matrix!")
172 # make the stress in voigt notation
173 # TODO: check the order!
174 value = np.array(
175 [
176 stress_ev_a3[0, 0],
177 stress_ev_a3[1, 1],
178 stress_ev_a3[2, 2],
179 stress_ev_a3[1, 2],
180 stress_ev_a3[0, 2],
181 stress_ev_a3[0, 1],
182 ]
183 )
184 elif name == "stress_equiv":
185 # Only store the size up to the max. periodic directions,
186 # let the atom parser decide how to transform the matrix
187 value = raw_value * GPa
188 elif name == "stress_1d":
189 value = raw_value * Hartree / Bohr
190 elif name == "stress_2d":
191 value = raw_value * Hartree / (Bohr**2)
192 elif name == "lattice":
193 value = raw_value * Bohr
195 # Non-frac coord
196 if value is not None:
197 block_dict[name] = value
198 # Finally, update the atomic positions
199 # TODO: should we keep a default?
200 if coord_dict != {}:
201 block_dict["atoms"] = coord_dict
203 return block_dict
206def _add_cell_info(static_steps, cell=None):
207 """Use the cell information to convert positions
208 if lattice exists in each step, use it to convert coord_frac
209 else use the external cell (for example from inpt file)
210 Args:
211 static_steps: raw list of steps
212 cell: external lattice information
213 """
214 new_steps = []
215 for step in static_steps:
216 new_step = step.copy()
217 if "lattice" in step:
218 lat = step["lattice"]
219 elif cell is not None:
220 lat = cell
221 else:
222 lat = None
224 if (lat is not None) and (step.get("atoms", None) is not None):
225 coord_frac = new_step["atoms"]["coord_frac"]
226 coord_cart = np.dot(coord_frac, lat)
227 new_step["atoms"]["coord"] = coord_cart
228 new_steps.append(new_step)
229 return new_steps
232@writer
233def _write_static(
234 fileobj,
235 data_dict,
236):
237 raise NotImplementedError("Writing static file from SPARC-X-API not supported!")