From 6e61ad6500d76a43f10e3b7e1d7c83fad159b8ae Mon Sep 17 00:00:00 2001
From: Malik Scheifinger <malik.scheifinger@kit.edu>
Date: Tue, 13 Feb 2024 09:11:36 +0100
Subject: [PATCH] Added complex sparse matrix instantiations.

---
 source/lac/sparse_matrix.inst.in   | 50 +++++++++++++---
 tests/lac/complex_sparse_matrix.cc | 91 ++++++++++++++++++++++++++++++
 2 files changed, 132 insertions(+), 9 deletions(-)
 create mode 100644 tests/lac/complex_sparse_matrix.cc

diff --git a/source/lac/sparse_matrix.inst.in b/source/lac/sparse_matrix.inst.in
index 343ed03ba6..037beafbed 100644
--- a/source/lac/sparse_matrix.inst.in
+++ b/source/lac/sparse_matrix.inst.in
@@ -112,15 +112,6 @@ for (S1, S2, S3 : REAL_SCALARS; V1, V2 : DEAL_II_VEC_TEMPLATES)
     template void SparseMatrix<S1>::Tvmult_add(V1<S2> &, const V2<S3> &) const;
   }
 
-for (S1 : REAL_SCALARS; S2, S3 : COMPLEX_SCALARS;
-     V1, V2 : DEAL_II_VEC_TEMPLATES)
-  {
-    template void SparseMatrix<S1>::vmult(V1<S2> &, const V2<S3> &) const;
-    template void SparseMatrix<S1>::Tvmult(V1<S2> &, const V2<S3> &) const;
-    template void SparseMatrix<S1>::vmult_add(V1<S2> &, const V2<S3> &) const;
-    template void SparseMatrix<S1>::Tvmult_add(V1<S2> &, const V2<S3> &) const;
-  }
-
 for (S1 : REAL_SCALARS)
   {
     template void SparseMatrix<S1>::vmult(
@@ -149,7 +140,48 @@ for (S1, S2, S3 : REAL_SCALARS)
                                            const bool) const;
   }
 
+// mixed instantiations
+
+for (S1 : REAL_SCALARS; S2, S3 : COMPLEX_SCALARS;
+     V1, V2 : DEAL_II_VEC_TEMPLATES)
+  {
+    template void SparseMatrix<S1>::vmult(V1<S2> &, const V2<S3> &) const;
+    template void SparseMatrix<S1>::Tvmult(V1<S2> &, const V2<S3> &) const;
+    template void SparseMatrix<S1>::vmult_add(V1<S2> &, const V2<S3> &) const;
+    template void SparseMatrix<S1>::Tvmult_add(V1<S2> &, const V2<S3> &) const;
+  }
 
+for (S1, S2 : COMPLEX_SCALARS; S3 : REAL_SCALARS;
+     V1, V2 : DEAL_II_VEC_TEMPLATES)
+  {
+    template void SparseMatrix<S1>::vmult(V1<S2> &, const V2<S3> &) const;
+    template void SparseMatrix<S1>::Tvmult(V1<S2> &, const V2<S3> &) const;
+    template void SparseMatrix<S1>::vmult_add(V1<S2> &, const V2<S3> &) const;
+    template void SparseMatrix<S1>::Tvmult_add(V1<S2> &, const V2<S3> &) const;
+  }
+
+for (S1 : COMPLEX_SCALARS; S2 : REAL_SCALARS)
+  {
+    template SparseMatrix<S1> &SparseMatrix<S1>::copy_from<S2>(
+      const SparseMatrix<S2> &);
+
+    template void SparseMatrix<S1>::copy_from<S2>(const FullMatrix<S2> &);
+
+    template void SparseMatrix<S1>::add<S2>(const S1, const SparseMatrix<S2> &);
+
+    template void SparseMatrix<S1>::add<S2>(const size_type,
+                                            const size_type,
+                                            const size_type *,
+                                            const S2 *,
+                                            const bool,
+                                            const bool);
+
+    template void SparseMatrix<S1>::set<S2>(const size_type,
+                                            const size_type,
+                                            const size_type *,
+                                            const S2 *,
+                                            const bool);
+  }
 
 // complex instantiations
 
diff --git a/tests/lac/complex_sparse_matrix.cc b/tests/lac/complex_sparse_matrix.cc
new file mode 100644
index 0000000000..36d7a7c6f7
--- /dev/null
+++ b/tests/lac/complex_sparse_matrix.cc
@@ -0,0 +1,91 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Test SparseMatrix<complex<number>> copy_from, add from SparseMatrix<number>,
+// and vmult, vmult_add, ... for mixed real and complex vectors
+
+#include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/lac/vector.h>
+
+#include <complex>
+
+#include "../tests.h"
+
+
+template <typename number>
+void
+check_operations()
+{
+  using namespace std::complex_literals;
+
+  SparsityPattern pattern(1, 1, 1);
+  pattern.add(0, 0);
+  pattern.compress();
+  SparseMatrix<number>               real_matrix(pattern);
+  SparseMatrix<std::complex<number>> cmpl_matrix(pattern);
+
+  real_matrix.set(0, 0, 1.0);
+  deallog << 'real_matrix ' << real_matrix(0, 0) << std::endl;
+
+
+  cmpl_matrix.copy_from(real_matrix);
+  Assert((cmpl_matrix(0, 0).real() == real_matrix(0, 0)) &&
+           (cmpl_matrix(0, 0).imag() == 0.0),
+         ExcInternalError());
+  deallog << 'cmpl_matrix ' << cmpl_matrix(0, 0) << std::endl;
+
+
+  cmpl_matrix.add(0, 0, 1i);
+  Assert((cmpl_matrix(0, 0).real() == 1.0) && (cmpl_matrix(0, 0).imag() == 1.0),
+         ExcInternalError());
+  deallog << 'cmpl_matrix ' << cmpl_matrix(0, 0) << std::endl;
+
+  Vector<number>               real_vec({2.});
+  Vector<std::complex<number>> cmpl_vec({0.});
+
+  cmpl_matrix.vmult(cmpl_vec, real_vec);
+  Assert((cmpl_vec(0).real() == 2.0) && (cmpl_vec(0).imag() == 2.0),
+         ExcInternalError());
+  deallog << 'cmpl_vec ' << cmpl_vec(0) << std::endl;
+
+  cmpl_matrix.vmult_add(cmpl_vec, real_vec);
+  Assert((cmpl_vec(0).real() == 4.0) && (cmpl_vec(0).imag() == 4.0),
+         ExcInternalError());
+  deallog << 'cmpl_vec ' << cmpl_vec(0) << std::endl;
+
+  cmpl_vec.reinit(1);
+
+  cmpl_matrix.Tvmult(cmpl_vec, real_vec);
+  Assert((cmpl_vec(0).real() == 2.0) && (cmpl_vec(0).imag() == 2.0),
+         ExcInternalError());
+  deallog << 'cmpl_vec ' << cmpl_vec(0) << std::endl;
+
+  cmpl_matrix.Tvmult_add(cmpl_vec, real_vec);
+  Assert((cmpl_vec(0).real() == 4.0) && (cmpl_vec(0).imag() == 4.0),
+         ExcInternalError());
+  deallog << 'cmpl_vec ' << cmpl_vec(0) << std::endl;
+}
+
+int
+main()
+{
+  initlog();
+
+  deallog << std::setprecision(0) << std::fixed;
+
+  check_operations<float>();
+  check_operations<double>();
+}
-- 
2.39.5