Notebook_Header_01.png

実際のバックエンドでVQCをテストする(2024年 Lab4)#

この Lab では、Qiskit Patterns のワークフローを使用して、理想的なバックエンド上でシンプルな変分量子分類器 (Variational Quantum Classifier, VQC) を構築し、学習させる方法を学びます。学習された VQC の性能が、ノイズに対してどのように影響を受けるかを体験します。最後に、VQC の深さを浅くする方法を学び、得られた回路を量子ハードウェア上で実行し、エラー抑制やエラー緩和の技術が結果にどのような影響を与えるかを確認します。

まず、この Lab で必要となるすべてのライブラリーをインポートすることから始めましょう:

Setup#

# ### Install Qiskit and relevant packages, if needed
# ### IMPORTANT: Make sure you are on 3.10 > python < 3.12
# %pip install qiskit[visualization]==1.0.2
# %pip install qiskit-ibm-runtime
# %pip install qiskit-aer
# %pip install graphviz
# %pip install qiskit-serverless -U
# %pip install qiskit-transpiler-service -U
# Import all in one cell

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import minimize

from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit.circuit.library import RealAmplitudes
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit.transpiler import InstructionProperties
from qiskit.visualization import plot_distribution
from qiskit.providers.fake_provider import GenericBackendV2
from qiskit.primitives import StatevectorEstimator

from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import (
    QiskitRuntimeService, 
    EstimatorV2 as Estimator,
    SamplerV2 as Sampler,
    EstimatorOptions
)
# Load your API token in .env
import os
from dotenv import load_dotenv

os.environ.pop('QXInstance', None)
os.environ.pop('QXToken', None)

load_dotenv()

instance = os.environ['QXInstance']
token = os.environ['QXToken']

Part 1: 理想的なバックエンドで VQC を設計する#

変分量子分類器 (VQC) は変分量子アルゴリズムの一種で、ある種の分類問題を解くことができるます。特に、VQC はハイブリッドアーキテクチャーを特徴としています。このアーキテクチャーは、パラメーター化された量子回路を初期状態に適用した形の ansatz に基づいています。出力はコスト関数の形で測定され、このコスト関数は回路のパラメーターに対して古典的に最適化され、最小値に収束します。

この Lab では、10羽の鳥のリストが与えられます:

Falcon, Hummingbird, Eagle, Osprey, Heron, Peacock, Parrot, Swan, Toucan, Cardinal.

(ハヤブサ、ハチドリ、ワシ、ミサゴ、サギ、クジャク、オウム、白鳥、オオハシ、ショウジョウコウカンチョウ)

リストの最初の5羽は IBM Quantum の鳥です。プロセッサー・タイプの名前がこれらの鳥にちなんで付けられているためです - こちら をご確認ください。残りの鳥は IBM Quantum 以外の鳥です。

データセットは bird_dataset.csv という辞書型で与えられます。それぞれの鳥は各行に対応し、状態ベクトル形式で与えられた5量子ビットの量子状態と関連付けられています。したがって、各鳥に対して \(2^5 = 32\) 個の係数が見つかります。

この問題では、IBM Quantum の鳥はエンタングルされた量子状態に関連付けられ、IBM Quantum ではない鳥はエンタングルされていない(積の)量子状態に関連付けられます。エンタングルされた量子状態とそうでない量子状態を区別することは、一般に簡単な問題ではないことがわかります。この Lab の目標は、この2つの鳥を区別できる VQC を作ることです。

でも心配しないでください! Qiskit Patterns** に基づいたワークフローを使って、アルゴリズムの全ステップをガイドします。

** Qiskit Patternsは、変分量子アルゴリズムだけでなく、どのようなタイプの問題にも適用できる一般的なワークフローであることに注意してください!

Step 1: 古典的な入力を量子問題にマップする#

Qiskit Patterns のこのステップでは、量子ハードウェア上で実行可能な量子回路と測定セットに問題を変換します。これは一般的に難しい問いであることに注意してください。

そのため、考慮すべき項目リストを以下に示します。

  • 古典的に難しい問題を選びます。エンタングル状態と非エンタングル状態の識別は、学習目的のために小規模な問題を選んでいるとはいえ、古典的に難しい問題です。

  • 結果を期待値やコスト関数に変換する必要があります。コスト関数は、問題の結果を最小化または最大化すべきものとして定義する問題固有の関数であり、化学、機械学習、金融、最適化など様々な用途の変分量子アルゴリズムに一般的に使用されます。

  • 可能な限り、使用するハードウェアの制限と利点を念頭に置くべきです。特に、バックエンドの接続性は重要です。多くのSWAPゲートを実装する必要がないように、IBM Quantum プロセッサーの heavy-hex トポロジーにできるだけ近い回路に問題をマッピングします。

とりあえず、まず古典的なデータを取得しましょう。この場合、それぞれの鳥の状態ベクトル係数です。

注意

このノートブックには同じ階層のdataフォルダ配下に追加のスクリプトが必要です。
birds_dataset_2024.csv, params_0_list_2024.npy, opt_params_shallow_VQC_2024.npy

# Define num_qubits, the number of qubits, for the rest of the Lab
num_qubits = 5

# Load the dictionary 
birds_dataset = pd.read_csv('./data/birds_dataset_2024.csv') 

# Check if the dataset is loaded correctly - coefficients should be complex numbers
for i in range(2**num_qubits):
    key = 'c%.0f' %i
    birds_dataset[key] = birds_dataset[key].astype(np.complex128)

# Print the dataset
birds_dataset
names c0 c1 c2 c3 c4 c5 c6 c7 c8 ... c22 c23 c24 c25 c26 c27 c28 c29 c30 c31
0 Falcon 0.707107+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.707107+0.000000j
1 Hummingbird 0.000000+0.000000j 0.707107+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.707107+0.000000j 0.000000+0.000000j
2 Eagle 0.000000+0.000000j 0.000000+0.000000j 0.707107+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.707107+0.000000j 0.000000+0.000000j 0.000000+0.000000j
3 Osprey 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.707107+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.707107+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j
4 Heron 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.707107+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.707107+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j
5 Peacock 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 1.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j
6 Parrot 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 1.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j
7 Swan 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 1.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j
8 Toucan 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 1.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j
9 Cardinal 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j ... 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.0+0.0j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j 0.000000+0.000000j

10 rows × 33 columns

Exercise 1: データセットの探索#

Exercise

各量子状態の係数を取得し、それらを複素数を含むリストのリスト[[], [], [], ...] に保存してください。

IBM Quantum の鳥の場合は1のラベルを、そうでない鳥の場合は0のラベルを付与したリスト[int, int, ...]を作成してください。

Hide code cell content
list_coefficients = birds_dataset.filter(like='c', axis=1).values.tolist() # Add your code here
list_labels = [1, 1, 1, 1, 1, 0, 0, 0, 0, 0] # Add your code here

データの準備ができたので、それを量子回路にマッピングします。この場合、Qiskit では initialization によって簡単に実装できる 振幅埋め込み(amplitude embedding) 技術を使います。\(32\) 個の係数のそれぞれが、5量子ビットのヒルベルト空間の基底状態の 1 つの振幅に対応するため、振幅埋め込みはこの種のデータセットに対して簡単なマッピングです。後述するように、 qc.initialize はあらゆる量子状態をマッピングできる一般的な関数です。しかし、この関数は非常に深い回路につながる可能性があります。まだデータセットに関する完全な情報がないため、今のところはこのマッピング手法を使うことにします。

より一般的な分類問題に対しては、他のマッピングを考えることができます。マッピングの選択は、与えられたデータセットの構造に依存します。一般的には、データはパラメータ化された量子ゲートの集合のパラメーターにマッピングされます。典型的な例は、ZZ Feature Map です。

振幅埋め込み回路を構築して、確認してみましょう:

index_bird = 2 # You can check different birds by changing the index
amplitudes = list_coefficients[index_bird]

# Build the amplitude embedding 
qc = QuantumCircuit(5)
qc.initialize(amplitudes, range(num_qubits))
qc.measure_all()

# Draw the amplitude embedding circuit
qc.draw(output="mpl")
_images/62ab9395486d8b61b2b55783a14f4925ce32cb7f00506de12111feef36b4c707.png

回路を描画する前に回路を分解(decompose)すると、その構造がよくわかります。上述したように、 qc.initialize() 関数は非常に深い回路になります。実は、任意の量子状態を効率的に量子回路にマッピングすることを目的とした状態準備の問題は、まだ未解決の研究課題です。

# Draw the decomposition of the amplitude embedding circuit
qc.decompose(reps=8).draw(output="mpl", fold=40)
_images/60e52bd6c31fac4ef2b15096e96b7e0b6df9ebc35cfe12325fe857af61211372.png

初期状態が準備できたら、VQCの変分ansatzを構築する必要があります。このansatzは通常、パラメーター化された回転ゲートとエンタングルレイヤーの集合からなり、所定の回数繰り返されます。回転のパラメーターを調整することで、ヒルベルト空間の特定の領域を探索することができます。

目的の解に収束する適切なansatzを選択することは簡単な問題ではありません。一般に、考慮すべきトレードオフがいくつかあります:

  • スピード: 探索空間を減らすことで、ゲートの数とansatzの深さを減らし、アルゴリズムを高速に実行します。

  • 精度: 探索空間の縮小は、問題の実際の解を除外してしまう危険性があり、最適解を下回るリスクがあります。

  • ノイズ: より深い回路はノイズの影響を受けるので、ansatzの接続性、ゲート、ゲートの fidelity を実験する必要があります。

変分量子回路のansatzの設計方法について詳しく知りたい場合は、IBM Learning Platform の Variational algorithm design コースを受講することを強くお勧めします。

Qiskit circuit library を確認すると、TwoLocalEfficientSU2RealAmplitudes など、ansatzes を構築するための多くのアーキテクチャがあることがわかります。この Exercise では、RealAmplitudes ansatzを選択します。

Exercise 2: ansatzの構築#


Qiskit circuit library の ```RealAmplitudes``` 回路を使用して以下のansatzを構築してください。

