restart; Digits:=14: with(LinearAlgebra): with(orthopoly): Laguerre basis set for central potentials. L(2,1,x);

LCgiIiQiIiJJInhHNiIhIiQqJEYlIiIjI0YkRik=

Radial orbitals (multiplied by r, leaving out Y_Lm): look almost like hydrogenic orbitals, but form an enumerably infinite set. uL:=(n,l,beta)->r^(l+1)*exp(-beta*r/2)*L(n-1,2*l+2,beta*r);

Zio2JUkibkc2IkkibEdGJUklYmV0YUdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlKigpSSJyR0YlLCY5JSIiIkYwRjBGMC1JJGV4cEdGJTYjLCQqJjkmRjBGLUYwIyEiIiIiI0YwLUkiTEdGJTYlLCY5JEYwRjhGMCwmRi9GOUY5RjBGNUYwRiVGJUYl

plot([seq(uL(n,0,1.5),n=1..4)],r=0..20,color=[red,blue,green,black]);

NiktJSdDVVJWRVNHNiQ3ZW83JCQiIiEhIiIkIiIhISIiNyQkIjJLTEwzRldZcyMhIz0kIjJiYEAjXDtgcEUhIz03JCQiMmttbTthKUdcYSEjPSQiMjkncFhKdzFKXyEjPTckJCIsRCJHJFI8KSEjNyQiMVgtM2pZKXlvKCEjPDckJCIyTExMJDN4JikqMyIhIzwkIjJ0eGgoM3JKLzUhIzw3JCQiLERjJ3lNOyEjNiQiMiZmXXgocFpoVyIhIzw3JCQiMm1tbW1UOih6QCEjPCQiMmt5U2MiMyk0Jj0hIzw3JCQiMUxMJDNGV1lzIyEjOyQiMmR3LSkpNHI1QSMhIzw3JCQiK0RKZHBLISM1JCIxVVMieiQ+YmVEISM7NyQkIjFubTt6Pl05USEjOyQiMXIpPkhkVGEnRyEjOzckJCIxTExMTDNWZlYhIzskIjEyT1QlKilbTzkkISM7NyQkIjFtbSJ6PjV4SSYhIzskIjJ4Jz1eIlIzWmMkISM8NyQkIixEYyopZkQnISM2JCIyWSEpUXkpKlJKIlIhIzw3JCQiMUxMM0Yqb1U/KCEjOyQiMWJyRSM9SnA+JSEjOzckJCIxbm1tIkhbRDopISM7JCIwZCYqeTE6TFUlISM6NyQkIi12JHBVJkc1ISM2JCIyRVM5QjsiZmJaISM8NyQkIjJMTEwkZTAkPUMiISM7JCIxMFR4JkdtSCpbISM7NyQkIjJMTCQzeGZdJkgiISM7JCIwKGVIc3UvLlwhIzo3JCQiMk1MTGVSIj1cOCEjOyQiMllIZi1ZOlohXCEjPDckJCIyT0wkZTlvJkdTIiEjOyQiMVFZU3UlPScpKlshIzs3JCQiMk5MTExCS2xYIiEjOyQiMU1IZCRvbWApWyEjOzckJCIyT0xMMzIkKVFjIiEjOyQiMTxsWFs0bVJbISM7NyQkIjJPTEwkM1JCcjshIzskIjEvSjlVPXpyWiEjOzckJCIyLStdUDk6XCk9ISM7JCIyMnVyJT4yJ1xlJSEjPDckJCIybW1tInpqZik0IyEjOyQiMUJMP3cnRylbViEjOzckJCIsdkIxbkgjISM1JCIyYTk/RiNIQi1UISM8NyQkIjJNTExlNDtbXCMhIzskIjI4QiNccCNHMyVRISM8NyQkIjFubTt6dCUqKnAjISM6JCIyOVNzKGVcIVJjJCEjPDckJCIsRG15XSFIISM1JCIydyczbTAqZXlHJCEjPDckJCIxbm07SGRBPEohIzokIjJOdT4pZVUrNEkhIzw3JCQiMlBMTGV6cyRITCEjOyQiMUMqSExhWDV1IyEjOzckJCIybm1tVF1SM2EkISM7JCIyKVJlKEdILXdbIyEjPDckJCIsREAxQnYkISM1JCIyKWU3PUcyYVxBISM8NyQkIjFNTGU5ZCMpcFIhIzokIjEsNyssJilwQD8hIzs3JCQiMW5tbTtfTSg9JSEjOiQiMSUpZTlcIXk5Ij0hIzs3JCQiMU1MTDN5X3FYISM6JCIyViVHJDQzc0xbIiEjPDckJCIyMCsrK2wrPismISM7JCIyRWY8d2llWTwiISM8NyQkIip2V11WJiEiKSQiMDhfJXoiKVxCIyohIzs3JCQiKk5mQyZlISIpJCIxOio9PTZOQUUoISM8NyQkIjFNTCRlejY6QichIzokIjFcI1sxVWMiPmUhIzw3JCQiMW5tbTs9QyNvJyEjOiQiMVAhNCdmLTxdVyEjPDckJCIxbm1tbSNwUzEoISM6JCIyOi0kPj1BKEhgJCEjPTckJCIsRE9EIzN2ISM1JCIyOG1HaUZON3AjISM9NyQkIjFtbW1tKHk4IXohIzokIjFSR11aWCIqM0AhIzw3JCQiMSwrXWkudEskKSEjOiQiMjE1P2cqSEs0OyEjPTckJCIsdjN6TXUpISM1JCIyOzY8PycpUjRDIiEjPTckJCIxb21tIkhfPzwqISM6JCIxUFlEXWNEUiUqISM9NyQkIjFubTt6aWhsJiohIzokIjExWzsjUjgiR3QhIz03JCQiMUxMTDMjRywqKiohIzokIjFuSWBpWVdtYiEjPTckJCIyTUwkZXp3NVY1ISM6JCIxWzpKUEdedlQhIz03JCQiLXYkUSNcIjMiISM1JCIxJjMtcTI5aUMkISM9NyQkIjJOTExlIipbSDciISM6JCIxRyM+NyUqNCpwQyEjPTckJCIqZHhkOyIhIigkIjJqbnJKK20nZj0hIz43JCQiLEQweHc/IiEiKiQiMmlfcUMrLXFTIiEjPjckJCIyLCtdaSZwQFs3ISM6JCIyKSpSNWVxb0gyIiEjPjckJCIsdmdIS0giISIqJCIxVmMlMyI0YEp6ISM+NyQkIjJsbW1tWnZPTCIhIzokIjAwaF8lKmYkUmchIz03JCQiMiwrK10yZ29QIiEjOiQiMSVwJnkpKlsjKjRYISM+NyQkIjJLTCRlUjwqZlQiISM6JCIyYjguJmYpUiVlTSEjPzckJCIyLSsrXSlIeGU5ISM6JCIyajZBaVUmKlxlIyEjPzckJCIybW07SCFvLSpcIiEjOiQiMigpUmlYZzdUJz4hIz83JCQiMiwrXTdrLjZhIiEjOiQiMjhQJj5BSHhzOSEjPzckJCIybW1tO1dUQWUiISM6JCIyWkN1SipIbTU2ISM/NyQkIi1EMSozYGkiISM1JCIxMHZUIil6c2YjKSEjPzckJCIyTkxMTCp6eW07ISM6JCIxUXUlPTkiKmU/JyEjPzckJCIyT0xMM04xIzQ8ISM6JCIxPiVHOSkqKnBIWSEjPzckJCIxbm0iSFl0N3YiISM5JCIxKClRVSFRPSxZJCEjPzckJCIyLSsrK3hHKip5IiEjOiQiMnlCMz0lUlhZRSEjQDckJCIxbm07OUBCTT0hIzkkIjFuQSYzKCpmXyU+ISM/NyQkIjJPTExMYmRRKD0hIzokIjImZTQnR3ZkaloiISNANyQkIi1ET2w1Oz4hIzUkIjItMlhATm0nKjQiISNANyQkIjItK11QP1dsJj4hIzokIjExMWYlNCQ+IkgpISNANyQkIiQrIyEiIiQiMTtsLjVrLz1oISNALSUmQ09MT1JHNiYlJFJHQkckIiM1ISIiJCIiISEiIiQiIiEhIiItJSdDVVJWRVNHNiQ3am83JCQiIiEhIiIkIiIhISIiNyQkIjJLTEwzRldZcyMhIz0kIjE2biZmMyNcKip5ISM8NyQkIjJrbW07YSlHXGEhIz0kIjJlIip6Uyo9Y0U6ISM8NyQkIixEIkckUjwpISM3JCIxdFknKT5dNTdBISM7NyQkIjJMTEwkM3gmKSozIiEjPCQiMjoyJypbInB3W0chIzw3JCQiLERjJ3lNOyEjNiQiMjtCbyhwOyNRKVIhIzw3JCQiMm1tbW1UOih6QCEjPCQiMFJeUnRdeCVcISM6NyQkIjFMTCQzRldZcyMhIzskIjFIJFEoKioqb2F2JiEjOzckJCIrREpkcEshIzUkIjFqU1o6KlwzVSchIzs3JCQiMW5tO3o+XTlRISM7JCIwQSc0TykqeWNwISM6NyQkIjFMTExMM1ZmViEjOyQiMSRvPWBwb19QKCEjOzckJCIxbW0iej41eEkmISM7JCIxLzs/NyRmZyZ5ISM7NyQkIixEYyopZkQnISM2JCIxdSNRYV9OdDEpISM7NyQkIjFtO3pXI0gsdCchIzskIjF3NGhuJT5tMykhIzs3JCQiMUxMM0Yqb1U/KCEjOyQiMTNWKDRUP2EwKSEjOzckJCIudiQ0JzMleXchIzgkIjF4c2RKWmh5eiEjOzckJCIxbm1tIkhbRDopISM7JCIxJFJYJ1s8dmd5ISM7NyQkIi12JHBVJkc1ISM2JCIxJlFfLW4heUhwISM7NyQkIjJMTEwkZTAkPUMiISM7JCIxbShvYUZZWGMmISM7NyQkIjJNTExlUiI9XDghIzskIjFQNF8vLDcpeSUhIzs3JCQiMk5MTExCS2xYIiEjOyQiMSNvLjgpKmVEKVIhIzs3JCQiMk9MTDMyJClRYyIhIzskIjFSSWQxcilmOyQhIzs3JCQiMk9MTCQzUkJyOyEjOyQiMmpTI0dKXT9gQiEjPDckJCIyb21UZ191IXk8ISM7JCIxRF9PeVYpKWY6ISM7NyQkIjItK11QOTpcKT0hIzskIjEsWG96SipbInohIzw3JCQiMk1MZTl3YjwqPiEjOyQiMVYmKT01TyYqSGIhIz03JCQiMm1tbSJ6amYpNCMhIzskITFPc007OW9KayEjPDckJCIsdkIxbkgjISM1JCEyTFF0UDlQZCM9ISM8NyQkIjJNTExlNDtbXCMhIzskITFtWWh0YXZdRyEjOzckJCIxbm07enQlKipwIyEjOiQhMVp4Qi4qPT11JCEjOzckJCIsRG15XSFIISM1JCEydlBPN1RjT1klISM8NyQkIjFubTtIZEE8SiEjOiQhMVUjKW9sYWdVXSEjOzckJCIyUExMZXpzJEhMISM7JCExLVg7bXEhZVkmISM7NyQkIjJubW1UXVIzYSQhIzskITE0IXlfI1FcXGQhIzs3JCQiLERAMUJ2JCEjNSQhMSJwJiopKkdFRyJmISM7NyQkIjFuO2FqZjFoUSEjOiQhMU0/JlFbb2MmZiEjOzckJCIxTUxlOWQjKXBSISM6JCEwNzU7YSFmdGYhIzo3JCQiLkRjWSZleVMhIzckITFDeGlpIT0pb2YhIzs3JCQiMW5tbTtfTSg9JSEjOiQhMXVRJSlvPFxWZiEjOzckJCIxLCtdN2wkKnlWISM6JCExSnchKiopZTFiZSEjOzckJCIxTUxMM3lfcVghIzokITFUWkB3UGQ+ZCEjOzckJCIyMCsrK2wrPismISM7JCEwL2dtdTckKkcmISM6NyQkIip2V11WJiEiKSQhMjBDJXBcJHBDdiUhIzw3JCQiKk5mQyZlISIpJCEyPHdqUCQpPm0+JSEjPDckJCIxTUwkZXo2OkInISM6JCExX2VPYlZkJHAkISM7NyQkIjFubW07PUMjbychIzokITFvVDcjKmZeREohIzs3JCQiMW5tbW0jcFMxKCEjOiQhMmp3P1hTI28kbyMhIzw3JCQiLERPRCMzdiEjNSQhMTZtdWcmKmVCQSEjOzckJCIxbW1tbSh5OCF6ISM6JCExYUBxSGIjbyc9ISM7NyQkIjEsK11pLnRLJCkhIzokITElKlsqKipHNihHOiEjOzckJCIsdjN6TXUpISM1JCEyJUhPRCx5QmI3ISM8NyQkIjFvbW0iSF8/PCohIzokITJcYXFuXyNbOjUhIzw3JCQiMW5tO3ppaGwmKiEjOiQhMSg+LGo0YWlKKSEjPDckJCIxTExMMyNHLCoqKiEjOiQhMSh6TXgoNFxybSEjPDckJCIyTUwkZXp3NVY1ISM6JCEyL0MvITQvaCFHJiEjPTckJCItdiRRI1wiMyIhIzUkITFoIVx0WnBBSCUhIzw3JCQiMk5MTGUiKltINyIhIzokITEpUSJvYDVTPk0hIzw3JCQiKmR4ZDsiISIoJCEyMjhLLmROU3AjISM9NyQkIixEMHh3PyIhIiokITJra21nVzBuNyMhIz03JCQiMiwrXWkmcEBbNyEjOiQhMic9VFYnKWYwKG8iISM9NyQkIix2Z0hLSCIhIiokITEmbyFvandrKzghIzw3JCQiMmxtbW1adk9MIiEjOiQhMiM0ZXFSNitGNSEjPTckJCIyLCsrXTJnb1AiISM6JCEwJlIlPiRlS2h6ISM8NyQkIjJLTCRlUjwqZlQiISM6JCExTUJdXjk6M2ohIz03JCQiMi0rK10pSHhlOSEjOiQhMVVhSjhPKTMpWyEjPTckJCIybW07SCFvLSpcIiEjOiQhMGA0ZlBfciNRISM8NyQkIjIsK103ay42YSIhIzokITFCZi4kcDZGJ0ghIz03JCQiMm1tbTtXVEFlIiEjOiQhMllaT2BpMUdJIyEjPjckJCItRDEqM2BpIiEjNSQhMmRdbiN5JioqZXciISM+NyQkIjJOTExMKnp5bTshIzokITFGJD5BKiozYU8iISM9NyQkIjJPTEwzTjEjNDwhIzokITIvcmpAKmUyWzUhIz43JCQiMW5tIkhZdDd2IiEjOSQhMWRFOkxaUV4hKSEjPjckJCIyLSsrK3hHKip5IiEjOiQhMU83JTMqKTQ6SichIz43JCQiMW5tOzlAQk09ISM5JCEyTGBhQFo0Jm9aISM/NyQkIjJPTExMYmRRKD0hIzokITE5RzJrKT1vcSQhIz43JCQiLURPbDU7PiEjNSQhMUNiIj4jcHJJRyEjPjckJCIyLStdUD9XbCY+ISM6JCEyS0ZyZTJ4WD0jISM/NyQkIiQrIyEiIiQhMiNmKTQyYHM9bCIhIz8tJSZDT0xPUkc2JiUkUkdCRyQiIiEhIiIkIiIhISIiJCIjNSEiIi0lJ0NVUlZFU0c2JDdlcDckJCIiISEiIiQiIiEhIiI3JCQiMktMTDNGV1lzIyEjPSQiMmkyOlp6KyRlOiEjPDckJCIya21tO2EpR1xhISM9JCIxJmUqcD1YTnBIISM7NyQkIixEIkckUjwpISM3JCIxKFFdbSd5WVRVISM7NyQkIjJMTEwkM3gmKSozIiEjPCQiMD0qW0lgZSNRJiEjOjckJCIybG07YThBQk8iISM8JCIxUDBtRUZJK2shIzs3JCQiLERjJ3lNOyEjNiQiMSUpSHhuLSk9SSghIzs3JCQiMk1MJGUqKTREMj4hIzwkIjE4PiVHPmJVNCkhIzs3JCQiMm1tbW1UOih6QCEjPCQiMCg0KCpwUjAleSkhIzo3JCQiK0RKZHBLISM1JCIyYnJjMiEqeVIxIiEjOzckJCIxTExMTDNWZlYhIzskIjI9JiozbC5JNjgiISM7NyQkIjJrO3pXbitsZiUhIzwkIjF0T1h2WzNMNiEjOjckJCIuRGNecU4kWyEjOCQiMiMpPVBLYlI2OCIhIzs3JCQiMUszeGMua3FdISM7JCIyKVEpUTt1YGI3IiEjOzckJCIxbW0iej41eEkmISM7JCIxVERQQmZkOzYhIzo3JCQiMUwkMy0pKVw9eSYhIzskIjJHZ2kyZSVRKjMiISM7NyQkIixEYyopZkQnISM2JCIyR0hePThXODAiISM7NyQkIjFMTDNGKm9VPyghIzskIjE0PiQ+OFsxXCohIzs3JCQiMW5tbSJIW0Q6KSEjOyQiMT12eiwqPjBAKSEjOzckJCIxTUxlOXcpKj0jKiEjOyQiMTwrSnNdbXlsISM7NyQkIi12JHBVJkc1ISM2JCIxMTMhUlokUVhbISM7NyQkIjJsbVRnaSc9TjYhIzskIjI7Qm9CYC9VNCQhIzw3JCQiMkxMTCRlMCQ9QyIhIzskIjF0eiFSYlsjKlEiISM7NyQkIjJNTExlUiI9XDghIzskITJ0KipIaGBbIT1CISM9NyQkIjJOTExMQktsWCIhIzskITI6PWx3cic+QTwhIzw3JCQiMk9MTDMyJClRYyIhIzskITJgIm8zJikpW3owJCEjPDckJCIyT0xMJDNSQnI7ISM7JCEyJmVNOztBQ0NVISM8NyQkIjItK11QOTpcKT0hIzskITEzR1hOYWo8ZyEjOzckJCIybW1tInpqZik0IyEjOyQhMS4kcCU9MSYpPXIhIzs3JCQiMk9MTCQzOGwoPiMhIzskITEwSyEzKHB0NHUhIzs3JCQiLHZCMW5IIyEjNSQhMV4pMypHZSFHZCghIzs3JCQiMk5MJDMtUEJZQiEjOyQhMSlbZ2E1WiQ0dyEjOzckJCIybG1tbTtoZFIjISM7JCExeHJyZGxhPHchIzs3JCQiLURKJylHWEMhIzYkITE9PEx6OG8pZighIzs3JCQiMk1MTGU0O1tcIyEjOyQhMXczajIiW1NiKCEjOzckJCIxbm07enQlKipwIyEjOiQhMSI+PS80cEw3KCEjOzckJCIsRG15XSFIISM1JCExYDsjKVFZWmxqISM7NyQkIjFubTtIZEE8SiEjOiQhMVIucEZhIjNMJiEjOzckJCIyUExMZXpzJEhMISM7JCEwKFtRbFIheTclISM6NyQkIjJubW1UXVIzYSQhIzskITJFd3lMJj1oT0chIzw3JCQiLERAMUJ2JCEjNSQhMm9UYTdxRGpeIiEjPDckJCIxTUxlOWQjKXBSISM6JCExZlBxLDVSND0hIzw3JCQiMW5tbTtfTSg9JSEjOiQiMjsxJ2Z2QmsqMyIhIzw3JCQiMSwrXTdsJCp5ViEjOiQiMD1pbjhYLjgjISM6NyQkIjFNTEwzeV9xWCEjOiQiMWlsMSZvNUAzJCEjOzckJCIxbm07SFVAJ3klISM6JCIwbG1pIikqXE5TISM6NyQkIjIwKysrbCs+KyYhIzskIjJtaC02JUc+ZFshIzw3JCQiKnZXXVYmISIpJCIxdi1IIWY5eDUnISM7NyQkIiowX1BrJiEiKSQiMSYqZlo1WTdHbCEjOzckJCIqTmZDJmUhIikkIjEoUSg+KFFEJlJvISM7NyQkIjFubSJIZCYpPi8nISM6JCIxQyM9PmpwXy4oISM7NyQkIjFNTCRlejY6QichIzokIjEjXHY7eC9jOighIzs3JCQiMW47LyxWPldqISM6JCIwVyRvRGBZJT4oISM6NyQkIi1EMW8ob1gnISM2JCIwRCdbOiEpKjNAKCEjOjckJCIxTSRlOUpmJnBsISM6JCIxOjA9KFIiWzFzISM7NyQkIjFubW07PUMjbychIzokIjE8RW5FWngjPSghIzs3JCQiMW9tbVRiOnRvISM6JCIxOicqUW5VJ0c1KCEjOzckJCIxbm1tbSNwUzEoISM6JCIxbzVVX3A4enAhIzs3JCQiLERPRCMzdiEjNSQiMUYsYWUlXChlbCEjOzckJCIxbW1tbSh5OCF6ISM6JCIxeVNyJXBxJXpnISM7NyQkIjEsK11pLnRLJCkhIzokIjEqXFA9Y3gwXCYhIzs3JCQiLHYzek11KSEjNSQiMS4pKVE7TzcyXCEjOzckJCIxb21tIkhfPzwqISM6JCIyTnVIZigqW19JJSEjPDckJCIxbm07emlobCYqISM6JCIyJD5yVlU+RnhQISM8NyQkIjFMTEwzI0csKioqISM6JCIxTTpEXiQ0dEMkISM7NyQkIjJNTCRlenc1VjUhIzokIjJZRnohZScpUltGISM8NyQkIi12JFEjXCIzIiEjNSQiMTFXMiIpPnlmQiEjOzckJCIyTkxMZSIqW0g3IiEjOiQiMF8lKXBoZ3opPiEjOjckJCIqZHhkOyIhIigkIjIjNGYoeXohM2E7ISM8NyQkIixEMHh3PyIhIiokIjIpenQjPmopXHQ4ISM8NyQkIjIsK11pJnBAWzchIzokIjJyIj42KlErOjkiISM8NyQkIix2Z0hLSCIhIiokIjFNdDYlUXlZQyohIzw3JCQiMmxtbW1adk9MIiEjOiQiMUJVYzVOYjl3ISM8NyQkIjIsKytdMmdvUCIhIzokIjJhYyUzSU1DamghIz03JCQiMktMJGVSPCpmVCIhIzokIjImM3l2KkguLjImISM9NyQkIjItKytdKUh4ZTkhIzokIjJPNy5RLiEzIjMlISM9NyQkIjJtbTtIIW8tKlwiISM6JCIxQV5xNEtdO0whIzw3JCQiMiwrXTdrLjZhIiEjOiQiMiYqeiJcO1toaEUhIz03JCQiMm1tbTtXVEFlIiEjOiQiMiQ+QSRbWE4uOSMhIz03JCQiLUQxKjNgaSIhIzUkIjFwcGBaWnQpcCIhIzw3JCQiMk5MTEwqenltOyEjOiQiMkoheSozU0NpTiIhIz03JCQiMk9MTDNOMSM0PCEjOiQiMlUyb18mKnBYMiIhIz03JCQiMW5tIkhZdDd2IiEjOSQiMWJDPWVWUTUmKSEjPTckJCIyLSsrK3hHKip5IiEjOiQiMSR5eipwMEhibyEjPTckJCIxbm07OUBCTT0hIzkkIjFDKilIUyczJ1FgISM9NyQkIjJPTExMYmRRKD0hIzokIjF5KVxKInlvZ1UhIz03JCQiLURPbDU7PiEjNSQiMlohZV1yQHpWTCEjPjckJCIyLStdUD9XbCY+ISM6JCIxJTRySDomM1pFISM9NyQkIiQrIyEiIiQiMiRwQXgkZmpjMCMhIz4tJSZDT0xPUkc2JiUkUkdCRyQiIiEhIiIkIiM1ISIiJCIiISEiIi0lJ0NVUlZFU0c2JDdecjckJCIiISEiIiQiIiEhIiI3JCQiMm1tO2E4QUJPIiEjPSQiMnglb2hWJGU1SyIhIzw3JCQiMktMTDNGV1lzIyEjPSQiMVJbYSpcUzpjIyEjOzckJCItRDFrJ3AzJSEjOCQiMTxCOSgqXEpDUCEjOzckJCIya21tO2EpR1xhISM9JCIxKip6VHFPPDdbISM7NyQkIixEIkckUjwpISM3JCIxO2gmKUd3JFJ4JyEjOzckJCIyTExMJDN4JikqMyIhIzwkIjFXTiYqeXNwbiUpISM7NyQkIjJsbTthOEFCTyIhIzwkIjEoKilIZ1EnPjgqKiEjOzckJCIsRGMneU07ISM2JCIyRyE0OGgyIkg2IiEjOzckJCIyTUwkZSopNEQyPiEjPCQiMi01dSFlITNMQCIhIzs3JCQiMm1tbW1UOih6QCEjPCQiMkE6UGUhKXlUSCIhIzs3JCQiMUxMJDNGV1lzIyEjOyQiMSlmNCpmeGEuOSEjOjckJCIrREpkcEshIzUkIjJadGg/Rm1EWCIhIzs3JCQiMk9MJDNfdi5VTiEjPCQiMmQoKWVqOCZwZDkhIzs3JCQiMW5tO3o+XTlRISM7JCIyJHBPNGNnW145ISM7NyQkIi1EMWsncDMlISM3JCIyayZ6QlRzLk45ISM7NyQkIjFMTExMM1ZmViEjOyQiMU08RSwsUTQ5ISM6NyQkIjFtbSJ6PjV4SSYhIzskIjIkXDJbdmhhaDchIzs3JCQiLERjKilmRCchIzYkIjIjW106djtqWzUhIzs3JCQiMW07elcjSCx0JyEjOyQiMWB3QXgiUXZFKiEjOzckJCIxTEwzRipvVT8oISM7JCIxJnlxOS1FYyl6ISM7NyQkIi52JDQnMyV5dyEjOCQiMUVOQWg2IWVtJyEjOzckJCIxbm1tIkhbRDopISM7JCIxaip5JFJXL0lgISM7NyQkIjEsXTdgencmbykhIzskIjEvOyJ5KT1SS1EhIzs3JCQiMU1MZTl3KSo9IyohIzskIjItTVVQPlc7TyMhIzw3JCQiMW47L3dzP18oKiEjOyQiMVItJ3AybS1QKiEjPDckJCItdiRwVSZHNSEjNiQhMk5NLDp0WWxEJSEjPTckJCIySyRlKilmWSc9MyIhIzskITFyM14wb2o4PCEjOzckJCIybG1UZ2knPU42ISM7JCExWEkjKmV1Jm8iSCEjOzckJCIvdj0jZjMmKT0iISM4JCEyTyRlImZMVXctJSEjPDckJCIyTExMJGUwJD1DIiEjOyQhMTY5X0cyVVNdISM7NyQkIjJNTExlUiI9XDghIzskITFkYzlMWm1vbiEjOzckJCIyTkxMTEJLbFgiISM7JCExaU4iUjg/VjIpISM7NyQkIjJPTEwzMiQpUWMiISM7JCExeCVwMydIUmwqKSEjOzckJCIyT0xMJDNSQnI7ISM7JCEwJj00bE1LayUqISM6NyQkIjE8L3dpUyV6cCIhIzokITFGNUV3TjpKJiohIzs3JCQiMixdKD08VWxDPCEjOyQhMU0jeSJbIWVpZCohIzs3JCQiMk5lOTtQazh2IiEjOyQhMCpwPkA4QisnKiEjOjckJCIyb21UZ191IXk8ISM7JCExQWtMNmxvLicqISM7NyQkIjItdm8vbyV5Lz0hIzskITAkZT5NR0QoZSohIzo3JCQiMk4kZSpbJFtcSj0hIzskITEmZiY0Ij10OmIqISM7NyQkIjE8SEsqKVw/ZT0hIzokITFmNWpTPUkoXCohIzs3JCQiMi0rXVA5OlwpPSEjOyQhMXQnKTNmQjVEJSohIzs3JCQiMk1MZTl3YjwqPiEjOyQhMSRROClwITQwKCopISM7NyQkIjJtbW0iempmKTQjISM7JCExJlFUXSZvSCRHKSEjOzckJCIsdkIxbkgjISM1JCExJFxTVykpMyRSbCEjOzckJCIyTUxMZTQ7W1wjISM7JCEyJylSeCJcWDctVyEjPDckJCIsdnQiUShmIyEjNSQhMU5PI1x6cXRAJCEjOzckJCIxbm07enQlKipwIyEjOiQhMic9KHArIz11Oj8hIzw3JCQiMk5MTDMtOEQhRyEjOyQhMSJma1U7LExAKSEjPDckJCIsRG15XSFIISM1JCIxTiMzTyoqKnBbTSEjPDckJCIxbm07SGRBPEohIzokIjEsNUhGQGUkZiMhIzs3JCQiMlBMTGV6cyRITCEjOyQiMmNrSU1HKnBGWCEjPDckJCIybm1tVF1SM2EkISM7JCIxQWQuSVdBcGchIzs3JCQiLERAMUJ2JCEjNSQiMURlNGVVYCE+KCEjOzckJCIxbjthamYxaFEhIzokIjE9UkkleXY4ZyghIzs3JCQiMU1MZTlkIylwUiEjOiQiMGQoUixYdSx6ISM6NyQkIjFuVDUhZjBVLSUhIzokIjE3SytOJz06LCkhIzs3JCQiLkRjWSZleVMhIzckIjE0dCdSJFEvJjQpISM7NyQkIjFNZTlUYCdIOCUhIzokIjEvaTZhIj5IOikhIzs3JCQiMW5tbTtfTSg9JSEjOiQiMCJ6XnQsIWU9KSEjOjckJCIxLF1pU0lDTlUhIzokIjEoNEd1PT5ZPikhIzs3JCQiMU1MZWszOSRHJSEjOiQiMWNGWTNfNSY9KSEjOzckJCIxbzthKW9RNUwlISM6JCIxWl1NLEZ6ZCIpISM7NyQkIjEsK103bCQqeVYhIzokIjFkNz4oeU9LNikhIzs3JCQiMW9tVGdAdHVXISM6JCIxWTEkbzItWih6ISM7NyQkIjFNTEwzeV9xWCEjOiQiMCdmYXBPRHV4ISM6NyQkIjFubTtIVUAneSUhIzokIjEiPS1gZScpPjcoISM7NyQkIjIwKysrbCs+KyYhIzskIjFUaWRWKm8jUWkhIzs3JCQiKUZaPV8hIigkIjF5Sl9DTjl2XiEjOzckJCIqdlddViYhIikkIjFwXlAleTMsKlIhIzs3JCQiKjBfUGsmISIpJCIyVzBKVm5ZInpGISM8NyQkIipOZkMmZSEiKSQiMmBRejNgai5hIiEjPDckJCIxbm0iSGQmKT4vJyEjOiQiMVdcPWAnZSMqPiUhIzw3JCQiMU1MJGV6NjpCJyEjOiQhMVVnVyVIMyVSbiEjPDckJCItRDFvKG9YJyEjNiQhMiMqPj9dOiVvOj4hIzw3JCQiMW5tbTs9QyNvJyEjOiQhMk03KlJua05ySSEjPDckJCIxb21tVGI6dG8hIzokITFqaCcpKSplOnFSISM7NyQkIjFubW1tI3BTMSghIzokITBbPDxMZHd5JSEjOjckJCIsRE9EIzN2ISM1JCExVnomeXBcUE4nISM7NyQkIjFtbW1tKHk4IXohIzokITEmKVJpLSFvTk0oISM7NyQkIjFNTGVrWDA8IikhIzokITEnSHhEJ1sqXHQoISM7NyQkIjEsK11pLnRLJCkhIzokITBmNlYiPTdFISkhIzo3JCQiK0RaNVEmKSEiKiQhMSg0Jio9REZsQCkhIzs3JCQiLHYzek11KSEjNSQhMUR4RGAzMUgkKSEjOzckJCIxTkxlKm9seCYqKSEjOiQhMTFsbl41OXIkKSEjOzckJCIxb21tIkhfPzwqISM6JCExKT5SPHNFV00pISM7NyQkIjFvbVQmR00pbyQqISM6JCExQGtgLmpPbSMpISM7NyQkIjFubTt6aWhsJiohIzokITEteiUpeiYpXFYiKSEjOzckJCIxTExMMyNHLCoqKiEjOiQhMShSd0E0SlJ2KCEjOzckJCIyTUwkZXp3NVY1ISM6JCExKilHJVI5bnRAKCEjOzckJCItdiRRI1wiMyIhIzUkITFIQFJZMSg+bychIzs3JCQiMk5MTGUiKltINyIhIzokITFCZk4lPSdSbmchIzs3JCQiKmR4ZDsiISIoJCExJ3lONydbc0FhISM7NyQkIixEMHh3PyIhIiokITJjNCl5SElLMFshIzw3JCQiMiwrXWkmcEBbNyEjOiQhMiU9JGUsK2ZkQiUhIzw3JCQiLHZnSEtIIiEiKiQhMUdraShlTXNrJCEjOzckJCIybG1tbVp2T0wiISM6JCExJnBWTys0UzskISM7NyQkIjIsKytdMmdvUCIhIzokITFqZCI0dz8oKXAjISM7NyQkIjJLTCRlUjwqZlQiISM6JCExb2cmKVFWZEFCISM7NyQkIjItKytdKUh4ZTkhIzokITJ3XS1hPmEkZj4hIzw3JCQiMm1tO0ghby0qXCIhIzokITIvInBQYVQpM20iISM8NyQkIjIsK103ay42YSIhIzokITIqKSpcKkcnZU8hUiIhIzw3JCQiMm1tbTtXVEFlIiEjOiQhMnhtJ1FMIVtKOyIhIzw3JCQiLUQxKjNgaSIhIzUkITF2Kj1ZKWVuMCcqISM8NyQkIjJOTExMKnp5bTshIzokITF2Kj5TZzloJnohIzw3JCQiMk9MTDNOMSM0PCEjOiQhMTwvKT1oTmlgJyEjPDckJCIxbm0iSFl0N3YiISM5JCExLyVvaCpwKillYCEjPDckJCIyLSsrK3hHKip5IiEjOiQhMSI0JjNXJz46WCUhIzw3JCQiMW5tOzlAQk09ISM5JCExKlsuM0Jvb2UkISM8NyQkIjJPTExMYmRRKD0hIzokITIlXHJGemlQW0ghIz03JCQiLURPbDU7PiEjNSQhMm1PKT5fKnBiUSMhIz03JCQiMi0rXVA/V2wmPiEjOiQhMmNFTShmJyl5VT4hIz03JCQiJCsjISIiJCExdiNcIil5JClSYiIhIzwtJSZDT0xPUkc2JiUkUkdCRyQiIiEhIiIkIiIhISIiJCIiISEiIi0lJVZJRVdHNiQ7JCIiISEiIiQiJCsjISIiJShERUZBVUxURy0lK0FYRVNMQUJFTFNHNictSSNtaUc2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkc2IjY1USJyNiIvJSdmYW1pbHlHUSE2Ii8lJXNpemVHUSMxMDYiLyUlYm9sZEdRJmZhbHNlNiIvJSdpdGFsaWNHUSV0cnVlNiIvJSp1bmRlcmxpbmVHUSZmYWxzZTYiLyUqc3Vic2NyaXB0R1EmZmFsc2U2Ii8lLHN1cGVyc2NyaXB0R1EmZmFsc2U2Ii8lK2ZvcmVncm91bmRHUShbMCwwLDBdNiIvJStiYWNrZ3JvdW5kR1EuWzI1NSwyNTUsMjU1XTYiLyUnb3BhcXVlR1EmZmFsc2U2Ii8lK2V4ZWN1dGFibGVHUSZmYWxzZTYiLyUpcmVhZG9ubHlHUSZmYWxzZTYiLyUpY29tcG9zZWRHUSZmYWxzZTYiLyUqY29udmVydGVkR1EmZmFsc2U2Ii8lK2ltc2VsZWN0ZWRHUSZmYWxzZTYiLyUscGxhY2Vob2xkZXJHUSZmYWxzZTYiLyU2c2VsZWN0aW9uLXBsYWNlaG9sZGVyR1EmZmFsc2U2Ii8lLG1hdGh2YXJpYW50R1EnaXRhbGljNiJRITYiLSUlRk9OVEc2JSUoREVGQVVMVEclKERFRkFVTFRHIiM1JStIT1JJWk9OVEFMRyUrSE9SSVpPTlRBTEctJSVST09URzYnLSUpQk9VTkRTX1hHNiMkIiRxJiEiIi0lKUJPVU5EU19ZRzYjJCIjXSEiIi0lLUJPVU5EU19XSURUSEc2IyQiJXFrISIiLSUuQk9VTkRTX0hFSUdIVEc2IyQiJT9QISIiLSUpQ0hJTERSRU5HNiI=

