cuMPS/python/main.py at main · LeonardEyer/cuMPS · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
import os
import sys

sys.path.append("/home/robinpr")
sys.path.append("/home/robinpr/tenpy-main")
from itertools import cycle
import numpy as np
import random as rng
from tenpy.networks.site import SpinHalfSite
from tenpy.networks.mps import MPS
from tenpy.networks.site import kron
from own_TeNPy import update_bond
import logging

rng = np.random.default_rng(seed=42)
# Set the logging level for TenPy to WARNING to suppress info messages
logging.getLogger('tenpy_git').setLevel(logging.WARNING)


class InitialState:
    def down_state(L):
        # print("Create state |00..1..0>")
        state = []
        for i in range(L):
            if i == L // 2:
                state.append('down');
            else:
                state.append('down');
        return state


def execute(psi_t_initial, t_final, parameter, Ntraj=1):
    L = parameter['L']
    p = parameter['p']
    theta = parameter['theta']
    phi = parameter['phi']
    cos_theta_2 = np.cos(theta / 2)
    sin_theta_2 = np.sin(theta / 2)
    # print(psi_t_initial)
    cos_phi_2 = np.cos(phi / 2)
    sin_phi_2 = np.sin(phi / 2)

    site = SpinHalfSite(conserve=None, sort_charge=False)
    op_o = np.array([[1., 0.], [0., 0.]])
    op_i = np.array([[0., 0.], [0., 1.]])
    op_Rx = np.array([[cos_theta_2, -1j * sin_theta_2], [-1j * sin_theta_2, cos_theta_2]])

    op_K1 = np.array([[0., -1j / 2 * sin_phi_2], [0., 0.]])

    op_K2 = np.array([[1., 0.], [0., cos_phi_2]])

    site.add_op('Rx', op_Rx)
    site.add_op('ii', op_i)
    site.add_op('oo', op_o)
    site.add_op('K1', op_K1)
    site.add_op('K2', op_K2)

    opR = kron(site.get_op('oo'), site.get_op('Id'), group=False) + kron(site.get_op('ii'), site.get_op('Rx'),
                                                                         group=False)
    opL = kron(site.get_op('Id'), site.get_op('oo'), group=False) + kron(site.get_op('Rx'), site.get_op('ii'),
                                                                         group=False)

    #with open(output_file, 'w') as f:
    for _ in range(Ntraj):
        psi_t = MPS.from_product_state([site] * L, psi_t_initial, 'finite')
        state_list = [np.sum(psi_t.expectation_value('ii')) / L]

        #pbar = tqdm.tqdm(range(t_final - 1))
        for t in range(t_final - 1):
            # Maximal optimierte Variante
            # prob = psi_t.expectation_value('K2')
            # decisions = rng.random(size=L) > prob
            # print(f"{np.mean(prob) = }", f"{np.sum(decisions) = }")
            #
            # # Batch-Weise Anwendung (keine Parallelität, aber vektorisiert)
            # for op, indices in [('K1', np.where(decisions)[0]),
            #                     ('K2', np.where(~decisions)[0])]:
            #     for i in indices:
            #         psi_t.apply_local_op(i, op, unitary=False)


            # 1. Definition der Update-Sequenz als Zyklus
            update_sequence = [
                (range(0, L - 1, 2), opR),  # Phase 1: Gerade Bonds mit opR
                (range(1, L - 1, 2), opL),  # Phase 2: Ungerade Bonds mit opL
                (range(1, L - 1, 2), opR),  # Phase 3: Ungerade Bonds mit opR
                (range(0, L - 1, 2), opL)  # Phase 4: Gerade Bonds mit opL
            ]

            # 2. Komprimierte Schleife mit paralleler Indexgenerierung
            for indices, op in update_sequence:
                for i in indices:
                    update_bond.SVD_based(psi_t, i, op, parameter['tebd_params']['trunc_params'])
            state_list.append(np.sum(psi_t.expectation_value('ii')) / L)
            #print(f"{t = }, max(chi) = {max(psi_t.chi)}, mean(chi) = {sum(psi_t.chi) / L:.2f}, mean(ii) = {state_list[-1]:.4f}")
            print(f"{t = }, max(chi) = {max(psi_t.chi)}, mean(chi) = {sum(psi_t.chi) / L:.2f}, mean(P1) = {state_list[-1]:.4f}")
        #print(np.array(state_list))

        #f.write(str(state_list) + '\n')
        print(np.array(state_list))
    return


if __name__ == '__main__':
    np.set_printoptions(precision=3)
    parameter = {
        'state_info': 'down_state',
        'theta': 0.14,
        'phi': 0.1222,
        'p': 1,
        'L': 30,
        'model_params': {
            'bc_MPS': 'finite',
            'omega': 6,
            'gamma': 1,
        },
        'tebd_params': {
            'N_steps': 2,
            'dt': 0.05,
            'order': 1,
            'trunc_params': {'chi_max': 1224, 'svd_min': 1.e-12}
        }
    }
    # param_theta = float(parameter['theta'])
    # param_p = str(parameter['p'])
    # base_output_file = f"long_output_theta_{param_theta:.3f}_phi0.185.txt"

    # Check if the file already exists and append a counter if necessary
    # output_file = base_output_file
    # counter = 1
    # while os.path.exists(output_file):
    #    output_file = f"long_output_theta_{param_theta:.3f}_{counter}_phi0.185.txt"
    #    counter += 1

    import time

    start = time.time()
    execute(InitialState.down_state(parameter['L']), t_final=31, parameter=parameter, Ntraj=1)
    duration = time.time() - start
    print(f"Duration: {duration:.6f} seconds")