<img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABQoAAAHDCAYAAABs9Sc3AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/H5lhTAAAACXBIWXMAAA9hAAAPYQGoP6dpAAB8MElEQVR4nO3dd3xUVf7/8fdMekKANEhIKKEXadJFQREFBLErLva+grAuC7u6NuyIqyisLCrurqhIERVBBKRIBKSj1IAhCaQBIQnpbWZ+f/gza74EmUlm5mYmr+fj4WOXmXvueU9yJ3POZ+6512Sz2WwCAAAAAAAA0KCZjQ4AAAAAAAAAwHgUCgEAAAAAAABQKAQAAAAAAABAoRAAAAAAAACAKBQCAAAAAAAAEIVCAAAAAAAAAKJQCAAAAAAAAEAUCgEAAAAAAACIQiEAAAAAAAAAUSgEAAAAAAAAIAqFAAAAAAAAAEShEAAAAAAAAIAoFAIAAAAAAAAQhUIAAAAAAAAAolAIAAAAAAAAQBQKAQAAAAAAAIhCIQAAAAAAAABRKAQAAAAAAAAgCoUAAAAAAAAARKEQAAAAAAAAgCgUAgAAAAAAABCFQgAAAAAAAACiUAgAAAAAAABAFAoBAAAAAAAAiEIhAAAAAAAAAFEoBAAAAAAAACAKhQAAAAAAAABEoRAAAAAAAACAKBQCAAAAAAAAEIVCAAAAAAAAAKJQCAAAAAAAAEAUCgEAAAAAAACIQiEAAAAAAAAAUSgEAAAAAAAAIAqFAAAAAAAAAEShEAAAAAAAAIAoFAIAAAAAAAAQhUIAAAAAAAAAolAIAAAAAAAAQBQKAQAAAAAAAIhCIQAAAAAAAABRKAQAAAAAAAAgCoUAAAAAAAAARKEQAAAAAAAAgCgUAgAAAAAAABCFQgAAAAAAAACiUAgAAAAAAABAFAoBAAAAAAAAiEIhAAAAAAAAAFEoBAAAAAAAACAKhQAAAAAAAABEoRAAAAAAAACAKBQCAAAAAAAAEIVCAAAAAAAAAKJQCAAAAAAAAEAUCgEAAAAAAACIQiEAAAAAAAAAUSgEAAAAAAAAIAqFAAAAAAAAAEShEAAAAAAAAIAoFAIAAAAAAAAQhUIAAAAAAAAAolAIAAAAAAAAQBQKAQAAAAAAAIhCIQAAAAAAAABRKAQAAAAAAAAgCoUAAAAAAAAARKEQAAAAAAAAgCgUAgAAAAAAABCFQgAAAAAAAACiUAgAAAAAAABAFAoBAAAAAAAASPI1OgAAeIIdO3bYvW12draWLVumG2+8UZGRkXa369evX22iAQAAAC7nyHhYqt2YmPEwYDzOKAQAJ8vOztb777+v7Oxso6MAAAAAhmBMDHgmCoUAAAAAAAAAKBQCAAAAAAAAoFAIAAAAAAAAQBQKAcDpQkNDNXLkSIWGhhodBQAAADAEY2LAM5lsNpvN6BAAUN85epe32uAubwAAAKivGA8DDQNnFAKAk5WVlenEiRMqKyszOgoAAABgCMbEgGeiUAgATpacnKybbrpJycnJRkcBAAAADMGYGPBMvkYHAGAcm82myhLP+obPNyhAJpPJ6BgAAAAAvABzIqA6CoVAA1ZZUqaP291hdAyHjE/6SH7BgUbHAAAAAOAFmBMB1bH0GAAAAAAAAACFQgAAAAAAAAAsPQYAp+vcubO2b99udAwAAADAMIyJAc/EGYUAAAAAAAAAKBQCgLOlpqbqvvvuU2pqqtFRAAAAAEMwJgY8E4VCAHCykpIS7d+/XyUlJUZHAQAAAAzBmBjwTBQKAQAAAAAAAHAzEwCOix7UTSOXTa/2WEVRifKPZSpp6SYdmv+1bBarQekAAAAAwLWYE8FbUSgEUGvHliUobf1uyWRSUFRTtb9lqPpPv0dNOsRq69R5RscDAAAAAJdiTgRvQ6EQQK2d2ZesY58lVP078T+rdUPCW+r4hyu1+9WFKjuTb2A648TExGj69OmKiYkxOgoAAABgiIYyJmZOBG/DNQoBOE1lSZlO7z4qk9msxq2bGx3HME2aNNGoUaPUpEkTo6MAAAAAhmioY2LmRPB0FAoBOFVom18+DMvyCg1OYpzc3FwtWbJEubm5RkcBAAAADNGQx8TMieDJWHoMoNZ8g/wVEB5adT2OTnddrYjubXV691HlH8s0Op5hTp48qZkzZ6p79+4KCwszOg4AAADgdg1lTMycCN6GQiHgJNnZ2Xrttde0bNkypaWlKSoqSjfeeKNefvllTZo0SR988IFmz56tiRMnGh3VaXpPG6fe08ZVeyxl5Q/a9sT7BiVCfVBUXKGc/DL5+pgV2TRQfn4N6+T1ykqrsvNKVV5hVXiTADUK9jM6ktvl5ZfpbGG5ggN9FdE0UGazyehIcCOr1aYzeaUqLq1Uk0b+ato4wOhIAAC4DHMieBsKhYAT7N27V6NGjVJWVpZCQkLUtWtXZWRk6O2331ZSUpJycnIkSb169TI2qJMlLlijlK+2yuznq7DOrXTRhOsVEhMhS1l51TZD5z4umU367uE3qh7zb9pI1298Uzuf/1DHliXUtGt4GKvVpnXbMvTOokNavvG4rFabJKlpqL/uvb6D/nhrF3Vo7d3Xp0lJL9C8pYf1/rIjys4tlSSZTNLIwXGaMK6LRg6Ok4+P9xZNS8sqtXh1sv756SFt33+66vG2caH6461ddO/1HRTRNNDAhHC10zklmv/5Ef1ryWGlZvxvqdUlvZrp0du66Oar4hXg72NgQgAAnI85EbyN985YADfJzs7Wtddeq6ysLE2ZMkWZmZnavXu3srKyNGPGDK1cuVI7duyQyWRSjx49jI7rVPnHspSZsE/p6/do/ztfat3dryqyVzsNmvFw1TZbn3hPzfp1Uvz1g6seG/jyAzq1/TAfiF6ioKhcYyau0dUPf6Mv1qdWFQklKa+gXG8uOKBOY5fqH//dJ5vN9jt78lzvLT2s9mOW6NX5P1UVCSXJZpNWfZ+mMRPXatgDq5SbX2ZgStf5+Xi+ut/0ue5+alO1IqEkHUsr0NQ3tqvtNYu17ocMgxLC1VYlnFDbaxbribd2VisSStKWvad0xxPfqdctnyslvcCghAAAuAZzIngbCoVAHU2aNElpaWmaOHGiXn/9dYWGhlY9N23aNPXs2VOVlZVq06aNGjdubGBS1zu9M1FJSzcp/vrBiurbSZJUnleoLVPmasBLDyioeZhajx6o6Eu6aetf5xmc1nWCg4M1YMAABQcHGx3F5crKLRo76Vut+j7td7ez2aS//GO7/vHf/W5K5j7vf5aoh57fLIvl94ugm3ZlacQj36iouMJNydzjRFahht63Uj8fz//d7fILK3TNhNXatJNr9XibtVvTNXbSWhUWV/7udoeTz2rofSuVebrYTckAAEZqSGPi32JOBE9HoRCog0OHDmnRokWKjIzUK6+8UuM2ffr0kST17Nmz2uPJyckaO3asQkNDFRYWprvuuktnzpxxeWZX+/HNpbJWWtR76m1Vj6Vv2KuUr7ZoyJxJGvjqg9oyZa7Kcr33DmCtWrXS7Nmz1apVK6OjuNysj/Zr4w77Cz9T39iuw8l5rgvkZukni/THFzfbvf2O/dl68b29rgtkgD++uEUZp+wr/JRXWHX73zaqosLq4lRwl9KySo3/20ZVXqBQ/qvjmUV67JWtLk4FAKgPGtKY+P9iTgRPRqEQqIOFCxfKarVq/PjxatSoUY3bBAUFSapeKCwoKNAVV1yhtLQ0LVy4UO+++64SEhI0ZswYWa2ePYEuSMlS8peb1WJIDzUb0KXq8Z3TP1RofLTS1+9R2rrdBiZ0PYvFosLCQlksFqOjuJTFYtXcxYcdbjd30SEXpDHGu0sT7S6Q/Or9ZUdUWvb7Z155imNp+fo64YRDbTJOFWv5xlQXJYK7LVmTrNO/WW5vjy82pCr9ZJGLEgEA6ouGMiauCXMieDIKhUAdrF+/XpJ0xRVXnHebtLRflmT+tlD47rvvKj09XV988YXGjBmjW265RZ988ol++OEHLV++3LWh3eCntz6T1VL9G7TKkjIVpp5S7qHjBiZzj6NHj2rYsGE6evSo0VFcau3WjHOuRWaP/yw/6hWFMpvNpveWJTrcLju3VJ+v845C2fxlR1Sby07OW+p4gRn107tLHX8PWCw2ffDFERekAQDUJw1lTHw+DX1OBM9lsnnrleUBN2jZsqXS0tK0Z8+eGu9oXFlZqZiYGGVnZyspKUlt27aV9L/C4oYNG6pt365dO11++eWaP3++w1n69u2rrKwsh9r42cx61trf4b5qa+Rn03Vi7S4d+Ffti6HTzdtVYXL/WZc333yz3dueOnVKCxcu1O23365mzZrZ3W7p0qW1iWaYwoCBOhsyqlZtm+e9IV/rWScnci+rApQZ/mSt2oYWb1Dj0o3ODWSAM41uU6l/V4fb+VjyFH32TRckgrtlNp0qq7nmM+p/T1DZjwovWuaCRAAAV3FkPCzVbkxsxHi4oc2JoqOjtXPnzlr3De/na3QAwJMVFf2ydKqkpKTG5xctWqTs7GyFhoYqPj6+6vGDBw/qlltuOWf7bt266eDBg7XKkpWVpfT0dIfa+Jt8pOa16s4wGZkZKre5f/nCr79re/x6PJSUlDjUztHfn+Eii6SQ2jU9eTJbKj/l3Dzu5tNICq9d04LCYhWc9LDfd01aV0j+jjezWD3weEfNmtRucUpJSTnHAAB4GEfGtVLtxsRGfDYwJwKqo1AI1EF0dLRyc3O1e/duDRo0qNpzmZmZmjp1qiSpR48eMplMVc/l5uaqadOm5+wvPDxciYmOL+P6NYuj/GxmycMuidgipoUhZxSGhNhfEft1IBQUFORQu9jYWIdzGanI30d5tWwb0yxUZpufM+O4nU1mZdgsksnH4bZNQsxq5GG/75rkBlhUm/vX+ppK1dwLXj+kLJXIIsfvZhkSaFVTjgEA8CiOjGul2o2JjRgPN7Q5UW3mjWhYKBQCdTB8+HAdOnRIM2bM0FVXXaWOHTtKknbs2KE777xT2dnZklTjsmRnq83p4xXFpfq43R0uSFOzb256ts77OHL0iPyCA52QxjE7duywe9vDhw9r4cKFGjVqlDp37mx3u1mzZtUimXHO5JUq7qpPVVrm2LeZIy6J1Tf/+tlFqdxr3LT1WvRNskNtzGaT9m3+j1pGO75cs75ZsyVNIx5Z7XC76Y+P1pMP/t0FieBuT83eqZfe+9HhdqsWPq/L+rzrgkQAAFdxZDws1W5MbMR4mDkRUB03MwHqYNq0aYqIiNCJEyfUrVs3de/eXR06dFD//v3Vtm1bDRs2TFL1G5lIUlhYmPLy8s7ZX05OjsLDa7mWEfVG+/bttXr1arVv397oKC4V0TRQt41o63C7R2/rcuGNPMSjtzr+WsZe3sorioSSNHxgrNq1DHWojZ+vWfff2MlFieBuD93cWWaz6cIb/ka3dk116cUetsYLAOCwhjImBrwNhUKgDuLi4pSQkKDRo0crMDBQKSkpCg8P17x587Ry5UodOfLLXR3/b6GwS5cuNV6L8ODBg+rSxXuKKA2Vr6+vwsLC5Ovr/SdtP/lATzUNtf8idUP7Rmv0kJYuTORel/WJ1rVDW9m9fUiQr559pLcLE7mX2WzSjD/1c6jNX+7uruYRQS5KBHdrFdNIk8fbf0Mbs9mkGY/3q3Y5DgCAd2pIY2LAm1AoBOqoS5cuWrFihQoKClRQUKBt27bpoYceUlFRkVJSUmQ2m3XRRRdVazNmzBh9//33SktLq3ps27ZtSkpK0rXXXuvulwAnS0tL05QpU6r9fr1VxzZNtGLO1XYVCwd0j9Lns4bLx8d7PnpMJpMWzrhcw/rHXHDbkCBffT5ruHp1jnBDMve56ap4vfP3S2RP3efBmzrpxcf6uD4U3Grmn/vrrmsvfLaIj49J7z17qUYPsb+4DgDwXA1pTAx4E++ZrQH1zIEDB2Sz2dShQwcFB1e/0PtDDz2kmJgYXXfddVqxYoWWLl2q22+/Xf3799d1111nUGI4S2FhoRISElRYWGh0FLcY3Lu5fvjoWt0xpp38/c79WImODNIzD/fW+vevUVjjAAMSulZIsJ9WzR2hlyf1Vcvocy/U7etj0q0j4rVlwbW6apB33rzhj7d10Zp5I3XlgBY1Pt+9Q5jmT79M854Z7PAyVdR/Pj5m/efFIfrX04PVtV3TGrcZcUmsvn13lO67oaN7wwEADNPQxsSAt+AcYMBF9u3bJ+ncZceS1LhxY61fv16TJ0/WuHHj5OvrqzFjxujNN9+U2Uz9Hp6nU3xTLXj5cr3xlwH6Yn2q/vKP7covqlB4E38dXz1OfjUUEL2Jv5+Pnnigp6be013fbE7T+L9tVH5RhZqG+ungFzcrJsrxu8J6muEDYzV8YKwSk/O0eku6npqzSwVFFYpsGqAfl97AUlMvZzKZ9PAtnfXQzZ20aVeWxk5aq/zCCjVp5Kedn16v9q0aGx0RAAAAdvDumRtgoN8rFEpSu3bttGLFChUWFiovL08fffSRoqKi3BkRcLqo8CA9eHNnhYb4SZKCAny9vkj4W76+Zo0Z2qrq9YcE+TWIIuFvdYpvqknju6nx//8ZBPj7UCRsQEwmk4b2jVFo8C+//0bBfhQJAQAAPAhnFAIucqFCoacK69pal7z+iPwaBako7bQ2PTZbYR1bavjHTyo/KUNrxr2g0jP58gny1+B/PKrIXu1ks9q0+5VPlLryB0lS36fvVJvrLlHOvmStv/c1g18RAAAAANjP3jnRpbMmKGZID5WdyZckZWz6STtfWCCJORHqLwqFgIusX7/e6Agucemsidr8+D+VcyBF7ccNU79n7tLPn25QflKGll81tWq7ix4ZK2t5hZZd8pgatWym0V+/oqwt+1WWW6idLyxQXuIJtRrV38BX4jpRUVGaPHkyZ4gCAACgwfLmMbG9cyJJOjB3uQ6+t/KcfXj7nAieq+GsBwNQZ+EXxauyuFQ5B1IkSUmLN6rl1X1l9j/3O4c21w1W4odrJEmFJ04pa8sBtRo1wJ1xDRMREaHx48crIsK77m4LAAAA2Mtbx8SOzIkAT8SRDMBujVo1U1jnVhq7dmbVY75BAQqODj9329hIFaadrvp34YlTComNdEtOo+Xn52v79u3q37+/Gjfm2lwAAABoeLx1TOzInEiSujxwjdqPG6ai9GztmbGwqsAI1FcUCgE45PSeo1p7+4tV/x63f76BaeqnjIwMPfnkk/rwww+9alAEAAAA2Mubx8T2zol2v7pQxSdzJZtNrUb11/CP/65llzymyuJSd0UFHMbSYwB2K0w9We2sQL9GQfINDlRxVs6526Znq1Hc/65H0qhlMxWlZ7slJwAAAAC4giNzouKsHMlmkyQdX7VdFQXFatK+hduyArVBoRCA3XIOpMhaYVHMkB6SpM73jFDyl5tlLa88Z9vUr7aq011XS/qlSBh9STcd/2a7W/MCAAAAgDM5MicKjvnfcuSoizsoIDxU+clZbssK1AZLjwE4ZNOEt3TprAnyf/VBFaRkadPEtxXWqdU52+1/50sNfvNR3bh1jmwWq7Y9+b7KcgoMSAwAAAAAzmPvnOjSWRMVFNVENotVlaXl2vjgP1RRUGxAYsB+FAoBOCTv8HGtGPnXC25XWVKm7x550w2J6p+AgAB16tRJAQEBRkcBAAAADOHNY2J750RrbnveDWkA52LpMYA6s1RUKiAsVGPXzlRgxIUvVNz36TvV/bEbVJZX6IZ07hcfH68FCxYoPj7e6CgAAACAIRramJg5EbwFZxQCqLPTOxO1pO8jdm+/84UF2vnCAhcmAgAAAAD3YU4Eb8EZhQDgZImJiRo8eLASExONjgIAAAAYgjEx4JkoFAKAk9lsNlVUVMhmsxkdBQAAADAEY2LAM7H0GGjAfIMCND7pI6NjOMQ3yPsuhgwAAADAGMyJgOooFAINmMlkkl9woNExAAAAAMAQzImA6lh6DAAAAAAAAIAzCgHA2dq0aaOFCxcqNjbW6CgAAACAIRgTA56JQiEAOFlgYKDatWtndAwAAADAMIyJAc/E0mMAcLLMzEy9+OKLyszMNDoKAAAAYAjGxIBnolAIAE529uxZLV++XGfPnjU6CgAAAGAIxsSAZ6JQCAAAAAAAAIBCIQAAAAAAAAAKhQAAAAAAAABEoRAAnM5sNqt3794ym/kTCwAAgIaJMTHgmXjHAoCTWa1W7dmzR1ar1egoAAAAgCEYEwOeiUIhAAAAAAAAAAqFAAAAAAAAACgUAgAAAAAAABCFQgBwutDQUI0cOVKhoaFGRwEAAAAMwZgY8Ey+RgcAAG8TGxur559/3ugYAAAAgGEYEwOeiTMKAcDJysrKdOLECZWVlRkdBQAAADAEY2LAM1EoBAAnS05O1k033aTk5GSjowAAAACGYEwMeCYKhQAAAAAAAAC4RiHQkNlsNlWWeNZSAN+gAJlMJqNjAAAAAPACzImA6igUAg1YZUmZPm53h9ExHDI+6SP5BQcaHQMAAACAF2BOBFTH0mMAAAAAAAAAnFEIAM7WuXNnbd++3egYAAAAgGEYEwOeiTMKAQAAAAAAAFAoBABnS01N1X333afU1FSjowAAAACGYEwMeCaWHgOAk5WUlGj//v0qKSkxOopbnckr1a6D2UpMOauCogpJUlFJpbbvO60eHcMUGMBHDgAAQEPRUMfEgKdj1gYAqLXktALNW3pYi745ppSMwnOezyso14Dxy+Xra9LFnSN1/40d9Ydr2qlRsJ8BaQEAAAAAv4dCIQDAYUkn8vX4a9u0YtNx2WwX3r6y0qbt+09r+/7TmvrGdj16Wxc983BvBQXyMQQAAAAA9QUzNAAOix7UTSOXTa/2WEVRifKPZSpp6SYdmv+1bBarQengSlarTf/89KD+Nmuniksra7WP/MIKvTr/Jy37NlX/fuEyXdKruZNTAgAAAK7FnAjeikIhgFo7tixBaet3SyaTgqKaqv0tQ9V/+j1q0iFWW6fOMzqeYWJiYjR9+nTFxMQYHcWpSssqddvUDVq+8bhT9nck9awuvXuF5jwxSI+O6+qUfQIAAKB+8NYx8f/FnAjehkIhgFo7sy9Zxz5LqPp34n9W64aEt9TxD1dq96sLVXYm38B0xmnSpIlGjRpldAynKiu3aOyktVq7NcOp+7XZpAkvb1V5hVV/uvMip+4bAAAAxvHGMXFNmBPB25iNDgDAe1SWlOn07qMymc1q3LrhLifNzc3VkiVLlJuba3QUp7nvmQSnFwl/6/GZ27RkTbLL9g8AAAD38sYxsT2YE8HTUSgE4FShbX75MCzLO/cOuA3FyZMnNXPmTJ08edLoKE6xePUxffJ1kkNtdiwcqxNrx2nHwrF2t3nkhc3Kyi52NB4AAADqIW8bEzuCORE8GUuPAdSab5C/AsJDq67H0emuqxXRva1O7z6q/GOZRseDE5w6U6IJL291uF10ZLDimoc41CbnbJn++OIWLXvzSplMJof7BAAAANyNORG8DYVCwEmys7P12muvadmyZUpLS1NUVJRuvPFGvfzyy5o0aZI++OADzZ49WxMnTjQ6qtP0njZOvaeNq/ZYysoftO2J9w1KBGd7/b/7lJ1b6rb+vlifqh9+OqVBPVmmAc9js9m0aVeWPvzqZ53I+uUMglbRjXTv9R11Sa9mFMAbgPIKi5Z9m6Jl61J1Jq9UgQE+6tYuTA/d3FntWzU2Op5b7D18Ru8vS9TPx/NVabEpOjJIf7imnUZcEisfHxYzAfA+zIngbSgUAk6wd+9ejRo1SllZWQoJCVHXrl2VkZGht99+W0lJScrJyZEk9erVy9igTpa4YI1Svtoqs5+vwjq30kUTrldITIQsZeVV2wyd+7hkNum7h9+oesy/aSNdv/FN7Xz+Qx1bllDTrlEPlJRWav7nR9ze7zuLDlEohMfZuCNTE17eooNJeec8N//zI+rRMVxzn7pEl/Ti2PZWcxcd0nNzd+tUTvUvV75OSNPM/+zTqEvj9P5zl6pFM8fOtvYUB5Ny9cBz32vrj6fOee7jlUmKjw3VP/7SXzdc2cb94QDAhZgTwdvwtR5QR9nZ2br22muVlZWlKVOmKDMzU7t371ZWVpZmzJihlStXaseOHTKZTOrRo4fRcZ0q/1iWMhP2KX39Hu1/50utu/tVRfZqp0EzHq7aZusT76lZv06Kv35w1WMDX35Ap7Yf9toPxODgYA0YMEDBwcFGR6mTJWuSlXO2zO39Ll6drNM5JW7vF6itrzYe19UPf1NjkfBXPx3J0bAHvtbqzWnuCwa3eXrOLj360pZzioS/ter7NA284yulZhS4MZl77DqYrUvuXFFjkfBXyekFuunP6/T+Z4luTAbASN4yJr4Q5kTwNhQKgTqaNGmS0tLSNHHiRL3++usKDQ2tem7atGnq2bOnKisr1aZNGzVu7N3Ljk7vTFTS0k2Kv36wovp2kiSV5xVqy5S5GvDSAwpqHqbWowcq+pJu2vrXeQandZ1WrVpp9uzZatWqldFR6mTV98YUNMorrNqwg+u5wDMkJufptqnrVVFpveC2ZeVW3TxlnZLTvK9Q1JB9sjJJL767165tT2QVafSENaq043jxFLn5ZRo9YY3OFpZfcFubTXr4hc1K2JXlhmQAjOYtY2JHMSeCp6NQCNTBoUOHtGjRIkVGRuqVV16pcZs+ffpIknr27Fn12K+Fxf79+ysgIMCrrlv145tLZa20qPfU26oeS9+wVylfbdGQOZM08NUHtWXKXJXleu8dwCwWiwoLC2WxWIyOUie7DmUb1/dB4/oGHPH2JwdVUmb/e72wuFL//PSgCxPBnWw2m16Z/6NDbQ4k5WllwgkXJXK//3x5VCfP2H8WuNVq0+v/3efCRADqC28ZE9cGcyJ4MgqFQB0sXLhQVqtV48ePV6NGjWrcJigoSFL1QuHPP/+szz77TNHR0erXr59bsrpLQUqWkr/crBZDeqjZgC5Vj++c/qFC46OVvn6P0tbtNjCh6x09elTDhg3T0aNHjY5Sa2cLynU0Nd+w/ncdPGNY34C9CorK9eFXjr/PP/jiiEpKK12QCO72/e6T2v9zrsPt3ll0yAVp3M9qtWnuYsdfy4pNJ3Q8k8kx4O28YUxcW8yJ4MkoFAJ1sH79eknSFVdccd5t0tJ+Wb7520LhkCFDlJmZqeXLl2v48OGuDWmAn976TFZL9W/QKkvKVJh6SrmHjhuYDPZKP1VkaP9pJ43tH7DHtn2nVVjseMEvN7+cs2a9xLc/ZNS6ndVqc3Ia9zueWVirL5WsVpvWb6/dzw4APAVzIngq7noM1EFqaqokqXXr1jU+X1lZqc2bN0uqXig0m51fo+/bt6+yshy75o+fzaxn1d/hvrK2HtB/Ym4+7/Nnj6brw7jbzvt8XXTs0FEVJvdf2+nmm8//ev+vU6d+uZj7qlWrtGvXLrvb3XDDDQ7ncpVyn2ipyR9rfG7HwrGKjvz9i1JHRwZV/e+JtePOu11WdrH63b78nMd/TkpRXFycA4nrl8ymf5bMTZSZlenRr6MuGsLPoMSvixR6/uP791x/0+0KrPDeMywawu9fkvKCR0qBgxxuZ7XaFNeqncy68HX96rMKn+ZSk0dr1fZPf35CT03Y5uREAFzJkfGwVLsxsRHj4YY2J4qOjtbOnTudnAjehEIhUAdFRb+c9VRSUvO1eRYtWqTs7GyFhoYqPj7epVmysrKUnp7uUBt/k4/U3EWBXCQjM0PlNvdf5+TX37U9fj0eSkpKHGrn6O/PpQJsUpOan4qODFZc8xC7duPrY7Z729+yVJbVr5+Ho0ItklmyWiye/TrqoiH8DBo1lUIvuFWNzpxKl4q99OciNYzfvyQ1PyMF1q5pZnqKJA8/q9C//LyfFRdyNvekzuZ68bEBeCFHxrVS7cbERnxmMCcCqqNQCNRBdHS0cnNztXv3bg0aVP2MgszMTE2dOlWS1KNHD5ffsCQ6OtrhNn42s+RhN15sEdPCkDMKQ0LsL3b9OhAKCgpyqF1sbKzDuVzFYgrW+c5PzcouvmD76Mgg+fqYVWmxKiv7/Be5P9++/M2liqpHPw9HZfr4yCrJ7OOjGA9+HXXREH4GFlO5smyVksmx4ZTJVq7ocKvMYd75c5Eaxu9fkkr98lSbK6r6VZ5Qs9gWTs/jbjaZlWXNl9Xc2IFGNslkUrNGRfIL9t5jA/BGjoxrpdqNiY0YDze0OVFt5o1oWEw2m83Dv8oEjDNp0iTNnj1bLVu21LfffquOHTtKknbs2KE777xTx44dU0VFhSZMmKA5c+bUuI/nnntO06dPlxFvxYriUn3c7g6391sX45M+kl9wLU/fqIMdO3bYvW1lZaUKCgoUGhoqX1/7Cwj17cY2rUd8quOZtbtW4Im14xTXPERpJ4vU8qpPHW4/eXw3zfrrwFr1XR/EDV+o9FPFim0WrLRvbzc6jiEays/gzic36qMVSQ61efCmTnr32UtdlKh+aCi/f4vFqvajlyglw7Ebc3z40hDdeW0HF6Vyr+lzd+u5uXscanN5vxhtmH+NixIBcBVHxsNS7cbERoyHmRMB1XEzE6AOpk2bpoiICJ04cULdunVT9+7d1aFDB/Xv319t27bVsGHDJFW/PiG8n6+vr8LCwhwqEtZHfbpGGth3hGF9A46YPL6bzGb7zxj39TFp4u1dXZgI7uTjY9af7ujmUJvYZsG65WrXXo7EnR68qZMaBTv2efe4gz8zAJ7JW8bEQENDoRCog7i4OCUkJGj06NEKDAxUSkqKwsPDNW/ePK1cuVJHjhyRRKGwoUlLS9OUKVOq7njtqa4cYMyyOLPZpMv7xRjSN+Covt2i9M7fL7FrW5NJev+5y9SjY7iLU8GdJo3vpjvHtLdr28aN/LRiztUKDPCeSXOLZiFa8vqV8vO1b1rxzMO9NfaKmm8CB8C7eMuYGGhoKBQCddSlSxetWLFCBQUFKigo0LZt2/TQQw+pqKhIKSkpMpvNuuiii4yOCTcqLCxUQkKCCgsdW4pW39wxur1Cgtw/mR0zpKVaRjdye79AbT18S2d99MpQhTcJOO82kWGBWjxzmO6+zjuWm+J/TCaT/v3CZZp2b3f5+51/aN2lbVNt/u8Y9ersfWdMj7w0Tqv/NeJ3b17VKNhXb04doOkTLnZjMgBG8pYxMdDQeM/XmUA9c+DAAdlsNnXs2FHBwcHnPL906VJJ0sGDB6v9u02bNurbt6/7ggLn0STUX3eMaa95Sw67td8J47q4tT/AGcaPbq8br2yjxauT9Z/lR5WwO0sWi00Bfma9P/0y3XJ1vAL8fYyOCRfx8TFrxuP9NeWu7vrgiyP67NsU7Tl0RharTYEBPlo552pd0T/G5Tc2M9IV/VsoedWtWrHpuN77LFFrtqSr0mKTv59Zb0wdoDvHtFfjRv5GxwQAABfAGYWAi+zbt0/S+Zcd33LLLbrlllu0ZMmSav8+301PACNMu7e7ggPd953Spb2ba/hA7oIJzxQU6Ku7r+ugDfOvUXREkKRfziS8Y0x7ioQNRLOIIP3t/p7asfA6RUf+cgxENAnQsAEtvLpI+CtfX7OuH9ZGK/85Qs3//3sgKixQE8Z1pUgIAICH4IxCwEUuVCjkhuPwBG3jGmvG4/302CtbXd5XUKCP/v3CEIduDAEAAAAAcB4KhYCLXKhQ6KnCurbWJa8/Ir9GQSpKO61Nj81WWMeWGv7xk8pPytCacS+o9Ey+2o8bpm4PjVaTDnHa+fyHOvjeyqp99H36TrW57hLl7EvW+ntfM/DVuEZUVJQmT56sqKgoo6M4xaO3ddEX61O1bluG3W2ysour/a89Zvypn9q3auxwPgAAANQ/3jYm/i1750SS1OnuEepy/yjZKi2yWW1aOfoJWcoqvH5OBM9FoRBwkfXr1xsdwSUunTVRmx//p3IOpKj9uGHq98xd+vnTDcpPytDyq6ZWbXfmpyRtfPgNdX/sxnP2sfOFBcpLPKFWo/q7M7rbREREaPz48UbHcBqz2aQl/ximK+7/Wj8m5tjVpt/tyx3qY/L4bpp4e9faxAMAAEA95G1j4t+yd07UckQ/tb3xMq0c/aQqCooVENFY1gqLJO+fE8FzcY1CAHYLvyhelcWlyjmQIklKWrxRLa/uK7P/ud855B5M1dmj6ZLV6uaUxsvPz9e3336r/Px8o6M4TVjjAK17b5T6dot0+r7/cnd3vTltQIO4fhcAAEBD4Y1jYsmxOdFFj47Vj28sUUXBL6tsys7ky9YA50fwLJxRCMBujVo1U1jnVhq7dmbVY75BAQqODjcwVf2TkZGhJ598Uh9++KEaN/aepbQRTQP13Qej9dScnZr10QHV9TKbEU0D9M8nL9GtI+IpEgIAAHgZbx0TOzInatohThE926nnn2+RT4CfkpZ8p0Pzv3ZnXMBhFAoBOOT0nqNae/uLVf8et3++gWngbsFBvnpj6kDdNDxef3xxs/YdzXV4HyaTdOuIeL3110FVd8UEAAAAPIW9cyKTr49CWzXTNzc8I/+mIRq57HkVpJ5U2re73BUVcBhLjwHYrTD1pEJi/7f01K9RkHyDA1WcZd916+A9Bvdurh+X3qCNH1yjW0fEy8/3wh8nUWGB+tv9PXTs61v16WvDKBICAADA4zgyJypKz9axz7+XzWpVWU6B0tftVlSfDu6MCziMMwoB2C3nQIqsFRbFDOmhzE0/qfM9I5T85WZZyyuNjgYDmEwmDe0bo6F9Y1RWbtG+oznadfCMjqSeVXFJpXx8TGoc4q+encLVp2uk2rUMZYkxAAAAPJojc6Jjn3+v2Ct6KWvzfvkE+iv6km7a988vDUgN2I9CIQCHbJrwli6dNUH+rz6ogpQsbZr4tsI6tTpnu/a3Xq7ef71d/k1D1Gpkf3V7ZKzW3f2qcvYnG5DavQICAtSpUycFBAQYHcVtAvx91LdblPp2izI6CgAAAOoBbx4T2zsnOjjvKw167WFd/92bstmk1JU/KHXFVgMSA/ajUAjAIXmHj2vFyL9ecLufF2/Uz4s3uj5QPRQfH68FCxYYHQMAAAAwjDePie2dE1nKKvT95DluSAQ4D9coBFBnlopKBYSFauzamQqMuPAdzfo+fae6P3aDyvIK3ZAOAAAAAFyLORG8BWcUAqiz0zsTtaTvI3Zvv/OFBdr5gnd+uyhJiYmJuu+++/TBBx+oU6dORscBAAAA3K6hjYmZE8FbcEYhADiZzWZTRUWFbDab0VEAAAAAQzAmBjwTZxQCDZhvUIDGJ31kdAyH+AZ538WQAQAAABiDORFQHYVCoAEzmUzyCw40OgYAAAAAGII5EVAdS48BAAAAAAAAcEYhADhbmzZttHDhQsXGxhodBQAAADAEY2LAM1EoBAAnCwwMVLt27YyOAQAAABiGMTHgmVh6DABOlpmZqRdffFGZmZlGRwEAAAAMwZgY8EwUCgHAyc6ePavly5fr7NmzRkcBAAAADMGYGPBMFAoBAAAAAAAAUCgEAAAAAAAAQKEQAAAAAAAAgCgUAoDThYeH6+6771Z4eLjRUQAAAABDMCYGPBOFQgBwMrPZLD8/P5nN/IkFAABAw8SYGPBMvGMBwMmys7P1/vvvKzs72+goAAAAgCEYEwOeiUIhAAAAAAAAAAqFAAAAAAAAACgUAgAAAAAAABCFQgBwutDQUI0cOVKhoaFGRwEAAAAMwZgY8Ey+RgcAAG8TGxur559/3ugYAAAAgGEYEwOeiTMKAcDJysrKdOLECZWVlRkdBQAAADAEY2LAM1EoBAAnS05O1k033aTk5GSjowAAAACGYEwMeCaWHgMNmM1mU2WJZ33D5xsUIJPJZHQMAAAAAF6AORFQHYVCoAGrLCnTx+3uMDqGQ8YnfSS/4ECjYwAAAADwAsyJgOpYegwAAAAAAACAQiEAAAAAAAAAlh4DgNN17txZ27dvNzoGAAAAYBjGxIBn4oxCAAAAAAAAABQKAcDZUlNTdd999yk1NdXoKDBAeYVFVqtNkmS12mSz2QxOBHcrLqmU5f8fA/z+AQANFWNiwDOx9BgAnKykpET79+9XSUmJ0VHgBgVF5Vq8Olnf7zmpXQezdfBYniyWX4pDmdklir7iE/XtFqm+XaN001Vt1KNjuMGJ4WxpWUVauCpJ2/ef1q6DZ5ScXlD1XMbpErUfvVh9ukZqQPco3T6qnWKigg1MCwCAezAmBjwThUIAAGrh5+P5mvXRfn341c8qKKo473anckr1dUKavk5I0/Pz9ujS3s01YVwX3TqircxmkxsTw9k27czUrI8PaPnG41XF4ZoknShQ0okCLV6drL/O2qEbhrXRn+7opkt6NXdjWgAAAODCKBQCcFj0oG4auWx6tccqikqUfyxTSUs36dD8r2WzWA1KB7iWxWLVGx/u19P/3K2ycovD7b/fc1Lf7zmpeUsTNX/6pWob19gFKeFKefll+vPr2/TvL4463Lay0qYla5K1ZE2yHr6ls2b+uZ9CQ/xdkBIAALgScyJ4KwqFAGrt2LIEpa3fLZlMCopqqva3DFX/6feoSYdYbZ06z+h4gNOlnyzSzVPW6YefTtd5Xxt3ZKr7jZ/rnb9foruv6+CEdHCH73dnady0DUo/VVznfc1bclirvj+hJa9fqf7do5yQDgAAuBtzIngbbmYCoNbO7EvWsc8SdGzpJh2Yu1wrRz+povRsdfzDlQqIaLhnScXExGj69OmKiYkxOgqcKDmtQJfevcIpRcJfFZdW6p6nN+mND/c5bZ9wnTVb0nTVw984pUj4q+OZRRr2wNfauCPTafsEAKA+aChjYuZE8DYUCgE4TWVJmU7vPiqT2azGrRvutbeaNGmiUaNGqUmTJkZHgZOcPFOi4Q+tUkpGoUv2P+X17Xpv6WGX7BvOsWXvSV3/p29VWub4cvMLKSqp1JiJa7Rjv/OK0AAAGK2hjomZE8HTUSgE4FShbX75MCzLc01BxRPk5uZqyZIlys3NNToKnMBms+mBZxN0LK3gwhvXwcRXtmr/0RyX9oHaOVtQrnHTNqik1PlFwl8VlVRq3LQNKio+/41xAADwJA15TMycCJ6MQiGAWvMN8ldAeKgCIhqraedWGvDyA4ro3landx9V/rGGu4zu5MmTmjlzpk6ePGl0FDjBgq9+1opNJxxqs2PhWJ1YO047Fo61u015hVX3PpOgykouel3f/OUf23Qiq8ihNrU5Bo6lFeiJt3c6Gg8AgHqpoYyJmRPB23AzE8CJsrOz9dprr2nZsmVKS0tTVFSUbrzxRr388suaNGmSPvjgA82ePVsTJ040OqpT9J42Tr2njav2WMrKH7TtifcNSgQ4V35huSbP+MHhdtGRwYprHuJwu50HsjV38SE99oduDreFa2zZe1LvLzvicLvaHgOzPzmoe8Z20MVdIx1uC6B+OHmmRFt/PKnC4ko1CvbVoJ7N1TwiyOhYblNZadWmXVnKyi6Rj49J7VqGqk/XSJlMJqOjAS7BnAjehkIh4CR79+7VqFGjlJWVpZCQEHXt2lUZGRl6++23lZSUpJycX5YU9urVy9igTpS4YI1Svtoqs5+vwjq30kUTrldITIQsZeVV2wyd+7hkNum7h9+oesy/aSNdv/FN7Xz+Qx1blmBEdMAuH634WXkF5Rfe0IlmLzyoCeO6ymxmQlUfzP7koNv7nPPpQX3w/BC39wugbnbsP603F+zX0rUpqvjN2eF+vmbdNLyNHr/zIq++w3nO2TL989ODmrfk8Dk3ferVOVyP3tpF917fUb6+LGqDd2FOBG/DX2nACbKzs3XttdcqKytLU6ZMUWZmpnbv3q2srCzNmDFDK1eu1I4dO2QymdSjRw+j4zpN/rEsZSbsU/r6Pdr/zpdad/eriuzVToNmPFy1zdYn3lOzfp0Uf/3gqscGvvyATm0/zAci6jWbzaZ3Fh1ye79HU/O1bluG2/vFubKyi/XZtylu73fhqmPKOVvm9n4B1N6/vziiQXd+pYWrjlUrEkpSRaVVn35zTJfc9ZXmL0s0KKFrJacVaOAdy/XMP3fXeGf4vYdz9NDzm3Xd5LVcixVehzkRvA2FQsAJJk2apLS0NE2cOFGvv/66QkNDq56bNm2aevbsqcrKSrVp00aNGzc2MKlrnd6ZqKSlmxR//WBF9e0kSSrPK9SWKXM14KUHFNQ8TK1HD1T0Jd209a/zDE7rOsHBwRowYICCg4ONjoI6+DExRweS8gzpe8FXPxvSL6pbsib5nAm/O5SWWfTZ2mS39wugdpZ9m6L7n02QxWL73e0sFpseeO57LV3jXe/vM3mluvqRb3Q0Nf+C236dkKY//G2jLBaux9sQNNQxMXMieDoKhUAdHTp0SIsWLVJkZKReeeWVGrfp06ePJKlnz55Vjy1dulQ33XSTWrdureDgYHXu3Fl///vfVVjo2XfG+vHNpbJWWtR76m1Vj6Vv2KuUr7ZoyJxJGvjqg9oyZa7Kcj37df6eVq1aafbs2WrVqpXRUVAH2/efNqzvHQeM6xv/Y+wxkG1Y3wDsV1lp1aRXt8r2+zXCaibN+EEVFd5TKPvHf/fr5+MXLhL+avnG4/pmc5oLE6G+aMhjYuZE8GQUCoE6WrhwoaxWq8aPH69GjRrVuE1Q0C8XsP5tofD111+Xj4+PXn75Za1atUp//OMfNXfuXI0cOVJWq+cOHgtSspT85Wa1GNJDzQZ0qXp85/QPFRofrfT1e5S2breBCV3PYrGosLBQFovF6Ciog10HjSvUJKacVUGRe6+NiHPtOnjGwL4pFAKe4Kvvjte41Pb3ZJ4u1pcbUl2UyL3Kyi16vxbLqY24tAfcryGPiZkTwZNRKATqaP369ZKkK6644rzbpKX98q3pbwuFX331lRYvXqzx48dr6NChmjx5subMmaPNmzfr+++/d21oF/vprc9ktVT/Bq2ypEyFqaeUe+i4gcnc4+jRoxo2bJiOHj1qdBTUwaFjeYb1bbP9UiyEcSwWqw4n5xnW/yED+wZgvw9reamI2rarb1ZvTtPp3FKH2636Pk2nzpS4IBHqk4Y+Jm7ocyJ4Lu56DNRRauov3wi3bt26xucrKyu1efNmSdULhVFR5971rm/fvpKk9PR0h3P07dtXWVlZDrXxs5n1rPo73FfW1gP6T8zN533+7NF0fRh323mfr4uOHTqqwuT+My5vvvn8r/f/OnXqlCRp1apV2rVrl93tbrjhBodzwXVONX5Y8m1R43M7Fo5VdOT5r7cTHRlU9b8n1o773X6ysovV7/bl5zw+avT1Cqj03DNOMpv+WTI3UWZWpuLi4oyO4zCr/GQLf+q8zzvrGDjf77+k1KLYuJYyyYH1jPWIp//+naGh/wwayus/1fhBydfx1/fNuh8UF3evCxK5V1FAXynkWofb2WxS94svlZ/lpAtSwVUcGQ9LtRsTGzEebmhzoujoaO3cudPJieBNKBQCdVRUVCRJKimp+VvRRYsWKTs7W6GhoYqPj//dfW3YsEGS1KVLl9/driZZWVkOFxj9TT5Sc4e7MlRGZobKbe5fvvDr79kevx4LJSUlDrWrTYEYLhRYft5PyejIYMU1D7ngLnx9zHZtV5Ps06elYg8+JkItklmyWiyeeWyb/KTw8z/tjmMgI92Dr+Hl6b9/Z2joP4OG8voDy2o1o6ooL/OOn0t4O6l2f+J06mSWVJbh3DxwKUfGtVLtxsRGvC+YEwHVUSgE6ig6Olq5ubnavXu3Bg0aVO25zMxMTZ06VZLUo0cPmUym8+4nPT1dTz/9tEaOHKlevXrVKoej/GxmycMuh9gipoUhZxSGhNg/Cv51IBQUFORQu9jYWIdzwXVO+9l0vqsEZmX//vWooiOD5OtjVqXFqqzs319adb59RUU2kb/Fc4+JTB8fWSWZfXwU44HHtk0mZdgsksmnxueddQycbz8mW5laeODP7Vee/vt3hob+M2gorz/Hp1i1WUAb6FusCC/4uZT6SbW6mqvNopioYJltnv8zaEgcGddKtRsTGzEebmhzotrMG9GwmGw2R+7RBeD/mjRpkmbPnq2WLVvq22+/VceOHSVJO3bs0J133qljx46poqJCEyZM0Jw5c2rcR2FhoS6//HJlZWVpx44diomJcUv2iuJSfdzuDrf05Szjkz6SX3Cg2/vdsWOH3dsePnxYd911lz788EN17tzZ7nb9+vWrTTS4yGOvbNWchQdr1fbE2nGKax6itJNFannVpw639/ExqWDrXQoK9Nzv8+KGL1T6qWLFNgtW2re3Gx2nVnre/Ll+OpJTq7Z1PQYG9ojS1o/G1qrv+sAbfv911dB/Bg3l9a/ZkqYRj6x2uN2qd0Zo5KWevyTbYrGq7TWLdTzTsTPNbh0Rr0Uzh7koFVzFkfGwVLsxsRHjYeZEQHXczASoo2nTpikiIkInTpxQt27d1L17d3Xo0EH9+/dX27ZtNWzYL4Og316f8LdKSkp07bXXKjk5WWvWrHFbkRCu0759e61evVrt27c3OgrqoE/XCMP67tYuzKOLhN6ib7dIw/ru09W4vgHYb/jAWLVv1dihNm3jQnX1Jd5xJp2Pj1kP32z/l6K/evQ2xy+zA8/DmBjwTBQKgTqKi4tTQkKCRo8ercDAQKWkpCg8PFzz5s3TypUrdeTIEUk1FworKip08803a+fOnVq1apW6du3q7vhwAV9fX4WFhcnXl0KPJxvQvZmBfZ97syO4n5G/B44BwDOYzSa9/9yl8vezb1rl52vW+89dKrP5/Jej8TSTx3dz6IuVB27sqCF9WPrYEDAmBjwThULACbp06aIVK1aooKBABQUF2rZtmx566CEVFRUpJSVFZrNZF110UbU2VqtV48eP17p16/Tll1+qf3/H77SF+iktLU1TpkxRWpoH34gA6tK2qfpfZEyx5p7rOhjSL6q75ep4BQXWfI1CVwoN8dMNV7Zxe78Aamdo3xh9MWu4QoJ+vxgSHOirz2ddqSv6t3BTMvcICfbTqndGaGCPC39m3n9DR819avDvXrcb3oMxMeCZKBQCLnTgwAHZbDZ16NBBwcHB1Z6bMGGClixZoscff1zBwcH64Ycfqv47ffq0QYnhDIWFhUpISFBhYaHRUVBHRiyN6tkpXIN6Gnc2I/4nrHGA/jCqndv7vXtsBzUK9nN7vwBqb9RlLbV/2Y2adm93RTQNqPacySRNvae79i+7UaOHtDIooWtFhgVqw/xr9O8XLqvx7MJrh7bSqndG6L3nLpWvL1PQhoIxMeCZ+CsNuNC+ffsk1bzseNWqVZKkV199VYMGDar238qVK92aE0DNbh0RrxbNgi+8oRM9fsdFnGlRj0wa382tSwR9fU2aOI5rdwGeqE1sqGY83l9pa8dp56fXKfL/FwxjIoP02p/7Kz4u1OCErhUY4Kt7ruuoHQuvU+Lym6tef3REoJbPvkojL43j8w0APACFQsCFfq9QmJKSIpvNVuN/99xzj5uTAqhJUKCv3n1msNv6Gz6whe4aywW/65MeHcM19Z7ubuvvqQd7qVN8U7f1B8D5AgN81adrpAL8f7l0QUMsjnVs06Tq9fv4MOUEAE/CVUUBF/q9QqGnCuvaWpe8/oj8GgWpKO20Nj02W2EdW2r4x08qPylDa8a9oNIz+br4iT+o9TUDZCmrkLXSot2vfqKMjT9Kkro+NEad7xmhyqJSLb9qqsGvCPh9o4e00t1jO+i/y4/a3SYru7ja/9qjUbCf3n/u0gY5oazvnvtjby3feFyHjuXZ3aY2x0CvzuF68oFeDqYDAADuZu+c6IoPpiq0ZbNq7dbf+5pOrNnJnAj1FoVCwIXWr19vdASnu3TWRG1+/J/KOZCi9uOGqd8zd+nnTzcoPymj2gfcyW2H9OObS2UpLVdY19Ya9fnzWtzrIVWWlOnguyuUsy9Z/Z+/x7gX4kJRUVGaPHmyoqK4a6m3ePtvA7XvaI52Hzpj1/b9bl/u0P59fEz66JWhat3Cu5elearAAF8t/ccwXXbPSuWcLbOrjaPHQLPwQC16bZj87LxzKgAA9Z03j4ntnRNtuG9m1f+P6NlOV33yd6Vv2CtJXj8ngudiNArAbuEXxauyuFQ5B1IkSUmLN6rl1X1l9j/3O4f09XtkKS2XJOUeOi6ZTAqMaOzOuIaJiIjQ+PHjFRERYXQUOEnjRv76Zu4I9eoc7vR9+/iYtOClobruitZO3zecp2u7MK2ZN1LhTQIuvLGDmoUHau27o9SxTROn7xsAAKN465jYkTnRb3W4fZiOfbZJ1opKN6QEao8zCgHYrVGrZgrr3Epj1/7vmzHfoAAFR/9+8aTDuCtUmHpShWkN427O+fn52r59u/r376/GjRtGcbQhiAoP0sb5o3XvM5v0+bpUp+yzeUSQ/vviEI0YHOeU/cG1+nSN1Pf/HaNx0zbopyM5TtnnxV0i9OlrV6hDa4qEAADv4q1j4trMiXwC/dX2+kv19fVPuyMiUCecUQjAIaf3HNXyq6ZW/VdR+PvX34q5tLt6TrlFGx95000JjZeRkaEnn3xSGRkZRkeBkzUJ9ddnb1ypT169vM5nlo0f3U4Hv7iJIqGH6dK2qXYsHKtnHu4tX9/aX0/Sz9esFyf20Q8fjaVICADwSt48JnZ0TtR6zECdPZahvMPH3ZQQqD3OKARgt8LUkwqJjaz6t1+jIPkGB6o4q+Yza5oP6qrBsx7VurteVX6S9w0Q0DCZTCbdfk07XTUoVvM/T9TcxYeVmlFoV1s/X7NuvqqNJozrqsG9m7s4KVzF389H0ydcrPGj2+mdRYf0n+VHdbag3K62YY39de/1HfXobV3UrqX3nF0BAEBD4eicSJI63n6lji70vuvXwztRKARgt5wDKbJWWBQzpIcyN/2kzveMUPKXm2UtP/c6G80HdtFlsx/T+nteU+5B5yzTBOqTyLBA/fW+nvrL3d21bluGvt9zUrsOZmtvYo7y8stVabEqKNBX7Vs2Vt9ukerTJUJjr2it5hFBRkeHk3Rs00Sz/jpQLz3WR199d1zb95/WroNndCApV4XFlTKbpZAgP3Vr11R9u0aqf/cojRnSSsFBDL8AAPBUjsyJJCm0TbQierbVurtfdXNSoHYYqQJwyKYJb+nSWRPk/+qDKkjJ0qaJbyusU6tzthv8j0fl4++nS9989H9tH5vN6fbwOj4+Zl19SZyuvoQlxA1VSLCfxo1qp3Gj2hkdBQAAuIG9cyLpl5uYpK7cporCEjenBGqHQiEAh+QdPq4VI/96we2WDX7MDWnqp4CAAHXq1EkBAc6/OyoAAADgCbx5TGzvnEiSdr/yiYvTAM7FzUwA1JmlolIBYaEau3amAiMufM2trg+N0cBXH1BpToEb0rlffHy8FixYoPj4eKOjAAAAAIZoaGNi5kTwFpxRCKDOTu9M1JK+j9i9/cF3V+jguytcmAgAAAAA3Ic5EbwFZxQCgJMlJiZq8ODBSkxMNDoKAAAAYAjGxIBnolAIAE5ms9lUUVEhm81mdBQAAADAEIyJAc/E0mOgAfMNCtD4pI+MjuEQ3yDvuxgyAAAAAGMwJwKqo1AINGAmk0l+wYFGxwAAAAAAQzAnAqpj6TEAAAAAAAAAzigEAGdr06aNFi5cqNjYWKOjAAAAAIZgTAx4JgqFAOBkgYGBateundExAAAAAMMwJgY8E0uPAcDJMjMz9eKLLyozM9PoKAAAAIAhGBMDnolCIQA42dmzZ7V8+XKdPXvW6CgAAACAIRgTA56JQiEAAAAAAAAACoUAAAAAAAAAKBQCAAAAAAAAEIVCAHA6s9ms3r17y2zmTywAAAAaJsbEgGfiHQsATma1WrVnzx5ZrVajowAAAACGYEwMeCYKhQAAAAAAAAAoFAIAAAAAAACgUAgAAAAAAABAFAoBwOlCQ0M1cuRIhYaGGh0FAAAAMARjYsAz+RodAAC8TWxsrJ5//nmjYwAAAACGYUwMeCbOKAQAJysrK9OJEydUVlZmdBQAAADAEIyJAc9EoRAAnCw5OVk33XSTkpOTjY4CAAAAGIIxMeCZWHoMNGA2m02VJZ71DZ9vUIBMJpPRMQAAAAB4AeZEQHUUCoEGrLKkTB+3u8PoGA4Zn/SR/IIDjY4BAAAAwAswJwKqY+kxAAAAAAAAAAqFAAAAAAAAAFh6DABO17lzZ23fvt3oGAAAAIBhGBMDnokzCgEAAAAAAABQKAQAZ0tNTdV9992n1NRUo6MAAAAAhmBMDHgmlh4DgJOVlJRo//79KikpMToK4BZWq02JKWe162C2dh7I1omsImXnlUqScvPLNW/JYfXpGqEeHcPl7+djcFrA+Ww2m1LSC7XzYLZ2HczWz8fzf/MeKNPbHx9Qn66R6t05QsFBDL+9UWWlVQeP5WnX/z8GMk4XVx0DeQXl+uDzI+rTNULd2oXJ15dzNdAwMCYGPBMjFQAAUCtn8kr17y+Oau7iQzqWVlDjNsWllXrkhc2SpIimAbrv+o565NbOahvX2J1RAZcoKq7QJ18n6Z3Fh7T3cE6N2xSXWjR5xg+SpJAgX40f3U6P3tZFPTtFuDMqXCTjVJHe+yxR736WqIxTxTVuU1RSqfufTZAktWgWrAdv7KQHb+qk2OYh7owKAIBdKBQCAACHVFRY9cr8H/Xy+z+qrNxid7szeWWa+Z99ev2/+3TXtR305rQBCmsc4MKkgGvYbDa991mipr25Q2cLyu1uV1RSqXeXJurdpYkadWmc5j0zWC2jG7kwKVyluKRST83Zqbc/OSiLxWZ3u4xTxZr+rz168b29mjiuq156rI9Cgv1cmBQAAMdQKATgsOhB3TRy2fRqj1UUlSj/WKaSlm7Soflfy2axGpQOgCvtO5Kju5/apD2Hz9R6Hzab9N/lR7X2h3S9+8xgjR7SyokJAdc6kVWo+59N0NqtGXXaz6rv03TRjcv0xl8G6L4bOspkMjkpIVxty96TuufpTTqaml/rfVgsNr318QGt2HRc/3lhiC69ONqJCQG4A3MieCsKhQBq7diyBKWt3y2ZTAqKaqr2twxV/+n3qEmHWG2dOs/oeIaJiYnR9OnTFRMTY3QUwKm+/SFd10/+VkUllU7ZX8apYo2ZuFZv/22gHvtDN6fsE3ClfUdydPUj3ygr2znX28ovrNADz32vfUdz9cbUATKbKRbWd4tXH9P4JzaqstL+swh/T9KJAl3xwNda8NJQjRvVzin7BOqLhjImZk4Eb0OhEECtndmXrGOfJVT9O/E/q3VDwlvq+IcrtfvVhSo7U/tv2j1ZkyZNNGrUKKNjAE61cUemxkxc69BSY3tNevUHmUwmTby9q9P3DTjLoWN5GvbgKmXnljp93299fEA2m02z/jqQMwvrsaVrknX7XzfKanVOkfBXlZU2/eFvG2U2m3TriLZO3TdgpIYyJmZOBG/DLbcAOE1lSZlO7z4qk9msxq2bGx3HMLm5uVqyZIlyc3ONjgI4xfHMQl032TVFwl899spWrd2a7rL9A3VRUFSuax5d7ZIi4a/e/uSg/rX4sMv2j7r5MfGM/vA35xcJf2WzSeOf2KjdB7Ndsn/ACA11TMycCJ6OQiEApwpt88uHYVleocFJjHPy5EnNnDlTJ0+eNDoKUGc2m00PPJeg/MIKh9rtWDhWJ9aO046FY+1uc/+zCcovtP/GEIC7THtjh1IyHPtcq817YOob25V8njuIwzgVFVbd8/QmVVQ6dq0xR4+Bykqb7n0mQeUVrvtSBnCnhjwmZk4ET0ahEHCS7OxsTZs2Te3bt1dgYKBatmypyZMnq6ioSPfff79MJpPmzJljdEyn8g3yV0B4qAIiGqtp51Ya8PIDiujeVqd3H1X+sUyj4wFwgvnLjtTqpg3RkcGKax6i6Mhgu9ucyCrS1De2O9wX4EobtmfoX0scP9OvNu+BopJKPfBcgmw215y1htp59YMftfdwjsPtanMM/HQkRy+/96PDfQEwDnMieBuuUQg4wd69ezVq1ChlZWUpJCREXbt2VUZGht5++20lJSUpJ+eXwWWvXr2MDepkvaeNU+9p46o9lrLyB2174n2DEgFwJovFqpfe2+vWPud/fkTPPtJbLZqFuLVf4HxefHevW/tbvz1TP/x0SoN6slytPigqrtA/Ptzv1j7f/Gi//nJPdzUK9nNrv6i70rJKLV6drHlLD2vf0VyVlFUqsmmgbhjWWn+8tYu6dww3OiJcgDkRvA1nFAJ1lJ2drWuvvVZZWVmaMmWKMjMztXv3bmVlZWnGjBlauXKlduzYIZPJpB49ehgd16kSF6zR6luna+34l7TzhQUqzSlQSEyELGX/Wzo4dO7jGjrvz9Xa+TdtpFv3vqe2N17m7sgAHPDN5jSHl1vWlcVi03ufJbq1T+B8Difnaf12958N8s6iQ27vEzVbuOqYzha495II+YUV+uTrJLf2ibrbtDNTra5epLuf2qQte0+poKhClZU2ZWWXaO7iw+px8+caN229SkorjY4KJ2NOBG9DoRCoo0mTJiktLU0TJ07U66+/rtDQ0Krnpk2bpp49e6qyslJt2rRR48aNDUzqfPnHspSZsE/p6/do/ztfat3dryqyVzsNmvFw1TZbn3hPzfp1Uvz1g6seG/jyAzq1/bCOLUuoabceLzg4WAMGDFBwsP1LjYD66N2lxhTs3v0s0WU3DAAcYVTRevHqZOWcLTOkb1T37mfG3GCGG9t4lu92Zuqqh7/R6Qvc8GjRN8kaO2ltg7kOZUMZEzMngrehUAjUwaFDh7Ro0SJFRkbqlVdeqXGbPn36SJJ69uxZ9VhCQoKGDx+umJgYBQQEKC4uTrfddpsOHfLsMwhO70xU0tJNir9+sKL6dpIklecVasuUuRrw0gMKah6m1qMHKvqSbtr613kGp3WdVq1aafbs2WrVqpXRUYBas9ls+m5XliF9Z5wq1jFu6IB64LudxrwHyius2vbTKUP6xv8UFldo5wFj7kK85/AZbu7kIQqLK3Tj4+tUXmHfzW6+/SFDL73bMK5D2VDHxMyJ4OkoFAJ1sHDhQlmtVo0fP16NGjWqcZugoCBJ1QuFubm56t69u95++22tWbNGM2bM0IEDBzRo0CClpaW5Jbur/PjmUlkrLeo99baqx9I37FXKV1s0ZM4kDXz1QW2ZMldlud57BzCLxaLCwkJZLA3j22J4p6QTBW5fbvdbOw+cNqxvQJLKyi366YjjN7Bwll2HjClQ4X/2HDojI+8rs+fwGeM6h90++TrJ4TOA/7XkUIM4q7Ahj4mZE8GTUSgE6mD9+vWSpCuuuOK82/xa+PttoXDs2LF68803dcstt2jo0KEaP368li1bprNnz+qzzz5zbWgXK0jJUvKXm9ViSA81G9Cl6vGd0z9UaHy00tfvUdq63QYmdL2jR49q2LBhOnr0qNFRgFozeoK6+xATZBjrYFKuKirtO0PIFXgPGI+/g7DH3MWOrwg6lVOqZd+mOD9MPdOQx8TMieDJuOsxUAepqamSpNatW9f4fGVlpTZv3iypeqGwJhEREZIkX9/avS379u2rrCzHlkj52cx6Vv1r1d/v+emtzxR//WD1nnqbVt/8nCSpsqRMhamnlHvoeJ323bFDR1WY3D9xu/nmm+3e9tSpX5aLrVq1Srt27bK73Q033OBwLsBVCgP6SSFjanxux8Kxio78/esNRUcGVf3vibXjzrtdVnax+t2+/JzH/znvv/rkjRsdSIz6JLPpnyVzE2VmZSouLs7oOLVS6ttWanx3jc+54z2w8pvvFLfkHvsD1zPecAzkB10hBV1e43PuOAaeff41/eOJ9fYHrme84Ri4EJtMygh7RjI5fv7Ng5Ne0F9KvnVBKtdxZDws1W5MbMR4uKHNiaKjo7Vz58469Q/vRqEQqIOioiJJUklJSY3PL1q0SNnZ2QoNDVV8fPw5z1ssFlmtVqWmpuqJJ55QdHS0br311lplycrKUnp6ukNt/E0+UvNa9LX1gP4Tc/6Bwtmj6fow7rbzPl8XGZkZKre5f/nCr79re/x6PJSUlDjUztHfH+BSEZ2kkJqfio4MVlzz8zz5f/j6mO3e9reKi0tVzHvCc4VaJLNktVg8929bozDpPPcgc8d7oLy80nN/dpJ3HAPNi6Sgmp9yxzFQUFisgiwP/dlJ3nEMXIjJXwqv3SK9wqJyFWZ61s/FkXGtVLsxsRHHCnMioDoKhUAdREdHKzc3V7t379agQYOqPZeZmampU6dKknr06CGTyXRO+6FDh1adcdi+fXutX79eUVFRtc7iKD+bWTJuVVWttIhpYcgZhSEh9g/wfx0IBQUFOdQuNjbW4VyAqxQFNFLeeZ7Lyi6+YPvoyCD5+phVabEqK7vmL1N+b1/BQf4K4z3hsTJ9fGSVZPbxUYyH/h7LfJvqfFcJdMd7wN/frCgP/dlJ3nEM5AcG63y3VXLHMRDaKEiNPfRnJ3nHMXAhNkkZtgrJ5Odw29AQH4/7/ToyrpVqNyY2Yjzc0OZEtZk3omGhUAjUwfDhw3Xo0CHNmDFDV111lTp27ChJ2rFjh+68805lZ/8yxejVq1eN7efPn6+8vDwlJydr5syZuvrqq7V58+Za3RmsNqePVxSX6uN2dzjczkhHjh6RX3Cg2/vdsWOH3dsePnxYCxcu1KhRo9S5c2e7282aNasWyQDXWLs1XVc//E2Nz9W0RO7/OrF2nOKahygru0Qtr/rU4f6fmvaInnhgrsPtUD/EDV+o9FPFiomOUdp+z7xJV3Jagdpes7jG59zxHrjj1pGaP/0lh9vVF95wDHy4/KjufmpTjc+54xh4+/Wndc91CxxuV194wzFgj1umrNPStSkOt9vw5Zvq0zXS+YFcyJHxsFS7MbER42HmREB13MwEqINp06YpIiJCJ06cULdu3dS9e3d16NBB/fv3V9u2bTVs2DBJ578+YadOnTRgwACNGzdO69atU0FBgV577TV3vgS3+uamZ3XgXxceWHu69u3ba/Xq1Wrfvr3RUYBau7hLhKH99+lqbP9Am9hGCmvsb1j/fQx+D0KGF3H6dPGsIlJD9ehtXS680f/R/6Iow48vd2BMXLOGMieC56JQCNRBXFycEhISNHr0aAUGBiolJUXh4eGaN2+eVq5cqSNHjki68I1MJKlp06Zq3769fv75Z1fHhov5+voqLCys1jemAeqDiKaBio8NNaz/hjCBQv1mMpnUt5txxyHvAeN1jm+i4EBjPssDA3zUpW1TQ/qGYy7vF6MhfexfymkySU8/3Mt1geoRxsSAZ6JQCNRRly5dtGLFChUUFKigoEDbtm3TQw89pKKiIqWkpMhsNuuiiy664H5OnTqlxMREtWvXzg2p4UppaWmaMmWK0tK8d5kNGoaxlzt+GQRnGNy7uSKaspwGxht7eWtD+m3RLFgXczaZ4Xx8zBoztKUhfY++rKV8fZmqeQKTyaTPZw1Xz07hdm3/9t8GacxQYz5f3Y0xMeCZ+PQBXOTAgQOy2Wzq0KGDgoODqz13xx136LnnntMXX3yhjRs36r333tPll18uX19fPf744wYlhrMUFhYqISFBhYWFRkcB6uSPt9p/jU1nevRWx5dxAa5w55j2Cgly/5kwD9/cWX5+DNPrg9osK3WGCeP4O+hJwpsEaNO/R+vhWzqf9yzU7h3C9PmsKzXx9q5uTmccxsSAZ2IEArjIvn37JNW87HjgwIH6+uuvde+992rUqFGaOXOmLrvsMu3du5dreACoNzrFN9WVA1q4tc+osEDddFUbt/YJnE+TUH/dMca9n8u+viY9cGNHt/aJ8xvSJ1rd2jV1a5+d45vo8n4xbu0Tdde4kb/+9fRgpX87TrOfGFRVMAwJ8lXCf0brx6U36PphbYwNCQB2oFAIuMjvFQonTpyo7du3Kzc3VyUlJTpy5IjmzZun1q2NWeIEAOfz+pT+8vU1ua2/mX/urwB/H7f1B1zIs4/0dutNTf52X0+1aBbitv7w+0wmk96cNtCtfc6aNlAmk/v+7sK5mjYO0MTbu1b93Wga6q9LL47mdwrAY3BVUcBFfq9Q6MnCurbWJa8/Ir9GQSpKO61Nj81WWMeWGv7xk8pPytCacS+o9Ey+ev/tdrW6uq9sFqskad+cL5T85WZJUt+n71Sb6y5Rzr5krb/Xe+/yDHiDXp0j9PcHemn6v/a4vK/RQ1rqrrGcVY36JSYqWG//bZDufPI7l/fVvUNYg7nJgSe5alCsHrq5k95dmujyvh64saNGDI5zeT8A6sbeOVFofLQumfmIApqEyCfQXye+3aWdzy+QbDZ1fWiMOt8zQpVFpVp+1VSjXxJQhUIh4CLr1683OoJLXDprojY//k/lHEhR+3HD1O+Zu/TzpxuUn5RR7QPuwDtfas+rCyVJwdHhun7TLGUk/KSynALtfGGB8hJPqNWo/ka9DJeKiorS5MmTFRUVZXQUwCmefLCnVm9J0w8/nba7TVZ2cbX/vZCYqGC9+8xgzrhAvTR+dDut2HRci75JtruNo++B4EBf/ffFIfL344za+uj1Kf313c4sJaactbuNo8dAh9aN9fqUAbXKB9RH3jwmtndO1O+Zu3R81XYdmv+1fAL8NGbVq8q8Yp/S1+/RwXdXKGdfsvo/f49xLwSoAYVCAHYLvyhelcWlyjmQIklKWrxR/Z69S8eWJZyzbXn+/wbFviGBMplMDaYAEBERofHjxxsdA3Aafz8frfznCF1x/9f66UiOXW363b7c7v1HNA3Qmn+NYLkl6i2TyaT/vjhUefnlWr0l3a42jrwHAgN89MVbw9WbOx3XW6Eh/lr77kgNuWelUjLsuzGDI8dA6xaN9O27o9Qk1H3L3AFX89YxsSNzIpvNJr/Gv9zY0ifQX2Y/X5WcynNjWsBxFAoB2K1Rq2YK69xKY9fOrHrMNyhAwdHhNW7f5f5r1PmeEQpuEaEtU+aq9Ey+u6IaKj8/X9u3b1f//v3VuHFjo+MAThHeJEAb5l+jsZPWavOek07bb6uYEH0zd6S6tG3qtH0CrhDg76Mv375Kdz75nZassf/MwgtpGuqvL94arqF9uXlFfdcyupG+/+8Yjfzjau3/Oddp++3arqm+mTtCLaMbOW2fQH3grWNiR+ZE25/5t6787xPqfNfV8m8Soh9nfaac/c77DAFcgZuZAHDI6T1HtfyqqVX/VRSefznNoflf6/PLJuvra/+uHpNuVEBYwxgAZ2Rk6Mknn1RGRobRUQCnCm8SoI3zr9HLk/rK36/uQ4j7b+ion5beSJEQHiPA30eLZl6hD56/TI0b+dV5f6OHtNT+ZTdSJPQgsc1DtP2TsZp6T3eZzXVbKWEySX+5u7t2LryOIiG8kjePie2dE3W+e6SSv/hei3s/pKX9/qi2N16mmCE93JwWcAyFQgB2K0w9qZDY/y2L8msUJN/gQBVn/f5SxNyDqSrOylH0Jd1cHRGAi/n6mvXEAz2169PrNPbyVrWaKA/q2UzfzB2h96dfxjI7eByTyaR7r++oA8tu0l3Xtq/VXbq7tmuqBS8P1Vezr1Jsc5bce5qgQF+99uf+2vzfMRo+sEWt9nHlgBba/OEYzZzSX0GBLPICPIkjc6LO947Qz4s3SpJKz+Qrfd1u5kSo9/hUAmC3nAMpslZYFDOkhzI3/aTO94xQ8pebZS2vPGfbJh3jdPZImiQptHVzhV8Ur7z//28Anu+iDuH68u2rdDyzUPOWHNbX35/Q/p9zVVlpq3H7Dq0b64p+MXrkls5chw1eIS46RP99aaj+8ZcB+vcXR/TZtyn68UiOSsssNW7fKiZEl/aO1kM3d9KQPtEN5rq93mxgz2Za++4oJSbn6V9LDmvt1nQdSj4rq/Xcv4Nms0md45voqoGxeuTWzuoc39T9gQE4hSNzooLUU4q9ord+/nS9fIMCFD34Ih34l/3XLwWMQKEQgEM2TXhLl86aIP9XH1RBSpY2TXxbYZ1anbNd36fuVKNWzWSrqJTVYtUPT87X2aP2XQAegOdoFdNIL03qq5cm9VVpWaV+OpKrE1mFKi23yN/PR1FhgerdOYIzB+G1IsMCNfXeHpp6bw9VVFh1KDlPPx/PV0lZpfx8zQprHKBencIVFR5kdFS4SKf4pnpz2kBJUlFxhfYm5ijzdLHKKiwK8PNRTFSwenUKV0hw3ZerA6gf7J0TfT95jga8dL+6PnCNzP6+OrF6p5K/2GxAYsB+FAoBOCTv8HGtGPnXC2637q5X3JCmfgoICFCnTp0UEBBgdBTArQIDfNW/e5T6d48yOgpgCD8/s3p0DFePjjXf5AveLyTYT4N7Nzc6BlAvePOY2N45Uc7+ZK267ik3JAKch2sUAqgzS0WlAsJCNXbtTAVGXPiOZn2fvlPdH7tBZXmFbkjnfvHx8VqwYIHi4+ONjgIAAAAYoqGNiR2dE3V9aIwGvvqASnMK3JAOsB9nFAKos9M7E7Wk7yN2b7/zhQXa+cICFyYCAAAAAPdxdE508N0VOvjuChcmAmqHMwoBwMkSExM1ePBgJSYmGh0FAAAAMARjYsAzUSgEACez2WyqqKiQzVbz3V8BAAAAb8eYGPBMLD0GGjDfoACNT/rI6BgO8Q3yvoshAwAAADAGcyKgOgqFQANmMpnkFxxodAwAAAAAMARzIqA6lh4DAAAAAAAA4IxCAHC2Nm3aaOHChYqNjTU6CgAAAGAIxsSAZ6JQCABOFhgYqHbt2hkdAwAAADAMY2LAM7H0GACcLDMzUy+++KIyMzONjgIAAAAYgjEx4JkoFAKAk509e1bLly/X2bNnjY4CAAAAGIIxMeCZKBQCAAAAAAAAoFAIAAAAAAAAgEIhAAAAAAAAAFEoBACnCw8P1913363w8HCjowAAAACGYEwMeCYKhQDgZGazWX5+fjKb+RMLAACAhokxMeCZeMcCgJNlZ2fr/fffV3Z2ttFRAAAAAEMwJgY8E4VCAAAAAAAAABQKAQAAAAAAAFAoBAAAAAAAACAKhQDgdKGhoRo5cqRCQ0ONjgIAAAAYgjEx4Jl8jQ4AAN4mNjZWzz//vNExAAAAAMMwJgY8E2cUAoCTlZWV6cSJEyorKzM6CgAAAGAIxsSAZ6JQCABOlpycrJtuuknJyclGRwEAAAAMwZgY8EwsPQYaMJvNpsoSz/qGzzcoQCaTyegYAAAAALwAcyKgOgqFQANWWVKmj9vdYXQMh4xP+kh+wYFGxwAAAADgBZgTAdWx9BgAAAAAAAAAhUIAAAAAAAAALD0GAKfr3Lmztm/fbnQMAAAAwDCMiQHPxBmFAAAAAAAAACgUAoCzpaam6r777lNqaqrRUQDArYqKK7T1x5MqKauUJJWWWbTvSI4qK60GJ4O7ZOeWav22jGrHQHJagWw2m8HJ4A42m00p6QX65vu0qmOgrNyi0zklBidzn/zCciXsyqr2HjicnCerteG9BxgTA56JpccA4GQlJSXav3+/SkoazqAYQMO170iO5i09rI07MnUo+Wy1yfCZs2XqcfPnCgr0Ua9OERozpKUeuLGTmkUEGZgYzmS12rRuW4bmf56orT+e0vHMomrPnzlbprbXLFZ4kwD17Rap20e11W0j2iookGmItygprdTi1clauCpJOw5kK+dsWbXns/PK1OzyT9QqJkSDejbTvdd11FWDYmU2mwxK7Fw2m0079mfr3aWHlbD7pI6knq32/JmzZepy3WdqFOyni7tE6MYrW+vusR3UtHGAQYndhzEx4Jn4hAYAAIDDvv0hXc//a48Sdp+84LYlpRZt/fGUtv54Ss/N3aNbro7XCxMvVtu4xm5IClewWm1677NE/ePDfTqamn/B7XPOlmnNlnSt2ZKuP8/cpgdv6qSnHuql0BB/N6SFKxQUleul937Ue58lnlMcrMnxzCIdz0zWom+S1a5lqKbc1V0P39LZowuGy75N0cvv/6hdB7MvuG1hcYU27crSpl1ZevLtXRo/up2en3CxoiOD3ZAUAOxHoRCAw6IHddPIZdOrPVZRVKL8Y5lKWrpJh+Z/LZuFZWYA4I3yC8v1l39s13ufJdaqfUWlVZ98naQvN6RqxuP99Mdbu3h0oaAhSjqRr/ueSdCmXVm1ap+bX67X/r1Pn35zTPOnX6bhA2OdnBCutn5bhu57NkGpGYW1ap90okCPvrRFn3ydpH+/METtW3nWlwanzpRowstbtHRtSq3aF5dW6r3PErV0bbJmPzFIf7imnUwm/g56GuZE8FYUCgHU2rFlCUpbv1symRQU1VTtbxmq/tPvUZMOsdo6dZ7R8QAATnYwKVfXTFhT6+LAbxWVVGriy1v11cbjWvrGlWoU7OeEhHC1z9el6I4nvlNxaWWd93U8s0hXPfSN/nZ/D708qS+FEg9gs9n09Jxdeum9H52yv+/3nFSPm5dpwUtDddNV8U7Zp6v98OMpjZ20VqdzS+u8r9z8ct3xxHf6OiFN/37hMvn7+TghIdyNORG8DTczAVBrZ/Yl69hnCTq2dJMOzF2ulaOfVFF6tjr+4UoFRHjWN8POFBMTo+nTpysmJsboKADgNPuO5GjIvSudUiT8rdVb0nX1w9+osLjCqfuF8y38Okk3T1nvlCLhb706/ydNfHkrNzyp52w2mybP+MFpRcJflZRadOvUDfp45c9O3a8rJOzK0vCHVjmlSPhbn3ydpBsfX6fyCotT92u0hjImZk4Eb0OhEIDTVJaU6fTuozKZzWrcurnRcQzTpEkTjRo1Sk2aNDE6CgA4RVpWka5+5BudybvwdchqY+uPp3Tj49/KwhKtemvt1nTd+ffvXHbn1ncWHdLz/9rjkn3DOV58d69mf3LQJfu2Wm26+6lNWr05zSX7d4YDP+dqzGNrVFTi3EL5r1ZuOqH7n03wqoJ5Qx0TMyeCp6NQCMCpQtv88mFYlufcM048SW5urpYsWaLc3FyjowBAndlsNj3wXIKysu2/a+WOhWN1Yu047Vg41u42a7dm6G0XFSFQN7n5Zbrn6U2yWOwvYNTmGHh+3l798OOp2kSEi23fd1rPzXWskOvoMWCx2HTvMwnKzXfNFxJ1UVFh1V1//075hfaf+Vyb98BHK5K06JtjtYlYLzXkMTFzIngyCoUAas03yF8B4aEKiGispp1bacDLDyiie1ud3n1U+ccyjY5nmJMnT2rmzJk6efLCdwIFgPru318c1eot6Q61iY4MVlzzEIfv5vnk2zt1JOWsQ23gen+euU0Zp4odalObY8BqteneZzaptMw1Z2yhdkrLKnXvM5scPpu0NsdA5uli/WnGD45GdLnX/v2Tdh8641Cb2v4dnPjKVp08Y/8XM/VZQxkTMyeCt6FQCDhJdna2pk2bpvbt2yswMFAtW7bU5MmTVVRUpPvvv18mk0lz5swxOqZT9Z42Trcf+Ldu3/+Brt/whrrcO1IpK3/Q+ntmGB0NAOAEJaWVmvrGdrf1V1pm0V9n7XBbf7iwXQez9Z8vj7qtv8PJZ/XOokNu6w8XNm/JYR1MynNbfx9+9bN27D/ttv4u5NSZEj0/z33L4s/klWn63N1u6w91x5wI3oZCIeAEe/fuVffu3TVz5kxlZWWpa9euqqio0Ntvv63bbrtNhw79MuDt1auXsUGdLHHBGq2+dbrWjn9JO19YoNKcAoXERMhSVl61zdC5j2vovD9Xa+fftJFu3fue2t54mbsjAwAcsGj1MeWcde8ywOUbj+tEFku16gsjinZzFx922bUQ4Rir1aZ/GnAM1Kdi8QdfHFF5hXuvn7pgRZLyC8svvCHqBeZE8DYUCoE6ys7O1rXXXqusrCxNmTJFmZmZ2r17t7KysjRjxgytXLlSO3bskMlkUo8ePYyO61T5x7KUmbBP6ev3aP87X2rd3a8qslc7DZrxcNU2W594T836dVL89YOrHhv48gM6tf2wji1LMCI2AMBORkzWrVab3l2a6PZ+ca7c/DJ98nWS2/v9+Xi+vv3BseXucI312zN0NDXf7f1++s0xnclz7p2Fa8NisepfSw67vd/C4gp9tKL+3wUav2BOBG9DoRCoo0mTJiktLU0TJ07U66+/rtDQ0Krnpk2bpp49e6qyslJt2rRR48aNDUzqeqd3Jipp6SbFXz9YUX07SZLK8wq1ZcpcDXjpAQU1D1Pr0QMVfUk3bf3rPIPTuk5wcLAGDBig4GDHrkkDAPVJ+ski7difbUjfn69LMaRfVLd6c5pKyyyG9P35ulRD+kV1Rv0eSsssWr3Z+GLx3sQcpWYYc4bz5+s9/z3QUMfEzIng6SgUAnVw6NAhLVq0SJGRkXrllVdq3KZPnz6SpJ49e553P6NGjZLJZNJzzz3niphu9eObS2WttKj31NuqHkvfsFcpX23RkDmTNPDVB7VlylyV5XrvsrJWrVpp9uzZatWqldFRAKDWdh00pkgoSYeSz6qo2P67i8I1dh107OYNTu37kHHHH/7HyL8D9eEYMPT1H8yWzebZS/Ab8piYORE8GYVCoA4WLlwoq9Wq8ePHq1GjRjVuExQUJOn8hcLFixdr7969rorodgUpWUr+crNaDOmhZgO6VD2+c/qHCo2PVvr6PUpb590XaLZYLCosLJTFYsxZGADgDEYWiaxWm348kmNY//iFkUWSn47kqryCz1EjVVRYDX0fGnn81YcMufnlSkn37CJSQx4TMyeCJ6NQCNTB+vXrJUlXXHHFebdJS0uTVHOhMD8/X3/605/0+uuvuyagQX566zNZLdW/QassKVNh6inlHjpuYDL3OHr0qIYNG6ajR913l0gAcLbk9AJD+z+WZmz/MPYYKCu3KPN0sWH9Q8o6U2zY0nOpfvwNMDqD0f3XVUMfEzf0ORE8l6/RAQBPlpr6y7VDWrduXePzlZWV2rx5s6SaC4V///vf1bFjR40fP1533HFHnbL07dtXWVlZDrXxs5n1rPo73FfW1gP6T8zN533+7NF0fRh323mfr4uOHTqqwuTeO89J0s03n//1/l+nTp2SJK1atUq7du2yu90NN9zgcC4AcJWckFukgItqfG7HwrGKjjz/NaeiI4Oq/vfE2nG/209WdrH63b78nMcfmzxFf3vE/r+hcL7Mpn+RzKE1PuesY+B8v39J6j/wUvlZjTuztaGrNIdLTSef93lXHwPpGScVFxfnQGLnOx16n+RX8zjfHX8Hb739DgVVHHEgsWs5Mh6WajcmNmI83NDmRNHR0dq5c6eTE8GbUCgE6qCoqEiSVFJSUuPzixYtUnZ2tkJDQxUfH1/tuZ07d+q9995zqJD0e7KyspSe7thFn/1NPlJzp3TvNhmZGSq3uf/b7V9/1/b49XgoKSlxqJ2jvz8AcKmWRVJAzU9FRwYrrnnIBXfh62O2a7ua5OVkKy+Pv4uGalQh+df8lDuOgVNZGVKF8ctPGyy/Uqnp+Z929TFgtVQYPzaKL5H8an7KHe+BnOxTUmH9+TvoyLhWqt2Y2IjfOXMioDoKhUAdREdHKzc3V7t379agQYOqPZeZmampU6dKknr06CGTyVT1nMVi0cMPP6yJEyeqW7duTsviKD+bWXL/yXl10iKmhSFnFIaE2D/A+3UgFBQU5FC72NhYh3MBgKvkBUrnm9ZlZf/+ktDoyCD5+phVabEqK7vmL9MutK/wpoEKCuHvopFOmct0vlvKOOsY+L39xDRvIrPtPNVquJzVFKDM33ne1ceAr6lczQ0eG53xt6j0PM+54+9gVHiw/JvUn7+DjoxrpdqNiY0YDze0OVFt5o1oWEw2T7+VEmCgSZMmafbs2WrZsqW+/fZbdezYUZK0Y8cO3XnnnTp27JgqKio0YcIEzZkzp6rdW2+9pZkzZ+rw4cNVN0ExmUx69tln3Xrn44riUn3crm5Lnt1tfNJH8gsOdHu/O3bssHvbw4cP66677tKHH36ozp07292uX79+tYkGAC4xf1miHnju+1q1PbF2nOKahyjtZJFaXvVp7faxZpziomt3Fg6c456nNum/y2t3bbG6HgNtWjRS8jeuWbIH+7W7ZnGtr5NX12PgzjHt9eHLQ2vVt7O89O5ePTWndqt/6vr6fX1Nyt9yl4IC68+5PY6Mh6XajYmNGA8zJwKq42YmQB1MmzZNEREROnHihLp166bu3burQ4cO6t+/v9q2bathw4ZJqn59wuzsbD399NN65plnVFlZqby8POXl5UmSSktLlZeXJ6vVw77SQjXt27fX6tWr1b59e6OjAECt9ekaaVjfzcIDFdv8/Nf+gnv06RphYN/GHX/4HyN/D0Yef//LYNzr79YurF4VCWuDMTHgmSgUAnUQFxenhIQEjR49WoGBgUpJSVF4eLjmzZunlStX6siRXy4+/NtCYVpamgoKCvTwww8rLCys6j9JmjFjhsLCwnT8OHfB8mS+vr4KCwuTr69nD+4ANGzd2oUpvIkxyz4vuzi62iU7YIwhfYxbnnbZxR52wTAvZewxYPzyyAE9ouTvZ8yUuT68/rpiTAx4JgqFQB116dJFK1asUEFBgQoKCrRt2zY99NBDKioqUkpKisxmsy666H93jWzfvr02bNhwzn+SdPfdd2vDhg1cN8LDpaWlacqUKUpLSzM6CgDUmp+fWfde18GQvh++xf7LNsB1enaKUP+Lotzeb2CAj+681phjD9WNH91OQYE+bu+3b7dIXVwPzioNaxygW66Ov/CGLvDQTZ0M6deZGBMDnolCIeAiBw4ckM1mU4cOHRQc/L/lU40aNdLll19+zn+S1KZNG11++eUKDOR6E56ssLBQCQkJKiwsNDoKANTJH2/r4vY+O7RurCsHtHB7v6jZhHHuPwb+cE07w85mRXVhjQM0/pp2bu/XiOPufB414O/gZRc3V/eO4W7v19kYEwOeiUIh4CL79u2TVH3ZMQAAnqRdy8a6fVRbt/b5xP09ZTaz7Li+uHVEvNq1DHVbf/5+Zk2566ILbwi3mXJ3d7cuv20bF6rbRrj3787vGdSzmS7vF+PWPp98oJdb+wOA36JQCLiIo4VCm83m1jseAwBgj7f+OlBRYe45033k4DjdY9ByZ9QsMMBX/35+iNx1ychnH+mtru3C3NMZ7NI5vqmmP3qx2/r74PnL6tVNPEwmk+ZPv1TBbsp0z3UdNPLSOLf0BQA1oVAIuIi3nlEY1rW1Rn/9iq7fNEtXffJ3BUQ0VvSgbrrj2Mcau3amAiMaV9u+SYdY3XHsY/V//p6qx7o+NEY3bpmtsWtnujk9AMBRUeFBmvvUJQ61ycouVtrJImVlF9vdpkmov959djA3MamHLusTrUl/6OZQm9ocA327RWravT0cjQc3+Mvd3R2+XmVtjoHH/tBVQ/u69+w9e7SNa6wZj/dzqE1tXn9ss2C9OXWAo/FgAHvnRKHx0br606c19tvXdd3GN9Vm7P8+T5kTob6qP1/VAF5m/fr1RkdwiUtnTdTmx/+pnAMpaj9umPo9c5d+/nSD8pMytPyqqdW2Nfn66JKZjyh11fZqjx98d4Vy9iVXKx56k6ioKE2ePFlRUe6/ADwAuMJNV8Vr+qMX69l3dtu1fb/blzu0/wB/H30xa7haRjeqTTy4wWt/7qcjqWe16nv7bkrg6DHQpkUjfTFruHx9OY+hPvL1NevzWVfq0rtXKjm9wK42jh4DIwfH6fUp/WsTzy0mjOuig0m5mrv4sF3bO/r6m4b6a+U/r1bTxt5zfU5vHhPbOye6dNZE/bxog45+sk4BEY117TczdGr7YRVn5Xj9nAiei09iAHYLvyhelcWlyjmQIklKWrxRLa/uK7N/zd859PrzLUr5aqsKjmW6MaXxIiIiNH78eEVERBgdBQCc5umHe7lk+WFIkK9WzLnK7dcAg2P8/Xz02RtXasyQlk7fd7uWodow/xrFNg9x+r7hPC2ahWjD/FHq0LrxhTd20DWXxemzN66Uv5/777BsL5PJpDlPXuKSm5tENA3Q2ndHqmcn7xo7euuY2JE5UXi31kpb98uXbGVn8pVzIEVtrnPsLH3A3SgUArBbo1bNFNa5lcaunamxa2fq2tUz5BsUoODoc+/KFtm7g6L6dNSh+V8bkNRY+fn5+vbbb5Wfn290FABwGpPJpGce6a2FMy532h1pe3YK1+YPx2j4wFin7A+uFRToq2VvDtczD/eWr69zlojfeGUbbfnwWrWJdd8NU1B7rVuEasuH1+qWq+Odsj8fH5OeeqiXvph1lYKD6v9iN7PZpDlPDtI7f79EIU7Ke2nv5tr28Vj17eZ9Z91565jYkTnRmZ+Oqd1NQ6raNevXSY3imrk7MuAQCoUAHHJ6z1Etv2pq1X8Vheded8UnyF8DX31AW/7yLwMSGi8jI0NPPvmkMjIyjI4CAE43blQ7HfziJt14ZZta78Pfz6zn/thb2z8Z63Vn0Hg7Pz+zpk+4WNs/Hquenc6dFNsrKixQn752hZa+MUzNIoKcmBCuFhkWqMWvD9Pi14epWXjtb3TUo2O4tn08Vi9M7CM/N95Vua5MJpP+eFsX7fvsRl05oEWt9xMS5KtZ0wbou3+PVruWzj9Lsz7w5jGxPXMiSfp+0hxFXtxBY9fOVP/p9yozYZ9sFoub0wKOqf9f2wCoNwpTTyokNrLq336NguQbHKjirJxq24W2jlaj2EiN/Ow5SZJ/4xDJbJJ/k0b6fvIcd0YGALhA84ggffbmldp/NEdzFx/Wh1/9rMLiigu2a92ikR65pbPuu74jxSEP17tLpHYvul7rt2fon58e0vKNx2W12i7Yrv9FUXr0ti66dUR8vbqzLRx3y9XxGjOkpZasSdY7iw5p277TF2xjNpt07dCWevS2Lho+MFZms+fevCg+LlTfvjdK2/ed1juLDunTb46prPzCBaBObZro0du66K5r23vV9QgbEnvnRJJUmHZaGx94verfV33yd2V896NbcgK1xaczALvlHEiRtcKimCE9lLnpJ3W+Z4SSv9wsa3llte3yDh/XpxfdX/XvXlNulX+TYG1/5j9uTgwAcKWLOoTrn3+/RK893k+7D53RroPZ2n3ojLLzSlVeYVVggI9axzRSn66R6tM1Qhe1D5OPj+ecOYTfZzabNHxgrIYPjNXpnBLtPJCtXYey9WNijgqKKmSx2hQc6KvO8U3Up2uk+naLVNs47zxzqqEKCvTVXWM76K6xHXQsLf+XY+Bgtg4nn1VxaaXMZpMah/ipR8dw9ekaqX7dIhUV7l1fEvTvHqX+3aP01l8HatfBX17/nsNnlJtfropKq4ICfNSuZWP16RqhPl0j1aVtU+7u7uHsnRNJUmBkE5WeyZdsNrW4vKeadIzTsc+/NyA1YD8KhQAcsmnCW7p01gT5v/qgClKytGni2wrr1MroWAAAA4UE++myPtG6rE+00VFgkKjwII26rKVGXeb8m53AM7SNa6y2cY1164i2RkcxRJNQfw0b0ELD6rAcGZ7D3jlRy6v7qvvE62WzWFV8Mlffjn9ZltJyAxID9qNQCMAheYePa8XIvzrUZu8/FrsoTf0UEBCgTp06KSCA5SQAAABomLx5TGzvnOjoJ+t09JN1bkgEOA9rPwDUmaWiUgFhoRq7dqYCIy68pKjrQ2M08NUHVJpT4IZ07hcfH68FCxYoPt45dwQEAAAAPE1DGxMzJ4K34IxCAHV2emeilvR9xO7tD767QgffXeHCRAAAAADgPsyJ4C04oxAAnCwxMVGDBw9WYmKi0VEAAAAAQzAmBjwThUIAcDKbzaaKigrZbDajowAAAACGYEwMeCaWHgMNmG9QgMYnfWR0DIf4BnnfxZABAAAAGIM5EVAdhUKgATOZTPILDjQ6BgAAAAAYgjkRUB1LjwEAAAAAAABwRiEAOFubNm20cOFCxcbGGh0FAAAAMARjYsAzUSgEACcLDAxUu3btjI4BAAAAGIYxMeCZWHoMAE6WmZmpF198UZmZmUZHAQAAAAzBmBjwTBQKAcDJzp49q+XLl+vs2bNGRwEAAAAMwZgY8EwUCgEAAAAAAABQKAQAAAAAAABAoRAAAAAAAACAJJPNZrMZHQIAvElmZqbmzZunhx9+WDExMUbHAQAAANyOMTHgmSgUAgAAAAAAAGDpMQAAAAAAAAAKhQAAAAAAAABEoRAAAAAAAACAKBQCAAAAAAAAEIVCAAAAAAAAAKJQCAAAAAAAAEAUCgEAAAAAAACIQiEAAAAAAAAAUSgEAAAAAAAAIAqFAAAAAAAAAEShEAAAAAAAAIAoFAIAAAAAAAAQhUIAAAAAAAAAolAIAAAAAAAAQBQKAQAAAAAAAIhCIQAAAAAAAABRKAQAAAAAAAAgCoUAAAAAAAAARKEQAAAAAAAAgCgUAgAAAAAAABCFQgAAAAAAAACiUAgAAAAAAABAFAoBAAAAAAAAiEIhAAAAAAAAAFEoBAAAAAAAACAKhQAAAAAAAABEoRAAAAAAAACAKBQCAAAAAAAAEIVCAAAAAAAAAKJQCAAAAAAAAEAUCgEAAAAAAACIQiEAAAAAAAAAUSgEAAAAAAAAIAqFAAAAAAAAAEShEAAAAAAAAIAoFAIAAAAAAAAQhUIAAAAAAAAAolAIAAAAAAAAQBQKAQAAAAAAAIhCIQAAAAAAAABRKAQAAAAAAAAgCoUAAAAAAAAARKEQAAAAAAAAgCgUAgAAAAAAABCFQgAAAAAAAACiUAgAAAAAAABAFAoBAAAAAAAAiEIhAAAAAAAAAFEoBAAAAAAAACAKhQAAAAAAAABEoRAAAAAAAACAKBQCAAAAAAAAEIVCAAAAAAAAAKJQCAAAAAAAAEDS/wOWI2u/tBHcVwAAAABJRU5ErkJggg=="/>
Hide code cell content
num_qubits = 5 # Add your code here
reps = 1 # Add your code here
entanglement = 'full' # Add your code here
ansatz = RealAmplitudes(num_qubits=num_qubits, entanglement=entanglement, reps=reps) # Add your code here
# Add code here to draw the ansatz circuit