int(uL(1,0,1.5)^2,r=0..infinity);

JCIvZiNmI2YjZiNmISM5

This basis is orthogonal, but not yet normalized. Orthogonality is demonstrated below. int(uL(1,0,1.5)*uL(2,0,1.5),r=0..infinity);

JCIiIUYj

# Define a basis with fixed parameter beta: phi:=n->uL(n,0,3/2);

Zio2I0kibkc2IkYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLUkjdUxHRiU2JTkkIiIhIyIiJCIiI0YlRiVGJQ==

OL:=(m,n)->int(phi(n)*phi(m),r=0..infinity);

Zio2JEkibUc2IkkibkdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLUkkaW50RzYkJSpwcm90ZWN0ZWRHSShfc3lzbGliR0YlNiQqJi1JJHBoaUdGJTYjOSUiIiItRjI2IzkkRjUvSSJyR0YlOyIiIUkpaW5maW5pdHlHRi1GJUYlRiU=

OL(2,2);

IyIjOyIiKg==

# a few orthogonality checks: OL(1,2);

IiIh

OL(2,4);

IiIh

# Now define the exponential potential and its matrix elements: V0:=-2;

ISIj

r0:=4.0; # to be consistent with the potential used in PH5000L5.mw where phase shifts are calculated.

JCIjUyEiIg==

