{"cells":[{"cell_type":"markdown","metadata":{},"source":["<br />\n","\n","<div style=\"text-align: center;\">\n","<font size=\"7\">数値計算試験問題</font>\n","</div>\n","<br />\n","<div style=\"text-align: right;\">\n","<font size=\"4\">2022/12/23 実施</font>\n","<br />\n","<font size=\"4\">cc by Shigeto R. Nishitani 2022 </font>\n","</div>\n","\n"]},{"attachments":{"image.png":{"image/png":"iVBORw0KGgoAAAANSUhEUgAAAXIAAAD8CAYAAABq6S8VAAAAAXNSR0IArs4c6QAAAERlWElmTU0AKgAAAAgAAYdpAAQAAAABAAAAGgAAAAAAA6ABAAMAAAABAAEAAKACAAQAAAABAAABcqADAAQAAAABAAAA/AAAAACR+ioEAAAvqUlEQVR4Ae2dB3xUVfbHfwmEZuiETgpIQDok9A6KdEFBRSyshT+6oIJYEF0FQdZCseC6rN3FAoKIVBFC7116DR2kKQlNQuZ/zpvMOsSUSTKTee/N73w+hzcz7777zvme4czNfbcAFBIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARLIEoEgD0vHS7kE0euiSaKxohQSIAESIAELEYgXW0tZyF6aSgIkQAIBQyA4YDyloyRAAiRgUwKedq0cFP/PizpE/y06STS19JcPVFGgQIGY8PDw1Odt8T45ORnBwfb9/aN/1v6aMn7Wjd+ePXvOiPVh2fHA00ReQSo/JlpadIHoINGlomlKdHS0Y/fu3Wmes/qHixcvRps2bazuRrr207900VjiBONniTClaWRQUNAGOZGt54+eNi01iav8Kvq9aCN9QyEBEiABEvA/AU8S+U1iZuEUU/V1B9FtKe95IAESIAES8DOBvB7cv4yU0Va4ipb/SnSevqGQAAmQAAn4n4AnifyAmFnX/6bSAhIgAW8QuHbtGo4ePYorV654o7pcraNo0aLYuXNnrt7T2zeTwSCoWLEiQkJCvFa1J4ncazdjRSRAAv4noEm8cOHCiIyMhDxg879BWbAgISHBsD0Ll5iqqMPhwNmzZ40f0qioKK/Z5kkfudduxopIgAT8T0Bb4iVLlrRcEvc/uZxboD+cyt7bfw0xkec8NqyBBCxHwGotccsBzsBgX7BnIs8AOE+RAAmQgBUIMJFbIUq0kQQClEDnzp3x22+/Zej9P/7xD/z8888ZlknvpE6g6tq1a3qnLfM5H3ZaJlQ0lAQCh4A+FFSdM2dOpk6PHDky0zJ2L8AWud0jTP9IwKQExo0bh1q1ahk6YcIExMfHo1q1anjwwQeNz44cOWKMrDlzRpcgAV577TU0aNAALVq0QJ8+ffD2228bn/fr1w/fffed8VpH4rzyyitGudq1a2PXrl3G52vXrkXTpk1Rv359NGvWDHZbQoQtciPM/IcEApPA008Dmzd71/d69QDJyxnKhg0b8Omnn2LNmjVGy7tx48Zo3bo19u7di88//xxNmjS54fp169Zh2rRpWLlypS7KZyTqmJiYG8q43pQqVQobN27EBx98YCT7jz76CNWrV8eyZcuQN29eoxvmxRdfNOpzXWP1IxO51SNI+0nAggSWL1+Onj174qabdNUP4M477zQSbURExF+SuJ5fsWIF7rjjDiOJ6xj4bt266cdpitalool++vTpxuvff/8dDz30kPFDoaNGdFKUnYSJ3E7RpC8kkEUCmbWcs1hdjou7EntOKsqfP79xeZ48eZCUpBuaAS+//DLatm2L77//3ujCsdsKpuwjN8LMf0iABHKTQMuWLTFjxgxcunQJFy9eNBKsfpaeNG/eHD/++KMxkSYxMRGzZs1Kr2ian2uLvEIFXY0b+Oyzz4yjnf5hIrdTNOkLCViEgD601IeUjRo1gvaPP/rooyhevHi61jds2BDdu3c3Hlh26tQJ+iBT113xVJ577jkMGzbMeNjpaqV7em3AltONJewqcXFxdnXN8Iv+WTu8nsRvx44dlnRS1llxXLhwwSEteIf0fzvkgakl/VCj04qB/GCsz+6PBvvIs0uO15EACeQqgf79+2Pbtm34448/jAeX2qqnOAkwkfObQAIkYAkCX331Fay++qGvQLOP3FdkWS8JkAAJ5BIBJvJcAs3bkAAJkICvCDCR+4os6yUBEiCBXCLARJ5LoHkbEiABEvAVASZyX5FlvSRAAukS0IWrMhNdSEsnDPladILQwIEDM7yNLner67xkVXQRL9eiX1m9NivlmcizQotlSYAEvELAk6SYnUR+/fp1r9iXupLsJvLU9fjqPRO5r8iyXhIggXQJhIaGGuc0Qeq6J7169TJWKOzbt6+xGuK7776L48ePG+uj6BopKj/99BPat29vrHzYu3dv6FR9FW31Pv/888bnU6dONep76qmnUE+WYdRlcnUJW5Vz586hR48eqFOnjrEw19atW43P3f/RZQB0pqkud3vrrbfi1KlTxtosH374IcaPH2/Uqasonj59GnfddRd0xqmqLuqlohsrd+jQATVr1jRmq8rcH/fqffaa48h9hpYVk4AFCPhrHVs3NJs2bcL27dtRvnx56JoqmhSffPJJ6HrlMlMVuiytdk+MGjUKM2fORNmyZfHGG28Y53V3IBXd0FiXrlXRpKtdMptlfd6lS5fi4YcfNiYS6TrlmqB1jZdFixYZ655rGXfRtc5Xr15tbEyty9+++eabGDt2LAYMGAD98Rk6dKhR/L777sPgwYONtdEPHz6M22+/HTt37sSIESOMz9Su2bNn4+OPP3av3mevmch9hpYVkwAJeEJA11upWLGiUVRb0fGywYQmVHfR5CrT2o3WbnBwsDG7UzeKcMk999zjemkcdeMJlVatWkGm9RvbxenSubqmuUq7du2M1rOec5ejR49C6zpx4oRxj6ioKPfT/3utW8upPS7RevQvBP3hcC2d26VLlwzXj3Fd640jE7k3KLIOErAqAROsY+tadlYRui89645Uuyhuu+02TJo0CboeeWpJvfxt6p3qU79Pfb3r/aBBgzBkyBBjgS7t9nn11Vddp244JicnGy133eTCDMI+cjNEgTaQAAn8hYAmbJ2Sr6I7BmmXy/79+433uvTtnj17jNdp/fPtt98aH2srXFdJVNVlcidPnmx8rklau2yKFClyw+Xuy93qTkUucbdFP9N+8Pfee8912ujG0Tf6F4AuJaAyd+5cnD9/3njt63+YyH1NmPWTAAlki4AuktWxY0fjgWdYWJixjrj2d+vDSu1Wce3HmVbl2lLW/nDt23b1U2vrWreY0+tfeOEFY0u51NdqGX2QqrsLaaJ3ie5IpJtSaNePPuzUh7Hr16836qpRo4bRL69ltR9eu1f0Yad2sYSHh7uqsN6Ry9jqQpXWFE+WQbWmZ06r6V/aS6haJaa6jG1mInt/OmSPz8yK+fW8t5exZYvcer+TtJgESIAEbiDAh5034OAbEiABqxPQ/u9AE7bIAy3i9JcEhID0K5CDnwj4gj0TuZ+CyduSgL8I6INAnYHoi4TiL5+scl9lruy9PWyRXStW+QbQThLwEgGdfKMTX3SaudXkypUrXk+Cuc1Ak7hrApS37s1E7i2SrIcELEIgJCQE6c1YNLsL2v+twwopNxJg18qNPPiOBEiABCxHICuJPI94t0l0VmZeFtYZV7IimUyjyqwoz5MACZAACeSQQFYS+VNyr50e3+/QIUBmZjGZe0yMBUmABEggWwQ87SPXpcm6iI4WHZLZnZKQF6dRDEGXHLjywDBcfPlLJFeMQP4alVG+ZxMUaNkQKFQos2p4ngRIgARIwAMCQR6U0SLfiY4R1WXHdEHerqKpRZrfUEU0Cse8i6a4Jgk9BEkoGXQO4Y5DKC3pXUU/P1S8Jn6NbYKgu5vgapUIyALAxjmz/6NLVboWxTe7rdmxj/5lh5p5rmH8zBOLrFoiG2hskGtis3qdp+U1aX+QUriNHDPtI4/R+QYujYhwJCU5HPHxDse8yWccn/We5fii0jDHiqDmjusIMsodK1rdsfeBEY7rR4/LMEtzC9fqMHd8MrOO8cuMkLnP2zl+klvXp+TZLB886SNvLrV2F40X/Ua0neh/RTMX7T4ZPVrWGAYipNF9+30l8dCULnjg8Ouof3E55n9yHJ81+gD7E8vi5i9fwfWK4dhRrw/Oz16Zed0sQQIkQAIkYBDwJJEPk5LaRx4peq/oItH7RTMWzdyyCDxkD760pGBBoNPfyqLfmsfR6GIcZo/fgxkVBqLCljko3rU5dldoh5PfLU/rUn5GAiRAAiTgRsCTRO5W3LOXCdHR0n6PTzeJp64lf355kvp0VfQ+Oh7H1x7DlKbjUfT4DpTt3RLbK3bA8ZnZ/osj9a34ngRIgARsRyCriXyxEEjrQafXwNzSMBR3r3waSbsPYHqztxF2bDPK3tEI6+s+gsT9p7x2H1ZEAiRAAnYhkNVEnmt+V4wuhDtXPIM/tu/D3BpDUWerDGGsGo01fcYj+Y+kXLODNyIBEiABsxMwbSJ3gatYowi6bH8TO7/9BduLNkPjb4Zgd6lmODRU9suLjARkR23OInXR4pEESCAQCZg+kbuCUvfuamh8Zg4W9f8GYQkHUW7sM1h2qBKuOcQFziJ1YeKRBEggAAlYJpFrbILzBKHdv+/B9bIVsBaN0BLLsUemHx2CbHB66RIwfHgAhpAukwAJBDoBSyVyV7DKnNqKFliBlWiC8jiBEjiHJWgFh7bMKSRAAiQQYAQsmcgRLi1wkWZYjSvIj/2ojNZYiuV5WuP3Y4kBFkK6SwIkEOgErJnIZbaoa9GtcjiFOtiKZcGt0ez6MpyOaoSd0z1fpDHQvwD0nwRIwPoErJnIdbaozhrV2aOy2FawHFt+8Rh2TFiAYtfPoNJdDbFowBTrR4cekAAJkIAHBDxdxtaDqnK5iCbzVNP/a4sJZ9pswqHWdxsPReev3Yh2K0cjpIDuiUEhARIgAXsSsGaLPINYlKpbAdVPxGF1vQG4fdMbWFehB84cuJDBFTxFAiRAAtYmYLtEruHIUzAfmmz6F9b2+wCNzs3F+epNsGv2fmtHitaTAAmQQDoEbJnIXb42+vRx7P/XApRKOoVS3ZpgzXguj+tiwyMJkIB9CNg6kWuYqg1oi6tLVuNSSDHUHdIOP/fnQ1D7fH3pCQmQgBKwfSJXJ8u2rIoSu1dhf7FY3PqfezC3zRtIvi57GFFIgARIwAYEAiKRa5xCI0uh2pGfsa7Kvei05AX8XPMpXLuabIMQ0gUSIIFAJxAwiVwDnTe0AGJ3T8bqZoPRYfd7WBnZB4lnrwb6d4D+kwAJWJxAQCVyjVVQnmA0WTEOa3q9hdYnp2BXVEecPfC7xcNI80mABAKZQMAlclewG08dio2Dv0TdhOU4WbMdTmw97TrFIwmQAAlYikDAJnKNUoNx92PXmBmofGUHLsa2QvyyI5YKHo0lARIgASUQ0IlcAdR+oQsOTfoJpZOOI2+bFtj94x79mEICJEACliEQ8IlcI1X9sZY4M3UxCuAyit3RCtu/+cUyAaShJEACJMBEnvIdqHxXfVyetxSO4Dwod18bbPl4Pb8dJEACJGAJAkzkbmGqdFt1JC9ehot5iiDy0fbY+N4Kt7N8SQIkQALmJMBEniou5VtURv41y3AuX1lUe7ID1r0Zl6oE35IACZCAuQgwkacRj9INKqLIpqU4USAKtZ7vjHWjf0qjFD8iARIgAXMQYCJPJw4la5RByS1xOFIwGrVf6o41r8xJpyQ/JgESIAH/EmAiz4B/8egwlP5lEQ4Wqon6I3tgTdEOaN2uHRAZCUyenMGVPEUCJEACuUeAiTwT1sWqlET51wdiL6qi/oXFWOtoCBw6BPTvz2SeCTueJgESyB0CTOQecC46fgQq4gj2IBoNsBGr0Qi4dAkYPtyDq1mEBEiABHxLgIncE76HD6MoEhCOQ9iF6ohxJXP5nEICJEAC/ibARO5JBMLDjVJFkIgIxGMnbjGS+drC7T25mmVIgARIwKcEmMg9wTt6NFCokFFSk3kkDhrJ3Ogzf36aJzWwDAmQAAn4jAATuSdo+/YFJk0CIiLgCApCkYiSiHhnCHaGNkT9N+/FuheYzD3ByDIkQAK+IcBE7ilXTebx8ViyaJFxLPpkP0TsmGck83pv3Iv1L073tCaWIwESIAGvEvAkkReQO64V3SK6XXSEKEUIFK1UBOHbJZnfFIu6Y+7BhpdnkAsJkAAJ5DoBTxK5bmrZTrSuaD3RjqJNRClCoFh4EVSSZL6rUAxqj7obG1+dSS4kQAIkkKsEPEnkDrEoMcWqEDmq6meUFALFI4qi4vb52F2oPmqN6IVNr80iGxIgARLINQJBHt4pj5TbIHqz6ETR50VTi0x1hCrCwsJipkyZkvq8Ld4nJiYiNDQ0TV8uHruMSo+8hGpXt2H2I2+jxP210yxn5g8z8s/MdntqG/3zlJQ5y9k5fm3bttUcG5sb5IvJTXRd11oZ3Sw6OtphV4mLi8vQtTN7zzm2F2jguIJ8js3/nJthWTOezMw/M9qcFZvoX1Zoma+sneMnOTXbu9l40rXinrN/kzeayLWfnJIGgZI3F0fpLQuwv0BNVHuhB7a+zSVw08DEj0iABLxIwJNEHib305a4SkHR20R36RtK2gRKRZdAqU0/42D+W1D12TuwdeyCtAvyUxIgARLwAgFPEnk5uY+2wreKrhPVrMSneQIhIyldvQRKbPwZ8fmroerQ7vhlwsKMivMcCZAACWSbQF4PrtQEXt+DciySikCZGiXhWP8zDse0Q5XB3bAtzyzUGtQuVSm+JQESIIGcEfCkRZ6zOwT41WVrlUKxDQtxJF8VVH6yK7a/r3/cUEiABEjAewSYyL3HMt2aytQKQ9F1C3E0X2VEDeqC7RMXp1uWJ0iABEggqwSYyLNKLJvly9YpjcJrF+FYvihEDpRk/sGSbNbEy0iABEjgRgJM5Dfy8Om7cnVLI3TNIhzPF4nIv3dmMvcpbVZOAoFDgIk8l2Ndrl4ZI5mfCIkwkvm29xfnsgW8HQmQgN0IMJH7IaKazG9aG2e0zCsP6oxf3uUDUD+EgbckAdsQYCL3Uyg1mRdeF2c8AK3yVBdsnSDrnFNIgARIIBsEmMizAc1bl+gD0KLrF+Fo/iqoOrgLtrzNGaDeYst6SCCQCDCR+znaZWqXlnHmi3A4fzSqPdsNm8bM87NFvD0JkIDVCDCRmyBipWuGoeSWRThYoAZqvHgH1r/KFRBMEBaaQAKWIcBEbpJQlapWEmW2LcS+QnVQZ8SdWNt9FBAZCQRLiPQ4ebJJLKUZJEACZiPARG6iiJSoUhwVti/A7vx10eDHV7HqkKxX5nAAhw7Jlh2yZweTuYmiRVNIwDwEmMjNEwvDkmKRxRBRKhHbUBONZM/r5WjmtPDSJWD4cJNZS3NIgATMQICJ3AxRSGVDkeO7URV7sVn2um6BlVgq/xpy+HCqknxLAiRAAtIDSwgmJBAejptwWdrk26RNHotW0i6PQxtAPqeQAAmQQGoCTOSpiZjh/ejRQKFCKIA/pE2+GavQGG2xGIsKdoEjWfrMKSRAAiTgRoCJ3A2GaV727QtMmgRERCBf0HU0qnQKK8r0RLtdHyCu4XNIvs5kbppY0RASMAEBJnITBCFNEzSZx8cDycnIc/ggmh79Dsvq/B3tNr6NJTUGIOnq9TQv44ckQAKBR4CJ3CIxD84bjBab3sOyVsPRds8krKncB1cu/GER62kmCZCALwkwkfuSrpfrDgoOQsslo7C8x9tofnwqtkZ2x4UTF718F1ZHAiRgNQJM5FaLmNjb4vtnsPqxjxFzfgHib74Vv+48a0EvaDIJkIC3CDCRe4tkLtfTZNLD2DJ8KqIvbcKFui0Rv+xILlvA25EACZiFABO5WSKRDTsajLoT+yfOR+mkYwhp0ww7vtuRjVp4CQmQgNUJMJFbPII1n2iNs9OXIgRJKNu7BdaOW25xj2g+CZBAVgkwkWeVmAnLR/Woi+TlK3EhXxjqPHMr4gZOM6GVNIkESMBXBJjIfUU2l+st2zQKJXetwP6iDdB6Ym/M7/KusXBiLpvB25EACfiBABO5H6D76paFo0qh6qGF2FipB26f8xQW1B6Ca1c4cchXvFkvCZiFABO5WSLhJTvyFS2ImANTsarhIHTYPh5rwnvj9xOyBC6FBEjAtgSYyG0Y2qC8edB07btY3ecdNDs9A4crt8GRdSdt6CldIgESUAJM5Db+HjT56klse20GKl/ZDjRpjC1fbLGxt3SNBAKXABO5zWNf56Xu+HWqDE8MSkKVh5pj8ZCZNveY7pFA4BFgIg+AmEf1ikG+zetwNPQWtBrfA/Pbv8mlcAMg7nQxcAgwkQdIrEvUKo/KR5ZgQ1Rv3L7oeSyLegAJv14OEO/pJgnYmwATub3je4N3+YoVQuy+b7Ci8yi0PPIVjkS0wKHlh28owzckQALWI8BEbr2Y5chiXQq3+ezh2DryB1S6sheFWsVi/bilOaqTF5MACfiXgCeJvJKYGCeqKzLJ8Ac8JUqxOIF6L3fD+flrkRhSHPWeaYeFDYbCER6B1u3aAZGRwOTJFveQ5pNA4BDwJJEnCY5nRGuINhH9e8prOVCsTCC8Q3WEHViLjcXbo/2msVh5pCISHYWAQ4eA/v2ZzK0cXNoeUAQ8SeQnhMjGFCoJctwpWiHlPQ8WJxBaoSgahu7CErSUX+nVOIWy2I8o4JLMBh0+3OLe0XwSCAwCQVl0M1LKLxWtJXpB1F2kCQdVhIWFxUyZMsX9nG1eJyYmIjQ01Db+qCPanRLkcMivdT1UwlEUwiX8IiFuHLQOSxYtspWvdoyfe4DonzsNa71u27btBrE41tdWa/bSG92Z2Y2io6MddpW4uDj7uRYRIWncWCzRcRxlHJLQjfcrQ1o6rp6/aCt/bRk/twjRPzcYFnspeXV9Zrk1vfOedK3otSGiusi1PgGbLkqxE4HRo4FC0jcuUk46V2pJe3xBUAc0vbYMR8o1wpG52+zkLX0hAdsR8CSRa/fLx6LaNz7OdgToENC3LzBpEhARAUdQEEIiKuK2Lx/Ein/MR+jVMyjVuSE2PvYvZ5udvEiABExHwJNE3lysfkC0nejmFO0sR4qdCGgyj4939onLUZN78xEdcHXtVmwu2gYNPnoCWyr3xOXDp+3kNX0hAVsQ8CSR6yaQ2iqvI1ovRefIkRIABMJjSyP21Gz82HYcqsfPxcXKtRH/AcMfAKGnixYi4Ekit5A7NNUXBELyB6PbosHY8OF6/BpUGpF/74JfWj4OR0KiL27HOkmABLJIgIk8i8ACuXiz/6uNUgfWYVrloai5/N84UaYuzkxfGshI6DsJmIIAE7kpwmAdI0pXyo87972FHwYvwZUrQIm72mBv16edE4is4wYtJQFbEWAit1U4c8cZGdiCnuNaImnDVnwX9gSqzn4Hp8rURsJMXZKHQgIkkNsEmMhzm7iN7hdd/yb0PPY+PnsoDgmJQSh8Rzsc7iSTe3/7zUZe0hUSMD8BJnLzx8jUFobIVLF+n7VBwvKt+KTUs6gw72P8Vr4GLn42lePOTR05GmcnAkzkdoqmH32p37wQ+h59Ex/2W4P9l8vhpr/djZMNuxpj0/1oFm9NAgFBgIk8IMKcO07mzy9rHH8qa/6sXoM3y41H6IYluHpzDSQMHwNcvZo7RvAuJBCABJjIAzDovnY5pnFeDD70ND4ZuhNzkzui8Osv4veI2nDMm+/rW7N+EghIAkzkARl23zutfedPvlUJ1XdMx7O15uLXUw4EdeqIhNtk8cwDB3xvAO9AAgFEgIk8gILtD1erVwfe2NIRy/+1DSMLvI6gn39CUvQtuPbsi0CC7lNCIQESyCkBJvKcEuT1mRIIlm/Z3wbkx4BDw/BSr9346vo9CHl7DK5ERMPxkSysef16pnWwAAmQQPoEmMjTZ8MzXiZQujQwYWoFRC75AvdFrcLG81EIeuxRXK3ZAFiwwMt3Y3UkEDgEmMgDJ9am8bRVK+Dz3U2wZuwKPFRwCo7tli6WDh2QdGtHWShZV0qmkAAJZIUAE3lWaLGs1wjow9DBQ4Lw5sHeGPPATgzBWFxYtA6oXx/Jfe8HDh702r1YEQnYnQATud0jbHL/ypQB/vNFfty7Zgj6NNyPMXgBf3w9DcnR1YBBg4CTJ03uAc0jAf8TYCL3fwxogRBo1AiYt7oYqnw7Bq3L78OkpIdxfeK/kFy5CvCijHA5d46cSIAE0iHARJ4OGH6c+wR0VcW77waW7KuAxLc+RMPQXfj6ck8kj/knkiOjgBEjgN9/z33DeEcSMDkBJnKTBygQzStQABg6VAayHLwZG4f8FzF5t+KHi7cCr77qTOivv84x6IH4xaDP6RJgIk8XDU/4m0DJksDYscCMfbUw88FpiA3agJ8SmgHDhzsT+j//CSRyuzl/x4n39z8BJnL/x4AWZEIgIgL49FPgi20NMKn7LDTCGiy8IJ3qw4YhOSIS0ITOWaKZUORpOxNgIrdzdG3mW40awPTpwIcbGuHdjnPQBKuwMMEtoWuXy4ULNvOa7pBA5gSYyDNnxBImI9BAJoL++CPwzuommHDbHKOF/nNCE2eXi7bQR43iQ1GTxYzm+JYAE7lv+bJ2HxJo3BiYPRv4YF0jTOw8G7FYh3kJLYCXX3Z2uegoF24758MIsGqzEGAiN0skaEe2CcTKXhY//AB8tCkWX/SaaTwUnZXQxhjlcj08EnjlFeD8+WzXzwtJwOwEmMjNHiHa5zGBevWAb74Bvt7dAD8+/D0a5t2EHxLaAyNHIqlSpJHQ8/KhqMc8WdA6BJjIrRMrWuohgapVZdr/f4CZh+th3QvT0Dx0C76/2MFI6LG970Pyy9JCZ5eLhzRZzAoEmMitECXamC0C5coBY2S70HnH6+D4hKnoVH4LZl29HcGjRuJKuUhcHv4aR7lkiywvMhsBJnKzRYT2eJ1A4cLAU08Bsw7XwXaZHfpQvS2Ye6UtCr7+DySWjsKpwTIO/eJFr9+XFZJAbhFgIs8t0ryP3wnkyQO0bn0Gn2+qg8iN32NEt/VYdq0JykwYhrMlbsaGhyfij8Q//G4nDSCBrBJgIs8qMZa3BQFZ9hyvzIxB7MnZmPz4cuwLjkbMpwNxsmg1fNt9MvbvTbaFn3QiMAgwkQdGnOllOgTCwoC+HzRHw8TFWDdqPq4VLo57frwfv0fH4vkGC/Dtt8DVq+lczI9JwCQEmMhNEgia4V8CwXmC0HB4B1Q5tx7n3puMysXO4Y1NHVDk3k5oW2YHnnySu9D5N0K8e0YEmMgzosNzgUcgOBglBt6HYid3I/nNt3HrTauw7EId3DJxIG6tfwZ16wLjxgGnTgUeGnpsXgKeJPJPxPxfRbeZ1w1aRgJeJpA/P4KffQYh8fuQ54kBGBD0IY4WrIr7f5+IZ5+5jgoVgM6dgcmTOeDFy+RZXTYIeJLIP5N6O2ajbl5CAtYnUKoU8P77CNq6FQWaxeDZQwNx8ZZYvH/fSuzYAdwv+0TrvqN9+8rwxlnAHxz0Yv2YW9ADTxL5UvGLGyZaMLg02YsEdA3dBQuAKVNQIOEMBnzZHAdvfQwrZ583kvi8eUC3bkDZssCjjwI//QRcu+bF+7MqEsiAgOyS6JFESilpb6BWBqX7yzlVhIWFxUyRL7wdJVF2pAkNDbWja4ZP9C/z0Oa5fBkRn3+OSlOn4lrRotg7aBCON2+LDRtLYNGi0lixohQuXcqLIkWuoUWLMzJ2/TTq1z+PkBBH5pXnsATjl0OAfry8bdu2G+T2sb40IVIq3+bpDaKjox12lbi4OLu6ZvhF/7IQ3o0bHY6YGIekZ4ejWzeH4/hx4+LLlx2OGTMcjr59HY7QUOfpokUdjvvvdzimT3c4EhOzcI8sFmX8sgjMRMUlv673NMemLudJ10rqa/ieBEhACeisotWrnRuLardLzZrAV1+hQH4H7rgD+O9/gdOnnZtg9OwJzJkD3HknoN3uev4TGUbA0S/8KnmDABO5NyiyjsAlkDcvMGSIc5B5tWrOp569egFnzxpMChQAunZ17jl68iSk60X6H6UDcvNm4JFHAF3Yq4lsbjR6NCDPU6VpH7go6Xn2CXiSyL+W6leJyrcUR0Xl60chARK4gYAm8eXLnRtB6z50deoAcXE3FAkJAdq2lS3q3gHi44FNmwDdxChZVgN46SUYY9QrVXIm+u+/58KMN8DjmwwJeJLI+0gN0m6AfA1RUfRjUQoJkEBqAroq1/PPO7tb9IF4+/bGPqJpDV8JkmEGuhGG7EqHtWuB48flP5b8z9LWuW6OoV0wJUsCbdo4fxs2bnQm/NS35HsSUAKeJHKSIgESyAoB3R16gwxA+NvfgNdfB9q1A06cyLAG7WJ5+GHgu++AM2ecjfmhQ517SA8bBsTEOMer33uvM+Fri55CAi4CTOQuEjySgDcJaItcm9g69VOb0/pgdOlSj+6QL5+zJa6bYmj3i/4GfPmlcybpkiXOcepRUUCVKsD//Z+zBa/975TAJcBEHrixp+e5QeC++4A1awAZb260zHWhliw+0dRJRjqDVIauG10w22QgsPaz16rlTOJ9pPNTW/Q6Z2n8+KpG10wmfwDkhue8Ry4SkEfuFBIgAZ8S0Iy7bp2zq+WZZ2RGhmTiDz8EtOmdRdG+dR3lqKorMiYlOVvt+lxVdcGCMpg501npzTcDLVs6tUULQN/r9RT7EWAit19M6ZEZCRQpAshMUGOYysiRwP79wLRpzkHlObBXRz82bOjU554DFi5cIY3/1tAuGO3J+eEH59BHvUXp0kDz5k5t1gzQrnxZG4xiAwJM5DYIIl2wCAFZItdI5NWrO1vnjRsDc+cC0dFecyBPHgdiY2Wet6g2/nVo465dwLJlkKUDgJUrAR3aqKJ/EGjXfdOmgJqiGhnJVrsBx2L/MJFbLGA01wYEtFNbn1Z27+5sHuuUT21W+0D0t0P7zlX1waiKPhhdtco5SlKP2sszYYLznLba1ZRGjf5s6etMVIq5CTCRmzs+tM6uBHTAuDaRO3RwzhLSZvJtt+WKt/rwVJcMUFXRVRp/+cX5TFafy+q4dv1tcT2TDQ93tvB1CKSqdsnoFnkU8xBgIjdPLGhJoBGoWtXZ19GpE9Cli3OM4T335DoFnXGqyVn18cedt79wwTlqcv16WclJVIfFT5/+p2kVZWqgdsvopCbXkd0yf/LJ7VdM5LlNnPcjAXcCOm5Qn0xqN4sOVdTmsY419LPos9k2bZzqMuW335xrxOiweE3sOsZ99uw/Z5zqNboygW6Hp1q7tnOIpI1XfXah8fuRidzvIaABAU9Ax5hrX4Ym8wcfdI4p7NfPdFiKFftrcr90ydkto4uAbdniVB3vLsv2/08qV3YmdB2FqapDJ/X5ri4oRvEOASZy73BkLSSQMwI33eTcK65HD+dcfR0grlsNmVwKFfpzxIvLVB0pc+iQczVHXdFx+3ZnstffKnVLRR/C6rh214NYPd5yi6zMVw1QFJSsEWAizxovliYB3xEoWNA58FtXzHrsMWeT1QTdLFl1WJO0DspR1XXXXaL7me7eDWOvU03uqjt3On+/XAley+oKkDpCU1UTu0t1w2tK2gSYyNPmwk9JwD8EtL9Bnyrqw0/tXilc+MZs6B+rvHJXHbeu/eaq7qIJft8+Z1LXMe+qmuA//fTGLhpFU758rPFwVZ8Tu1Rb9joSJ5BnrTKRu3+j+JoEzEBAM9aMGc7hiHff7ew/1yVxbSqa4F1dLO4u6vBHXTNGW/F79jh11aorMlQy1Jix6t6K1y4eTei6kJiq9su7jjp8MhurIbibYvrXTOSmDxENDEgC2hLXTuXWrZ0tct1aSGfpBJBoC7t8eafqhhwqixdvk9E0bYy+9sOHgb17na15bdGramtesV296iyv/2pXjw6XdHX3uI6RkTBmsuo9dCl5KwsTuZWjR9vtTaBECeCnn5yzP7t1c07F1CxEga4xo61u1dtvvxGIPmzVjTp0OZuDB52qr+PjnTj1nLtoXdovHxHxp2or3qV6Tlv8ZhYmcjNHh7aRgI4z1yamLoii/eY6G7R4cXLJgICrBa6tcP2DJrVcueIcVaMjazS5q+pr1YULnT8C+mPgLrpbkyZ0Va1X1fVaH8Lqe38meyZy92jxNQmYkYAO33D1md91FzBvnv07fX0YB30E4RoJk9ZtdE6Wttq160aT+5EjTtX3qvpbeu7cX6/Ucfaa1FW1u8Z1dHUP6W+yPpTVmbTeFiZybxNlfSTgCwLatPzkE+CBB5yrX+nrQB6m4QvGKXVqonV1s+h67mmJToQ6ehQ4dsx5dL3W96o6tFIXJ7t+/a9X6yJkmtRdqsldX+dEmMhzQo/XkkBuEtAx5fp0T9cz19WrBg7MzbvzXm4EtBtFZ6dmtAKxJvHTp52JXUffqGpL3/VajzrMUhO+/hWQE2Eizwk9XksCuU3glVeci5wMHuxc0CS9JmNu28X7/YWAjoTR1rZqRqL98efP52yPkeCMbsBzJEACJiOgT/J0J2YdrtGrl/PvepOZSHOyRkBDqg9TcyJM5Dmhx2tJwB8EdJEtXb9cO2r14af7oGl/2MN7+p0AE7nfQ0ADSCAbBHQqpC4zqLtADBuWjQp4iZ0IMJHbKZr0JbAI6OJa+sBz/HjnwuCB5T29dSPARO4Ggy9JwHIE3nrL+dCzXz/n8AjLOUCDvUGAidwbFFkHCfiLgM5u+fZbZ3+5Dk9Ma+Cyv2zjfXONABN5rqHmjUjARwR0muLEibqiFMK//tpHN2G1ZibARG7m6NA2EvCUwEMPAbLkbaQ+ANU91ygBRYCJPKDCTWdtS0Cn60urPEmXv9X+ct2tgRIwBJjIAybUdNT2BGQRj91Dhji3uh892vbu0sE/CTCR/8mCr0jA8gTOtmjhXFhLE/mGDZb3hw54RoCJ3DNOLEUC1iHwzjtAmTKA9puzi8U6ccuBpUzkOYDHS0nAlAR044l//9u5lurYsaY0kUZ5l4Cnibyj3Ha36D7RF7xrAmsjARLwOoGuXZ3rsOiStwcOeL16VmguAp4kclmMETJIFZ1EZYEH9Ek5yoFCAiRgWgLaxaK7JDzxBKBb0lNsS8CTRN5IvNeWuP6s65imb0TvEKWQAAmYmYDuNaYPPefPB6ZMMbOltC2HBGTwaaYiix5Du1YeTSkpe02hsWjq7Un6y2eqKrVEtxmv7PePbNSEM/Zz638e0b//obDkC8bPkmEzjJYpupCJAL6RXlLtR25VayJ/3+19Wi/Xp/WhTT6zs28aIvpn7S8q42fd+GU7dp50rRwTLpXc2FSU1/oZhQRIgARIwAQEPEnk68TOqqJRovlE7xWdKUohARIgARIwAQEdkZKZyNagkK27MVl0kOh/RaeJZiYbMitg4fN29k3DQv8s/OVk/CwdPLv/37N0cGg8CZAACZAACZAACZAACZAACZAACZAACaQikNm0/fxSXvagMiYTrZFjpKiVJDP/+okzp0U3p6hrnL28Nb18Ihb+KrotHUt1fsG7ojoRbKtoA1ErSWb+tRFnfhd1xe4fVnJObNVRZHGiO0S3iz4lmlqsGkNPfGsjzlo1fgXE9rWiW0Q1diNEU0uu5U59SLpftLKojmRRo2qIussT8ubDlA90pIsmdauIJ/71E2cyG09vVn9biWGanNNL5J3l3FxRTQZNRPWH2EqSmX9txJlZVnIola3l5L3rx7WwvN4jmvr/n1Vj6IlvbcRfq8ZP/0+FiqqEiOr/Lf0/5i5Zzp2eDD90v4HrtSfT9nUa/+cpF3wnx/ai6oQVxBP/rOBHejYulRPn0jspn2vsvhDVBTpWixYT1f9gVpHM/LOKH+nZeUJObEw5mSDHnaIyH/8GsWoMPfHtBkct9kb/TyWm2KyJXDX1QjhZzp3ZTeT6pTmSYowejoqm/iK5l0mS8/qnUElRK4i77WpvWv7p53eJateD/lC5T5qSt5YWT/23spNNxXj9S1L/8qhpYUcixfb6oqn/arJDDNPzTdyFleOnf/Frt96vogtEM4qdR7kzu4lc7h3w8qMQiBStI6rBcP31IS8pJiegrdkI0bqi74nOELWi6J/o00SfFr1gRQcysDkj36wev+vidz1RnSWvf/3XEs2RZDeRH5O7urdA1SD9zF3cy+SVE0VFz7oXMPFrd9vVzLT8U1+upvjwkRxjUl7b4eCJ/1b2U5Oe68/bOfI6RLSUxRxSmzWJTxadnobtVo5hZr7ZIX4ast9E40R1YIW7uMfOp7lTKz8gGiXqetiZ+s/Tv8s594edU+S9VcQT/9z7jHuKY9qXbCWJFGO3pWNwF/nc/WHn2nTKmfnjSDEuPf/KyjnX8xptER12ey8vTS9q+xeiEzKw1Kox9MQ3K8cvTGJWLCVuBeW4TLRrynvXIVdzZ2e5qz4t19Erw1MsGCnH7imvdZjNVNF9opoIKotaSTLzb4w4o8OHtJ9Vf1Wri1pFvhZD9aHSNdGjoo+IDkhRORhJbqIcNba/iMaKWkky82+gOOOKnf4AN7OSc2JrC1F9QKbPZzanqH5f7RBDT3yzcvzqSJw2iWrstKHhGvpqp9wpblFIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgARIgAR8SeD/AZ8IeEHSk87QAAAAAElFTkSuQmCC"}},"cell_type":"markdown","metadata":{},"source":["# 1 簡単な行列計算(25点)\n","\n","関数\n","$$\n","f(x) = \\frac{4}{1+x^2}\n","$$\n","を多項式\n","$$\n","F(x) = a_0 + a_1 x + a_2 x^2 + a_3 x^3 + a_4 x^4\n","$$\n","で補間することを試みる．\n","与関数のx=[0, 0.25, 0.5, 0.75, 1.0]\n","での値から直接逆行列から多項式補間する手法を試す．\n","連立方程式の係数行列$\\mathbf{A}$(ヴァンデルモンド行列と呼ばれる)およびデータベクトル$\\mathbf{y}$は\n","```\n","A = [[1.         0.         0.         0.         0.        ]\n"," [1.         0.25       0.0625     0.015625   0.00390625]\n"," [1.         0.5        0.25       0.125      0.0625    ]\n"," [1.         0.75       0.5625     0.421875   0.31640625]\n"," [1.         1.         1.         1.         1.        ]]\n","y = [4.         3.76470588 3.2        2.56       2.        ]\n","```\n","となる．ヴァンデルモンド行列$\\mathbf{A}$の逆行列をデータベクトル$\\mathbf{y}$に掛けることで，\n","係数の値を求めよ．\n","\n","以下は$\\mathbf{A}$，$\\mathbf{y}$を求めるコードと，与関数，補間関数のプロットである．\n","```python\n","import scipy.linalg as linalg   # SciPy Linear Algebra Library\n","import numpy as np\n","\n","nn = 5 # x_i number\n","\n","def func(x):\n","    return 4.0/(1+x**2)\n","xx = []\n","yy = []\n","for x in np.linspace(0,1,nn,endpoint=True):\n","    xx.append(x)\n","    yy.append(func(x))\n","\n","print(xx)\n","print(yy)\n","\n","a_matrix = []\n","y_vector = []\n","for i in range(nn):\n","    for j in range(nn):\n","        a_matrix.append(xx[i]**j)\n","    y_vector.append(yy[i])\n","        \n","A = np.array(a_matrix).reshape(nn,nn)\n","y = np.array(y_vector)\n","print(A)\n","print(y)\n","\n","inv_A = linalg.inv(A)\n","print(np.dot(inv_A,y))\n","```\n","\n","![image.png](attachment:image.png)"]},{"attachments":{"image.png":{"image/png":"iVBORw0KGgoAAAANSUhEUgAAAYoAAAEGCAYAAAB7DNKzAAAAAXNSR0IArs4c6QAAAERlWElmTU0AKgAAAAgAAYdpAAQAAAABAAAAGgAAAAAAA6ABAAMAAAABAAEAAKACAAQAAAABAAABiqADAAQAAAABAAABBgAAAABgCvtYAAA5sUlEQVR4Ae2dB3hUZfbGX3qREgihaOgBhAQIBEKHIChlV3QVV1DBSlFYRd1VFAurICp/K6KCgGBlsSOKIJjQe4cgCIhIpJdgQGrmf85MRidhmEwyd2buvXnP87y5/Su/7zKH7/vuPRegkQAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkICVCBSyUmH9LWtkZKSjVq1a/p6e7byTJ0/isssuy7bPqhusi/lajm1ivjbREtmlXQKtx5o1aw4LjihztpJxpbpWkpoYExPjyK8lJyfn91LTXce6mK5JHGwT87WJlsgu7RJoPeT3c7W3n+PC3nZaeN/XUvaB5cuXt3AVWHQSIAESMBcBuzkKc9FlaUiABEjABgToKGzQiKwCCZAACQSTQNFgJs60SYAESMBMBM6dO4e9e/fi9OnT2Yqlw9Vbt27Nts+KG/7Wo2TJkoiOjkaxYsX8qiYdhV+YeBIJkIAdCKiTKFu2LPSpyEKF/nro8/fff3fut3od/amHzN3jyJEjTodZu3Ztv6pshaEnfVZ1mugd0a1+1YonkQAJkIAXAtqTkMfnszkJL6fZepc6SGWQs1flq9LhchRTpFAHRZtzFK67bG8T7RANzzp2gyw/FQ0Q9craxwUJkAAJ5IuAZ08iXwnY4KK8MgiXo5gqrNUpeFoR2Rgv6iFqJOqbtYyW5a8itQuuRXD+zhm9GinPHELm+czgZMBUSYAESMCCBMI1R7FQWNXKwStRtrUnsStr/3RZXifaK1JnsV7ky7ENlOMq59hbSkqKrubJprx1FjPSbsL3Easx7Km9iEqMyNP1Zjs5IyMD+eFgtnpoeexSF7vUw6ptopO9Oo6f0y5cuOB1f87zzL6dl3ro0JMVfh9qCXTPoafesj3JoyH6yfobIp2jeFf0lsivOYqEhIR8vf6ZeSHTMbrHDEdEoWOOkjjleKFHsuPcH+fylZYZLgr0LU0z1MFdBrvUxS710HaxYl1SU1Pdt1S25YkTJ7JtW3XjUvWQp70uqpI3FvIba9k3s09K4e8U3Sv6UOTLnCE80tPTfZ1zyWOFChdC20eikLr2DHpU24BHZyehVcWfsGGGTpvQSIAESCBwArt378aVV16JW2+9FQ0bNkTv3r1x6tQpzJ8/H82aNUPjxo1x11134cyZM1i1ahVuuEGnaYGvvvoKpUqVwtmzZ50T0XXq1HHu37lzJ7p37w75DzK6deuGH3/80bn/jjvuwODBg9GqVSs88sgjzn35/ROuoSdv5U2TndU9Duhwk+4LuVWLr4LP9lbGZ/9ZhiGvxKDFzREYPi4FT3zTBiXKlQh5eZghCZBAEAgMGyYD2jqiDZSSoScU0WnSAC0+Hnj11VwT2bZtGyZPnox27do5ncLLL7+MCRMmOJ1F/fr10b9/f7z11lsYOnSoFNFVxkWLFiEuLs7pPM6fP+90AJrRwIED8fbbb6NevXr44YcfcN999zmXekwfB166dKlULbC6+Rrz13xCaasks3oifbC3uKiPaKYoL2ZYrCftXfR+qQ22/lQUt9ZdgVGLkxAfJdAnbMpLeXguCZAACVxEoHr16k4noQduu+02p4PQdxrUSajdfvvtWLhwIYoWLYq6des6XwZcuXIlHnroIed+dRodOnRwzt2pI7jpppsQL05qmDi/ffv2OdPQP7o/UCeh6YSrR/Gx5J0kqiTSyeqnRZNFQ0VzROr+9BHaLaK8mA49XZvfoSdvGVWsWwFTd7RHX3kiauDTVdF+cG38660FGP1dAspULePtEu4jARKwAgGP//n/EeIX7nI+nhoREeF8Cc4bto4dO2L27NnOt6i7du0KHVLSSeuxY8ciMzMTeq2715HzhTujPpkQrh6FPvpaTVRMpENM6iTUvhWpS60rGi3KqxnWo8iZcbcRLbB5T3kMabwI4zZ0QFz0ccwdsybnadwmARIggVwJ7NmzB8uWLXOe99FHH6FFixbQuYsdO3Y4973//vvo1KmTc117Dq+KU2vTpg2ioqKcDkWHrnQYqly5ctCeyCeffOI8V2assWHDBue6kX/C5SiMrINnWgFNZnsm5G297OVlMW5jJywcvxkli5xFt8cTcGe9RTj283Fvp3MfCZAACXgl0KBBA4wfP945mX3s2DE8+OCDePfdd51DRTqZXbhwYedEtF6sk9EHDhyA9izUmjRp4pzwdvdKPvzwQ+d8R9OmTZGYmOic9HaeaOCfcA09GViFbElpj+JreVZ6QLa9Bm+0v68J1t9yGs/0TMGLy9rju5gjGP/wctzwYmuDc2JyJEACdiSgcw8ffPBBtqp16dIF69aty7ZPN/RJJ30Cym0TJ050rzqX2qP47rvvnOueQ09Tp07Ndl4gG3brUQTCIk/XlowoieeWJmHVRztQtcQx3Di2NXpHL8P+jRqZhEYCJEAC9iFgN0cR1KEnb83erO+VWHm4Lp67JgWz0pqhUXwxTBuwGI5Mh7fTuY8ESKCAE9DItZs3e75rbH4gdnMUQZvM9tWUxUoXw2NzkrD+m9/QqMyvuGNSe/SovAa/LNnr6zIeIwESIAFLELCbowgr9Ct71sHCo3EY13sBFh+5ErHtI/DGTQsYZDCsrcLMSYAEAiVgN0cR8qGnnA1QuGhhDP2kE7YsPo72kT/iX592QseKm7FttjvWYc4ruE0CJEAC5iZgN0cRlqEnb01cs100Zh9MwNR7FiM1ozqa9rwcY7ql4Nypc95O5z4SIAESMC0BuzkKU4HWMCC3v9MeWzeex7VXrMPjcyXIYKUdWPfxj6YqJwtDAiRAAr4I0FH4omPQsSpxUfhkbxsJMrgc+85URMtbYvB42xScPn7aoByYDAmQgBUJ6JvUGobDSNOAgUab3RxF2OcofDWQvpCXuqME+tdbhjHLJMhgld+w+M2Nvi7hMRIgAZsR0FAd+ma2RojVMBzPPvssWrZs6Xzj+umnNeydy9577z3nPn3jul+/fs6dGufp008/dZ+CMmVc8eb0A0Qa6uPmm29Go0b6gVBjraixyYU9tZC8mR1ILSvUjsCU7R3Q9/k1GPhkZXQY0gRD3l6AMd81h4YIoZEACYSGgEeUcQmyVyqUUcbx008/Ydq0aZAPDTl/+DUyrPYuevXq5YwOGxkZiVGjRjlDhFeqVAlHjx7NFcratWuxfPlyZ3iPXE/O4wl261HksfrhO/3q4QnY9GsF3N90Ad7cJEEGa6RDv9lNIwESsD+BmjVronXr1pg7d65T+sGi5s2bOz86pE5EvyuhIcLVSahVrFgxVyga50lf5guG2a1HEQxGQUtTw5S/tr4TbpZvXNx9f2l0f6IF+k9ZjFfmxkLDm9NIgASCR8Ajyrh8L/sPlC0buh69O/y39iIee+wxDBo0KFtFx40bl23bvaExotxzGrrUr925zZ2me9vIJXsURtLMZ1ptBzXGugNXYES7FHy0qxUa1juPTx92hSDOZ5K8jARIwAIE9NOlU6ZMcX6ASIublpaGgwcP4qqrrnKGDj9y5IizFu6hJ+0xrFmzxrlv5syZkG9hO9eD/cdujuJaATbRyA8XBbsB3OlrkEH9it6q6bsQXfIwbnq5DW68Qp6SWn/AfQqXJEACNiNwzTXX4JZbbnF+a0LDi+v3szUCbGxsLEaMGOH8JoVOZuuX7dQGDBiABQsWQPfp9yyC2YuwGeqLqyMfGZceXf4sOTk5fxcaeNW5P845nu+e7CiBPxwRhY45pty50JF5ITPPOZihLnku9CUusEtd7FIPbSYr1iU1NdXrHSaTyl73W21nXurhjYX8mnqdKLVbj+Jir2HBPUVLFsWjs5Ow8bt9aFz2F9z1bgd0i1qL3YsZZNCCzckik4DlCdBRmLgJ63erjZQjjfFmn4VYdrQ+4jpE4PUbF+DC2QsmLjWLRgIkYDcCdBQmb1ENMnjvxx2xZekJdIzaigc+74QOkanYOmunyUvO4pGAOQnIcJI5CxbCUuWVAR1FCBsnkKxqtLkC3+xvgfcGLca2k1cg/tpojL6aQQYDYcprCx6BkiVLQp8kyusPpZ1Iad2VgbLw1/gehb+kTHCeBhns93Z7dPvXIdzfYw2emJeEGZHbMOWdTCTc1tAEJWQRSMDcBKKjo7F3714cOnQoW0FPnz6dpx/ObBebaMPfeqiTUBb+mhUcRR2pzAhReVFvfytm5/Mqx0Zh+p4o9H18Be59oSYS+0Xh3+NSMHJ2K5SqWMrOVWfdSCAgAsWKFUPt2rUvSkNjJenb0Va3YNUj2ENPUwT8QVHOD8R2l33bRDtEw0W+TL/4c7evEwrqseuea4XUXaVwV4MleHFlEppW3Y8Fr60vqDhYbxIggSARCLajmCrlVqfgaUVkY7yoh0jDHPbNWjaW5awcqizbNB8EImqWxzs/dsS8F9figqMwkobF497YhTix94SPq3iIBEiABPwnUMj/U/N9Zi25Uh1AXFYKbWQ5UtQta/uxrOWYrOWlFp/KAV9DTwPluApVqlRJmD59uq7m2TIyMv4M3Zvni8N8wZljZ/HZ8DOYvL0nqhXej4f+mYzmg/wfhwxz8X1mb+V28ayYXeqhdWJdPFvWHOuBtknnzp01PkiLcNSmlmS62SNj/bGf5LGtgdbf8NjOuRopO94W6fOgbqeS8xz3tjOER0xMjEzs58+s+LZpzpoun7TJEVtiuzwD6HDcWmux49CPh3OeYrltO7SLQrdLPVgXc/4TCvT+kh9Sy76ZfUQKP1hUV5Rbr8M038yWsobNWt0dh7WHa+KBpl9ixu6WaNjQgen3L4Ujk8+Ph61RmDEJWJhAsOcovKFJk53VPQ7o2IjuM8KcPQorBgU0ovKeaRQvUxzXvxqBNZ/uRu3SB9B3XFtcf8VKpK3e53ka10mABEggVwLhcBSrpFT1RPqMWnFRH9FMkRHGHkUOio1vrC/hP67E//09Bd/vb4xGLUvjnf4L2bvIwYmbJEAClyYQbEfxsWS9TNRApBHt9DHX86KhojmiraIZoi0iI4w9Ci8UixQvgoe/liCD8w6hecQuDHy/I7pErsfOH37xcjZ3kQAJkEB2AsF2FH0lu2qiYiIdYposUvtWVF+k8w6jRUYZexQ+SMZ0qYn5h5piwq0LseZ4HTTuEoWXr0thkEEfzHiIBEgACLajCDVj9ihyIa5BBgd+IEEGV55Cl8qb8fDMJLStuBWbv/gplyt5mARIoKASsJujYI/Czzs5umU1zNzXEh8NXYpdp6qi+Q018d/OKTib8dc3eP1MiqeRAAnYnIDdHAV7FHm4YTXIoD4NlSozRDfVXIWRKUlIqPQLVk1LzUMqPJUESMDuBOzmKNijyMcdG9WwEj7c3Q5fP7kSx86VQes7GuDfLVJw6vCpfKTGS0iABOxGwG6Owm7tE9L6/P2ZRGz5uTQGNFyCl9YkoUm1Q0h5lUEGQ9oIzIwETEjAbo6CQ08B3mTla5TH26kdkfyKy0F0fjAegxouRPqe9ABT5uUkQAJWJWA3R8GhJ4PuRI1Cu3GffOdChqAm/dgOjWqfcg5NGZQ8kyEBErAQAbs5CguhN39RS1cqjbGrkrB86jZEFv8dvUYl4pZaS3Bo62HzF54lJAESMIwAHYVhKO2bUMvbG2H1oVrOx2c//UWCDMYWwkdDljAMiH2bnDUjgWwE7OYoOEeRrXmN29Agg0/9kIR1X+5BTOl9uPXNduhVbRX2rtpnXCZMiQRIwJQE7OYoOEcR5Nss9roYLDna0Bn6Y/7BODRKvMwZEiTzfGaQc2byJEAC4SJgN0cRLo4FKl8NMvjgl0nYnHwYLSvswOCPJMhg1AbsmM8ggwXqRmBlCwwBOooC09TGV7ROUg3MO9xMwpYvwloNMti1sjOc+fnTGiCYRgIkYBcCdBR2ackw1UPDgNwzrQNSV53CNVU34j/fSJDBStuw6bPtYSoRsyUBEjCaAB2F0UQLaHpXtKiGL9MSnZ9c3X2qMpr3ro2nO6XgzIkzBZQIq00C9iFgN0fBp57CeG9q7+Lm19pi67Yi6FN7BZ5ZmITmlX/F8kmbw1gqZk0CJBAoAbs5Cj71FOgdYcD1kfUq4v1d7fHNyFU4cb402g5ohIcSUnDy4EkDUmcSJEACoSZgN0cRan7MzweBnk+3xJbdZTA4djFeWZuExlccwfyxa31cwUMkQAJmJEBHYcZWsVGZykWXw5ubO2LB6xtQtNAFdH2kOQZcuRDHf2GQQRs1M6ticwJ0FDZvYLNUr+O/mmLD/qp4JDEFU7ZJkME6f+Crx1eYpXgsBwmQgA8CdBQ+4PCQsQRKVSyFF1YkYcW0bYgqno7rx7TCzTWW4sDmQ8ZmxNRIgAQMJWAVR3G91Pod0f9E1xhKgImFnECL/hJk8EgdPNslBV/+moBGTYrgg3sZZDDkDcEMScBPAqFwFFOkLAdFOZ+R7C77tol2iIaLfNmXcnCAaLDoZl8n8pg1CBQrXQxPzEvC+q/3osFlaej3djv8repq7FmWZo0KsJQkUIAIhMJRTBWe6hQ8rYhsjBf1EDUS9c1aNpblrByqLNtue0JW9DqaTQg0/HtdLDrSCK/dsAALDjVCbFuZ/O6zAAwyaJMGZjVsQaBQiGpRS/JRBxCXlV8bWY4UdcvafixrOSZrmXOh5Xxe9L1oXs6DWdsDZalClSpVEqZPn561O2+LjIwMlClTJm8XmfRsq9Xl2Kbf8foTUUg50RptLluDB5/6FVGJEU66VqvLpW4Ju9RD68e6XKqVw7c/0Dbp3LnzGil9i3DVoJZk7Dn01Fu2J3kUpp+sv+GxnXP1ftmhFXhbpMNPPi0hIcGRX0tOTs7vpaa7zop1ybyQ6Xj37kWOiELHHCXwh+P57smOc3+cc1ixLt5uCLvUQ+vGunhr4fDuC7RN5Id1tbcf11AMPXnLN6/7XpcLEkTqJNRZXMoYwuNSZCyyX8OA3DGpPbauP4u/Xb4ew79LQquKPyHth6MWqQGLSQL2IxAuR6EzltU9cEbLOmcxPYAU9NWqTSrjs7TW+PTfy5B2JhK3P9sLI9ql4PTx0wUdDetPAiEnEC5HsUpqWk9UW1Rc1Ec0UxSoMdZToARNdv2NY9sgdXsx3Hh5Mp5bmoRmVdKwdMImk5WSxSEBexMIhaP4WBAuEzUQ7RXdLTovGiqaI9oqmiHaIgrUOPQUKEETXl+xbgXc+2ExfDdqNf64UALtB8fi/qYLkLE/w4SlZZFIwH4EQuEo9NHXaqJiIh1imixS+1ZUX1RXNFpkhLFHYQRFk6bRbUQLbN4bgaFNFuGNjR0QF30cc8foMw40EiCBYBIIhaMIZvlzps0eRU4iNtsuU7UMXt/QCYve3IySRc6i2+MJuLPeIhz7+bjNasrqkIB5CNjNUbBHYZ57K6glaXdvE6w/cDkeb5uC93e0QaOYM/j8keVBzZOJk0BBJWA3R1FQ27FA1rtkREmMXpKE1dN3olqJo7hxbGv0jl6G/Rs1YgyNBEjAKAJ2cxQcejLqzrBQOvE3N8CKwzEY0y0Fs9KaoVF8MUy9ZzEcmQ4L1YJFJQHzErCbo+DQk3nvtaCWTIMM6st5G779DbFl9+DOye3RvfIa7F6sD9rRSIAEAiFgN0cRCAteawMCDXrUwYIjjfHGTQuw9EgDxHWIwLjeDDJog6ZlFcJIwG6OgkNPYbyZzJJ14aKFMWRGJ2xenI72kT/i/s86oWPFzfjx211mKSLLQQKWImA3R8GhJ0vdfsEtbM120Zh9MAHTBixGakZ1NP3bFXjumhScO3UuuBkzdRKwGQG7OQqbNQ+rEygBDTLYf6IEGdx4HtdFr8WI75OQWGkn1n6oAQFoJEAC/hCgo/CHEs+xPIEqcVGY8Wsb57sW+89UQOJt9fBYmxT8cYxBBi3fuKxA0AnYzVFwjiLot4y1M/jHC62RuqMEbq+3DM8vT0J8lX1Y/OZGa1eKpSeBIBOwm6PgHEWQbxg7JF+hdgQmb++A719Yi7OOougwpInEj1qA33/73Q7VYx1IwHACdnMUhgNigvYl0PWR5tj0awUMi0/Bm5s6ILbGCcx+1usHvuwLgTUjAT8I0FH4AYmn2JeABhl8ZV0SlkzYgjJFTqPnUy3Qv+4SHPmJX9Szb6uzZnklQEeRV2I835YE2gxsjHWHovFk+2R8vCsRDRtkYsaDyxgGxJatzUrllYDdHAUns/N6B/D8PwmUKFcCzyzqjDUzdqFGqYO4+dU2+McVK/Hb2v1/nsMVEiiIBOzmKDiZXRDvYoPr3OSmBlh+pD5e7JmCOfuboFFCKUy+cxF7FwZzZnLWIZAXR9FWqnWLqL+HrFNTlpQE8kCgaMmi+M83Sdg4Zz+alv8Z90ztgKuj1uHnhb/mIRWeSgL2IOCvo3hfqvt/ovailllqIUsaCdiaQL1raiP5cBO81XchVh6NQVyninjthgW4cPaCrevNypGAJ4Ginhs+1tUpNBIxwL8PSDxkTwIaZHDwRx3xtwd+w6Br0zDsi074X+QmTPqwFBr1irFnpVkrEvAg4G+PYrNcU9XjOq6SQIEjUL3V5fhmfwt8cO8SbD95OZpdVx2jujLIYIG7EQpghf11FJWETapojmimh2SVRgIFh4AGGbz1zXZI3ZSJf9RYgyfnJ6FF5C6seV//edBIwJ4E/B16GhnG6jeUvB8QqbOaL3pLRCOBsBKoHBuF6b9Eoe+Ilbj3+Rpo1b8SHn4jBSNnt0KpiqXCWjZmTgJGE/C3R7FAMvam3MozRU7QL93r0JWndZeNbaIdouGeB7ysazzowaJ/itp5Oc5dJBA2AteNTkTqrlK4s8EyvLgyCU2r7seC19aHrTzMmASCQcBfR3FHPjOfKtepU/C0IrIxXtRDpBPkfbOWjWU5K4cqy7ZaL9E3om91g0YCZiIQUbM83vmxA+aNXYcLjsJIGhaPwbELkf7rCTMVk2UhgXwTKOTHlU/JOfVFt/lxrrdTaslOdQBxWQfbyHKkqFvW9mNZyzFZS18LdRZ/u8QJA2W/ClWqVEmYPn36JU7zvTsjIwNlypTxfZJFjrIuoW+os8fO4LPhZzFpe09ULXwAw29bgsZ3Rv1ZELbJnyhMtWKXdgm0Hp07d14jDdMir40zUS74WORvz8Nb+rVkp+fQU2/ZnuRxYj9Zf8NjO+dqkux4XTRBNETky5whPGJiYhz5teTk5PxearrrWJfwNcmKyZsccSW2yePkDkffmkscB7cedhaGbRK+NvGVs13aJdB6yI/ram8/sLk5AB0WelaU6e3iEO1LkXzuFw0S6ZCVL2MID190eCxkBBLvisOaQzXx36RkfPpLCzSUQdYPhzLIYMgagBkZSiA3R6H/Q58hqmtgrmmSVnWP9KJlXfcZYc4eRXp6uhFpMQ0SCIhA8bIl8FRyZ6z7fDfqlU7DbePb4LkbS2Dvqn0BpcuLSSDUBHJzFClSoD6iDwws2CpJq56otqi4SNPXdzOMMPYojKDINAwlEPuP+lh8pBFe6fUDlh5vitjE0pjYj0EGDYXMxIJKIDdHoZnr/MKN+SyFzm8sEzUQ7RXdLTovGiqaI9JHX7XHskVkhLFHYQRFpmE4gSIlimLYV1fh/ZfnoUXETgz6oAO6RK7HruRfDM+LCdqLQHKyDOnImM6pU+Grlz+OQkv3Wz6LqHMc1UTFRDrENFmkpo+51hfpkNZokVHGHoVRJJlOUAhUbFYe8w7HY+JtC7H6eF00vqoSXrs+mUEGg0LbHok+9BCwa5e8eLYtfPXJzVE8klW0cbLUJ49yKnwl954zexTeuXCviQgUKlIYA97viNSVJ5FUOVV6Gp3RIXILfpy53USlZFFI4C8CuTkKHRpS00em1niRHjOTsUdhptZgWXwSiG5ZDbP2tcD79y3DtpPRiL+uBsZ0nY/zp876vI4HSSDUBHKL9aQ/vGrTXAvT/9UexbV86sn07cQCZhHQIIP6NNTVQ45gaPcNeHx+F3wamYop71xA09sakxMJmIJAbj0KdyF1PkFfvpsr+sFDsmoqY4/CVM3BwvhLoEqjSHyypxU+Hb4aaWcroUW/KzGi1TycPhrGGUx/C8/zbE/AX0fxiZBYJ3pC9B8PySqNBEjAKAI3jmmB1J0lcVv9VXhuZVc0qXoAKa/oPz0aCYSPgL+OQh9p1fDeK0WecxXhK7n3nHXoaSKHnrzD4V5rEKhYqxze3dYW817aICERCqPzQ81wz5WLcHT3CWtUgKW0HYHcHEVFqbFKh3Q0zpI+6urep0uzGYeezNYiLE++CXR5qCk2/haFRxOTMXVbGzSsewb/e2SNRI/Kd5K8kATyRSA3R6G9B33i6XbRv0VLsrZ1n4pGAiQQRAKlK5XG8ys6Y/UH21CjxAH0GZuAa6PXYs+6I0HMlUmTQHYCuTkKDbNRRyQhzZwRXjfIUr/KMk4UK6KRAAmEgED8rbFYfqQ+Xu4+F8m/NUCj5iXwWr/VuHCe3YsQ4C/wWeTmKNyA9PFY/STp6yJ1Euo4zPjILOcopGFo9iRQpFRxPDj7GmyZk4YO5Tdh2Act0DZqO7bMY5BBe7a4eWrlr6OIkyLfI5KoI04NkKXuM5txjsJsLcLyGE6g1jX18e3hRHzUbzZ2Ha+IhKsr4KXeS3HhXDi/BmB4NZmgiQj46yjWSplbe5S7laxzjsIDCFdJIJQEChUtgr7v9cCWFSfRI2o1/v1ZW3SO2oRd838OZTGYlxcC42TM5cknvRyw8C5/HUWC1HGpaHeWNCJsS9Em0UYRjQRIIAwEKifWwuf722HaXQuwIb0WmnSNwsQb58BxTp9op4WDwP3ymbVRo8KRc/DyzC2Ehzvn7u4VLkmABMxFQMOA9J/cCUmDD+Cunvsx6PNu+LLSEkz6rAIu76rTiTQSCIyAvz2KXyQbXwqsFLyaBEggYAI1WlbB3P1NMO6utUg50QxxV1fF9BtmAGfOBJw2EyjYBPx1FFahxKeerNJSLGdQCBQuUghDJzfH+uVn0KDSUfT94p+4ufIPODKHU4pBAV5AErWbo+BTTwXkxmU1fROo36oCFu2Lwejbt+GLE10Q1/0KzPjbNDgyTvq+kEdJwAsBuzkKL1XkLhIomASKygzk41MbYNXis6ha6Txu/vZ2tK+0FcvH6WfraSTgPwE6Cv9Z8UwSsCSBpu3KYPX+6pj0yHbsOl8Dbe5viT51VmL3xhOWrA8LHXoCdBShZ84cSSDkBIoUAe5+oT5+2lcWT7adj5k/x+HKpsXx6A3bkZ4e8uIwQ4sRoKOwWIOxuCQQCIEyUaXwzJIu2P7NDvSpMBdjv4hBTOUTGD/mBM6dCyRlXmtnAlZxFJdJI+hjG3+3c2OwbiQQKgLRPZtg6oEeWH3fu4g7txZDHy+HJjXTMetrB8OYh6oRLJRPsB3FFGFxULQ5BxN9gW+baIdoeI5j3jYflZ3yQDiNBEjAMALFiqH5+Lvxw+Yq+KrBI8jctx/X9iqEPr1O4vRpw3JhQjYgEGxHMVUY5XyrW0ZLMV7UQ6SvjfbNWjaW5awcqizbV4tSRepwaCRAAgYTKNSoIXptGYPNr8zD6GIjMWPWZejWOA3HjxaMIIMbJQhRRob+LNEuRaDQpQ4YuL+WpKUOIC4rzTayHCnqlrX9WNZyTNYy52K07NChJ3Uqf4j+IfJ2Bw+U/SpUqVIlYfr06bqaZ8vIyECZMmXyfJ0ZL2BdzNcqZm+Tkvv3Y9PjmzDk56dRr+RuPP/8RpRtWskrSLPXxWuhvezs3DkJ9esfw4QJ+rmdwE3TU0tOTnEuA/0zYEACduwoi4kTV6NevQyfyQXaJp07d14jGbTwmUmQDtaSdDd7pN1b1id5bPeT9Tc8ti+1eocc8GuOIiEhwZFfS05Ozu+lpruOdTFdkzgs0SaZmY7v//2dowxOOGrgF0fqQ+84HOfOXQTTEnW5qNQX74BzVubi/fndY3R68fH68VuHY+3a3EsUaJvIb+xqbz/AwR568pZnfvdNlQtn5XIxQ3jkAoiHSSBXAoUKoevYblg49wzOlCiL9i//A8ti7wE2GPM/7lzz5wmmIxAOR5EmFKp7kIiWdd1HIwESMBGBZldXwtLUCoisVhxdtr+Jr5s/DTzxBDjTbaJGClFRwuEoNH5APVFtUXFRH9FMkRHGWE9GUGQaJJBFoE4dYMmGsohrVhzXOz7HpNH7gWbN5Os0+nkaWkEhEGxH8bGA1I8cNRDtFd0tOi8aKpoj2irSx163iIwwDj0ZQZFpkIAHgago4IeFRXFNt8IYINOLz+67B4527RHz+uv6uJDHmVy1K4FgO4q+Aq6aqJhIh5gmi9S+FdUX1RWNFhll7FEYRZLpkIAHAX0QcKb0+/v3B55Kfxj3xS1E1S9kR1wcMHeux5lctSMBiS9pK9MexbXpDF5jq0ZlZcxBQN7Pw9SpwOWXQx6bbY9ZlQ4j5uA2RHb7FZH1UxDZsxUqVS+FyEhkU3WZkSxVyhx1YCnyR8BujkJ7FF+XL19+QP5w8CoSIAFfBOSBKIyRN57qy3jA5Mnncd7REqnba+PIdoeoGC54ubi2zEbu3AnotYGaPii6RQaqtSNDCx2BYA89ha4mzIkESCBkBO68Exg1ajMWLSmM1EOVcWDdPpyLT8QxRGBHtyFY8e0RfCsDzA88APz8s7xItdmYos2eDTRu7HIWxqTIVPwhYDdHwclsf1qd55CA0QTi41Fo1UpEPP8Y6qZMRuKt9dDjwFQ8OEy6AGI//GBMhup81E6ccC35NzQE7OYoOJkdmvuGuZDAxQT0k3qPSvxOfTEvNhaQbkfNQd1Rt+Y5wxzF/PkXZ8s9wSdgN0cRfGLMgQRIwDeBBvI0/IIFEvpTYn/K+xZXpX2AFHnL+/xZbyHafCfleTQtDfjxR889XA8VAbs5Cg49herOYT4k4ItAYflpue8+52TCVbEHcOJ0CaxrIc+YbNVXp/Jn7E3kj5sRV9nNUXDoyYi7gmmQgFEEatRA5+9kOEps/o6agMxl4LnnkJ/P6c2bZ1ShmE5eCdjNUeS1/jyfBEggyASqVC3kfJz1h0RxGNdfD4wYAbRsCaxd63fO+lisOoqKFf2+hCcaSICOwkCYTIoESMA7gauuAhavLIEz7/0P+OIL4MABIDFRvm8pH7j84w/vF3ns1bmJffuALl08dnI1ZATs5ig4RxGyW4cZkYD/BNRRqD9YsUKu0V5Faipwxx3ACy+4hqMWLfKZmHvYiY7CJ6agHbSbo+AcRdBuFSZMAvkn0KkToPPbf05IV6ggny+bBHz/PXD2LNCxIzBkCPD7714zUUehkWz1LW9a6AnYzVGEniBzJAESyJVARATQvLmXF++6dnW9tj1sGPDWW673L/T1aw87fx5ISQH0VFp4CNBRhIc7cyWBAkdAh42WLwdOnsxR9csuA155xfWNi7JlgZ49XWFqjxxxnrh6tetNbA475eAWwk06ihDCZlYkUJAJ6DyF9g4WL74EhdatXU9CPfkk8PHHQKNGwCefYP48VxgQvd5t+hQULXQE7OYoOJkdunuHOZFAngi0aycfpinmZfjJM5USJYBnngHWrJEPJlcH/vlPzHttM5rFnUWlSp4ncj2UBOzmKDiZHcq7h3mRQB4I6AiTdhr+nND2dW2TJs5xqlOjXsbSw/Xlm91va1xzgF0JX9SCdsxujiJooJgwCZBA4AR0nkHfszt2zI+0JMjg4pYP4ixKoGv9PcA997jeu/DjUp5iLAE6CmN5MjUSIAEfBHSeQTsFGjPQH9PHYnW4qv3SF11PRW3b5rrsf/Li3oUL/iTBcwwgQEdhAEQmQQIk4B+BVq1cn0X19/sU6ijatgUuKys/VYMHA+9MdGX0+mviPdq7XtzzL2ueFQABOooA4PFSEiCBvBEoXhzo0MG/eYrDh4H163O8PxFV2ZXhU08DP/0ENGsGPPus66W9vBWFZ+eBAB1FHmDxVBIggcAJ6DyFRvDYv993WsnJrmEqr+9PdOvmSuSGG4CnngJatABWrfKdII/mm4AVHEWS1G6RSB57gK7TSIAELEzA/T6EOgJfpk9H6ft3GmjWbYUKuddkWVl6F/q+xVdfAfpynj5S9cgjwKlTHidx1QgCwXYUU6SQB0U5P63eXfbprNQO0XCRL5OpL2SISor2+jqRx0iABMxPQEeLypfP5X0KqYbOTyQlAfqFVZ/Wq5fzA0m4+25g7FigaVP/Z8t9JsyDbgLBdhRTJSN1Cp5WRDbkG4noIZJXL9E3a9lYlrNySP7L4OxN6LmPiv4ropEACViYQBH5BVAH4GtCe/duYOfOHPMTvuqswaQmykS3dkMy5ZOrmsG997pif/i6jsf8IpCbr/YrER8nLZRjtXIcT5Rt7Unsyto/XZbXicaI/p61z9tCn7wu4e1A1r6BslRh7969EkQsRVfzbBkZGfm+Ns+ZBfkC1iXIgPORPNvEBa169StkxKgepk9fjqpVT19E8ptvqsq+K1Gu3Er59/jXUNKGDRVkf1N5F2OtBJ09cdF1GqK2sHyru/aUKYgWx3Hms8+w/cEHcbRNm4vP/XNPknMtv78Zfybz54qx6WVkJEjKZbFagl6lp+vgyqXNyvdXLamW59BTb9me5FHVfrL+hsd2zlWZrcIEkTw4jSSRL7tWDk6MiYlx5NeSk5Pze6nprmNdTNckDraJq002b9a3KRyOyZO9t1GfPg5HtWoOR2Zm9uNz57quW7w4+36vW8uXOxyxsa4Lbr3V4Th0yOtprilzr4fytdPo9OLjXVVYuzb34gR6f8nvp4RgvNgKX7zLdHs+lxINEt0sShH5Mobw8EWHx0jAJAQ03p/ORXsbftKRIx1B0qedsk1e57Xs+tKGvgb+9NPAjBlAw4aQLozrdzyvaRXw88PhKNKEeXUP7tGyrvuMMGePIj093Yi0mAYJkECQCKgDuOoql6PQ/y972mYZfzh0yKDPnhYvDowc6QoyqF896itTovqFvTSjfnI8S27f9XA4Cn3YuZ5IWg3SiugjmikywtijMIIi0yCBEBBQR6HfwXZH5XBnqU87qXl9f8J1KO9/G8uzMsuWAf/3f66v6mmX5p132Lvwk2SwHcXHUg5pHTQQ6aOt8vwazouGiuaItoqkT4gtIiOMPQojKDINEggBAbcjyBlNVrcbyC+GRhk31PRxq4cfBjZudH1ub+BAg72RoaU1VWLBdhTSz0M1UTGRDjFNFql9K6ovqisaLTLK2KMwiiTTIYEgE9CRoJo1s89T6OezNWCg24kEpQgxMa5JkAkTXENS7kwYZNBN4qJlsB3FRRkGeQd7FEEGzORJwCgC7nkKfUNbJ7DVVqxwfSo16N/HlsdooT2KLR6DGfoI7aZNroLwbzYCdnMU7FFka15ukIC5Ceg8hX6bYsMGVzl12El/w5OSQlTuaB3oyLKffwYSElyT32fOuPdyKQTs5ijYo+BtTQIWIqCOQs09T6ET2fpbXUHfqwu1bZUp05tukvgP/3UVQrs3NCcBuzkK9ih4Y5OAhQhcfrm8f32la57i999dQ0++hp0Ceq8iNy76Ue4PPwS+lp+R48cBHYp66CHXWFhu19r8uN0chc2bi9UjAfsR0F7FwoWuXsV5eSYyqBPZ/uD7u0QS0jjogwYBr7wC6Pe7vb0Z6E9aNjnHbo6CQ082uTFZjYJDQB3FyZMS7G2MhIguCbRrZ4K6lyvn+vRqSopr0kS914ABrp6GCYoX6iLYzVFw6CnUdxDzI4EACSQluUJ1rFzp+rqpOovcLOfb3Lmdn+/jnTq5Ztr/8x9AAg0iNlZeDzbq/eB8lyrkF9rNUYQcIDMkARIIjEBkJBAf70oj7MNO3qpSujTw4ouuCRQt7HUS7LqPBJQ4eNDb2bbcR0dhy2ZlpUjAWgTcTz/5msgORo3y1DPRz61KqG/nN7q/+EK+ptPINfmdp0SCUYvgp2k3R8E5iuDfM8yBBAwnMFSC+jz7rCuyhuGJG5mgBhl84glg3TqJWFcPuO02+YqOTH7/+quRuZguLbs5Cs5RmO4WY4FIIHcCtWq5fn/1ZTtLmPYmFi8GXn1VPn6Q4pq7eOutv14xt0Ql/C+kVZrF/xrxTBIgARIIBQENMvjAA/JZNomLrt++uO8+oHPnUOQc8jzoKEKOnBmSAAnYioBGN5w7V0KeSsxTdywSraC+FGITo6OwSUOyGiRAAmEkoK+M33WX60U9dzFat87uONz7Lbiko7Bgo7HIJEACJiWgMUncphPc+qTUk08CFg8yaDdHwaee3DcplyRAAuEloGFAbrkFGDUKaNbM9YW98JYo37nbzVHwqad83wq8kARIwFAC+nLetGnA7NmuGCUam2TYMEsGGbSbozC0nZkYCZCAuQgENXpssKravbvrySh9Kuq114C4OMD9YfBg5WlwunQUBgNlciRAAsEnYLmXocuWBd54wxUmV1/au/pq4O67LRNkkI4i+Pc0cyABEjApgZA7nA4dXE9CDR/uGpZq2BDQcCAmNzoKkzcQi0cCJGAzAhoeV2Oqa7jcqlWBG24A/vlP4MAB01bUCo5CyzhaNE50u2lJsmAkQAIkkBcCzZu7nMVo+Xn76itAexfvvQeEvJuTe6GD7SgkgDsOiuQd92wmszvYJtohkj6YT7tOjkaLzon2+jyTB0mABEjASgSKFQMef9w1HKWO4nb5v3DPnsCePaaqRbAdxVSprToFT5MAKRgv6iFqJOqbtWwsy1k5VFm2G4iWiuTjtbhXRCMBEiABexHQD4cvWgS8/rprqR9IGi8/k5mZpqhn0SCXYqGkXytHHomyrT2JXVn7p8tSew0yaAeJ13uRaS/ibNbeCxcd/WvHQFlVYe/evRLQMUVX82wZGRn5vjbPmQX5AtYlyIDzkTzbJB/QPC5Zvz5CtuIlyvc6XLiQ7nEkf6uu3+Ek58X5/c24OOcA0mvcGCUnTUL9l15CRYm9fnzCBJxKT5EsKsqnMFYjPT3j4uw89lj5/qol9fAceuot25M86tZP1t/w2M65Wlp2TBbpHMWQnAe9bSckJDjya8nJyfm91HTXsS6maxIH2ySwNpk3TwfwHY6FCwNLx331hQuu9DRNo8w1yRBgapmZDse77zocERGO+ELrnHVeu+JsrokGen/J7+lqb7+pwR568pZnXvedkgvkgWP8S6RDVr6MITx80eExEiABaxDQNwvvuAPYuhUoV85V5v79XR9MCkMNwuEo0qSe1T3qqhPVuo9GAiRAAiTgSUAfn61dx7Xn0CGgZUtgxAjg9GnPs4K+Hg5HsUpqVU9UW1Rc1Ec0U2SEMdaTERSZBgmQgPkIfPYZ0E9G6p97TqZp4oElS0JWxmA7io+lJstE+uSSTkrrENJ50VDRHJH0qzBDtEVkhHHoyQiKTIMESMB8BMqXB959V3455adTexT6lve/ZET+99+DXtZgOwp99LWaqJhIh5h0UlrtW1F9UV3RaJFRxh6FUSSZDgmYkIA7KKAJ30kLHa1rrnEFGZSnopyP0GqQQXUeQbRgO4ogFt1r0uxReMXCnSRAAt4IWNbhlCnz1zsXpeXBUI1QK5PfRU+c8FbNgPfZzVGwRxHwLcEESIAELENAv3Eh75Q4J7g/+ACJ+qTU0qWGF99ujsJwQEyQBEiABExNQIMM6lf05IW8jLoymh8TY3hx7eYoOPRk+C3CBEmABCxBQJ6E2jh2LFBZIx8Za3ZzFBx6Mvb+YGokQAIkALs5CjYpCZAACZCAwQSKGpxeuJPToadr09MDDxYW7oowfxIgARIwCwG79Sg49GSWO4vlIAESsA0BuzkK2zQMK0ICJEACZiFAR2GWlmA5SIAESMCkBDhHYdKGYbFIgARIwCwE7Naj4ByFWe4sloMEgkjAsqE3AmASzjrL1zFsaRK4Hb/ks2aV5LrD+bzWbJexLmZrEYBtYr420RLZpV0CrUdNYRFlziYyV6lWm6s4AZWGdQkIX1AuZpsEBWvAidqlXYJSD7sNPQV8tzABEiABEiCB7AToKLLz4BYJkAAJkEAOAkVybHPTRWCNjUCwLuZrTLaJ+dpES2SXdrFLPcx5l7BUJEACJEACJEACJEACJEACJEACJEACJEACJBA8AvLhWWwT7RAND142QU95t+SwSbReFJTH5STdYNkUSfigaLNHBhVl/XvRT1nLCh7HzLzqrS4jpcBpIm0bVU+R2a26FDBZlCraInpApGbFdrlUXUZKfazWLvJpO6wUbRBpu/xXpFZbtEKkv2P/ExUX0QwioBP7O0V1RApW4TcSWdF2S6H1xRsrWkcpdHORp6N4UbbdjluXL4isYN7qMlIK/m8rFN6jjNVkXdtEraxou0j/bVixXS5Vl5FSH6u1i74wXUakVkykzqG1aIaoj0jtbdG9zrUA/vDx2L/gJcqqeuBdorOi6aLrRLTQElgo2R3NkaW2w7Ssfbq8Psdxs256q4tZy+qrXPvk4NqsE36X5VbRFSIrtsul6pJVPUstHFLajKwSq6NQ6b6rRJ+K1Az590JH4YKpf/XG//WvTezN2uexyzKrerPMFeljcgMtU+pLF7SKHNJ/4Gr7RbptZRsqhd8o0qEpqwyjuXnXkpVmIv3fq9XbxbMuUh1YsV10JESHMA+KdHhWR0WOi86L1Az5HaOjcMG029/2UiEdKughGiLSIRC7mDpBlVXtLSl4XVG8SJ3fSyKrmA5zfCYaJjqRo9BWa5ecdbFqu1yQdtB7KVqkoyJXigw3Ooq/kOpElk50uU3B6z4rmrvc+r+ML0R6A1nZDkjhdWxZTZdaL6ua1kX/cWeK3hFZpW10WEOdxIeiz0VqVm2XS9XFiu3iaglXLyJZNtqIIkTuT0gY8jtGR+HGDKyS1Xqi2iKdzNbJoJkiq9llUuCyWYXW9WtEnhPDWYcstdB2uD2rxLr8ylKlz15Yt8PTvf8QWaFtCkk5J4u2il4Wuc2K7XKpulixXTTKqzoFtVKiq0XaRuoweovUrP7vxVULk/3tKeXRJzp2ikaYrGz+FqeOnLghS/rInNXq8bGUWYdkzol0fPVuUaRovugn0TxRRZEVzFtd3peC66PLG0X6Q+v5AyWbprT2UiqHSMu8Pkv6b8WK7XKpulixXZpIG6wTabvofzieEqnpb8BK0Q7RJ6ISIhoJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJkAAJGEBgmKRR2oB0mAQJkAAJkIBNCeyWelWyad1YLRIgARIggTwS0LfgvxHpy4768tPTIo1CrC/X6ZuyavqW/DKRRmLVF6E03pDabtGLIj1XX5SKEdFIgARIgARsRuBGqc87HnUqL+u7Re4ehS4XitShqD0qcr9Ju1vW3W/Q95f1WSIaCZAACZCAzQjUl/rsFr0g6iBS2y1yO4q/y/phkTsERqqsa/wktd0iDbugpgHrjjjX+IcETELAHWHQJMVhMUjAsgQ0RlhzUU/RKJHGpvI0DUb3vaiv506PdY2l5DbPdfc+LkmABEiABCxO4HIpv37DWE17D1+KdM5BoxGraaTPPSL3/IMOQWkvRG23aLiuiN0m+tq5xj8kYBIC7FGYpCFYDMsTaCw1GCvS70ycE90raiP6TvSbqLPoDtHHohIitSdE2hNRqyDSKKBnRJfqdcghGgmQAAmQQEEksFsq7Z7LKIj1Z51NTqCwycvH4pEACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZAACZBAHgj8P+cxTHjIxtqqAAAAAElFTkSuQmCC"}},"cell_type":"markdown","metadata":{},"source":["# 2 誤差・精度(25点)\n","\n","最初に紛れ込んだ丸め誤差が次第に拡大されて，最後には真の値と全く違う値となる\n","「不安定な」例として「Numeriacl Recipes in C」で紹介されている漸化式の誤差を検討する．\n","\n","次のような，いわゆる黄金比\n","$$\n","\\phi = \\frac{\\sqrt{5} - 1}{2} = 0.61803399\n","$$\n","の累乗計算を考える．素直に累乗(power)で求めた場合と，\n","\\begin{align}\n","\\phi^{n+1} &= \\phi^{n-1} - \\phi^n \\\\\n","\\phi^0 &= 1.0 \\\\\n","\\phi^1 &= 0.61803398\n","\\end{align}\n","という漸化式(recurrence formula)で求めた場合とで，数値を%20.15fで出力して，\n","それらの誤差をn=30程度までで議論せよ．\n","\n","以下は，それぞれのリスト(phi_power, phi_recur)を片対数でプロットした結果である．\n","\n","「Numeriacl Recipes in C, C言語による数値計算のレシピ」, W.H.Press他(技術評論社, 1993), p.44.\n","\n","``` python\n","import matplotlib.pyplot as plt\n","\n","phi1 = 0.61803399\n","phi_recur = [1]\n","phi_recur.append(phi1)\n","phi_power = [1]\n","phi_power.append(phi1)\n","step = [0, 1]\n","\n","for i in range(2,31):\n","    # ...\n","    #ここを考える\n","    # ...\n","\n","plt.plot(step, phi_power, color = 'r', label=\"power\")\n","plt.plot(step, phi_recur, color = 'b', label=\"recur\")\n","plt.legend()\n","plt.xlabel('step')\n","plt.ylabel('phi^n')\n","plt.yscale('log')\n","plt.grid()\n","plt.show()\n","```\n","\n","![image.png](attachment:image.png)"]},{"cell_type":"markdown","metadata":{},"source":["# 3 数値積分，解，収束性(25点)\n","\n","数値積分に際して，関数値の計算にコストがかかる場合，限られた点での値から高精度で計算するNewton-Cotesの公式を使う．\n","4区間5点の閉区間(closed)公式(Boole's formula)は次のとおりである．\n","\\begin{align}\n","h = & \\frac{b-a}{4} \\\\\n","\\int_a^b f(x) dx = & \\frac{2}{45}h\\left\\{7f(a)+32f(a+h)+12f(a+2h)\\right. \\\\\n","& \\left. +32f(a+3h)+7f(a+4h)\\right\\}\n","\\end{align}\n","\n","これに従って求めた\n","\\begin{align}\n","\\int_0^1 \\frac{4}{1+x^2} dx & \\,(= \\pi)\n","\\end{align}\n","の値は以下の通りである．\n","\n","数値積分の中点則を用いて求めた場合，同じ程度の精度を達成するには何点が必要となるか．\n","``` python\n","import numpy as np\n","def func(x):\n","    return 4.0/(1+x**2)\n","\n","n = 4    \n","h = 1.0/n\n","\n","boole = 2/45*h*(7*func(0)+32*func(h)+12*func(2*h)+32*func(3*h)+7*func(4*h))\n","print(\"%26s : %15.10f\" % (\"from Boole's formula\",boole) )\n","print(\"%26s : %15.10f\" % (\"actual Error\", boole-np.pi))\n","print(\"%26s : %15.10f\" % (\"Estimated Error f^(6)(0.5)\", h**7*8/945*1311.8)) # 1311.8=$f^6(0.5)$\n","```\n","\n","```\n","      from Boole's formula :    3.1421176471\n","              actual Error :    0.0005249935\n","Estimated Error f^(6)(0.5) :    0.0006778067\n","```"]},{"cell_type":"markdown","metadata":{},"source":["ちなみに問１で求めた補間多項式を積分した値は，3.142183333となりBoole公式で求めた値と大体一致する．"]},{"cell_type":"markdown","metadata":{},"source":["# 4 微分方程式 (25点)\n","\n","\n","Verlet法による小惑星軌道のシミュレーションを次のような条件で行った．\n","```python\n","def force(pos):\n","    x=pos[0]\n","    y=pos[1]\n","    L=(x*x+y*y)**(3/2)\n","    return [-x/L,-y/L]\n","\n","def Verlet(r0,rh):\n","    f=force(r0)\n","    x=2*r0[0]-rh[0]+h**2/m*f[0]\n","    y=2*r0[1]-rh[1]+h**2/m*f[1]\n","    return [x,y]\n","\n","h=0.1\n","dx, dy=0.0, h\n","m=0.2\n","xx=[3.0, 3.0-dx]\n","yy =[0.0,-dy]\n","```\n","\n","1. この小惑星軌道をplotせよ．\n","1. この小惑星の公転周期が規格化した単位でいかほどになるか答えよ．\n","1. 異なる初期条件を使って軌跡を表示し，その飛翔体の振る舞いを解説せよ．\n"]}],"metadata":{"kernelspec":{"display_name":"Python 3","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.8.5"},"toc":{"base_numbering":1,"nav_menu":{},"number_sections":true,"sideBar":true,"skip_h1_title":false,"title_cell":"Table of Contents","title_sidebar":"Contents","toc_cell":false,"toc_position":{},"toc_section_display":true,"toc_window_display":true},"vscode":{"interpreter":{"hash":"f3f87633aac09da3bda522f97956bee375b5501d1579e6458804e567301cb62a"}}},"nbformat":4,"nbformat_minor":4}