Step 2: 問題の最適化#

Step 3: Qiskit primitivesを使って実行#

Lab 2 で紹介した Qiskit Patterns の Step 2 では、量子ハードウェアに回路を効率的にマッピングするためのトランスパイル処理を扱います。
Transpilation には、回路の各論理量子ビットをデバイス内の実際の物理量子ビットにマッピングすることと、個々の量子ゲートを量子ハードウェア上で実行される演算シーケンスに変換することが含まれます。場合によっては、heavy-hex パターンに従うデバイスのトポロジーに適応させるために、SWAPゲートを追加する必要があります。

そして Step 3 では、Qiskit primitives を使って最終的な回路を実行します。今回は、VQC のコスト関数を構築するためにEstimator primitive を使用します。Estimator primitive は、複雑な代数プロセスを用いて、オブザーバブルを既知の固有基底を持つ他のオブザーバブルの組み合わせに分解することで、実際の量子デバイス上の期待値を推定します。

コスト関数は、データセットの各鳥の回路の出力に対するオブザーバブル \(\hat{O}\) の期待値で定義されます:

\(C(\theta) = \sum_{i\in \text{birds}} (\langle \psi_i(\theta) \vert \hat{O} \vert \psi_i(\theta) \rangle - L_i),\)

ここで、 \(\psi_i(\theta)\) は鳥 \(i\) の回路の出力状態、\(L_i\) は同じ鳥のラベルです。オブザーバブルは \(\hat{O} = ZZZZZ\) であり、 \(\theta\) は ansatz のパラメータベクトルです。

