# This script's purpose is to demonstrate and explain the following functionalities of the SciPAL SymMps toolkit
# 1. create a spin-1/2 lattice with L=100 Spins
# 1a. create a fully polarized state |\psi> = |\uparrow \cdots \uparrow>
# 1b. create a random near-vacuum state by repeated application of \sum_{i}c_i S^+_i |\downarrow \cdots \downarrow> with random weights c_i
# 2. apply a local perturbation to a state flipping a spin at lattice site i=50
#
# Dependencies:
#   none
#
# Input arguments:
#   path to compiled SciPAL-SymMps toolkit
#   type of created state (full|vacuum)

# Force OMP-/MKL thread numbers to be one otherwise thread spawning will slow down simulations
export OMP_NUM_THREADS=1
export MKL_NUM_THREADS=1

# clear log folder to avoid unnecessary blow-up of output files (each executable will append its output onto an existing output file)
rm log/*.output

# remove already existing lattice files to avoid loading of lattice data which we may not want to have
rm lattice.*

# We exit script if an executable returns with error code
set -eu

# just define a shortcut for the passed location of the toolbox ($1 must only define a path to the folder of the tool-kit containing bin (compiled binaries) and link (symbolic links with some reasonable names), e.g., in any singularity container you may pass /opt/symmps/gnu or /opt/symmps/mkl
EXEC_DIR=$1

# specify type of created state
# TYPE=full   : a fully polarized state |\uparrow \cdots \uparrow> 
# TYPE=vacuum : a random near-vacuum state by repeated application of \sum_{i}c_i S^+_i |\downarrow \cdots \downarrow> with random weights c_i
TYPE=$2

# create a lattice file containing defining a subspace of tensor-product Hilbert space with the following shortoptions:
# 1. -L defines number of lattice sites, i.e., copies of local Hilbert spaces
# 2. -F defines fundamental datatype of the Hilbert space to be double
# 3. -b defines used basis set, here we choose 'custom' basis set which contains the Spin-Basis, to get more information on any basis set <basis_set> run sym-create-lattice -b <basis_set> -G
# 4. -d defines dimension of the local Hilbert space, we want to describe S=1/2 degrees of freedom thus d=2
# 5. -g defines local symmetry generator which is used to construct the global symmetry generator in this case S^z=\sum_{i=1}^{L}S^z_i, note that only the local operator needs to be specified and since they all share the same representation in each individual local Hilbert space the index 'i' is dropped
# 6. -t defines target-symmetry sector for the global symmetry generator, in this case S^z=0
# 7. -o defines output filename for the created lattice (there is no particular naming convention, however you may want to encode the fundamental datatype in the filename)
case $TYPE in
	full)
		${EXEC_DIR}/link/sym-create-lattice -L 100 -F d -b custom -d 2 -g Sz -t 50.0 -o lattice.real
		;;
	vacuum)
                ${EXEC_DIR}/link/sym-create-lattice -L 100 -F d -b custom -d 2 -g Sz -t 0.0 -o lattice.real
		;;
esac

# create a state with the following shortoptions:
# 1. -i defines lattice input file containing Hilbert space and local basis
# 2. -t defines type of created mps, choose 'full' means that all sites are initialized in eigenstate of the local symmetry generators with maximum eigenvalue, i.e., s_i=1/2
# 3. -o defines dest state filename
${EXEC_DIR}/link/sym-create-mps -i lattice.real -t $TYPE -o state.symmps

# apply local perturbation by flipping spin at site i=50 with the following shortoptions:
# 1. -i defines lattice input file containing Hulbert space and local basis
# 2. -S defines input file for source state
# 3. -D defines output file for target state
# 4. -O defines MPO identifier, in case you want to apply an MPO from passed lattice insert MPO-identifier here, we want to apply only a local perturbation which can be done by specifying the local operator and the lattice on which it should act with syntax <local_operator>@<site>, the MPO is then created automatically
# 5. -m defines application method, here we choose variational application method which first performs a warmup sweep and then variationally minimizes the l2-distance to the target state <MPO>|\psi>
# 6. -n if type is variational this defines the max. number of sweeps during variational optimization towards the target state
# 7. -e if type is variational this defines the escape precision which is given by the change in the l2-norm || <MPO>|\psi> - |\phi> || between the two last sweeps where |\phi> is the optimized state
# 8. -d defines the maximum discarded weight for the target state
# 9. -w defines the maximum bond dimension of the target state
${EXEC_DIR}/link/sym-apply-mpo -i lattice.real -S state.symmps -D excited.symmps -O 'Sminus@50' -m variational -n 10 -e '1e-12' -d 1e-12 -w 100