Implem dot prod with SpinOrbit'

RemiHelleboid · Sep 15, 2023 · 55a6634 · 55a6634
1 parent 4b58536
commit 55a6634
Show file tree

Hide file tree

Showing 4 changed files with 34 additions and 4 deletions.
diff --git a/src/EPP/Hamiltonian.cpp b/src/EPP/Hamiltonian.cpp
@@ -9,6 +9,7 @@
 
 #include "Constants.hpp"
 #include "NonLocalFunctional.hpp"
+#include "SpinOrbitFunctional.hpp"
 #include "Vector3D.h"
 
 namespace EmpiricalPseudopotential {
@@ -41,13 +42,13 @@ void Hamiltonian::SetConstantNonDiagonalMatrix() {
     }
 }
 
-void Hamiltonian::SetMatrix(const Vector3D<double>& k, bool add_non_local_correction) {
+void Hamiltonian::SetMatrix(const Vector3D<double>& k, bool add_non_local_correction, bool enable_soc) {
     const unsigned int basisSize       = static_cast<unsigned int>(m_basisVectors.size());
     const double       latticeConstant = m_material.get_lattice_constant_meter();
     constexpr double   one_eight       = 1.0 / 8.0;
     Vector3D<double>   tau{one_eight * latticeConstant, one_eight * latticeConstant, one_eight * latticeConstant};
     const double       fourier_factor = 2.0 * M_PI / latticeConstant;
-    matrix = m_constant_non_diagonal_matrix;
+    matrix                            = m_constant_non_diagonal_matrix;
     // diagonal elements
     const double diag_factor = pow(Constants::h_bar, 2) / (2.0 * Constants::m0 * Constants::q);
     for (unsigned int i = 0; i < basisSize; ++i) {
@@ -71,7 +72,17 @@ void Hamiltonian::SetMatrix(const Vector3D<double>& k, bool add_non_local_correc
             }
         }
     }
-} 
+
+    if (enable_soc) {
+        std::size_t      matrix_size = matrix.rows();
+        std::size_t new_matrix_size = 2 * matrix_size;
+        Eigen::MatrixXcd ZerosMat = Eigen::MatrixXcd::Zero(matrix_size, matrix_size);
+        Eigen::MatrixXcd new_matrix(new_matrix_size, new_matrix_size);
+        new_matrix << matrix, ZerosMat, ZerosMat, matrix;
+        matrix = new_matrix;
+
+    }
+}
 
 void Hamiltonian::Diagonalize(bool keep_eigenvectors) {
     solver.compute(matrix, keep_eigenvectors ? Eigen::ComputeEigenvectors : Eigen::EigenvaluesOnly);

diff --git a/src/EPP/Hamiltonian.h b/src/EPP/Hamiltonian.h
@@ -12,7 +12,7 @@ class Hamiltonian {
     Hamiltonian(const Material& material, const std::vector<Vector3D<int>>& basisVectors);
     void SetConstantNonDiagonalMatrix();
 
-    void SetMatrix(const Vector3D<double>& k, bool add_non_local_correction = false);
+    void SetMatrix(const Vector3D<double>& k, bool add_non_local_correction = false, bool enable_soc = false);
     void Diagonalize(bool keep_eigenvectors = false);
 
     const Eigen::VectorXd&  get_eigenvalues() const { return solver.eigenvalues(); }

diff --git a/src/EPP/SpinOrbitFunctional.cpp b/src/EPP/SpinOrbitFunctional.cpp
@@ -77,4 +77,15 @@ double SpinOrbitCorrection::compute_lambda_antisym(const Vector3D<double>& K, co
     return lambda_antisym;
 }
 
+Eigen::Matrix<std::complex<double>, 2, 2> SpinOrbitCorrection::compute_pauli_state_dot_product(Vector3D<double> myVect) {
+    using namespace std::complex_literals;
+    std::complex<double>  a00 = myVect.Z;
+    std::complex<double>  a01 = myVect.X - myVect.Y * 1i;
+    std::complex<double>  a10 = myVect.X + myVect.Y * 1i;
+    std::complex<double>  a11 = -myVect.Z;
+    Eigen::Matrix<std::complex<double>, 2, 2> res_matrix;
+    res_matrix << a00, a10, a01, a11;
+    return res_matrix;
+}
+
 }  // namespace EmpiricalPseudopotential
diff --git a/src/EPP/SpinOrbitFunctional.hpp b/src/EPP/SpinOrbitFunctional.hpp
@@ -16,6 +16,9 @@
 #include <system_error>
 #include <vector>
 
+#include <Eigen/Dense>
+#include <complex>
+
 #include "Constants.hpp"
 #include "Material.h"
 #include "SpinOrbitParameters.hpp"
@@ -55,6 +58,11 @@ class SpinOrbitCorrection {
     double compute_lambda_sym(const Vector3D<double>& K, const Vector3D<double>& Kp) const;
     double compute_lambda_antisym(const Vector3D<double>& K, const Vector3D<double>& Kp) const;
 
+    static Eigen::Matrix<std::complex<double>, 2, 2> compute_pauli_state_dot_product(Vector3D<double> a);
+
+
 };
 
+
+
 }  // namespace EmpiricalPseudopotential