VQC の訓練に成功した場合、つまりコスト関数を最小化する最適なパラメーターセット \(\theta^{opt}\) に到達した場合、VQC はIBM Quantumの鳥の場合は期待値\(\langle ZZZZZ \rangle = 1\)、非IBM Quantumの鳥の場合は期待値 \(\langle ZZZZZ \rangle = 0\) を出力します。

# Define the observable 
obs = SparsePauliOp("ZZZZZ")
# Define the estimator and pass manager
estimator = StatevectorEstimator() #To train we use StatevectorEstimator to get the exact simulation
pm = generate_preset_pass_manager(backend=AerSimulator(), optimization_level=3, seed_transpiler=0) 
# Define the cost function
def cost_func(params, list_coefficients, list_labels, ansatz, obs, estimator, pm, callback_dict):

    """Return cost function for optimization

    Parameters:
        params (ndarray): Array of ansatz parameters
        list_coefficients (list): List of arrays of complex coefficients
        list_labels (list): List of labels
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        obs (SparsePauliOp): Observable
        estimator (EstimatorV2): Statevector estimator primitive instance
        pm (PassManager): Pass manager
        callback_dict (dict): Dictionary to store callback information

    Returns:
        float: Cost function estimate
    """

    cost = 0
    for amplitudes,label in zip(list_coefficients, list_labels):
        qc = QuantumCircuit(num_qubits)
        # Amplitude embedding
        qc.initialize(amplitudes)
        # Compose initial state + ansatz
        classifier = qc.compose(ansatz) 
        # Transpile classifier
        transpiled_classifier = pm.run(classifier)
        # Transpile observable
        transpiled_obs = obs.apply_layout(layout=transpiled_classifier.layout)
        # Run estimator
        pub = (transpiled_classifier, transpiled_obs, params)
        job = estimator.run([pub])
        # Get result
        result = job.result()[0].data.evs
        # Compute cost function (cumulative)
        cost += np.abs(result - label)

    callback_dict["iters"] += 1
    callback_dict["prev_vector"] = params
    callback_dict["cost_history"].append(cost)

    # Print the iterations to screen on a single line
    print(
        "Iters. done: {} [Current cost: {}]".format(callback_dict["iters"], cost),
        end="\r",
        flush=True,
    )
    
    return cost

