text:
using Rubin
using Tests
using Elliptic
using HypergeometricFunctions
using Polylogarithms
using SpecialFunctions

@test integrate((c+d*x)^4*cos(a+b*x), x) == :(b^-1*(c+d*x)^4*sin(a+b*x)+24*b^-5*d^4*sin(a+b*x)+-24*b^-4*d^3*(c+d*x)*cos(a+b*x)+-12*b^-3*d^2*(c+d*x)^2*sin(a+b*x)+4*d*b^-2*(c+d*x)^3*cos(a+b*x))
@test integrate((c+d*x)^3*cos(a+b*x), x) == :(b^-1*(c+d*x)^3*sin(a+b*x)+-6*b^-4*d^3*cos(a+b*x)+-6*b^-3*d^2*(c+d*x)*sin(a+b*x)+3*d*b^-2*(c+d*x)^2*cos(a+b*x))
@test integrate((c+d*x)^2*cos(a+b*x), x) == :(b^-1*(c+d*x)^2*sin(a+b*x)+-2*b^-3*d^2*sin(a+b*x)+2*d*b^-2*(c+d*x)*cos(a+b*x))
@test integrate((c+d*x)*cos(a+b*x), x) == :(d*b^-2*cos(a+b*x)+b^-1*(c+d*x)*sin(a+b*x))
@test integrate((c+d*x)^-1*cos(a+b*x), x) == :(d^-1*cos(a+-1*b*c*d^-1)*cosintegral(b*x+b*c*d^-1)+-1*d^-1*sin(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1))
@test integrate((c+d*x)^-2*cos(a+b*x), x) == :(-1*d^-1*(c+d*x)^-1*cos(a+b*x)+-1*b*d^-2*cos(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1)+-1*b*d^-2*cosintegral(b*x+b*c*d^-1)*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^-3*cos(a+b*x), x) == :(-1//2*d^-1*(c+d*x)^-2*cos(a+b*x)+(1/2)*b*d^-2*(c+d*x)^-1*sin(a+b*x)+(1/2)*b^2*d^-3*sin(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1)+-1//2*b^2*d^-3*cos(a+-1*b*c*d^-1)*cosintegral(b*x+b*c*d^-1))
@test integrate((c+d*x)^-4*cos(a+b*x), x) == :(-1//3*d^-1*(c+d*x)^-3*cos(a+b*x)+1//6*b*d^-2*(c+d*x)^-2*sin(a+b*x)+1//6*b^2*d^-3*(c+d*x)^-1*cos(a+b*x)+1//6*b^3*d^-4*cos(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1)+1//6*b^3*d^-4*cosintegral(b*x+b*c*d^-1)*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^4*cos(a+b*x)^2, x) == :(1//10*d^-1*(c+d*x)^5+-1//2*d*b^-2*(c+d*x)^3+3//4*x*b^-4*d^4+d*b^-2*(c+d*x)^3*cos(a+b*x)^2+(1/2)*b^-1*(c+d*x)^4*cos(a+b*x)*sin(a+b*x)+-3//2*b^-4*d^3*cos(a+b*x)^2*(c+d*x)+3//4*b^-5*d^4*cos(a+b*x)*sin(a+b*x)+-3//2*b^-3*d^2*(c+d*x)^2*cos(a+b*x)*sin(a+b*x))
@test integrate((c+d*x)^3*cos(a+b*x)^2, x) == :(1//8*d^-1*(c+d*x)^4+-3//8*b^-4*d^3*cos(a+b*x)^2+-3//8*b^-2*d^3*x^2+(1/2)*b^-1*(c+d*x)^3*cos(a+b*x)*sin(a+b*x)+-3//4*c*x*b^-2*d^2+3//4*d*b^-2*(c+d*x)^2*cos(a+b*x)^2+-3//4*b^-3*d^2*(c+d*x)*cos(a+b*x)*sin(a+b*x))
@test integrate((c+d*x)^2*cos(a+b*x)^2, x) == :(1//6*d^-1*(c+d*x)^3+-1//4*x*b^-2*d^2+(1/2)*d*b^-2*cos(a+b*x)^2*(c+d*x)+(1/2)*b^-1*(c+d*x)^2*cos(a+b*x)*sin(a+b*x)+-1//4*b^-3*d^2*cos(a+b*x)*sin(a+b*x))
@test integrate(cos(a+b*x)^2*(c+d*x), x) == :((1/2)*c*x+1//4*d*x^2+1//4*d*b^-2*cos(a+b*x)^2+(1/2)*b^-1*(c+d*x)*cos(a+b*x)*sin(a+b*x))
@test integrate((c+d*x)^-1*cos(a+b*x)^2, x) == :((1/2)*d^-1*log(c+d*x)+(1/2)*d^-1*cos(2a+-2*b*c*d^-1)*cosintegral(2*b*x+2*b*c*d^-1)+-1//2*d^-1*sin(2a+-2*b*c*d^-1)*sinintegral(2*b*x+2*b*c*d^-1))
@test integrate((c+d*x)^-2*cos(a+b*x)^2, x) == :(-1*d^-1*(c+d*x)^-1*cos(a+b*x)^2+-1*b*d^-2*cos(2a+-2*b*c*d^-1)*sinintegral(2*b*x+2*b*c*d^-1)+-1*b*d^-2*cosintegral(2*b*x+2*b*c*d^-1)*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^-3*cos(a+b*x)^2, x) == :(-1//2*d^-1*(c+d*x)^-2*cos(a+b*x)^2+b^2*d^-3*sin(2a+-2*b*c*d^-1)*sinintegral(2*b*x+2*b*c*d^-1)+-1*b^2*d^-3*cos(2a+-2*b*c*d^-1)*cosintegral(2*b*x+2*b*c*d^-1)+b*d^-2*(c+d*x)^-1*cos(a+b*x)*sin(a+b*x))
@test integrate((c+d*x)^4*cos(a+b*x)^3, x) == :(-8//81*b^-5*d^4*sin(a+b*x)^3+2//3*b^-1*(c+d*x)^4*sin(a+b*x)+488//27*b^-5*d^4*sin(a+b*x)+-160//9*b^-4*d^3*(c+d*x)*cos(a+b*x)+-80//9*b^-3*d^2*(c+d*x)^2*sin(a+b*x)+-8//27*b^-4*d^3*cos(a+b*x)^3*(c+d*x)+1//3*b^-1*(c+d*x)^4*cos(a+b*x)^2*sin(a+b*x)+4//9*d*b^-2*(c+d*x)^3*cos(a+b*x)^3+8//3*d*b^-2*(c+d*x)^3*cos(a+b*x)+-4//9*b^-3*d^2*(c+d*x)^2*cos(a+b*x)^2*sin(a+b*x))
@test integrate((c+d*x)^3*cos(a+b*x)^3, x) == :(-40//9*b^-4*d^3*cos(a+b*x)+-2//27*b^-4*d^3*cos(a+b*x)^3+2//3*b^-1*(c+d*x)^3*sin(a+b*x)+2*d*b^-2*(c+d*x)^2*cos(a+b*x)+-40//9*b^-3*d^2*(c+d*x)*sin(a+b*x)+1//3*d*b^-2*(c+d*x)^2*cos(a+b*x)^3+1//3*b^-1*(c+d*x)^3*cos(a+b*x)^2*sin(a+b*x)+-2//9*b^-3*d^2*cos(a+b*x)^2*(c+d*x)*sin(a+b*x))
@test integrate((c+d*x)^2*cos(a+b*x)^3, x) == :(-14//9*b^-3*d^2*sin(a+b*x)+2//3*b^-1*(c+d*x)^2*sin(a+b*x)+2//27*b^-3*d^2*sin(a+b*x)^3+1//3*b^-1*(c+d*x)^2*cos(a+b*x)^2*sin(a+b*x)+2//9*d*b^-2*cos(a+b*x)^3*(c+d*x)+4//3*d*b^-2*(c+d*x)*cos(a+b*x))
@test integrate(cos(a+b*x)^3*(c+d*x), x) == :(1//3*b^-1*(2c+2*d*x)*sin(a+b*x)+1//9*d*b^-2*cos(a+b*x)^3+2//3*d*b^-2*cos(a+b*x)+1//3*b^-1*cos(a+b*x)^2*(c+d*x)*sin(a+b*x))
@test integrate((c+d*x)^-1*cos(a+b*x)^3, x) == :(-3//4*d^-1*sin(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1)+-1//4*d^-1*sin(3a+-3*b*c*d^-1)*sinintegral(3*b*x+3*b*c*d^-1)+1//4*d^-1*cos(3a+-3*b*c*d^-1)*cosintegral(3*b*x+3*b*c*d^-1)+3//4*d^-1*cos(a+-1*b*c*d^-1)*cosintegral(b*x+b*c*d^-1))
@test integrate((c+d*x)^-2*cos(a+b*x)^3, x) == :(-1*d^-1*(c+d*x)^-1*cos(a+b*x)^3+-3//4*b*d^-2*cos(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1)+-3//4*b*d^-2*cos(3a+-3*b*c*d^-1)*sinintegral(3*b*x+3*b*c*d^-1)+-3//4*b*d^-2*cosintegral(b*x+b*c*d^-1)*sin(a+-1*b*c*d^-1)+-3//4*b*d^-2*cosintegral(3*b*x+3*b*c*d^-1)*sin(3a+-3*b*c*d^-1))
@test integrate((c+d*x)^-3*cos(a+b*x)^3, x) == :(-1//2*d^-1*(c+d*x)^-2*cos(a+b*x)^3+-9//8*b^2*d^-3*cos(3a+-3*b*c*d^-1)*cosintegral(3*b*x+3*b*c*d^-1)+-3//8*b^2*d^-3*cos(a+-1*b*c*d^-1)*cosintegral(b*x+b*c*d^-1)+3//8*b^2*d^-3*sin(a+-1*b*c*d^-1)*sinintegral(b*x+b*c*d^-1)+9//8*b^2*d^-3*sin(3a+-3*b*c*d^-1)*sinintegral(3*b*x+3*b*c*d^-1)+3//2*b*d^-2*(c+d*x)^-1*cos(a+b*x)^2*sin(a+b*x))
@test integrate(x^3*cos(a+b*x)^4, x) == :(3//32*x^4+-45//128*b^-4*cos(a+b*x)^2+-45//128*b^-2*x^2+-3//128*b^-4*cos(a+b*x)^4+3//16*b^-2*x^2*cos(a+b*x)^4+9//16*b^-2*x^2*cos(a+b*x)^2+-45//64*x*b^-3*cos(a+b*x)*sin(a+b*x)+-3//32*x*b^-3*cos(a+b*x)^3*sin(a+b*x)+1//4*b^-1*x^3*cos(a+b*x)^3*sin(a+b*x)+3//8*b^-1*x^3*cos(a+b*x)*sin(a+b*x))
@test integrate(x^2*cos(a+b*x)^4, x) == :(1//8*x^3+-15//64*x*b^-2+-15//64*b^-3*cos(a+b*x)*sin(a+b*x)+-1//32*b^-3*cos(a+b*x)^3*sin(a+b*x)+1//8*x*b^-2*cos(a+b*x)^4+3//8*x*b^-2*cos(a+b*x)^2+1//4*b^-1*x^2*cos(a+b*x)^3*sin(a+b*x)+3//8*b^-1*x^2*cos(a+b*x)*sin(a+b*x))
@test integrate(x*cos(a+b*x)^4, x) == :(3//16*x^2+1//16*b^-2*cos(a+b*x)^4+3//16*b^-2*cos(a+b*x)^2+1//4*x*b^-1*cos(a+b*x)^3*sin(a+b*x)+3//8*x*b^-1*cos(a+b*x)*sin(a+b*x))
@test integrate(x^-1*cos(a+b*x)^4, x) == :(3//8*log(x)+(1/2)*cos(2a)*cosintegral(2*b*x)+-1//2*sin(2a)*sinintegral(2*b*x)+-1//8*sin(4a)*sinintegral(4*b*x)+1//8*cos(4a)*cosintegral(4*b*x))
@test integrate(x^-2*cos(a+b*x)^4, x) == :(-1*x^-1*cos(a+b*x)^4+-1*b*cos(2a)*sinintegral(2*b*x)+-1*b*cosintegral(2*b*x)*sin(2a)+-1//2*b*cos(4a)*sinintegral(4*b*x)+-1//2*b*cosintegral(4*b*x)*sin(4a))
@test integrate(x^-3*cos(a+b*x)^4, x) == :(-1//2*x^-2*cos(a+b*x)^4+b^2*sin(2a)*sinintegral(2*b*x)+b^2*sin(4a)*sinintegral(4*b*x)+-1*b^2*cos(2a)*cosintegral(2*b*x)+-1*b^2*cos(4a)*cosintegral(4*b*x)+2*b*x^-1*cos(a+b*x)^3*sin(a+b*x))
@test integrate((c+d*x)^3*sec(a+b*x), x) == :(-6*im*b^-4*d^3*Polylogarithms.polylog(4,-1*im*exp(im*(a+b*x)))+-6*b^-3*d^2*(c+d*x)*Polylogarithms.polylog(3,-1*im*exp(im*(a+b*x)))+-2*im*b^-1*(c+d*x)^3*arctan(exp(im*(a+b*x)))+6*im*b^-4*d^3*Polylogarithms.polylog(4,im*exp(im*(a+b*x)))+6*b^-3*d^2*(c+d*x)*Polylogarithms.polylog(3,im*exp(im*(a+b*x)))+-3*im*d*b^-2*(c+d*x)^2*Polylogarithms.polylog(2,im*exp(im*(a+b*x)))+3*im*d*b^-2*(c+d*x)^2*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x))))
@test integrate((c+d*x)^2*sec(a+b*x), x) == :(-2*b^-3*d^2*Polylogarithms.polylog(3,-1*im*exp(im*(a+b*x)))+2*b^-3*d^2*Polylogarithms.polylog(3,im*exp(im*(a+b*x)))+-2*im*b^-1*(c+d*x)^2*arctan(exp(im*(a+b*x)))+-2*im*d*b^-2*(c+d*x)*Polylogarithms.polylog(2,im*exp(im*(a+b*x)))+2*im*d*b^-2*(c+d*x)*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x))))
@test integrate((c+d*x)*sec(a+b*x), x) == :(im*d*b^-2*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x)))+-1*im*d*b^-2*Polylogarithms.polylog(2,im*exp(im*(a+b*x)))+-2*im*b^-1*(c+d*x)*arctan(exp(im*(a+b*x))))
@test integrate((c+d*x)^3*sec(a+b*x)^2, x) == :(b^-1*(c+d*x)^3*tan(a+b*x)+-1*im*b^-1*(c+d*x)^3+3//2*b^-4*d^3*Polylogarithms.polylog(3,-1*exp(2*im*(a+b*x)))+3*d*b^-2*(c+d*x)^2*log(1+exp(2*im*(a+b*x)))+-3*im*b^-3*d^2*(c+d*x)*Polylogarithms.polylog(2,-1*exp(2*im*(a+b*x))))
@test integrate((c+d*x)^2*sec(a+b*x)^2, x) == :(b^-1*(c+d*x)^2*tan(a+b*x)+-1*im*b^-1*(c+d*x)^2+-1*im*b^-3*d^2*Polylogarithms.polylog(2,-1*exp(2*im*(a+b*x)))+2*d*b^-2*(c+d*x)*log(1+exp(2*im*(a+b*x))))
@test integrate(sec(a+b*x)^2*(c+d*x), x) == :(d*b^-2*log(cos(a+b*x))+b^-1*(c+d*x)*tan(a+b*x))
@test integrate((c+d*x)^3*sec(a+b*x)^3, x) == :((1/2)*b^-1*(c+d*x)^3*sec(a+b*x)*tan(a+b*x)+-1*im*b^-1*(c+d*x)^3*arctan(exp(im*(a+b*x)))+-3*im*b^-4*d^3*Polylogarithms.polylog(2,im*exp(im*(a+b*x)))+-3*im*b^-4*d^3*Polylogarithms.polylog(4,-1*im*exp(im*(a+b*x)))+-3*b^-3*d^2*(c+d*x)*Polylogarithms.polylog(3,-1*im*exp(im*(a+b*x)))+3*im*b^-4*d^3*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x)))+3*im*b^-4*d^3*Polylogarithms.polylog(4,im*exp(im*(a+b*x)))+3*b^-3*d^2*(c+d*x)*Polylogarithms.polylog(3,im*exp(im*(a+b*x)))+-3//2*d*b^-2*(c+d*x)^2*sec(a+b*x)+-6*im*b^-3*d^2*(c+d*x)*arctan(exp(im*(a+b*x)))+-3//2*im*d*b^-2*(c+d*x)^2*Polylogarithms.polylog(2,im*exp(im*(a+b*x)))+3//2*im*d*b^-2*(c+d*x)^2*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x))))
@test integrate((c+d*x)^2*sec(a+b*x)^3, x) == :(b^-3*d^2*Polylogarithms.polylog(3,im*exp(im*(a+b*x)))+b^-3*d^2*arctanh(sin(a+b*x))+-1*b^-3*d^2*Polylogarithms.polylog(3,-1*im*exp(im*(a+b*x)))+(1/2)*b^-1*(c+d*x)^2*sec(a+b*x)*tan(a+b*x)+-1*im*b^-1*(c+d*x)^2*arctan(exp(im*(a+b*x)))+-1*d*b^-2*(c+d*x)*sec(a+b*x)+im*d*b^-2*(c+d*x)*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x)))+-1*im*d*b^-2*(c+d*x)*Polylogarithms.polylog(2,im*exp(im*(a+b*x))))
@test integrate(sec(a+b*x)^3*(c+d*x), x) == :(-1//2*d*b^-2*sec(a+b*x)+(1/2)*im*d*b^-2*Polylogarithms.polylog(2,-1*im*exp(im*(a+b*x)))+(1/2)*b^-1*(c+d*x)*sec(a+b*x)*tan(a+b*x)+-1*im*b^-1*(c+d*x)*arctan(exp(im*(a+b*x)))+-1//2*im*d*b^-2*Polylogarithms.polylog(2,im*exp(im*(a+b*x))))
@test integrate((c+d*x)^5//2*cos(a+b*x), x) == :(b^-1*(c+d*x)^5//2*sin(a+b*x)+-15//4*b^-3*d^2*(c+d*x)^(1/2)*sin(a+b*x)+5//2*d*b^-2*(c+d*x)^3//2*cos(a+b*x)+15//8*2^(1/2)*pi^(1/2)*b^-7//2*d^5//2*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+15//8*2^(1/2)*pi^(1/2)*b^-7//2*d^5//2*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^3//2*cos(a+b*x), x) == :(b^-1*(c+d*x)^3//2*sin(a+b*x)+3//2*d*b^-2*(c+d*x)^(1/2)*cos(a+b*x)+-3//4*2^(1/2)*pi^(1/2)*b^-5//2*d^3//2*cos(a+-1*b*c*d^-1)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+3//4*2^(1/2)*pi^(1/2)*b^-5//2*d^3//2*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^(1/2)*cos(a+b*x), x) == :(b^-1*(c+d*x)^(1/2)*sin(a+b*x)+-1//2*2^(1/2)*pi^(1/2)*b^-3//2*d^(1/2)*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1//2*2^(1/2)*pi^(1/2)*b^-3//2*d^(1/2)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^-1//2*cos(a+b*x), x) == :(2^(1/2)*pi^(1/2)*b^-1//2*d^-1//2*cos(a+-1*b*c*d^-1)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1*2^(1/2)*pi^(1/2)*b^-1//2*d^-1//2*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^-3//2*cos(a+b*x), x) == :(-2*d^-1*(c+d*x)^-1//2*cos(a+b*x)+-2*2^(1/2)*pi^(1/2)*b^(1/2)*d^-3//2*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-2*2^(1/2)*pi^(1/2)*b^(1/2)*d^-3//2*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^-5//2*cos(a+b*x), x) == :(-2//3*d^-1*(c+d*x)^-3//2*cos(a+b*x)+4//3*b*d^-2*(c+d*x)^-1//2*sin(a+b*x)+-4//3*2^(1/2)*pi^(1/2)*b^3//2*d^-5//2*cos(a+-1*b*c*d^-1)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+4//3*2^(1/2)*pi^(1/2)*b^3//2*d^-5//2*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^-7//2*cos(a+b*x), x) == :(-2//5*d^-1*(c+d*x)^-5//2*cos(a+b*x)+4//15*b*d^-2*(c+d*x)^-3//2*sin(a+b*x)+8//15*b^2*d^-3*(c+d*x)^-1//2*cos(a+b*x)+8//15*2^(1/2)*pi^(1/2)*b^5//2*d^-7//2*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+8//15*2^(1/2)*pi^(1/2)*b^5//2*d^-7//2*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^5//2*cos(a+b*x)^2, x) == :(1//7*d^-1*(c+d*x)^7//2+-5//16*d*b^-2*(c+d*x)^3//2+(1/2)*b^-1*(c+d*x)^5//2*cos(a+b*x)*sin(a+b*x)+-15//64*b^-3*d^2*(c+d*x)^(1/2)*sin(2a+2*b*x)+5//8*d*b^-2*(c+d*x)^3//2*cos(a+b*x)^2+15//128*pi^(1/2)*b^-7//2*d^5//2*cos(2a+-2*b*c*d^-1)*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+15//128*pi^(1/2)*b^-7//2*d^5//2*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^3//2*cos(a+b*x)^2, x) == :(1//5*d^-1*(c+d*x)^5//2+-3//16*d*b^-2*(c+d*x)^(1/2)+(1/2)*b^-1*(c+d*x)^3//2*cos(a+b*x)*sin(a+b*x)+3//8*d*b^-2*(c+d*x)^(1/2)*cos(a+b*x)^2+-3//32*pi^(1/2)*b^-5//2*d^3//2*cos(2a+-2*b*c*d^-1)*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+3//32*pi^(1/2)*b^-5//2*d^3//2*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^(1/2)*cos(a+b*x)^2, x) == :(1//3*d^-1*(c+d*x)^3//2+1//4*b^-1*(c+d*x)^(1/2)*sin(2a+2*b*x)+-1//8*pi^(1/2)*b^-3//2*d^(1/2)*cos(2a+-2*b*c*d^-1)*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1//8*pi^(1/2)*b^-3//2*d^(1/2)*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^-1//2*cos(a+b*x)^2, x) == :(d^-1*(c+d*x)^(1/2)+(1/2)*pi^(1/2)*b^-1//2*d^-1//2*cos(2a+-2*b*c*d^-1)*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1//2*pi^(1/2)*b^-1//2*d^-1//2*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^-3//2*cos(a+b*x)^2, x) == :(-2*d^-1*(c+d*x)^-1//2*cos(a+b*x)^2+-2*pi^(1/2)*b^(1/2)*d^-3//2*cos(2a+-2*b*c*d^-1)*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-2*pi^(1/2)*b^(1/2)*d^-3//2*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^-5//2*cos(a+b*x)^2, x) == :(-2//3*d^-1*(c+d*x)^-3//2*cos(a+b*x)^2+-8//3*pi^(1/2)*b^3//2*d^-5//2*cos(2a+-2*b*c*d^-1)*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+8//3*b*d^-2*(c+d*x)^-1//2*cos(a+b*x)*sin(a+b*x)+8//3*pi^(1/2)*b^3//2*d^-5//2*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^-7//2*cos(a+b*x)^2, x) == :(-16//15*b^2*d^-3*(c+d*x)^-1//2+-2//5*d^-1*(c+d*x)^-5//2*cos(a+b*x)^2+32//15*b^2*d^-3*(c+d*x)^-1//2*cos(a+b*x)^2+8//15*b*d^-2*(c+d*x)^-3//2*cos(a+b*x)*sin(a+b*x)+32//15*pi^(1/2)*b^5//2*d^-7//2*cos(2a+-2*b*c*d^-1)*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+32//15*pi^(1/2)*b^5//2*d^-7//2*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1))
@test integrate((c+d*x)^-9//2*cos(a+b*x)^2, x) == :(-16//105*b^2*d^-3*(c+d*x)^-3//2+-2//7*d^-1*(c+d*x)^-7//2*cos(a+b*x)^2+32//105*b^2*d^-3*(c+d*x)^-3//2*cos(a+b*x)^2+-128//105*pi^(1/2)*b^7//2*d^-9//2*fresnels(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(2a+-2*b*c*d^-1)+-128//105*b^3*d^-4*(c+d*x)^-1//2*cos(a+b*x)*sin(a+b*x)+8//35*b*d^-2*(c+d*x)^-5//2*cos(a+b*x)*sin(a+b*x)+128//105*pi^(1/2)*b^7//2*d^-9//2*cos(2a+-2*b*c*d^-1)*fresnelc(2*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2)))
@test integrate((c+d*x)^5//2*cos(a+b*x)^3, x) == :(2//3*b^-1*(c+d*x)^5//2*sin(a+b*x)+-45//16*b^-3*d^2*(c+d*x)^(1/2)*sin(a+b*x)+-5//144*b^-3*d^2*(c+d*x)^(1/2)*sin(3a+3*b*x)+1//3*b^-1*(c+d*x)^5//2*cos(a+b*x)^2*sin(a+b*x)+5//3*d*b^-2*(c+d*x)^3//2*cos(a+b*x)+5//18*d*b^-2*(c+d*x)^3//2*cos(a+b*x)^3+5//864*6^(1/2)*pi^(1/2)*b^-7//2*d^5//2*cos(3a+-3*b*c*d^-1)*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+5//864*6^(1/2)*pi^(1/2)*b^-7//2*d^5//2*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1)+45//32*2^(1/2)*pi^(1/2)*b^-7//2*d^5//2*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+45//32*2^(1/2)*pi^(1/2)*b^-7//2*d^5//2*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^3//2*cos(a+b*x)^3, x) == :(2//3*b^-1*(c+d*x)^3//2*sin(a+b*x)+d*b^-2*(c+d*x)^(1/2)*cos(a+b*x)+1//3*b^-1*(c+d*x)^3//2*cos(a+b*x)^2*sin(a+b*x)+1//6*d*b^-2*(c+d*x)^(1/2)*cos(a+b*x)^3+-9//16*2^(1/2)*pi^(1/2)*b^-5//2*d^3//2*cos(a+-1*b*c*d^-1)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1//144*6^(1/2)*pi^(1/2)*b^-5//2*d^3//2*cos(3a+-3*b*c*d^-1)*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+1//144*6^(1/2)*pi^(1/2)*b^-5//2*d^3//2*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1)+9//16*2^(1/2)*pi^(1/2)*b^-5//2*d^3//2*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1))
@test integrate((c+d*x)^(1/2)*cos(a+b*x)^3, x) == :(1//12*b^-1*(c+d*x)^(1/2)*sin(3a+3*b*x)+3//4*b^-1*(c+d*x)^(1/2)*sin(a+b*x)+-3//8*2^(1/2)*pi^(1/2)*b^-3//2*d^(1/2)*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-3//8*2^(1/2)*pi^(1/2)*b^-3//2*d^(1/2)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1)+-1//72*6^(1/2)*pi^(1/2)*b^-3//2*d^(1/2)*cos(3a+-3*b*c*d^-1)*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1//72*6^(1/2)*pi^(1/2)*b^-3//2*d^(1/2)*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1))
@test integrate((c+d*x)^-1//2*cos(a+b*x)^3, x) == :(-3//4*2^(1/2)*pi^(1/2)*b^-1//2*d^-1//2*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1)+-1//12*6^(1/2)*pi^(1/2)*b^-1//2*d^-1//2*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1)+1//12*6^(1/2)*pi^(1/2)*b^-1//2*d^-1//2*cos(3a+-3*b*c*d^-1)*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+3//4*2^(1/2)*pi^(1/2)*b^-1//2*d^-1//2*cos(a+-1*b*c*d^-1)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2)))
@test integrate((c+d*x)^-3//2*cos(a+b*x)^3, x) == :(-2*d^-1*(c+d*x)^-1//2*cos(a+b*x)^3+-3//2*2^(1/2)*pi^(1/2)*b^(1/2)*d^-3//2*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-3//2*2^(1/2)*pi^(1/2)*b^(1/2)*d^-3//2*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1)+-1//2*6^(1/2)*pi^(1/2)*b^(1/2)*d^-3//2*cos(3a+-3*b*c*d^-1)*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1//2*6^(1/2)*pi^(1/2)*b^(1/2)*d^-3//2*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1))
@test integrate((c+d*x)^-5//2*cos(a+b*x)^3, x) == :(-2//3*d^-1*(c+d*x)^-3//2*cos(a+b*x)^3+4*b*d^-2*(c+d*x)^-1//2*cos(a+b*x)^2*sin(a+b*x)+2^(1/2)*pi^(1/2)*b^3//2*d^-5//2*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1)+6^(1/2)*pi^(1/2)*b^3//2*d^-5//2*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1)+-1*2^(1/2)*pi^(1/2)*b^3//2*d^-5//2*cos(a+-1*b*c*d^-1)*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+-1*6^(1/2)*pi^(1/2)*b^3//2*d^-5//2*cos(3a+-3*b*c*d^-1)*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2)))
@test integrate((c+d*x)^-7//2*cos(a+b*x)^3, x) == :(-2//5*d^-1*(c+d*x)^-5//2*cos(a+b*x)^3+-16//5*b^2*d^-3*(c+d*x)^-1//2*cos(a+b*x)+24//5*b^2*d^-3*(c+d*x)^-1//2*cos(a+b*x)^3+4//5*b*d^-2*(c+d*x)^-3//2*cos(a+b*x)^2*sin(a+b*x)+2//5*2^(1/2)*pi^(1/2)*b^5//2*d^-7//2*cos(a+-1*b*c*d^-1)*fresnels(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+2//5*2^(1/2)*pi^(1/2)*b^5//2*d^-7//2*fresnelc(2^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(a+-1*b*c*d^-1)+6//5*6^(1/2)*pi^(1/2)*b^5//2*d^-7//2*cos(3a+-3*b*c*d^-1)*fresnels(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))+6//5*6^(1/2)*pi^(1/2)*b^5//2*d^-7//2*fresnelc(6^(1/2)*pi^-1//2*b^(1/2)*d^-1//2*(c+d*x)^(1/2))*sin(3a+-3*b*c*d^-1))
@test integrate(x^3//2*cos(x), x) == :(x^3//2*sin(x)+3//2*x^(1/2)*cos(x)+-3//4*2^(1/2)*pi^(1/2)*fresnelc(2^(1/2)*pi^-1//2*x^(1/2)))
@test integrate(x^(1/2)*cos(x), x) == :(x^(1/2)*sin(x)+-1//2*2^(1/2)*pi^(1/2)*fresnels(2^(1/2)*pi^-1//2*x^(1/2)))
@test integrate(x^-1//2*cos(x), x) == :(2^(1/2)*pi^(1/2)*fresnelc(2^(1/2)*pi^-1//2*x^(1/2)))
@test integrate(x^-3//2*cos(x), x) == :(-2*x^-1//2*cos(x)+-2*2^(1/2)*pi^(1/2)*fresnels(2^(1/2)*pi^-1//2*x^(1/2)))
@test integrate((c+d*x)^4//3*cos(a+b*x), x) == :(b^-1*(c+d*x)^4//3*sin(a+b*x)+4//3*d*b^-2*(c+d*x)^1//3*cos(a+b*x)+-2//9*im*b^-3*d^2*(im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+2//9*im*b^-3*d^2*(-1*im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^2//3*cos(a+b*x), x) == :(b^-1*(c+d*x)^2//3*sin(a+b*x)+1//3*d*b^-2*(im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+1//3*d*b^-2*(-1*im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^1//3*cos(a+b*x), x) == :(b^-1*(c+d*x)^1//3*sin(a+b*x)+1//6*d*b^-2*(im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+1//6*d*b^-2*(-1*im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^-1//3*cos(a+b*x), x) == :((1/2)*im*b^-1*(im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+-1//2*im*b^-1*(-1*im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^-2//3*cos(a+b*x), x) == :((1/2)*im*b^-1*(im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+-1//2*im*b^-1*(-1*im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^-4//3*cos(a+b*x), x) == :(-3*d^-1*(c+d*x)^-1//3*cos(a+b*x)+3//2*d^-1*(im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+3//2*d^-1*(-1*im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^-5//3*cos(a+b*x), x) == :(-3//2*d^-1*(c+d*x)^-2//3*cos(a+b*x)+3//4*d^-1*(im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+3//4*d^-1*(-1*im*b*d^-1*(c+d*x))^2//3*(c+d*x)^-2//3*SpecialFunctions.gamma(1//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^-7//3*cos(a+b*x), x) == :(-3//4*d^-1*(c+d*x)^-4//3*cos(a+b*x)+9//4*b*d^-2*(c+d*x)^-1//3*sin(a+b*x)+-9//8*im*b*d^-2*(im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+9//8*im*b*d^-2*(-1*im*b*d^-1*(c+d*x))^1//3*(c+d*x)^-1//3*SpecialFunctions.gamma(2//3,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate(cos(a+b*x)^(1/2), x) == :(2*b^-1*Elliptic.E((1/2)*a+(1/2)*b*x,2))
@test integrate(cos(a+b*x)^3//2, x) == :(2//3*b^-1*Elliptic.F((1/2)*a+(1/2)*b*x,2)+2//3*b^-1*cos(a+b*x)^(1/2)*sin(a+b*x))
@test integrate(x*cos(a+b*x)^3//2+-1//3*x*cos(a+b*x)^-1//2, x) == :(4//9*b^-2*cos(a+b*x)^3//2+2//3*x*b^-1*cos(a+b*x)^(1/2)*sin(a+b*x))
@test integrate(cos(a+b*x)^-1//2, x) == :(2*b^-1*Elliptic.F((1/2)*a+(1/2)*b*x,2))
@test integrate(cos(a+b*x)^-3//2, x) == :(-2*b^-1*Elliptic.E((1/2)*a+(1/2)*b*x,2)+2*b^-1*cos(a+b*x)^-1//2*sin(a+b*x))
@test integrate(x*cos(a+b*x)^(1/2)+x*cos(a+b*x)^-3//2, x) == :(4*b^-2*cos(a+b*x)^(1/2)+2*x*b^-1*cos(a+b*x)^-1//2*sin(a+b*x))
@test integrate(x*cos(x)^(1/2)+x*cos(x)^-3//2, x) == :(4*cos(x)^(1/2)+2*x*cos(x)^-1//2*sin(x))
@test integrate(x*cos(x)^-5//2+-1//3*x*cos(x)^-1//2, x) == :(-4//3*cos(x)^-1//2+2//3*x*cos(x)^-3//2*sin(x))
@test integrate(x*cos(x)^-7//2+3//5*x*cos(x)^(1/2), x) == :(-4//15*cos(x)^-3//2+12//5*cos(x)^(1/2)+2//5*x*cos(x)^-5//2*sin(x)+6//5*x*cos(x)^-1//2*sin(x))
@test integrate(x^2*cos(x)^(1/2)+x^2*cos(x)^-3//2, x) == :(-16*Elliptic.E((1/2)*x,2)+8*x*cos(x)^(1/2)+2*x^2*cos(x)^-1//2*sin(x))
@test integrate(x*sec(x)^-3//2+-1//3*x*sec(x)^(1/2), x) == :(4//9*sec(x)^-3//2+2//3*x*sec(x)^-1//2*sin(x))
@test integrate(x*sec(x)^-5//2+-3//5*x*sec(x)^-1//2, x) == :(4//25*sec(x)^-5//2+2//5*x*sec(x)^-3//2*sin(x))
@test integrate(x*sec(x)^-7//2+-5//21*x*sec(x)^(1/2), x) == :(4//49*sec(x)^-7//2+20//63*sec(x)^-3//2+2//7*x*sec(x)^-5//2*sin(x)+10//21*x*sec(x)^-1//2*sin(x))
@test integrate(x^2*sec(x)^-3//2+-1//3*x^2*sec(x)^(1/2), x) == :(-16//27*sec(x)^-1//2*sin(x)+8//9*x*sec(x)^-3//2+-16//27*cos(x)^(1/2)*sec(x)^(1/2)*Elliptic.F((1/2)*x,2)+2//3*x^2*sec(x)^-1//2*sin(x))
@test integrate((c+d*x)^m*cos(a+b*x)^3, x) == :(-3//8*im*b^-1*(-1*im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1))+3//8*im*b^-1*(im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+-1//8*im*3^(-1+-1m)*b^-1*(-1*im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,-3*im*b*d^-1*(c+d*x))*exp(3*im*(a+-1*b*c*d^-1))+1//8*im*3^(-1+-1m)*b^-1*(im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,3*im*b*d^-1*(c+d*x))*exp(-3*im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^m*cos(a+b*x)^2, x) == :((1/2)*d^-1*(1+m)^-1*(c+d*x)^(1+m)+im*2^(-3+-1m)*b^-1*(im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,2*im*b*d^-1*(c+d*x))*exp(-2*im*(a+-1*b*c*d^-1))+-1*im*2^(-3+-1m)*b^-1*(-1*im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,-2*im*b*d^-1*(c+d*x))*exp(2*im*(a+-1*b*c*d^-1)))
@test integrate((c+d*x)^m*cos(a+b*x), x) == :((1/2)*im*b^-1*(im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,im*b*d^-1*(c+d*x))*exp(-1*im*(a+-1*b*c*d^-1))+-1//2*im*b^-1*(-1*im*b*d^-1*(c+d*x))^(-1m)*(c+d*x)^m*SpecialFunctions.gamma(1+m,-1*im*b*d^-1*(c+d*x))*exp(im*(a+-1*b*c*d^-1)))
@test integrate(x^(3+m)*cos(a+b*x), x) == :(-1//2*b^-4*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(4+m,im*b*x)*exp(-1*im*a)+-1//2*b^-4*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(4+m,-1*im*b*x)*exp(im*a))
@test integrate(x^(2+m)*cos(a+b*x), x) == :((1/2)*im*b^-3*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(3+m,-1*im*b*x)*exp(im*a)+-1//2*im*b^-3*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(3+m,im*b*x)*exp(-1*im*a))
@test integrate(x^(1+m)*cos(a+b*x), x) == :((1/2)*b^-2*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(2+m,im*b*x)*exp(-1*im*a)+(1/2)*b^-2*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(2+m,-1*im*b*x)*exp(im*a))
@test integrate(x^m*cos(a+b*x), x) == :((1/2)*im*b^-1*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(1+m,im*b*x)*exp(-1*im*a)+-1//2*im*b^-1*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(1+m,-1*im*b*x)*exp(im*a))
@test integrate(x^(-1+m)*cos(a+b*x), x) == :(-1//2*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(m,im*b*x)*exp(-1*im*a)+-1//2*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(m,-1*im*b*x)*exp(im*a))
@test integrate(x^(-2+m)*cos(a+b*x), x) == :((1/2)*im*b*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(-1+m,-1*im*b*x)*exp(im*a)+-1//2*im*b*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(-1+m,im*b*x)*exp(-1*im*a))
@test integrate(x^(-3+m)*cos(a+b*x), x) == :((1/2)*b^2*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(-2+m,im*b*x)*exp(-1*im*a)+(1/2)*b^2*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(-2+m,-1*im*b*x)*exp(im*a))
@test integrate(x^(3+m)*cos(a+b*x)^2, x) == :(x^(4+m)*(8+2m)^-1+-1*2^(-6+-1m)*b^-4*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(4+m,2*im*b*x)*exp(-2*im*a)+-1*2^(-6+-1m)*b^-4*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(4+m,-2*im*b*x)*exp(2*im*a))
@test integrate(x^(2+m)*cos(a+b*x)^2, x) == :(x^(3+m)*(6+2m)^-1+im*2^(-5+-1m)*b^-3*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(3+m,-2*im*b*x)*exp(2*im*a)+-1*im*2^(-5+-1m)*b^-3*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(3+m,2*im*b*x)*exp(-2*im*a))
@test integrate(x^(1+m)*cos(a+b*x)^2, x) == :(x^(2+m)*(4+2m)^-1+2^(-4+-1m)*b^-2*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(2+m,2*im*b*x)*exp(-2*im*a)+2^(-4+-1m)*b^-2*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(2+m,-2*im*b*x)*exp(2*im*a))
@test integrate(x^m*cos(a+b*x)^2, x) == :(x^(1+m)*(2+2m)^-1+im*2^(-3+-1m)*b^-1*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(1+m,2*im*b*x)*exp(-2*im*a)+-1*im*2^(-3+-1m)*b^-1*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(1+m,-2*im*b*x)*exp(2*im*a))
@test integrate(x^(-1+m)*cos(a+b*x)^2, x) == :((1/2)*m^-1*x^m+-1*2^(-2+-1m)*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(m,2*im*b*x)*exp(-2*im*a)+-1*2^(-2+-1m)*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(m,-2*im*b*x)*exp(2*im*a))
@test integrate(x^(-2+m)*cos(a+b*x)^2, x) == :(-1*x^(-1+m)*(2+-2m)^-1+im*b*2^(-1+-1m)*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(-1+m,-2*im*b*x)*exp(2*im*a)+-1*im*b*2^(-1+-1m)*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(-1+m,2*im*b*x)*exp(-2*im*a))
@test integrate(x^(-3+m)*cos(a+b*x)^2, x) == :(-1*x^(-2+m)*(4+-2m)^-1+2^(-1m)*b^2*x^m*(im*b*x)^(-1m)*SpecialFunctions.gamma(-2+m,2*im*b*x)*exp(-2*im*a)+2^(-1m)*b^2*x^m*(-1*im*b*x)^(-1m)*SpecialFunctions.gamma(-2+m,-2*im*b*x)*exp(2*im*a))
@test integrate((c+d*x)^3*(a+a*cos(e+f*x)), x) == :(1//4*a*d^-1*(c+d*x)^4+a*f^-1*(c+d*x)^3*sin(e+f*x)+-6*a*d^3*f^-4*cos(e+f*x)+-6*a*d^2*f^-3*(c+d*x)*sin(e+f*x)+3*a*d*f^-2*(c+d*x)^2*cos(e+f*x))
@test integrate((c+d*x)^2*(a+a*cos(e+f*x)), x) == :(1//3*a*d^-1*(c+d*x)^3+a*f^-1*(c+d*x)^2*sin(e+f*x)+-2*a*d^2*f^-3*sin(e+f*x)+2*a*d*f^-2*(c+d*x)*cos(e+f*x))
@test integrate((a+a*cos(e+f*x))*(c+d*x), x) == :((1/2)*a*d^-1*(c+d*x)^2+a*d*f^-2*cos(e+f*x)+a*f^-1*(c+d*x)*sin(e+f*x))
@test integrate((c+d*x)^-1*(a+a*cos(e+f*x)), x) == :(a*d^-1*log(c+d*x)+a*d^-1*cos(-1*e+c*f*d^-1)*cosintegral(f*x+c*f*d^-1)+a*d^-1*sin(-1*e+c*f*d^-1)*sinintegral(f*x+c*f*d^-1))
@test integrate((c+d*x)^-2*(a+a*cos(e+f*x)), x) == :(-1*a*d^-1*(c+d*x)^-1+-1*a*d^-1*(c+d*x)^-1*cos(e+f*x)+a*f*d^-2*cosintegral(f*x+c*f*d^-1)*sin(-1*e+c*f*d^-1)+-1*a*f*d^-2*cos(-1*e+c*f*d^-1)*sinintegral(f*x+c*f*d^-1))
@test integrate((a+a*cos(e+f*x))^2*(c+d*x)^3, x) == :(3//8*a^2*d^-1*(c+d*x)^4+-12*a^2*d^3*f^-4*cos(e+f*x)+2*a^2*f^-1*(c+d*x)^3*sin(e+f*x)+-3//8*a^2*d^3*f^-4*cos(e+f*x)^2+-3//8*a^2*d^3*f^-2*x^2+(1/2)*a^2*f^-1*(c+d*x)^3*cos(e+f*x)*sin(e+f*x)+-12*a^2*d^2*f^-3*(c+d*x)*sin(e+f*x)+6*d*a^2*f^-2*(c+d*x)^2*cos(e+f*x)+-3//4*c*x*a^2*d^2*f^-2+3//4*d*a^2*f^-2*(c+d*x)^2*cos(e+f*x)^2+-3//4*a^2*d^2*f^-3*(c+d*x)*cos(e+f*x)*sin(e+f*x))
@test integrate((a+a*cos(e+f*x))^2*(c+d*x)^2, x) == :((1/2)*a^2*d^-1*(c+d*x)^3+-4*a^2*d^2*f^-3*sin(e+f*x)+2*a^2*f^-1*(c+d*x)^2*sin(e+f*x)+-1//4*x*a^2*d^2*f^-2+(1/2)*d*a^2*f^-2*cos(e+f*x)^2*(c+d*x)+(1/2)*a^2*f^-1*(c+d*x)^2*cos(e+f*x)*sin(e+f*x)+4*d*a^2*f^-2*(c+d*x)*cos(e+f*x)+-1//4*a^2*d^2*f^-3*cos(e+f*x)*sin(e+f*x))
@test integrate((a+a*cos(e+f*x))^2*(c+d*x), x) == :((1/2)*c*x*a^2+(1/2)*a^2*d^-1*(c+d*x)^2+1//4*d*a^2*x^2+2*d*a^2*f^-2*cos(e+f*x)+2*a^2*f^-1*(c+d*x)*sin(e+f*x)+1//4*d*a^2*f^-2*cos(e+f*x)^2+(1/2)*a^2*f^-1*(c+d*x)*cos(e+f*x)*sin(e+f*x))
@test integrate((a+a*cos(e+f*x))^2*(c+d*x)^-1, x) == :(3//2*a^2*d^-1*log(c+d*x)+(1/2)*a^2*d^-1*cos(-2*e+2*c*f*d^-1)*cosintegral(2*f*x+2*c*f*d^-1)+(1/2)*a^2*d^-1*sin(-2*e+2*c*f*d^-1)*sinintegral(2*f*x+2*c*f*d^-1)+2*a^2*d^-1*cos(-1*e+c*f*d^-1)*cosintegral(f*x+c*f*d^-1)+2*a^2*d^-1*sin(-1*e+c*f*d^-1)*sinintegral(f*x+c*f*d^-1))
@test integrate((a+a*cos(e+f*x))^2*(c+d*x)^-2, x) == :(-4*a^2*d^-1*(c+d*x)^-1*cos((1/2)*e+(1/2)*f*x)^4+f*a^2*d^-2*cosintegral(2*f*x+2*c*f*d^-1)*sin(-2*e+2*c*f*d^-1)+-1*f*a^2*d^-2*cos(-2*e+2*c*f*d^-1)*sinintegral(2*f*x+2*c*f*d^-1)+-2*f*a^2*d^-2*cos(-1*e+c*f*d^-1)*sinintegral(f*x+c*f*d^-1)+2*f*a^2*d^-2*cosintegral(f*x+c*f*d^-1)*sin(-1*e+c*f*d^-1))
@test integrate((a+a*cos(e+f*x))^-1*(c+d*x)^3, x) == :(a^-1*f^-1*(c+d*x)^3*tan((1/2)*e+(1/2)*f*x)+-1*im*a^-1*f^-1*(c+d*x)^3+12*a^-1*d^3*f^-4*Polylogarithms.polylog(3,-1*exp(im*(e+f*x)))+6*d*a^-1*f^-2*(c+d*x)^2*log(1+exp(im*(e+f*x)))+-12*im*a^-1*d^2*f^-3*(c+d*x)*Polylogarithms.polylog(2,-1*exp(im*(e+f*x))))
@test integrate((a+a*cos(e+f*x))^-1*(c+d*x)^2, x) == :(a^-1*f^-1*(c+d*x)^2*tan((1/2)*e+(1/2)*f*x)+-1*im*a^-1*f^-1*(c+d*x)^2+-4*im*a^-1*d^2*f^-3*Polylogarithms.polylog(2,-1*exp(im*(e+f*x)))+4*d*a^-1*f^-2*(c+d*x)*log(1+exp(im*(e+f*x))))
@test integrate((a+a*cos(e+f*x))^-1*(c+d*x), x) == :(a^-1*f^-1*(c+d*x)*tan((1/2)*e+(1/2)*f*x)+2*d*a^-1*f^-2*log(cos((1/2)*e+(1/2)*f*x)))
@test integrate((a+a*cos(e+f*x))^-2*(c+d*x)^3, x) == :(4*a^-2*d^3*f^-4*Polylogarithms.polylog(3,-1*exp(im*(e+f*x)))+4*a^-2*d^3*f^-4*log(cos((1/2)*e+(1/2)*f*x))+-1//3*im*a^-2*f^-1*(c+d*x)^3+1//3*a^-2*f^-1*(c+d*x)^3*tan((1/2)*e+(1/2)*f*x)+2*d*a^-2*f^-2*(c+d*x)^2*log(1+exp(im*(e+f*x)))+2*a^-2*d^2*f^-3*(c+d*x)*tan((1/2)*e+(1/2)*f*x)+-1//2*d*a^-2*f^-2*(c+d*x)^2*sec((1/2)*e+(1/2)*f*x)^2+1//6*a^-2*f^-1*(c+d*x)^3*sec((1/2)*e+(1/2)*f*x)^2*tan((1/2)*e+(1/2)*f*x)+-4*im*a^-2*d^2*f^-3*(c+d*x)*Polylogarithms.polylog(2,-1*exp(im*(e+f*x))))
@test integrate((a+a*cos(e+f*x))^-2*(c+d*x)^2, x) == :(-1//3*im*a^-2*f^-1*(c+d*x)^2+1//3*a^-2*f^-1*(c+d*x)^2*tan((1/2)*e+(1/2)*f*x)+2//3*a^-2*d^2*f^-3*tan((1/2)*e+(1/2)*f*x)+-4//3*im*a^-2*d^2*f^-3*Polylogarithms.polylog(2,-1*exp(im*(e+f*x)))+-1//3*d*a^-2*f^-2*sec((1/2)*e+(1/2)*f*x)^2*(c+d*x)+1//6*a^-2*f^-1*(c+d*x)^2*sec((1/2)*e+(1/2)*f*x)^2*tan((1/2)*e+(1/2)*f*x)+4//3*d*a^-2*f^-2*(c+d*x)*log(1+exp(im*(e+f*x))))
@test integrate((a+a*cos(e+f*x))^-2*(c+d*x), x) == :(-1//6*d*a^-2*f^-2*sec((1/2)*e+(1/2)*f*x)^2+1//3*a^-2*f^-1*(c+d*x)*tan((1/2)*e+(1/2)*f*x)+2//3*d*a^-2*f^-2*log(cos((1/2)*e+(1/2)*f*x))+1//6*a^-2*f^-1*sec((1/2)*e+(1/2)*f*x)^2*(c+d*x)*tan((1/2)*e+(1/2)*f*x))
@test integrate((a+-1*a*cos(e+f*x))^-1*(c+d*x)^3, x) == :(-1*im*a^-1*f^-1*(c+d*x)^3+-1*a^-1*f^-1*(c+d*x)^3*cot((1/2)*e+(1/2)*f*x)+12*a^-1*d^3*f^-4*Polylogarithms.polylog(3,exp(im*(e+f*x)))+6*d*a^-1*f^-2*(c+d*x)^2*log(1+-1*exp(im*(e+f*x)))+-12*im*a^-1*d^2*f^-3*(c+d*x)*Polylogarithms.polylog(2,exp(im*(e+f*x))))
@test integrate((a+-1*a*cos(e+f*x))^-1*(c+d*x)^2, x) == :(-1*im*a^-1*f^-1*(c+d*x)^2+-1*a^-1*f^-1*(c+d*x)^2*cot((1/2)*e+(1/2)*f*x)+-4*im*a^-1*d^2*f^-3*Polylogarithms.polylog(2,exp(im*(e+f*x)))+4*d*a^-1*f^-2*(c+d*x)*log(1+-1*exp(im*(e+f*x))))
@test integrate((a+-1*a*cos(e+f*x))^-1*(c+d*x), x) == :(-1*a^-1*f^-1*(c+d*x)*cot((1/2)*e+(1/2)*f*x)+2*d*a^-1*f^-2*log(sin((1/2)*e+(1/2)*f*x)))
@test integrate(x^3*(a+a*cos(c+d*x))^(1/2), x) == :(-96*d^-4*(a+a*cos(c+d*x))^(1/2)+12*d^-2*x^2*(a+a*cos(c+d*x))^(1/2)+-48*x*d^-3*(a+a*cos(c+d*x))^(1/2)*tan((1/2)*c+(1/2)*d*x)+2*d^-1*x^3*(a+a*cos(c+d*x))^(1/2)*tan((1/2)*c+(1/2)*d*x))
@test integrate(x^2*(a+a*cos(c+d*x))^(1/2), x) == :(-16*d^-3*(a+a*cos(c+d*x))^(1/2)*tan((1/2)*c+(1/2)*d*x)+8*x*d^-2*(a+a*cos(c+d*x))^(1/2)+2*d^-1*x^2*(a+a*cos(c+d*x))^(1/2)*tan((1/2)*c+(1/2)*d*x))
@test integrate(x*(a+a*cos(c+d*x))^(1/2), x) == :(4*d^-2*(a+a*cos(c+d*x))^(1/2)+2*x*d^-1*(a+a*cos(c+d*x))^(1/2)*tan((1/2)*c+(1/2)*d*x))
@test integrate((a+a*cos(c+d*x))^(1/2), x) == :(2*a*d^-1*(a+a*cos(c+d*x))^-1//2*sin(c+d*x))
@test integrate(x^-1*(a+a*cos(c+d*x))^(1/2), x) == :((a+a*cos(c+d*x))^(1/2)*cos((1/2)*c)*cosintegral((1/2)*d*x)*sec((1/2)*c+(1/2)*d*x)+-1*(a+a*cos(c+d*x))^(1/2)*sec((1/2)*c+(1/2)*d*x)*sin((1/2)*c)*sinintegral((1/2)*d*x))
@test integrate(x^-2*(a+a*cos(c+d*x))^(1/2), x) == :(-1*x^-1*(a+a*cos(c+d*x))^(1/2)+-1//2*d*(a+a*cos(c+d*x))^(1/2)*cos((1/2)*c)*sec((1/2)*c+(1/2)*d*x)*sinintegral((1/2)*d*x)+-1//2*d*(a+a*cos(c+d*x))^(1/2)*cosintegral((1/2)*d*x)*sec((1/2)*c+(1/2)*d*x)*sin((1/2)*c))
@test integrate(x^-3*(a+a*cos(c+d*x))^(1/2), x) == :(-1//2*x^-2*(a+a*cos(c+d*x))^(1/2)+1//4*d*x^-1*(a+a*cos(c+d*x))^(1/2)*tan((1/2)*c+(1/2)*d*x)+-1//8*d^2*(a+a*cos(c+d*x))^(1/2)*cos((1/2)*c)*cosintegral((1/2)*d*x)*sec((1/2)*c+(1/2)*d*x)+1//8*d^2*(a+a*cos(c+d*x))^(1/2)*sec((1/2)*c+(1/2)*d*x)*sin((1/2)*c)*sinintegral((1/2)*d*x))
@test integrate(x^3*(a+a*cos(x))^(1/2), x) == :(-96*(a+a*cos(x))^(1/2)+12*x^2*(a+a*cos(x))^(1/2)+-48*x*(a+a*cos(x))^(1/2)*tan((1/2)*x)+2*x^3*(a+a*cos(x))^(1/2)*tan((1/2)*x))
@test integrate(x^2*(a+a*cos(x))^(1/2), x) == :(-16*(a+a*cos(x))^(1/2)*tan((1/2)*x)+8*x*(a+a*cos(x))^(1/2)+2*x^2*(a+a*cos(x))^(1/2)*tan((1/2)*x))
@test integrate(x*(a+a*cos(x))^(1/2), x) == :(4*(a+a*cos(x))^(1/2)+2*x*(a+a*cos(x))^(1/2)*tan((1/2)*x))
@test integrate((a+a*cos(x))^(1/2), x) == :(2*a*(a+a*cos(x))^-1//2*sin(x))
@test integrate(x^-1*(a+a*cos(x))^(1/2), x) == :((a+a*cos(x))^(1/2)*cosintegral((1/2)*x)*sec((1/2)*x))
@test integrate(x^-2*(a+a*cos(x))^(1/2), x) == :(-1*x^-1*(a+a*cos(x))^(1/2)+-1//2*(a+a*cos(x))^(1/2)*sec((1/2)*x)*sinintegral((1/2)*x))
@test integrate(x^-3*(a+a*cos(x))^(1/2), x) == :(-1//2*x^-2*(a+a*cos(x))^(1/2)+-1//8*(a+a*cos(x))^(1/2)*cosintegral((1/2)*x)*sec((1/2)*x)+1//4*x^-1*(a+a*cos(x))^(1/2)*tan((1/2)*x))
@test integrate(x^3*(a+-1*a*cos(x))^(1/2), x) == :(-96*(a+-1*a*cos(x))^(1/2)+12*x^2*(a+-1*a*cos(x))^(1/2)+-2*x^3*(a+-1*a*cos(x))^(1/2)*cot((1/2)*x)+48*x*(a+-1*a*cos(x))^(1/2)*cot((1/2)*x))
@test integrate(x^2*(a+-1*a*cos(x))^(1/2), x) == :(8*x*(a+-1*a*cos(x))^(1/2)+16*(a+-1*a*cos(x))^(1/2)*cot((1/2)*x)+-2*x^2*(a+-1*a*cos(x))^(1/2)*cot((1/2)*x))
@test integrate(x*(a+-1*a*cos(x))^(1/2), x) == :(4*(a+-1*a*cos(x))^(1/2)+-2*x*(a+-1*a*cos(x))^(1/2)*cot((1/2)*x))
@test integrate((a+-1*a*cos(x))^(1/2), x) == :(-2*a*(a+-1*a*cos(x))^-1//2*sin(x))
@test integrate(x^-1*(a+-1*a*cos(x))^(1/2), x) == :((a+-1*a*cos(x))^(1/2)*csc((1/2)*x)*sinintegral((1/2)*x))
@test integrate(x^-2*(a+-1*a*cos(x))^(1/2), x) == :(-1*x^-1*(a+-1*a*cos(x))^(1/2)+(1/2)*(a+-1*a*cos(x))^(1/2)*cosintegral((1/2)*x)*csc((1/2)*x))
@test integrate(x^-3*(a+-1*a*cos(x))^(1/2), x) == :(-1//2*x^-2*(a+-1*a*cos(x))^(1/2)+-1//4*x^-1*(a+-1*a*cos(x))^(1/2)*cot((1/2)*x)+-1//8*(a+-1*a*cos(x))^(1/2)*csc((1/2)*x)*sinintegral((1/2)*x))
@test integrate(x^3*(a+a*cos(x))^3//2, x) == :(-1280//9*a*(a+a*cos(x))^(1/2)+16*a*x^2*(a+a*cos(x))^(1/2)+-64//27*a*(a+a*cos(x))^(1/2)*cos((1/2)*x)^2+-640//9*a*x*(a+a*cos(x))^(1/2)*tan((1/2)*x)+8//3*a*x^2*(a+a*cos(x))^(1/2)*cos((1/2)*x)^2+8//3*a*x^3*(a+a*cos(x))^(1/2)*tan((1/2)*x)+-32//9*a*x*(a+a*cos(x))^(1/2)*cos((1/2)*x)*sin((1/2)*x)+4//3*a*x^3*(a+a*cos(x))^(1/2)*cos((1/2)*x)*sin((1/2)*x))
@test integrate(x^2*(a+a*cos(x))^3//2, x) == :(-224//9*a*(a+a*cos(x))^(1/2)*tan((1/2)*x)+32//3*a*x*(a+a*cos(x))^(1/2)+8//3*a*x^2*(a+a*cos(x))^(1/2)*tan((1/2)*x)+16//9*a*x*(a+a*cos(x))^(1/2)*cos((1/2)*x)^2+32//27*a*(a+a*cos(x))^(1/2)*sin((1/2)*x)^2*tan((1/2)*x)+4//3*a*x^2*(a+a*cos(x))^(1/2)*cos((1/2)*x)*sin((1/2)*x))
@test integrate(x*(a+a*cos(x))^3//2, x) == :(16//3*a*(a+a*cos(x))^(1/2)+8//9*a*(a+a*cos(x))^(1/2)*cos((1/2)*x)^2+8//3*a*x*(a+a*cos(x))^(1/2)*tan((1/2)*x)+4//3*a*x*(a+a*cos(x))^(1/2)*cos((1/2)*x)*sin((1/2)*x))
@test integrate(x^-1*(a+a*cos(x))^3//2, x) == :((1/2)*a*(a+a*cos(x))^(1/2)*cosintegral(3//2*x)*sec((1/2)*x)+3//2*a*(a+a*cos(x))^(1/2)*cosintegral((1/2)*x)*sec((1/2)*x))
@test integrate(x^-2*(a+a*cos(x))^3//2, x) == :(-2*a*x^-1*(a+a*cos(x))^(1/2)*cos((1/2)*x)^2+-3//4*a*(a+a*cos(x))^(1/2)*sec((1/2)*x)*sinintegral((1/2)*x)+-3//4*a*(a+a*cos(x))^(1/2)*sec((1/2)*x)*sinintegral(3//2*x))
@test integrate(x^-3*(a+a*cos(x))^3//2, x) == :(-1*a*x^-2*(a+a*cos(x))^(1/2)*cos((1/2)*x)^2+-9//16*a*(a+a*cos(x))^(1/2)*cosintegral(3//2*x)*sec((1/2)*x)+-3//16*a*(a+a*cos(x))^(1/2)*cosintegral((1/2)*x)*sec((1/2)*x)+3//2*a*x^-1*(a+a*cos(x))^(1/2)*cos((1/2)*x)*sin((1/2)*x))
@test integrate(x^3*(a+a*cos(c+d*x))^-1//2, x) == :(-96*im*d^-4*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(4,-1*im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+-48*x*d^-3*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(3,-1*im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+48*x*d^-3*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(3,im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+96*im*d^-4*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(4,im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+-12*im*d^-2*x^2*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+-4*im*d^-1*x^3*(a+a*cos(c+d*x))^-1//2*arctan(exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+12*im*d^-2*x^2*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x))
@test integrate(x^2*(a+a*cos(c+d*x))^-1//2, x) == :(-16*d^-3*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(3,-1*im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+16*d^-3*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(3,im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+-8*im*x*d^-2*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+-4*im*d^-1*x^2*(a+a*cos(c+d*x))^-1//2*arctan(exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+8*im*x*d^-2*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x))
@test integrate(x*(a+a*cos(c+d*x))^-1//2, x) == :(-4*im*d^-2*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+4*im*d^-2*(a+a*cos(c+d*x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x)+-4*im*x*d^-1*(a+a*cos(c+d*x))^-1//2*arctan(exp((1/2)*im*(c+d*x)))*cos((1/2)*c+(1/2)*d*x))
@test integrate((a+a*cos(c+d*x))^-1//2, x) == :(2^(1/2)*a^-1//2*d^-1*arctanh((1/2)*2^(1/2)*a^(1/2)*(a+a*cos(c+d*x))^-1//2*sin(c+d*x)))
@test integrate(x^3*(a+-1*a*cos(x))^-1//2, x) == :(-96*im*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(4,-1*exp((1/2)*im*x))*sin((1/2)*x)+-48*x*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(3,-1*exp((1/2)*im*x))*sin((1/2)*x)+-4*x^3*(a+-1*a*cos(x))^-1//2*arctanh(exp((1/2)*im*x))*sin((1/2)*x)+48*x*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(3,exp((1/2)*im*x))*sin((1/2)*x)+96*im*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(4,exp((1/2)*im*x))*sin((1/2)*x)+-12*im*x^2*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(2,exp((1/2)*im*x))*sin((1/2)*x)+12*im*x^2*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*exp((1/2)*im*x))*sin((1/2)*x))
@test integrate(x^2*(a+-1*a*cos(x))^-1//2, x) == :(-16*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(3,-1*exp((1/2)*im*x))*sin((1/2)*x)+16*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(3,exp((1/2)*im*x))*sin((1/2)*x)+-4*x^2*(a+-1*a*cos(x))^-1//2*arctanh(exp((1/2)*im*x))*sin((1/2)*x)+-8*im*x*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(2,exp((1/2)*im*x))*sin((1/2)*x)+8*im*x*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*exp((1/2)*im*x))*sin((1/2)*x))
@test integrate(x*(a+-1*a*cos(x))^-1//2, x) == :(-4*im*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(2,exp((1/2)*im*x))*sin((1/2)*x)+-4*x*(a+-1*a*cos(x))^-1//2*arctanh(exp((1/2)*im*x))*sin((1/2)*x)+4*im*(a+-1*a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*exp((1/2)*im*x))*sin((1/2)*x))
@test integrate((a+-1*a*cos(x))^-1//2, x) == :(-1*2^(1/2)*a^-1//2*arctanh((1/2)*2^(1/2)*a^(1/2)*(a+-1*a*cos(x))^-1//2*sin(x)))
@test integrate(x^3*(a+a*cos(x))^-3//2, x) == :(-3*a^-1*x^2*(a+a*cos(x))^-1//2+(1/2)*a^-1*x^3*(a+a*cos(x))^-1//2*tan((1/2)*x)+-24*im*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*x))*cos((1/2)*x)+-24*im*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(4,-1*im*exp((1/2)*im*x))*cos((1/2)*x)+-12*x*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(3,-1*im*exp((1/2)*im*x))*cos((1/2)*x)+12*x*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(3,im*exp((1/2)*im*x))*cos((1/2)*x)+24*im*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*x))*cos((1/2)*x)+24*im*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(4,im*exp((1/2)*im*x))*cos((1/2)*x)+-1*im*a^-1*x^3*(a+a*cos(x))^-1//2*arctan(exp((1/2)*im*x))*cos((1/2)*x)+-24*im*x*a^-1*(a+a*cos(x))^-1//2*arctan(exp((1/2)*im*x))*cos((1/2)*x)+-3*im*a^-1*x^2*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*x))*cos((1/2)*x)+3*im*a^-1*x^2*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*x))*cos((1/2)*x))
@test integrate(x^2*(a+a*cos(x))^-3//2, x) == :(-2*x*a^-1*(a+a*cos(x))^-1//2+(1/2)*a^-1*x^2*(a+a*cos(x))^-1//2*tan((1/2)*x)+-4*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(3,-1*im*exp((1/2)*im*x))*cos((1/2)*x)+4*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(3,im*exp((1/2)*im*x))*cos((1/2)*x)+4*a^-1*(a+a*cos(x))^-1//2*arctanh(sin((1/2)*x))*cos((1/2)*x)+-1*im*a^-1*x^2*(a+a*cos(x))^-1//2*arctan(exp((1/2)*im*x))*cos((1/2)*x)+-2*im*x*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*x))*cos((1/2)*x)+2*im*x*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*x))*cos((1/2)*x))
@test integrate(x*(a+a*cos(x))^-3//2, x) == :(-1*a^-1*(a+a*cos(x))^-1//2+(1/2)*x*a^-1*(a+a*cos(x))^-1//2*tan((1/2)*x)+im*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,-1*im*exp((1/2)*im*x))*cos((1/2)*x)+-1*im*a^-1*(a+a*cos(x))^-1//2*Polylogarithms.polylog(2,im*exp((1/2)*im*x))*cos((1/2)*x)+-1*im*x*a^-1*(a+a*cos(x))^-1//2*arctan(exp((1/2)*im*x))*cos((1/2)*x))
@test integrate(x^3*(a+b*cos(x))^-1, x) == :(-6*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(4,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+6*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(4,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+-3*x^2*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(2,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+3*x^2*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(2,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+im*x^3*(a^2+-1*b^2)^-1//2*log(1+b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+-1*im*x^3*(a^2+-1*b^2)^-1//2*log(1+b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+-6*im*x*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(3,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*x))+6*im*x*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(3,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*x)))
@test integrate(x^2*(a+b*cos(c+d*x))^-1, x) == :(-2*im*d^-3*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(3,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-2*x*d^-2*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(2,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+2*im*d^-3*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(3,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+2*x*d^-2*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(2,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+im*d^-1*x^2*(a^2+-1*b^2)^-1//2*log(1+b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-1*im*d^-1*x^2*(a^2+-1*b^2)^-1//2*log(1+b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x))))
@test integrate(x*(a+b*cos(c+d*x))^-1, x) == :(d^-2*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(2,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-1*d^-2*(a^2+-1*b^2)^-1//2*Polylogarithms.polylog(2,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+im*x*d^-1*(a^2+-1*b^2)^-1//2*log(1+b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-1*im*x*d^-1*(a^2+-1*b^2)^-1//2*log(1+b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x))))
@test integrate((a+b*cos(c+d*x))^-2*(e+f*x), x) == :(-1*f*d^-2*(a^2+-1*b^2)^-1*log(a+b*cos(c+d*x))+a*f*d^-2*(a^2+-1*b^2)^-3//2*Polylogarithms.polylog(2,-1*b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-1*a*f*d^-2*(a^2+-1*b^2)^-3//2*Polylogarithms.polylog(2,-1*b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+im*a*d^-1*(a^2+-1*b^2)^-3//2*(e+f*x)*log(1+b*(a+(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-1*im*a*d^-1*(a^2+-1*b^2)^-3//2*(e+f*x)*log(1+b*(a+-1*(a^2+-1*b^2)^(1/2))^-1*exp(im*(c+d*x)))+-1*b*d^-1*(a+b*cos(c+d*x))^-1*(a^2+-1*b^2)^-1*(e+f*x)*sin(c+d*x))

max_stars_repo_path: test/4 Trig functions/4.2 Cosine/4.2.10 (c+d x)^m (a+b cos)^n.jl


text:
test_extend_from_vk(name, ex) = test_ex(FromVk(StructDefinition{false}(struct_by_name(name))), :(from_vk(T::Type{$(VulkanGen.remove_vk_prefix(name))}, x::$name) = $ex))

@testset "Overloads" begin
    @testset "`from_vk`" begin
        test_extend_from_vk(:VkLayerProperties, :(
            T(from_vk(String, x.layerName), from_vk(VersionNumber, x.specVersion), from_vk(VersionNumber, x.implementationVersion), from_vk(String, x.description))
        ))

        test_extend_from_vk(:VkQueueFamilyProperties, :(
            T(x.queueFlags, x.queueCount, x.timestampValidBits, from_vk(_Extent3D, x.minImageTransferGranularity))
        ))

        test_extend_from_vk(:VkPhysicalDeviceMemoryProperties, :(
            T(x.memoryTypeCount, from_vk.(MemoryType, x.memoryTypes), x.memoryHeapCount, from_vk.(MemoryHeap, x.memoryHeaps))
        ))

        test_extend_from_vk(:VkDisplayPlaneCapabilities2KHR, :(
            T(x.sType, x.pNext, from_vk(DisplayPlaneCapabilitiesKHR, x.capabilities))
        ))

        test_extend_from_vk(:VkDrmFormatModifierPropertiesListEXT, :(
            T(x.sType, x.pNext, unsafe_wrap(Vector{DrmFormatModifierPropertiesEXT}, x.pDrmFormatModifierProperties, x.drmFormatModifierCount; own = true))
        ))

        test_extend_from_vk(:VkPhysicalDeviceGroupProperties, :(
            T(x.sType, x.pNext, x.physicalDeviceCount, PhysicalDevice.(x.physicalDevices), from_vk(Bool, x.subsetAllocation))
        ))
    end

    @testset "Parent navigation" begin
        test_ex(Parent(HandleDefinition(handle_by_name(:VkDevice))), :(
            parent(device::Device) = device.physical_device
        ))

        test_ex(Parent(HandleDefinition(handle_by_name(:VkSurfaceKHR))), :(
            parent(surface::SurfaceKHR) = surface.instance
        ))
    end
end

ext: jl


text:
module password

export check, check_strict, is_candidate, is_candidate_strict

using StatsBase

function check(first::Int, last::Int)
    return count([is_candidate(pass) for pass in first:last])
end

function check_strict(first::Int, last::Int)
    return count([is_candidate(pass) && is_candidate_strict(pass) for pass in first:last])
end

function is_candidate(password::String)
    is_candidate(parse(password, Int))
end

function is_candidate(password::Int)
    pass = digits(password)
    # password must be six characters
    size = length(pass)
    # password must be ascending, ie sorted smallest to biggest
    ascending = issorted(pass, rev=true)
    # password must use a digit twice
    repeated = length(sort(unique(pass))) < size

    return size == 6 && ascending && repeated
end

function is_candidate_strict(password::String)
    is_candidate_strict(parse(password, Int))
end

function is_candidate_strict(password::Int)
    pass = digits(password)
    # password must use a digit exactly twice
    return 2 in values(countmap(pass))
end

end


text:
struct SUInt128{V} <: Unsigned
    function SUInt128{V}() where V
        !(typeof(V) === UInt128) && throw(ArgumentError("SUInt128 only supports static UInt128 storage, got $(typeof(V))"))
        return new{V}()
    end
end

struct SUInt64{V} <: Unsigned
    function SUInt64{V}() where V
        !(typeof(V) === UInt64) && throw(ArgumentError("SUInt64 only supports static UInt64 storage, got $(typeof(V))"))
        new{V}()
    end
end

const SUInt{V} = SUInt64{V}

struct SUInt32{V} <: Unsigned
    function SUInt32{V}() where V
        !(typeof(V) === UInt32) && throw(ArgumentError("SUInt32 only supports static UInt32 storage, got $(typeof(V))"))
        new{V}()
    end
end

struct SUInt16{V} <: Unsigned
    function SUInt16{V}() where V
        !(typeof(V) === UInt16) && throw(ArgumentError("SUInt16 only supports static UInt16 storage, got $(typeof(V))"))
        new{V}()
    end
end

struct SUInt8{V} <: Unsigned
    function SUInt8{V}() where V
        !(typeof(V) === UInt8) && throw(ArgumentError("SUInt8 only supports static UInt8 storage, got $(typeof(V))"))
        new{V}()
    end
end

const SUnsigned{V} = Union{SUInt8{V},SUInt16{V},SUInt32{V},SUInt64{V},SUInt128{V}}

SUnsigned(x::UInt8) = SUInt8(x)
SUnsigned(x::UInt16) = SUInt16(x)
SUnsigned(x::UInt32) = SUInt32(x)
SUnsigned(x::UInt64) = SUInt64(x)
SUnsigned(x::UInt128) = SUInt128(x)

Base.widen(::Type{<:SUInt8{X}}) where X = SUInt16
Base.widen(::Type{<:SUInt16{X}}) where X = SUInt32
Base.widen(::Type{<:SUInt32{X}}) where X= SUInt64
Base.widen(::Type{<:SUInt64{X}}) where X = SUInt128
Base.widen(::Type{<:SUInt128{X}}) where X = SUInt128

Base.unsigned(::Type{T}) where {T<:SUnsigned} = T
Base.unsigned(x::SUnsigned) = x

Base.Signed(x::SUnsigned) = signed(x)

function SUnsigned(val::Val{V}) where V
    if V isa UInt8
        SUInt8{V}()
    elseif V isa UInt16
        SUInt16{V}()
    elseif V isa UInt32
        SUInt32{V}()
    elseif V isa UInt128
        SUInt128{V}()
    else
        SUInt64{V}()
    end
end

const AbstractUInt128 = Union{SUInt128,UInt128}
const AbstractUInt64  = Union{SUInt64,UInt64}
const AbstractUInt32  = Union{SUInt32,UInt32}
const AbstractUInt16  = Union{SUInt16,UInt16}
const AbstractUInt8   = Union{SUInt8,UInt8}

SUI2BUI = Dict(SUInt128 => UInt128,
               SUInt16 => UInt16,
               SUInt32 => UInt32,
               SUInt64 => UInt64,
               SUInt8 => UInt8)

BUI2SUI = Dict(UInt128 => SUInt128,
               UInt16 => SUInt16,
               UInt32 => SUInt32,
               UInt64 => SUInt64,
               UInt8 => SUInt8)

for (ST,BT) in SUI2BUI
    defbasics(ST, BT)
    defbool(ST, BT)
end

SUInt128One, SUInt128Zero = defmath(SUInt128, UInt128)

SUInt64One, SUInt64Zero = defmath(SUInt64, UInt64)
SUInt64ZeroType = typeof(SUInt64Zero)
SUInt64OneType = typeof(SUInt64One)

SUInt32One, SUInt32Zero = defmath(SUInt32, UInt32)
SUInt32ZeroType = typeof(SUInt32Zero)
SUInt32OneType = typeof(SUInt32One)

SUInt16One, SUInt16Zero = defmath(SUInt16, UInt16)
SUInt8One, SUInt8Zero = defmath(SUInt8, UInt8)


text:
SRC_DIR = joinpath(Base.source_dir(), "..", "..", "src")
OUT_DIR = joinpath(Base.source_dir(), "..", "..", "output")
mkpath(OUT_DIR)

#include(joinpath(SRC_DIR, "BridgeSDEInference.jl"))
#using Main.BridgeSDEInference
include(joinpath(SRC_DIR, "BridgeSDEInference_for_tests.jl"))

using StaticArrays
using Distributions
using Random
# Let's generate the data
# -----------------------
using Bridge
DIR = "auxiliary"
include(joinpath(SRC_DIR, DIR, "data_simulation_fns.jl"))
include(joinpath(SRC_DIR, DIR, "utility_functions.jl"))
Random.seed!(4)
#θ₀ = [10.0, 28.0, 8.0/3.0, 3.0, 3.0, 3.0]
θˣ = [10.0, 28.0, 8.0/3.0, 3.0]
Pˣ = LorenzCV(θˣ...)

x0, dt, T = ℝ{3}(1.5, -1.5, 25.0), 1/5000, 4.0
tt = 0.0:dt:T
XX, _ = simulate_segment(ℝ{3}(0.0, 0.0, 0.0), x0, Pˣ, tt)

θ_init = [5.0, 15.0, 6.0, 8.0]
Pˣ = LorenzCV(θ_init...)

skip = 1000

Σdiagel = 10^0
Σ = SMatrix{2,2}(1.0I)*Σdiagel
L = @SMatrix[0.0 1.0 0.0;
             0.0 0.0 1.0]

obs_time, obs_vals = XX.tt[1:skip:end], [rand(Gaussian(L*x, Σ)) for x in XX.yy[1:skip:end]]

aux_flag = Val{(false,true,true)}()
P̃ = [LorenzCVAux(θ_init..., t₀, u, T, v, aux_flag, x0[3]) for (t₀, T, u, v)
     in zip(obs_time[1:end-1], obs_time[2:end], obs_vals[1:end-1], obs_vals[2:end])]

model_setup = DiffusionSetup(Pˣ, P̃, PartObs())
set_observations!(model_setup, [L for _ in P̃], [Σ for _ in P̃], obs_vals, obs_time) # uses default fpt
set_imputation_grid!(model_setup, 1/2000)
set_x0_prior!(model_setup,
              GsnStartingPt(x0, @SMatrix [20.0 0.0 0.0;
                                          0.0 20.0 0.0;
                                          0.0 0.0 400.0]),
              x0)
set_auxiliary!(model_setup; skip_for_save=10^0,
               adaptive_prop=Adaptation(x0, [0.7, 0.4, 0.2, 0.2, 0.2],
                                        [500, 500, 500, 500, 500], 1))
initialise!(eltype(x0), model_setup, Vern7(), false, NoChangePt(100))

mcmc_setup = MCMCSetup(
      Imputation(NoBlocking(), 0.96, Vern7()),
      ParamUpdate(ConjugateUpdt(), [1,2,3], θ_init, nothing,
                  MvNormal(fill(0.0, 3), diagm(0=>fill(1000.0, 3))),
                  UpdtAuxiliary(Vern7(), check_if_recompute_ODEs(P̃, [1,2,3]))),
      ParamUpdate(MetropolisHastingsUpdt(), 4, θ_init,
                  UniformRandomWalk(0.5, true), ImproperPosPrior(),
                  UpdtAuxiliary(Vern7(), check_if_recompute_ODEs(P̃, 4))))

schedule = MCMCSchedule(4*10^3, [[1,2]],
                        (save=1*10^3, verbose=10^2, warm_up=100,
                         readjust=(x->x%100==0), fuse=(x->false)))

Random.seed!(4)
out = mcmc(mcmc_setup, schedule, model_setup)
#out, elapsed = @timeit mcmc(mcmc_setup, schedule, model_setup)
#display(out.accpt_tracker)
include(joinpath(SRC_DIR, DIR, "plotting_fns.jl"))
plot_acceptance([out[2].updates[3]])
plot_acceptance([out[2].updates[1]], [500, 500, 500, 500, 500])
plot_chains(out[2]; truth=[10.0, 28.0, 8.0/3.0, 3.0],
            ylims=[nothing, (25,30), (2,5), (0,10)])
plot_paths(out[1], out[2], schedule; obs=(times=obs_time[2:end],
                     vals=[[v[1] for v in obs_vals[2:end]],
                           [v[2] for v in obs_vals[2:end]]], indices=[2,3]))


text:
function augment_states(m::HetDSGE{T}, TTT::Matrix{T}, TTT_jump::Matrix{T}, RRR::Matrix{T},
                        CCC::Vector{T}, C_eqn::Vector{T}) where {T<:AbstractFloat}
    endo     = m.endogenous_states #augment_model_states(m.endogenous_states, n_model_states(m))
    endo_new = m.endogenous_states_augmented
    exo      = m.exogenous_shocks

    n_endo = n_model_states(m)
    n_exo  = n_shocks_exogenous(m)

    @assert (n_endo, n_endo) == size(TTT)
    @assert (n_endo, n_exo) == size(RRR)
    @assert n_endo == length(CCC)

    # Initialize augmented matrices
    n_new_states = length(endo_new)
    n_new_eqs = n_new_states
    TTT_aug = zeros(n_endo + n_new_eqs, n_endo + n_new_states)
    TTT_aug[1:n_endo, 1:n_endo] = TTT
    RRR_aug = [RRR; zeros(n_new_eqs, n_exo)]
    CCC_aug = [CCC; zeros(n_new_eqs)]

    ### TTT modifications

    # Track Lags
    TTT_aug[endo_new[:i_t1], first(endo[:i′_t])] = 1.0
    TTT_aug[endo_new[:c_t1], endo_new[:c_t]] = 1.0 #1:get_setting(m, :n_backward_looking_states)] = C_eqn
    TTT_aug[endo_new[:c_t], 1:get_setting(m, :n_backward_looking_states)] = -C_eqn  #endo_new[:c_t]] = 1.0
   #= TTT_aug[endo_new[:i_t1], endo[:i_t]] = 1.0
    TTT_aug[endo_new[:w_t1], endo[:w_t]] = 1.0
    TTT_aug[endo_new[:π_t1], endo[:π_t]] = 1.0
    TTT_aug[endo_new[:L_t1], endo[:L_t]]  = 1.0=#

    ## We construct state for expected inflation using the fact
    ##
    ##   E_t[s_{t+1}] = CCC + TTT*s_{t} = CCC + TTT * (CCC + TTT*s_{t-1} + RRR)
    ##                = (CCC+TTT*CCC) + (TTT^2)s_{t-1} + TTT*RRR
    ##
    ## So to construct the state for E_t[p_{t+1}], we need to make all these
    ## objects, and just index out the row relevant to pi_t

   #= T2 = TTT^2
    TR = TTT*RRR
    CTC = CCC+TTT*CCC

    TTT_aug[endo_new[:Et_π_t],:] = [T2[endo[:π_t],:]; zeros(n_new_states)]

    RRR_aug[endo_new[:Et_π_t],:] = TR[endo[:π_t],:]
    CCC_aug[endo_new[:Et_π_t],:] = CTC[endo[:π_t],:]=#

    return TTT_aug, RRR_aug, CCC_aug
end


text:
module InterfacesTests

include("../test_interfaces.jl")

nparts = 4
prun(test_interfaces,sequential,nparts)

nparts = (2,2)
prun(test_interfaces,sequential,nparts)

end # module


text:
n = rand(40:100)
A = matrixdepot("gilbert", n)
@test issym(A)

B = matrixdepot("gilbert", n, 0.2)
C = matrixdepot("gilbert", Int, n, 0.5)
@test nnz(B) <= nnz(C)
println("'gilbert' passed test ...")


text:
@compat (::Type{Sphere})(x...) = HyperSphere(x...)

@compat (::Type{Circle})(x...) = HyperSphere(x...)

widths{N, T}(c::HyperSphere{N, T}) = Vec{N, T}(radius(c)*2)
radius(c::HyperSphere) = c.r
origin(c::HyperSphere) = c.center
minimum{N, T}(c::HyperSphere{N, T}) = Vec{N, T}(origin(c)) - Vec{N, T}(radius(c))
maximum{N, T}(c::HyperSphere{N, T}) = Vec{N, T}(origin(c)) + Vec{N, T}(radius(c))


text:
using ProjectiveDictionaryPairLearning
using MAT, Base.Test

###### Accelerating the normalization with testing
exampleMat = [1.0 2 3; 4 5 6; 7 8 9; 3 4 5]

normcol_equal(matin) = matin ./ repmat(sqrt(sum(matin .^ 2, 1) + eps()), size(matin, 1), 1)

tic()
exampleResult = normcol_equal(exampleMat)
normcol_equalTime = toq()

# tic()
# normalize!(exampleMat, 2 ,1)
# normalizeTime = toq()

# @test_approx_eq(exampleMat, exampleResult)
# @test normalizeTime < normcol_equalTime

# println("normcol_equalTime : $normcol_equalTime s")
# println("normalizeTime     : $normalizeTime s")

###### Accelerating the Class-specific reconstruction error calculation
# load TtData, DictMat, EncoderMat, DictSize
# prepared using
    # matwrite("class-specific.mat", @compat Dict(
    #     "TtData" => TtData,
    #     "DictMat" => DictMat,
    #     "EncoderMat" => EncoderMat,
    #     "DictSize" => DictSize
    # ))
matFileName = "class-specific.mat"
data = matread(joinpath(dirname(@__FILE__), matFileName))
TtData, DictMat, EncoderMat, DictSize = data["TtData"], data["DictMat"], data["EncoderMat"], data["DictSize"]

# prepare
ClassNum    = length(DictMat)
PredictCoef = EncoderMat * TtData

ErrorMATLAB = zeros(ClassNum, size(TtData, 2))
ErrorSubtract = zeros(ClassNum, size(TtData, 2))

tic()
# Class-specific reconstruction error calculation
for i = 1:ClassNum
    ErrorMATLAB[i, :]= sum((DictMat[i] * PredictCoef[(i - 1) * DictSize + 1 : i * DictSize, :]- TtData) .^ 2, 1)
end
ErrorMATLABTime = toq()

tic()
# Class-specific reconstruction error calculation
for i=1:ClassNum
    @inbounds tempMat = DictMat[i] * PredictCoef[(i-1)*DictSize+1:i*DictSize, :]
    tempMat -= TtData
    @inbounds ErrorSubtract[i,:] = sum(abs2(tempMat), 1)
end
ErrorSubtractTime = toq()

@test_approx_eq(ErrorMATLAB, ErrorSubtract)
@test ErrorSubtractTime < ErrorMATLABTime

println("Finished performance tests!")


text:
using MultidimensionalSpectroscopy, PyPlot, QuantumOptics, LinearAlgebra, FFTW, Colors,
        Printf, DelimitedFiles, Random, Combinatorics, Distributions

# make sure to set script directory as pwd()
cd(@__DIR__)

pygui(true)

calc_2d = true

cmp = create_colormap("bright");

# https://qudev.phys.ethz.ch/static/content/science/papers/Fink2008.pdf

ωr = 4
Na = 1               # number of cavity fock states
A  = FockBasis(Na)
a  = destroy(A)
at = create(A)

ωa = 4.1
Nσ = 2              # number of levels of matter
Σ  = NLevelBasis(Nσ)
sm  = transition(Σ, 1, 2)
sp = transition(Σ, 2, 1)
sz = sp * sm

H_r = (ωr * at * a)
H_a = (ωa * sz)

"""
generate H for Nm TLSs
H_am = H_a ⊗ one(H_a) ⊗ one(H_a) ⊗ ... + one(H_a) ⊗ H_a ⊗ one(H_a) ⊗ ... +
        one(H_a) ⊗ one(H_a) ⊗ H_a ⊗ ... + ...
"""
Nm = 2

# vector with molecule indices
R = collect(1:Nm)

for i in R
    if i == 1
        H_am = H_a
        smm = sm
    else
        global H_am = H_am ⊗ one(H_a)
        global smm = smm ⊗ one(H_a)
    end
end

spp = dagger(smm)

disorder = zeros(Nm)

for j in R
    println(circshift(R,j-1))
    if j == 1
        disorder[j] = 1
        global H_amb = H_am
    else
        disorder[j] = 1 + (rand(Float64, 1)[1] - 0.5) / 20
        global H_amb = H_amb + disorder[j] * permutesystems(H_am,circshift(R,j-1))
    end
end

for i in R
    if i == 1
        println(circshift(R,i-1))
        global H_x = at ⊗ one(H_amb) * (one(A) ⊗ smm) + (one(A) ⊗ spp) * (a ⊗ one(H_amb))

    else
        global H_x = H_x + at ⊗ one(H_amb) * (one(A) ⊗ permutesystems(smm,circshift(R,i-1))) +
                    one(A) ⊗ permutesystems(spp,circshift(R,i-1)) * (a ⊗ one(H_amb))
    end
end

#H_x = at ⊗ one(H_am) * (one(A) ⊗ sm1) + one(A) ⊗ sp1 * (a ⊗ one(H_am)) +
#        at ⊗ one(H_am) * (one(A) ⊗ sm2) + one(A) ⊗ sp2 * (a ⊗ one(H_am)) +
#          at ⊗ one(H_am) * (one(A) ⊗ sm3) + one(A) ⊗ sp3 * (a ⊗ one(H_am))

#H_x = at ⊗ one(Σ) * (one(A) ⊗ sm) + one(A) ⊗ sp * (a ⊗ one(Σ))

H_TC = H_r ⊗ one(H_amb) + (one(A) ⊗ H_amb)  + 0.4 * H_x

figure(figsize=(5,4));
imshow(real(H_TC.data-H_x.data + 20*H_x.data),origin="lower",cmap="Blues")
title("H_TC with off-diagonal elements * 20")
#pcolormesh(real(H_x.data))

figure(figsize=(5,4));
imshow(real(H_TC.data),cmap="Greys")

title("H_TC")
colorbar()

num_of_exc = 0
ground = tensor(fockstate(A,0),repeat([nlevelstate(Σ,1)],Nm-num_of_exc)...)

E = eigvals(dense(H_TC).data)
xn = 0
figure(figsize=(5,4));
title("Energy levels (Eigenvalues) of H_TC")
hlines([0, ωr], -.65, -.35)
hlines([0],  .35,  .65)
hlines([ωa .* disorder[:]],  .35,  .65)
for i in 1:length(E)
    println(i)
    if i > 1 && E[i] - E[i-1] < 0.1
        hlines(E[i],-.15,.15)
        global xn += .02
    else
        hlines(E[i],-.15,.15)
        global xn = 0
    end
end
xlim([-.8, .8])

for i in R
    if i == 1
        Psi0 = fockstate(A,0) ⊗ nlevelstate(Σ,1)
    else
        global Psi0 =  Psi0 ⊗ nlevelstate(Σ,1)
    end
end

# GS density matrix
rho0 = dm(Psi0)

# excite density matrix
rho1 = (at ⊗ one(smm)) * rho0 * dagger((at ⊗ one(smm)))

# debug: view rho0 after "excitation"
# imshow(real(rho0.data),origin="lower",cmap="Blues")

tlist = [0:0.1:100;]

# [relaxation within matter, relaxtion with light field]
J = [one(at) ⊗ (smm), a ⊗ one(smm)]
Γ = [.6, .06]

# master equation time evolution
tout, rhot = timeevolution.master(tlist,rho1,H_TC,J;rates=Γ)

# visualize evolution of density matrix
# cmds.view_dm_evo(rhot,5)

# calculate population in light and matter
nc = real(expect((at ⊗ one(smm)) * (a ⊗ one(smm)), rhot))
na = real(expect(one(at) ⊗ (spp * smm), rhot))

# plot population
figure(figsize=(10,4));
plot(tout,nc)
plot(tout,na)

# calculate correlation function
corr = timecorrelations.correlation(tlist, rho1, H_TC, J, (at ⊗ one(smm)),
            (a ⊗ one(smm)); rates=Γ)

# calculate spectrum from correlation
ω, spec = timecorrelations.correlation2spectrum(tlist, corr; normalize_spec=true)

# plot both
figure(figsize=(10,4));
ax1 = subplot(211);
plot(tlist,corr)
ax2 = subplot(212);
plot(ω, spec)

# parameters for 2D
F = 0.6 * J
H = H_TC

# interaction light with cavity and TLSs
#μ12 = at ⊗ (sm*sp) + a ⊗ (sm*sp) + (a*at) ⊗ sp + (a*at) ⊗ sm
# interact light only with cavity
μ12 = at ⊗ (smm*spp) + a ⊗ (smm*spp)
# interact with both (?) for higher states
#μ23 = at ⊗ (spp*smm) + a ⊗ (spp*smm) + (at*a) ⊗ spp + (at*a) ⊗ smm
# or only the atom ?
#μ23 = (at*a) ⊗ spp + (at*a) ⊗ smm
# or only the atom ?
μ23 = at ⊗ (spp*smm) + a ⊗ (spp*smm)

# GS density matrix (reset it you idiot)
rho0 = dm(Psi0)

if calc_2d
    ## calculate (complex) 3rd order corr function (with T=0)
    zp = 11 # zeropad up to 2^zp

    ## calculate 2D spectra at
    T = [0] #fs

    spectra2d = Array{out2d}(undef, length(T))
    for i = 1:length(T)
        spectra2d[i] = make2Dspectra(tlist,[rho0, rho0],[H, H],[F, F],[μ12, μ12],[μ23, μ23],T[i],
                                            "lindblad";debug=true,zp=zp)
    end

    ## crop 2D data and increase dw
    spectra2d = [crop2d(spectra2d[i],1;w_max=20,step=2) for i = 1:length(T)]

    ## assign ω-axis from output
    ω = spectra2d[1].ω
end

if calc_2d
    ## plot 2D spectra for each(?) T
    # what to plot
    rep="absorptive"
    scal="lin"

    ## make  subplot layout
    nplots = length(T); ncols = ceil(sqrt(nplots)); nrows = ceil(nplots / ncols)

    ## create figure
    figure(figsize=(ncols*3,nrows*3));
    suptitle(rep * " 2D spectrum (" * scal * ". scaling)")
    for i = 1:nplots
        subplot(nrows,ncols,i,aspect="equal")
        plot2d(ω,round.(spectra2d[i].full2d,digits=1);repr=rep,scaling=scal)
        title("2D spectrum at $(T[i]) fs")
        #xlim([0, E₁+E₂]); ylim([0, E₁+E₂]);
    end
end

# The end
#########################################################################

#H_TC = H_r ⊗ one(Σ) ⊗ one(Σ) + one(A) ⊗ H_a ⊗ one(Σ)

"""
b = NLevelBasis(2)
H = diagonaloperator(b, [0, 1])
j12 = transition(b, 1, 2) # decay from 2nd into 1st level

Nb = 3
b_mb = ManyBodyBasis(b, bosonstates(b, [1, 1, 1]))
H_mb = manybodyoperator(b_mb, H)
j12_mb = manybodyoperator(b_mb, j12)

wc = 1
N = 1               # number of cavity fock states
b_fock = FockBasis(N)
a = destroy(b_fock) ⊗ one(b_mb)
at = dagger(a)

wa = 1
sm = one(b_fock) ⊗ j12_mb
sp = dagger(sm)
println(dense(sm))

Psi0 = fockstate(b_fock,0) ⊗ basisstate(b_mb,[Nb,0])    # start with an excited cavity
rho0 = dm(Psi0)

rho0 = rho0 * a
g = 0.24

### this should be on the cloud

H_full = wc * at * a + wa * sp * sm + g * (at * sm + a * sp)

Γ = [0.01]
J_mb = [sm]

tlist = [0:0.1:50;]

tout, rhot = timeevolution.master(tlist,rho0,H_full,J_mb;rates=Γ)

nc = real(expect(at * a, rhot))
na = real(expect(sp * sm, rhot))

figure(figsize=(6,3))
plot(tout,nc)
"""

max_stars_repo_name: lbmnky/CMDS

max_issues_repo_name: lbmnky/CMDS

max_forks_repo_name: lbmnky/CMDS


text:
@testset "Input and output" begin

@test begin
    x = Dict{Symbol, Vector}()
    x[:A] = [1,2,3]
    x[:B] = [2.0,4.0,6.0]
    table = readtable(Val{(:A,:B)}, x)
    table == @Table(A=[1,2,3], B=[2.0,4.0,6.0])
end

@test begin
    fileio = open("data1.csv")
    table = readtable(Table{(:A,:B,:C,:D),Tuple{Vector{Int},Vector{Float64},Vector{Compat.ASCIIString},Vector{Bool}}}, fileio, header = true)
    close(fileio)
    table == @Table(A=[1,2,3], B=[2.0,4.0,6.0], C=["A","AB","ABC"], D=[true,false,true])
end

@test begin
    stream = IOBuffer()
    t = @Table(A=[1,2,3], B=[2.0,4.0,6.0])
    writetable(stream, t)
    @compat String(stream) == "1,2.0\n2,4.0\n3,6.0\n"
end

end


text:
@testset "Krylov Trust Region" begin

@testset "Toy test problem 1" begin
    # Test on actual optimization problems.
    function f(x::Vector)
        (x[1] - 5.0)^4
    end

    function fg!(g, x)
        g[1] = 4.0 * (x[1] - 5.0)^3
        f(x)
    end

    function fg2!(_f, g, x)
        g  === nothing || (g[1] = 4.0 * (x[1] - 5.0)^3)
        _f === nothing || return f(x)
    end

    function hv!(Hv, x, v)
        Hv[1] = 12.0 * (x[1] - 5.0)^2 * v[1]
    end

    d = Optim.TwiceDifferentiableHV(f, fg!, hv!, [0.0])
    d2 = Optim.TwiceDifferentiableHV(NLSolversBase.only_fg_and_hv!(fg2!, hv!), [0.0])

    result = Optim.optimize(d, [0.0], Optim.KrylovTrustRegion())
    @test norm(Optim.minimizer(result) - [5.0]) < 0.01
    result = Optim.optimize(d2, [0.0], Optim.KrylovTrustRegion())
    @test norm(Optim.minimizer(result) - [5.0]) < 0.01
end

@testset "Toy test problem 2" begin
    eta = 0.9

    function f2(x::Vector)
        0.5 * (x[1]^2 + eta * x[2]^2)
    end

    function fg2!(storage::Vector, x::Vector)
        storage[:] = [x[1], eta * x[2]]
        f2(x)
    end

    function hv2!(Hv::Vector, x::Vector, v::Vector)
        Hv[:] = [1.0 0.0; 0.0 eta] * v
    end

    d2 = Optim.TwiceDifferentiableHV(f2, fg2!, hv2!, Float64[127, 921])

    result = Optim.optimize(d2, Float64[127, 921], Optim.KrylovTrustRegion())
    @test result.g_converged
    @test norm(Optim.minimizer(result) - [0.0, 0.0]) < 0.01
end

@testset "Stock test problems" begin
    for (name, prob) in MultivariateProblems.UnconstrainedProblems.examples
      if prob.istwicedifferentiable
            hv!(storage::Vector, x::Vector, v::Vector) = begin
                n = length(x)
                H = Matrix{Float64}(undef, n, n)
                MVP.hessian(prob)(H, x)
                storage .= H * v
            end
            fg!(g::Vector, x::Vector) = begin
                MVP.gradient(prob)(g,x)
                MVP.objective(prob)(x)
            end
        ddf = Optim.TwiceDifferentiableHV(MVP.objective(prob), fg!, hv!, prob.initial_x)
        result = Optim.optimize(ddf, prob.initial_x, Optim.KrylovTrustRegion())
        @test norm(Optim.minimizer(result) - prob.solutions) < 1e-2
      end
    end
end

end


text:
module ParticleGridMapper
using StaticArrays
using OctreeBH
#using LinearAlgebra  #norm()
using .Threads

export map_particle_to_2Dgrid_loopP,
map_particle_to_2Dgrid_loopP_thread,
map_particle_to_3Dgrid_loopP_thread,
map_particle_to_2Dgrid_loopP_noCar,
map_particle_to_3Dgrid_NGP,
map_particle_to_3Dgrid_NGP_thread,
buildtree, DataP2G, get_AMRgrid_volumes, get_AMRgrid, get_AMRfield,
map_particle_to_AMRgrid_SPH!,
map_particle_to_AMRgrid_MFM!,
map_particle_to_AMRgrid_NGP!,
map_particle_to_AMRgrid_loopP_recursive!,
project_AMRgrid_to_image,
get_max_tree_depth,
set_max_depth_AMR!,
balance_all_level!

import OctreeBH: AbstractData, assign_additional_node_data!

include("kernel.jl")
include("cartesian_grid.jl")
include("mapAMRgrid.jl")

end


text:
#Ost - "E6" - восток, Nord - "N8" - север, West - "W4" - запад, Sud - "S2" - юг
function main_function9!(r)
    steps = 1
    side = Nord
    while ismarker(r)==false
        for _ in 1:2
            special_move!(r,side,steps)
            side=counterclockwise(side)
        end
        steps += 1
    end
end
function special_move!(r, side, num::Int)
    for _ in 1:num
        if ismarker(r)
            return nothing
        end
        move!(r,side)
    end
end
function counterclockwise(side)
    return HorizonSide(mod(Int(side)+1,4))
end


text:
### A Pluto.jl notebook ###
# v0.17.5

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el)
        el
    end
end

# ╔═╡ 08980845-f030-4a64-a4b2-ab027b3a2721
   begin  using Pkg
       Pkg.activate(".testenv")
	   Pkg.add("Revise")
       using Revise
       Pkg.add(["Triangulate","ExtendableGrids"])
	   Pkg.add("GridVisualize")
	   Pkg.add("PlutoUI")
	   Pkg.develop("PlutoVista")
   end

# ╔═╡ d6c0fb79-4129-444a-978a-bd2222b53df6
begin
    using PlutoUI
    using PlutoVista
    using GridVisualize
    using Printf
    using Triangulate
	using ExtendableGrids
end

# ╔═╡ 93ca4fd0-8f61-4174-b459-55f5395c0f56
md"""
# Test Notebook for [PlutoVista](https://github.com/j-fu/PlutoVista.jl) and [vtk.js](https://kitware.github.io/vtk-js/index.html)
"""

# ╔═╡ 7c06fcf0-8c98-49f7-add8-435f57a9c9da
function maketriangulation(maxarea)
    triin=Triangulate.TriangulateIO()
    triin.pointlist=Matrix{Cdouble}([-1.0 -1.0 ; 1.0 -1.0 ; 1.0 2 ; -1.0 1.0]')
    triin.segmentlist=Matrix{Cint}([1 2 ; 2 3 ; 3 4 ; 4 1; 4 2 ]')
    triin.segmentmarkerlist=Vector{Int32}([1, 2, 3, 4, 5])
    triin.regionlist=Matrix{Cdouble}([-0.9 0.9; -0.9 0.9; 1 2 ; 2.0*maxarea maxarea])
    (triout, vorout)=triangulate("paADQ", triin)
    triout.pointlist, triout.trianglelist,Int.(vec(triout.triangleattributelist)),triout.segmentlist,triout.segmentmarkerlist
end

# ╔═╡ 60dcfcf5-391e-418f-8e7c-3a0fe94f1e0d
p=PlutoVTKPlot(resolution=(300,300),tickfontsize=10)

# ╔═╡ db2823d9-aa6d-4be3-af5c-873c072cfd2b
md"""
Change grid resolution: $(@bind resolution Slider(5:200))
"""

# ╔═╡ 890710fe-dac0-4256-b1ba-79776f1ea7e5
(pts,tris,markers,edges,edgemarkers)=maketriangulation(1/resolution^2)

# ╔═╡ b8a976e3-7fef-4527-ae6a-4da31c93a04f
func=0.5.+0.5*[sin(10*pts[1,i])*cos(10*pts[2,i]) for i=1:size(pts,2)]

# ╔═╡ 401b36bd-fa8f-4a9c-9556-bbc82c3ddbca
 md"""
Change time: $(@bind time Slider(0:0.1:10,show_value=true))
"""

# ╔═╡ 6fd4a1ee-7a4a-405b-8e1f-5819eababe10
ft=0.5*[sin(10*pts[1,i]-time)*cos(10*pts[2,i]-time) for i=1:size(pts,2)];

# ╔═╡ e76f8a6a-ab91-454a-b200-cfc8b57eb331
triplot!(p,pts,tris,ft)

# ╔═╡ bce0cfe7-4112-4bb8-aac6-43885f3746a9
md"""Number of gridpoints: $(size(pts,2)) """

# ╔═╡ 81046dcd-3cfb-4133-943f-61b9b3cdb183
let
	p=PlutoVTKPlot(resolution=(300,300),axisfontsize=20)
	tricontour!(p,pts,tris,func;cmap=:spring,levels=(0.1:0.2:1),limits=(0,1),xlabel="a",aspect=1.5)
end

# ╔═╡ 7019ce3f-f2db-4581-8bd9-64f76231a62a
let
	p=PlutoVTKPlot(resolution=(300,300),legendfontsize=15)
	trimesh!(p,pts,tris;
		markers=markers,
		edges=edges,edgemarkers=edgemarkers
	)
end

# ╔═╡ 2e3546f6-eb47-4693-aa00-902570fab7b5
function grid3d(;n=15)
    X=collect(0:1/n:1)
    g=simplexgrid(X,X,X)
end

# ╔═╡ bd0a59a2-564d-42bd-ab6f-a50b26f9241f
function func3d(;n=15)
    g=grid3d(n=n)
    g, map((x,y,z)->sinpi(2*x)*sinpi(3.5*y)*sinpi(1.5*z),g)
end

# ╔═╡ 368b8cf5-fabd-4b84-b33c-b15c4452393b
	g,f=func3d(;n=20)

# ╔═╡ f64729e4-d2b4-40d3-acbb-1395dbe0337d
p3d=PlutoVTKPlot(resolution=(300,300))

# ╔═╡ 7be35f33-5f7a-4765-8bb6-1487e209efc8
md"""
f: $(@bind flevel Slider(0:0.01:1,show_value=true,default=0.45))

x: $(@bind xplane Slider(0:0.01:1,show_value=true,default=0.45))
y: $(@bind yplane Slider(0:0.01:1,show_value=true,default=0.45))
z: $(@bind zplane Slider(0:0.01:1,show_value=true,default=0.45))
"""

# ╔═╡ 3681ef5b-c794-44da-9fe7-cedcd68b426c
tetcontour!(p3d,g[Coordinates],g[CellNodes],f;levels=[flevel],
	faces=g[BFaceNodes],
	facemarkers=g[BFaceRegions],
	xplanes=[xplane],yplanes=[yplane],zplanes=[zplane],outlinealpha=0.1, levelalpha=0.25)

# ╔═╡ 0f440c27-7ff1-4db5-b4eb-8ce1e9018ef1
g[BFaceNodes]

# ╔═╡ 606f6837-f3b7-4a52-b9c3-034799c7bf93
g3x=let
	n=10
X=collect(0:1/n:1)
simplexgrid(X,X,X)
end

# ╔═╡ ecb3bb5e-6ae5-4d6e-9834-d52ce977b3fc
p3dx=PlutoVTKPlot(resolution=(300,300))

# ╔═╡ 90ff6ffc-84dc-45fd-8d09-9eb916397630
 md"""
x: $(@bind gxplane Slider(0:0.01:1,show_value=true,default=0.45))
y: $(@bind gyplane Slider(0:0.01:1,show_value=true,default=0.45))
z: $(@bind gzplane Slider(0:0.01:1,show_value=true,default=0.45))
"""

# ╔═╡ ae5707d9-41b7-4937-924a-fb54b83c31db
alpha=0.15

# ╔═╡ 519d106f-3f6f-4db1-b4e4-c5e7ef176857
tetmesh!(p3dx,g3x[Coordinates],g3x[CellNodes];
		markers=g3x[CellRegions],
	faces=g3x[BFaceNodes],
	facemarkers=g3x[BFaceRegions],
	xplanes=gxplane,yplanes=gyplane,zplanes=gzplane,outlinealpha=0)

# ╔═╡ 4263e897-d878-4fab-acae-a6c4dae37c5e
qp=PlutoVTKPlot(resolution=(500,500))

# ╔═╡ 5dfb1dde-06ed-483a-bac3-6a21a7f98856
@bind x0 Slider(1:0.01:9)

# ╔═╡ 6d7547c6-5667-463f-99dc-2f3fbaeb4c4d
let
	n=50
	X=0:10/n:10
	g=simplexgrid(X,X)
	pts=g[Coordinates]
    f(x,y)=sin(x-x0)*cos(y)+0.05*(x-x0)*y
    fx(x,y)=-cos(x-x0)*cos(y)-0.05*y
    fy(x,y)=sin(x-x0)*sin(y)-0.05*(x-x0)
    v=map(f,g)
    ∇v=0.5*vcat(map(fx,g)',map(fy,g)')

	tricontour!(qp,pts,g[CellNodes],v,colormap=:summer,levels=5)
	quiver2d!(qp,pts,∇v)
end

# ╔═╡ 23c75823-e69f-4ca0-a3ca-783612ba9c3c
html"""<hr>"""

# ╔═╡ Cell order:
# ╟─93ca4fd0-8f61-4174-b459-55f5395c0f56
# ╠═d6c0fb79-4129-444a-978a-bd2222b53df6
# ╠═7c06fcf0-8c98-49f7-add8-435f57a9c9da
# ╠═890710fe-dac0-4256-b1ba-79776f1ea7e5
# ╠═b8a976e3-7fef-4527-ae6a-4da31c93a04f
# ╠═60dcfcf5-391e-418f-8e7c-3a0fe94f1e0d
# ╟─db2823d9-aa6d-4be3-af5c-873c072cfd2b
# ╟─401b36bd-fa8f-4a9c-9556-bbc82c3ddbca
# ╠═6fd4a1ee-7a4a-405b-8e1f-5819eababe10
# ╠═e76f8a6a-ab91-454a-b200-cfc8b57eb331
# ╟─bce0cfe7-4112-4bb8-aac6-43885f3746a9
# ╠═81046dcd-3cfb-4133-943f-61b9b3cdb183
# ╠═7019ce3f-f2db-4581-8bd9-64f76231a62a
# ╠═2e3546f6-eb47-4693-aa00-902570fab7b5
# ╠═bd0a59a2-564d-42bd-ab6f-a50b26f9241f
# ╠═368b8cf5-fabd-4b84-b33c-b15c4452393b
# ╠═f64729e4-d2b4-40d3-acbb-1395dbe0337d
# ╟─7be35f33-5f7a-4765-8bb6-1487e209efc8
# ╠═3681ef5b-c794-44da-9fe7-cedcd68b426c
# ╠═0f440c27-7ff1-4db5-b4eb-8ce1e9018ef1
# ╠═606f6837-f3b7-4a52-b9c3-034799c7bf93
# ╠═ecb3bb5e-6ae5-4d6e-9834-d52ce977b3fc
# ╠═90ff6ffc-84dc-45fd-8d09-9eb916397630
# ╠═ae5707d9-41b7-4937-924a-fb54b83c31db
# ╠═519d106f-3f6f-4db1-b4e4-c5e7ef176857
# ╠═4263e897-d878-4fab-acae-a6c4dae37c5e
# ╠═5dfb1dde-06ed-483a-bac3-6a21a7f98856
# ╠═6d7547c6-5667-463f-99dc-2f3fbaeb4c4d
# ╟─23c75823-e69f-4ca0-a3ca-783612ba9c3c
# ╠═08980845-f030-4a64-a4b2-ab027b3a2721


text:
module MatrixBlobs

using Blobs
import Blobs: @logmsg, load, save, flush
using Base.Random: UUID
import Base: serialize, deserialize, getindex, setindex!, size, append!, flush, *

export DenseMatBlobs, SparseMatBlobs, size, getindex, setindex!, serialize, deserialize, save, load, flush, *

if isless(Base.VERSION, v"0.5.0-")
    view = sub
end

const BYTES_128MB = 128 * 1024 * 1024

def_cache() = max(BYTES_128MB, floor(Int, Base.Sys.free_memory()/2/nworkers()))

def_sz{T}(::Type{T}, blk_sz::Int=BYTES_128MB) = floor(Int, blk_sz / sizeof(T))

relidx(range::UnitRange{Int}, i1::Int, i2::Int) = i1, (i2 - first(range) + 1)

type SparseMatBlobs{Tv,Ti} <: AbstractMatrix{Tv}
    metadir::AbstractString
    sz::Tuple{Int,Int}
    splits::Vector{Pair{UnitRange{Int},Base.Random.UUID}}
    coll::BlobCollection{SparseMatrixCSC{Tv,Ti}}
end

type DenseMatBlobs{T} <: AbstractMatrix{T}
    metadir::AbstractString
    sz::Tuple{Int,Int}
    splits::Vector{Pair{UnitRange{Int},Base.Random.UUID}}        # keep a mapping of index ranges to blobs
    coll::BlobCollection{Matrix{T}}
end

typealias MatBlobs Union{SparseMatBlobs,DenseMatBlobs}

size{T<:MatBlobs}(dm::T) = dm.sz

split_ranges{T<:MatBlobs}(dm::T) = [p.first for p in dm.splits]

function splitidx{T<:MatBlobs}(dm::T, splitdim_idx::Int)
    splits = UnitRange{Int}[x.first for x in dm.splits]
    for splitnum in 1:length(splits)
        if splitdim_idx in splits[splitnum]
            return splitnum
        end
    end
    throw(BoundsError("MatBlobs $(dm.sz)", splitdim_idx))
end

function serialize(s::SerializationState, sm::MatBlobs)
    Serializer.serialize_type(s, typeof(sm))
    serialize(s, sm.metadir)
    serialize(s, sm.sz)
    serialize(s, sm.splits)

    coll = sm.coll
    serialize(s, coll.id)
    serialize(s, coll.mutability)
    serialize(s, coll.reader)
    serialize(s, coll.blobs)
    serialize(s, coll.maxcache)
end

function matblob(metadir::AbstractString; maxcache::Int=def_cache())
    @logmsg("reading back matrix from $metadir")
    open(joinpath(metadir, "meta"), "r") do io
        mat = deserialize(SerializationState(io))
        mat.metadir = metadir
        max_cached!(mat.coll, maxcache)
        mat
    end
end

function flush(dm::MatBlobs, wrkrs::Vector{Int}=Int[]; callback::Bool=true)
    isempty(wrkrs) ? flush(dm.coll; callback=callback) : flush(dm.coll, wrkrs; callback=callback)
    nothing
end

##
# SparseMatBlobs specific functions
sersz{Tv,Ti}(sp::SparseMatrixCSC{Tv,Ti}) = (sizeof(Int64)*3 + sizeof(sp.colptr) + sizeof(sp.rowval) + sizeof(sp.nzval))

function load{Tv<:Real,Ti<:Integer}(meta::FileMeta, reader::FileBlobIO{SparseMatrixCSC{Tv,Ti}})
    open(meta.filename, "r+") do fhandle
        seek(fhandle, meta.offset)
        header = Array(Int64, 3)
        pos1 = position(fhandle)
        header = read!(fhandle, header)
        m = header[1]
        n = header[2]
        nz = header[3]

        pos1 += sizeof(header)
        colptr = reader.use_mmap ? blobmmap(fhandle, Vector{Ti}, (n+1,), pos1) : read!(fhandle, Array(Ti, n+1))

        pos1 += sizeof(colptr)
        rowval = reader.use_mmap ? blobmmap(fhandle, Vector{Ti}, (nz,), pos1) : read!(fhandle, Array(Ti, nz))

        pos1 += sizeof(rowval)
        nzval = reader.use_mmap ? blobmmap(fhandle, Vector{Tv}, (nz,), pos1) : read!(fhandle, Array(Tv, nz))
        return SparseMatrixCSC{Tv,Ti}(m, n, colptr, rowval, nzval)
    end
end

function save{Tv<:Real,Ti<:Integer}(spm::SparseMatrixCSC{Tv,Ti}, meta::FileMeta, writer::FileBlobIO{SparseMatrixCSC{Tv,Ti}})
    if writer.use_mmap && ismmapped(spm.colptr) && ismmapped(spm.rowval) && ismmapped(spm.nzval)
        syncmmapped(spm.colptr)
        syncmmapped(spm.rowval)
        syncmmapped(spm.nzval)
    else
        header = Int64[spm.m, spm.n, length(spm.nzval)]

        touch(meta.filename)
        open(meta.filename, "r+") do fhandle
            seek(fhandle, meta.offset)
            write(fhandle, header)
            write(fhandle, spm.colptr)
            write(fhandle, spm.rowval)
            write(fhandle, spm.nzval)
        end
    end
    nothing
end

function load(sm::SparseMatBlobs, col::Int)
    splitnum = splitidx(sm, col)
    p = sm.splits[splitnum]
    range = p.first
    bid = p.second
    data = load(sm.coll, bid)
    #@logmsg("loaded split idx:$(col) from splitnum $splitnum ($bid) with range: $range, sz: $(size(data))")
    data, range
end

getindex{Tv}(sm::SparseMatBlobs{Tv}, i::Int) = getindex(sm, ind2sub(size(sm), i)...)
function getindex{Tv}(sm::SparseMatBlobs{Tv}, i1::Int, i2::Int)
    part, range = load(sm, i2)
    part[relidx(range, i1, i2)...]
end

function getindex{Tv}(sm::SparseMatBlobs{Tv}, ::Colon, i2::Int)
    part, range = load(sm, i2)
    reli1, reli2 = relidx(range, 1, i2)
    part[:,reli2]
end

function deserialize{Tv,Ti}(s::SerializationState, ::Type{SparseMatBlobs{Tv,Ti}})
    metadir = deserialize(s)
    sz = deserialize(s)
    splits = deserialize(s)

    coll_id = deserialize(s)
    coll_mut = deserialize(s)
    coll_reader = deserialize(s)
    coll_blobs = deserialize(s)
    coll_maxcache = deserialize(s)

    coll = BlobCollection(SparseMatrixCSC{Tv,Ti}, coll_mut, coll_reader; maxcache=coll_maxcache, strategy=maxmem, id=coll_id)
    coll.blobs = coll_blobs
    SparseMatBlobs{Tv,Ti}(metadir, sz, splits, coll)
end

function SparseMatBlobs{Tv,Ti}(::Type{Tv}, ::Type{Ti}, metadir::AbstractString; maxcache::Int=def_cache())
    T = SparseMatrixCSC{Tv,Ti}
    mut = Mutable(BYTES_128MB, FileBlobIO(T, true))
    coll = BlobCollection(T, mut, FileBlobIO(T, true); maxcache=maxcache, strategy=maxmem)
    SparseMatBlobs{Tv,Ti}(metadir, (0,0), Pair[], coll)
end

function append!{Tv,Ti}(sp::SparseMatBlobs{Tv,Ti}, S::SparseMatrixCSC{Tv,Ti})
    m,n = size(S)
    if isempty(sp.splits)
        sp.sz = (m, n)
        idxrange = 1:n
    else
        (sp.sz[1] == m) || throw(BoundsError("SparseMatBlobs $(sp.sz)", (m,n)))
        old_n = sp.sz[2]
        idxrange = (old_n+1):(old_n+n)
        sp.sz = (m, old_n+n)
    end

    fname = joinpath(sp.metadir, string(length(sp.splits)+1))
    meta = FileMeta(fname, 0, sersz(S))

    blob = append!(sp.coll, SparseMatrixCSC{Tv,Ti}, meta, StrongLocality(myid()), Nullable(S))
    push!(sp.splits, idxrange => blob.id)
    @logmsg("appending blob $(blob.id) of size: $(size(S)) for idxrange: $idxrange, sersz: $(meta.size)")
    blob
end

function save(sp::SparseMatBlobs, wrkrs::Vector{Int}=Int[])
    isempty(wrkrs) ? save(sp.coll) : save(sp.coll, wrkrs)
    save(sp.coll)
    open(joinpath(sp.metadir, "meta"), "w") do io
        serialize(SerializationState(io), sp)
    end
    nothing
end

SparseMatBlobs(metadir::AbstractString; maxcache::Int=def_cache()) = matblob(metadir; maxcache=maxcache)

##
# DenseMatBlobs specific functions
sersz{T}(d::Matrix{T}) = (sizeof(Int64)*2 + sizeof(d))

function load{T}(dm::DenseMatBlobs{T}, splitdim_idx::Int)
    splitnum = splitidx(dm, splitdim_idx)
    p = dm.splits[splitnum]
    range = p.first
    bid = p.second
    load(dm.coll, bid), range
end

getindex{T}(dm::DenseMatBlobs{T}, i::Int) = getindex(dm, ind2sub(size(dm), i)...)

function getindex{T}(dm::DenseMatBlobs{T}, i1::Int, i2::Int)
    part, range = load(dm, i2)
    part[relidx(range, i1, i2)...]
end

function getindex{T}(dm::DenseMatBlobs{T}, ::Colon, i2::Int)
    part, range = load(dm, i2)
    reli1, reli2 = relidx(range, 1, i2)
    @inbounds v = part[:,reli2]
    v
end

function getindex{T}(dm::DenseMatBlobs{T}, ::Colon, idxs::Vector{Int})
    res = Array(T, dm.sz[1], length(idxs))
    sp = sortperm(idxs)
    sidxs = idxs[sp]
    remain = length(idxs)
    startpos = 1
    endpos = 0
    for spair in dm.splits
        (remain > 0) || break
        r,b = spair.first, spair.second
        endpos = searchsortedlast(sidxs, last(r))
        (endpos < startpos) && continue
        A = load(dm.coll, b)
        ir = startpos:endpos
        fidxs = sidxs[ir]
        @inbounds res[:, sp[ir]] = A[:, fidxs-start(r)+1]
        remain -= length(ir)
        startpos = endpos+1
    end
    res
end

setindex!{T}(dm::DenseMatBlobs{T}, v::T, i::Int) = setindex!(dm, v, ind2sub(size(dm), i)...)

function setindex!{T}(dm::DenseMatBlobs{T}, v::T, i1::Int, i2::Int)
    part, range = load(dm, i2)
    part[relidx(range, i1, i2)...] = v
end

function setindex!{T}(dm::DenseMatBlobs{T}, v, ::Colon, i2::Int)
    part, range = load(dm, i2)
    reli1, reli2 = relidx(range, 1, i2)
    @inbounds part[:,reli2] = v
    v
end

function *{T1,T2}(A::Matrix{T1}, B::DenseMatBlobs{T2})
    m,n = size(B)
    (size(A, 2) == m) || throw(DimensionMismatch("A has dimensions $(size(A)) but B has dimensions $(size(B))"))
    Rm = size(A,1)
    res = Array(promote_type(T1,T2), Rm, n)
    for idx in 1:length(B.splits)
        p = B.splits[idx]
        part, r = load(B, first(p.first))
        Base.A_mul_B!(view(res, 1:Rm, r), A, part)
    end
    res
end

function deserialize{T}(s::SerializationState, ::Type{DenseMatBlobs{T}})
    metadir = deserialize(s)
    sz = deserialize(s)
    splits = deserialize(s)

    coll_id = deserialize(s)
    coll_mut = deserialize(s)
    coll_reader = deserialize(s)
    coll_blobs = deserialize(s)
    coll_maxcache = deserialize(s)

    coll = BlobCollection(Matrix{T}, coll_mut, coll_reader; maxcache=coll_maxcache, strategy=maxmem, id=coll_id)
    coll.blobs = coll_blobs
    DenseMatBlobs{T}(metadir, sz, splits, coll)
end

function DenseMatBlobs{Tv}(::Type{Tv}, metadir::AbstractString; maxcache::Int=def_cache())
    T = Matrix{Tv}
    io = FileBlobIO(Array{Tv}, true)
    mut = Mutable(BYTES_128MB, io)
    coll = BlobCollection(T, mut, io; maxcache=maxcache, strategy=maxmem)
    DenseMatBlobs{Tv}(metadir, (0,0), Pair[], coll)
end

function append!{Tv}(dm::DenseMatBlobs{Tv}, M::Matrix{Tv})
    m,n = size(M)
    if dm.sz[1] > 0
        (dm.sz[1] == m) || throw(BoundsError("DenseMatBlobs with size $(size(dm))", (m,n)))
    else
        dm.sz = (m,0)
    end
    if isempty(dm.splits)
        dm.sz = (m, n)
        idxrange = 1:n
    else
        old_split_dim = dm.sz[2]
        new_split_dim = old_split_dim + n
        idxrange = (old_split_dim+1):new_split_dim
        dm.sz = (dm.sz[1], new_split_dim)
    end

    fname = joinpath(dm.metadir, string(length(dm.splits)+1))
    meta = FileMeta(fname, 0, sersz(M))

    blob = append!(dm.coll, Matrix{Tv}, meta, StrongLocality(myid()), Nullable(M))
    push!(dm.splits, idxrange => blob.id)
    @logmsg("appending blob $(blob.id) of size: $(size(M)) for idxrange: $idxrange, sersz: $(meta.size)")
    blob
end

DenseMatBlobs(metadir::AbstractString; maxcache::Int=def_cache()) = matblob(metadir; maxcache=maxcache)

function save(dm::DenseMatBlobs, wrkrs::Vector{Int}=Int[])
    isempty(wrkrs) ? save(dm.coll) : save(dm.coll, wrkrs)
    open(joinpath(dm.metadir, "meta"), "w") do io
        serialize(SerializationState(io), dm)
    end
    nothing
end

end # module


text:
# guard against some future changes silently making simple vector graphics be
# rasterized if they are using features unsupported by the SVG spec
function svg_isnt_rasterized(x)
    mktempdir() do path
        path = joinpath(path, "test.svg")
        save(path, x)
        # this is rough but an easy way to catch rasterization,
        # if an image element is present in the svg
        return !occursin("<image id=", read(path, String))
    end
end

@testset "SVG rasterization" begin
    @test svg_isnt_rasterized(Scene())
    @test svg_isnt_rasterized(begin f = Figure(); Axis(f[1, 1]); f end)
    @test svg_isnt_rasterized(scatter(1:3))
    @test svg_isnt_rasterized(lines(1:3))
    @test svg_isnt_rasterized(heatmap(rand(5, 5)))
    @test !svg_isnt_rasterized(image(rand(5, 5)))
end


text:
addprocs() #Can be turned on for faster calculations
using Base.Test
using ThreePhotons

function test_file(filename)
  print_with_color(:green, "Testing: $filename\n")
  include(filename)
end

#General test parameters
K = 15
L = 12
N = 16
KMAX = 38
LMAX = 25

test_file("test_structure.jl")
# test_file("test_cubic.jl")
# test_file("test_spherical_harmonics.jl")
# test_file("test_datagen.jl")
# test_file("test_correlations.jl")
test_file("test_cuda.jl")
# test_file("test_determination.jl")
# test_file("test_phases.jl")
# test_file("test_data_processing.jl")
# test_file("test_utilities.jl")

print_with_color(:green, "Testing complete!\n")


text:
using FriendlyDoLoop
using Test

@testset "FriendlyDoLoop" begin

  include("FoosTests.jl")

  include("BarsTests.jl")

end