Generic Programming – Opening New Frontiers
To solve a limitation of Quasar, in which __kernel__
functions in some circumstances needed to be duplicated for different container types (e.g. vec[int8]
, vec[scalar]
, vec[cscalar]
), there is now finally support for generic programming.
Consider the following program that extracts the diagonal elements of a matrix and that is supposed to deal with arguments of either type mat
or type cmat
:
function y : vec = diag(x : mat)
assert(size(x,0)==size(x,1))
N = size(x,0)
y = zeros(N)
parallel_do(size(y), __kernel__
(x:mat, y:vec, pos:int) -> y[pos] = x[pos,pos])
end
function y : cvec = diag(x : cmat)
assert(size(x,0)==size(x,1))
N = size(x,0)
y = czeros(N)
parallel_do(size(y), __kernel__
(x:cmat, y:cvec, pos : int) -> y[pos] = x[pos,pos])
end
Although function overloading here greatly solves part of the problem (at least from the user’s perspective), there is still duplication of the function diag
. In general, we would like to specify functions that can “work” irrespective of their underlying type.
The solution is to use Generic Programming. In Quasar, this is fairly easy to do:
function y = diag[T](x : mat[T])
assert(size(x,0)==size(x,1))
N = size(x,0)
y = vec[T](N)
parallel_do(size(y), __kernel__
(pos) -> y[pos] = x[pos,pos])
end
As you can see, the types of the function signature have simply be omitted. The same holds for the __kernel__
function.
In this example, the type parameter T
is required because it is needed for the construction of vector y
(through the vec[T]
constructor). If T==scalar
, vec[T]
reduces to zeros
, while ifT==cscalar
, vec[T]
reduces to czeros
(complex-valued zero matrix). In case the type parameter is not required, it can be dropped, as in the following example:
function [] = copy_mat(x, y)
assert(size(x)==size(y))
parallel_do(size(y), __kernel__
(pos) -> y[pos] = x[pos])
end
Remarkably, this is still a generic function in Quasar; no special syntax is needed here.
Note that in previous versions of Quasar, all kernel function parameters needed to be explicitly typed. This is now no longer the case: the compiler will deduce the parameter types by calls to diag
and by applying the internal type inference mechanism. The same holds for the __device__
functions.
When calling diag
with two different types of parameters (for example once with x:mat
and a second time with x:cmat
), the compiler will make two generic instantiations of diag
. Internally, the compiler may either:
- Keep the generic definition (type erasion)
function y = diag(x)
- Make two instances of
diag
(reification):function y : vec = diag(x : mat) function y : cvec = diag(x : cmat)
The compiler will combine these two techniques in a transparent way, such that:
- For kernel-functions explicit code is generated for the specific data types,
- For less performance-critical host code type erasion is used (to avoid code duplication).
The selection of the code to run is made at compile-time, so correspondingly the Spectroscope Debugger needs special support for this.
Of course, when calling the diag
function with a variable of type that cannot be determined at compile-time, a compiler error is generated:
The type of the arguments ('op') needs to be fully defined for this function call!
This is then similar to the original handling of kernel functions.
Extensions
There are several extensions possible to fine-tune the behavior of the generic code.
Type classes
Type classes allow the type range of the input parameters to be narrowed. For example:
function y = diag(x : [mat|cmat])
This construction only allows variables of the type mat
and cmat
to be passed to the function. This is useful when it is already known in advance which types are relevant (in this case a real-valued or complex-valued matrix).
Equivalently, type class aliases can be defined. The type:
type AllInt : [int|int8|int16|int32|uint8|uint32|uint64]
groups all integer types that exist in Quasar. Type classes are also useful for defining reductions:
type RealNumber: [scalar|cube|AllInt|cube[AllInt]]
type ComplexNumber: [cscalar|ccube]
reduction (x : RealNumber) -> real(x) = x
Without type classes, the reduction would need to be written 4 times, one for each element.
Type parameters
Levels of genericity
There are three levels of genericity (for which generic instances can be constructed):
- Type constraints: a type constraint binds the type of an input argument of the function.
- Value constraints: gives an explicit value to the value of an input argument
- Logic predicates that are not type or value constraints
As an example, consider the following generic function:
function y = __device__ soft_thresholding(x, T)
if abs(x)>=T
y = (abs(x) - T) * (x / abs(x))
else
y = 0
endif
end
reduction x : scalar -> abs(x) = x where x >= 0
Now, we can make a specialization of this function to a specific type:
soft_thresholding_real = $specialize(soft_thresholding,
type(x,"scalar") && type(T, "scalar"))
But also for a fixed threshold:
soft_thresholding_T = $specialize(soft_thresholding,T==10)
We can even go one step further and specify that x>0
:
soft_thresholding_P = $specialize(soft_thresholding,x>0)
Everything combined, we get:
soft_thresholding_E = $specialize(soft_thresholding,
type(x,"scalar") && type(T,"scalar") && T==10 && x>0)
Based on this knowledge (and the above reduction), the compiler will then generate the following function:
function y = __device__ soft_thresholding_E(x : scalar, T : scalar)
if x >= 10
y = x - 10
else
y = 0
endif
end
There are two ways of performing this type of specialization:
- Using the
$specialize
function. Note that this approach is only recommended for testing. - Alternatively, the specializations can be performed automatically, using the
assert
function from the calling function:function [] = __kernel__ denoising(x : mat, y : mat) assert(x[pos]>0) y[pos] = soft_thresholding(x[pos], 10) end