Step 4: 古典的な後処理#

Qiskit Pattern の最後のステップでは、結果に対して古典的な後処理を行います。変分量子アルゴリズムの場合、このステージでは、定義したコスト関数を最小化するための ansatz パラメーターの古典最適化が含まれます。これはVQCを 訓練する 部分です。

訓練を開始するには、パラメーターの初期セットを定義する必要があります。この初期点の選択は非常に重要なステップです!初期点を誤ると、トレーニングのプロセスが深刻な問題になる可能性があります。特に、オプティマイザが局所的最小値にはまり込んでしまうと、最適解に収束できなくなります。

ここでは、params_0_list.npy にある VQC の初期パラメーターの5つの設定を使用します。

# Intialize the lists to store the results from different runs
cost_history_list = []
res_list = []

# Retrieve the initial parameters
params_0_list = np.load("./data/params_0_list_2024.npy")

for it, params_0 in enumerate(params_0_list):

    print('Iteration number: ', it)
    
    # Initialize a callback dictionary
    callback_dict = {
        "prev_vector": None,
        "iters": 0,
        "cost_history": [],
    }

    # Minimize the cost function using scipy
    res = minimize(
        cost_func, 
        params_0, 
        args=(list_coefficients, list_labels, ansatz, obs, estimator, pm, callback_dict), 
        method="cobyla", # Classical optimizer
        options={'maxiter': 200}) # Maximum number of iterations
    
    # Print the results after convergence
    print(res)

    # Save the results from different runs
    res_list.append(res)
    cost_history_list.append(callback_dict["cost_history"])
