Is there a straightforward way to do so? I've used the naive approach of defining the p_state as an array with the coefficients but it does not work.
Here's my code:
Python: Select all
def H_model(D, spin, lmbd, p):
# Define model parameters
model_params = dict(
L=2, # Number of sites
S=spin, # Spin value (e.g., 0.5 for spin-1/2)
Jz=-4, # Coupling in the z-direction
hx=lmbd*2, # Transverse field in the x-direction
Jzp=-p*4, # Next-nearest-neighbor coupling in the z-direction
Jx=-p*lmbd*4, # Next-nearest-neighbor coupling in the x-direction
bc_MPS='infinite', # Boundary conditions for MPS
conserve=None # No conservation laws
)
# Initialize the model
M = SpinChainNNN(model_params)
# Define the |+> state as a list of local states
# For spin-1/2, |+> corresponds to the state vector [1/sqrt(2), 1/sqrt(2)]
plus_state_vector = np.array([1.0/np.sqrt(2), 1.0/np.sqrt(2)])
# Initialize the MPS in the |+> state for all sites
# The p_state argument must be a list of length L, where each element is the state vector for that site
p_state = np.array([plus_state_vector, plus_state_vector])
psi = MPS.from_product_state(
M.lat.mps_sites(), # Use the lattice sites from the model
p_state, # Pass the list of state vectors
bc=model_params['bc_MPS'] # Use the same boundary conditions as the model
)
return M, psi