Improve robustness of simultaneous phase and boundary condition handling

[weinbe2]

Currently, the t_boundary member of the GaugeField class is serving double-duty.

• In the case of applying/removing staggered phases, it keeps track of what boundary conditions the field had when it was created. This manifests itself as not updating GaugeField's t_boundary member when the apply/remove staggered phases routine is called, even if the routine is being used to add/remove anti-periodic boundary conditions.
• In all other cases (namely accessors), t_boundary is used to maintain the current state of the gauge field.

In many cases, this is harmless, because the minus signs cancel in local closed loops---so it's safe for most gauge observables, forces, etc. This doesn't apply in the dslash routines either, because phases are kept on the fields appropriately.

This does manifest itself for gauge observables that wrap around the lattice, though, such as the temporal Polyakov loop. Ultimately, this is the motivation for improving the robustness of phase/boundary condition tracking.

[kostrzewa]

We have a different but related situation when using QUDA from within tmLQCD. Because tmLQCD always uses twisted boundary conditions (always in time, optionally also in the spatial dimensions if one wishes to use momenta shifted using twisted boundary conditions), we bake in phases before sending the gauge field off to QUDA. This means that we can't use gauge compression, which is not a big deal, but does affect performance a little.

When dealing with the gauge force we do not remove the phases.

In some sense, I believe that these phase factors should be properties of the Dirac operators rather than of the gauge field, such that compression and twisted boundary conditions could be used independently (at the cost of an extra complex multiply per direction).

[weinbe2]

@kostrzewa it would take a good amount of hacking in the interface and deeper in the code, but we could proof-of-concept compiling the twisted mass/clover operators with recon-13 and recon-9 (normally used for the long links in the improved staggered operator), which ultimately just enables gauge compression on U(3) fields because the thirteenth component is the argument to the U(1) phase. If you pass in the fields with phases already applied, you can "lie" and say the fields have no phases/boundary conditions and everything should still work. (With a billion asterisks due to the necessary hacking.)

The headaches documented here are very specific to how applying/removing staggered signs and boundary conditions are handled internally. If you can pass in the fields with signs already applied and don't need to change it (or only need to change it very infrequently -> you can manage the data motion hit), a lot of this won't apply.

We can sync up offline if you want to dig into this deeper.