Extended gauge field routines do not handle boundary conditions correctly

[maddyscientist]

The extended gauge copy and ghost extraction routines do not presently handle the time boundary conditions correctly. When extended, the boundary condition lies on the the R[3]-1 on the first time slice of nodes and on X[3]-R[3]-1 on the last time slice of nodes (where X[] is the local extended dimensions and R[] is the extended radii dimensions).

I think to fix this, we need a modified variant of the timeBoundary functor used by the gauge field order structs. This requires a couple of additions to the GaugeField class such that any gauge field that is extended knows it is extended and it knows what the depth of the extended region is (the R[] array).

  /**
     timeBoundary variant for extended gauge field
     @param idx extended field linear index
     @param X the extended gauge field dimensions
     @param R the radii dimenions of the extended region
     @param tBoundary the boundary condition
     @param isFirstTimeSlice if we're on the first time slice of nodes
     @param isLastTimeSlide if we're on the last time slice of nodes
  */
  template <typename Float>
    __device__ __host__ inline Float timeBoundary(int idx, const int X[QUDA_MAX_DIM], const int R[QUDA_MAX_DIM], 
                          QudaTboundary tBoundary, bool isFirstTimeSlice, bool isLastTimeSlice) {
    if ( idx >= (R[3]-1)*X[0]*X[1]*X[2]/2 && idx < R[3]*X[0]*X[1]*X[2]/2 ) {
      // the boundary condition is on the R[3]-1 time slice
      return isFirstTimeSlice ? tBoundary : 1.0;
    } else if ( idx >= (X[3]-R[3]-1)*X[0]*X[1]*X[2]/2 && idx < (X[3]-R[3])*X[0]*X[1]*X[2]/2 ) { 
      // the boundary condition lies on the X[3]-R[3]-1 time slice
      return isLastTimeSlice ? tBoundary : 1.0;
    } else {
      return 1.0;
    }
  }

[AlexVaq]

Actually adding the R[] array to the

Alex

El 5/11/2014, a las 9:17, mikeaclark notifications@github.com escribió:

The extended gauge copy and ghost extraction routines do not presently
handle the time boundary conditions correctly. When extended, the boundary
condition lies on the the R[3]-1 on the first time slice of nodes and on
X[3]-R[3]-1 on the last time slice of nodes (where X[] is the local
extended dimensions and R[] is the extended radii dimensions).

I think to fix this, we need a modified variant of the timeBoundary functor
used by the gauge field order structs. This requires a couple of additions
to the GaugeField class such that any gauge field that is extended knows it
is extended and it knows what the depth of the extended region is (the R[]
array).

/**
timeBoundary variant for extended gauge field
@param idx extended field linear index
@param X the extended gauge field dimensions
@param R the radii dimenions of the extended region
@param tBoundary the boundary condition
@param isFirstTimeSlice if we're on the first time slice of nodes
@param isLastTimeSlide if we're on the last time slice of nodes
_/
template
device host inline Float timeBoundary(int idx, const int
X[QUDA_MAX_DIM], const int R[QUDA_MAX_DIM],
QudaTboundary tBoundary, bool
isFirstTimeSlice, bool isLastTimeSlice) {
if ( idx >= (R[3]-1)_X[0]_X[1]_X[2]/2 && idx < R[3]_X[0]_X[1]_X[2]/2 ) {
// the boundary condition is on the R[3]-1 time slice
return isFirstTimeSlice ? tBoundary : 1.0;
} else if ( idx >= (X[3]-R[3]-1)_X[0]_X[1]_X[2]/2 && idx <
(X[3]-R[3])_X[0]_X[1]*X[2]/2 ) {
// the boundary condition lies on the X[3]-R[3]-1 time slice
return isLastTimeSlice ? tBoundary : 1.0;
} else {
return 1.0;
}
}

—
Reply to this email directly or view it on GitHub
#173.

[AlexVaq]

Actually adding the R[] array to the GaugeField class might be enough. One
can just set it to zero when the fields are not extended, and use the
extended formulae for computing the boundaries, the sizes, etc. everywhere.
The final result should remain correct for extended and non-extended gauge
fields.

Alex

El 5/11/2014, a las 9:17, mikeaclark notifications@github.com escribió:

The extended gauge copy and ghost extraction routines do not presently
handle the time boundary conditions correctly. When extended, the boundary
condition lies on the the R[3]-1 on the first time slice of nodes and on
X[3]-R[3]-1 on the last time slice of nodes (where X[] is the local
extended dimensions and R[] is the extended radii dimensions).

I think to fix this, we need a modified variant of the timeBoundary functor
used by the gauge field order structs. This requires a couple of additions
to the GaugeField class such that any gauge field that is extended knows it
is extended and it knows what the depth of the extended region is (the R[]
array).

/**
timeBoundary variant for extended gauge field
@param idx extended field linear index
@param X the extended gauge field dimensions
@param R the radii dimenions of the extended region
@param tBoundary the boundary condition
@param isFirstTimeSlice if we're on the first time slice of nodes
@param isLastTimeSlide if we're on the last time slice of nodes
_/
template
device host inline Float timeBoundary(int idx, const int
X[QUDA_MAX_DIM], const int R[QUDA_MAX_DIM],
QudaTboundary tBoundary, bool
isFirstTimeSlice, bool isLastTimeSlice) {
if ( idx >= (R[3]-1)_X[0]_X[1]_X[2]/2 && idx < R[3]_X[0]_X[1]_X[2]/2 ) {
// the boundary condition is on the R[3]-1 time slice
return isFirstTimeSlice ? tBoundary : 1.0;
} else if ( idx >= (X[3]-R[3]-1)_X[0]_X[1]_X[2]/2 && idx <
(X[3]-R[3])_X[0]_X[1]*X[2]/2 ) {
// the boundary condition lies on the X[3]-R[3]-1 time slice
return isLastTimeSlice ? tBoundary : 1.0;
} else {
return 1.0;
}
}

—
Reply to this email directly or view it on GitHub
#173.

[maddyscientist]

I believe this issue is now fixed (as of commit a26d948).