Iteration number:  0
 message: Maximum number of function evaluations has been exceeded.
 success: False
  status: 2
     fun: 4.000011467451044
       x: [ 3.517e+00 -8.492e-04  5.588e-04  3.655e+00  4.717e+00
            1.571e+00  1.804e-03  6.281e+00  2.629e+00  3.145e+00]
    nfev: 200
   maxcv: 0.0
Iteration number:  1
 message: Maximum number of function evaluations has been exceeded.
 success: False
  status: 2
     fun: 2.9486038843201667e-08
       x: [-1.309e-07 -1.571e+00  1.571e+00 -1.104e+00  1.934e+00
            1.571e+00  1.571e+00 -2.869e-05  1.571e+00 -2.305e-05]
    nfev: 200
   maxcv: 0.0
Iteration number:  2
 message: Optimization terminated successfully.59]]
 success: True
  status: 1
     fun: 4.113279744917959
       x: [ 4.050e+00  3.146e+00  3.213e+00  3.439e+00  1.043e+00
            2.321e+00  3.105e+00  3.063e+00  3.192e+00  2.586e+00]
    nfev: 145
   maxcv: 0.0
Iteration number:  3
 message: Maximum number of function evaluations has been exceeded.
 success: False
  status: 2
     fun: 4.4506198016960825
       x: [ 1.637e+00  9.279e-01  2.235e-01  6.615e-01  2.320e+00
            8.706e-01  1.651e+00  6.076e-01  1.083e+00 -4.224e-01]
    nfev: 200
   maxcv: 0.0
