tnc/builders/
peps.rs

1use itertools::iproduct;
2
3use crate::tensornetwork::tensor::Tensor;
4
5/// Generate the initial state for an `length` x `depth` PEPS.
6///
7/// # Arguments
8/// * `length` - length of the PEPS
9/// * `depth` - depth of the PEPS
10/// * `physical_dim` - physical dimension of the PEPs
11/// * `virtual_dim` - virtual bond dimension between lattice sites in the PEPs
12///
13/// # Returns
14/// [`Tensor`] of initial state of PEPS
15fn peps_init(length: usize, depth: usize, physical_dim: u64, virtual_dim: u64) -> Tensor {
16    let mut pep = Tensor::default();
17
18    let physical_up = length * depth;
19    let virtual_vertical = (length - 1) * depth;
20    let virtual_horizontal = (depth - 1) * length;
21    let total_edges = physical_up + virtual_vertical + virtual_horizontal;
22
23    let mut tensors = vec![Tensor::default(); length * depth];
24
25    // Consider the corners
26    tensors[0] = Tensor::new(
27        vec![0, physical_up, physical_up + virtual_vertical],
28        vec![physical_dim, virtual_dim, virtual_dim],
29    );
30    tensors[length - 1] = Tensor::new(
31        vec![
32            length - 1,
33            physical_up + length - 2,
34            physical_up + virtual_vertical + length - 1,
35        ],
36        vec![physical_dim, virtual_dim, virtual_dim],
37    );
38    tensors[length * (depth - 1)] = Tensor::new(
39        vec![
40            length * (depth - 1),
41            physical_up + (length - 1) * (depth - 1),
42            physical_up + virtual_vertical + length * (depth - 2),
43        ],
44        vec![physical_dim, virtual_dim, virtual_dim],
45    );
46    tensors[length * depth - 1] = Tensor::new(
47        vec![
48            length * depth - 1,
49            physical_up + (length - 1) * depth - 1,
50            physical_up + virtual_vertical + length * (depth - 1) - 1,
51        ],
52        vec![physical_dim, virtual_dim, virtual_dim],
53    );
54
55    // Consider the horizontal edges
56    for j in 1..(length - 1) {
57        tensors[j] = Tensor::new(
58            vec![
59                j,
60                physical_up + j - 1,
61                physical_up + j,
62                physical_up + virtual_vertical + j,
63            ],
64            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
65        );
66
67        tensors[physical_up - j - 1] = Tensor::new(
68            vec![
69                physical_up - j - 1,
70                physical_up + virtual_vertical - j - 1,
71                physical_up + virtual_vertical - j,
72                total_edges - j - 1,
73            ],
74            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
75        );
76    }
77
78    // Consider the vertical edges
79    for i in 1..(depth - 1) {
80        tensors[i * length] = Tensor::new(
81            vec![
82                i * length,
83                physical_up + i * (length - 1),
84                physical_up + virtual_vertical + (i - 1) * length,
85                physical_up + virtual_vertical + i * length,
86            ],
87            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
88        );
89
90        tensors[(i + 1) * length - 1] = Tensor::new(
91            vec![
92                (i + 1) * length - 1,
93                physical_up + (i + 1) * (length - 1) - 1,
94                physical_up + virtual_vertical + i * length - 1,
95                physical_up + virtual_vertical + (i + 1) * length - 1,
96            ],
97            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
98        );
99    }
100
101    // Consider the remaining bulk not on the edges
102    for (i, j) in iproduct!(1..(depth - 1), 1..(length - 1)) {
103        let index = i * length + j;
104        tensors[index] = Tensor::new(
105            vec![
106                index,
107                physical_up + i * (length - 1) + j - 1,
108                physical_up + i * (length - 1) + j,
109                physical_up + virtual_vertical + (i - 1) * length + j,
110                physical_up + virtual_vertical + i * length + j,
111            ],
112            vec![
113                physical_dim,
114                virtual_dim,
115                virtual_dim,
116                virtual_dim,
117                virtual_dim,
118            ],
119        );
120    }
121
122    pep.push_tensors(tensors);
123    pep
124}
125
126/// Generate an intermediate PEP-operator an n x n PEPs with shared bond dimension of `dimension`.
127///
128/// # Arguments
129/// * `peps` - mutable [`Tensor`] representing initial PEPS and successive PEPOs
130/// * `length` - length of the PEPS
131/// * `depth` - depth of the PEPS
132/// * `layer` - layer of PEPO to be applied
133/// * `physical_dim` - physical dimension of the PEPs
134/// * `virtual_dim` - virtual bond dimension between lattice sites in the PEPs
135///
136/// # Returns
137/// Updated PEPS [`Tensor`] with `layer`th PEPO applied.
138fn pepo(
139    mut peps: Tensor,
140    length: usize,
141    depth: usize,
142    layer: usize,
143    physical_dim: u64,
144    virtual_dim: u64,
145) -> Tensor {
146    let physical_up = length * depth;
147    let virtual_vertical = (length - 1) * depth;
148    let virtual_horizontal = (depth - 1) * length;
149    let total_edges = physical_up + virtual_vertical + virtual_horizontal;
150    let last = total_edges * layer;
151    let start = total_edges * (layer + 1);
152
153    let mut tensors = vec![Tensor::default(); length * depth];
154
155    // Consider the corners
156    tensors[0] = Tensor::new(
157        vec![
158            last,
159            start,
160            start + physical_up,
161            start + physical_up + virtual_vertical,
162        ],
163        vec![physical_dim, physical_dim, virtual_dim, virtual_dim],
164    );
165    tensors[length - 1] = Tensor::new(
166        vec![
167            last + length - 1,
168            start + length - 1,
169            start + physical_up + length - 2,
170            start + physical_up + virtual_vertical + length - 1,
171        ],
172        vec![physical_dim, physical_dim, virtual_dim, virtual_dim],
173    );
174    tensors[length * (depth - 1)] = Tensor::new(
175        vec![
176            last + length * (depth - 1),
177            start + length * (depth - 1),
178            start + physical_up + (length - 1) * (depth - 1),
179            start + physical_up + virtual_vertical + length * (depth - 2),
180        ],
181        vec![physical_dim, physical_dim, virtual_dim, virtual_dim],
182    );
183    tensors[length * depth - 1] = Tensor::new(
184        vec![
185            last + length * depth - 1,
186            start + length * depth - 1,
187            start + physical_up + (length - 1) * depth - 1,
188            start + physical_up + virtual_vertical + length * (depth - 1) - 1,
189        ],
190        vec![physical_dim, physical_dim, virtual_dim, virtual_dim],
191    );
192
193    // Consider the horizontal edges
194    for j in 1..(length - 1) {
195        tensors[j] = Tensor::new(
196            vec![
197                last + j,
198                start + j,
199                start + physical_up + j - 1,
200                start + physical_up + j,
201                start + physical_up + virtual_vertical + j,
202            ],
203            vec![
204                physical_dim,
205                physical_dim,
206                virtual_dim,
207                virtual_dim,
208                virtual_dim,
209            ],
210        );
211
212        tensors[physical_up - j - 1] = Tensor::new(
213            vec![
214                last + physical_up - j - 1,
215                start + physical_up - j - 1,
216                start + physical_up + virtual_vertical - j - 1,
217                start + physical_up + virtual_vertical - j,
218                start + total_edges - j - 1,
219            ],
220            vec![
221                physical_dim,
222                physical_dim,
223                virtual_dim,
224                virtual_dim,
225                virtual_dim,
226            ],
227        );
228    }
229
230    // Consider the vertical edges
231    for i in 1..(depth - 1) {
232        tensors[i * length] = Tensor::new(
233            vec![
234                last + i * length,
235                start + i * length,
236                start + physical_up + i * (length - 1),
237                start + physical_up + virtual_vertical + (i - 1) * length,
238                start + physical_up + virtual_vertical + i * length,
239            ],
240            vec![
241                physical_dim,
242                physical_dim,
243                virtual_dim,
244                virtual_dim,
245                virtual_dim,
246            ],
247        );
248
249        tensors[(i + 1) * length - 1] = Tensor::new(
250            vec![
251                last + (i + 1) * length - 1,
252                start + (i + 1) * length - 1,
253                start + physical_up + (i + 1) * (length - 1) - 1,
254                start + physical_up + virtual_vertical + i * length - 1,
255                start + physical_up + virtual_vertical + (i + 1) * length - 1,
256            ],
257            vec![
258                physical_dim,
259                physical_dim,
260                virtual_dim,
261                virtual_dim,
262                virtual_dim,
263            ],
264        );
265    }
266
267    // Consider the remaining bulk not on the edges
268    for (i, j) in iproduct!(1..(depth - 1), 1..(length - 1)) {
269        let index = i * length + j;
270        tensors[index] = Tensor::new(
271            vec![
272                last + index,
273                start + index,
274                start + physical_up + i * (length - 1) + j - 1,
275                start + physical_up + i * (length - 1) + j,
276                start + physical_up + virtual_vertical + (i - 1) * length + j,
277                start + physical_up + virtual_vertical + i * length + j,
278            ],
279            vec![
280                physical_dim,
281                physical_dim,
282                virtual_dim,
283                virtual_dim,
284                virtual_dim,
285                virtual_dim,
286            ],
287        );
288    }
289    peps.push_tensors(tensors);
290    peps
291}
292
293/// Applies the final state in a PEPS.
294///
295/// # Arguments
296/// * `peps` - mutable [`Tensor`] representing initial PEPS and successive PEPOs
297/// * `length` - length of the PEPS
298/// * `depth` - depth of the PEPS
299/// * `physical_dim` - physical dimension of the PEPs
300/// * `virtual_dim` - virtual bond dimension between lattice sites in the PEPs
301///
302/// # Returns
303/// Updated PEPS [`Tensor`] with final PEPS applied
304fn peps_final(
305    mut peps: Tensor,
306    length: usize,
307    depth: usize,
308    physical_dim: u64,
309    virtual_dim: u64,
310    layers: usize,
311) -> Tensor {
312    let physical_up = length * depth;
313    let virtual_vertical = (length - 1) * depth;
314    let virtual_horizontal = (depth - 1) * length;
315    let mut total_edges = physical_up + virtual_vertical + virtual_horizontal;
316    let last = total_edges * layers;
317    let start = total_edges * (layers + 1);
318    total_edges -= physical_up;
319
320    let mut tensors = vec![Tensor::default(); length * depth];
321
322    // Consider the corners
323    tensors[0] = Tensor::new(
324        vec![last, start, start + virtual_vertical],
325        vec![physical_dim, virtual_dim, virtual_dim],
326    );
327    tensors[length - 1] = Tensor::new(
328        vec![
329            last + length - 1,
330            start + length - 2,
331            start + virtual_vertical + length - 1,
332        ],
333        vec![physical_dim, virtual_dim, virtual_dim],
334    );
335    tensors[length * (depth - 1)] = Tensor::new(
336        vec![
337            last + length * (depth - 1),
338            start + (length - 1) * (depth - 1),
339            start + virtual_vertical + length * (depth - 2),
340        ],
341        vec![physical_dim, virtual_dim, virtual_dim],
342    );
343    tensors[length * depth - 1] = Tensor::new(
344        vec![
345            last + length * depth - 1,
346            start + (length - 1) * depth - 1,
347            start + virtual_vertical + length * (depth - 1) - 1,
348        ],
349        vec![physical_dim, virtual_dim, virtual_dim],
350    );
351
352    // Consider the horizontal edges
353    for j in 1..(length - 1) {
354        tensors[j] = Tensor::new(
355            vec![
356                last + j,
357                start + j - 1,
358                start + j,
359                start + virtual_vertical + j,
360            ],
361            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
362        );
363
364        tensors[physical_up - j - 1] = Tensor::new(
365            vec![
366                last + physical_up - j - 1,
367                start + virtual_vertical - j - 1,
368                start + virtual_vertical - j,
369                start + total_edges - j - 1,
370            ],
371            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
372        );
373    }
374
375    // Consider the vertical edges
376    for i in 1..(depth - 1) {
377        tensors[i * length] = Tensor::new(
378            vec![
379                last + i * length,
380                start + i * (length - 1),
381                start + virtual_vertical + (i - 1) * length,
382                start + virtual_vertical + i * length,
383            ],
384            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
385        );
386
387        tensors[(i + 1) * length - 1] = Tensor::new(
388            vec![
389                last + (i + 1) * length - 1,
390                start + (i + 1) * (length - 1) - 1,
391                start + virtual_vertical + i * length - 1,
392                start + virtual_vertical + (i + 1) * length - 1,
393            ],
394            vec![physical_dim, virtual_dim, virtual_dim, virtual_dim],
395        );
396    }
397
398    // Consider the remaining bulk not on the edges
399    for (i, j) in iproduct!(1..(depth - 1), 1..(length - 1)) {
400        let index = i * length + j;
401        tensors[index] = Tensor::new(
402            vec![
403                last + index,
404                start + i * (length - 1) + j - 1,
405                start + i * (length - 1) + j,
406                start + virtual_vertical + (i - 1) * length + j,
407                start + virtual_vertical + i * length + j,
408            ],
409            vec![
410                physical_dim,
411                virtual_dim,
412                virtual_dim,
413                virtual_dim,
414                virtual_dim,
415            ],
416        );
417    }
418    peps.push_tensors(tensors);
419    peps
420}
421
422/// Generates the structure for a PEPS with `length` x `depth` dimensions and with `layers` layers.
423///
424/// The `EdgeIndex` in the PEPS `bond_dims` for each layer is ordered as `physical bonds to previous layer,
425/// physical bonds to next layer, vertical virtual bonds, horizontal virtual bonds`.
426///
427/// Each new layer other than the final layer adds:
428/// * p = length * depth physical bonds
429/// * vv = (length - 1) * depth vertical virtual bonds
430/// * vh = (depth - 1) * length horizontal virtual bonds
431///
432/// Bonds for the `k`th layer, where k is not the initial or final layer, then run from:
433/// * physical bonds to previous layer: (k-1) * (p+vv+vh): (kp) + (k-1) * (vv+vh)
434/// * physical bonds to next layer: k * (p+vv+vh)
435///
436/// # Arguments
437/// * `length` - length of the PEPS
438/// * `depth` - depth of the PEPS
439/// * `physical_dim` - physical dimension of the PEPs
440/// * `virtual_dim` - virtual bond dimension between lattice sites in the PEPs
441/// * `layers` - number of operator layers in PEPS, 0 layers returns an inner product of two states.
442///
443/// # Returns
444/// [`Tensor`] representing a PEPS.
445///
446/// # Panics
447/// Panics if `length < 2` or `depth < 2`.
448#[must_use]
449pub fn peps(
450    length: usize,
451    depth: usize,
452    physical_dim: u64,
453    virtual_dim: u64,
454    layers: usize,
455) -> Tensor {
456    assert!(length > 1, "PEPS should have length greater than 1");
457    assert!(depth > 1, "PEPS should have depth greater than 1");
458    let mut new_peps = peps_init(length, depth, physical_dim, virtual_dim);
459    for layer in 0..layers {
460        new_peps = pepo(new_peps, length, depth, layer, physical_dim, virtual_dim);
461    }
462    peps_final(new_peps, length, depth, physical_dim, virtual_dim, layers)
463}
464
465#[cfg(test)]
466mod tests {
467    use super::*;
468
469    use std::iter::zip;
470
471    use rustc_hash::FxHashMap;
472
473    use crate::tensornetwork::tensor::Tensor;
474
475    #[test]
476    fn test_pep_init() {
477        let length = 3;
478        let depth = 3;
479        let physical_dim = 4;
480        let virtual_dim = 10;
481
482        let bond_dims = FxHashMap::from_iter([
483            (0, 4),
484            (1, 4),
485            (2, 4),
486            (3, 4),
487            (4, 4),
488            (5, 4),
489            (6, 4),
490            (7, 4),
491            (8, 4),
492            (9, 10),
493            (10, 10),
494            (11, 10),
495            (12, 10),
496            (13, 10),
497            (14, 10),
498            (15, 10),
499            (16, 10),
500            (17, 10),
501            (18, 10),
502            (19, 10),
503            (20, 10),
504        ]);
505        let tensors = vec![
506            Tensor::new_from_map(vec![0, 9, 15], &bond_dims),
507            Tensor::new_from_map(vec![1, 9, 10, 16], &bond_dims),
508            Tensor::new_from_map(vec![2, 10, 17], &bond_dims),
509            Tensor::new_from_map(vec![3, 11, 15, 18], &bond_dims),
510            Tensor::new_from_map(vec![4, 11, 12, 16, 19], &bond_dims),
511            Tensor::new_from_map(vec![5, 12, 17, 20], &bond_dims),
512            Tensor::new_from_map(vec![6, 13, 18], &bond_dims),
513            Tensor::new_from_map(vec![7, 13, 14, 19], &bond_dims),
514            Tensor::new_from_map(vec![8, 14, 20], &bond_dims),
515        ];
516        let ref_tensor = Tensor::new_composite(tensors);
517
518        let new_peps = peps_init(length, depth, physical_dim, virtual_dim);
519        for (t1, t2) in zip(new_peps.tensors().iter(), ref_tensor.tensors().iter()) {
520            assert_eq!(t1.legs(), t2.legs());
521            assert_eq!(t1.bond_dims(), t2.bond_dims());
522        }
523    }
524
525    #[test]
526    fn test_pepo() {
527        let length = 3;
528        let depth = 3;
529        let physical_dim = 4;
530        let virtual_dim = 10;
531        let layers = 1;
532
533        let bond_dims = FxHashMap::from_iter([
534            (0, 4),
535            (1, 4),
536            (2, 4),
537            (3, 4),
538            (4, 4),
539            (5, 4),
540            (6, 4),
541            (7, 4),
542            (8, 4),
543            (9, 10),
544            (10, 10),
545            (11, 10),
546            (12, 10),
547            (13, 10),
548            (14, 10),
549            (15, 10),
550            (16, 10),
551            (17, 10),
552            (18, 10),
553            (19, 10),
554            (20, 10),
555            (21, 4),
556            (22, 4),
557            (23, 4),
558            (24, 4),
559            (25, 4),
560            (26, 4),
561            (27, 4),
562            (28, 4),
563            (29, 4),
564            (30, 10),
565            (31, 10),
566            (32, 10),
567            (33, 10),
568            (34, 10),
569            (35, 10),
570            (36, 10),
571            (37, 10),
572            (38, 10),
573            (39, 10),
574            (40, 10),
575            (41, 10),
576        ]);
577        let tensors = vec![
578            Tensor::new_from_map(vec![0, 9, 15], &bond_dims),
579            Tensor::new_from_map(vec![1, 9, 10, 16], &bond_dims),
580            Tensor::new_from_map(vec![2, 10, 17], &bond_dims),
581            Tensor::new_from_map(vec![3, 11, 15, 18], &bond_dims),
582            Tensor::new_from_map(vec![4, 11, 12, 16, 19], &bond_dims),
583            Tensor::new_from_map(vec![5, 12, 17, 20], &bond_dims),
584            Tensor::new_from_map(vec![6, 13, 18], &bond_dims),
585            Tensor::new_from_map(vec![7, 13, 14, 19], &bond_dims),
586            Tensor::new_from_map(vec![8, 14, 20], &bond_dims),
587            Tensor::new_from_map(vec![0, 21, 30, 36], &bond_dims),
588            Tensor::new_from_map(vec![1, 22, 30, 31, 37], &bond_dims),
589            Tensor::new_from_map(vec![2, 23, 31, 38], &bond_dims),
590            Tensor::new_from_map(vec![3, 24, 32, 36, 39], &bond_dims),
591            Tensor::new_from_map(vec![4, 25, 32, 33, 37, 40], &bond_dims),
592            Tensor::new_from_map(vec![5, 26, 33, 38, 41], &bond_dims),
593            Tensor::new_from_map(vec![6, 27, 34, 39], &bond_dims),
594            Tensor::new_from_map(vec![7, 28, 34, 35, 40], &bond_dims),
595            Tensor::new_from_map(vec![8, 29, 35, 41], &bond_dims),
596        ];
597        let ref_tensor = Tensor::new_composite(tensors);
598
599        let mut new_peps = peps_init(length, depth, physical_dim, virtual_dim);
600        for layer in 0..layers {
601            new_peps = pepo(new_peps, length, depth, layer, physical_dim, virtual_dim);
602        }
603
604        for (t1, t2) in zip(new_peps.tensors().iter(), ref_tensor.tensors().iter()) {
605            assert_eq!(t1.legs(), t2.legs());
606            assert_eq!(t1.bond_dims(), t2.bond_dims());
607        }
608    }
609
610    #[test]
611    fn test_peps_final() {
612        let length = 3;
613        let depth = 3;
614        let physical_dim = 4;
615        let virtual_dim = 10;
616        let layers = 1;
617
618        let bond_dims = FxHashMap::from_iter([
619            (0, 4),
620            (1, 4),
621            (2, 4),
622            (3, 4),
623            (4, 4),
624            (5, 4),
625            (6, 4),
626            (7, 4),
627            (8, 4),
628            (9, 10),
629            (10, 10),
630            (11, 10),
631            (12, 10),
632            (13, 10),
633            (14, 10),
634            (15, 10),
635            (16, 10),
636            (17, 10),
637            (18, 10),
638            (19, 10),
639            (20, 10),
640            (21, 4),
641            (22, 4),
642            (23, 4),
643            (24, 4),
644            (25, 4),
645            (26, 4),
646            (27, 4),
647            (28, 4),
648            (29, 4),
649            (30, 10),
650            (31, 10),
651            (32, 10),
652            (33, 10),
653            (34, 10),
654            (35, 10),
655            (36, 10),
656            (37, 10),
657            (38, 10),
658            (39, 10),
659            (40, 10),
660            (41, 10),
661            (42, 10),
662            (43, 10),
663            (44, 10),
664            (45, 10),
665            (46, 10),
666            (47, 10),
667            (48, 10),
668            (49, 10),
669            (50, 10),
670            (51, 10),
671            (52, 10),
672            (53, 10),
673        ]);
674        let tensors = vec![
675            Tensor::new_from_map(vec![0, 9, 15], &bond_dims),
676            Tensor::new_from_map(vec![1, 9, 10, 16], &bond_dims),
677            Tensor::new_from_map(vec![2, 10, 17], &bond_dims),
678            Tensor::new_from_map(vec![3, 11, 15, 18], &bond_dims),
679            Tensor::new_from_map(vec![4, 11, 12, 16, 19], &bond_dims),
680            Tensor::new_from_map(vec![5, 12, 17, 20], &bond_dims),
681            Tensor::new_from_map(vec![6, 13, 18], &bond_dims),
682            Tensor::new_from_map(vec![7, 13, 14, 19], &bond_dims),
683            Tensor::new_from_map(vec![8, 14, 20], &bond_dims),
684            Tensor::new_from_map(vec![0, 21, 30, 36], &bond_dims),
685            Tensor::new_from_map(vec![1, 22, 30, 31, 37], &bond_dims),
686            Tensor::new_from_map(vec![2, 23, 31, 38], &bond_dims),
687            Tensor::new_from_map(vec![3, 24, 32, 36, 39], &bond_dims),
688            Tensor::new_from_map(vec![4, 25, 32, 33, 37, 40], &bond_dims),
689            Tensor::new_from_map(vec![5, 26, 33, 38, 41], &bond_dims),
690            Tensor::new_from_map(vec![6, 27, 34, 39], &bond_dims),
691            Tensor::new_from_map(vec![7, 28, 34, 35, 40], &bond_dims),
692            Tensor::new_from_map(vec![8, 29, 35, 41], &bond_dims),
693            Tensor::new_from_map(vec![21, 42, 48], &bond_dims),
694            Tensor::new_from_map(vec![22, 42, 43, 49], &bond_dims),
695            Tensor::new_from_map(vec![23, 43, 50], &bond_dims),
696            Tensor::new_from_map(vec![24, 44, 48, 51], &bond_dims),
697            Tensor::new_from_map(vec![25, 44, 45, 49, 52], &bond_dims),
698            Tensor::new_from_map(vec![26, 45, 50, 53], &bond_dims),
699            Tensor::new_from_map(vec![27, 46, 51], &bond_dims),
700            Tensor::new_from_map(vec![28, 46, 47, 52], &bond_dims),
701            Tensor::new_from_map(vec![29, 47, 53], &bond_dims),
702        ];
703        let ref_tensor = Tensor::new_composite(tensors);
704
705        let mut new_peps = peps_init(length, depth, physical_dim, virtual_dim);
706        for layer in 0..layers {
707            new_peps = pepo(new_peps, length, depth, layer, physical_dim, virtual_dim);
708        }
709        let new_peps = peps_final(new_peps, length, depth, physical_dim, virtual_dim, layers);
710        for (t1, t2) in zip(new_peps.tensors().iter(), ref_tensor.tensors().iter()) {
711            assert_eq!(t1.legs(), t2.legs());
712            assert_eq!(t1.bond_dims(), t2.bond_dims());
713        }
714    }
715
716    #[test]
717    fn test_peps() {
718        let length = 3;
719        let depth = 3;
720        let physical_dim = 4;
721        let virtual_dim = 10;
722        let layers = 1;
723
724        let bond_dims = FxHashMap::from_iter([
725            (0, 4),
726            (1, 4),
727            (2, 4),
728            (3, 4),
729            (4, 4),
730            (5, 4),
731            (6, 4),
732            (7, 4),
733            (8, 4),
734            (9, 10),
735            (10, 10),
736            (11, 10),
737            (12, 10),
738            (13, 10),
739            (14, 10),
740            (15, 10),
741            (16, 10),
742            (17, 10),
743            (18, 10),
744            (19, 10),
745            (20, 10),
746            (21, 4),
747            (22, 4),
748            (23, 4),
749            (24, 4),
750            (25, 4),
751            (26, 4),
752            (27, 4),
753            (28, 4),
754            (29, 4),
755            (30, 10),
756            (31, 10),
757            (32, 10),
758            (33, 10),
759            (34, 10),
760            (35, 10),
761            (36, 10),
762            (37, 10),
763            (38, 10),
764            (39, 10),
765            (40, 10),
766            (41, 10),
767            (42, 10),
768            (43, 10),
769            (44, 10),
770            (45, 10),
771            (46, 10),
772            (47, 10),
773            (48, 10),
774            (49, 10),
775            (50, 10),
776            (51, 10),
777            (52, 10),
778            (53, 10),
779        ]);
780        let tensors = vec![
781            Tensor::new_from_map(vec![0, 9, 15], &bond_dims),
782            Tensor::new_from_map(vec![1, 9, 10, 16], &bond_dims),
783            Tensor::new_from_map(vec![2, 10, 17], &bond_dims),
784            Tensor::new_from_map(vec![3, 11, 15, 18], &bond_dims),
785            Tensor::new_from_map(vec![4, 11, 12, 16, 19], &bond_dims),
786            Tensor::new_from_map(vec![5, 12, 17, 20], &bond_dims),
787            Tensor::new_from_map(vec![6, 13, 18], &bond_dims),
788            Tensor::new_from_map(vec![7, 13, 14, 19], &bond_dims),
789            Tensor::new_from_map(vec![8, 14, 20], &bond_dims),
790            Tensor::new_from_map(vec![0, 21, 30, 36], &bond_dims),
791            Tensor::new_from_map(vec![1, 22, 30, 31, 37], &bond_dims),
792            Tensor::new_from_map(vec![2, 23, 31, 38], &bond_dims),
793            Tensor::new_from_map(vec![3, 24, 32, 36, 39], &bond_dims),
794            Tensor::new_from_map(vec![4, 25, 32, 33, 37, 40], &bond_dims),
795            Tensor::new_from_map(vec![5, 26, 33, 38, 41], &bond_dims),
796            Tensor::new_from_map(vec![6, 27, 34, 39], &bond_dims),
797            Tensor::new_from_map(vec![7, 28, 34, 35, 40], &bond_dims),
798            Tensor::new_from_map(vec![8, 29, 35, 41], &bond_dims),
799            Tensor::new_from_map(vec![21, 42, 48], &bond_dims),
800            Tensor::new_from_map(vec![22, 42, 43, 49], &bond_dims),
801            Tensor::new_from_map(vec![23, 43, 50], &bond_dims),
802            Tensor::new_from_map(vec![24, 44, 48, 51], &bond_dims),
803            Tensor::new_from_map(vec![25, 44, 45, 49, 52], &bond_dims),
804            Tensor::new_from_map(vec![26, 45, 50, 53], &bond_dims),
805            Tensor::new_from_map(vec![27, 46, 51], &bond_dims),
806            Tensor::new_from_map(vec![28, 46, 47, 52], &bond_dims),
807            Tensor::new_from_map(vec![29, 47, 53], &bond_dims),
808        ];
809        let ref_tensor = Tensor::new_composite(tensors);
810
811        let new_peps = peps(length, depth, physical_dim, virtual_dim, layers);
812        for (t1, t2) in zip(new_peps.tensors().iter(), ref_tensor.tensors().iter()) {
813            assert_eq!(t1.legs(), t2.legs());
814            assert_eq!(t1.bond_dims(), t2.bond_dims());
815        }
816    }
817
818    #[test]
819    fn test_inner_product() {
820        let length = 2;
821        let depth = 2;
822        let physical_dim = 4;
823        let virtual_dim = 10;
824        let layers = 0;
825
826        let bond_dims = FxHashMap::from_iter([
827            (0, 4),
828            (1, 4),
829            (2, 4),
830            (3, 4),
831            (4, 10),
832            (5, 10),
833            (6, 10),
834            (7, 10),
835            (8, 10),
836            (9, 10),
837            (10, 10),
838            (11, 10),
839        ]);
840        let tensors = vec![
841            Tensor::new_from_map(vec![0, 4, 6], &bond_dims),
842            Tensor::new_from_map(vec![1, 4, 7], &bond_dims),
843            Tensor::new_from_map(vec![2, 5, 6], &bond_dims),
844            Tensor::new_from_map(vec![3, 5, 7], &bond_dims),
845            Tensor::new_from_map(vec![0, 8, 10], &bond_dims),
846            Tensor::new_from_map(vec![1, 8, 11], &bond_dims),
847            Tensor::new_from_map(vec![2, 9, 10], &bond_dims),
848            Tensor::new_from_map(vec![3, 9, 11], &bond_dims),
849        ];
850        let ref_tensor = Tensor::new_composite(tensors);
851
852        let mut new_peps = peps_init(length, depth, physical_dim, virtual_dim);
853        for layer in 0..layers {
854            new_peps = pepo(new_peps, length, depth, layer, physical_dim, virtual_dim);
855        }
856        let new_peps = peps_final(new_peps, length, depth, physical_dim, virtual_dim, layers);
857        for (t1, t2) in zip(new_peps.tensors().iter(), ref_tensor.tensors().iter()) {
858            assert_eq!(t1.legs(), t2.legs());
859            assert_eq!(t1.bond_dims(), t2.bond_dims());
860        }
861    }
862
863    #[test]
864    #[should_panic(expected = "PEPS should have length greater than 1")]
865    fn test_mps() {
866        let length = 1;
867        let depth = 2;
868        let physical_dim = 4;
869        let virtual_dim = 10;
870        let layers = 0;
871
872        let bond_dims = FxHashMap::from_iter([(0, 4), (1, 4), (2, 10)]);
873        let tensors = vec![
874            Tensor::new_from_map(vec![0, 2], &bond_dims),
875            Tensor::new_from_map(vec![1, 2], &bond_dims),
876        ];
877        let ref_tensor = Tensor::new_composite(tensors);
878
879        let new_peps = peps(length, depth, physical_dim, virtual_dim, layers);
880        for (t1, t2) in zip(new_peps.tensors().iter(), ref_tensor.tensors().iter()) {
881            assert_eq!(t1.legs(), t2.legs());
882            assert_eq!(t1.bond_dims(), t2.bond_dims());
883        }
884    }
885}