V:=V0*exp(-r/r0);

LCQtSSRleHBHNiQlKnByb3RlY3RlZEdJKF9zeXNsaWJHNiI2IywkSSJyR0YoJCEvKysrKysrRCEjOSEiIw==

VME:=(m,n)->int(phi(n)*V*phi(m),r=0..infinity);

Zio2JEkibUc2IkkibkdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLUkkaW50RzYkJSpwcm90ZWN0ZWRHSShfc3lzbGliR0YlNiQqKC1JJHBoaUdGJTYjOSUiIiJJIlZHRiVGNS1GMjYjOSRGNS9JInJHRiU7IiIhSSlpbmZpbml0eUdGLUYlRiVGJQ==

VME(1,1);

JCEvP0peb2JqdSEjOQ==

Had we defined the basis with a non-float (using 15/10 instead of 1.5) the numbers would be computed symbolically exactly. VME(2,1);

JCEvMyo+QW4nKT4kISM5

VME(3,1);

JCEvIilvMCNcIVIiKiEjOg==

VME(1,6);

JCEvdidSVGdiSyohIzw=

TME:=(m,n)->-1/2*int(phi(n)*diff(phi(m),r$2),r=0..infinity);

Zio2JEkibUc2IkkibkdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLCQtSSRpbnRHNiQlKnByb3RlY3RlZEdJKF9zeXNsaWJHRiU2JComLUkkcGhpR0YlNiM5JSIiIi1JJWRpZmZHRi42JC1GMzYjOSQtSSIkR0YuNiRJInJHRiUiIiNGNi9GQDsiIiFJKWluZmluaXR5R0YuIyEiIkZBRiVGJUYl

