from sympt import *
import sympy as sp
import numpy as np
import matplotlib.pyplot as plt
from tqdm import tqdm

cmap = plt.get_cmap('inferno')
plt.rc('text', usetex = True)
plt.rc('font', family = 'serif', weight = 'normal')
plt.rc('font', size = 18)
color_list = [cmap(i) for i in np.linspace(0, 1, 16)]
fntsize = 20

SymPT Mask

def letter(num):
    # Define the patterns for each letter
    patterns = {
        0: np.array([
            [1, 1, 1],
            [1, 0, 0],
            [1, 1, 1],
            [0, 0, 1],
            [1, 1, 1]
        ]),  # Capital "S"
        1: np.array([
            [0, 0, 0, 0],
            [1, 0, 0, 1],
            [0, 1, 1, 1],
            [0, 0, 0, 1],
            [1, 1, 1, 1]
        ]),  # Lowercase "y"
       2: np.array([
            [0, 0, 0, 0, 0],
            [1, 0, 0, 0, 1],
            [1, 1, 0, 1, 1],
            [1, 0, 1, 0, 1],
            [1, 0, 0, 0, 1]
        ]), # Lowercase "m"
        3: np.array([
            [1, 1, 1],
            [1, 0, 1],
            [1, 1, 1],
            [1, 0, 0],
            [1, 0, 0]
        ]),  # Capital "P"
        4: np.array([
            [1, 1, 1],
            [0, 1, 0],
            [0, 1, 0],
            [0, 1, 0],
            [0, 1, 0]
        ]),  # Capital "T"
        }

    # Return the correct pattern or raise an error if the input is invalid
    if num in patterns:
        return patterns[num]
    else:
        raise ValueError("Input must be a number from 0 to 3.")
space = np.zeros((5,1))
zero_b = np.zeros((5,5)) 
last_pad = np.zeros((12,37))

SymPT = np.block(
    [
        [zero_b, zero_b, zero_b, zero_b, zero_b, zero_b, space, space, space, space, space, space, space],
        [zero_b, zero_b, zero_b, letter(0),space,letter(1),space,letter(2),space,letter(3),space,letter(4)], 
        [zero_b, zero_b, zero_b, zero_b, zero_b, zero_b, space, space, space, space, space, space, space],
        [zero_b, zero_b, zero_b, zero_b, zero_b, zero_b, space, space, space, space, space, space, space],
        [zero_b, zero_b, zero_b, zero_b, zero_b, zero_b, space, space, space, space, space, space, space],
        [last_pad]
    ]
)
SymPT = SymPT + SymPT.T
mask = Block(SymPT)
plt.imshow(SymPT)
SymPT.shape

(37, 37)

sympt mask

Hamiltonian

I = SymPT.shape[0]
J = SymPT.shape[1]
omegas = [RDSymbol(f"\\omega_{i}", order = 0, real = True) for i in range(I)]
lambdas = []

H = sp.zeros(I, J)
for i in range(I):
    H[i, i] = omegas[i]
    for j in range(i+1, J):  # Loop only over upper triangular part (including diagonal)
        lambdas.append(RDSymbol(f"\\lambda_{i}_{j}", order = 1, real = True))
        H[i, j] = lambdas[-1]
        H[j, i] = lambdas[-1]  # Ensure symmetry

def substitution_rules(I, J):
    """Returns 2 dictionaries containing substitution rules for H"""
    omega_values = np.sort([np.random.uniform(-10, 10) for _ in range(I)])  
    substitution_dict = {omegas[i]: omega_values[i] for i in range(I)}
    substitution_omegas = {omegas[i]: omega_values[i] for i in range(I)}

    lambda_index = 0
    for i in range(I):
        for j in range(i+1, J):  # Iterate only over the upper triangular part
            bound = 0.1**(lambdas[lambda_index].order) * abs(omega_values[i] - omega_values[j])            
            lambda_value = bound             
            substitution_dict[lambdas[lambda_index]] = lambda_value 
            lambda_index += 1

    return substitution_dict
values = substitution_rules(I,J)

Eff_Frame = EffectiveFrame(H, symbol_values = values)

[32mCreating the EffectiveFrame object with matrix form.[0m
[32mThe EffectiveFrame object has been initialized successfully.[0m

# Calculate the effective model using mask up to the second order
Eff_Frame.solve(max_order = 3, method = "ACE", mask=mask)
# Obtaining the result in the matrix form
result = Eff_Frame.get_H(return_form = "matrix")

[32mSubstituting the symbol values in the Hamiltonian and perturbative interactions.[0m
Computing the effective Hamiltonian: 100%|██████████| 3/3 [00:21<00:00,  7.29s/it]
[32mThe Hamiltonian has been solved successfully. Please use the get_H method to get the result in the desired form.[0m
Converting to matrix form: 100%|██████████| 4/4 [00:00<00:00, 83.49it/s]

Plotting

np_H = np.array(H.xreplace(values), dtype = float)

fig, ax = plt.subplots(dpi = 100)
mapp = ax.imshow(np_H, cmap = "inferno")
fig.colorbar(mapp, ax = ax)
#plt.savefig('ACE_non_masked_hamiltonian.pdf', bbox_inches = 'tight')

Total Hamiltonian

np_result = np.array(result).astype(float)
np_result[np.abs(np_result) == 0] = None

fig, ax = plt.subplots(dpi = 100)
mapp = ax.imshow(np_result, cmap = "inferno")
fig.colorbar(mapp)
#plt.savefig('ACE_masked_hamiltonian.pdf',  bbox_inches='tight')

Effective Hamiltonian

License

SymPT is licensed under the MIT License. See the LICENSE file for details.

Citation

If you use SymPT in your research, please cite the following paper:

BibTeX Entry:

@misc{diotallevi2024symptcomprehensivetoolautomating,
      title={SymPT: a comprehensive tool for automating effective Hamiltonian derivations}, 
      author={Giovanni Francesco Diotallevi and Leander Reascos and Mónica Benito},
      year={2024},
      eprint={2412.10240},
      archivePrefix={arXiv},
      primaryClass={quant-ph},
      url={https://arxiv.org/abs/2412.10240}, 
}

APA Citation:

Diotallevi, G. F., Reascos, L., & Benito, M. (2024). SymPT: a comprehensive tool for automating effective Hamiltonian derivations. arXiv preprint arXiv:2412.10240.

IEEE Citation:

G. F. Diotallevi, L. Reascos, and M. Benito, "SymPT: a comprehensive tool for automating effective Hamiltonian derivations," arXiv preprint arXiv:2412.10240, 2024.

References

SymPT: a comprehensive tool for automating effective Hamiltonian derivations