The usage of contractor#

Overview#

In this tutorial, we will demonstrate how to utilize different types of TensorNetwork contractors for the circuit simulation to achieve a better space-time tradeoff. The customization of the contractor is the main highlight for the TensorCircuit package since a better contractor can make better use of the power of the TensorNetwok simulation engine. [WIP]

Setup#

[1]:
import tensorcircuit as tc
import numpy as np
import cotengra as ctg
import opt_einsum as oem

K = tc.set_backend("tensorflow")

Testbed System#

We provide tensor networks for two circuits, and test the contraction efficiency for these two systems, the first system is small while the second one is large.

[2]:
# get state for small system
def small_tn():
    n = 10
    d = 4
    param = K.ones([2 * d, n])
    c = tc.Circuit(n)
    c = tc.templates.blocks.example_block(c, param, nlayers=d)
    return c.state()
[3]:
# get expectation for extra large system
def large_tn():
    n = 60
    d = 8
    param = K.ones([2 * d, n])
    c = tc.Circuit(n)
    c = tc.templates.blocks.example_block(c, param, nlayers=d, is_split=True)
    # the two qubit gate is split and truncated with SVD decomposition
    return c.expectation([tc.gates.z(), [n // 2]], reuse=False)

Opt-einsum#

There are several contractor optimizers provided by opt-einsum and shipped with the TensorNetwork package. Since TensorCircuit is built on top of TensorNetwork, we can use these simple contractor optimizers. Though for any moderate system, only greedy optimizer works, other optimizers come with exponential scaling and fail in circuit simulation scenarios.

We always set contraction_info=True (default False) for the contractor system, which will print contraction information summary including contraction size, flops, and writes. For the definition of these metrics, also refer to cotengra docs.

[4]:
# if we set nothing, the default optimizer is greedy, i.e.
tc.set_contractor("greedy", debug_level=2, contraction_info=True)
# We set debug_level=2 to not really run the contraction computation
# i.e. by set debug_level>0, only contraction information and the return shape is correct, the result is wrong
[4]:
functools.partial(<function custom at 0x7fd5a0a3d430>, optimizer=<function contraction_info_decorator.<locals>.new_algorithm at 0x7fd588e281f0>, memory_limit=None, debug_level=2)
[5]:
small_tn()
------ contraction cost summary ------
log10[FLOPs]:  5.132  log2[SIZE]:  11  log2[WRITE]:  13.083
[5]:
<tf.Tensor: shape=(1024,), dtype=complex64, numpy=
array([0.+0.j, 0.+0.j, 0.+0.j, ..., 0.+0.j, 0.+0.j, 0.+0.j],
      dtype=complex64)>
[6]:
large_tn()
------ contraction cost summary ------
log10[FLOPs]:  17.766  log2[SIZE]:  44  log2[WRITE]:  49.636
[6]:
<tf.Tensor: shape=(), dtype=complex64, numpy=0j>
[7]:
# we can use more fancy contractor in opt-einsum but not in tensornetwork
# custom_stateful is used for contraction path finder which has a life cycle of one-time path solver
tc.set_contractor(
    "custom_stateful",
    optimizer=oem.RandomGreedy,
    max_time=60,
    max_repeats=128,
    minimize="size",
    debug_level=2,
    contraction_info=True,
)
[7]:
functools.partial(<function custom_stateful at 0x7fd5a0a3d4c0>, optimizer=<class 'opt_einsum.path_random.RandomGreedy'>, opt_conf=None, contraction_info=True, debug_level=2, max_time=60, max_repeats=128, minimize='size')
[8]:
small_tn()
------ contraction cost summary ------
log10[FLOPs]:  4.925  log2[SIZE]:  10  log2[WRITE]:  12.531
[8]:
<tf.Tensor: shape=(1024,), dtype=complex64, numpy=
array([0.+0.j, 0.+0.j, 0.+0.j, ..., 0.+0.j, 0.+0.j, 0.+0.j],
      dtype=complex64)>
[9]:
large_tn()
------ contraction cost summary ------
log10[FLOPs]:  11.199  log2[SIZE]:  26  log2[WRITE]:  28.183
[9]:
<tf.Tensor: shape=(), dtype=complex64, numpy=0j>

Cotengra#

for more advanced contractors, we ask help for the sota contractor optimizer on the market: cotengra

[10]:
opt = ctg.ReusableHyperOptimizer(
    methods=["greedy", "kahypar"],
    parallel=True,
    minimize="write",
    max_time=120,
    max_repeats=1024,
    progbar=True,
)
tc.set_contractor(
    "custom", optimizer=opt, preprocessing=True, contraction_info=True, debug_level=2
)
## for more setup on cotengra optimizers, see the reference in
## https://cotengra.readthedocs.io/en/latest/advanced.html
## preprocessing=True will merge all single-qubit gates into neighboring two-qubit gates
[10]:
functools.partial(<function custom at 0x7fd5a0a3d430>, optimizer=<function contraction_info_decorator.<locals>.new_algorithm at 0x7fd588e28ee0>, memory_limit=None, debug_level=2, preprocessing=True)
[11]:
small_tn()
log2[SIZE]: 10.00 log10[FLOPs]: 4.90: 100%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆ| 1024/1024 [00:28<00:00, 35.45it/s]
------ contraction cost summary ------
log10[FLOPs]:  4.900  log2[SIZE]:  10  log2[WRITE]:  12.255

[11]:
<tf.Tensor: shape=(1024,), dtype=complex64, numpy=
array([0.+0.j, 0.+0.j, 0.+0.j, ..., 0.+0.j, 0.+0.j, 0.+0.j],
      dtype=complex64)>
[12]:
large_tn()
log2[SIZE]: 20.00 log10[FLOPs]: 9.50:   4%|β–ˆβ–Š                                         | 43/1024 [02:09<49:22,  3.02s/it]
------ contraction cost summary ------
log10[FLOPs]:  9.501  log2[SIZE]:  20  log2[WRITE]:  24.090
[12]:
<tf.Tensor: shape=(), dtype=complex64, numpy=0j>

We can also apply subtree reconf after cotengra find the path, which in general further (and usually greatly) improve flops and write for the contraction. Indeed, the subtree reconf postprocessing is in general more important than increasing the search time or repeats for the optimizer.

[13]:
opt = ctg.ReusableHyperOptimizer(
    minimize="combo",
    max_repeats=1024,
    max_time=240,
    progbar=True,
)


def opt_reconf(inputs, output, size, **kws):
    tree = opt.search(inputs, output, size)
    tree_r = tree.subtree_reconfigure_forest(
        progbar=True, num_trees=10, num_restarts=20, subtree_weight_what=("size",)
    )
    return tree_r.get_path()


tc.set_contractor(
    "custom",
    optimizer=opt_reconf,
    contraction_info=True,
    preprocessing=True,
    debug_level=2,
)
[13]:
functools.partial(<function custom at 0x7fd5a0a3d430>, optimizer=<function contraction_info_decorator.<locals>.new_algorithm at 0x7fd58c832700>, memory_limit=None, debug_level=2, preprocessing=True)
[14]:
small_tn()
log2[SIZE]: 10.00 log10[FLOPs]: 4.87: 100%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆ| 1024/1024 [02:29<00:00,  6.86it/s]
log2[SIZE]: 10.00 log10[FLOPs]: 4.86: 100%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆ| 20/20 [00:31<00:00,  1.57s/it]
------ contraction cost summary ------
log10[FLOPs]:  4.859  log2[SIZE]:  10  log2[WRITE]:  12.583

[14]:
<tf.Tensor: shape=(1024,), dtype=complex64, numpy=
array([0.+0.j, 0.+0.j, 0.+0.j, ..., 0.+0.j, 0.+0.j, 0.+0.j],
      dtype=complex64)>
[15]:
large_tn()
log2[SIZE]: 21.00 log10[FLOPs]: 9.62:   9%|β–ˆβ–ˆβ–ˆβ–‰                                       | 93/1024 [04:04<40:49,  2.63s/it]
log2[SIZE]: 17.00 log10[FLOPs]: 8.66: 100%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆ| 20/20 [03:08<00:00,  9.44s/it]
------ contraction cost summary ------
log10[FLOPs]:  8.657  log2[SIZE]:  17  log2[WRITE]:  25.035
[15]:
<tf.Tensor: shape=(), dtype=complex64, numpy=0j>

We can also extract the tensornetwork directly for the circuit or observable calculations and we can do the contraction by ourselves using any method we like. Moreover, all these contractors or our customized external contractions can still be compatible with jit, automatic differentiation. Specifically, the contraction path solver, though taking some time overhead, is only evaluated once due to the jit infrastructure (note for demonstration usage, we don’t decorate our contraction function with K.jit here).