TME(1,1);

IyIiIiIiJw==

TME(1,3);

IyIiIiIiJA==

TME(1,2),TME(2,1);

NiQjIiIiIiIkRiM=

TME(10,10);

IyIkOigiIic=

TME(10,11);

IyIkJlEiIiQ=

The kinetic energy is not diagonal, but symmetric! Nb:=10;

IiM1

VM:=Matrix(Nb,Nb,shape=symmetric): HM:=Matrix(Nb,Nb,shape=symmetric): OM:=Matrix(Nb,Nb,shape=symmetric): for m from 1 to Nb do: for n from 1 to m do: VM[m,n]:=VME(m,n): OM[m,n]:=OL(m,n): od: od: #print(map(evalf,VM)); for m from 1 to Nb do: for n from 1 to m do: HM[m,n]:=evalf(VM[m,n])+evalf(TME(m,n)): od: od: E10:=(convert(sort(Re(Eigenvalues(HM,OM))),list));

NywkITMqKio0L3ZzeCk0NSEjPCQhM3c+dFBqTCM0J1whIz0kITMsXUdXXktQTkBGKCQhMyIqNHAmbzQqKltJJyEjPiQhM0hJNnMjSDoncCkpISNAJCIzKikqXGEkeWhDSFhGLSQiMygpelBGLE1LXjhGKCQiMykpNDZAeCN6SkMlRigkIjMtIT11XE0iR0U+RiUkIjMkKio0KCpbeSIqeT4iISM7

Nb:=20;

IiM/

VM:=Matrix(Nb,Nb,shape=symmetric): HM:=Matrix(Nb,Nb,shape=symmetric): OM:=Matrix(Nb,Nb,shape=symmetric): for m from 1 to Nb do: for n from 1 to m do: VM[m,n]:=VME(m,n): OM[m,n]:=OL(m,n): od: od: #print(map(evalf,VM)); for m from 1 to Nb do: for n from 1 to m do: HM[m,n]:=evalf(VM[m,n])+evalf(TME(m,n)): od: od: E20:=(convert(sort(Re(Eigenvalues(HM,OM))),list));

