4 template <
typename Float,
int writeX,
int writeY,
int writeZ,
int writeW,
int writeV,
typename SpinorX,
5 typename SpinorY,
typename SpinorZ,
typename SpinorW,
typename SpinorV,
typename Functor>
6 void genericBlas(SpinorX &
X, SpinorY &Y, SpinorZ &
Z, SpinorW &W, SpinorV &
V, Functor f)
10 for (
int x = 0; x < X.VolumeCB(); x++) {
11 for (
int s = 0;
s < X.Nspin();
s++) {
12 for (
int c = 0; c < X.Ncolor(); c++) {
13 complex<Float> X_(
X(
parity, x,
s, c));
14 complex<Float> Y_ = Y(
parity, x,
s, c);
15 complex<Float> Z_ =
Z(
parity, x,
s, c);
16 complex<Float> W_ = W(
parity, x,
s, c);
17 complex<Float> V_ =
V(
parity, x,
s, c);
18 f(X_, Y_, Z_, W_, V_);
20 if (writeY) Y(
parity, x,
s, c) = Y_;
22 if (writeW) W(
parity, x,
s, c) = W_;
30 template <
typename Float,
typename yFloat,
int nSpin,
int nColor,
QudaFieldOrder order,
int writeX,
int writeY,
31 int writeZ,
int writeW,
int writeV,
typename Functor>
33 ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v, Functor f)
35 colorspinor::FieldOrderCB<Float, nSpin, nColor, 1, order>
X(x),
Z(z), W(w);
36 colorspinor::FieldOrderCB<yFloat, nSpin, nColor, 1, order> Y(y),
V(v);
37 genericBlas<yFloat, writeX, writeY, writeZ, writeW, writeV>(
X, Y,
Z, W,
V, f);
40 template <
typename Float,
typename yFloat,
int nSpin,
QudaFieldOrder order,
int writeX,
int writeY,
int writeZ,
41 int writeW,
int writeV,
typename Functor>
43 ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v, Functor f)
45 if (x.Ncolor() == 3) {
46 genericBlas<Float, yFloat, nSpin, 3, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
47 }
else if (x.Ncolor() == 4) {
48 genericBlas<Float, yFloat, nSpin, 4, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
49 }
else if (x.Ncolor() == 6) {
50 genericBlas<Float, yFloat, nSpin, 6, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
51 }
else if (x.Ncolor() == 8) {
52 genericBlas<Float, yFloat, nSpin, 8, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
53 }
else if (x.Ncolor() == 12) {
54 genericBlas<Float, yFloat, nSpin, 12, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
55 }
else if (x.Ncolor() == 16) {
56 genericBlas<Float, yFloat, nSpin, 16, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
57 }
else if (x.Ncolor() == 20) {
58 genericBlas<Float, yFloat, nSpin, 20, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
59 }
else if (x.Ncolor() == 24) {
60 genericBlas<Float, yFloat, nSpin, 24, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
61 }
else if (x.Ncolor() == 32) {
62 genericBlas<Float, yFloat, nSpin, 32, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
64 errorQuda(
"nColor = %d not implemented", x.Ncolor());
68 template <
typename Float,
typename yFloat,
QudaFieldOrder order,
int writeX,
int writeY,
int writeZ,
int writeW,
69 int writeV,
typename Functor>
71 ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v, Functor f)
74 genericBlas<Float, yFloat, 4, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
75 }
else if (x.Nspin() == 2) {
76 genericBlas<Float, yFloat, 2, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
77 #ifdef GPU_STAGGERED_DIRAC 78 }
else if (x.Nspin() == 1) {
79 genericBlas<Float, yFloat, 1, order, writeX, writeY, writeZ, writeW, writeV, Functor>(x, y, z, w, v, f);
82 errorQuda(
"nSpin = %d not implemented", x.Nspin());
86 template <
typename Float,
typename yFloat,
int writeX,
int writeY,
int writeZ,
int writeW,
int writeV,
typename Functor>
88 ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v, Functor f)
91 genericBlas<Float, yFloat, QUDA_SPACE_SPIN_COLOR_FIELD_ORDER, writeX, writeY, writeZ, writeW, writeV, Functor>(
enum QudaFieldOrder_s QudaFieldOrder
void genericBlas(SpinorX &X, SpinorY &Y, SpinorZ &Z, SpinorW &W, SpinorV &V, Functor f)