シュレディンガー方程式を解いて、量子系の時間発展の様子を計算します。 量子系のハミルトニアンをとすると、系の波動関数の時間発展はシュレディンガー方程式に従います。
このシュレディンガー方程式を時間領域で差分化して、数値的に波動関数の時間発展を求めます。 時刻での波動関数と時刻での波動関数との間には
という関係が成り立ちます。 ハミルトニアンの指数関数を求めてもいいのですが(行列の指数関数を求めたい)、 ここではこの関係式を、
と近似し、からを 逐次的に求めてゆきます。ComplexMatrix::solve(ComplexColumnVector& b) メソッドを用いてこの方程式を解いています。 この近似法では確率保存則が保たれることが知られています。 http://www-cms.phys.s.u-tokyo.ac.jp/~naoki/CIPINTRO/CIP/wave.html のページの解説を参考にしました。詳しくはそちらを参照してください。
#ifndef _QuantumDynamics_h #define _QuantumDynamics_h #include <iostream> #include <vector> #include <octave/config.h> #include <octave/CColVector.h> #include <octave/CMatrix.h> #include <CDiagMatrix.h> using namespace std; //================================================== struct Record { double time; ComplexColumnVector phi; Record(double itime, ComplexColumnVector iphi) : time(itime), phi(iphi) {}; }; typedef vector<Record> Data; //================================================== class QuantumDynamics { ComplexMatrix H; // Hamiltonian ComplexColumnVector init_phi; // initival state double dt; double t_max; Data data; int skip; public: QuantumDynamics(const ComplexMatrix& iH, const ComplexColumnVector& iinit_phi, double idt, double it_max, int iskip = 1) : H(iH), init_phi(iinit_phi), dt(idt), t_max(it_max), skip(iskip) {}; void calc(); Data getData(); }; void QuantumDynamics::calc() { const Complex IU(0.0, 1.0); ComplexMatrix Unit(ComplexDiagMatrix(H.rows(),H.cols(),1.0)); ComplexMatrix lhMatrix = Unit + 0.5 * dt * IU * H; ComplexMatrix rhMatrix = Unit - 0.5 * dt * IU * H; ComplexColumnVector phi = init_phi; ComplexColumnVector rhs; int count = 0; double t = 0.0; while (t < t_max) { if (count % skip == 0) { data.push_back(Record(t, phi)); } rhs = rhMatrix * phi; phi = lhMatrix.solve(rhs); t += dt; count++; } } Data QuantumDynamics::getData() { return data; } #endif
#include "QuantumDynamics.h" using namespace std; void print_probabilities(const Data& data) { for (int i = 0; i < data.size(); i++) { cout << i << " " <<data[i].time << " "; for (int j = 0; j < data[i].phi.rows(); j++) { cout << norm(data[i].phi(j)) << " "; } cout << endl; } } int main() { // Hamiltonian double del = 0.0; double eps = 1.0; ComplexMatrix H(2,2); H(0,0) = 0.0; H(0,1) = eps; H(1,0) = eps; H(1,1) = del; // initial state ComplexColumnVector init_phi(2); init_phi(0) = 1.0; init_phi(1) = 0.0; double dt = 1.e-2; // calculation of dynamics QuantumDynamics qDyn(H, init_phi, dt, 20.0, 1); qDyn.calc(); print_probabilities(qDyn.getData()); }
サンプルプログラムでは、エネルギー準位が縮退した2つの状態間にトンネルがあるような系の時間発展を求めています。 下に各状態の占有確率の時間発展を示します。 自由粒子の運動を求めるような、系の状態数が大きく、ハミルトニアンがスパースになる場合には、ここで用いたComplexMatrix::solve()を使った方法は効率的でないと思われます。 反復解法などを使うべきでしょう。
&ref(): File not found: "wikilogo1.png" at page "量子系のダイナミクス"; |