NzYkITMoKlI1KlJbeSk0NSEjPCQhMykqZmo4ZngvaFwhIz0kITMuSUFNUHlyTkBGKCQhMzNTaSo0KVFxM2ohIz4kITM1XW11XF0lUnEkISM/JCIzJypmWyQpWyx1eVxGMCQiMyoqZnMlXE5nc24iRi0kIjNIcVsvZz5aJ1wkRi0kIjMkKjRLWTF1O3BnRi0kIjNwPytmJUhDLWcqRi0kIjMtNVxkJ2U9PVciRigkIjMxK1MzYVk2MUBGKCQiMyIqNEUkKSo0QkYvJEYoJCIzN2cpZSU9SEc9V0YoJCIzXUltJkdjaSp6bEYoJCIzLiFSdmE/LChcNUYlJCIzJypwOl1YMjcwPkYlJCIzQyEpUS84Y0BtU0YlJCIzM11IYDU/VCEzIiEjOyQiMy1dejJPMnlMW0ZN

# We see that the first 4 negative-E eigenvalues are confirmed in their lower significant digits. # Problem 10.8 in Quantum Mechanics by Ghatak and Lokanathan, QC174.12 G485 2004 asks to transform the SE for the s-states of an exp. potential to the Bessel equation. # The exact s-state eigenvalues are obtained from: g:=sqrt(r0^2*abs(8*V0));

JCIvKysrKysrOyEjNw==

plot(BesselJ(nu,g),nu=0..30);

NiYtJSdDVVJWRVNHNiQ3XXM3JCQiIiEhIiIkIS9qJClSMioqWzwhIzk3JCQiLERjJ3lNOyEjNiQhLyczXS8iW1w5ISM5NyQkIitESmRwSyEjNSQhLy9CbkRRZTUhIzk3JCQiLUQxaydwMyUhIzckIS9OYm5HbW0kKSEjOjckJCIsdm9mViFcISM2JCEvRiI+VSlbQmchIzo3JCQiLXZvSHZAZCEjNyQhLyR6XyYzYSNmJCEjOjckJCIqRFkiUmwhIiokIS8uYTc5I0c2IiEjOjckJCIudm9IbDoneiEjOCQiLzdqJm96RD8kISM6NyQkIi12VlYpUlEqISM3JCIvQTx5eEpZdCEjOjckJCIvRDFSLmshMyIhIzgkIi9ELTh5Mzg2ISM5NyQkIi12VkEpR0EiISM2JCIvdGlkWC9SOSEjOTckJCIvdj1VIltHUSIhIzgkIi88NWBjI2ZzIiEjOTckJCIuRDEvOUdhIiEjNyQiLiNSJj5sWSI+ISM4NyQkIi9RJSkqKXB6QTshIzgkIi8vYGo+S3A+ISM5NyQkIi9EMVIqekZxIiEjOCQiL1wqUUc7bio+ISM5NyQkIi89bjg5eFU8ISM4JCIvKFxDZFMsKyMhIzk3JCQiLzdHKSlHdyN5IiEjOCQiLyIqcHZzdCcqPiEjOTckJCIvMSpHT2FGIz0hIzgkIi89TTI3YicpPiEjOTckJCIsdiRldWk9ISM1JCIvWnRTKlwncD4hIzk3JCQiLXYkNHNQLSMhIzYkIi8nM1A0KHpOPSEjOTckJCIqTil6JT0jISIpJCIvIiopPmFGVWciISM5NyQkIi1EMVkjZU0jISM2JCItZiU0JTQqRyIhIzc3JCQiLEQnMyZvXSMhIzUkIi8pNC9xI3AlMyohIzo3JCQiL0QxKnk2cm0jISM4JCIvLSEpKipbMWBbISM6NyQkIi5EY3JzdCNHISM3JCIvW3JRNXUwUyEjOzckJCIvdj1VT2ooKUghIzgkIS9sY3VYKXovJSEjOjckJCItdm9YKnk5JCEjNiQhLyQzdSNIaHYjKSEjOjckJCIsRCdwWidIJCEjNSQhL1R4VkNTIz0iISM5NyQkIi1EYyRmXVckISM2JCEvPWE4M0soWyIhIzk3JCQiKnZUT2YkISIpJCEvYmlkbSg0dCIhIzk3JCQiLXZWVEFVUCEjNiQhLyU0Rz5bWyE+ISM5NyQkIit2JFsiPlEhIiokIS9uXCplOWInPiEjOTckJCItRDFFMicqUSEjNiQhL0RKIWY+YisjISM5NyQkIix2JG8qSChSISM1JCEvbzheR2xDPyEjOTckJCItdm81IypcUyEjNiQhL0tPYlgjSC0jISM5NyQkIi1ESiZwUD8lISM2JCEvTndpOip6Jj4hIzk3JCQiLXYkKnpoZFYhIzYkIS8lPikqKT0mXCI9ISM5NyQkIi12JEhHbl4lISM2JCEuV2o5WURmIiEjODckJCItdiRmUWVuJSEjNiQhL0hfOnY3MDghIzk3JCQiLXYkKilbXCRbISM2JCEvWCcpcEssZScqISM6NyQkIi12JD5mUypcISM2JCEvVSE9RGFJKmUhIzo3JCQiK0QjZkU6JiEiKiQhLz4xKVJfaCM+ISM6NyQkIi1EYyNmN0omISM2JCIvPngkZlZqNCMhIzo3JCQiLHZHZilwYSEjNSQiLiNbeEJKP2chIzk3JCQiLXY9JGYlR2MhIzYkIi93TSlSODdxKiEjOjckJCIvREpYKilmInomISM4JCIvQFhMUXI0OCEjOTckJCIudj1kUVomZiEjNyQiL2ZlanB4KGYiISM5NyQkIi92Vik+eXk2JyEjOCQiL3ApUWcnb0M9ISM5NyQkIitEeSwiRychIiokIi8kbzNOKVwkKT4hIzk3JCQiL3YkNGNrR04nISM4JCIvUCcqUW1tST8hIzk3JCQiLnZvSDZaVSchIzckIi9GPVhBZ2o/ISM5NyQkIi9EIkcuZWxcJyEjOCQiL2whZU9FQTMjISM5NyQkIi12b1pTb2whIzYkIi9APUU8YyczIyEjOTckJCIuPm44R1ZnJyEjNyQiLz0hNEgwTTMjISM5NyQkIi92by86RFNtISM4JCIuQG9xSG4yIyEjODckJCIuY0UoWzx3bSEjNyQiLyNmWilcY20/ISM5NyQkIi5EMUMpNDduISM3JCIvendIbyVIMCMhIzk3JCQiL0Rjd1wlUnknISM4JCIvVk5lMV86PyEjOTckJCIsRHIiemJvISM1JCIvWSFlO1hbJz4hIzk3JCQiLkQiR2xiPHEhIzckIi90S0t2JFshPSEjOTckJCItdlY4S3pyISM2JCIvOWwlcGxwZSIhIzk3JCQiLnYkZmgzVHQhIzckIi9iJzNgJWU+OCEjOTckJCIrdjQmR10oISIqJCIvaztxVUU3NSEjOTckJCIsRCwhR2x3ISM1JCIvKCoqUlIzKFJuISM6NyQkIiowNHgjeSEiKSQiL3J4J2ZxJm9KISM6NyQkIix2M1EsKnohIzUkIS8iUXgib2MvWyEjOzckJCIrRHJjXyIpISIqJCEvQnM2diIqKjQlISM6NyQkIil3NDQkKSEiKCQhL0VEaGUvbHUhIzo3JCQiK3YhR2NZKSEiKiQhL1lZcyQ+RDEiISM5NyQkIipiZUBpKSEiKSQhL20rdDRrXTghIzk3JCQiK0QhKm95KCkhIiokIS9wKHlaVVhnIiEjOTckJCIudiRmJHlIMSohIzckIS8oKjNwND1pPiEjOTckJCItdiRwbnNNKiEjNiQhLzJRRyk0LDwjISM5NyQkIi5jRWR6PFYqISM3JCEvKFFjIltqLEEhIzk3JCQiL0RjXjlIOyYqISM4JCEvbnpybkA+QSEjOTckJCIucC9MLjNnKiEjNyQhL1RqXiRwSEEjISM5NyQkIi52JDRfSiZvKiEjNyQhLz0iPTwkNDhBISM5NyQkIi92PW4qUVYmKSohIzgkIS4iR0A4bWBAISM4NyQkIixERk9CKyIhIiokIS94J1xjbFMvIyEjOTckJCItREpMKDQuIiEjNSQhL0RONDAwYzwhIzk3JCQiKlI1J2Y1ISIoJCEvJFI7Nmk8TyIhIzk3JCQiL1FmVmlFdzUhIzckIS5AQG47XjQiISM4NyQkIi92PSg0QUg0IiEjNyQhLyVSWiQpPWEzKSEjOjckJCIvN3ldemQ0NiEjNyQhLik9SzNEIjMmISM5NyQkIi52ViFRQkU2ISM2JCEvSkdvQzApKj4hIzo3JCQiLzcuZHEoNDkiISM3JCIvTCwpejddXighIzs3JCQiL3ZvNC5zYjYhIzckIi8pKVtXdHAiWyQhIzo3JCQiL1FNaU5ZcTYhIzckIi95cmJlXGVoISM6NyQkIis6bz8mPSIhIikkIi8kNFI0QjV2KSEjOjckJCIvUSUpKlwjUSw3ISM3JCIvJVEjeikqZlk2ISM5NyQkIi92byU9ZXZAIiEjNyQiL2FYWlNgLDkhIzk3JCQiLzdgcFF0TDchIzckIi9LPkx0SlA7ISM5NyQkIi52VmI0Klw3ISM2JCIvdV9oIikqPSY9ISM5NyQkIi12M2RyIUciISM1JCIvSFxELSRvPiMhIzk3JCQiLkRKJz1fNjghIzYkIi87PS44MGFDISM5NyQkIi92VmBbbVY4ISM3JCIvKj0wKHkvSEUhIzk3JCQiLXZWeSFlUCIhIzUkIi88J3AiKXlTciMhIzk3JCQiL2Neb3I9JFEiISM3JCIvQCE9VUw9cyMhIzk3JCQiLzdHJFxtMFIiISM3JCIvN0skKkhaREYhIzk3JCQiL28vPWUlelIiISM3JCIvUSJlKj05REYhIzk3JCQiL0QiRzlEYFMiISM3JCItOlkqKik0cyMhIzc3JCQiL1FNI3okMz85ISM3JCIvT3NfXSU9cSMhIzk3JCQiLnY9V1VbViIhIzYkIi9iKXkpXEpwRSEjOTckJCIvRGNPM29tOSEjNyQiL2A8KHl1KGZEISM5NyQkIi1ESiM+JilcIiEjNSQiL1dxdydIdFMjISM5NyQkIi52JD46bWs6ISM2JCIuayUpZl52KyMhIzg3JCQiLkRjZFFBaSIhIzYkIi9cPTh0NEc7ISM5NyQkIi12dExVJW8iISM1JCIvRCQpKUdjJVI3ISM5NyQkIitiam1bPCEiKSQiL281SlJBUyopISM6NyQkIi12eWJeNj0hIzUkIi9ASHQoWzpDJyEjOjckJCIudlZWREIoPSEjNiQiLzdtI3A7M0UlISM6NyQkIi1ENlclKVI+ISM1JCIuUl1NdycqbyMhIzk3JCQiKzpLXis/ISIpJCIvPT9xJCk9RTwhIzo3JCQiLEQ2IUhsPyEiKiQiL0pURnAhUi8iISM6NyQkIi52JDR3KVI3IyEjNiQiL3pEJHoxbVgnISM7NyQkIix2WmYiKT0jISIqJCIvb2swRmtBUCEjOzckJCIudlY/UyZbQSEjNiQiLypIZTxMbzsjISM7NyQkIi52PVliO0ojISM2JCIvLExQSTAuNyEjOzckJCIsRDtpTFAjISIqJCIvMS5oKnBTaSchIzw3JCQiLnYkZkwnelYjISM2JCIvRVpvYD9xTSEjPDckJCIqKj49K0QhIigkIi80bSlcKWZDPSEjPDckJCItREUmNFFjIyEjNSQiLyJcKUhFU3EjKiEjPTckJCIudlY+NXBpIyEjNiQiL19iLDVoWlkhIz03JCQiK2JKKltvIyEiKSQiL3cjMydcI1xVIyEjPTckJCItRHIiWzh2IyEjNSQiL2ZucT9pSDYhIz03JCQiKkwneTVHISIoJCIvaDV0QFU4YyEjPjckJCIudlYhKWZUKEchIzYkIi9FajVdXz9FISM+NyQkIi5EY0k7WyRIISM2JCIvInp3eDNcQyIhIz43JCQiJCskISIiJCIvd0lrJ1FfXSYhIz8tJSZDT0xPUkc2JiUkUkdCRyQiIzUhIiIkIiIhISIiJCIiISEiIi0lJVZJRVdHNiQ7JCIiISEiIiQiJCskISIiJShERUZBVUxURy0lK0FYRVNMQUJFTFNHNictSSNtaUc2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkc2IjY1UScmIzk1Nzs2Ii8lJ2ZhbWlseUdRITYiLyUlc2l6ZUdRIzEwNiIvJSVib2xkR1EmZmFsc2U2Ii8lJ2l0YWxpY0dRJmZhbHNlNiIvJSp1bmRlcmxpbmVHUSZmYWxzZTYiLyUqc3Vic2NyaXB0R1EmZmFsc2U2Ii8lLHN1cGVyc2NyaXB0R1EmZmFsc2U2Ii8lK2ZvcmVncm91bmRHUShbMCwwLDBdNiIvJStiYWNrZ3JvdW5kR1EuWzI1NSwyNTUsMjU1XTYiLyUnb3BhcXVlR1EmZmFsc2U2Ii8lK2V4ZWN1dGFibGVHUSZmYWxzZTYiLyUpcmVhZG9ubHlHUSZmYWxzZTYiLyUpY29tcG9zZWRHUSZmYWxzZTYiLyUqY29udmVydGVkR1EmZmFsc2U2Ii8lK2ltc2VsZWN0ZWRHUSZmYWxzZTYiLyUscGxhY2Vob2xkZXJHUSZmYWxzZTYiLyU2c2VsZWN0aW9uLXBsYWNlaG9sZGVyR1EmZmFsc2U2Ii8lLG1hdGh2YXJpYW50R1Enbm9ybWFsNiJRITYiLSUlRk9OVEc2JSUoREVGQVVMVEclKERFRkFVTFRHIiM1JStIT1JJWk9OVEFMRyUrSE9SSVpPTlRBTEctJSVST09URzYnLSUpQk9VTkRTX1hHNiMkIiQhZiEiIi0lKUJPVU5EU19ZRzYjJCIjXSEiIi0lLUJPVU5EU19XSURUSEc2IyQiJSE0KiEiIi0lLkJPVU5EU19IRUlHSFRHNiMkIiUhKVEhIiItJSlDSElMRFJFTkc2Ig==