Iteration number:  4
 message: Maximum number of function evaluations has been exceeded.
 success: False
  status: 2
     fun: 1.6973675975627892e-06
       x: [ 3.142e+00  1.571e+00  1.571e+00  1.660e+00  1.510e+00
            1.571e+00  1.571e+00  1.292e-04  1.571e+00 -5.931e-04]
    nfev: 200
   maxcv: 0.0

最適化されたパラメーターのリストを保存し、各初期点のコスト収束履歴をプロットします。

fig, ax = plt.subplots(1, 1, figsize=(7,5))
ax.set_title('Cost history')
ax.set_ylabel('Cost')
ax.set_xlabel('Iterations')

# Add your code here
y = cost_history_list[4]
x = range(len(y))
plt.plot(x, y)
[<matplotlib.lines.Line2D at 0x7ff260f8c750>]
_images/f6d279e2283415fe0a3de51d5a19ff1a33a8d393ce2a5ada9eddd0c7463717da.png

Exercise 3: 最適なパラメータの検証#

どの初期パラメーターセットが最もよく収束するかをチェックする必要があります。そのために、Lab を通して使う2つの関数を定義してください。

1. 関数 `test_VQC` :データセットの各鳥に最適なパラメーターを持つ回路を適用し、コスト関数の収束値を出力する
2. 関数 `compute performance`:各最適パラメーターセットに対するトータルパフォーマンス(P)を出力します。


$P = 100 - 100\cdot\sum_{i\in \text{birds}} \frac{\vert (\langle \psi_i(\theta_{opt}) \vert \hat{O} \vert \psi_i(\theta_{opt})\rangle - L_i\vert}{2^5}$.
Hide code cell content
def test_VQC(list_coefficients, list_labels, ansatz, obs, opt_params, estimator, pm):
    
    """Return the performance of the classifier

    Parameters:
        list_coefficients (list): List of arrays of complex coefficients
        list_labels (list): List of labels
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        obs (SparsePauliOp): Observable
        opt_params (ndarray): Array of optimized parameters
        estimator (EstimatorV2): Statevector estimator
        pm (PassManager): Pass manager for transpilation 
        
    Returns:
        list: List of test results
    """
    ### Write your code below here ###
    results_test = []
    for coefficients, labels in zip(list_coefficients, list_labels):
        num_qubits = 5
        qc = QuantumCircuit(num_qubits)
        # Amplitude embedding
        qc.initialize(coefficients)
        # Initial state + ansats
        classifier = qc.compose(ansatz.decompose(reps=6)) 
        # Transpile classifier
        transpiled_classifier = pm.run(classifier) 
        # Run estimator
        pub = (transpiled_classifier, obs, opt_params)
        job = estimator.run([pub])
        # Get result
        result = np.abs(job.result()[0].data.evs)
        results_test.append(result)
        
    ### Don't change any code past this line ###                                 
    return results_test
Hide code cell content
def compute_performance(result_list, list_labels):
    
    """Return the performance of the classifier

    Parameters:
        result_list (list): List of results
        list_labels (list): List of labels
        
    Returns:
        float: Performance of the classifier
    """
    ### Write your code below here ###
    performance = 100
    for result, label in zip(result_list, list_labels):
        performance -= 100 * np.abs(result - label) / 2**num_qubits
        
    ### Don't change any code past this line ###                                    
    return performance
fig, ax = plt.subplots(1, 1, figsize=(7,5))
ax.set_title('Cost on test data')
ax.set_ylabel('Cost')
ax.set_xlabel('State index')

ax.plot(list_labels, 'k-', linewidth=3, alpha=0.6, label='Labels')

best_result_index = 0
best_performance = 0

for index in range(len(res_list)):
    opt_params = res_list[index].x
    results_test = test_VQC(list_coefficients, list_labels, ansatz, obs, opt_params, estimator, pm)
    performance = compute_performance(results_test, list_labels)
    print(f"Performance for trial {index}: {performance}")
    ax.plot(results_test, 'o--', label='Predictions trial '+str(index))

    if performance > best_performance:
        best_result_index = index

ax.legend()
Performance for trial 0: 93.97555178184847
Performance for trial 1: 99.99999990785612
Performance for trial 2: 90.85954362890645
Performance for trial 3: 86.57582851919992
Performance for trial 4: 99.99999469572624
<matplotlib.legend.Legend at 0x7ff2615dac90>
_images/3fd5df13447f4c3deca0f81d941ec7093defaa3d4dca4bdd51f40ddbf8bfa451.png
best_result_index
4

Part 2: ノイズがVQCの性能に与える影響#

VQCの訓練が完了したので、制御された方法でノイズを追加したときに、VQCが2つの鳥のセットを区別できるかどうかをテストしたいと思います。Lab のこのパートでは、ローカルテストに焦点を当てます。つまり、ローカルマシン上で量子ノイズをシミュレートします。このため、カスタムの fake backend を作成します。

その前に、量子ノイズの種類について少し説明しましょう:

  1. デコヒーレンス: 量子ビットは初期化された後、特にアイドル状態が続くと、時間とともに量子状態の情報を失います。通常、デコヒーレンスにはエラー率を用いず、T1とT2の時間を用います。これは、量子ビットが情報を失うまでの時間です。

  2. 測定エラー: 量子ビットを測定するとエラーが発生する可能性があり、つまり、0の代わりに1が測定されたり、その逆が起こったりします。これは古典的なチャンネルと同様に機能します。

  3. ゲートエラー: ゲートは完全ではないので、適用時にエラーを引き起こす可能性がわずかにあります。これは特にCX、CZ、ECRゲートのような2量子ビット・ゲートに当てはまり、通常1量子ビット・ゲートの約10倍のエラー率を持ちます。

  4. クロストーク・エラー: ある量子ビットにゲートを適用すると、他の量子ビット、特に隣接する量子ビットも影響を受ける可能性があります。これは、これらの量子ビットがアイドル状態であっても同様です。幸いなことに、最新のHeronデバイスでは、この問題は少なくなっていますが、それでも注意すべきことです。

IBM Quantum systems page では、量子デバイスのコヒーレンス時間やエラー率を直接見ることができます。以下はその概要です:

ノイズについてもっと知りたい方は、Qiskit video で詳しく説明しています。

さて、ノイズについて学んだところで、Qiskit を使ってカスタム fake backend を作成し、ノイズを作ってみたいと思います。このために、 GenericBackendV2 class クラスを使用します。このクラスでは、ノイズレベルやシミュレートしたいバックエンドの接続性を指定することができます。ここでは、 GenericBackendV2 のデフォルトの全対全接続性を仮定します。

fake_backend = GenericBackendV2(
    num_qubits=5,
    basis_gates=["id", "rz", "sx", "x", "cx"]
  )

汎用バックエンドをそのまま使ってもあまり意味がありません。なぜなら、最終的には、量子回路を動作させたい特定のバックエンドとできるだけ同じような動作をするバックエンドモックを作りたいからです。この Lab では、先ほど作成した汎用的なバックエンドモックのプロパティを変更し、ニーズに合わせる方法を紹介します。

そのために update_error_rate というヘルパー関数を作ります。

def update_error_rate(backend, error_rates):

    """Updates the error rates of the backend

    Parameters: 
        backend (BackendV2): Backend to update
        error_rates (dict): Dictionary of error rates

    Returns:    
        None
    """

    default_duration=1e-8
    if "default_duration" in error_rates:
        default_duration = error_rates["default_duration"]

    # Update the 1-qubit gate properties 
    for i in range(backend.num_qubits):
        qarg = (i,)
        if "rz_error" in error_rates:
            backend.target.update_instruction_properties('rz', qarg, InstructionProperties(error=error_rates["rz_error"], duration=default_duration))
        if "x_error" in error_rates:
            backend.target.update_instruction_properties('x', qarg, InstructionProperties(error=error_rates["x_error"], duration=default_duration))
        if "sx_error" in error_rates:
            backend.target.update_instruction_properties('sx', qarg, InstructionProperties(error=error_rates["sx_error"], duration=default_duration)) 
        if "measure_error" in error_rates:
            backend.target.update_instruction_properties('measure', qarg, InstructionProperties(error=error_rates["measure_error"], duration=default_duration))            
    
    # Update the 2-qubit gate properties (CX gate) for all edges in the chosen coupling map
    if "cx_error" in error_rates:    
        for edge in backend.coupling_map:
            backend.target.update_instruction_properties('cx', tuple(edge), InstructionProperties(error=error_rates["cx_error"], duration=default_duration))

ここで、この関数を使用してすべてのゲートのエラー率を \(1e-8\) に変更し、ゲートの持続時間も同じ数に変更します。

error_rates = {
  "default_duration": 1e-8,
  "rz_error": 1e-8,
  "x_error": 1e-8,
  "sx_error": 1e-8,
  "measure_error": 1e-8,
  "cx_error": 1e-8
}

update_error_rate(fake_backend, error_rates)

Exercise 4: カスタムFake backendの構築#

Exercise

カスタム fake backend を定義し、エラーレートを変えてもう一度 VQC をテストし、Exercise 3 での結果を検証しましょう。

  • fake backend や本物のbackendで作業するときは、必ず回路とオブザーバブルの両方をトランスパイルする必要があります。この例は、最近のQiskit SDK v1.0 のブログNew primitives design セクションにあります。

  • RZ ゲートと CX ゲートのエラーレートを変えて VQC をテストします。いずれの場合も update_error_rate 関数を使用する必要があります。

  • Exercise 3 での最適なパラメーターを使用してください。

  • 以前に作成した関数を使用して、各エラー・レートの合計パフォーマンス(P)を計算し、各鳥のラベルと比較した最終コストをプロットしてください。

fig, ax = plt.subplots(1, 1, figsize=(7,5))
ax.set_title('Cost on test data')
ax.set_ylabel('Cost')
ax.set_xlabel('State index')
ax.plot(list_labels, 'k-', linewidth=3, alpha=0.6, label='Labels')

error_rate_list = [1e-1, 1e-2, 1e-3, 1e-4]

fake_backend = GenericBackendV2(
    num_qubits=5,
    basis_gates=["id", "rz", "sx", "x", "cx"]
  )

for error_rate_value in error_rate_list:

    update_error_rate # Add your code here
    
    estimator = # Add your code here
    pm = # Add your code here

    opt_params = # Add your code here
    results_test = test_VQC(list_coefficients, list_labels, ansatz, obs, opt_params, estimator, pm)

    print(f"Performance for run {index}: {compute_performance(results_test, list_labels)}")
    ax.plot(results_test, 'o--', label='Predictions error rate '+str(error_rate_value))
    
ax.legend()
Hide code cell content
fig, ax = plt.subplots(1, 1, figsize=(7,5))
ax.set_title('Cost on test data')
ax.set_ylabel('Cost')
ax.set_xlabel('State index')
ax.plot(list_labels, 'k-', linewidth=3, alpha=0.6, label='Labels')

error_rate_list = [1e-1, 1e-2, 1e-3, 1e-4]

fake_backend = GenericBackendV2(
    num_qubits=5,
    basis_gates=["id", "rz", "sx", "x", "cx"]
  )

for error_rate_value in error_rate_list:

    error_rates["rz_error"] = error_rate_value
    error_rates["cx_error"] = error_rate_value
    update_error_rate(fake_backend, error_rates) # Add your code here
    
    estimator = Estimator(fake_backend) # Add your code here
    pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend) # Add your code here

    opt_params = res_list[best_result_index].x # Add your code here
    results_test = test_VQC(list_coefficients, list_labels, ansatz, obs, opt_params, estimator, pm)

    print(f"Performance for run {index}: {compute_performance(results_test, list_labels)}")
    ax.plot(results_test, 'o--', label='Predictions error rate '+str(error_rate_value))
    
ax.legend()
Performance for run 4: 84.381103515625
Performance for run 4: 86.67449951171875
Performance for run 4: 97.03216552734375
Performance for run 4: 99.530029296875
<matplotlib.legend.Legend at 0x7ff2610344d0>
_images/1d69a8ea1d6de99340c5871993059e61eda6febcb10e9c541b28470027948e41.png
解説

RZCXのエラーレートをupdate_error_rate関数を使って小さくしていくほど予測精度が改善していることがわかります。

Part 3: 量子ハードウェアでVQCをテストする#

このパートでは、VQC を実際のハードウェア上で テスト します!

そのために、いくつか変更が必要です。現状の主な問題は、これまで考えてきた回路(振幅埋め込み + ansatz)がハードウェアの制限を考慮していないことです。特に、振幅埋め込みとして利用してきた qc.initialize は任意の状態ベクトル状態を回路にマッピングできる一般的な関数ではありますが、効率の面では芳しくありません。

データセットの構造に合わせたカスタムマッピングを作成し、回路の深さを削減できるはずです。さらに、これまではすべての量子ビットのペアに対して接続を持つ RealAmplitudes ansatz を検討してきましたが、これは heavy-hex 接続を持つ実際のデバイスにトランスパイルすると多数の SWAP ゲートを必要としてしまいます。

このことを明らかにするために、まず、これまで使ってきた VQC を接続性の限られたデバイスにトランスパイルするとどうなるかをチェックしましょう。
AerSimulator を使用して、実際のバックエンドデータからシミュレータを作成し、デバイスとして使用します。

ヒント

AerSimulator.from_backend(backend) を使用すると、IBM Quantum backend を模倣するようにシミュレータが自動的に設定されます。これは、指定したバックエンドの現在のノイズモデル・基本ゲート・カップリングマップを実装することを意味します。

IBM Quantum backend を模倣したシミュレータを得るもう1つの方法は、 qiskit_ibm_runtime.fake_provider.FakeBackend() オブジェクトを使うことです。これも似たような構成のシミュレーターを自動的に作成しますが、ノイズモデルは現状サポートしていません。

Lab のこのパートでは、アクセスできる実際の IBM Quantum backend を選択することをお勧めします。後で本物のバックエンドでコードを実行するためのガイドとして、このシミュレーションを利用するからです。デバイス毎に特性が大きく変わる可能性があるため、テストはデバイスによってはあまりうまくいかない可能性があることを考慮してください。

我々は ibm_sherbrooke を使用しました。IBM Quantum platform の Systems タブ でアクセスできるシステムのリストを見つけるか、以下のコードを実行してください:

    from qiskit_ibm_runtime import QiskitRuntimeService

    service = QiskitRuntimeService()
    # list all available systems
    service.backends()
    # get the least busy system
    backend = service.least_busy()
# Choose a real backend
service = QiskitRuntimeService()
backend = service.backend("ibm_sherbrooke")

# Define a fake backend with the same properties as the real backend
fake_backend = AerSimulator.from_backend(backend)

まず、 qc.initialize を利用してトランスパイルした後の2量子ビット回路の深さを確認してみましょう。

index_bird = 0 #you can check different birds by changing the index
qc = QuantumCircuit(num_qubits)
qc.initialize(list_coefficients[index_bird])
pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend)
transpiled_qc = pm.run(qc)

print('Depth of two-qubit gates: ', transpiled_qc.depth(lambda x: len(x.qubits) == 2))
transpiled_qc.draw(output="mpl", idle_wires=False, fold=40)
Depth of two-qubit gates:  39
_images/226a7f9a51c64df551882ab7937c47ac968bb62af2c3c016db4cb17882456a78.png

この回路は実際の量子ハードウェアで動かすには深すぎます!データセットの状態をマッピングする別の方法を探してみましょう。

Exercise 5: データセットに合わせたカスタムマッピング#


状態をより効率的にマッピングする関数を作成します。そのためには、データセットをよく見て、どのように構築されたかを理解する必要があります。ハードウェア上でコードを実行したいときはいつでも、データの構造と問題を賢く利用する必要があります!

- **IBM Quantum の鳥**は、辞書の最初の5つのエントリーに対応し、インデックスは $0$、$1$、$2$、$3$、$4$ です。この状態は次のようにして作成されます:
    - まず、関数  `generate_GHZ` を使ってGHZ 状態を生成することから始めます。
    - 次に、インデックスを5桁の2進数表現に変更します。例えば、鳥番号 $2$ は2進数 $00010$ に対応します。つまり、最後から 2 番目の量子ビットにビット反転を追加する必要があります。
    - これにより最終的なGHZ状態は、鳥番号 $2$ の場合は $\vert 00010 \rangle + \vert 11101 \rangle$ となり、鳥番号 $3$ の場合は $\vert 00011 \rangle + \vert 11100 \rangle$ になります。

- **IBM Quantum でない鳥**は、辞書の後半の5つのエントリーに対応し、インデックスは$5$、$6$、$7$、$8$、$9$ です。ここで、上記の2番目のステップで示したように、状態は鳥のインデックスの2進表現に対応する積状態です。
Hide code cell content
def amplitude_embedding(num_qubits, bird_index):
    """Create amplitude embedding circuit

    Parameters: 
        num_qubits (int): Number of qubits for the ansatz
        bird_index (int): Data index of the bird

    Returns:
        qc (QuantumCircuit): Quantum circuit with amplitude embedding of the bird
    """
    def generate_GHZ(qc):
        #GHZ state
        qc.h(0)
        for i, j in zip(range(num_qubits-1), range(1,num_qubits)):
            qc.cx(i, j)

            
    ### Write your code below here ###
    qc = QuantumCircuit(num_qubits)
    bird_index_ibmq = (0, 1, 2, 3, 4)
    if bird_index in bird_index_ibmq:
        # generating a GHZ state
        generate_GHZ(qc)
        qc.barrier()
    # change the index to its binary rep with 5 digits
    bird_index_bin = format(bird_index, f'0{num_qubits}b')
    # add a bit flip                                                
    for i in range(len(bird_index_bin)):
        if bird_index_bin[i] == '1':
            qc.x(num_qubits-i-1)

        
    ### Don't change any code past this line ###
    return qc

下のセルを使って、正しい状態を作り出せているか確認できます:

index_bird = 0 # You can check different birds by changing the index

# Build the amplitude embedding 
qc = amplitude_embedding(num_qubits, index_bird)
qc.measure_all()

# Define the backend and the pass manager
aer_sim = AerSimulator()
pm = generate_preset_pass_manager(backend=aer_sim, optimization_level=3)
isa_circuit = pm.run(qc)

# Define the sampler with the number of shots
sampler = Sampler(backend=aer_sim)
result = sampler.run([isa_circuit]).result()
samp_dist = result[0].data.meas.get_counts()
plot_distribution(samp_dist, figsize=(15, 5))
_images/00b8e2d665431c35e6782ba3a9c9d4ce995aca1306515e316281a333641b8ac9.png

さて、新しい振幅埋め込み回路の深さをチェックしてみましょう:

index_bird = 0 #You can check different birds by changing the index
qc = amplitude_embedding(num_qubits, index_bird)
pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend)
transpiled_qc = pm.run(qc)

print('Depth of two-qubit gates: ', transpiled_qc.depth(lambda x: len(x.qubits) == 2))
transpiled_qc.draw(output="mpl", fold=False, idle_wires=False)
Depth of two-qubit gates:  4
_images/b3e317faa5ec0a2d6b6b959b43222ae3583d28c9d78e7261f0ed0b3360e33694.png

