netket.optimizer.SR#

class netket.optimizer.SR[source]#

Bases: AbstractLinearPreconditioner

Stochastic Reconfiguration or Natural Gradient preconditioner for the gradient.

Constructs the structure holding the parameters for using the Stochastic Reconfiguration/Natural gradient method.

This preconditioner changes the gradient \(\nabla_i E\) such that the preconditioned gradient \(\Delta_j\) solves the system of equations

\[(S_{i,j} + \delta_{i,j}(\epsilon_1 S_{i,i} + \epsilon_2)) \Delta_{j} = \nabla_i E\]

Where \(S\) is the Quantum Geometric Tensor (or Fisher Information Matrix), preconditioned according to the diagonal scale \(\epsilon_1\) (diag_scale) and the diagonal shift \(\epsilon_2\) (diag_shift). The default regularisation takes \(\epsilon_1=0\) and \(\epsilon_2=0.01\).

Depending on the arguments, an implementation is chosen. For details on all possible kwargs check the specific SR implementations in the documentation.

You can also construct one of those structures directly.

Warning

NetKet also has an experimental implementation of the SR preconditioner using the kernel trick, also known as MinSR. This implementation relies on inverting the \(T = X^T X\) matrix, where \(X\) is the Jacobian of wavefunction and is therefore much more efficient than the standard SR for very large numbers of parameters.

Look at netket.experimental.driver.VMC_SRt for more details.

Inheritance
Inheritance diagram of netket.optimizer.SR
__init__(qgt=None, solver=<function cg>, *, diag_shift=0.01, diag_scale=None, solver_restart=False, **kwargs)[source]#

Constructs the structure holding the parameters for using the Stochastic Reconfiguration/Natural gradient method.

Depending on the arguments, an implementation is chosen. For details on all possible kwargs check the specific SR implementations in the documentation.

You can also construct one of those structures directly.

Parameters:
Attributes
diag_shift: Union[Any, Callable[[Union[Array, ndarray, bool, number, float, int]], Union[Array, ndarray, bool, number, float, int]]]#

Diagonal shift added to the S matrix. Can be a Scalar value, an optax schedule or a Callable function.

diag_scale: Union[Any, Callable[[Union[Array, ndarray, bool, number, float, int]], Union[Array, ndarray, bool, number, float, int]], None]#

Diagonal shift added to the S matrix. Can be a Scalar value, an optax schedule or a Callable function.

qgt_constructor: Callable#

The Quantum Geometric Tensor type or a constructor.

qgt_kwargs: dict#

The keyword arguments to be passed to the Geometric Tensor constructor.

solver: SolverT#

Function used to solve the linear system.

solver_restart: bool#

If False uses the last solution of the linear system as a starting point for the solution of the next.

x0: PyTree | None#

Solution of the last linear system solved.

info: Any#

Additional information returned by the solver when solving the last linear system.

Methods
__call__(vstate, gradient, step=None, *args, **kwargs)[source]#

Call self as a function.

Return type:

Any

Parameters:
lhs_constructor(vstate, step=None)[source]#

This method constructs the left-hand side (LHS) operator for the linear system.

Parameters:
replace(**kwargs)[source]#

Replace the values of the fields of the object with the values of the keyword arguments. If the object is a dataclass, dataclasses.replace will be used. Otherwise, a new object will be created with the same type as the original object.

Return type:

TypeVar(P, bound= Pytree)

Parameters:
  • self (P)

  • kwargs (Any)