QUDA  v1.1.0
A library for QCD on GPUs
copy_gauge_inc.cu
Go to the documentation of this file.
1 #include <gauge_field_order.h>
2 #include <copy_gauge_helper.cuh>
3 #include <instantiate.h>
4 
5 namespace quda {
6 
7  template <typename FloatOut, typename FloatIn, int length, typename InOrder>
8  void copyGauge(const InOrder &inOrder, const GaugeField &out, const GaugeField &in,
9  QudaFieldLocation location, FloatOut *Out, FloatOut **outGhost, int type) {
10  if (out.isNative()) {
11  // this overrides the check that the texture maps to the gauge
12  // pointer - this is safe here since it only occurs when running
13  // the copier on the host when we will not be using texture
14  // reads
15  const bool override = true;
16  if (out.Reconstruct() == QUDA_RECONSTRUCT_NO) {
17  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_NO>::type G;
18  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
19  } else if (out.Reconstruct() == QUDA_RECONSTRUCT_12) {
20 #if QUDA_RECONSTRUCT & 2
21  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_12>::type G;
22  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
23 #else
24  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-12", QUDA_RECONSTRUCT);
25 #endif
26  } else if (out.Reconstruct() == QUDA_RECONSTRUCT_8) {
27 #if QUDA_RECONSTRUCT & 1
28  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_8>::type G;
29  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
30 #else
31  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-8", QUDA_RECONSTRUCT);
32 #endif
33 #ifdef GPU_STAGGERED_DIRAC
34  } else if (out.Reconstruct() == QUDA_RECONSTRUCT_13) {
35 #if QUDA_RECONSTRUCT & 2
36  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_13>::type G;
37  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
38 #else
39  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-13", QUDA_RECONSTRUCT);
40 #endif
41  } else if (out.Reconstruct() == QUDA_RECONSTRUCT_9) {
42 #if QUDA_RECONSTRUCT & 1
43  if (out.StaggeredPhase() == QUDA_STAGGERED_PHASE_MILC) {
44  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_9, 18, QUDA_STAGGERED_PHASE_MILC>::type G;
45  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
46  } else if (out.StaggeredPhase() == QUDA_STAGGERED_PHASE_TIFR) {
47 #ifdef BUILD_TIFR_INTERFACE
48  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_9, 18, QUDA_STAGGERED_PHASE_TIFR>::type G;
49  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
50 #else
51  errorQuda("TIFR interface has not been built so TIFR phase type not enabled\n");
52 #endif
53  } else if (out.StaggeredPhase() == QUDA_STAGGERED_PHASE_NO) {
54  typedef typename gauge_mapper<FloatOut, QUDA_RECONSTRUCT_9>::type G;
55  copyGauge<FloatOut, FloatIn, length>(G(out, Out, outGhost, override), inOrder, out, in, location, type);
56  } else {
57  errorQuda("Staggered phase type %d not supported", out.StaggeredPhase());
58  }
59 #else
60  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-9", QUDA_RECONSTRUCT);
61 #endif
62 #endif // GPU_STAGGERED_DIRAC
63  } else {
64  errorQuda("Reconstruction %d and order %d not supported", out.Reconstruct(), out.Order());
65  }
66  } else if (out.Order() == QUDA_QDP_GAUGE_ORDER) {
67 
68 #ifdef BUILD_QDP_INTERFACE
69  copyGauge<FloatOut,FloatIn,length>
70  (QDPOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
71 #else
72  errorQuda("QDP interface has not been built\n");
73 #endif
74 
75  } else if (out.Order() == QUDA_QDPJIT_GAUGE_ORDER) {
76 
77 #ifdef BUILD_QDPJIT_INTERFACE
78  copyGauge<FloatOut,FloatIn,length>
79  (QDPJITOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
80 #else
81  errorQuda("QDPJIT interface has not been built\n");
82 #endif
83 
84  } else if (out.Order() == QUDA_CPS_WILSON_GAUGE_ORDER) {
85 
86 #ifdef BUILD_CPS_INTERFACE
87  copyGauge<FloatOut,FloatIn,length>
88  (CPSOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
89 #else
90  errorQuda("CPS interface has not been built\n");
91 #endif
92 
93  } else if (out.Order() == QUDA_MILC_GAUGE_ORDER) {
94 
95 #ifdef BUILD_MILC_INTERFACE
96  copyGauge<FloatOut,FloatIn,length>
97  (MILCOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
98 #else
99  errorQuda("MILC interface has not been built\n");
100 #endif
101 
102  } else if (out.Order() == QUDA_MILC_SITE_GAUGE_ORDER) {
103 
104 #ifdef BUILD_MILC_INTERFACE
105  copyGauge<FloatOut,FloatIn,length>
106  (MILCSiteOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
107 #else
108  errorQuda("MILC interface has not been built\n");
109 #endif
110 
111  } else if (out.Order() == QUDA_BQCD_GAUGE_ORDER) {
112 
113 #ifdef BUILD_BQCD_INTERFACE
114  copyGauge<FloatOut,FloatIn,length>
115  (BQCDOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
116 #else
117  errorQuda("BQCD interface has not been built\n");
118 #endif
119 
120  } else if (out.Order() == QUDA_TIFR_GAUGE_ORDER) {
121 
122 #ifdef BUILD_TIFR_INTERFACE
123  copyGauge<FloatOut,FloatIn,length>
124  (TIFROrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
125 #else
126  errorQuda("TIFR interface has not been built\n");
127 #endif
128 
129  } else if (out.Order() == QUDA_TIFR_PADDED_GAUGE_ORDER) {
130 
131 #ifdef BUILD_TIFR_INTERFACE
132  copyGauge<FloatOut,FloatIn,length>
133  (TIFRPaddedOrder<FloatOut,length>(out, Out, outGhost), inOrder, out, in, location, type);
134 #else
135  errorQuda("TIFR interface has not been built\n");
136 #endif
137 
138  } else {
139  errorQuda("Gauge field %d order not supported", out.Order());
140  }
141 
142  }
143 
144  template <typename FloatOut, typename FloatIn, int length>
145  void copyGauge(GaugeField &out, const GaugeField &in, QudaFieldLocation location, FloatOut *Out, FloatIn *In,
146  FloatOut **outGhost, FloatIn **inGhost, int type)
147  {
148 
149  // reconstruction only supported on FloatN fields currently
150  if (in.isNative()) {
151  // this overrides the check that the texture maps to the gauge
152  // pointer - this is safe here since it only occurs when running
153  // the copier on the host when we will not be using texture
154  // reads
155  const bool override = true;
156  if (in.Reconstruct() == QUDA_RECONSTRUCT_NO) {
157  typedef typename gauge_mapper<FloatIn, QUDA_RECONSTRUCT_NO>::type G;
158  copyGauge<FloatOut, FloatIn, length>(G(in, In, inGhost, override), out, in, location, Out, outGhost, type);
159  } else if (in.Reconstruct() == QUDA_RECONSTRUCT_12) {
160 #if QUDA_RECONSTRUCT & 2
161  typedef typename gauge_mapper<FloatIn,QUDA_RECONSTRUCT_12>::type G;
162  copyGauge<FloatOut,FloatIn,length> (G(in,In,inGhost,override), out, in, location, Out, outGhost, type);
163 #else
164  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-12", QUDA_RECONSTRUCT);
165 #endif
166  } else if (in.Reconstruct() == QUDA_RECONSTRUCT_8) {
167 #if QUDA_RECONSTRUCT & 1
168  typedef typename gauge_mapper<FloatIn,QUDA_RECONSTRUCT_8>::type G;
169  copyGauge<FloatOut,FloatIn,length> (G(in,In,inGhost,override), out, in, location, Out, outGhost, type);
170 #else
171  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-8", QUDA_RECONSTRUCT);
172 #endif
173 #ifdef GPU_STAGGERED_DIRAC
174  } else if (in.Reconstruct() == QUDA_RECONSTRUCT_13) {
175 #if QUDA_RECONSTRUCT & 2
176  typedef typename gauge_mapper<FloatIn,QUDA_RECONSTRUCT_13>::type G;
177  copyGauge<FloatOut,FloatIn,length> (G(in,In,inGhost,override), out, in, location, Out, outGhost, type);
178 #else
179  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-13", QUDA_RECONSTRUCT);
180 #endif
181  } else if (in.Reconstruct() == QUDA_RECONSTRUCT_9) {
182 #if QUDA_RECONSTRUCT & 1
183  if (in.StaggeredPhase() == QUDA_STAGGERED_PHASE_MILC) {
184  typedef typename gauge_mapper<FloatIn, QUDA_RECONSTRUCT_9, 18, QUDA_STAGGERED_PHASE_MILC>::type G;
185  copyGauge<FloatOut, FloatIn, length>(G(in, In, inGhost, override), out, in, location, Out, outGhost, type);
186  } else if (in.StaggeredPhase() == QUDA_STAGGERED_PHASE_NO) {
187  typedef typename gauge_mapper<FloatIn, QUDA_RECONSTRUCT_9>::type G;
188  copyGauge<FloatOut, FloatIn, length>(G(in, In, inGhost, override), out, in, location, Out, outGhost, type);
189  } else {
190  errorQuda("Staggered phase type %d not supported", in.StaggeredPhase());
191  }
192 #else
193  errorQuda("QUDA_RECONSTRUCT=%d does not enable reconstruct-9", QUDA_RECONSTRUCT);
194 #endif
195 #endif // GPU_STAGGERED_DIRAC
196  } else {
197  errorQuda("Reconstruction %d and order %d not supported", in.Reconstruct(), in.Order());
198  }
199  } else if (in.Order() == QUDA_QDP_GAUGE_ORDER) {
200 
201 #ifdef BUILD_QDP_INTERFACE
202  copyGauge<FloatOut,FloatIn,length>(QDPOrder<FloatIn,length>(in, In, inGhost),
203  out, in, location, Out, outGhost, type);
204 #else
205  errorQuda("QDP interface has not been built\n");
206 #endif
207 
208  } else if (in.Order() == QUDA_QDPJIT_GAUGE_ORDER) {
209 
210 #ifdef BUILD_QDPJIT_INTERFACE
211  copyGauge<FloatOut,FloatIn,length>(QDPJITOrder<FloatIn,length>(in, In, inGhost),
212  out, in, location, Out, outGhost, type);
213 #else
214  errorQuda("QDPJIT interface has not been built\n");
215 #endif
216 
217  } else if (in.Order() == QUDA_CPS_WILSON_GAUGE_ORDER) {
218 
219 #ifdef BUILD_CPS_INTERFACE
220  copyGauge<FloatOut,FloatIn,length>(CPSOrder<FloatIn,length>(in, In, inGhost),
221  out, in, location, Out, outGhost, type);
222 #else
223  errorQuda("CPS interface has not been built\n");
224 #endif
225 
226  } else if (in.Order() == QUDA_MILC_GAUGE_ORDER) {
227 
228 #ifdef BUILD_MILC_INTERFACE
229  copyGauge<FloatOut,FloatIn,length>(MILCOrder<FloatIn,length>(in, In, inGhost),
230  out, in, location, Out, outGhost, type);
231 #else
232  errorQuda("MILC interface has not been built\n");
233 #endif
234 
235  } else if (in.Order() == QUDA_MILC_SITE_GAUGE_ORDER) {
236 
237 #ifdef BUILD_MILC_INTERFACE
238  copyGauge<FloatOut,FloatIn,length>(MILCSiteOrder<FloatIn,length>(in, In, inGhost),
239  out, in, location, Out, outGhost, type);
240 #else
241  errorQuda("MILC interface has not been built\n");
242 #endif
243 
244  } else if (in.Order() == QUDA_BQCD_GAUGE_ORDER) {
245 
246 #ifdef BUILD_BQCD_INTERFACE
247  copyGauge<FloatOut,FloatIn,length>(BQCDOrder<FloatIn,length>(in, In, inGhost),
248  out, in, location, Out, outGhost, type);
249 #else
250  errorQuda("BQCD interface has not been built\n");
251 #endif
252 
253  } else if (in.Order() == QUDA_TIFR_GAUGE_ORDER) {
254 
255 #ifdef BUILD_TIFR_INTERFACE
256  copyGauge<FloatOut,FloatIn,length>(TIFROrder<FloatIn,length>(in, In, inGhost),
257  out, in, location, Out, outGhost, type);
258 #else
259  errorQuda("TIFR interface has not been built\n");
260 #endif
261 
262  } else if (in.Order() == QUDA_TIFR_PADDED_GAUGE_ORDER) {
263 
264 #ifdef BUILD_TIFR_INTERFACE
265  copyGauge<FloatOut,FloatIn,length>(TIFRPaddedOrder<FloatIn,length>(in, In, inGhost),
266  out, in, location, Out, outGhost, type);
267 #else
268  errorQuda("TIFR interface has not been built\n");
269 #endif
270 
271  } else {
272  errorQuda("Gauge field order %d not supported", in.Order());
273  }
274  }
275 
276  void checkMomOrder(const GaugeField &u);
277 
278  template <typename FloatOut, typename FloatIn, int length, typename Out, typename In, typename Arg>
279  void copyMom(Arg &arg, const GaugeField &out, const GaugeField &in, QudaFieldLocation location) {
280  CopyGauge<FloatOut,FloatIn,length, Arg> momCopier(arg, out, in, location);
281  momCopier.apply(0);
282  }
283 
284  template <typename FloatOut, typename FloatIn> struct GaugeCopy {
285  GaugeCopy(GaugeField &out, const GaugeField &in, QudaFieldLocation location, void *Out_, void *In_,
286  void **outGhost_, void **inGhost_, int type)
287  {
288  FloatOut *Out = reinterpret_cast<FloatOut*>(Out_);
289  FloatIn *In = reinterpret_cast<FloatIn*>(In_);
290  FloatOut **outGhost = reinterpret_cast<FloatOut**>(outGhost_);
291  FloatIn **inGhost = reinterpret_cast<FloatIn**>(inGhost_);
292 
293  if (in.Ncolor() != 3 && out.Ncolor() != 3) {
294  errorQuda("Unsupported number of colors; out.Nc=%d, in.Nc=%d", out.Ncolor(), in.Ncolor());
295  }
296 
297  if (out.Geometry() != in.Geometry()) {
298  errorQuda("Field geometries %d %d do not match", out.Geometry(), in.Geometry());
299  }
300 
301  if (in.LinkType() != QUDA_ASQTAD_MOM_LINKS && out.LinkType() != QUDA_ASQTAD_MOM_LINKS) {
302  // we are doing gauge field packing
303  copyGauge<FloatOut,FloatIn,18>(out, in, location, Out, In, outGhost, inGhost, type);
304  } else {
305  if (out.Geometry() != QUDA_VECTOR_GEOMETRY) errorQuda("Unsupported geometry %d", out.Geometry());
306 
307  checkMomOrder(in);
308  checkMomOrder(out);
309 
310  // this overrides the check that the texture maps to the gauge
311  // pointer - this is safe here since it only occurs when running
312  // the copier on the host when we will not be using texture
313  // reads
314  const bool override = true;
315 
316  // momentum only currently supported on MILC (10), TIFR (18) and Float2 (10) fields currently
317  if (out.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
318  if (in.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
319  typedef FloatNOrder<FloatOut,10,2,10> momOut;
320  typedef FloatNOrder<FloatIn,10,2,10> momIn;
321  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out, 0, override), momIn(in, In, 0, override), in);
322  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
323  } else if (in.Order() == QUDA_QDP_GAUGE_ORDER) {
324 #ifdef BUILD_QDP_INTERFACE
325  typedef FloatNOrder<FloatOut,10,2,10> momOut;
326  typedef QDPOrder<FloatIn,10> momIn;
327  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out, 0, override), momIn(in, In), in);
328  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
329 #else
330  errorQuda("QDP interface has not been built\n");
331 #endif
332  } else if (in.Order() == QUDA_MILC_GAUGE_ORDER) {
333 #ifdef BUILD_MILC_INTERFACE
334  typedef FloatNOrder<FloatOut,10,2,10> momOut;
335  typedef MILCOrder<FloatIn,10> momIn;
336  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out, 0, override), momIn(in, In), in);
337  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
338 #else
339  errorQuda("MILC interface has not been built\n");
340 #endif
341  } else if (in.Order() == QUDA_MILC_SITE_GAUGE_ORDER) {
342 #ifdef BUILD_MILC_INTERFACE
343  typedef FloatNOrder<FloatOut,10,2,10> momOut;
344  typedef MILCSiteOrder<FloatIn,10> momIn;
345  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out, 0, override), momIn(in, In), in);
346  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
347 #else
348  errorQuda("MILC interface has not been built\n");
349 #endif
350  } else if (in.Order() == QUDA_TIFR_GAUGE_ORDER) {
351 #ifdef BUILD_TIFR_INTERFACE
352  typedef FloatNOrder<FloatOut,18,2,11> momOut;
353  typedef TIFROrder<FloatIn,18> momIn;
354  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out, 0, override), momIn(in, In), in);
355  copyMom<FloatOut,FloatIn,18,momOut,momIn>(arg,out,in,location);
356 #else
357  errorQuda("TIFR interface has not been built\n");
358 #endif
359  } else if (in.Order() == QUDA_TIFR_PADDED_GAUGE_ORDER) {
360 #ifdef BUILD_TIFR_INTERFACE
361  typedef FloatNOrder<FloatOut,18,2,11> momOut;
362  typedef TIFRPaddedOrder<FloatIn,18> momIn;
363  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out, 0, override), momIn(in, In), in);
364  copyMom<FloatOut,FloatIn,18,momOut,momIn>(arg,out,in,location);
365 #else
366  errorQuda("TIFR interface has not been built\n");
367 #endif
368  } else {
369  errorQuda("Gauge field orders %d not supported", in.Order());
370  }
371  } else if (out.Order() == QUDA_QDP_GAUGE_ORDER) {
372  typedef QDPOrder<FloatOut,10> momOut;
373 #ifdef BUILD_QDP_INTERFACE
374  if (in.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
375  typedef FloatNOrder<FloatIn,10,2,10> momIn;
376  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In, 0, override), in);
377  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
378  } else if (in.Order() == QUDA_MILC_GAUGE_ORDER) {
379  typedef MILCOrder<FloatIn,10> momIn;
380  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In), in);
381  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
382  } else {
383  errorQuda("Gauge field orders %d not supported", in.Order());
384  }
385 #else
386  errorQuda("QDP interface has not been built\n");
387 #endif
388  } else if (out.Order() == QUDA_MILC_GAUGE_ORDER) {
389  typedef MILCOrder<FloatOut,10> momOut;
390 #ifdef BUILD_MILC_INTERFACE
391  if (in.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
392  typedef FloatNOrder<FloatIn,10,2,10> momIn;
393  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In, 0, override), in);
394  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
395  } else if (in.Order() == QUDA_MILC_GAUGE_ORDER) {
396  typedef MILCOrder<FloatIn,10> momIn;
397  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In), in);
398  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
399  } else if (in.Order() == QUDA_QDP_GAUGE_ORDER) {
400  typedef QDPOrder<FloatIn,10> momIn;
401  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In), in);
402  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
403  } else {
404  errorQuda("Gauge field orders %d not supported", in.Order());
405  }
406 #else
407  errorQuda("MILC interface has not been built\n");
408 #endif
409  } else if (out.Order() == QUDA_MILC_SITE_GAUGE_ORDER) {
410  typedef MILCSiteOrder<FloatOut,10> momOut;
411 #ifdef BUILD_MILC_INTERFACE
412  if (in.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
413  typedef FloatNOrder<FloatIn,10,2,10> momIn;
414  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In, 0, override), in);
415  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
416  } else if (in.Order() == QUDA_MILC_SITE_GAUGE_ORDER) {
417  typedef MILCSiteOrder<FloatIn,10> momIn;
418  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In), in);
419  copyMom<FloatOut,FloatIn,10,momOut,momIn>(arg,out,in,location);
420  } else {
421  errorQuda("Gauge field orders %d not supported", in.Order());
422  }
423 #else
424  errorQuda("MILC interface has not been built\n");
425 #endif
426  } else if (out.Order() == QUDA_TIFR_GAUGE_ORDER) {
427  typedef TIFROrder<FloatOut,18> momOut;
428 #ifdef BUILD_TIFR_INTERFACE
429  if (in.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
430  // FIX ME - 11 is a misnomer to avoid confusion in template instantiation
431  typedef FloatNOrder<FloatIn,18,2,11> momIn;
432  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In, 0, override), in);
433  copyMom<FloatOut,FloatIn,18,momOut,momIn>(arg,out,in,location);
434  } else if (in.Order() == QUDA_TIFR_GAUGE_ORDER) {
435  typedef TIFROrder<FloatIn,18> momIn;
436  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In), in);
437  copyMom<FloatOut,FloatIn,18,momOut,momIn>(arg,out,in,location);
438  } else {
439  errorQuda("Gauge field orders %d not supported", in.Order());
440  }
441 #else
442  errorQuda("TIFR interface has not been built\n");
443 #endif
444  } else if (out.Order() == QUDA_TIFR_PADDED_GAUGE_ORDER) {
445  typedef TIFRPaddedOrder<FloatOut,18> momOut;
446 #ifdef BUILD_TIFR_INTERFACE
447  if (in.Order() == QUDA_FLOAT2_GAUGE_ORDER) {
448  // FIX ME - 11 is a misnomer to avoid confusion in template instantiation
449  typedef FloatNOrder<FloatIn,18,2,11> momIn;
450  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In, 0, override), in);
451  copyMom<FloatOut,FloatIn,18,momOut,momIn>(arg,out,in,location);
452  } else if (in.Order() == QUDA_TIFR_PADDED_GAUGE_ORDER) {
453  typedef TIFRPaddedOrder<FloatIn,18> momIn;
454  CopyGaugeArg<momOut,momIn> arg(momOut(out, Out), momIn(in, In), in);
455  copyMom<FloatOut,FloatIn,18,momOut,momIn>(arg,out,in,location);
456  } else {
457  errorQuda("Gauge field orders %d not supported", in.Order());
458  }
459 #else
460  errorQuda("TIFR interface has not been built\n");
461 #endif
462  } else {
463  errorQuda("Gauge field orders %d not supported", out.Order());
464  }
465  }
466  }
467  };
468 
469  template <typename FloatIn>
470  void copyGenericGauge(GaugeField &out, const GaugeField &in, QudaFieldLocation location, void *Out, void *In,
471  void **ghostOut, void **ghostIn, int type)
472  {
473  instantiatePrecision2<GaugeCopy, FloatIn>(out, in, location, Out, In, ghostOut, ghostIn, type);
474  }
475 
476 } // namespace quda