では、完全な接続性を使って、トランスパイルされた RealAmplitudes ansatzの深さをチェックしてみましょう。

old_ansatz = RealAmplitudes(num_qubits, reps=1, entanglement='full', insert_barriers=True)
pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend)
transpiled_ansatz = pm.run(old_ansatz)

print('Depth of two-qubit gates: ', transpiled_ansatz.depth(lambda x: len(x.qubits) == 2))
transpiled_ansatz.draw(output="mpl", idle_wires=False, fold=40)
Depth of two-qubit gates:  16
_images/281d36ad0689ccb10acc046228a3b144525365fa4cedba2abbe536eb89615836.png

接続性を pairwise に変更し、回路の深さを再度確認します。深さはどのくらい減少しましたか?

ansatz = RealAmplitudes(num_qubits, reps=1, entanglement='pairwise', insert_barriers=True) # Add your code here
pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend) # Add your code here
transpiled_ansatz = pm.run(ansatz) # Add your code here

print('Depth of two-qubit gates: ', transpiled_ansatz.depth(lambda x: len(x.qubits) == 2))
transpiled_ansatz.draw(output="mpl", fold=False, idle_wires=False)
Depth of two-qubit gates:  2
_images/d6ec05b8b491c014eb5c37565ee62e65fe819d96861f0d2b34b0da70db3c6f8c.png

ここで、振幅埋め込みと ansatz を組み合わせてトランスパイルすることで、新旧VQC回路の2量子ビットのゲート深さの合計を比較することができます。

old_mapping = QuantumCircuit(num_qubits)
old_mapping.initialize(list_coefficients[index_bird])
old_classifier = old_mapping.compose(old_ansatz)

new_mapping = amplitude_embedding(num_qubits, index_bird)
new_classifier = new_mapping.compose(ansatz)

pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend)
old_transpiled_classifier = pm.run(old_classifier)
new_transpiled_classifier = pm.run(new_classifier)

print('Old depth of two-qubit gates: ', old_transpiled_classifier.depth(lambda x: len(x.qubits) == 2))
print('Current depth of two-qubit gates: ', new_transpiled_classifier.depth(lambda x: len(x.qubits) == 2))
Old depth of two-qubit gates:  66
Current depth of two-qubit gates:  6

新しいansatzによって、深さは10分の1削減できました!これで、量子ハードウェア上でVQCをテストする準備が整いました。

Exercise 6: 結果取得の準備#

更新された回路で新しい `test_shallow_VQC` 関数を作成してください。この関数は、新しい振幅埋め込みのステップと新しい ansatz を実装する必要があります。ここでは、各鳥をインデックスで直接マッピングしているので、係数のリストはもう必要ないことに注意してください。
```python
def test_shallow_VQC(list_labels, ansatz, obs, opt_params, estimator, pm):
    
    """Return the performance of the classifier

    Parameters:
        list_labels (list): List of labels
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        obs (SparsePauliOp): Observable
        opt_params (ndarray): Array of optimized parameters
        estimator (EstimatorV2): Statevector estimator
        pm (PassManager): Pass manager for transpilation 
        
    Returns:
        results_test (list): List of test results
    """
    ### Write your code below here ###



        
    ### Don't change any code past this line ###                       
    return results_test
```
Hide code cell content
def test_shallow_VQC(list_labels, ansatz, obs, opt_params, estimator, pm):
    
    """Return the performance of the classifier
    Parameters:
        list_labels (list): List of labels
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        obs (SparsePauliOp): Observable
        opt_params (ndarray): Array of optimized parameters
        estimator (EstimatorV2): Statevector estimator
        pm (PassManager): Pass manager for transpilation 
        
    Returns:
        results_test (list): List of test results
    """
    ### Write your code below here ###
    results_test = []
    for index_bird in range(len(list_labels)):
        amp_embed = amplitude_embedding(ansatz.num_qubits, index_bird)
        new_classifier = amp_embed.compose(ansatz)
        # transpile classifier and observable
        new_transpiled_classifier = pm.run(new_classifier)
        transpiled_obs = obs.apply_layout(layout=new_transpiled_classifier.layout)
        # Run estimator
        pub = (new_transpiled_classifier, transpiled_obs, opt_params)
        job = estimator.run([pub])
        # Get result
        # result = np.abs(job.result()[0].data.evs.item())
        result = np.abs(job.result()[0].data.evs)
        results_test.append(result)


        
    ### Don't change any code past this line ###                       
    return results_test

fake backend で shallow VQC (浅いVQC) をテストしてみましょう!新しいVQCはすでにトレーニング済みなので、 opt_params_shallow_VQC.npy で与えられた最適なパラメーターを使用できます。

estimator = Estimator(backend=fake_backend)
estimator.options.default_shots = 5000
pm = generate_preset_pass_manager(optimization_level=3, backend=fake_backend)

opt_params = np.load('./data/opt_params_shallow_VQC_2024.npy') # Load optimal parameters
results_test = test_shallow_VQC(list_labels, ansatz, obs, opt_params, estimator, pm)

print(f"Performance: {compute_performance(results_test, list_labels)}")

fig, ax = plt.subplots(1, 1, figsize=(7,5))
ax.set_title('Cost on test data')
ax.set_ylabel('Cost')
ax.set_xlabel('State index')
ax.plot(list_labels, 'k-', linewidth=3, alpha=0.6, label='Labels')
ax.plot(results_test, 'o--', label='Fake Backend')
ax.legend()
/home/dmurata/workspace/quantum/ibm-quantum-challenge-textbook/.venv/lib/python3.11/site-packages/qiskit_ibm_runtime/fake_provider/local_service.py:243: UserWarning: Options {'default_shots': 5000} have no effect in local testing mode.
  warnings.warn(f"Options {options_copy} have no effect in local testing mode.")
Performance: 97.265625
<matplotlib.legend.Legend at 0x7ff20629e010>
_images/5b7ac3e9499989e015ded02616bbeabe2d2f3567e72a0b7d46d90607ca829380.png

このシミュレーションは、前に選んだ実際のデバイスで得られる結果と似たものになるはずであることに注意してください。ノイズの影響で性能が低下していることがわかりますが、それでも信頼できる結果が得られています。

実機でのテスト#

この Lab の最後のパートでは、実際の量子ハードウェア上でコードを実行します。

Fake backend を使うことで、ゲートのノイズをシミュレートすることができますが、量子ビットのデコヒーレンスは考慮されていないことに注意してください。この特殊なケースでは、デバイスのほとんどの量子ビットが長期間使用されないため、デコヒーレンスに対処する必要があります。このため、ダイナミック・デカップリング(DD)のようなエラー抑制技術を使用することが非常に重要です。さらに、resilience_level を変化させることで、異なるエラー緩和技術を実装します:

  • resilience_level = 0: エラー緩和なし

  • resilience_level = 1: 測定トワリング+トワリング読み出しエラー抑制(TREX) を実装。

  • resilience_level = 2: resilience_level_1 + ゼロノイズ外挿 (ZNE) + ゲート・トワリング を実装。

これらの各手法について詳しく説明しましょう。

Dynamical Decoupling (DD)#

上のセクションで述べたように、量子ビットはデコヒーレンスによって時間とともに情報を失い、さらにクロストークによって他の量子ビットに適用される演算の影響を受ける可能性があります。これらの影響は、量子ビットがすでに量子状態で初期化されているにもかかわらず、計算中に長期間使用されなかった場合に特に顕著に現れます。

このような効果を排除するために、 Dynamical Decoupling が利用されます。これは、アイドル状態の量子ビットの状態を反転させるパルスシーケンス( dynamical decoupling sequences として知られています)を追加することで、デコヒーレンスによるノイズの影響を打ち消し、抑制します。

下の2つの図では、アイドル状態の量子ビットにXゲートが適用されているのがわかります。偶数個のXゲート(パルスシーケンス)を適用していますが、XX=Iであるため、論理的には量子回路の計算結果に対して影響を与えることなく、物理的な量子ビットに対してノイズの影響を抑制することが出来ます。

詳細は こちら をご覧ください。

Twirled Readout Error eXtinction (TREX)#

TREXは、パウリオブザーバブルの期待値測定時における測定エラーを緩和する手法です。このチャネルは、測定直前にXゲート、直後にを介して量子ビットをランダムに反転させることで、測定時の読み出しエラーを対角化することができます。通常のエラーは異なる状態間で相互に影響し合っていますが、対角化することで異なる状態に対してエラーを独立させることができ補正が容易になります。
また、測定後に古典ビット反転を挟んで、論理的な計算結果に影響を与えないようにします。

対角化されたノイズの寄与は、ゼロ状態で初期化されたランダム回路(キャリブレーション回路)によって学習され、補正されます。このキャリブレーション回路の実行が追加で必要となるためにTREXは全体のオーバーヘッドが増加するという特性もあります。

詳細は こちら をご覧ください。

Zero Noise Extrapolation (ZNE)#

ZNEは、Estimator primitive で使用できるエラー緩和手法です。これには2つの異なるフェーズがあります:

  • 回路内のノイズ増幅:異なるノイズレベルで期待値を計算

  • 外挿:各ノイズレベルでの計算結果を使用して、ノイズがない(ゼロノイズ)場合の期待値を外挿により推定

最初のステップでノイズを増幅する方法いくつかあります。例として、下図のように回路内のゲートを複数回繰り返すゲートフォールディングがあります。

Qiskit runtime でエラー緩和を設定する方法についての詳細は、このチュートリアルドキュメント に記載されています。

Estimator primitiveを使用する場合、 EstimatorOptions() オブジェクトを使用してエラー緩和とエラー抑制技術のオプションを設定することができます。この表は最もよく使われるものを示しています:

Options

Sub-options

Sub-sub-options

Choices

Default

default_shots

4096

optimization_level

0/1

1

resilience_level

0/1/2

1

dynamical_decoupling

enable

True/False

False

sequence_type

'XX'/'XpXm'/'XY4'

'XX'

extra_slack_distribution

'middle'/'edges'

'middle'

scheduling_method

'asap'/'alap'

'alap'

resilience

measure_mitigation

True/False

True

measure_noise_learning

num_randomizations

32

shots_per_randomization

'auto'

zne_mitigation

True/False

False

zne

noise_factors

(1, 3, 5)

extrapolator

'exponential'/
'linear'/
'double_exponential'/
'polynomial_degree_(1 <= k <= 7)'

('exponential', 'linear')

twirling

enable_gates

True/False

False

enable_measure

True/False

True

num_randomizations

'auto'

shots_per_randomization

'auto'

strategy

'active'/
'active-circuit'/
'active-accum'/
'all'

'active-accum'

まず、VQC をテストする IBM Quantum backend を選択することから始めましょう:

service = QiskitRuntimeService()
backend = service.backend("ibm_sherbrooke")

Exercise 7: 実機での結果取得の準備#

IBM Quantum backend で shallow VQC (浅い VQC) をテストする関数を作成してください。この関数は `test_shallow_VQC` に似ていますが、各鳥の pub `pub = (transpiled_classifier, transpiled_obs, opt_params)` を含む pubs のリストを作成する必要があります。そして、pubs のリストを使用して Estimator primitiveを呼び出し、job ID を表示します。

関数の準備ができたら、2つの異なるケースについて `EstimatorOptions()` を定義する必要があります:

- `options_0`: ダイナミック・デカップリングもエラー緩和手法も実装しません。
- `options_1`: シーケンスタイプ `XpXm` と測定トワリング + TREX を用いてダイナミック・デカップリングを実装します。

両ケースとも、ショット数を `5000` 、 `optimization_level=0` に設定してください。

```python
def test_shallow_VQC_QPU(list_labels, ansatz, obs, opt_params, options, backend):
    
    """Return the performance of the classifier

    Parameters:
        list_labels (list): List of labels 
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        obs (SparsePauliOp): Observable
        opt_params (ndarray): Array of optimized parameters
        options (EstimatorOptions): Estimator options
        backend (service.backend): Backend to run the job 
        
    Returns:
        job_id (str): Job ID
    """

    estimator = Estimator(backend=backend, options=options)
    pm = generate_preset_pass_manager(optimization_level=3, backend=backend)

    pubs = []
    for bird, label in enumerate(list_labels):
        ### Write your code below here ###



        
        ### Don't change any code past this line ###
        pub = (transpiled_classifier, transpiled_obs, opt_params)
        pubs.append(pub)

    job = estimator.run(pubs)
    job_id = job.job_id()
    print(f"Job ID: {job_id}")
    print(f"Status: {job.status()}")
                                           
    return job_id
```
Hide code cell content
def test_shallow_VQC_QPU(list_labels, anstaz, obs, opt_params, options, backend):
    
    """Return the performance of the classifier

    Parameters:
        list_labels (list): List of labels 
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        obs (SparsePauliOp): Observable
        opt_params (ndarray): Array of optimized parameters
        options (EstimatorOptions): Estimator options
        backend (service.backend): Backend to run the job 
        
    Returns:
        job_id (str): Job ID
    """

    estimator = Estimator(backend=backend, options=options)
    pm = generate_preset_pass_manager(optimization_level=3, backend=backend)

    pubs = []
    for bird, label in enumerate(list_labels):
        ### Write your code below here ###
        amp_embed = amplitude_embedding(anstaz.num_qubits, bird)
        new_classifier = amp_embed.compose(anstaz)
        transpiled_classifier = pm.run(new_classifier)
        transpiled_obs = obs.apply_layout(layout=transpiled_classifier.layout)
       
        ### Don't change any code past this line ###
        pub = (transpiled_classifier, transpiled_obs, opt_params)
        pubs.append(pub)

    job = estimator.run(pubs)
    job_id = job.job_id()
    print(f"Job ID: {job_id}")
    print(f"Status: {job.status()}")
                                           
    return job_id
## No DD, no TREX (no ZNE)
options_0 = EstimatorOptions() #Add your code here
options_0.optimization_level = 0
options_0.default_shots = 5000
options_0.resilience_level = 0
options_0.dynamical_decoupling.enable = False
options_0.dynamical_decoupling.sequence_type = "XpXm" 
options_0.twirling.enable_measure = False

## DD + TREX (no ZNE)
options_1 = EstimatorOptions() #Add your code here
options_1.optimization_level = 0
options_1.default_shots = 5000
options_1.resilience_level = 1
options_1.dynamical_decoupling.enable = True
options_1.dynamical_decoupling.sequence_type = "XpXm" 

おめでとうございます。これで Lab は終了です。 システムでの実行はオプションで、Grader はありません。以下の事前実行結果を確認するか、自分でセルを実行して結果を確認してください。

警告

実機でコードを実行する前に、すべてが正しいことを確認してください!

各jobの実行時間は、量子デバイス上では30秒から50秒程度です。そのため、Open Plan (ユーティリティ・スケールのシステムで1ヶ月の実行時間は10分まで) でも実行できるはずです。もしキューが長すぎる場合は、job を実行したままにしておいて、後日 job_id を下記の retrieve_job 関数に渡すことで job を取得することができます。

def retrieve_job(job_id):

    ''' Retrieve results from job_id '''

    job = service.job(job_id)

    results_test = []
    errors_test = []
    for result in job.result():
        results_test.append(abs(result.data.evs))
        errors_test.append(abs(result.data.stds))

    return results_test, errors_test
## No DD, no TREX (no ZNE)
job_id_0 = test_shallow_VQC_QPU(list_labels, ansatz, obs, opt_params, options_0, backend)

## DD + TREX (no ZNE)
job_id_1 = test_shallow_VQC_QPU(list_labels, ansatz, obs, opt_params, options_1, backend)
Job ID: csstyjjzsqjg008t92hg
Status: QUEUED
Job ID: csstyktah4b0008cw6v0
Status: QUEUED

DD + TREX + ZNEのジョブを送信することもできます!これにより、さらに良い結果が得られます。ZNEはオーバーヘッドが大きいので、job の実行時間が長くなることに注意してください。

job が終了したら、以下のコードを使って結果をプロットできます:

results_test_0_DD, errors_test_0_DD = retrieve_job(job_id_0) #(Add job_id 0 here)
results_test_1_DD, errors_test_1_DD = retrieve_job(job_id_1) #(Add job_id 1 here)

fig, ax = plt.subplots(1, 1, figsize=(7,5))
ax.set_title('Cost on test data')
ax.set_ylabel('Cost')
ax.set_xlabel('State index')

print(f"Performance for no DD + no TREX: {compute_performance(results_test_0_DD, list_labels):.3f}")
print(f"Performance for DD + TREX: {compute_performance(results_test_1_DD, list_labels):.3f}")
ax.errorbar(range(10), results_test_0_DD, fmt='-o', yerr=errors_test_0_DD, color='tab:orange', label='ibm_sherbrooke no EM')
ax.errorbar(range(10), results_test_1_DD, fmt='-o',  yerr=errors_test_1_DD, color='tab:blue', label='ibm_sherbrooke TREX + DD')

ax.plot(list_labels, 'k-', label='Labels')

ax.legend() 
Performance for no DD + no TREX: 93.506
Performance for DD + TREX: 96.869
<matplotlib.legend.Legend at 0x7ff201479d50>
_images/9dd418459a8f238dc86077a87f4883bfa7f05e198b7da13a4943114917f41988.png

ノイズレベルが異なるため、使用するデバイスによって結果が異なる可能性があります。これは ibm_kyiv で得られた結果の一例です。さまざまなエラー緩和とエラー抑制のテクニックを組み合わせることで、結果が向上していることがわかります。

  • エラー緩和なしでのパフォーマンス: 86.824

  • DD + TREXを使用した場合のパフォーマンス: 92.979

  • DD + TREX + ZNEを使用した場合のパフォーマンス: 98.004

Qiskit Patterns ワークフローを使用した VQC の構築と訓練に関するこの Lab の完了、おめでとうございます!この Exercise を通して、あなたは理想的なバックエンド上で VQC を構築し、訓練する実践経験を得ただけでなく、ノイズが VQC 性能に与える重大な影響についても探求しました。

回路の深さを減らし、実際の量子ハードウェアで VQC を実行することで、現実の量子コンピューティングの課題と解決策を包括的に理解することができました。さらに、実際のバックエンドで VQC をテストすることで、量子エラーレートと計算の忠実度に関する実用的な側面について貴重な洞察を得ることができました。

この Lab と Challenge によって、量子プログラミングのスキルが強化され、量子エラーレートと計算の忠実度の複雑なバランスに対する理解が深まったことを願っています。これからも実験と学習を続け、量子技術の可能性の限界に挑戦してください。量子コンピューティングをマスターする旅は、チャレンジングであると同時にエキサイティングでもある。

さらなる学習とリソースについては、IBM Quantum Learning platform をご覧ください。

Happy coding with Qiskit!

追加情報#

Created by: Joana Fraxanet, Marcel Pfaffhauser, Junye Huang

Advised by: Vishal Sharathchandra Bajpe, Pedro Rivero, Blake Johnson, Paul Nation

Translated by: Kifumi Numata

Version: 1.2.0