covariant-kernel-learning-line.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from qulacs import QuantumState,QuantumCircuit\n",
    "from qulacs import ParametricQuantumCircuit\n",
    "from qulacs.gate import RY,CZ,merge\n",
    "from qulacs.state import inner_product\n",
    "import numpy as np \n",
    "import random\n",
    "\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "attachments": {
    "image.png": {
     "image/png": "iVBORw0KGgoAAAANSUhEUgAAAc0AAABMCAYAAAAP3uzUAAAAAXNSR0IArs4c6QAAAARnQU1BAACxjwv8YQUAAAAJcEhZcwAADsMAAA7DAcdvqGQAABrGSURBVHhe7d0HfBTFHgfw35X0hJBKek8IMQiCQYr0UKSJ/YmI2B6ggOUpDwuiD3vF91BUbIgIShMpihp6NVJDCCSEFEJ6J+WSXHmzdxNzIuDldvdu726+fu5zM5Ng9vZ257+zO0WWkJCgA8MwDMMwf0tO3xmGYRiG+RssaDIMwzCMiVjQZBiGYRgTsaDJMAzDMCZiQZNhGIZhTMSCJsMwDMOYiAVNhmEYhjERC5oMwzAMYyIWNBmGYRjGRCxoMgzDMIyJbDJo6nRs5j/GsbFzwHLYvrY+mUxGU9Znk0FTSjuQYazB+Bxg54O42P61PilduLDbswxj41hLyLrY/ncskl/lREb+C/BJhp93PLzcg+Hm4gOF3BkabRtaWutwsakE1fXnUFF7EmpNC/1XTGdwJ/3lrqZ1Tj7QeiZB5xYFOAcCSi/yhZDrLI0KaKsBVEVQNGVD1phD/wUjlNDgSMREJSC4Wzi6dvWFm6u7/ntqampAVU0liksKcDbvNKqqy+m/YITA1Te+XXrCxysOHq7BcHXuqG9a2+rQqCpFbUMuqupPkjJW3wjJ29sb3RN6ICIiEoEBgfD09IJcLoeqRYWamhoUF19AXl4uzuZat76RZNB0dvJEQsQExASPQETQjeSgVdCfXN35skM4V7IdOYVbyMFdQUuZztC6x0MXeBPgOxxyryRaenW61lpoa/YAlb9CWfkTLWU6q1+fwRg8YBRSyHuAXxAtvbr882dx6PAu7Nn/M7KyT9BSpjOclZ6ICh6PiICRCPEfROobJf3J1ZVUp6Ow7Ffkl25FUwu7eDFHbEwcRo4YhcE3DkViomn1TV1dHQ4c3Iudu3YgbfvPtNRyJBU0vT3C0TthOq6NmwI516Lh4VTeBhzLWU5aoFm0hLkaje8Q0ry5F3LunQddSwW0xd9AWfQ5icCkRcr8rcnjp2LyuCkID42mJebJyDqC7zevwM597MLlUpe7m+LlFoakyOlIjLyHd31z9sL3yMz/EtUX/1zfcH+T3b79q/43DMQ/7roHgwYOpiXmqayqxNp132LlN19BpWqmpeKSTNAc2PMJpPSYSXPCOZ6zEvsz3kGrupGWMMa07rHQRj8JhX8qLRGGrrUauoL/QlG8ipYwlxo2aCzuv+cxhIVE0RJhHD+Zjs9XvoeTWUdpCXOpPvGPo1fsLJoTTlbBShzJeZfUNw20hDEWHRWDWTNmY/hwYesb7vbtx8s+wLr139ES8Vg9aAb5XocRfV9CgE93WiK8i81l2PH7i8gr2U5LGI46ZCoU8QtoThyaqu1QZC+ErJXdvmqnVDrh8VkLcdPI22iJOJavXoKvVn9AcwwnsGsf9E9aCL8uibREeA2qMhzMegnny9JoCcO54/Z/YN5Tz9KcOHbv2YnX31iEikrxHs8p/Pz8XqRpi0uKvg2TBi+Fh5s/LRGHi5MnukdOgFqtQknVEVrq2NTxL0ER9SjNiUfuHg1twDhSk5yCrOUCLXVckWGxeHXBx+jfdygtEU/v5H6IjkjA/t+2Q6vV0FLHFR92G1L7LIW7i7j1DfeMNCZ4PNQaFcprWX3DmT/veTz8oPB3Ei8VGRmFUaljceZMFkpKimmpsKwWNLlnlyP6WvZPRwRxD/ldcL78AC1xTOqk96EImkxz4pORSkTW7RZoSeCUN+fRUsfTPS4Zry9cJvjt2KuJDI/FtSR47j+UhtZWx+3tyT27HHiNZeubUP9BUMpdUVy1n5Y4pjdeewfjx02iOfF5eHjq/95pEjgLCwtoqXCsEjR7xk7B8D7P05xlhQb0hQxyFFUcoiWORZ30HhQBN9GcBclkkAVOgLbuMOSqIlroOKLC4/DGi8vg6xNASyynW0AIeiT0wi87Nzpkp5Tu4XdjQJK4jyGupJsPqW9kCpRWd9Q3jjRZwiuL3kDqyDE0ZzncPh4zehyOHT+iH6oiJIsHzajg4bip/9s0Zx1hgf3Q0FSGitpTtMQxqGPmQxF8J81Zie8Q6Cp/hkxdTwvsn7OzC1594ROEBkfQEssLCgzVj/3cc+AXWuIYwvyHYVjvd2jOOoJ8U9CkKkdVfSYtcQyPzf0Xbpl8O81Zx8ABg7FjVxrq64Wrb/j1s+4kd1c/pKa8QnPWlZryMvy9xesMIDXqwElQhN9Pc9Yjc/aFNmERzTmGx2cuRFy09Y+1kUMm4K5bHqQ5++fm7IdBPaVR3wxKXgS/LqaNQ7QHY8eMx9Qp99Gc9fj4+OCZfwt7l8GiQXNwr2fgQQKnVAzpLW5PLqnQKbtCHvsczVmfwmcA1GEP05x9Gzl0AsaMuIXmrO+f9z2F2CjxeqpLSUrifNE7/XRGvx7SOQfF1KWLN558fB7NWV+/lP6Ydu8DNMefxYJmVNAwJEZOpDlpCO92A5Jj7qI5+6WJmktaeF1pThoUMU9C52LarDe2xviZ1f1THqMp6Zh+j/S2SWih/kMRG2K5ziemCPK5Ht0j7qY5+zXj4Uf0LTwpeWTmHHTrJkx9Y7Gg2TdRmi2L63vMoCn7xE1eoAi9h+YkRCaHJuwhmrEv7Z1t7pz8AIK7henTUjIwZThSrhtEc/apZ8w/aUparo227/omKjIad94hvQsDhUKBaVOFeTxlkaAZ0e1GhAVeT3PS4u0RiqRo6z6sFpM2eApNSY8i7F7onKRzu15ok8dL8GKFmjROutvGV4jfIAT7SrO+8XQLRnzYHTRnf26/Tbp37rhg7uvrS3Pms0jQ7BHF/5mOjHzW0TOBh14iTe1F5EXe758NJIXQX+AhKepWmrIvOpkTFCHCnKAvPA9kf2J4HXqSFgpA081y40UtaejAMfqhHnxoqw5i9eLpmD6tB0ZO7oFR0ybgP19tRq4AM0JyrU1Ljhe1pLhQ/ufz4AeAOa8D89+nr7eBB0kjMVGAEUPxYfZ5ka5UKnHzJGGf3xdsnIWU/tciZcBs/FTOf6jO+Jv437IXPWgqFa7oHjGe5sw38SEgljRKzv0CfLUAWLsRuOhFKqf7gGSeFw/c2E0fT/urQDQBY8g37EJz5nEhlcQnbwBTI8gBfI7sc1ouFFnAWJqyL0MG8Rubpq3/GS8/ej++bR6NBR9kIe37LGxb8iaSip7FI/PfwmkB5sIf1L9j/k97GTvITV4SHTyO5swz+l9k3/QCWrKBb0mwfP0xYNsBwDsJGPcgwPeGe7euveHtwW9yfikaPiwVrq5uNMeftnwz3lr8O+64VbhAPHLkKJoyn+hzz8aEpGLijfzmv/QjMfeu/sCFNGDjTlpIOPUhlTnZn+rTwIqVtNBMu46+gmM5X+nTMtdwEmxc9Wmb01YLXZth3kV14htQ8GjJefcElpNKIonsioPkYmVeHbDlTsP+vuFd+ktCODgMspYSKBRK/TRYtqqk5AKamw0rLWxefVi/Bqa5crfchTlf98TLXzyPPkaHoqboa8yc/QrC5+7BCyP49Qw9kfk7nnjuXkNG2RUyZ8tPvCAInRo6OtNUeOBIpPb5UJ82l5xcII4fQuqErwHj0X2jnwb6kIh5hpRvSKeFZvrt9Gv6VVE4WlfyP5ULF2wsitQ3clrfLFzwMiaMF6rzVT22LRiGd+VvYfWtRzB6ZgEWbfwAYwP5h6uJk8egtLSE5jpP9KA5uNd89OnO7wFs6lwgnrQqd78CXDo8ePyTQISSXAm+SVqhtMwcOee3YesB8ocIl7jn4BI+TZ+2NU2nnoa67Ad9WnvDDnIBYP4twqgU4MO7gSOrgOdJJRFEGiZiBE1N1r+gJFeVnE3fb0NQULA+bWsGDu6LtrY2XJPYG/99nf/qLrrWFhLILr1TcBT/mzwFO29dhXXTetMy8429oxfZ5lbIu/SBZ1/bXJGmrTINzRmP6NMp3f+N5GjhhhcYGzgbGBJPLmg2AGuMLt7NkV/6M3Ycm6NPq2OfgSJsuj5tazRZ88i5u1Gf3rj+R4SEhOrTfDUceBljnyvFf1Z/gMElb2KggEHz+RfmY9vPW2mu80S/PRvocw1Nmc+LBEzdxb8GTE69ikR+8nO+1WygbzJNkQOh0borg/OhaTijf9c5+/MKmJx8EijHkYsSLmCKSefRg6bIhU9eLk3Zlrz8c/qAyYmJEmYig78GTLKvGstRTd79Pfjddm/XPumCttFw3NgizcWObffz7jiPhdaFNgZrBJjYx5+7jdPOhusbeWO2/p3rYCNUwNS1HsNHb36LwY8/ixECBMlLde/O7/wUPWj6dImhKfO5X+VOKb0bBk+eF91cL1qF3Fmf1jYV6t9tka7ZsO06V9u5zSnz6DhGhJ4n0lKMtzs8VLx9X3t4M3YjDql9Oy40+OCm1tPTNELbUmVI2xid0eo5XUR4VsjdruU6BiWHAeUHgV8FWHWK60XL9ffgyJvP699tkUxlmBA9PEy4KSIzv1mA9S5z8dB4AXp5XkZkBL/zU9SgKZcp4OFqmRk5lOY/PvqDp1ugIUHv0dsaHbfQttZwFaF1oZ/FBsicu9GUYSV2W1RV1RFw/P06Po+QdKpd+OLDX5Bw6yLcLFBc9vPtOE50bbYZNLWthvNVJpPDg+u5JpCJzxp6zs7jOgYlAXm/AGsFvIPt3n6O2uhaszpNE2TkxfH3F6a+0RSswtsfKfDwgvsQq6CFAgsMMDrmzVjAQNSg6aT0oCnxqQ3fHS9OSk/9u05jo6uu0wNYT+lFEzZA0bGtTU0CjKewAuPtdncTft/rVGex7t05yE79BK9N7Q3DPRH+uGWU/kBamzaJu1gkuHUshbTpVUPP2fbes6HDgIfmk5YK/Tlf7fUNbLW+Map0vTz573ud5gI2LH4ZLXc/hynXCHWE/5WnZ8f5aU6vcXFbmnJhLhWartK93o0+Z2g4Znjng1vCR89GF+zV6Tq2+4/PYguMtlWr1dKUbdHqOrabG68mJJ3qKJY9cxvSwpbhnQeGoKuAZy3XY7md8fFjU/7Y9+JVZ0fXAT9lksAcDAwRbASEYd/LjI4dW2J8vMgV/Oubsm0v4/3jU/HkP1MgzBP7y+NmB+JD1KCp0bTSFD8XL5IDi1wcXK5LkQ+5wOE6CZnfgbiDVke3VyHmVyYimRNNEFobWnDYaFudncS7whST8XYLudiztnYXFj91N45etxzvT7sBHgKfsX/aVplt7ntuEg+ORitMfXMlBXlkf5F3H2H6u5ALRNoasNH6RibvqG/4H/NV2PvLbqiaVuGR4dcaJjSgr4EzVpAveTcWTOopyCQHLTy3VdSg2apuhFpj6FHIx9Fc8gW5AnHDaAHFjdP0I8G04QK/4SbtmlU1+neZUlqTm5tK7uxNU+QzqGtpSvp0asN+53hxXaVtkPF219Vz/Vv505SswQtzZyC732oSMIW7JWusrr7jOJE7dRw/tkROFyNQa5rIi3/gHDwDGHWZK/TIaHJxRN5rBOqr1tJm2Pc6hW3ud5nR8VJXV0dT5vLD7e9lIv3gib+89n98L/ljQ7DohwykH1jCe9hJbS2/ulHUoMmpb+R/hFVtAYpIazJkEDDiRoC7ex5IAubNYwDXJhJUfzT8Hh/cydbUYuiEInO10dU3ZEpyINOOV80CndmWoCqiCfK9Btrmvg8I6Oj8U1ZRTFPma877DHOfeAFVqeIFTE6F0bbKXcTpwCQ24+1uUPE77rkL8Z7xQB9ST48fRb5XWn7dbcBYEkh19cDJDbSQB41WTeobQwcgnattjkvmHqto6YQYJSX8j3lL4butogfNqjrDOB6+Nn0K5FYBMeRAnrYIuP1mcnVPAumu5eQgFuDCvqquY6yUws12p7iSe8Tq32VNZ/XvfH21uGPO2d13kn1OynwSO8rOkJ/znkmzsWNsZkSEUN0sLMt4uwsK+e/7hnMH9VPlZa+/GzdNNsw9+6fXLTPxqwDHff55w76Xc2N66RAIWyNz69j3tRf5jXlsOwKs/obsF9IYSSAX5Q/SuWfHDDAMp9z4X256Cf5qGzq2U+dmu7Ngwd1Q33DjlG1F7jl+56foMwJxS4LdeO1TNCddx3NWYufR/+jTzlFz4ORLmrQ2qOX8l1BXGJremuvWQN7lWn1ayrTHp0NRe0CfXvrBZ3B2ts1na0/Neww1NdUICYrAio+20VLpqq6txB3TB+vTCu8UuMZK/zy9HE3DaaiyF+rTydEPIaX70/q0lJ0uXIUDp17Up9WRj0LmO0SftjW6ouVQVhhm1/n806/RM1n69c0jcx5Gevohmus80YNmoE9P3D1qLc1J15Z9c3H2gvQrus5QRz0BReRMmpMobStke3uTs89Ge25ewRdLtiAijP/EHmJK27UJr74nnRX2heDX5RpMGrie5qRrx/HHkV8iwHMlCZk1YzYeuF+a65i2a21txZDhN0CjMb++Ef32bHlNBmouSnuGHe75Ql7JDpqzH/LqXTQlXZqqX+0uYHL2p2+nKek6+DvPCVQlqKo+E/USn9GLG55UVG5/9c3efbtpSrp279nBK2ByRA+anOzCTTQlTWcKfhC9u7o1yOuPQNtwmuakSVa2habsy8495k8IbQkNjfXYue8nmrMvucWGBQukKvfCRqg1AqztJjEZJ0/gzJksmpOmbT/zb91bJGieypP27dmM3O9oyv7oSqT72bSNuVBwLU070j7DSM65LKQf3atPS9Hmbd/a7EQSf+dskbTrm+yiNTRlfzZsXEdT0sN1Vtq5i/8dIIsEzfqmYn1HGynKvZCG0moh+sNJk7J4JXQqaXYH1xUZ1hO0J8ZzWa79QZqfj7s9tX7z1zRnfxpUJThVIM3PV1i+HeW1h2nO/qxb/51kF134ZtUKmuLHIkGTk561FG1q4WZKEcpvp/gtWGsLtAX8FgEXg7buOJSl9tvC5/x+dB92H5Be57IV332IqmrbnCTcVCdyl0KtkV59c+zsUpqyX5998TFNScfJzAx8L1Ar2GJBs1FVgX0Zb9GcNKRnfYLympM0Z7+UpWuhqZJYp4+8d2jCvn22YjFp2alpzvq428YrvrX/C8Xm1kqkn3mT5qQh49wyVNWfoDn79cOm77Fnr7Q6IS75cDFN8WexoMk5nrMC2eel0fngQsVh7M9wjIqbIz+7CLq2epqzLk3+EijqzB8nZUuKivOx+KOXaM76lixbRFP273Th1zgnkWEdZbXH8Hv22zRn/9565zXU1/OdWk8Yyz5disOHhVtJX+Hn52cYYWshhaV7EBU8HB6ufrTE8hpV1di8bxZUrdL4Ui1Bpq6HtikXssAJtMQ6NOVboCQB3JHknDsFD3dPJHXnuVI6T+8tXYi9B9NozjFcqNyNsMDhcHexzLq+l9PcWoO0wzP/mGvWETQ0XER+fh7GjL6JlljHL7/+hDfffo3mhGHxoKnRtqGo7ACiQ1Ph4mT5ybm556qb9s5ERa2ha7Q566nZKnlzHjQtZZD7j6AllqWt2Qdl5iM051i455tBgaGIi+5BSyzr85XvY90Py2nOcWh1bSip2o/wwJFWqW+4Oa3TjsxC9cVMWuI4CgrzUVZagqFDhtMSyzp4aD+e/vfjNCcciwdNjqq1Vt/iDAsYAHdXH1oqvuaWWmwiLcziSkNT3ZxVu22dvOEUCZylJHCOpCWWoanaAcXJWZDBPoc5mGLfoTQE+HVDfOzlFrkTz7IV72LV2k9ozvG0tNWRFuceBPkNgJuzLy0Vn4r8XS5gltX8RkscB1e3cg2SM9mnUVZm+cC5e89O/bSWYgyrskrQ5DS31CCn6Ef4domHj5f4ExaXVJ7A5n0z/9Txx5Famca4wKm9SPZD14GQKd1pqXg0RV9DeWaeQwfMdgfSDTPB9E7up38XEze05K3/PYuNW7+hJY7F+PxuaatBfslWeHvGwdtD/PqmvC4D2w/PRGV9hj7vaHWN8eflAmdmZgb69esPdzfx65s1a1dj4UvPijYO2WpBk6PWNONM4SZyVaJFeGB/Wiq8I2e+xI8H5+qfLTAG8uZ8cmZvhsY1GHKPOFoqLF1rLbQ5C6E877itnMs5fjId2WdPontcMrp4ibN269GMQ3jhtUdx+Nh+WsKotSrklWzW1zfBfuLVN5n5y7Hz2BzS0mT1TbvzRYXYtm0rgoKCERNtWBlFaNyanq+/uQhfLv+MlohD9AnbTeXvnYB+SbMRHz6GlvBXULpPPw6zuPJ3WsJcjjpgPGQRMyH3TKAl/GmKvoKiYAlkasfpbNVZcrkcU++chXtunwGlsmMVfD7KK0uwcs3H+hl/mCvz8YxHr7jZiA4aS0v4K6rcjxO5H6KsRriemvYodeQYPHj/w4iLE66+Wf3tSnz6+UcCLIb99yQTNNsF+fZGcuydSIycDIVcQUs750zhj8jMW4PzZftoCWMKTeBE6IJuh8LHvKtwXUsltGUbIC9ZDbnRwtLM1XGtzYlj/4HRw29GWIh5tw6zsk/gp7T1LFh2kr93LySG34WYkJtJfaOkpZ3DDWvJKVqD4ipW33TGmNHjMGniZPRLMa++qayqxNatP2DdBrLvLTgLkSSCJnf/+9JOOUqFC6KDRyC8W38E+vSCf9eEywZR7l9V1+WirPo4CssPoqBkp0MNJeGj/WH9pbSu4dD6DgW6Xg+ZZ0/I3cLoT/5Mp2k0TAhffxSy2gNQVEt3rlVbkdzjOvTtPYi890VcdOIVb99WVJXqJyo4kZmO9CN7kC/AwteOTCF3RnjgCAT79id1TW/4eiZAfoX6prYhF5W1x1FSfQhFFTv0HY0Y84WGhmFA/0Ho2+d69Ei8Rp+/nMbGBuSczcGJjGP47beDOPSbYQ1eS5NcS/NqPFwD4ersTQ5mZ1Lhq8nBWo+GplJyIF/5I1wpMDCm08ldoHMOIDVLF3LEyElUbYasrYa8qulvMObijk2uw8KVjlEvT2909faFm6ub/lhuam5CdU0FmlVN9DcYsbi7BMLFyVDfaLVqtKnr0ai6en3D8Ofi7AI/f394enrpH2G0qFSoratBTY00nhHbVNBkGHvELuwYxnZYdBq9zrr0lu2lWEUjHey7MB/bd4yj+7u6XkokHTT/rjKxpR1t79h3wTCMuWzpwlHSQZNhGIZhpIQFTYZhGIYxEQuaDMMwDGMiFjQZhmEYxkQsaDIMwzCMiVjQZBiGYRgTsaDJMAzDMCZiQZNhGIZhTMSCJsMwDMOYiAVNhmEYhjERC5oMwzAMYyIWNBmGYRjGRCxoMgzDMIyJWNBkGIZhGJMA/weMsJTxYeMMsQAAAABJRU5ErkJggg=="
    }
   },
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Prepare training/test data\n",
    "\n",
    "![image.png](attachment:image.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "copied from Noro-san's code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_qubits = 5\n",
    "\n",
    "# 用いるスタビライザー群の要素を決定する\n",
    "# 生成元 XZIII, ZXZII, IZXZI, IIZXZ, IIIZX \n",
    "S = np.zeros((n_qubits, 2*n_qubits))\n",
    "S[0] = np.array([1,0, 0,1, 0,0, 0,0, 0,0]) # XZIII (X->1,0 Z->0,1 I->0,0としてならべる)\n",
    "S[1] = np.array([0,1, 1,0, 0,1, 0,0, 0,0]) # ZXZII\n",
    "S[2] = np.array([0,0, 0,1, 1,0, 0,1, 0,0]) # IZXZI\n",
    "S[3] = np.array([0,0, 0,0, 0,1, 1,0, 0,1]) # IIZXZ\n",
    "S[4] = np.array([0,0, 0,0, 0,0, 0,1, 1,0]) # IIIZX\n",
    "\n",
    "def x_generator(c, select_bin_i):\n",
    "    SS = np.copy(S)\n",
    "    \n",
    "    #select = np.random.randint(0,2,n_qubits) # i番目の生成元を使う(1)か使わない(0)かを決める行列\n",
    "    se = np.zeros(n_qubits)\n",
    "    for j in range(n_qubits):\n",
    "        se[j] = int(select_bin_i[j])\n",
    "        \n",
    "    for i in range(n_qubits):\n",
    "        SS[i] = SS[i] * se[i]\n",
    "    op = np.sum(SS, axis=0) # それぞれのビットにそれぞれのゲート(X, Z)が何回掛かるかを示す行列\n",
    "\n",
    "    # \\theta_{+-}を計算する\n",
    "    theta = np.zeros(2*n_qubits)\n",
    "\n",
    "    for i in range(n_qubits):\n",
    "        theta[2*i] = c[2*i] + np.pi/2 * op[2*i]\n",
    "        theta[2*i+1] = (-1)**op[2*i] * c[2*i+1] + np.pi/2 * op[2*i+1]\n",
    "\n",
    "    for i in range(2*n_qubits): # 0～\\piの値にする\n",
    "        theta[i] = theta[i] % np.pi\n",
    "    \n",
    "    return theta"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [],
   "source": [
    "random.seed(0)\n",
    "n_data_train = 20\n",
    "epsilon = 0.0 # エラーの分散\n",
    "\n",
    "l = [i for i in range(2**n_qubits)] # 0 = 2**n_qubits-1 のリスト\n",
    "select = random.sample(l, n_data_train)\n",
    "select_bin = [format(select[i], f'0{n_qubits}b') for i in range(n_data_train)]\n",
    "\n",
    "# c_{+-}を決める\n",
    "c_plus = np.pi/2 * np.random.rand(2*n_qubits) - np.pi/4 # - \\pi/4 ～ \\pi/4 のランダムな数\n",
    "c_minus = np.pi/2 * np.random.rand(2*n_qubits) - np.pi/4\n",
    "\n",
    "x_plus = []\n",
    "x_minus = []\n",
    "#x_plus = x_generator(c_plus, select_bin[0]).reshape((1,-1))\n",
    "#x_minus = x_generator(c_minus, select_bin[0]).reshape((1,-1))\n",
    "for i in range(n_data_train):\n",
    "    xx = x_generator(c_plus, select_bin[i]).reshape((1,-1))\n",
    "    e = np.random.normal(0, np.sqrt(epsilon), (1, 2*n_qubits))\n",
    "    x_plus.append(xx + e)\n",
    "    xx = x_generator(c_minus, select_bin[i]).reshape((1,-1))\n",
    "    e = np.random.normal(0, np.sqrt(epsilon), (1, 2*n_qubits))\n",
    "    x_minus.append(xx + e)\n",
    "x_plus = np.concatenate(x_plus)\n",
    "x_minus = np.concatenate(x_minus)\n",
    "x_train = np.concatenate([x_plus,x_minus])\n",
    "y_plus = np.ones(n_data_train)\n",
    "y_minus = -np.ones(n_data_train)\n",
    "y_train = np.concatenate([y_plus, y_minus]).reshape(1,-1)\n",
    "#print(x_train)\n",
    "#print(y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_train = y_train.T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(40, 10)"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x_train.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create helper function to get fiducial states, kernel matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [],
   "source": [
    "def create_featured_gate(n_qubits,Lambda,edges,x):\n",
    "    circuit = ParametricQuantumCircuit(n_qubits)\n",
    "    for i in range(n_qubits):\n",
    "        circuit.add_parametric_RY_gate(i,-Lambda[i])\n",
    "    for control, target  in edges:\n",
    "        circuit.add_CZ_gate(control,target)\n",
    "\n",
    "    for i in range(n_qubits):\n",
    "        circuit.add_RZ_gate(i,x[2*i])\n",
    "        circuit.add_RX_gate(i,x[2*i+1])\n",
    "\n",
    "    return circuit\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [],
   "source": [
    "edges= [[0,1],[1,2],[2,3],[3,4]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_kernel_value(x1,x2,Lambda):\n",
    "\n",
    "    state1 = QuantumState(n_qubits)\n",
    "    state2 = QuantumState(n_qubits)\n",
    "    state1.set_zero_state()\n",
    "    state1.set_zero_state()\n",
    "\n",
    "    circuit1 = create_featured_gate(n_qubits,Lambda,edges,x1)\n",
    "    circuit2 = create_featured_gate(n_qubits,Lambda,edges,x2)\n",
    "    \n",
    "    circuit1.update_quantum_state(state1)\n",
    "\n",
    "    circuit2.update_quantum_state(state2)\n",
    "\n",
    "    return abs(inner_product(state1,state2))**2\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "def construct_kernel_matrix(Lambda,x):\n",
    "\n",
    "    kernel_matrix = np.identity(x.shape[0])\n",
    "    for i in range(x.shape[0]):\n",
    "        for j in range(i):\n",
    "            kernel_matrix[i,j] = get_kernel_value(x[i],x[j],Lambda)\n",
    "            kernel_matrix[j,i] = kernel_matrix[i,j].conjugate()\n",
    "    \n",
    "    return kernel_matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define helper function for cvxopt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [],
   "source": [
    "from cvxopt import matrix, solvers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cvxopt_solver(Kernel,y,C):\n",
    "    \n",
    "    P = matrix(np.outer(y,y) * Kernel )\n",
    "    q = matrix(- np.ones(y.shape[0]))\n",
    "    n = Kernel.shape[0] \n",
    "\n",
    "    G = matrix(np.vstack(( - np.eye(n), np.eye(n)  ))) \n",
    "    h = matrix(np.vstack((np.zeros((n, 1)), np.full((n, 1),C) )))\n",
    "    A = matrix(y, (1,n))\n",
    "    b = matrix(np.zeros(1), (1, 1))\n",
    "\n",
    "    solvers.options[\"maxiters\"] = 10000\n",
    "    solvers.options[\"show_progress\"] = False\n",
    "\n",
    "    sol = solvers.qp(P,q,G,h,A,b)\n",
    "    return sol "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Kernel alignment programme"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [],
   "source": [
    "C = 1 \n",
    "P = 20 # the number of steps for SPSA\n",
    "### SPSA parameters. \n",
    "c = 0.1\n",
    "gamma = 0.101\n",
    "a = 0.05\n",
    "A = 1\n",
    "sigma = 0.602\n",
    "\n",
    "lam0 = .1\n",
    "lam_history = []\n",
    "F_history = []\n",
    "K_history = []\n",
    "alpha_history = []\n",
    "\n",
    "\n",
    "np.random.seed(1111)\n",
    "\n",
    "lam = np.full(n_qubits, lam0)\n",
    "for i in range(P):\n",
    "\n",
    "    c_i =  c / (i + 1)** gamma\n",
    "    a_i = a / ( i + 1 + A)** sigma\n",
    "\n",
    "    Delta = 2*np.random.randint(0, 2, 1) - 1 \n",
    "    lam_plus = lam + c_i * Delta\n",
    "    lam_minus = lam - c_i * Delta\n",
    "\n",
    "    K = construct_kernel_matrix(lam, x_train)\n",
    "    K_plus = construct_kernel_matrix(lam_plus, x_train)\n",
    "    K_minus = construct_kernel_matrix(lam_minus, x_train)\n",
    "\n",
    "    sol = cvxopt_solver(K,y_train,C)\n",
    "    sol_plus = cvxopt_solver(K_plus,y_train,C)\n",
    "    sol_minus = cvxopt_solver(K_minus,y_train,C)\n",
    "    alpha_plus = np.array(sol_plus['x']).reshape(-1)\n",
    "    alpha_minus = np.array(sol_minus['x']).reshape(-1)\n",
    "    alpha = sol['x']\n",
    "\n",
    "    F_plus = - sol_plus['primal objective']\n",
    "    F_minus = - sol_minus['primal objective']\n",
    "    F = - sol['primal objective']\n",
    "    lam = lam - (F_plus - F_minus)* Delta* a_i / (2 *c_i)\n",
    "    lam_history.append(lam[0])\n",
    "    F_history.append(F)\n",
    "    K_history.append(K)\n",
    "    alpha_history.append(alpha)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.2353992352355068\n",
      "[ 9.91e-01]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 6.81e-01]\n",
      "[ 1.00e+00]\n",
      "[ 3.55e-01]\n",
      "[ 1.00e+00]\n",
      "[ 6.78e-01]\n",
      "[ 9.54e-01]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 9.45e-01]\n",
      "[ 1.00e+00]\n",
      "[ 7.88e-01]\n",
      "[ 1.00e+00]\n",
      "[ 9.32e-01]\n",
      "[ 6.43e-01]\n",
      "[ 9.04e-01]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 6.64e-01]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 7.02e-01]\n",
      "[ 7.52e-01]\n",
      "[ 6.76e-01]\n",
      "[ 6.75e-01]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "[ 9.10e-01]\n",
      "[ 9.69e-01]\n",
      "[ 1.00e+00]\n",
      "[ 9.20e-01]\n",
      "[ 1.00e+00]\n",
      "[ 6.03e-01]\n",
      "[ 1.00e+00]\n",
      "[ 1.00e+00]\n",
      "\n",
      "0.9911110833964032\n"
     ]
    }
   ],
   "source": [
    "print(lam_history[-1])\n",
    "print(alpha_history[-1])\n",
    "a = alpha_history[-1]\n",
    "print(alpha[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAefUlEQVR4nO3deZhU9Z3v8fcXel/otWjopgEbAUUU1BaXUUOuTqLGNatm1Ynj6BNzk8zNHc3NxEwmMzfX5EniZJNxvMYsc8VJ4ijxEp0kdxSTgAIKCihby9I00Cu97/29f1Q1Nr3Q1VDd1X3q83qeeuqc8/tV1dfD8cPhd36nytwdERGZ+qbFuwAREYkNBbqISEAo0EVEAkKBLiISEAp0EZGASIrXBxcWFvr8+fPj9fEiIlPS5s2ba909NFxb3AJ9/vz5bNq0KV4fLyIyJZnZ/pHaRh1yMbPHzKzazLadpM9KM9tiZtvN7MVTLVRERE5dNGPojwPXjNRoZrnAj4Ab3f0c4EMxqUxERMZk1EB393VA/Um6fBR4yt0PRPpXx6g2EREZg1jMclkE5JnZC2a22cw+OVJHM7vLzDaZ2aaampoYfLSIiPSLRaAnARcC7wPeC3zFzBYN19HdH3H3cncvD4WGvUgrIiKnKBazXCqBWndvBVrNbB2wDNgVg/cWEZEoxeIM/RngCjNLMrMM4GLgzRi8r4iIjMGoZ+hm9gSwEig0s0rgq0AygLuvcvc3zew54HWgD3jU3Uec4tivoqaVj/zz+hO2XX/ebD5x6Xzau3q5/cevDHnNBy+cw4fKS6lv7eKen28e0v7xS+Zxw7Jiqo6184Untwxp/8sryrh6SRF7a1r4H0+9MaT9s/9lIZcvLGR7VSN//+sdQ9r/5prFXDgvn8376/nmczuHtD9wwxLOKc7hD7tr+f7/2z2k/X++/1wWhLL43Y6j/MtLFUPav/uR5RTnpvPrrVX8fMPQqaYPf/xC8jNT+MWmg/xyc+WQ9sfvWEF6ynR+tn4fz75+eEj7k391KQCPrNvL79888dp1WvJ0fvIXKwD43u9388c9tSe052WksOoTFwLw4HNv8er+hhPaZ+ek8dCt5wPwtV9vZ0dV0wntZaFMvvH+8wD40lOvU1HTekL7kuIZfPWGcwD4/OrXONzYcUL7BfPyuO+aswC4+2ebaWjrOqH9z84s5L9etRCATz32Ch3dvSe0X3X2TO66cgHAkOMOdOzp2Ju6x95Aowa6u98WRZ9vAd8arZ+IiIwfi9cPXJSXl7vuFBURGRsz2+zu5cO16cu5REQCQoEuIhIQCnQRkYBQoIuIBIQCXUQkIBToIiIBoUAXEQkIBbqISEAo0EVEAkKBLiISEAp0EZGAUKCLiASEAl1EJCAU6CIiAaFAFxEJCAW6iEhAKNBFRAJCgS4iEhAKdBGRgFCgi4gEhAJdRCQgFOgiIgGhQBcRCQgFuohIQIwa6Gb2mJlVm9m2UfpdZGa9ZvbB2JUnIiLRiuYM/XHgmpN1MLPpwIPA8zGoSURETsGoge7u64D6Ubp9FvgVUB2LokREZOxOewzdzEqAW4BVUfS9y8w2mdmmmpqa0/1oEREZIBYXRR8C7nP33tE6uvsj7l7u7uWhUCgGHy0iIv2SYvAe5cBqMwMoBK4zsx53fzoG7y0iIlE67UB39zP6l83sceBZhbmIyMQbNdDN7AlgJVBoZpXAV4FkAHcfddxcREQmxqiB7u63Rftm7n77aVUjIiKnTHeKiogEhAJdRCQgFOgiIgGhQBcRCQgFuohIQCjQRUQCQoEuIhIQCnQRkYBQoIuIBIQCXUQkIBToIiIBoUAXEQkIBbqISEAo0EVEAkKBLiISEAp0EZGAUKCLiASEAl1EJCAU6CIiAaFAFxEJiFF/JFpERMams6eX+tYu6lq6qG8NP+pau6hv7aS+tYsrFoa47tzZMf9cBbqIyCg6ut8J6NrWzkhQd4ZDumVgYIcfLZ09w77P9GlGXkYKpfkZ41KnAl1EEk5fn3OsvZu6lk5qW7qoi4R0XUsnta1d1DaHw7quJby9eYSATp5u5GemUJCZSkFWCvMKMiLrKeRnpoaXs1KOb5uRlsy0aTZu/10KdBEJhL4+p76ti9qWTmqbu6hp6aC2Obxe09xJTSS8a5o7aWjrorfPh7yHGeRnpFCYFQ7oc+fkUpCZQmFWCgVZqRREArogM5X8rBSyU5MwG7+AHqtRA93MHgOuB6rdfekw7R8D7oustgD3uPvWmFYpIgnJ3Wnu7KG6qZPq5o5wMDcPDejalvDY9HAhnZI0jVBWKoVZKZTkprFsTs7xwC7ISqUwMxLWWSnkZaQwfRzPoMdbNGfojwM/AH46QvvbwLvcvcHMrgUeAS6OTXkiEkS9fU5dayfVTeFgrokEdnUkrKub3wnwju6+Ia9PmT6NwqwUCrNTKc5J47ySHAqzU8LBnZ16/LkwK5UZaZPrLHo8jRro7r7OzOafpP1PA1Y3AHNiUJeITEHuTkNbN0ebOgY8Oocs17Z0MszJNDPSkpg5I41QVioXzM1jZnYqoexUZmanHV8OZaeSk56cMCE9FrEeQ/808JuRGs3sLuAugLlz58b4o0VkPHV093KksYOqxnaqmzo5EgnpwctdvUPPqPMzU5iZncqsnDTOnp1N0YyBAf3Oclry9Dj8lwVHzALdzN5NONAvH6mPuz9CeEiG8vLyYf5+FpF4aOvq4XBjB0caOzjc2MHhY+0cbnpn/UhjOw1t3UNel5kynaKcNIqy0yifl3d8eVZOGkUzImfWM1JJTVJQT4SYBLqZnQc8Clzr7nWxeE8RiY3OnvCZ9aGGdg4da6fqWAdHmiLPjR0cbmynqWPotLz8zBRmzUijOCeNC+bmUpybzqwZ4bAOB3YaWamaKDeZnPafhpnNBZ4CPuHuu06/JBGJlrvT2N5NZUM7VcfCj/7QPhRZrmnuHPK6wqwUZuekM7cgg4vL8pmVk8bsnDRmzUinODcc1hr+mHqimbb4BLASKDSzSuCrQDKAu68CHgAKgB9FLlL0uHv5eBUskkjcnbrWLg7Wt3GwoZ3KhjYO1g8M7nbaunpPeE1q0jRKctMpzk3n3YtDlORmUJybdnzb7Nw0DYEEVDSzXG4bpf1O4M6YVSSSYJo6usOBXd8f2G1UNrRzsCH8PDiw8zNTKMlN58xQFlcuDFGcm8acvHBYF+emU5CZohkgCUoDYCLjrLfPqTrWzv66NvbVtbK/rpWD9eHAPljfNmT8Ois1iTl56cwryOTyM0OU5qdTmpdBaX4Gc/LSydS4tYxAR4ZIDHT39nGooT0S2G0nPB+sb6O7951JXSlJ05iTFw7p8+fmHg/r8HO65ljLKVOgi0Spp7ePgw3tVNS0sK+ujf11rcefKxvaT7jtPCNlOvMKMllclM17lsxiXkEG8woymF+QyawZaeP6BU2SuBToIoM0dXRTUdPK3uoW9taEHxU1reyraz3hTDs7LYn5BZmcW5LDDecVhwO7MJN5BRmEslJ1li0TToEuCamvzzl0rD0S2K2R0A4vD5zmlzTNmFeQwYJQFledXcSCUCZloSzKCjPJzdDQiEwuCnQJNPdwcO880szOo83h5yPNvF3bSmfPO7eo56QnsyCUycpFIRbMzGJBKIuyUCZz8zNInq5fapSpQYEugVHf2hUJ7Kbj4b3raMsJvx5TnJPGolnZXLGwkLJQOLgXhDLJ11Q/CQAFukw5bV097D7acuJZ99HmE4ZKcjOSWVyUzQcuKGHRrGzOmpXNwqJsZqQlx7FykfGlQJdJrbmjm+1VTWw71Mi2Q428caiRitpWPHJtMi15GouKslm5KMTiWdnhR1E2oWxdlJTEo0CXSaOxvZvtVf3B3cT2SHj3mzUjjaUlOdywrJizZ89gcVE2pfkZU/oXZkRiSYEucXGsrYtth5rYVhU+6952qJH9dW3H24tzwuF9y/klLJ2Tw9LiHELZqXGsWGTyU6DLuOvtc3YdbWbT/gY276vn1QPHOFD/TnjPyUtnaXEOHy4vZWlJDkuLZ1CQpfAWGSsFusRcW1cPWw4cY9P+Bjbtb+C1/Q00R2aahLJTuXBuHreuKOXckvCZd15mSpwrFgkGBbqctqNNHWza18DGffVs3t/AjsNNx2+DX1yUzQ3Liymfl0f5vHxK89N1sVJknCjQZUzcnZ1Hm9m4Lzx8sml/A5UN7UB4xsmyObnc/a4yyuflc8HcPHIyNE1QZKIo0GVUtS2d/GF3Let21bBudy21LeH53qHsVMrn5XH7ZfMpn5/POcUzdFelSBwp0GWIrp4+Xj3QEAnwGrYdagIgLyOZKxaGuGJhIRefUaDhE5FJRoEuAOyrbWXd7hrW7apl/d5aWrt6SZpmXDA3jy++ZxFXLgqxtDhHX/sqMokp0BNUS2cPf9pTezzE+6cRluanc/P5JVy5KMRlCwrI1q3yIlOGAj2B1DR3svaNw6x94zCb9zfQ0+dkpEzn0rICPn35GVy5KMT8ggwNo4hMUQr0gGts7+b57Uf49dYq/rinlj6HRUVZ3HlFGVcuKuTCeXn6BXiRgFCgB1B7Vy+/f+soa7ZU8cLOGrp6+yjNT+eelQu4cVkJi2dlx7tEERkHCvSA6O7t46XdNazZUsVvdxyltauXUHYqH7tkLjcuK2Z5aa6GUkQCToE+hfX1OS+/Xc+arVX8ZtthjrV1k5OezA3LirlxWTEXlxXomwhFEsiogW5mjwHXA9XuvnSYdgP+CbgOaANud/dXY12ohLk72w418cyWQzz7+mGONHWQnjydP19SxI3LirlyUYiUJN3cI5KIojlDfxz4AfDTEdqvBRZGHhcDD0eeJYb6+pzfv1XNwy/s4dUDx0iZPo13LQ7x5WVnc9XZM8lI0T+2RBLdqCng7uvMbP5JutwE/NTdHdhgZrlmNtvdD8eqyETW3dvHM1uq+OcX97K7uoU5eel87cZzuHl5ib4nRUROEIvTuhLg4ID1ysi2IYFuZncBdwHMnTs3Bh8dXG1dPax+5SCPvlRBVWMHZ83K5p9uXc77zp1Nkr4vRUSGEYtAH+6qmw/X0d0fAR4BKC8vH7ZPoqtv7eInf9rHT9bv41hbNyvOyOcfbzmXlYtDmqUiIicVi0CvBEoHrM8BqmLwvgnl0LF2Hn2pgtWvHKS9u5erzy7inpVlXDgvP96licgUEYtAXwPca2arCV8MbdT4efR2HW1m1Yt7WbMl/HfgjcuLuftdC1hUpJt/RGRsopm2+ASwEig0s0rgq0AygLuvAtYSnrK4h/C0xTvGq9gg2by/nodf2Mvv3qwmPXk6n7h0HndeUUZJbnq8SxORKSqaWS63jdLuwGdiVlHAvXWkiQee3s4r++rJzUjmc1ct5FOXzSdfv6spIqdJk5cnSF+f89gf3+abz+1kRnoSD1y/hFtXlGr+uIjEjNJkAhxt6uCLv9jKS7trufrsmfyvD5xHYVZqvMsSkYBRoI+z57Yd4UtPvU57dy//eMtSPrpirqYfisi4UKCPk9bOHr7+7A5WbzzI0pIZPPSR8zlzZla8yxKRAFOgj4MtB4/x+dWvsb++jXtWLuALVy/SF2aJyLhToMdQb5/z8At7+O7vdlOUncoTf3kJl5QVxLssEUkQCvQYOVjfxl//2xY27mvghmXF/MPNS8lJ15dnicjEUaDHwNOvHeIrT28D4LsfWcbNy0t04VNEJpwC/TQ0tnfzlae3sWZrFRfNz+M7H15OaX5GvMsSkQSlQD9FL1fU8df/tpUjTR188T2LuGflmfq5NxGJKwX6GPX09vGd3+7i4Rf3Mi8/g1/dcxnLS3PjXZaIiAJ9rL792108/MJebr2olK9cv4TMVO1CEZkclEZj8NLuGla9uJfbVpTyjfefF+9yREROoLtdolTT3MkXntzKmaEsHrj+nHiXIyIyhM7Qo9DX5/y3X2yluaObn9+5gvSU6fEuSURkCJ2hR+HRP1SwblcNf3v9Es6aNSPe5YiIDEuBPoqtB4/xzed28t5zivj4xXPjXY6IyIgU6CfR3NHNZ594jZnZqTz4gfN096eITGoaQx+Bu/O3T2+jsqGNJ//qUnIz9BNxIjK56Qx9BL/cXMkzW6r4/NWLuGh+frzLEREZlQJ9GHtrWnjgme1cUpbPZ959ZrzLERGJigJ9kM6eXj77f14jLXkaD33kfH0/i4hMGRpDH+Qba99ix+EmHv1kObNy0uJdjohI1HSGPsDvdhzl8T/t4/bL5nP1kqJ4lyMiMiYK9IgjjR38919uZcnsGXzpurPiXY6IyJhFFehmdo2Z7TSzPWZ2/zDtOWb2azPbambbzeyO2Jc6fnr7nM+tfo3Onj6+/9HzSU3Srf0iMvWMGuhmNh34IXAtsAS4zcyWDOr2GWCHuy8DVgLfNrMpM3H7h/+5h5ffrudrN57DglBWvMsRETkl0ZyhrwD2uHuFu3cBq4GbBvVxINvCt1JmAfVAT0wrHScb99Xz0O92cdPyYj544Zx4lyMicsqiCfQS4OCA9crItoF+AJwNVAFvAJ9z977Bb2Rmd5nZJjPbVFNTc4olx86xti4+98RrlOZn8A83L9Wt/SIypUUT6MOlnA9afy+wBSgGlgM/MLMhX0vo7o+4e7m7l4dCoTGWGlvuzn2/ep3q5k6+d+v5ZKclx7UeEZHTFU2gVwKlA9bnED4TH+gO4CkP2wO8DUzqqSI/f/kAz28/yt9cs5hl+k1QEQmAaAJ9I7DQzM6IXOi8FVgzqM8B4CoAMysCFgMVsSw0lt460sTXn93BlYtC3Hl5WbzLERGJiVHvFHX3HjO7F3gemA485u7bzezuSPsq4OvA42b2BuEhmvvcvXYc6z5l7V3hW/tnpCXz7Q8tY5pu7ReRgIjq1n93XwusHbRt1YDlKuA9sS1tfPzfNw6zu7qFH99xEaHs1HiXIyISMwl3p+j6vXXkZ6bwroXxvSgrIhJrCRXo7s6GijouKcvXUIuIBE5CBXplQzuHjrVzSVlBvEsREYm5hAr09XvrALhUgS4iAZRYgV5RR2FWCmfO1Pe1iEjwJEyg94+fX1xWoFv8RSSQEibQ99e1cbixQ+PnIhJYCRPoGyo0fi4iwZYwgb6+oo5QdioLQpnxLkVEZFwkRKC/M/9c4+ciElwJEehv17ZytKmTS8ry412KiMi4SYhA31BRD2j8XESCLSECfX1FHTOzUzmjUOPnIhJcgQ/0/vHzSxdo/FxEgi3wgb63ppWa5k7NPxeRwAt8oGv+uYgkisAH+vqKOmbnpDGvICPepYiIjKtAB7q787Lmn4tIggh0oO+pbqG2pUvDLSKSEAId6P3j57ogKiKJINCBvr6ijpLcdErz0+NdiojIuAtsoIfnn9dzcVm+xs9FJCEENtB3HW2hvlXj5yKSOAIb6Bo/F5FEE1Wgm9k1ZrbTzPaY2f0j9FlpZlvMbLuZvRjbMsdu/d465uSlU5qv+ecikhiSRutgZtOBHwJ/DlQCG81sjbvvGNAnF/gRcI27HzCzmeNUb1T6+pyX367jqrOL4lmGiMiEiuYMfQWwx90r3L0LWA3cNKjPR4Gn3P0AgLtXx7bMsdl5tJmGtm6Nn4tIQokm0EuAgwPWKyPbBloE5JnZC2a22cw+OdwbmdldZrbJzDbV1NScWsVRWL83Mn6+QIEuIokjmkAfbs6fD1pPAi4E3ge8F/iKmS0a8iL3R9y93N3LQ6HQmIuN1oaKOubmZ1CSq/nnIpI4ogn0SqB0wPocoGqYPs+5e6u71wLrgGWxKXFswuPn9fq5ORFJONEE+kZgoZmdYWYpwK3AmkF9ngGuMLMkM8sALgbejG2p0XnzSBON7d1cquEWEUkwo85ycfceM7sXeB6YDjzm7tvN7O5I+yp3f9PMngNeB/qAR91923gWPpLj4+e6ICoiCWbUQAdw97XA2kHbVg1a/xbwrdiVdmo2VNQzvyCD2TkaPxeRxBKoO0V7I/PPdXYuIokoUIH+5uEmmjt6NH4uIgkpUIGu8XMRSWSBCvQNFXWUFWZSNCMt3qWIiEy4wAR6T28fr7xdz8U6OxeRBBWYQN9xuInmTo2fi0jiCkygvzN+rjtERSQxBSbQN1TUsSCUycxsjZ+LSGIKRKD39PaxcV+DhltEJKEFItC3VTXR0tmj6YoiktACEeiafy4iEpBA31BRx8KZWRRmpca7FBGRuJnygd7d28fGffUaPxeRhDflA/2NQ420dfVquEVEEt6UD/T+8fOLz9D8cxFJbFM+0DdU1LG4KJsCjZ+LSIKb0oHe1dPHJs0/FxEBpnigv155jPbuXt3uLyLCFA/0DRX94+c6QxcRmdKBvr6ijrNmZZOXmRLvUkRE4m7KBnpnTy+b92v8XESk35QN9K0HG+no7tP8cxGRiCkb6Bsq6jDT/HMRkX5TNtDX763j7FkzyM3Q+LmICEQZ6GZ2jZntNLM9Znb/SfpdZGa9ZvbB2JU4VEd3L68e0Pi5iMhAowa6mU0HfghcCywBbjOzJSP0exB4PtZFDrbl4DE6ezR+LiIyUDRn6CuAPe5e4e5dwGrgpmH6fRb4FVAdw/qGtaGijmkGKzR+LiJyXDSBXgIcHLBeGdl2nJmVALcAq072RmZ2l5ltMrNNNTU1Y631uPV76zinOIec9ORTfg8RkaCJJtBtmG0+aP0h4D537z3ZG7n7I+5e7u7loVAoyhJP1NHdy2sHj+l2fxGRQZKi6FMJlA5YnwNUDepTDqw2M4BC4Doz63H3p2NR5ECvHmigq6dPF0RFRAaJJtA3AgvN7AzgEHAr8NGBHdz9jP5lM3sceHY8whwgefo03r04RPl8naGLiAw0aqC7e4+Z3Ut49sp04DF3325md0faTzpuHmsXzc/nx3esmMiPFBGZEqI5Q8fd1wJrB20bNsjd/fbTL0tERMZqyt4pKiIiJ1Kgi4gEhAJdRCQgFOgiIgGhQBcRCQgFuohIQCjQRUQCwtwHfy3LBH2wWQ2w/xRfXgjUxrCcWJvs9cHkr1H1nR7Vd3omc33z3H3YL8OKW6CfDjPb5O7l8a5jJJO9Ppj8Naq+06P6Ts9kr28kGnIREQkIBbqISEBM1UB/JN4FjGKy1weTv0bVd3pU3+mZ7PUNa0qOoYuIyFBT9QxdREQGUaCLiATEpA50M7vGzHaa2R4zu3+YdjOz70XaXzezCyawtlIz+08ze9PMtpvZ54bps9LMGs1sS+TxwETVF/n8fWb2RuSzNw3THs/9t3jAftliZk1m9vlBfSZ8/5nZY2ZWbWbbBmzLN7PfmtnuyHPeCK896fE6jvV9y8zeivwZ/ruZ5Y7w2pMeD+NY39+Z2aEBf47XjfDaeO2/JwfUts/Mtozw2nHff6fN3Sflg/CvI+0FyoAUYCuwZFCf64DfEP4h60uAlyewvtnABZHlbGDXMPWtJPxzfPHah/uAwpO0x23/DfNnfYTwDRNx3X/AlcAFwLYB274J3B9Zvh94cIT/hpMer+NY33uApMjyg8PVF83xMI71/R3wxSiOgbjsv0Ht3wYeiNf+O93HZD5DXwHscfcKd+8CVgM3DepzE/BTD9sA5JrZ7Ikozt0Pu/urkeVm4E2gZCI+O4bitv8GuQrY6+6neudwzLj7OqB+0OabgJ9Eln8C3DzMS6M5XselPnf/D3fviaxuIPxD7nExwv6LRtz2Xz8L/8r9h4EnYv25E2UyB3oJcHDAeiVDAzOaPuPOzOYD5wMvD9N8qZltNbPfmNk5E1sZDvyHmW02s7uGaZ8U+4/wD4+P9D9RPPdfvyJ3Pwzhv8iBmcP0mSz78i8I/6trOKMdD+Pp3siQ0GMjDFlNhv13BXDU3XeP0B7P/ReVyRzoNsy2wXMso+kzrswsC/gV8Hl3bxrU/CrhYYRlwPeBpyeyNuDP3P0C4FrgM2Z25aD2ybD/UoAbgV8M0xzv/TcWk2FffhnoAf51hC6jHQ/j5WFgAbAcOEx4WGOwuO8/4DZOfnYer/0Xtckc6JVA6YD1OUDVKfQZN2aWTDjM/9Xdnxrc7u5N7t4SWV4LJJtZ4UTV5+5Vkedq4N8J/7N2oLjuv4hrgVfd/ejghnjvvwGO9g9FRZ6rh+kT72PxU8D1wMc8MuA7WBTHw7hw96Pu3uvufcC/jPC58d5/ScD7gSdH6hOv/TcWkznQNwILzeyMyFncrcCaQX3WAJ+MzNa4BGjs/6fxeIuMt/1v4E13/84IfWZF+mFmKwjv77oJqi/TzLL7lwlfONs2qFvc9t8AI54VxXP/DbIG+FRk+VPAM8P0ieZ4HRdmdg1wH3Cju7eN0Cea42G86ht4XeaWET43bvsv4mrgLXevHK4xnvtvTOJ9VfZkD8KzMHYRvvr95ci2u4G7I8sG/DDS/gZQPoG1XU74n4SvA1sij+sG1XcvsJ3wFfsNwGUTWF9Z5HO3RmqYVPsv8vkZhAM6Z8C2uO4/wn+5HAa6CZ81fhooAH4P7I4850f6FgNrT3a8TlB9ewiPP/cfh6sG1zfS8TBB9f0scny9TjikZ0+m/RfZ/nj/cTeg74Tvv9N96NZ/EZGAmMxDLiIiMgYKdBGRgFCgi4gEhAJdRCQgFOgiIgGhQBcRCQgFuohIQPx/JL7vYhl0adYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "x = np.arange(0,P,1)\n",
    "y = np.array(lam_history)\n",
    "plt.plot(x,y)\n",
    "plt.axhline(np.pi/2,linestyle='--')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## トレーニングデータでの確認"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {},
   "outputs": [],
   "source": [
    "def predict(alpha, y_train, x_train, z, Lambda):\n",
    "    kernel_value = np.zeros(len(x_train))\n",
    "    for i in range(len(x_train)):\n",
    "        kernel_value[i] = y_train[i] * alpha[i] * get_kernel_value(x_train[i], z, Lambda)\n",
    "    return np.sign(np.sum(kernel_value))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0]\n",
      "[ 1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.\n",
      "  1.  1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1.\n",
      " -1. -1. -1. -1.]\n",
      "[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
      " 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
      "0.0\n"
     ]
    }
   ],
   "source": [
    "predict_train = []\n",
    "lam_vector = np.full(n_qubits, lam_history[-1])\n",
    "for i in range(len(x_train)):\n",
    "    predict_train.append(predict(alpha_history[-1], y_train, x_train, x_train[i], lam_vector))\n",
    "print(predict_train)\n",
    "print(y_train.reshape(-1))\n",
    "print(predict_train - y_train.reshape(-1))\n",
    "print(np.sum(np.abs(predict_train - y_train.reshape(-1)))/2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## テストデータでの確認"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "テストデータ生成"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_data_test = 20\n",
    "epsilon = 0.0 # エラーの分散\n",
    "\n",
    "l = [i for i in range(2**n_qubits)] # 0 = 2**n_qubits-1 のリスト\n",
    "select = random.sample(l, n_data_test)\n",
    "select_bin = [format(select[i], f'0{n_qubits}b') for i in range(n_data_test)]\n",
    "\n",
    "x_plus = []\n",
    "x_minus = []\n",
    "for i in range(n_data_test):\n",
    "    xx = x_generator(c_plus, select_bin[i]).reshape((1,-1))\n",
    "    e = np.random.normal(0, np.sqrt(epsilon), (1, 2*n_qubits))\n",
    "    x_plus.append(xx + e)\n",
    "    xx = x_generator(c_minus, select_bin[i]).reshape((1,-1))\n",
    "    e = np.random.normal(0, np.sqrt(epsilon), (1, 2*n_qubits))\n",
    "    x_minus.append(xx + e)\n",
    "x_plus = np.concatenate(x_plus)\n",
    "x_minus = np.concatenate(x_minus)\n",
    "x_test = np.concatenate([x_plus,x_minus])\n",
    "y_plus = np.ones(n_data_test)\n",
    "y_minus = -np.ones(n_data_test)\n",
    "y_test = np.concatenate([y_plus, y_minus]).reshape(1,-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1.0, -1.0, -1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, 1.0, -1.0, -1.0, 1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, 1.0]\n",
      "[ 1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.\n",
      "  1.  1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1. -1.\n",
      " -1. -1. -1. -1.]\n",
      "[ 0. -2. -2.  0.  0.  0.  0.  0. -2.  0.  0.  0.  0.  0.  0.  0.  0.  0.\n",
      "  0.  0.  0.  2.  0.  0.  0.  2.  0.  0.  2.  0.  0.  0.  0.  0.  0.  0.\n",
      "  0.  0.  0.  2.]\n",
      "7.0\n"
     ]
    }
   ],
   "source": [
    "predict_test = []\n",
    "lam_vector = np.full(n_qubits, lam_history[-1])\n",
    "for i in range(len(x_test)):\n",
    "    predict_test.append(predict(alpha_history[-1], y_train, x_train, x_test[i], lam_vector))\n",
    "print(predict_test)\n",
    "print(y_test.reshape(-1))\n",
    "print(predict_test - y_test.reshape(-1))\n",
    "print(np.sum(np.abs(predict_test - y_test.reshape(-1)))/2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Under construction..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "class KernelMatrix:\n",
    "\n",
    "    def __init__(self,n_qubit,Lambda,edges) -> None:\n",
    "        self.n_qubit = n_qubit\n",
    "        self.Lambda = Lambda\n",
    "        self.edges = edges\n",
    "        self.state1 = QuantumState(n_qubit)\n",
    "        self.state2 = QuantumState(n_qubit)\n",
    "\n",
    "    def make_fiducial_state(self,Lambda):\n",
    "        circuit = ParametricQuantumCircuit(self.n_qubit)\n",
    "        for i in range(self.n_qubit):\n",
    "            circuit.add_parametric_RY_gate(Lambda)\n",
    "        for control, target in self.edges:\n",
    "            circuit.add_CZ_gate(control,target)\n",
    "        return circuit\n",
    "\n",
    "    def create_featured_gate(self,theta):\n",
    "        circuit = ParametricQuantumCircuit(self.n_qubit)\n",
    "        for i in range(self.n_qubit):\n",
    "            circuit.add_parametric_RY_gate(i,self.Lambda)\n",
    "        for control, target in self.edges:\n",
    "            circuit.add_CZ_gate(control,target)\n",
    "\n",
    "        for i in range(self.n_qubit):\n",
    "            circuit.add_RZ_gate(i,theta[i,0])\n",
    "            circuit.add_RX_gate(i,theta[i,1])\n",
    "\n",
    "        return circuit\n",
    "        \n",
    "    def featured_gate(self,theta):\n",
    "        circuit = QuantumCircuit(self.n_qubit)\n",
    "        for i in range(self.n_qubit):\n",
    "            circuit.add_RZ_gate(i,theta[i,0])\n",
    "            circuit.add_RX_gate(i,theta[i,1])\n",
    "\n",
    "        return circuit\n",
    "\n",
    "    def get_kernel_value(self,x1,x2):\n",
    "        self.state1.set_zero_state()\n",
    "        self.state1.set_zero_state()\n",
    "        fiducial_gate = self.make_fiducial_state()\n",
    "\n",
    "\n",
    "        featured_gate1 = self.featured_gate(x1)\n",
    "        featured_gate2 = self.featured_gate(x2)\n",
    "        \n",
    "        circuit1 = merge(fiducial_gate,featured_gate1)\n",
    "        circuit1.update_quantum_state(state1)\n",
    "        circuit2 = merge(fiducial_gate,featured_gate2)\n",
    "        circuit2.update_quantum_state(state2)\n",
    "\n",
    "        return inner_product(state1,state2)\n",
    "    \n",
    "    def construct_kernel_matrix(self,x_data):\n",
    "        \"\"\"\n",
    "        Args:\n",
    "            x_data : MxNxD array of data where N is the number of samples,\n",
    "             N the number of qubits, D the feature dimension\n",
    "        \"\"\"\n",
    "        kernel_matrix = np.identity(x_data.shape[0])\n",
    "        for i in range(x_data.shape[0]):\n",
    "            for j in range(i):\n",
    "                kernel_matrix[i,j] = self.get_kernel_value(x_data[i],x_data[j])\n",
    "                kernel_matrix[j,i] = kernel_matrix[i,j].conjugate()\n",
    "        \n",
    "        return kernel_matrix\n",
    "\n",
    "    def spsa(self,x_data,y_label,n_step,initial_Lambda = 0.01):\n",
    "        Lambda = initial_Lambda\n",
    "        for i in n_step:\n",
    "            Delta = 2 * np.random.rand() - 1\n",
    "\n",
    "            c_i = c / ( i + 1)**gamma\n",
    "            Lambda__pm_i = Lambda + np.array([Delta*c_i, - Delta*c_i])\n",
    "            self.construct_kernel_matrix()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([1, 2])"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.array([0,1 ]) + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       " *** Quantum State ***\n",
       " * Qubit Count : 2\n",
       " * Dimension   : 4\n",
       " * State vector : \n",
       "(1,0)\n",
       "(0,0)\n",
       "(0,0)\n",
       "(0,0)"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "KernelMatrix(2,np.pi,[[0,1]]).state1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "state1 = QuantumState(n)\n",
    "state2 = QuantumState(n)\n",
    "def get_kernel_value(x1,x2):\n",
    "    state1.set_zero_state()\n",
    "    create_featured_gate()"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "4d6cf19fbe203395831accd42fb3bdc3b710eb7fbe04f788a536db11ba6ce36d"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}