Nothing
###########################################################################
# Test Script
#
# Author: nmorris
###############################################################################
#library(gmatrix,lib.loc="/home/nmorris/.RLIBS")
gtest =function() {
tn=100
tsn=5
scaler=10.13
testscalergpu=list(
as.gvector(scaler,type="d"),
as.gvector(scaler,type="s"),
as.gvector(scaler,type="i"),
as.gvector(scaler,type="l")
)
testscalercpu=list(
as.numeric(scaler),
as.numeric(scaler),
as.integer(scaler),
as.logical(scaler))
vec=((1:tn)-1)*(5/tn)
testvecgpu=list(
as.gvector(vec,type="d"),
as.gvector(vec,type="s"),
as.gvector(vec,type="i"),
as.gvector(vec,type="l")
)
testveccpu=list(
as.numeric(vec),
as.numeric(vec),
as.integer(vec),
as.logical(vec))
mat=matrix(0:(tn*tn-1),tn,tn)*(5/tn)
testmatgpu=list(
as.gmatrix(mat,type="d"),
as.gmatrix(mat,type="s"),
as.gmatrix(mat,type="i"),
as.gmatrix(mat,type="l"))
testmatcpu=list(
matrix(as.numeric(mat),tn,tn),
matrix(as.numeric(mat),tn,tn),
matrix(as.integer(mat),tn,tn),
matrix(as.logical(mat),tn,tn)
)
vec=5*((1:tsn)-1)/tsn
testsvecgpu=list(
as.gvector(vec,type="d"),
as.gvector(vec,type="s"),
as.gvector(vec,type="i"),
as.gvector(vec,type="l")
)
testsveccpu=list(
as.numeric(vec),
as.numeric(vec),
as.integer(vec),
as.logical(vec))
mat=matrix(0:(tsn*tsn-1),tsn,tsn)*(5/tn)
testsmatgpu=list(
as.gmatrix(mat,type="d"),
as.gmatrix(mat,type="s"),
as.gmatrix(mat,type="i"),
as.gmatrix(mat,type="l"))
testsmatcpu=list(
matrix(as.numeric(mat),tn,tn),
matrix(as.numeric(mat),tn,tn),
matrix(as.integer(mat),tn,tn),
matrix(as.logical(mat),tn,tn)
)
.type_name= function(type) {
if(type==0L)
return("double")
else if(type==1L)
return("single")
else if(type==2L)
return("integer")
else if(type==3L)
return("logical")
else
stop("Invalid Type")
}
.type_name2= function(type) {
if(type==0L)
return("double")
else if(type==1L)
return("single")
else if(type==2L)
return("integer")
else if(type==3L)
return("logical")
else
stop("Invalid Type")
}
warningslist=list()
###################################
# Test Matrix Multiplication
###################################
cputype=function(gpu,type) {
ret=type
if(!gpu)
if(type==1)
ret=0
return(ret)
}
checkmm=function(x1,y1,mat1,mat2,gpu1,gpu2, type1=cputype(gpu1,i-1),type2=cputype(gpu2,j-1)) {
yexpr=substitute(y1, list(i=i,j=j))
xexpr=substitute(x1, list(i=i,j=j))
tmpnm=deparse(expr=xexpr)
mywarnings=character(0)
x=tryCatch(eval(xexpr), error=function(e) {
mywarnings<<-conditionMessage(e)
return(-99999L)})
if(as.logical(x[1]!=-99999L)) {
y=eval(yexpr)
if(class(x) %in% c("gmatrix","gvector")) {
maxtype=max(type1,type2)
mintype=min(type1,type2)
if(maxtype>=1L & mintype==1L)
out_type=1L
else
out_type=0L
if(!gpu1 & type1==0L)
if(type2==1L)
out_type=1L
if(!gpu2 & type2==0L)
if(type1==1L)
out_type=1L
if(x@type!=out_type) {
mywarnings=c(mywarnings,paste("Incorrect type returned:",.type_name(x@type)))
# browser()
}
x=h(x)
}
if(class(y) %in% c("gmatrix","gvector")) {
y=h(y)
}
tmp=sum((y-x)^2)/sum(x^2)
if(!is.numeric(tmp)||any(is.na(tmp)))
mywarnings=c(mywarnings,paste(mywarnings, "result missing or not numeric."))
else if(tmp>10^(-6))
mywarnings=c(mywarnings,"not equal")
}
if(length(mywarnings)>0) {
c1=ifelse(gpu1,paste(ifelse(mat1,"gmatrix","gvector"),"of type",.type_name(type1)),
paste(ifelse(mat1,"matrix","vector"),"of type",.type_name2(type1)))
c2=ifelse(gpu2,paste(ifelse(mat2,"gmatrix","gvector"),"of type",.type_name(type2)),
paste(ifelse(mat2,"matrix","vector"),"of type",.type_name2(type2)))
mywarnings=list(mywarnings)
names(mywarnings)=paste("applying",opname," a",c1, "and a",c2,".")
warningslist<<-c(warningslist,mywarnings)
}
return(mywarnings)
}
#warningslist=list()
cat("Checking matrix multiplication, crossprod and tcrossprod... \n")
multOps=list(
list("matrix multiplication", function(a,b) a %*% b),
list("crossprod", function(a,b) crossprod(a, b)),
list("tcrossprod", function(a,b) tcrossprod(a , b))
)
for(op in multOps) {
# browser()
#suppressWarnings(setGeneric("%op%",op[[2]],where=globalenv()))
#print(op)
opname=op[[1]]
opf=op[[2]]
for(i in 1:4){
for(j in 1:4){
checkmm(opf(testmatgpu[[i]], testmatgpu[[j]]), opf(testmatcpu[[i]], testmatcpu[[j]]),mat1=TRUE,mat2=TRUE,gpu1=TRUE,gpu2=TRUE)
checkmm(opf(testmatcpu[[i]], testmatgpu[[j]]), opf(testmatcpu[[i]], testmatcpu[[j]]),mat1=TRUE,mat2=TRUE,gpu1=FALSE,gpu2=TRUE)
checkmm(opf(testmatgpu[[i]], testmatcpu[[j]]), opf(testmatcpu[[i]], testmatcpu[[j]]),mat1=TRUE,mat2=TRUE,gpu1=TRUE,gpu2=FALSE)
if(opname!="tcrossprod"){
checkmm(opf(testmatgpu[[i]], testvecgpu[[j]]), opf(testmatcpu[[i]], testveccpu[[j]]),mat1=TRUE,mat2=FALSE,gpu1=TRUE,gpu2=TRUE)
checkmm(opf(testmatgpu[[i]], testveccpu[[j]]), opf(testmatcpu[[i]], testveccpu[[j]]),mat1=TRUE,mat2=FALSE,gpu1=TRUE,gpu2=FALSE)
checkmm(opf(testmatcpu[[i]], testvecgpu[[j]]), opf(testmatcpu[[i]], testveccpu[[j]]),mat1=TRUE,mat2=FALSE,gpu1=FALSE,gpu2=TRUE)
}
checkmm(opf(testvecgpu[[i]], testmatgpu[[j]]), opf(testveccpu[[i]], testmatcpu[[j]]),mat1=FALSE,mat2=TRUE,gpu1=TRUE,gpu2=TRUE)
checkmm(opf(testvecgpu[[i]], testmatcpu[[j]]), opf(testveccpu[[i]], testmatcpu[[j]]),mat1=FALSE,mat2=TRUE,gpu1=TRUE,gpu2=FALSE)
checkmm(opf(testveccpu[[i]], testmatgpu[[j]]), opf(testveccpu[[i]], testmatcpu[[j]]),mat1=FALSE,mat2=TRUE,gpu1=FALSE,gpu2=TRUE)
checkmm(opf(testvecgpu[[i]], testvecgpu[[j]]), opf(testveccpu[[i]], testveccpu[[j]]),mat1=FALSE,mat2=FALSE,gpu1=TRUE,gpu2=TRUE)
checkmm(opf(testveccpu[[i]], testvecgpu[[j]]), opf(testveccpu[[i]], testveccpu[[j]]),mat1=FALSE,mat2=FALSE,gpu1=FALSE,gpu2=TRUE)
checkmm(opf(testvecgpu[[i]], testveccpu[[j]]), opf(testveccpu[[i]], testveccpu[[j]]),mat1=FALSE,mat2=FALSE,gpu1=TRUE,gpu2=FALSE)
}
}
}
#rm("%op%", envir=globalenv() )
#warningslist
cat("Checking outer product and kronecker product... \n")
if(any(h(testvecgpu[[1]] %o% testvecgpu[[1]])!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(any(h(testveccpu[[1]] %o% testvecgpu[[1]])!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(any(h(testvecgpu[[1]] %o% testveccpu[[1]])!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(any(h(testvecgpu[[1]] %o% testveccpu[[1]])!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(type(testvecgpu[[2]] %o% testveccpu[[1]])!= "single")
warningslist<-c(warningslist,paste("%o% not working correctly (should return single)"))
if(type(testveccpu[[1]] %o% testvecgpu[[2]])!= "single")
warningslist<-c(warningslist,paste("%o% not working correctly (should return single)"))
if(any(h(gouter(testvecgpu[[1]] , testvecgpu[[1]]))!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(any(h(gouter(testveccpu[[1]] , testvecgpu[[1]]))!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(any(h(gouter(testvecgpu[[1]] , testveccpu[[1]]))!=testveccpu[[1]] %o% testveccpu[[1]] ))
warningslist<-c(warningslist,paste("%o% not working correctly"))
if(any(h(gkroneckerProd(testmatgpu[[1]][1:10,1:10] , testmatgpu[[1]][1:10,1:10]))!=testmatcpu[[1]][1:10,1:10] %x% testmatcpu[[1]][1:10,1:10] ))
warningslist<-c(warningslist,paste("gkroneckerProd not working correctly"))
if(any(h(gkroneckerProd(testmatcpu[[1]][1:10,1:10] , testmatgpu[[1]][1:10,1:10]))!=testmatcpu[[1]][1:10,1:10] %x% testmatcpu[[1]][1:10,1:10] ))
warningslist<-c(warningslist,paste("gkroneckerProd not working correctly"))
if(any(h(testmatgpu[[1]][1:10,1:10] %x% testmatcpu[[1]][1:10,1:10])!=testmatcpu[[1]][1:10,1:10] %x% testmatcpu[[1]][1:10,1:10] ))
warningslist<-c(warningslist,paste("%x% not working correctly"))
if(any(h(testvecgpu[[1]] %x% testveccpu[[1]])!=testveccpu[[1]] %x% testveccpu[[1]]))
warningslist<-c(warningslist,paste("%x% not working correctly"))
######################################################
# Test Exchangeable + Nonexchangeable Binary elementwise operations #
######################################################
exOps=list(
c("mult","*", "min", "2"), #name, operator, output type, maximum output type
c("add" ,"+", "min", "2"),
c("eq" ,"==", "3", "2"),
c("ne" ,"!=", "3", "2" ),
c("and" ,"&", "3", "2"),
c("or" ,"|", "3", "2" ),
c("sub","-", "min", "2" ),
c("div","/", "min", "1"), #always returns at least a single
c("pow","^", "min", "1" ),#always returns at least a single
# c("mod","%%", "min", "2"), #cuda is not numericaly accurate on this...Todo: add note in
c("gt" ,">", "3", "2"),
c("lt" ,"<", "3", "2"),
c("gte" ,">=", "3", "2" ),
c("lte" ,"<=", "3", "2" )
)
checkbinobs=function(x1,y1,mat1,mat2,gpu1,gpu2, type1=i-1,type2=j-1, out_type=outt, op=binaryop, maxtype=as.integer(maxouttype)) {
if(out_type=="min") {
out_type=min(type1,type2)
if(out_type>maxtype) {
if(maxtype==2)
out_type=maxtype
else
out_type=0L
}
if(gpu1 && !gpu2) {
if(type1==1)
out_type=1
else if(type2==1)
out_type=0
}
if(!gpu1 && gpu2) {
if(type2==1)
out_type=1
else if(type1==1)
out_type=0
}
# if(gpu1 && gpu2)
# if(type1==1 || type2==1)
# out_type=0
} else
out_type=as.numeric(out_type)
yexpr=substitute(y1, list(i=i,j=j))
xexpr=substitute(x1, list(i=i,j=j))
mywarnings=character(0)
x=tryCatch(eval(xexpr), error=function(e) {
mywarnings<<-conditionMessage(e)
return(-99999L)
}
)
#browser()
if(as.logical(x[1]!=-99999L)) {
y=eval(yexpr)
# if(x@type<2)
if(x@type!=out_type) {
# browser()
mywarnings=c(mywarnings,paste("Incorrect type returned:",.type_name(x@type)))
}
x=h(x)
#browser()
delta=max(abs(x-y)/abs(x))
#tmp=sum((y-x)^2)/sum(x^2)
if(!is.finite(delta) || max(abs(x))<10^-6) {
if(any(is.finite(y)!=is.finite(x)) & exOp[2]!="^") { #Todo: add note about "^"
#browser()
c(mywarnings,paste("not equal: rel error not defined"))
}else {
# browser()
good=is.finite(x) & is.finite(y)
tmpx=x[good];tmpy=y[good]
delta=(tmpx-tmpy)/tmpx
delta=delta[is.finite(delta) & abs(tmpx)>10^-4]
if(length(delta)>0)
if(max(delta)>ifelse(max(type1,type2)==1| out_type==1,10^-4,10^-5)) {
#browser()
mywarnings=c(mywarnings,paste("not equal: rel error not defined"))
}
}
} else if(delta>ifelse(max(type1,type2)==1 | out_type==1,10^-4,10^-9)) {
mywarnings=c(mywarnings,paste("not equal: rel error = ", sum((y-x)^2)/sum(x^2)))
}
}
tmpf=function(mat,gpu,type) {
ap=ifelse(gpu,"g","")
if(mat==0)
ret= paste("scaler " , ap, "vector",sep="")
else if(mat==1)
ret= paste(ap, "vector",sep="")
else
ret= paste(ap, "matrix",sep="")
return(paste(ret,"of type",.type_name2(type)))
}
if(length(mywarnings)>0) {
c1=tmpf(mat1,gpu1,type1)
c2=tmpf(mat2,gpu2,type2)
mywarnings=list(mywarnings)
#browser()
names(mywarnings)=paste("Operation " , op, " failed for: (" ,c1,") ", op," (" ,c2,").",sep="")
warningslist<<-c(warningslist,mywarnings)
}
return(mywarnings)
}
cat("Checking Binary Operations... ")
for(exOp in exOps) {
#exOp=exOps[[8]]
binaryop=exOp[2]
maxouttype=exOp[4]
outt=exOp[3]
# suppressWarnings(
# setGeneric("%op%",
# function(a,b)
# return(eval(parse(text=paste("a",binaryop,"b")))
# ),
# where=globalenv()
# )
# )
opf=function(a,b)
return(eval(parse(text=paste("a",binaryop,"b"))))
cat(binaryop," ")
for(i in 1:4){
for(j in 1:4){
# cat(i,j,"\n")
checkbinobs(opf(testmatgpu[[i]], testmatgpu[[j]]), opf(testmatcpu[[i]], testmatcpu[[j]]),mat1=2,mat2=2,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testmatcpu[[i]], testmatgpu[[j]]), opf(testmatcpu[[i]], testmatcpu[[j]]),mat1=2,mat2=2,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testmatgpu[[i]], testmatcpu[[j]]), opf(testmatcpu[[i]], testmatcpu[[j]]),mat1=2,mat2=2,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testmatgpu[[i]], testvecgpu[[j]]), opf(testmatcpu[[i]], testveccpu[[j]]),mat1=2,mat2=1,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testmatgpu[[i]], testveccpu[[j]]), opf(testmatcpu[[i]], testveccpu[[j]]),mat1=2,mat2=1,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testmatcpu[[i]], testvecgpu[[j]]), opf(testmatcpu[[i]], testveccpu[[j]]),mat1=2,mat2=1,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testmatgpu[[j]]), opf(testveccpu[[i]], testmatcpu[[j]]),mat1=1,mat2=2,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testmatcpu[[j]]), opf(testveccpu[[i]], testmatcpu[[j]]),mat1=1,mat2=2,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testveccpu[[i]], testmatgpu[[j]]), opf(testveccpu[[i]], testmatcpu[[j]]),mat1=1,mat2=2,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testvecgpu[[j]]), opf(testveccpu[[i]], testveccpu[[j]]),mat1=1,mat2=1,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testveccpu[[i]], testvecgpu[[j]]), opf(testveccpu[[i]], testveccpu[[j]]),mat1=1,mat2=1,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testveccpu[[j]]), opf(testveccpu[[i]], testveccpu[[j]]),mat1=1,mat2=1,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testvecgpu[[i]], testscalergpu[[j]]), opf(testveccpu[[i]], testscalercpu[[j]]),mat1=1,mat2=0,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testveccpu[[i]], testscalergpu[[j]]), opf(testveccpu[[i]], testscalercpu[[j]]),mat1=1,mat2=0,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testscalercpu[[j]]), opf(testveccpu[[i]], testscalercpu[[j]]),mat1=1,mat2=0,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testscalergpu[[i]], testvecgpu[[j]]), opf(testscalercpu[[i]], testveccpu[[j]]),mat1=0,mat2=1,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testscalergpu[[i]], testveccpu[[j]]), opf(testscalercpu[[i]], testveccpu[[j]]) ,mat1=0,mat2=1,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testscalercpu[[i]], testvecgpu[[j]]), opf(testscalercpu[[i]], testveccpu[[j]]) ,mat1=0,mat2=1,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testmatgpu[[i]], testscalergpu[[j]]), opf(testmatcpu[[i]], testscalercpu[[j]]),mat1=2,mat2=0,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testmatcpu[[i]], testscalergpu[[j]]), opf(testmatcpu[[i]], testscalercpu[[j]]),mat1=2,mat2=0,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testmatgpu[[i]], testscalercpu[[j]]), opf(testmatcpu[[i]], testscalercpu[[j]]),mat1=2,mat2=0,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testscalergpu[[i]], testmatgpu[[j]]), opf(testscalercpu[[i]], testmatcpu[[j]]),mat1=0,mat2=2,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testscalergpu[[i]], testmatcpu[[j]]), opf(testscalercpu[[i]], testmatcpu[[j]]) ,mat1=0,mat2=2,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testscalercpu[[i]], testmatgpu[[j]]), opf(testscalercpu[[i]], testmatcpu[[j]]) ,mat1=0,mat2=2,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testsvecgpu[[j]]), opf(testveccpu[[i]], testsveccpu[[j]]),mat1=1,mat2=0,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testveccpu[[i]], testsvecgpu[[j]]), opf(testveccpu[[i]], testsveccpu[[j]]),mat1=1,mat2=0,gpu1=FALSE,gpu2=TRUE)
checkbinobs(opf(testvecgpu[[i]], testsveccpu[[j]]), opf(testveccpu[[i]], testsveccpu[[j]]),mat1=1,mat2=0,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testsvecgpu[[i]], testvecgpu[[j]]), opf(testsveccpu[[i]], testveccpu[[j]]),mat1=0,mat2=1,gpu1=TRUE,gpu2=TRUE)
checkbinobs(opf(testsvecgpu[[i]], testveccpu[[j]]), opf(testsveccpu[[i]], testveccpu[[j]]) ,mat1=0,mat2=1,gpu1=TRUE,gpu2=FALSE)
checkbinobs(opf(testsveccpu[[i]], testvecgpu[[j]]), opf(testsveccpu[[i]], testveccpu[[j]]) ,mat1=0,mat2=1,gpu1=FALSE,gpu2=TRUE)
}
}
}
#########
unOps= list(
c('sqrt','sqrt', "x@type"),
c('exp','exp', "x@type"),
c('expm1','expm1', "x@type"),
c('log','log', "x@type"),
c('log2','log2', 'x@type'),
c('log10','log10', 'x@type'),
c('log1p','log1p', 'x@type'),
c('sin','sin', 'x@type'),
c('cos','cos', 'x@type'),
c('tan','tan', 'x@type'),
c('asin','asin', 'x@type'),
c('acos','acos', 'x@type'),
c('atan','atan', 'x@type'),
c('sinh','sinh', 'x@type'),
c('cosh','cosh', 'x@type'),
c('tanh','tanh', 'x@type'),
c('asinh','asinh',' x@type'),
c('acosh','acosh', 'x@type'),
c('atanh','atanh', 'x@type'),
c('abs','abs', 'x@type'),
c('lgamma','lgamma', 'x@type'),
c('gamma','gamma', 'x@type'),
c('sign','sign', 'x@type'),
c('round','round', '2L'),
c('ceiling','ceil', '2L'),
c('floor','floor', '2L'),
c('is.na','isna', '3L'),
c('is.nan','isnan', '3L'),
c('is.finite','isfinite', '3L'),
c('is.infinite','isinfinite', '3L'),
c("!","!",""),
c("-","-",""),
c("+","+",""))
checkunobs=function(x1,y, op=unaryop) {
mywarnings=character()
x=tryCatch(eval(x1), error=function(e) {
mywarnings<<-conditionMessage(e)
return(-99999L)
}
)
if(as.logical(x[1]!=-99999L)) {
x=h(x)
y=suppressWarnings(y)
tmp=sum((y-x)^2)/sum(y^2)
if(!is.finite(tmp) || sum(x^2)<10^-3) {
if(any(is.nan(y)!=is.nan(x)))
c(mywarnings,paste("not equal: rel error not defined"))
else {
good=is.finite(x) & is.finite(y)
if(sum((y[good]-x[good])^2)>10^-4)
mywarnings=c(mywarnings,paste("not equal: rel error not defined"))
}
} else if(tmp>10^-4)
mywarnings=c(mywarnings,paste("not equal: rel error = ", sum((y-x)^2)/sum(x^2)))
}
if(length(mywarnings)>0) {
mywarnings=list(mywarnings)
names(mywarnings)=paste("Operation " , unaryop, " failed.",sep="")
warningslist<<-c(warningslist,mywarnings)
}
return(mywarnings)
}
cat("\nChecking Unary Operations/special functions... ")
for(exOp in unOps) {
unaryop=exOp[1]
callexpr=parse(text=paste(unaryop,"(x)"))
cat(unaryop," ")
myop = function(x)
return(eval(callexpr))
for(i in 1:4) {
checkunobs(myop(testvecgpu[[i]]), myop(testveccpu[[i]]))
}
}
cat("\nChecking ifelse and which... \n")
if(any(h(ifelse(gmatrix(c(T,F),4,4),c(1,2),c(1,2,3)))!=ifelse(matrix(c(T,F),4,4),c(1,2),c(1,2,3))))
warningslist<-c(warningslist,"ifelse not working correctly")
if(any(h(which(g.rep(c(T,F),100)))!=which(rep(c(T,F),100))))
warningslist<-c(warningslist,"which or g.rep not working correctly")
cat("Checking sort and order... \n")
for(i in 0:2) {
if(any(h(sort(as.gvector(c(1,4,3,2),type=i)))!=c(1,2,3,4)))
warningslist<-c(warningslist,paste("sort not working correctly:",i))
if(any(h(sort(as.gvector(c(1,4,3,2)),decreasing=TRUE))!=c(4,3,2,1)))
warningslist<-c(warningslist,paste("sort decreasing=TRUE not working correctly:",i))
}
cat("Checking max, min, sum and col/row sums/means... \n")
for(i in 0:2) {
if(h(max(as.gvector(c(1,4,3,2),type=i)))!=4)
warningslist<-c(warningslist,paste("max not working correctly for type:",i))
if(h(min(as.gvector(c(1,4,3,2),type=i)))!=1)
warningslist<-c(warningslist,paste("min not working correctly for type:",i))
if(h(sum(as.gvector(c(1,4,3,2),type=i)))!=10)
warningslist<-c(warningslist,paste("sum not working correctly for type:",i))
if(h(mean(as.gvector(c(1,4,3,2),type=i)))!=2.5)
warningslist<-c(warningslist,paste("mean not working correctly for type:",i))
if(h(mean(as.gvector(c(1,4,3,2),type=i)))!=2.5)
warningslist<-c(warningslist,paste("mean not working correctly for type:",i))
checkunobs(rowMeans(testmatgpu[[i+1]]), rowMeans(testmatcpu[[i+1]]), "rowMeans" )
checkunobs(rowSums( testmatgpu[[i+1]]), rowSums(testmatcpu[[i+1]]), "rowSums" )
checkunobs(colMeans(testmatgpu[[i+1]]), colMeans(testmatcpu[[i+1]]), "colMeans" )
checkunobs(colSums( testmatgpu[[i+1]]), colSums(testmatcpu[[i+1]]), "colSums" )
}
cat("Checking gsumby... \n")
if(any(h(gsumby(1:10, c(1,5),c(4,10)))!=c(10,45)))
warningslist<-c(warningslist,paste("gsumby not working correctly"))
cat("Checking transpose... \n")
for(i in 1:4){
if(any(h(t(testmatgpu[[i]])!=t(testmatcpu[[i]]))))
warningslist<-c(warningslist,paste("Error in transpose type:",i-1))
if(any(h(t(testvecgpu[[i]])!=t(testveccpu[[i]]))))
warningslist<-c(warningslist,paste("Error in transpose type:",i-1))
}
cat("Checking logRowSums... \n")
tmp=cbind(rep(7,10000),1,1,1)
gtmp=g(log(tmp))
if(any( abs(h(gRowLogSums(gtmp))-log(10))>10^-9))
warningslist<-c(warningslist, "Error in 'gRowLogSums.'")
cat("Checking indexing and diag functions... \n")
for(i in 1:4) {
for(j in 1:4)
testvecgpu[[i]][2:5]=testvecgpu[[j]][7:10,]
if(any(h(testvecgpu[[i]][2:5]!=testvecgpu[[j]][7:10])))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
testmatgpu[[i]][2:5,]=testmatgpu[[j]][7:10,]
if(any(h(testmatgpu[[i]][2:5,]!=testmatgpu[[j]][7:10,])))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
testmatgpu[[i]][2:5,]=testmatgpu[[j]][7:10,]
if(any(h(testmatgpu[[i]][2:5,]!=testmatgpu[[j]][7:10,])))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
testmatgpu[[i]][,2:5]=testmatgpu[[j]][,7:10]
if(any(h(testmatgpu[[i]][,2:5]!=testmatgpu[[j]][,7:10])))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
testmatgpu[[i]][2:5,2:5]=testmatgpu[[j]][8:11,7:10]
if(any(h(testmatgpu[[i]][2:5,2:5]!=testmatgpu[[j]][8:11,7:10])))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
diag(testmatgpu[[i]])=diag(testmatgpu[[j]])
if(any(h(diag(testmatgpu[[i]])!= convertType(diag(testmatgpu[[j]]) , as.integer(i-1)))))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
diag(testmatgpu[[i]])=diag(testmatcpu[[j]])
if(any(h(diag(testmatgpu[[i]])!=convertType(diag(testmatgpu[[j]]), as.integer(i-1)))))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
diag(testmatgpu[[i]])=(testscalercpu[[j]])
if(any(h(diag(testmatgpu[[i]])!=convertType(testscalergpu[[j]] , as.integer(i-1))) ))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
diag(testmatgpu[[i]])=(testscalergpu[[j]])
if(any(h(diag(testmatgpu[[i]])!=convertType(testscalergpu[[j]] , as.integer(i-1))) ))
warningslist<-c(warningslist,paste("Indexing problem for gvector"))
}
dosimcont =function(qdistfun=qnorm, rdistfun=rnorm, nsims=100000, cell_ct=100) {
boundaries= qdistfun((0:cell_ct)/cell_ct)
y=1:(cell_ct+1)
myfun=approxfun(boundaries, y, method = "constant")
mysims=as.numeric(rdistfun(nsims))
mytbl=as.vector(table(myfun(mysims)))
expected=nsims/cell_ct
chisq = sum( (mytbl - expected)^2/expected)
pval= pchisq(chisq,cell_ct-1, lower.tail=FALSE)
return(pval)
}
dosimint =function(qdistfun, pdistfun, rdistfun, nsims=100000, max1=.99, max2=NULL) { #integer random vars
if(is.null(max2)) {
max2= qdistfun(max1)
}
boundaries=.5:(max2-.5)
b_ct=length(boundaries)
probs=pdistfun(boundaries)
probs=c(probs[1],probs[2:b_ct]-probs[1:(b_ct-1)],1-probs[b_ct])
mysims=as.numeric(rdistfun(nsims))
mysims=ifelse(mysims>=max2,max2,mysims)
#browser()
mytbl=as.vector(table(factor(mysims, levels=0:max2)))
expected=nsims*probs
#deal with small expected cell counts so assymptotics works
bad=expected<5
if(sum(bad)>0) {
prosb=c(probs[!bad], sum(probs[bad]))
expected=c(expected[!bad], sum(expected[bad]))
mytbl=c(mytbl[!bad], sum(mytbl[bad]))
if(expected[length(expected)]<=.05) {
bad=which(min(expected)==expected)[1]
prosb=c(probs[-bad], sum(probs[bad]))
expected=c(expected[!bad], sum(expected[bad]))
mytbl=c(mytbl[-bad], sum(mytbl[bad]))
}}
#browser()
#calculate statistic
chisq = sum( (mytbl - expected)^2/expected)
pval= pchisq(chisq,length(expected)-1, lower.tail=FALSE)
return(pval)
}
mycheck=function(x, distfun) {
mywarnings=character()
ntests=length(x)
x=if(any(x<.005/ntests))
warningslist<<-c(warningslist,paste("Test of the goodness of fit of simulations from function",distfun, "yeilded the following p-values:",paste(as.vector(x), collapse=", ")))
}
cat("Testing Random Generators...\n")
#normal
mycheck(sapply(c(.5,10), function(mean) {
tmp = lapply(c(.01, .5, 1, 5, 10), function(sd) {
#cat(mean,sd,"\n");
dosimcont(
qdistfun=function(q) qnorm(q,mean,sd=sd),
rdistfun=function(n) grnorm(n,mean,sd=sd)
)
})
return(unlist(tmp))
}
), "grnorm")
#gamma
mycheck(sapply(c(.5,10), function(rate) {
tmp = lapply(c(.01, .5, 1, 5, 10), function(shape) {
#cat(rate,shape,"\n");
dosimcont(
qdistfun=function(q) qgamma(q,shape,rate=rate),
rdistfun=function(n) grgamma(n,shape,rate=rate)
)
})
return(unlist(tmp))
}
), "grgamma")
#beta
mycheck(sapply(c(.25, .5, 1, 5, 10), function(shape1) {
tmp = lapply(c(.5, 1, 5, 10), function(shape2) {
# cat(shape1,shape2,"\n");
dosimcont(
qdistfun=function(q) qbeta(q,shape1,shape2),
rdistfun=function(n) grbeta(n,shape1,shape2),cell_ct=25
)
})
return(unlist(tmp))
}
), "grbeta")
#uniform
mycheck(sapply(c(.01, .5, 1, 5, 10), function(min) {
tmp = lapply(c(.01, .5, 1, 5, 10)+10, function(max) {
dosimcont(
qdistfun=function(q) qunif(q,min,max),
rdistfun=function(n) grunif(n,min,max)
)
})
return(unlist(tmp))
}
), "grunif")
#binom
mycheck(sapply(c(10),function(size) {
tmp = lapply(c(.01, .25, .5, .8), function(p) {
dosimint(
qdistfun=function(q) qbinom(q,size,p),
pdistfun=function(q) pbinom(q,size,p),
rdistfun=function(n) grbinom(n,size,p),max2=size
)})
return(unlist(tmp))
}
), "grbinom")
#pois
mycheck(sapply(c(.01, .25, .5, .8), function(lamb) {
dosimint(
qdistfun=function(q) qpois(q,lamb),
pdistfun=function(q) ppois(q,lamb),
rdistfun=function(n) grpois(n,lamb),max1=.999999
)
}), "grpois")
#todo: check rsample
# tmp=cbind(rep(7,10000),1,1,1)
# gtmp=g(log(tmp))
# gRowLogSums(gtmp)
# table(h(rsample(gtmp)))
cat("Checking distribution functions...\n")
#Distribution Checks
checkd=function(a,b,distfun) {
tmp=abs(as.numeric(a)-b)/b
tmp=tmp[abs(b)>10^-20 & is.finite(b)]
if(any(tmp > 10^-8))
warningslist<<-c(warningslist,paste("Problem with function",distfun))
}
checkd( gdnorm ((-10):10, mean=c(1,2), sd=c(1,2,3,4)), dnorm((-10):10, mean=c(1,2), sd=c(1,2,3,4)), "gdnorm")
checkd( gdgamma((1:20)/5, shape=c(1,2), rate=c(1,2,3,4)) , dgamma((1:20)/5, shape=c(1,2), rate=c(1,2,3,4)),"gdgamma")
checkd( gdunif((1:20)/5, min=c(1,2), max=3:6) , dunif((1:20)/5, min=c(1,2), max=3:6), "gdunif")
checkd( gdbeta((1:20)/5, shape1=c(1,2), shape2=c(1,2,3,4)) , dbeta((1:20)/5, shape1=c(1,2), shape2=c(1,2,3,4)), "gdbeta")
checkd( gdbinom(1:20, size=c(20,30), prob=c(1,2,3,4)/5) , dbinom(1:20, size=c(20,30), prob=c(1,2,3,4)/5), "gdbinom")
checkd( gdpois(1:20, c(5,10)) , dpois(1:20, c(5,10)), "gdbinom")
checkd( gdnorm ((-10):10, mean=c(1,2), sd=c(1,2,3,4), log=TRUE), dnorm((-10):10, mean=c(1,2), sd=c(1,2,3,4), log=TRUE), "gdnorm")
checkd( gdgamma((1:20)/5, shape=c(1,2), rate=c(1,2,3,4), log=TRUE) , dgamma((1:20)/5, shape=c(1,2), rate=c(1,2,3,4), log=TRUE),"gdgamma")
checkd( gdunif((1:20)/5, min=c(1,2), max=3:6, log=TRUE) , dunif((1:20)/5, min=c(1,2), max=3:6, log=TRUE), "gdunif")
checkd( gdbeta((1:20)/5, shape1=c(1,2), shape2=c(1,2,3,4), log=TRUE) , dbeta((1:20)/5, shape1=c(1,2), shape2=c(1,2,3,4), log=TRUE), "gdbeta")
checkd( gdbinom(1:20, size=c(20,30), prob=c(1,2,3,4)/5, log=TRUE) , dbinom(1:20, size=c(20,30), prob=c(1,2,3,4)/5, log=TRUE), "gdbinom")
checkd( gdpois(1:20, c(5,10), log=TRUE) , dpois(1:20, c(5,10), log=TRUE), "gdbinom")
checkd( gqnorm(log((0:10)/10), mean=c(1,2), sd=c(1,2,3,4), log.p=TRUE,warn=FALSE),
qnorm(log((0:10)/10), mean=c(1,2), sd=c(1,2,3,4), log.p=TRUE), "gqnorm")
checkd( gqnorm((0:10)/10, mean=c(1,2), sd=c(1,2,3,4),warn=FALSE),
qnorm((0:10)/10, mean=c(1,2), sd=c(1,2,3,4)), "gqnorm")
checkd( gpnorm(((-10):10)/2, mean=c(1,2), sd=c(1,2,3,4), log.p=TRUE,warn=FALSE),
pnorm(((-10):10)/2, mean=c(1,2), sd=c(1,2,3,4), log.p=TRUE), "gpnorm")
checkd( gpnorm(((-10):10)/2, mean=c(1,2), sd=c(1,2,3,4),warn=FALSE),
pnorm(((-10):10)/2, mean=c(1,2), sd=c(1,2,3,4)), "gpnorm")
if(.Call("cudaVersion")>=7000L) {
cat("Checking solver functions...\n")
checkD=function(x1,y1, funnm) {
yexpr=substitute(y1)
xexpr=substitute(x1)
tmpnm=deparse(expr=xexpr)
mywarnings=character(0)
x=tryCatch(eval(xexpr), error=function(e) {
mywarnings<<-conditionMessage(e)
return(-99999L)})
if(as.logical(x[1]!=-99999L)) {
y=eval(yexpr)
tmp=sum((y-x)^2)/sum(x^2)
if(!is.numeric(tmp)||any(is.na(tmp)))
mywarnings=c(mywarnings,paste(funnm, ":result missing or not numeric."))
else if(tmp>10^(-6))
mywarnings=c(mywarnings,paste(funnm, ":result not equal to host value."))
}
warningslist<<-c(warningslist,mywarnings)
}
#SVD
hm=matrix(rnorm(10^2*2),20,10)
gdm=g(hm)
gsm=g(hm, type="s")
checkD(svd(hm)$d, h(svd(gdm)@S), "svd double")
checkD(svd(hm)$d, h(svd(gsm)@S), "svd single")
#QR
hb=matrix(rnorm(4),20,2)
gdb=g(hb)
gsb=g(hb, type="s")
hcoefqr = qr.coef(qr(hm) ,hb)
checkD(hcoefqr, h(gqr.coef(qr(gdm),gdb)) , "qr/gqr.coef double")
checkD(hcoefqr, h(gqr.coef(qr(gsm),gsb)) , "qr/gqr.coef single")
checkD(hcoefqr, h(gqr.coef(qr(gdm),hb)) , "qr/gqr.coef double / host")
#QR
hm =t(hm) %*% hm
gdm=g(hm)
gsm=g(hm, type="s")
checkD(chol(hm), h(chol(gdm)), "chol double")
checkD(chol(hm), h(chol(gsm)), "chol single")
}
if(length(warningslist)==0)
cat("No errors or warnings\n")
else {
cat("\n\nThe following errors/warnings were detected:\n")
#browser()
for(i in 1:length(warningslist)) {
if(is.null(names(warningslist[i])))
cat("Error:\n")
else if(nchar(names(warningslist[i]))!=0){
cat(names(warningslist[i]),":\n")
} else{
cat("Error:\n")
}
for(err in warningslist[i])
cat("***",err,"\n")
}
}
invisible(TRUE)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.