Skip to content

[WIP] Add Dynamics

Vicentini Filippo requested to merge PhilipVinc/dyn into master

Created by: PhilipVinc

This PR support dynamic as shown in the following example:

import netket as nk
import numpy as np
from netket.dynamics import TimeEvolution

import netket as nk

# 1D Lattice
g = nk.graph.Hypercube(length=3, n_dim=1, pbc=True)

# Hilbert space of spins on the graph
hi = nk.hilbert.Spin(s=1 / 2, N=g.n_nodes)

# Ising spin hamiltonian
ha = nk.operator.Ising(h=1.0, hilbert=hi, graph=g)
Sx = sum([nk.operator.spin.sigmax(hi, i) for i in range(g.n_nodes)])

# RBM Spin Machine
ma = nk.machine.RbmSpin(alpha=1, hilbert=hi)
ma.init_random_parameters(seed=1234, sigma=0.01)

# Metropolis Local Sampling
sa = nk.sampler.MetropolisLocal(ma, n_chains=32)

# Optimizer
op = nk.optimizer.Sgd(ma, learning_rate=0.1)

# Stochastic Reconfiguration
sr = nk.optimizer.SR(ma, diag_shift=0.1)

# Create the optimization driver
gs = nk.Vmc(ha, sampler=sa, optimizer=op, n_samples=1000, sr=sr)

# Run the optimization for 300 iterations
gs.run(n_iter=300, out="test")

# Now that we have the ground state, change the shift to a small value
sr = nk.optimizer.SR(ma, diag_shift=0.001, sparse_maxiter=1000)

tevo = TimeEvolution(
    ha, sampler=sa, sr=sr, n_samples=1000,
)
tevo.solver("rk45", dt=0.01, adaptive=False)

obs = {"Sx": Sx}

tevo.run(2.0, out="test_tevo", obs=obs)

➜ python Examples/Ising1d/tevo.py
100%|██████████| 300/300 [00:37<00:00,  7.95it/s, Energy=-12.7835-0.0022j ± 0.0042 [σ²=0.0175, R̂=1.0001]]]
100%|█| 2.0000000000000013/2.0 [03:16<00:00, 98.22s/it, Energy=-12.7793+0.0013j ± 0.0037 [σ²=0.0177, R̂=1.00]]]

the way this is implemented is that there is a new driver TimeEvolution which when constructed builds internally the driver for the problem you asked (VMC or StadyState, in this case). There is also a method that performs the time-evolution itself.

Changes to abstract_variational_driver are required to correctly handle floating-point steps in the progress bar and intervals for saving observables.

--

I wrap the solvers from scipy, and do some magic in order to disable their adaptivity which normally can't be disabled. I also wrote an Euler solver myself for scipy, because sometimes that's just useful and simple.

I still need to write some tests and there are a few rough edges here and there. I have been using this code for quite a while and would like to upstream it, so I'd like to hear some feedback about it before taking it to the last stretch.

Merge request reports