Obtaining the final wave function after RBM

[lynkeryn]

My code defines a 1D Heisenberg model and uses RBM to find the ground state.
How do I print the resulting wavefunction (amplitudes of the basis states)?

code:

# Import netket library
import netket as nk

# # Import Json, this will be needed to load log files
import json

# # Helper libraries
import numpy as np
import matplotlib.pyplot as plt
import time
plt.figure(figsize=(15,12))
niter = 150
L = 8
for i in range(6):
  #build state
  g = nk.graph.Hypercube(length=L, n_dim=1, pbc=True)
  hi = nk.hilbert.Spin(s=0.5, total_sz=0, N=g.n_nodes)
  ha = nk.operator.Heisenberg(hilbert=hi, graph=g)
  # RBM ansatz with alpha=1
  ma = nk.models.RBM(alpha=1)
  # Build the sampler
  sa = nk.sampler.MetropolisExchange(hilbert=hi,graph=g)

  # Optimizer
  op = nk.optimizer.Sgd(learning_rate=0.01+0.01*i)

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

  # The variational state
  vs = nk.vqs.MCState(sa, ma, n_samples=1000)

  # The ground-state optimization loop
  gs = nk.VMC(
      hamiltonian=ha,
      optimizer=op,
      preconditioner=sr,
      variational_state=vs)

  start = time.time()
  gs.run(out='RBM', n_iter=niter)
  end = time.time()

  print('### RBM calculation')
  print('Has',vs.n_parameters,'parameters')
  print('The RBM calculation took',end-start,'seconds')


  # import the data from log file
  data=json.load(open("RBM.log"))
    
  # Extract the relevant information
  iters_RBM = data["Energy"]["iters"]
  energy_RBM = data["Energy"]["Mean"]
      
  plt.subplot(2,3,i+1)

  plt.plot(iters_RBM, energy_RBM, color='red', label='Energy (RBM)')
  plt.plot(iters_RBM, [eigvals]*len(iters_RBM),  '--', label = f'Exact = {round(eigvals[0],3)}',)
  plt.title(f'learning_rate={round(0.01+0.01*i,2)}')
  plt.ylabel('Energy')
  plt.xlabel('Iteration')
  # plt.axis([0,iters_RBM[-1],exact_gs_energy-0.03,exact_gs_energy+0.2])
  # plt.axhline(y=exact_gs_energy, xmin=0,
  #                 xmax=iters_RBM[-1], linewidth=2, color='k', label='Exact')
  plt.legend()
plt.show()

[PhilipVinc]

vs.to_array() unpacks the NQS into a dense vector, though this will only work on relatively small Hilbert spaces where the full vector can fit in memory.

[gcarleo]

because you fixed total_sz=0 , the reduced size of the hilbert space is what you are seeing

[PhilipVinc]

(To understand what that means, try looking at all the basis elements in the reduced basis, that is, look at the output of hi.all_states())

[lynkeryn]

Thank you all!

[lynkeryn]

I have one more question: why, having removed the restriction on the total spin, the eigenvalue of the Hamiltonian did not change for me? As a result, it turns out that the model began to predict differently and converges to a different value, while its own value (which I strive for) remains the same @PhilipVinc @gcarleo

[PhilipVinc]

There are several reasons.
One is for example that your sampling is way off (look at the value of the $\hat{R}$ in the expectation value of the Energy, which if larger than 1.05 suggests that your sampling is wrong so everything won't work).

The reason for the sampling being off is because you are initialising your Markov chains with random states from the now unconstrained Hilbert space. So you have some chains with magnetization=0 , some with $m=2$, some with $m=4$... etc, but your sampler is an exchange sampler which preserves the magnetisation, so a chain initialised in the sector with $m=2$ will be stuck there, your gradients being biased and you will never be able to converge and even if you did your energy estimation would be wrong.

You could use a sampler that does not preserve the magnetisation and solve the optimisation in the full space, but that would be harder than solving it in the constrained space...