# Five bound states exist in this case: the roots wrt nu are to be found: nu0:=fsolve(BesselJ(nu,g),nu=4..5);

LUknZnNvbHZlRzYkJSpwcm90ZWN0ZWRHSShfc3lzbGliRzYiNiUtSShCZXNzZWxKR0YkNiRJI251R0YnJCIvKysrKysrOyEjN0YsOyIiJSIiJg==

E0:=-nu0^2/(8.)/r0^2; # this is the energy for given nu.

LCQqJC1JJ2Zzb2x2ZUc2JCUqcHJvdGVjdGVkR0koX3N5c2xpYkc2IjYlLUkoQmVzc2VsSkdGJjYkSSNudUdGKSQiLysrKysrKzshIzdGLjsiIiUiIiYiIiMkIS8rKysrXTd5ISM7

nu1:=fsolve(BesselJ(nu,g),nu=1.5..2);

LUknZnNvbHZlRzYkJSpwcm90ZWN0ZWRHSShfc3lzbGliRzYiNiUtSShCZXNzZWxKR0YkNiRJI251R0YnJCIvKysrKysrOyEjN0YsOyQiIzohIiIiIiM=

E1:=-evalf(nu1^2/(8.)/r0^2);

LCQqJC1JJ2Zzb2x2ZUc2JCUqcHJvdGVjdGVkR0koX3N5c2xpYkc2IjYlLUkoQmVzc2VsSkdGJjYkSSNudUdGKSQiLysrKysrKzshIzdGLjskIiM6ISIiIiIjRjYkIS8rKysrXTd5ISM7

nu2:=fsolve(BesselJ(nu,g),nu=5..6);

JCIvPE45eFxHXyEjOA==

E2:=-evalf(nu2^2/(8.)/r0^2);

JCEvVHdSeXJOQCEjOQ==

nu3:=fsolve(BesselJ(nu,g),nu=2..4);

JCIvQzpmI3o7JUchIzg=

E3:=-evalf(nu3^2/(8.)/r0^2);

JCEvWGE6UnEzaiEjOg==

nu4:=fsolve(BesselJ(nu,g),nu=0.5..1.5);

JCIvPzVsYjwvcCEjOQ==

E4:=-evalf(nu4^2/(8.)/r0^2);

JCEvPkUkUk1TcyQhIzs=

# We can study systematically for each of the 5 eigenvalues how well the numerical calculation is doing. # We optimize the calculation wrt the basis parameter. phi:=n->uL(n,0,1.0);

Zio2I0kibkc2IkYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLUkjdUxHRiU2JTkkIiIhJCIjNSEiIkYlRiVGJQ==

Diag:=proc(Nb,beta,N_E) local phi,VME,TME,OL,VM,HM,OM,m,n,res; phi:=n->uL(n,0,beta); VME:=(m,n)->int(phi(n)*V*phi(m),r=0..infinity); TME:=(m,n)->-1/2*int(phi(n)*diff(phi(m),r$2),r=0..infinity); OL:=(m,n)->int(phi(n)*phi(m),r=0..infinity); VM:=Matrix(Nb,Nb,shape=symmetric): HM:=Matrix(Nb,Nb,shape=symmetric): OM:=Matrix(Nb,Nb,shape=symmetric): for m from 1 to Nb do: for n from 1 to m do: VM[m,n]:=VME(m,n): OM[m,n]:=OL(m,n): od: od: for m from 1 to Nb do: for n from 1 to m do: HM[m,n]:=evalf(VM[m,n])+evalf(TME(m,n)): od: od: res:=(convert(sort(Re(Eigenvalues(HM,OM))),list)); [seq(res[j],j=1..N_E)]; end: Diag(5,1,5);

NyckITMhKipHKFItP08vNSEjPCQhMz5nd0JxYlpMWiEjPSQhMyYqcC4kZiswIVw8RigkITNHIVtGJD5ONmhQISM+JCIzOCFSUyhRMHdeXEYo

Diag(5,1/2,5);

NyckITNxUlZZIj1WaCEqKiEjPSQhMyIpSE82LiZlWkUlRiUkITMxU3M7ckopPkEiRiUkITMpKnB1YzMoUkZUIiEjPiQiMyUpcCJ5QCYpPSEzRyEjPw==

# A smaller beta-value results in worse (higher) energies for all states, except the last one (which is still coming out as unbound). Note this 5-by-5 calculation has only 5 eigenvalues! j:=4; # pick the state for which the analysis is carried out:

IiIl

LL5:=[]: for i from 1 to 20 do: beta:=3/10*i; LL5:=[op(LL5),[beta,Diag(5,beta,5)[j]]]: od: PL1:=plot(LL5): beta:='beta': PL0:=plot(E||(j-1),beta=0..2,color=green): plots[display](PL0,PL1,view=[0..2,-0.1..0.2]);

