API reference

GeometryOptimization.GeoOptDefaultCallbackType

Callback producing a convergence table summarising the geometry optimisation convergence. If always_show_header=true the header is shown in each iteration. This is helpful if the calculator produces output as well.

source
SciMLBase.OptimizationProblemMethod
Optimization.OptimizationProblem(system, calculator geoopt_state; kwargs...)

Turn system, calculator and geoopt_state::GeometryOptimizationState into a SciML-compatible OptimizationProblem. Note that the system is not updated automatically and that internally atomic units are used.

source
GeometryOptimization.clamp_atomsMethod
clamp_atoms(
    system,
    clamped_indexes::Union{Nothing, AbstractVector{<:Integer}}
) -> AtomsBase.FlexibleSystem{_A, _B, AtomsBase.PeriodicCell{_A1, _B1}} where {_A, _B, _A1, _B1}

Clamp given atoms in the system. Clamped atoms are fixed and their positions will not be optimized. The atoms to be clamped should be given as a list of indices corresponding to their positions in the system.

Experimental

This is this a very experimental interface and will likely change in the future.

source
GeometryOptimization.convert_to_updatableMethod
convert_to_updatable(system) -> Any

Convert the input system to a kind of system we can work with, i.e. one where we have the particles keyword argument for setting new particles and positions.

source
GeometryOptimization.minimize_energy!Function
minimize_energy!(
    system,
    calculator;
    ...
) -> NamedTuple{(:system, :converged, :energy, :forces, :virial, :state, :stats, :alg, :optimres), <:Tuple{AtomsBase.FlexibleSystem{_A, _B, AtomsBase.PeriodicCell{_A1, _B1}} where {_A, _B, _A1, _B1}, Any, Any, Any, Any, Any, Any, Union{Optim.ConjugateGradient{Float64, Nothing, Optim.var"#33#35", LineSearches.InitialHagerZhang{Float64}, LineSearches.BackTracking{Float64, Int64}}, Optim.Fminbox{Optim.ConjugateGradient{Float64, Nothing, Optim.var"#33#35", LineSearches.InitialHagerZhang{Float64}, LineSearches.BackTracking{Float64, Int64}}, Float64, Optim.var"#51#53"}}, Union{SciMLBase.OptimizationSolution{_A, _B, var"#s182", OptimizationBase.OptimizationCache{F, RC, LB, UB, LC, UC, S, O, D, P, C, M}, Optim.ConjugateGradient{Float64, Nothing, Optim.var"#33#35", LineSearches.InitialHagerZhang{Float64}, LineSearches.BackTracking{Float64, Int64}}, _C, Optim.MultivariateOptimizationResults{O1, Tx, Tc, Tf, M1, Tls, Tsb}, SciMLBase.OptimizationStats} where {_A, _B, T, N, var"#s182"<:AbstractArray{T, N}, F, RC, LB, UB, LC, UC, S, O, D, P, C, M, _C, O1, Tx, Tc, Tf, M1, Tls, Tsb}, SciMLBase.OptimizationSolution{_A, _B, _C, OptimizationBase.OptimizationCache{F, RC, LB, UB, LC, UC, S, O, D, P, C, M}, Optim.Fminbox{Optim.ConjugateGradient{Float64, Nothing, Optim.var"#33#35", LineSearches.InitialHagerZhang{Float64}, LineSearches.BackTracking{Float64, Int64}}, Float64, Optim.var"#51#53"}} where {_A, _B, _C, F, RC, LB, UB, LC, UC, S, O, D, P, C, M}}}}
minimize_energy!(
    system,
    calculator,
    solver;
    kwargs...
) -> NamedTuple{(:system, :converged, :energy, :forces, :virial, :state, :stats, :alg, :optimres), <:Tuple{AtomsBase.FlexibleSystem{_A, _B, AtomsBase.PeriodicCell{_A1, _B1}} where {_A, _B, _A1, _B1}, Vararg{Any, 8}}}

Minimise the energy of a system using the specified calculator. For now only optimises atomic positions. Returns a named tuple including the optimised system as first entry. Under the hood this constructs an Optimization.OptimizationProblem and uses Optimization.jl to solve it using the passed solver.

Typical arguments passed as solver are GeometryOptimization.Autoselect() (the default), GeometryOptimization.OptimLBFGS(), GeometryOptimization.OptimCG(), GeometryOptimization.OptimSD(). These automatically choose some heuristics for setting up the solvers, which we found to work well in practice. Beyond that any other solver compatible with Optimization.jl can also be employed here.

Keyword arguments:

  • maxiters: Maximal number of iterations
  • maxtime: Maximal allowed runtime (in seconds)
  • tol_energy: Tolerance in the energy to stop the minimisation (all tol_* need to be satisfied)
  • tol_forces: Tolerance in the force to stop the minimisation (all tol_* need to be satisfied)
  • tol_virial: Tolerance in the virial to stop the minimisation (all tol_* need to be satisfied)
  • maxstep: Maximal step size (in AU or length units) to be taken in a single optimisation step (not supported for all solvers)
  • verbosity: Printing level. The idea is that 0 is silent, 1 displays the optimisation progress and ≥ 2 starts displaying things from the calculator as well (e.g SCF iterations).
  • callback: A custom callback, which obtains the pair (optimization_state, geoopt_state) and is expected to return false (continue iterating) or true (halt iterations). Note that specifying this overwrites the default printing callback. The calculation thus becomes silent unless a GeoOptDefaultCallback is included in the callback.
  • kwargs: All other keyword arguments are passed to the call to solve. Note, that if special kwargs should be passed to the Optimization.OptimizationProblem the user needs to setup the problem manually (e.g. OptimizationProblem(system, calculator))
source
GeometryOptimization.update_not_clamped_positionsMethod
update_not_clamped_positions(
    system,
    positions::AbstractVector{<:Union{Unitful.Quantity{T, 𝐋, U}, Unitful.Level{L, S, Unitful.Quantity{T, 𝐋, U}} where {L, S}} where {T, U}}
) -> AtomsBase.FlexibleSystem{_A, _B, AtomsBase.PeriodicCell{_A1, _B1}} where {_A, _B, _A1, _B1}

Creates a new system based on $system$ where the non clamped positions are updated to the ones provided (in the order in which they appear in the system).

source
GeometryOptimization.update_positionsMethod
update_positions(
    system,
    positions::AbstractVector{<:AbstractVector{<:Union{Unitful.Quantity{T, 𝐋, U}, Unitful.Level{L, S, Unitful.Quantity{T, 𝐋, U}} where {L, S}} where {T, U}}}
) -> AtomsBase.FlexibleSystem{_A, _B, AtomsBase.PeriodicCell{_A1, _B1}} where {_A, _B, _A1, _B1}

Creates a new system based on $system$ but with atoms positions updated to the ones provided.

source