FIXME: atomic units aren't there yet. But I need them now. Must fix t…

…his later.
cosurgi · Jul 7, 2015 · 9209a1e · 9209a1e
1 parent 943f987
commit 9209a1e
Show file tree

Hide file tree

Showing 2 changed files with 8 additions and 7 deletions.
diff --git a/pkg/qm/SchrodingerPropagator.cpp b/pkg/qm/SchrodingerPropagator.cpp
@@ -168,7 +168,7 @@ Real SchrodingerKosloffPropagator::eMax()
 		int rank = psiGlobal->psiGlobalTable.rank();
 		Real Ekin(0);
 		for(int dim=0 ; dim<rank ; dim++)
-			Ekin += std::pow(psiGlobal->kMax(dim)* 1/* FIXME: must be `hbar` here */,2)/(2 /*FIXME: must be mass here psi->m */);
+			Ekin += std::pow(psiGlobal->kMax(dim)* 1/* FIXME: must be `hbar` here */,2)/(2 *FIXMEatomowe_MASS/*FIXME: must be mass here psi->m */);
 		ret=std::max(ret, Ekin );
 	}
 	// FIXME                                                                                ↓ ?  bez sensu, że w obu to się nazywa psiGlobalTable ....
@@ -180,7 +180,7 @@ Real SchrodingerKosloffPropagator::eMax()
 void SchrodingerKosloffPropagator::calc_Hnorm_psi(const NDimTable<Complexr>& psi_0,NDimTable<Complexr>& psi_1,
 	/*FIXME - remove*/QMStateDiscrete* psi)
 {
-	Real mass(1); // FIXME - this shouldn't be here
+	Real mass(FIXMEatomowe_MASS); // FIXME - this shouldn't be here
 	Real dt=scene->dt;
 
 	Real R   = calcKosloffR(dt); // FIXME -  that's duplicate here, depends on dt !!
@@ -212,7 +212,7 @@ void SchrodingerKosloffPropagator::calc_Hnorm_psi(const NDimTable<Complexr>& psi
 	// previous loop was:   };
 
 	if(Vpsi.rank() != 0)
-		Vpsi    .multMult(psi_0,dt/(hbar*R));// ψᵥ: ψᵥ=(dt V ψ₀)/(ℏ R)
+		Vpsi    .multMult(psi_0,dt/(FIXMEatomowe_hbar*R));// ψᵥ: ψᵥ=(dt V ψ₀)/(ℏ R)
 
 //? NDimTable<Complexr> psi_0c(psi_0);
 //? psi_0c.shiftByHalf();
@@ -223,7 +223,7 @@ void SchrodingerKosloffPropagator::calc_Hnorm_psi(const NDimTable<Complexr>& psi
 //? psi_1.shiftByHalf();
 	psi_1   .IFFT();                     // ψ₁: ψ₁=       ℱ⁻¹(-k²ℱ(ψ₀))
 //? psi_1.shiftByHalf();
-	psi_1 *= dt*hbar/(R*2*mass);         // ψ₁: ψ₁=(dt ℏ² ℱ⁻¹(-k²ℱ(ψ₀)) )/(ℏ R 2 m)
+	psi_1 *= dt*FIXMEatomowe_hbar/(R*2*mass);         // ψ₁: ψ₁=(dt ℏ² ℱ⁻¹(-k²ℱ(ψ₀)) )/(ℏ R 2 m)
 	psi_1   .mult2Add(psi_0,(1+G/R));    // ψ₁: ψ₁=(dt ℏ² ℱ⁻¹(-k²ℱ(ψ₀)) )/(ℏ R 2 m) + (1+G/R)ψ₀
 
 	if(Vpsi.rank() != 0)

diff --git a/pkg/qm/SchrodingerPropagator.hpp b/pkg/qm/SchrodingerPropagator.hpp
@@ -72,8 +72,8 @@ class SchrodingerKosloffPropagator: public GlobalEngine
 		virtual void action();
 		Real eMin();
 		Real eMax();
-		Real calcKosloffR(Real dt) { return dt*(eMax() - eMin())/(2*hbar);}; // calculate R parameter in Kosloff method
-		Real calcKosloffG(Real dt) { return dt*eMin()/(2*hbar);};            // calculate G parameter in Kosloff method
+		Real calcKosloffR(Real dt) { return dt*(eMax() - eMin())/(2*FIXMEatomowe_hbar);}; // calculate R parameter in Kosloff method
+		Real calcKosloffG(Real dt) { return dt*eMin()/(2*FIXMEatomowe_hbar);};            // calculate G parameter in Kosloff method
 		// FIXME: all ak can be precalculated, only recalculate if scene->dt changes
 		// FIXME: same with get_full_potentialInteractionGlobal_psiGlobalTable() - it can precalculate, and recalculate only upon dirty is set.
 		Complexr calcAK(int k,Real R) { return std::pow(Mathr::I,k)*(2.0 - Real(k==0))*(boost::math::cyl_bessel_j(k,R));};
@@ -91,7 +91,8 @@ operator  Û=exp(-iℍ̂t/ħ), and is following: ψ=Ûψ. The wavefunction ψ
 here in SchrodingerKosloffPropagator it is calculated using Tal-Ezer and Kosloff approach \
 found in [TalEzer1984]_"
 			, // attributes, public variables
-			((Real    ,hbar,1               ,,"Planck's constant $h$ divided by $2\\pi$"))
+			((Real    ,FIXMEatomowe_hbar,1               ,,"Planck's constant $h$ divided by $2\\pi$"))
+			((Real    ,FIXMEatomowe_MASS,1               ,,"FIXME - should use mass of the particle"))
 			((int     ,steps     ,-1     ,,"Override automatic selection of number of steps in Chebyshev expansion."))
 			((bool    ,virialCheck,false ,,"Check energies using virial theorem (Coulomb potential ONLY - FIXME!!!!!!!!)."))
 			, // constructor