NictJSdDVVJWRVNHNiQ3UzckJCIiISEiIiQhL1hhOlJxM2ohIzo3JCQiMUxMTEwzVmZWISM8JCEvWGE6UnEzaiEjOjckJCIxbm1tIkhbRDopISM8JCEvWGE6UnEzaiEjOjckJCIyTExMJGUwJD1DIiEjPCQhL1hhOlJxM2ohIzo3JCQiMkxMTCQzUkJyOyEjPCQhL1hhOlJxM2ohIzo3JCQiMmxtbSJ6amYpNCMhIzwkIS9YYTpScTNqISM6NyQkIjJLTExlNDtbXCMhIzwkIS9YYTpScTNqISM6NyQkIjIoKioqKlxpJ3ldIUghIzwkIS9YYTpScTNqISM6NyQkIjFMTCRlenMkSEwhIzskIS9YYTpScTNqISM6NyQkIjImKioqKlw3aUlfUCEjPCQhL1hhOlJxM2ohIzo3JCQiMW5tbTtfTSg9JSEjOyQhL1hhOlJxM2ohIzo3JCQiMk1MTCQzeV9xWCEjPCQhL1hhOlJxM2ohIzo3JCQiKmwrPismISIqJCEvWGE6UnEzaiEjOjckJCIqdlddViYhIiokIS9YYTpScTNqISM6NyQkIipOZkMmZSEiKiQhL1hhOlJxM2ohIzo3JCQiMUxMJGV6NjpCJyEjOyQhL1hhOlJxM2ohIzo3JCQiMW1tbTs9QyNvJyEjOyQhL1hhOlJxM2ohIzo3JCQiMW1tbW0jcFMxKCEjOyQhL1hhOlJxM2ohIzo3JCQiLERPRCMzdiEjNiQhL1hhOlJxM2ohIzo3JCQiMW1tbW0oeTgheiEjOyQhL1hhOlJxM2ohIzo3JCQiLERPSUZMKSEjNiQhL1hhOlJxM2ohIzo3JCQiLHYzek11KSEjNiQhL1hhOlJxM2ohIzo3JCQiMW5tbSJIXz88KiEjOyQhL1hhOlJxM2ohIzo3JCQiMW5tO3ppaGwmKiEjOyQhL1hhOlJxM2ohIzo3JCQiMUxMTDMjRywqKiohIzskIS9YYTpScTNqISM6NyQkIjJLTCRlenc1VjUhIzskIS9YYTpScTNqISM6NyQkIi12JFEjXCIzIiEjNiQhL1hhOlJxM2ohIzo3JCQiMktMTGUiKltINyIhIzskIS9YYTpScTNqISM6NyQkIipkeGQ7IiEiKSQhL1hhOlJxM2ohIzo3JCQiMikqKioqXF9xbjI3ISM7JCEvWGE6UnEzaiEjOjckJCIyKSoqKlxpJnBAWzchIzskIS9YYTpScTNqISM6NyQkIjIpKioqKlwyJ0hLSCIhIzskIS9YYTpScTNqISM6NyQkIjJsbW1tWnZPTCIhIzskIS9YYTpScTNqISM6NyQkIit2KydvUCIhIiokIS9YYTpScTNqISM6NyQkIjJLTCRlUjwqZlQiISM7JCEvWGE6UnEzaiEjOjckJCIrJilIeGU5ISIqJCEvWGE6UnEzaiEjOjckJCIybG07SCFvLSpcIiEjOyQhL1hhOlJxM2ohIzo3JCQiMioqKipcN2suNmEiISM7JCEvWGE6UnEzaiEjOjckJCIybW1tO1dUQWUiISM7JCEvWGE6UnEzaiEjOjckJCItRDEqM2BpIiEjNiQhL1hhOlJxM2ohIzo3JCQiMk1MTEwqenltOyEjOyQhL1hhOlJxM2ohIzo3JCQiMkxMTDNOMSM0PCEjOyQhL1hhOlJxM2ohIzo3JCQiMm1tO0hZdDd2IiEjOyQhL1hhOlJxM2ohIzo3JCQiKnhHKip5IiEiKSQhL1hhOlJxM2ohIzo3JCQiMm1tbVQ2S1UkPSEjOyQhL1hhOlJxM2ohIzo3JCQiMkxMTExiZFEoPSEjOyQhL1hhOlJxM2ohIzo3JCQiLURPbDU7PiEjNiQhL1hhOlJxM2ohIzo3JCQiLXYuVWFjPiEjNiQhL1hhOlJxM2ohIzo3JCQiIz8hIiIkIS9YYTpScTNqISM6LSUmQ09MT1JHNiYlJFJHQkckIiIhISIiJCIjNSEiIiQiIiEhIiItJSdDVVJWRVNHNiQ3NjckJCIiJCEiIiQiMEpzW0dFdV0nISM9NyQkIiInISIiJCEwdkJBR1knKlIjISM7NyQkIiIqISIiJCEwUCxvMiNHSkghIzs3JCQiIzchIiIkIS8memFZcyRIYSEjOjckJCIjOiEiIiQhMEs0WUBpcSVSISM7NyQkIiM9ISIiJCIwLDdjUjdQPiIhIzo3JCQiI0AhIiIkIjAtNys5KmVvVSEjOjckJCIjQyEiIiQiMCY0aFUmKnpJJSkhIzo3JCQiI0YhIiIkIjBMMyVbJz4mXDghIzk3JCQiI0khIiIkIjBUc0lsU3ckPiEjOTckJCIjTCEiIiQiMFZHPW85SGcjISM5NyQkIiNPISIiJCIwTDdGTG5FTSQhIzk3JCQiI1IhIiIkIjAlZShScSo9YlQhIzk3JCQiI1UhIiIkIi8lKip5YkIkUl0hIzg3JCQiI1ghIiIkIjB1VjxgWlUqZiEjOTckJCIjWyEiIiQiMFJSUi5lJD5xISM5NyQkIiNeISIiJCIvQSdRaCk+OSIpISM4NyQkIiNhISIiJCIweWoxKmVUeSMqISM5NyQkIiNkISIiJCIwWlNzM3Q2MCIhIzg3JCQiI2chIiIkIjB5M0AnPlIiPSIhIzgtJSZDT0xPUkc2JiUkUkdCRyQiIzUhIiIkIiIhISIiJCIiISEiIi0lJVZJRVdHNiQ7JCIiISEiIiQiIz8hIiI7JCEiIiEiIiQiIiMhIiItJStBWEVTTEFCRUxTRzYnLUkjbWlHNiMvSSttb2R1bGVuYW1lRzYiSSxUeXBlc2V0dGluZ0dJKF9zeXNsaWJHNiI2NVEnJiM5NDY7NiIvJSdmYW1pbHlHUShERUZBVUxUNiIvJSVzaXplR1EjMTA2Ii8lJWJvbGRHUSZmYWxzZTYiLyUnaXRhbGljR1EmZmFsc2U2Ii8lKnVuZGVybGluZUdRJmZhbHNlNiIvJSpzdWJzY3JpcHRHUSZmYWxzZTYiLyUsc3VwZXJzY3JpcHRHUSZmYWxzZTYiLyUrZm9yZWdyb3VuZEdRKFswLDAsMF02Ii8lK2JhY2tncm91bmRHUS5bMjU1LDI1NSwyNTVdNiIvJSdvcGFxdWVHUSZmYWxzZTYiLyUrZXhlY3V0YWJsZUdRJmZhbHNlNiIvJSlyZWFkb25seUdRJmZhbHNlNiIvJSljb21wb3NlZEdRJmZhbHNlNiIvJSpjb252ZXJ0ZWRHUSZmYWxzZTYiLyUraW1zZWxlY3RlZEdRJmZhbHNlNiIvJSxwbGFjZWhvbGRlckdRJmZhbHNlNiIvJTZzZWxlY3Rpb24tcGxhY2Vob2xkZXJHUSZmYWxzZTYiLyUsbWF0aHZhcmlhbnRHUSdub3JtYWw2IlEhNiItJSVGT05URzYlJShERUZBVUxURyUoREVGQVVMVEciIzUlK0hPUklaT05UQUxHJStIT1JJWk9OVEFMRy0lJVJPT1RHNictJSlCT1VORFNfWEc2IyQiJHEmISIiLSUpQk9VTkRTX1lHNiMkIiQ/IiEiIi0lLUJPVU5EU19XSURUSEc2IyQiJTVwISIiLSUuQk9VTkRTX0hFSUdIVEc2IyQiJV1QISIiLSUpQ0hJTERSRU5HNiI=

j;#

IiIl

#LL5; # Exercise: find the minimum value from the 5-by-5 calculation by optimizing wrt beta; how big is the relative error? # Try 10-by-10 calculations: LL10:=[]: for i from 1 to 5 do: beta:=1+1/10*i; LL10:=[op(LL10),[beta,Diag(10,beta,5)[j]]]: od: PL2:=plot(LL10,color=blue,style=point,symbol=cross): plots[display](PL0,PL1,PL2,view=[0..2,-0.065..-0.035]);

