Coverage for /builds/ericyuan00000/ase/ase/optimize/sciopt.py: 68.22%
107 statements
« prev ^ index » next coverage.py v7.5.3, created at 2025-06-18 01:20 +0000
1# fmt: off
3from typing import IO, Optional, Union
5import numpy as np
6import scipy.optimize as opt
8from ase import Atoms
9from ase.optimize.optimize import Optimizer
12class Converged(Exception):
13 pass
16class OptimizerConvergenceError(Exception):
17 pass
20class SciPyOptimizer(Optimizer):
21 """General interface for SciPy optimizers
23 Only the call to the optimizer is still needed
24 """
26 def __init__(
27 self,
28 atoms: Atoms,
29 logfile: Union[IO, str] = '-',
30 trajectory: Optional[str] = None,
31 callback_always: bool = False,
32 alpha: float = 70.0,
33 **kwargs,
34 ):
35 """Initialize object
37 Parameters
38 ----------
39 atoms: :class:`~ase.Atoms`
40 The Atoms object to relax.
42 trajectory: str
43 Trajectory file used to store optimisation path.
45 logfile: file object or str
46 If *logfile* is a string, a file with that name will be opened.
47 Use '-' for stdout.
49 callback_always: bool
50 Should the callback be run after each force call (also in the
51 linesearch)
53 alpha: float
54 Initial guess for the Hessian (curvature of energy surface). A
55 conservative value of 70.0 is the default, but number of needed
56 steps to converge might be less if a lower value is used. However,
57 a lower value also means risk of instability.
59 kwargs : dict, optional
60 Extra arguments passed to
61 :class:`~ase.optimize.optimize.Optimizer`.
63 """
64 restart = None
65 Optimizer.__init__(self, atoms, restart, logfile, trajectory, **kwargs)
66 self.force_calls = 0
67 self.callback_always = callback_always
68 self.H0 = alpha
69 self.max_steps = 0
71 def x0(self):
72 """Return x0 in a way SciPy can use
74 This class is mostly usable for subclasses wanting to redefine the
75 parameters (and the objective function)"""
76 return self.optimizable.get_positions().reshape(-1)
78 def f(self, x):
79 """Objective function for use of the optimizers"""
80 self.optimizable.set_positions(x.reshape(-1, 3))
81 # Scale the problem as SciPy uses I as initial Hessian.
82 return self.optimizable.get_potential_energy() / self.H0
84 def fprime(self, x):
85 """Gradient of the objective function for use of the optimizers"""
86 self.optimizable.set_positions(x.reshape(-1, 3))
87 self.force_calls += 1
89 if self.callback_always:
90 self.callback(x)
92 # Remember that forces are minus the gradient!
93 # Scale the problem as SciPy uses I as initial Hessian.
94 return - self.optimizable.get_forces().reshape(-1) / self.H0
96 def callback(self, x):
97 """Callback function to be run after each iteration by SciPy
99 This should also be called once before optimization starts, as SciPy
100 optimizers only calls it after each iteration, while ase optimizers
101 call something similar before as well.
103 :meth:`callback`() can raise a :exc:`Converged` exception to signal the
104 optimisation is complete. This will be silently ignored by
105 :meth:`run`().
106 """
107 if self.nsteps < self.max_steps:
108 self.nsteps += 1
109 f = self.optimizable.get_forces()
110 self.log(f)
111 self.call_observers()
112 if self.converged(f):
113 raise Converged
115 def run(self, fmax=0.05, steps=100000000):
116 self.fmax = fmax
118 try:
119 # As SciPy does not log the zeroth iteration, we do that manually
120 if self.nsteps == 0:
121 self.log()
122 self.call_observers()
124 self.max_steps = steps + self.nsteps
126 # Scale the problem as SciPy uses I as initial Hessian.
127 self.call_fmin(fmax / self.H0, steps)
128 except Converged:
129 pass
130 return self.converged()
132 def dump(self, data):
133 pass
135 def load(self):
136 pass
138 def call_fmin(self, fmax, steps):
139 raise NotImplementedError
142class SciPyFminCG(SciPyOptimizer):
143 """Non-linear (Polak-Ribiere) conjugate gradient algorithm"""
145 def call_fmin(self, fmax, steps):
146 output = opt.fmin_cg(self.f,
147 self.x0(),
148 fprime=self.fprime,
149 # args=(),
150 gtol=fmax * 0.1, # Should never be reached
151 norm=np.inf,
152 # epsilon=
153 maxiter=steps,
154 full_output=1,
155 disp=0,
156 # retall=0,
157 callback=self.callback)
158 warnflag = output[-1]
159 if warnflag == 2:
160 raise OptimizerConvergenceError(
161 'Warning: Desired error not necessarily achieved '
162 'due to precision loss')
165class SciPyFminBFGS(SciPyOptimizer):
166 """Quasi-Newton method (Broydon-Fletcher-Goldfarb-Shanno)"""
168 def call_fmin(self, fmax, steps):
169 output = opt.fmin_bfgs(self.f,
170 self.x0(),
171 fprime=self.fprime,
172 # args=(),
173 gtol=fmax * 0.1, # Should never be reached
174 norm=np.inf,
175 # epsilon=1.4901161193847656e-08,
176 maxiter=steps,
177 full_output=1,
178 disp=0,
179 # retall=0,
180 callback=self.callback)
181 warnflag = output[-1]
182 if warnflag == 2:
183 raise OptimizerConvergenceError(
184 'Warning: Desired error not necessarily achieved '
185 'due to precision loss')
188class SciPyGradientlessOptimizer(Optimizer):
189 """General interface for gradient less SciPy optimizers
191 Only the call to the optimizer is still needed
193 Note: If you redefine x0() and f(), you don't even need an atoms object.
194 Redefining these also allows you to specify an arbitrary objective
195 function.
197 XXX: This is still a work in progress
198 """
200 def __init__(
201 self,
202 atoms: Atoms,
203 logfile: Union[IO, str] = '-',
204 trajectory: Optional[str] = None,
205 callback_always: bool = False,
206 **kwargs,
207 ):
208 """Initialize object
210 Parameters
211 ----------
212 atoms: :class:`~ase.Atoms`
213 The Atoms object to relax.
215 trajectory: str
216 Trajectory file used to store optimisation path.
218 logfile: file object or str
219 If *logfile* is a string, a file with that name will be opened.
220 Use '-' for stdout.
222 callback_always: bool
223 Should the callback be run after each force call (also in the
224 linesearch)
226 alpha: float
227 Initial guess for the Hessian (curvature of energy surface). A
228 conservative value of 70.0 is the default, but number of needed
229 steps to converge might be less if a lower value is used. However,
230 a lower value also means risk of instability.
232 kwargs : dict, optional
233 Extra arguments passed to
234 :class:`~ase.optimize.optimize.Optimizer`.
236 """
237 restart = None
238 Optimizer.__init__(self, atoms, restart, logfile, trajectory, **kwargs)
239 self.function_calls = 0
240 self.callback_always = callback_always
242 def x0(self):
243 """Return x0 in a way SciPy can use
245 This class is mostly usable for subclasses wanting to redefine the
246 parameters (and the objective function)"""
247 return self.optimizable.get_positions().reshape(-1)
249 def f(self, x):
250 """Objective function for use of the optimizers"""
251 self.optimizable.set_positions(x.reshape(-1, 3))
252 self.function_calls += 1
253 # Scale the problem as SciPy uses I as initial Hessian.
254 return self.optimizable.get_potential_energy()
256 def callback(self, x):
257 """Callback function to be run after each iteration by SciPy
259 This should also be called once before optimization starts, as SciPy
260 optimizers only calls it after each iteration, while ase optimizers
261 call something similar before as well.
262 """
263 # We can't assume that forces are available!
264 # f = self.optimizable.get_forces()
265 # self.log(f)
266 self.call_observers()
267 # if self.converged(f):
268 # raise Converged
269 self.nsteps += 1
271 def run(self, ftol=0.01, xtol=0.01, steps=100000000):
272 self.xtol = xtol
273 self.ftol = ftol
274 # As SciPy does not log the zeroth iteration, we do that manually
275 self.callback(None)
276 try:
277 # Scale the problem as SciPy uses I as initial Hessian.
278 self.call_fmin(xtol, ftol, steps)
279 except Converged:
280 pass
281 return self.converged()
283 def dump(self, data):
284 pass
286 def load(self):
287 pass
289 def call_fmin(self, xtol, ftol, steps):
290 raise NotImplementedError
293class SciPyFmin(SciPyGradientlessOptimizer):
294 """Nelder-Mead Simplex algorithm
296 Uses only function calls.
298 XXX: This is still a work in progress
299 """
301 def call_fmin(self, xtol, ftol, steps):
302 opt.fmin(self.f,
303 self.x0(),
304 # args=(),
305 xtol=xtol,
306 ftol=ftol,
307 maxiter=steps,
308 # maxfun=None,
309 # full_output=1,
310 disp=0,
311 # retall=0,
312 callback=self.callback)
315class SciPyFminPowell(SciPyGradientlessOptimizer):
316 """Powell's (modified) level set method
318 Uses only function calls.
320 XXX: This is still a work in progress
321 """
323 def __init__(self, *args, **kwargs):
324 """Parameters:
326 direc: float
327 How much to change x to initially. Defaults to 0.04.
328 """
329 direc = kwargs.pop('direc', None)
330 SciPyGradientlessOptimizer.__init__(self, *args, **kwargs)
332 if direc is None:
333 self.direc = np.eye(len(self.x0()), dtype=float) * 0.04
334 else:
335 self.direc = np.eye(len(self.x0()), dtype=float) * direc
337 def call_fmin(self, xtol, ftol, steps):
338 opt.fmin_powell(self.f,
339 self.x0(),
340 # args=(),
341 xtol=xtol,
342 ftol=ftol,
343 maxiter=steps,
344 # maxfun=None,
345 # full_output=1,
346 disp=0,
347 # retall=0,
348 callback=self.callback,
349 direc=self.direc)