fix QR decomp for complex-valued matrix (#47)

Co-authored-by: Kazuki Komatsu <komatsu.kazuki.op@tut.jp>

source/kaleidic/lubeck2.d

@@ -701,7 +701,12 @@ struct QRResult(T)
         auto aScope = aSliced.lightScope.canonical;
         auto tauSliced = tau.as!T.rcslice;
         auto tauScope = tauSliced.lightScope;
-        orgqr!T(aScope, tauScope, work);
+
+        static if(is(T == double) || is(T == float))
+            orgqr!T(aScope, tauScope, work);
+        else
+            ungqr!T(aScope, tauScope, work);
+
         Q = aScope.transposed.as!T.rcslice;
     }
 }
@@ -752,6 +757,31 @@ unittest
     assert(equal!approxEqual(mtimes(res.Q, res.R), data));
 }
 
+pure nothrow
+unittest
+{
+    import mir.ndslice;
+    import mir.math;
+
+    auto data = mininitRcslice!(Complex!double)(3, 3);
+    data[] = [[12, -51,   4],
+              [ 6, 167, -68],
+              [-4,  24, -41]];
+
+    auto res = qr(data);
+    auto q = mininitRcslice!(Complex!double)(3, 3);
+    q[] = [[-6.0/7.0,   69.0/175.0, 58.0/175.0],
+           [-3.0/7.0, -158.0/175.0, -6.0/175.0],
+           [ 2.0/7.0,   -6.0/35.0 , 33.0/35.0 ]];
+    auto aE = function (Complex!double x, Complex!double y) => approxEqual(x, y, 0.00005, 0.00005);
+    auto r = mininitRcslice!(Complex!double)(3, 3);
+    r[] = [[-14,  -21,  14],
+           [  0, -175,  70],
+           [  0,    0, -35]];
+    assert(equal!approxEqual(mtimes(res.Q, res.R), data));
+}
+
+
 @safe pure @nogc
 EigenResult!(realType!T) eigen(T, SliceKind kind)(
     auto ref Slice!(const(T)*,2, kind) a,