NigtJSdDVVJWRVNHNiQ3UzckJCIiISEiIiQhL1hhOlJxM2ohIzo3JCQiMUxMTEwzVmZWISM8JCEvWGE6UnEzaiEjOjckJCIxbm1tIkhbRDopISM8JCEvWGE6UnEzaiEjOjckJCIyTExMJGUwJD1DIiEjPCQhL1hhOlJxM2ohIzo3JCQiMkxMTCQzUkJyOyEjPCQhL1hhOlJxM2ohIzo3JCQiMmxtbSJ6amYpNCMhIzwkIS9YYTpScTNqISM6NyQkIjJLTExlNDtbXCMhIzwkIS9YYTpScTNqISM6NyQkIjIoKioqKlxpJ3ldIUghIzwkIS9YYTpScTNqISM6NyQkIjFMTCRlenMkSEwhIzskIS9YYTpScTNqISM6NyQkIjImKioqKlw3aUlfUCEjPCQhL1hhOlJxM2ohIzo3JCQiMW5tbTtfTSg9JSEjOyQhL1hhOlJxM2ohIzo3JCQiMk1MTCQzeV9xWCEjPCQhL1hhOlJxM2ohIzo3JCQiKmwrPismISIqJCEvWGE6UnEzaiEjOjckJCIqdlddViYhIiokIS9YYTpScTNqISM6NyQkIipOZkMmZSEiKiQhL1hhOlJxM2ohIzo3JCQiMUxMJGV6NjpCJyEjOyQhL1hhOlJxM2ohIzo3JCQiMW1tbTs9QyNvJyEjOyQhL1hhOlJxM2ohIzo3JCQiMW1tbW0jcFMxKCEjOyQhL1hhOlJxM2ohIzo3JCQiLERPRCMzdiEjNiQhL1hhOlJxM2ohIzo3JCQiMW1tbW0oeTgheiEjOyQhL1hhOlJxM2ohIzo3JCQiLERPSUZMKSEjNiQhL1hhOlJxM2ohIzo3JCQiLHYzek11KSEjNiQhL1hhOlJxM2ohIzo3JCQiMW5tbSJIXz88KiEjOyQhL1hhOlJxM2ohIzo3JCQiMW5tO3ppaGwmKiEjOyQhL1hhOlJxM2ohIzo3JCQiMUxMTDMjRywqKiohIzskIS9YYTpScTNqISM6NyQkIjJLTCRlenc1VjUhIzskIS9YYTpScTNqISM6NyQkIi12JFEjXCIzIiEjNiQhL1hhOlJxM2ohIzo3JCQiMktMTGUiKltINyIhIzskIS9YYTpScTNqISM6NyQkIipkeGQ7IiEiKSQhL1hhOlJxM2ohIzo3JCQiMikqKioqXF9xbjI3ISM7JCEvWGE6UnEzaiEjOjckJCIyKSoqKlxpJnBAWzchIzskIS9YYTpScTNqISM6NyQkIjIpKioqKlwyJ0hLSCIhIzskIS9YYTpScTNqISM6NyQkIjJsbW1tWnZPTCIhIzskIS9YYTpScTNqISM6NyQkIit2KydvUCIhIiokIS9YYTpScTNqISM6NyQkIjJLTCRlUjwqZlQiISM7JCEvWGE6UnEzaiEjOjckJCIrJilIeGU5ISIqJCEvWGE6UnEzaiEjOjckJCIybG07SCFvLSpcIiEjOyQhL1hhOlJxM2ohIzo3JCQiMioqKipcN2suNmEiISM7JCEvWGE6UnEzaiEjOjckJCIybW1tO1dUQWUiISM7JCEvWGE6UnEzaiEjOjckJCItRDEqM2BpIiEjNiQhL1hhOlJxM2ohIzo3JCQiMk1MTEwqenltOyEjOyQhL1hhOlJxM2ohIzo3JCQiMkxMTDNOMSM0PCEjOyQhL1hhOlJxM2ohIzo3JCQiMm1tO0hZdDd2IiEjOyQhL1hhOlJxM2ohIzo3JCQiKnhHKip5IiEiKSQhL1hhOlJxM2ohIzo3JCQiMm1tbVQ2S1UkPSEjOyQhL1hhOlJxM2ohIzo3JCQiMkxMTExiZFEoPSEjOyQhL1hhOlJxM2ohIzo3JCQiLURPbDU7PiEjNiQhL1hhOlJxM2ohIzo3JCQiLXYuVWFjPiEjNiQhL1hhOlJxM2ohIzo3JCQiIz8hIiIkIS9YYTpScTNqISM6LSUmQ09MT1JHNiYlJFJHQkckIiIhISIiJCIjNSEiIiQiIiEhIiItJSdDVVJWRVNHNiQ3NjckJCIiJCEiIiQiMEpzW0dFdV0nISM9NyQkIiInISIiJCEwdkJBR1knKlIjISM7NyQkIiIqISIiJCEwUCxvMiNHSkghIzs3JCQiIzchIiIkIS8memFZcyRIYSEjOjckJCIjOiEiIiQhMEs0WUBpcSVSISM7NyQkIiM9ISIiJCIwLDdjUjdQPiIhIzo3JCQiI0AhIiIkIjAtNys5KmVvVSEjOjckJCIjQyEiIiQiMCY0aFUmKnpJJSkhIzo3JCQiI0YhIiIkIjBMMyVbJz4mXDghIzk3JCQiI0khIiIkIjBUc0lsU3ckPiEjOTckJCIjTCEiIiQiMFZHPW85SGcjISM5NyQkIiNPISIiJCIwTDdGTG5FTSQhIzk3JCQiI1IhIiIkIjAlZShScSo9YlQhIzk3JCQiI1UhIiIkIi8lKip5YkIkUl0hIzg3JCQiI1ghIiIkIjB1VjxgWlUqZiEjOTckJCIjWyEiIiQiMFJSUi5lJD5xISM5NyQkIiNeISIiJCIvQSdRaCk+OSIpISM4NyQkIiNhISIiJCIweWoxKmVUeSMqISM5NyQkIiNkISIiJCIwWlNzM3Q2MCIhIzg3JCQiI2chIiIkIjB5M0AnPlIiPSIhIzgtJSZDT0xPUkc2JiUkUkdCRyQiIzUhIiIkIiIhISIiJCIiISEiIi0lJ0NVUlZFU0c2JjcnNyQkIiM2ISIiJCEvdy1PJXklemkhIzo3JCQiIzchIiIkITBMel5lYkFJJyEjOzckJCIjOCEiIiQhMEpaMCFRbixqISM7NyQkIiM5ISIiJCEwS1ZaZCkpPkknISM7NyQkIiM6ISIiJCEwInAmbzQqKltJJyEjOy0lJkNPTE9SRzYmJSRSR0JHJCIiISEiIiQiIiEhIiIkIiM1ISIiLSUnU1lNQk9MRzYjJSZDUk9TU0ctJSZTVFlMRUc2IyUmUE9JTlRHLSUlVklFV0c2JDskIiIhISIiJCIjPyEiIjskISNsISIkJCEjTiEiJC0lK0FYRVNMQUJFTFNHNictSSZtZnJhY0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkc2IjYqLUkjbW5HNiMvSSttb2R1bGVuYW1lRzYiSSxUeXBlc2V0dGluZ0dJKF9zeXNsaWJHNiI2NVEiMzYiLyUnZmFtaWx5R1EoREVGQVVMVDYiLyUlc2l6ZUdRIzEwNiIvJSVib2xkR1EmZmFsc2U2Ii8lJ2l0YWxpY0dRJmZhbHNlNiIvJSp1bmRlcmxpbmVHUSZmYWxzZTYiLyUqc3Vic2NyaXB0R1EmZmFsc2U2Ii8lLHN1cGVyc2NyaXB0R1EmZmFsc2U2Ii8lK2ZvcmVncm91bmRHUShbMCwwLDBdNiIvJStiYWNrZ3JvdW5kR1EuWzI1NSwyNTUsMjU1XTYiLyUnb3BhcXVlR1EmZmFsc2U2Ii8lK2V4ZWN1dGFibGVHUSZmYWxzZTYiLyUpcmVhZG9ubHlHUSZmYWxzZTYiLyUpY29tcG9zZWRHUSZmYWxzZTYiLyUqY29udmVydGVkR1EmZmFsc2U2Ii8lK2ltc2VsZWN0ZWRHUSZmYWxzZTYiLyUscGxhY2Vob2xkZXJHUSZmYWxzZTYiLyU2c2VsZWN0aW9uLXBsYWNlaG9sZGVyR1EmZmFsc2U2Ii8lLG1hdGh2YXJpYW50R1Enbm9ybWFsNiItSSNtbkc2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkc2IjY1USIyNiIvJSdmYW1pbHlHUShERUZBVUxUNiIvJSVzaXplR1EjMTA2Ii8lJWJvbGRHUSZmYWxzZTYiLyUnaXRhbGljR1EmZmFsc2U2Ii8lKnVuZGVybGluZUdRJmZhbHNlNiIvJSpzdWJzY3JpcHRHUSZmYWxzZTYiLyUsc3VwZXJzY3JpcHRHUSZmYWxzZTYiLyUrZm9yZWdyb3VuZEdRKFswLDAsMF02Ii8lK2JhY2tncm91bmRHUS5bMjU1LDI1NSwyNTVdNiIvJSdvcGFxdWVHUSZmYWxzZTYiLyUrZXhlY3V0YWJsZUdRJmZhbHNlNiIvJSlyZWFkb25seUdRJmZhbHNlNiIvJSljb21wb3NlZEdRJmZhbHNlNiIvJSpjb252ZXJ0ZWRHUSZmYWxzZTYiLyUraW1zZWxlY3RlZEdRJmZhbHNlNiIvJSxwbGFjZWhvbGRlckdRJmZhbHNlNiIvJTZzZWxlY3Rpb24tcGxhY2Vob2xkZXJHUSZmYWxzZTYiLyUsbWF0aHZhcmlhbnRHUSdub3JtYWw2Ii8lLmxpbmV0aGlja25lc3NHUSIxNiIvJStkZW5vbWFsaWduR1EnY2VudGVyNiIvJSludW1hbGlnbkdRJ2NlbnRlcjYiLyUpYmV2ZWxsZWRHUSZmYWxzZTYiLyUrZm9yZWdyb3VuZEdRKFswLDAsMF02Ii8lK2JhY2tncm91bmRHUS5bMjU1LDI1NSwyNTVdNiJRITYiLSUlRk9OVEc2JSUoREVGQVVMVEclKERFRkFVTFRHIiM1JStIT1JJWk9OVEFMRyUrSE9SSVpPTlRBTEctJSVST09URzYnLSUpQk9VTkRTX1hHNiMkIiRJKCEiIi0lKUJPVU5EU19ZRzYjJCIkPychIiItJS1CT1VORFNfV0lEVEhHNiMkIiVJJykhIiItJS5CT1VORFNfSEVJR0hURzYjJCIlXUshIiItJSlDSElMRFJFTkc2Ig==

# we see that the 10-by-10 diagonalization is less sensitive to the basis parameter beta. # Exercise: track the relative error (compared to the exact answer) as a function of beta. # Now ask the questions: is there a bound state in the L=1 sector (p-wave)? Diag:=proc(Nb,beta,L,N_E) local phi,VME,TME,OL,VM,HM,OM,m,n,res; phi:=n->uL(n,L,beta); VME:=(m,n)->int(phi(n)*(V+L*(L+1)/2/r^2)*phi(m),r=0..infinity); TME:=(m,n)->-1/2*int(phi(n)*diff(phi(m),r$2),r=0..infinity); OL:=(m,n)->int(phi(n)*phi(m),r=0..infinity); VM:=Matrix(Nb,Nb,shape=symmetric): HM:=Matrix(Nb,Nb,shape=symmetric): OM:=Matrix(Nb,Nb,shape=symmetric): for m from 1 to Nb do: for n from 1 to m do: VM[m,n]:=VME(m,n): OM[m,n]:=OL(m,n): od: od: for m from 1 to Nb do: for n from 1 to m do: HM[m,n]:=evalf(VM[m,n])+evalf(TME(m,n)): od: od: res:=(convert(sort(Re(Eigenvalues(HM,OM))),list)); [seq(res[j],j=1..N_E)]; end: Diag(5,1,1,1);

NyMkITNsUkhZKHB6ISlmJyEjPQ==

Diag(5,2,1,1);

JCEzKikqXGlIVlhObCchIz0=

Diag(5,3,1,1);

JCEzJSpSb104LS9hbSEjPQ==

Diag(5,4,1,3);

NyUkITNUU0FwaUU0Vm0hIz0kITM1NXg6YTlkVj5GJSQiM19IaSp6OyxALSlGJQ==

Diag(10,4,1,3);

NyUkITNNZ0w8RUYvYW0hIz0kITNGU1ApKXpqNSY0JEYlJCEzTj5meTFLWiZlKSEjPg==

Diag(20,4,1,3);

NyUkITNwPnZgdUYvYW0hIz0kITN5PmY3LF8zKDQkRiUkITMrK3RaQjM7NTZGJQ==

Diag(30,4,1,4);

NyYkITNwUkFkdUYvYW0hIz0kITMkKkgibzg/JjMoNCRGJSQhMzAhKTRGPTlANTZGJSQhMy4hWyVvUDhUIXAiISM+

Diag(30,5,1,5);

NyckITNvejg1dUYvYW0hIz0kITNDPz8+LF8zKDQkRiUkITMwUyNlJmUlNC02IkYlJCEzKikqZTMqNCdIJiplIiEjPiQiM0tdWi1fNTJCTkYs

Diag(30,3,1,5);

NyckITNJSTRrdUYvYW0hIz0kITMyPyNlNz8mMyg0JEYlJCEzJSpmXCM9VTYtNiJGJSQhMykqKkg9ckVtLXEiISM+JCIzRCFlIW88XXJRbiEjPw==

Diag(30,2,1,5);

NyckITNCNU5ddUYvYW0hIz0kITMjKSoqPkosXzMoNCRGJSQhMyw/dyU9VTYtNiJGJSQhMzs1cXB6Q1crPCEjPiQiMy0/aUAoSCpcN0QhIz8=

# For L>0 we do not have an exact result to compare against. Thus, we need to be diligent when determining the best results from matrix diagonalization.