#	differential equation solver;
#	solve differential equation and fractional differential equation for unknown y with independent variable x;
#	dsolve(eq,y)
#	there are 3 way to input derivative y: y', y(1,x), or ds(y,x,1)
#	there are 3 way to input second order derivative y: y'', y(2,x), or ds(y,x,2)
#	e.g. dsolve( y' = x*y, y)
#	dsolve( y(0.5,x)=y^2 )
#	dsolve(a_,y_,x_,q_) = If(isfree(a,y),d(a,x,-q), If(isfree(a/y,y),mittag(q,d(a/y,x,-q)), If(isfree(a,x), solve(d(1/a, y,-q)=x^q/q!,y) )));



dsolve(a_ and b_, y_,x_,p_) := dsolve(a, y,x,p) and dsolve(b, y,x,p);
dsolve(a_ and b_ and c_ and d_, y_,x_,p_) := dsolve(a, y,x,p) and dsolve(b, y,x,p) and dsolve(c, y,x,p) and dsolve(d, y,x,p);

dsolve(a_*b_*c_,y_,x_,1):=if(hasnot(a,y) and has(b,c,y) and hasnot(b,c,x), int(1/b/c,y)=int(a,x)+C_1, 
	if(has(a,y) and hasnot(a,x) and has(b,x) and hasnot(b,c,y), int(1/a,y)=int(b*c,x)+C_1 ));

dsolve(c_+y_^m_, y_,x_,p_) := if(p>1 and has(c,x),block(s:=solve(y^m+c,y), if(hasnot(d(s,x,round(p-0.6)),x), s)));
dsolve(c_+b_*y_^m_, y_,x_,p_) := if(p>1 and has(c,x),block(s:=solve(y^m+c/b,y), if(hasnot(d(s,x,round(p-0.6)),x), s)));
dsolve(a_+y_, y_,x_,q_) := if(hasnot(a,y), gsolution(1,y,x,q)+If(isfree(a,x),-a,psolution(1,a,y,x,q)) );
dsolve(a_+b_*y_, y_,x_,q_) := if(hasnot(a,y), gsolution(b,y,x,q)+psolution(b,a,y,x,q));
dsolve(a_*y_+y_, y_,x_,q_) := If(isfree(a,y), C_1*mittag(q,d(a+1,x,-q)*q!) );

dsolve(b_*y_^n_, y_,x_,q_) := If(isfree(b,y),If(isfree(b,x), 
	if(n==-1 and q==2, -sqrt(-b)*exp(-inverseerf(sqrt(2)/sqrt(pi)*x+C_1)^2) and sqrt(-b)*exp(-inverseerf(sqrt(2)/sqrt(pi)*x+C_1)^2),
	block(f:=(1/b*fallingfactorial(q/(1-n),q))^(1/(n-1)),g:=(C_1+x)^(q/(1-n)),
	if(iseven(q), f*(C_1-x)^(q/(1-n)) and f*g, 
	if(isodd(n), -f*g and f*g, f*g )))), 
	if(n==-1, block(dd:=d(b,x,q),if(dd<>0 and hasnot(dd,x), -b/sqrt(dd) and b/sqrt(dd) )),
	block(g:=int(simplify(1/b),x,q),pp:=simplify(b^2*g^n),if(pp==1,g,if(pp==-1,-g ))) 
	)));
#dsolve(b_*y_^n_, y_,x_,2) := If(isfree(b,n,y),If(isfree(b,x),((2n+2)/b)^(1/(n-1))/((n-1)*(C_1-x))^(2/(n-1)) and ((2n+2)/b)^(1/(n-1))/((n-1)*(C_1+x))^(2/(n-1)), 
	(n/fallingfactorial(-n,2)*int(b,x,2))^(1/(2-n)) ));
#dsolve(b_/y_, y_,x_,p_) := block(dd:=d(b,x,p),if(dd<>0 and hasnot(dd,x), -b/sqrt(dd) and b/sqrt(dd) ));
#dsolve(b_/y_, y_,x_,2) := if(isconstant(b), -sqrt(-b)*exp(-inverseerf(sqrt(2)/sqrt(pi)*x+C_1)^2) and sqrt(-b)*exp(-inverseerf(sqrt(2)/sqrt(pi)*x+C_1)^2) );
dsolve(b_*x_*y_^n_, y_,x_,q_) := If(isconstant(b,n), block(k:=(1/b*fallingfactorial((1+q)/(1-n),q))^(1/(n-1))*(x)^((1+q)/(1-n)), if(isodd(n),-k and k,k)) );
dsolve(b_*x_^m_*y_^n_, y_,x_,q_) := If(isconstant(b,n), block(k:=(1/b*fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*(x)^((m+q)/(1-n)), if(isodd(n),-k and k,k)) );
#dsolve(b_*(c_+y_)^n_, y_,x_,q_) := if(isconstant(b,c),block(p:=(1/b*fallingfactorial(q/(1-n),q))^(1/(n-1)),g:=(C_1+x)^(q/(1-n)),
	if(iseven(q),p*(C_1-x)^(q/(1-n))-c and p*g-c, if(isodd(n),-p*g-c and p*g-c, p*g-c ))));
#dsolve(b_/(c_+y_), y_,x_,p_) := if(isconstant(c),block(dd:=d(b,x,p),if(dd<>0 and hasnot(dd,x), -b/sqrt(dd)-c and b/sqrt(dd)-c,
	dsolve(b/y,y,x,p)-c )));
dsolve(b_*(c_+y_), y_,x_,p_) := if(hasnot(b,y) and d(c,x,p)==0,  block(f:=dsolve(b*y,y,x,p), if(-left(f)==right(f), left(f)-c and right(f)-c,f-c  )));
dsolve(b*(c_+a_*y_), y_,x_,p_) := if(hasnot(b,y) and isconstant(a) and d(c,x,p)==0,  block(f:=dsolve(b*y,y,x,p), 
	if(-left(f)==right(f), left(f)/sqrt(a)-c/a and right(f)/sqrt(a)-c/a,f/sqrt(a)-c/a  )));
dsolve(b_*(c_+y_)^n_, y_,x_,p_) := if(hasnot(b,y) and d(c,x,p)==0,  block(f:=dsolve(b*y^n,y,x,p), if(-left(f)==right(f), left(f)-c and right(f)-c,f-c  )));
dsolve(b_*(c_+a_*y_)^n_, y_,x_,p_) := if(hasnot(b,y) and isconstant(a) and d(c,x,p)==0,  block(f:=dsolve(b*a^n*y^n,y,x,p), 
	if(-left(f)==right(f), left(f)-c/a and right(f)-c/a,f-c/a  )));
#dsolve(b_/(c_+a_*y_), y_,x_,p_) := if(isconstant(a,c),block(dd:=d(b,x,p),if(dd<>0 and hasnot(dd,x), -b/sqrt(a*dd)-c/a and b/sqrt(a*dd)-c/a ) ));
#dsolve(b_*x_^m_*(f_+y_)^n_, y_,x_,q_) := If(isconstant(b,f), block(k:=(1/b*fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*(x)^((m+q)/(1-n)), if(isodd(n),-k-f and k-f,k-f)) );
dsolve(b_*(c_+x_)^m_*y_^n_, y_,x_,q_) := If(isconstant(b,c,n), block(k:=(1/b*fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*(c+x)^((m+q)/(1-n)), if(isodd(n),-k and k,k)) );
dsolve(b_*(c_+x_)^m_*(f_+y_)^n_, y_,x_,q_) := If(isconstant(b,c,f,n), block(k:=(1/b*fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*(c+x)^((m+q)/(1-n)), if(isodd(n),-k-f and k-f,k-f)) );
dsolve(b_*y_, y_,x_,n_) := gsolution(b,y,x,n);
dsolve(b_*exp(y_), y_,x_,n_) := if(isconstant(b), if(iseven(n),log((-1)^n*n!*C_2^n/b/(C_1-C_2*x)^n) and log((-1)^n*n!*C_2^n/b/(C_1+C_2*x)^n), log((-1)^n*n!/b/(C_1+x)^n)));
dsolve(b_*exp(c_*y_), y_,x_,n_) := if(isconstant(b,c), if(iseven(n),log((-1)^n*n!*C_2^n/b/c/(C_1-C_2*x)^n) and log((-1)^n*n!*C_2^n/b/c/(C_1+C_2*x)^n), log((-1)^n*n!/b/c/(C_1+x)^n))/c );
dsolve(b_*exp(c_+y_), y_,x_,n_) := if(isconstant(b),if(d(c,x,n)==0,if(iseven(n),log((-1)^n*n!/b*C_2^n/(C_1-C_2*x)^n)-c and log((-1)^n*n!/b*C_2^n/(C_1+C_2*x)^n)-c, log((-1)^n*n!/b/(C_1+x)^n)-c),
	if(d(c,x,2)==0,-log(C_1-b/d(c,x)*exp(c)) )));

dsolve(x_^m_*y_^n_, y_,x_,q_) := if(isconstant(m,n),block(p:=(fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*x^((m+q)/(1-n)),if(isodd(n),-p and p, p )));
dsolve(x_*y_^n_, y_,x_,q_) := if(isconstant(n),block(p:=(fallingfactorial((1+q)/(1-n),q))^(1/(n-1))*x^((1+q)/(1-n)),if(isodd(n),-p and p, p )));
#dsolve(x_*(c_+y_)^n_, y_,x_,q_) := if(isconstant(c),block(p:=(fallingfactorial((1+q)/(1-n),q))^(1/(n-1))*x^((1+q)/(1-n)),if(isodd(n),-p*c and p-c, p-c )));
#dsolve(x_^m_*(c_+y_)^n_, y_,x_,q_) := if(isconstant(c),block(p:=(fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*x^((m+q)/(1-n)),if(isodd(n),-p-c and p-c, p-c )));
dsolve((b_+x_)^m_*y_^n_, y_,x_,q_) := If(isfree(b,x), block(p:=(fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*(b+x)^((m+q)/(1-n)),if(isodd(n),-p and p, p )));
dsolve((b_+x_)*(c_+y_)^n_, y_,x_,q_) := If(isfree(b,x) and d(c,x,q)==0, block(p:=(fallingfactorial((1+q)/(1-n),q))^(1/(n-1))*(b+x)^((1+q)/(1-n)),if(isodd(n),-p-c and p-c, p-c )));
dsolve((b_+a_*x_)*y_^n_, y_,x_,q_) := If(isconstant(a,b), block(k:=(fallingfactorial((1+q)/(1-n),q)/a)^(1/(n-1))*(b/a+x)^((1+q)/(1-n)),if(isodd(n),-k and k,k)) );
dsolve((b_+a_*x_)*(c_+y_)^n_, y_,x_,q_) := If(isfree(a,b,x) and d(c,x,q)==0, block(k:=(fallingfactorial((1+q)/(1-n),q)/a)^(1/(n-1))*(b/a+x)^((1+q)/(1-n)),if(isodd(n),-k-c and k-c,k-c)) );
dsolve((b_+x_)*y_^n_, y_,x_,q_) := If(isfree(b,x), block(p:=(fallingfactorial((1+q)/(1-n),q))^(1/(n-1))*(b+x)^((1+q)/(1-n)),if(isodd(n),-p and p, p )));
dsolve((b_+c_*x_)*y_^n_, y_,x_,q_) := If(isfree(b,c,x), block(p:=(1/c*fallingfactorial((1+q)/(1-n),q))^(1/(n-1))*(b/c+x)^((1+q)/(1-n)),if(isodd(n),-p and p, p )));
dsolve((b_+x_)^m_*(c_+y_)^n_, y_,x_,q_) := If(isfree(b,x) and d(c,x,q)==0, block(p:=(fallingfactorial((m+q)/(1-n),q))^(1/(n-1))*(b+x)^((m+q)/(1-n)),if(isodd(n),-p-c and p-c, p-c )));
dsolve(exp(xx_)*y_^n_, y_,x_,q_) := if(hasnot(xx,y) and d(xx,x,2)==0, block(f:=(d(xx,x)/(1-n))^(q/(n-1))*exp(xx/(1-n)), if(isodd(n), -f and f, f )));

#dsolve(1/b_, y_,x_,n_) := if(has(b,y),if(has(b,x), dsolve(replace(b,x,zz),zz,y,n)=x,C_1+int(b,y,n)=x^n/n! ));
#dsolve(1/b_, y_,x_,n_) := if( has(b,y), dsolve(replace(b,x,zz),zz,y,n)=x);
dsolve((c_+y_)^n_, y_,x_,q_) := if(d(c,x,q)==0,block(p:=fallingfactorial(q/(1-n),q)^(1/(n-1)),g:=(C_1+x)^(q/(1-n)),
	if(iseven(q),p*(C_1-x)^(q/(1-n))-c and p*g-c, 
	if(isodd(n),-p*g-c and p*g-c, p*g-c ))));
dsolve((c_+b_*y_)^n_, y_,x_,q_) := if(d(c,x,q)==0,block(p:=fallingfactorial(q/(1-n),q)^(1/(n-1))*b^(n/(1-n)),g:=(C_1+x)^(q/(1-n)),
	if(iseven(q),p*(C_1-x)^(q/(1-n))-c/b and p*g-c/b, 
	if(isodd(n),-p*g-c/b and p*g-c/b, p*g-c/b ))));
#dsolve((c_+y_)^n_, y_,x_,q_) := if(d(c,x,q)==0, block(f:=dsolve(y^n,y,x,q), if(-left(f)==right(f), left(f)-c and right(f)-c, f-c )));
#dsolve(y_^2, y_,x_,p_) :=  (2p)!/p!/2*(-1/(C_1+x))^p;
dsolve(y_^n_, y_,x_,q_) := if(isconstant(n),block(f:=fallingfactorial(q/(1-n),q)^(1/(n-1)),g:=(C_1+x)^(q/(1-n)),
	if(iseven(q),f*(C_1-x)^(q/(1-n)) and f*g, 
	if(isodd(n),-f*g and f*g, f*g ))));
dsolve(y_, y_,x_,n_) := gsolution(1,y,x,n);
dsolve(exp(b_*y_), y_,x_,n_) := if(isconstant(b),if(iseven(n),log((-1)^n*n!/b*C_2^n/(C_1-C_2*x)^n) and log((-1)^n*n!/b*C_2^n/(C_1+C_2*x)^n), log((-1)^n*n!/b/(C_1+x)^n) )/b);
dsolve(exp(y_), y_,x_,n_) := if(iseven(n),log((-1)^n*n!*C_2^n/(C_1-C_2*x)^n) and log((-1)^n*n!*C_2^n/(C_1+C_2*x)^n), log((-1)^n*n!/(C_1+x)^n));
dsolve(exp(c_+y_), y_,x_,n_) := if(d(c,x,n)==0,if(iseven(n),log((-1)^n*n!*C_2^n/(C_1-C_2*x)^n)-c and log((-1)^n*n!*C_2^n/(C_1+C_2*x)^n)-c, log((-1)^n*n!/(C_1+x)^n)-c),
	if(d(c,x,n+1)==0,-log(C_1-exp(c)/d(c,x,n)) ));


dsolve(exp(y_)*x_^n_, y_,x_,2) := if(n<0, log((n+2)*x^(-n-2)));
dsolve(x_*y_, y_,x_,2) := C_1*Ai(x)+C_2*Bi(x);
dsolve(1/(c_+y_), y_,x_,2) := if(d(c,x,2)==0, exp(-inverseerf(i* sqrt(2)/sqrt(pi)*(C_1 + x))^2)-c and exp(-inverseerf(-i* sqrt(2)/sqrt(pi)*(C_1 + x))^2)-c );
dsolve(1/y_, y_,x_,2) := -exp(-inverseerf(i* sqrt(2)/sqrt(pi)*x+C_1)^2) and exp(-inverseerf(i* sqrt(2)/sqrt(pi)*x+C_1)^2);

dsolve(a_+b_*y_, y_,x_,2) := if(isfree(a,b,y), if(a==b, gsolution(b,y,x,2)-1, gsolution(b,y,x,2)+psolution(b,a,y,x,2) ));
dsolve(a_+y_, y_,x_,2) := if(isfree(a,y),if(isfree(a,x), C_1*exp(x)+C_2*exp(-x)-a, C_1*exp(x)+C_2*exp(-x)+psolution(1,a,y,x,2) ));
dsolve(x_+b_*y_, y_,x_,2) := if(isfree(b,x), gsolution(b,y,x,2)-x/b, gsolution(b,y,x,2)+psolution(b,x,y,x,2) );
dsolve(x_+y_, y_,x_,2) := C_1*exp(x)+C_2*exp(-x)-x;
dsolve(a_*y_^3+b_*y_, y_,x_,2) := sqrt(2b/a)*csch(C_1+sqrt(b)*x);
dsolve(a_*y_^3+y_, y_,x_,2) := sqrt(2/a)*csch(C_1+x);
dsolve(y_^3-y_, y_,x_,2) := sqrt(2)*sec(C_1+x);
dsolve(y_-y_^3, y_,x_,2) := sqrt(2)*sech(C_1+x);
dsolve(-y_-y_^3, y_,x_,2) := sqrt(-2)*sec(C_1+x);
dsolve(y_^3+y_, y_,x_,2) := sqrt(2)*csch(C_1+x);


dsolve(a_+y_^2, y_,x_,1) := if(isconstant(a), tan(sqrt(a)*x+C_1)*sqrt(a),if(hasnot(a,y),block(u:=gsolution(-a,y,x,2), -d(u,x)/u) ));
dsolve(a_+c_*y_^2, y_,x_,1) := if(isconstant(a,c), tan(sgn(c)*sqrt(a*c)*x+C_1)*sqrt(a/c),if(hasnot(a,y),block(u:=gsolution(d(c,x)/c,-a*c,y,x,2), -d(u,x)/u/c) ));
dsolve(c_+x_*y_^2, y_,x_,1) := if(isconstant(c), -Ai(-cbrt(c)*x)/AiPrime(-cbrt(c)*x)*abs(c)^(2/3) and -Bi(-cbrt(c)*x)/BiPrime(-cbrt(c)*x)*abs(c)^(2/3) );
dsolve(c_+b_*x_*y_^2, y_,x_,1) := if(isconstant(b,c), -Ai(-cbrt(b*c)*x)/AiPrime(-cbrt(b*c)*x)*abs(b*c)^(2/3)/b and -Bi(-cbrt(b*c)*x)/BiPrime(-cbrt(b*c)*x)*abs(b*c)^(2/3)/b );
dsolve(b_+y_^0.5,y_,x_,1):= if(hasnot(b,y),b^2 *(-sgn(b)*(W(exp(-x/(2 b)+ C_1)) + 1))^2);
#dsolve(b_-y_^0.5,y_,x_,1):= if(hasnot(b,y),b^2 *(W(-exp(-(x + C_1)/(2 b) - 1)/b) + 1)^2);
dsolve(b_+a_*y_^0.5,y_,x_,1):= if(hasnot(b,y),(-b/a*(W(C_1*exp(-x*a^2/(2 b)))+1))^2);
dsolve(a_+b_*y_,y_,x_,1):=If(isfree(a,b,y), If(isfree(a,b,x), C_1*exp(b*x)-a/b,	C_1*exp(integrate(b,x))+psolution(b,a,y,x,1) ));
dsolve(a_+x_/y_,y_,x_,1):=replace(int(1/(zz-1/zz-replace(a,x/y,zz)),zz),zz,y/x)-log(x)+C_1=0;
#dsolve(a_+y_/x_,y_,x_,1):=replace(int(1/replace(a,y/x,zz),zz),zz,y/x)-log(x)+C_1=0;
dsolve(x_+y_/x_,y_,x_,1):= C_1*x+x*x;
dsolve(x_+x_/y_,y_,x_,1):= -W(C_1*exp(-x*x/2))-1;
dsolve(a_+a_*y_,y_,x_,1):=If(isfree(a,y), C_1*exp(integrate(a,x))-1);
#dsolve(y_+b_*y_^n_,y_,x_,1):=If(isfree(b,y), block(p:=(C_1+(1-n)*integrate(b*exp((n-1)*x),x)),exp(x)*p^(1/(1-n)) ));
dsolve(y_+b_*y_^n_,y_,x_,1):=If(isfree(b,y), exp(x)*(C_1+(1-n)*integrate(b*exp((n-1)*x),x))^(1/(1-n)) );
dsolve(a_*y_+y_^n_,y_,x_,1):=If(isfree(a,y), block(dsolve:=integrate(a,x),p:=(C_1+(1-n)*integrate(exp((n-1)*dsolve),x)),exp(dsolve)*p^(1/(1-n)) ));
dsolve(a_*y_+b_*y_^n_,y_,x_,1):=If(isfree(a,b,y), block(dsolve:=integrate(a,x),p:=(C_1+(1-n)*integrate(b*exp((n-1)*dsolve),x)),exp(dsolve)*p^(1/(1-n)) ));
dsolve(y_+y_^n_,y_,x_,1):= exp(x)*(C_1-exp((n-1)*x))^(1/(1-n));
#dsolve(a_*y_+y_^n_,y_,x_,1):=If(isfree(a,y), exp(integrate(a,x))*(C_1+(1-n)*integrate(exp((n-1)*integrate(a,x)),x))^(1/(1-n)) );
dsolve(x_+b_/y_, y_,x_,1) := if(isconstant(b),(2^(1/3)* x *Ai((x^2 - 2 y)/(2 *2^(1/3) *b^(2/3))) - 2 b^(1/3)* AiPrime((x^2 - 2 y)/(2 *2^(1/3)* b^(2/3))))/(2^(1/3) *x *Bi((x^2 - 2y )/(2 *2^(1/3) *b^(2/3))) - 2 b^(1/3) *BiPrime((x^2 - 2 *y)/(2 *2^(1/3) *b^(2/3)))) + C_1 = 0);
dsolve(b_+b_/y_, y_,x_,1) :=  -W(-exp(-int(b,x) + C_1 - 1)) - 1;
dsolve(1+x_/y_, y_,x_,1) := ((5 + sqrt(5))* log(-2 y/x + sqrt(5) + 1) - (sqrt(5) - 5) *log(2 y/x + sqrt(5) - 1)) = C_1 -10 log(x);
dsolve(x_+y_, y_,x_,1) := C_1*exp(x)-x-1;
dsolve(a_/x_+y_/x_,y_,x_,1):= if(isconstant(a),C_1*x-a, if(isfunction(a) and argument(a,1)==y(1,x),C_1*x-replace(a,y(1,x),C_1), if(hasnot(a,y), C_1*x+int(a/x^2,x)*x )));
dsolve(-a_/x_+1/x_*y_,y_,x_,1):= if(isconstant(a),C_1*x+a, if(isfunction(a) and argument(a,1)==y(1,x),C_1*x+replace(a,y(1,x),C_1), if(hasnot(a,y), C_1*x-int(a/x^2,x)*x )));
#dsolve(b_*x_+y_^2, y_,x_,1) := if(isconstant(b), cbrt(b)* AiPrime(cbrt(-b) *x)/Ai(cbrt(-b) *x) and  cbrt(b) *BiPrime(cbrt(-b) *x)/Bi(cbrt(-b) *x) );
dsolve(c_*x_^n_+b_*y_^2, y_,x_,1) := if(isconstant(b,c),-sqrt(b*c)*besselJ((-n-1)/(n+2),2/(n+2)*sqrt(b*c)*x^(n/2+1))/besselJ(1/(n+2),2/(n+2)*sqrt(b*c)*x^(n/2+1))*x^(n/2) );
dsolve(a_+a_*y_^2, y_,x_,1) := tan(C_1+int(a,x));
dsolve(x_^n_+y_^2, y_,x_,1) := -besselJ((-n-1)/(n+2),2/(n+2)*x^(n/2+1))/besselJ(1/(n+2),2/(n+2)*x^(n/2+1))*x^(n/2) and besselJ((n+1)/(n+2),2/(n+2)*x^(n/2+1))/besselJ(-1/(n+2),2/(n+2)*x^(n/2+1))*x^(n/2);
dsolve(x_+y_^2, y_,x_,1) := AiPrime(-x)/Ai(-x) and BiPrime(-x)/Bi(-x);
dsolve(-x_+y_^2, y_,x_,1) := -AiPrime(x)/Ai(x) and -BiPrime(x)/Bi(x);
dsolve(exp(x_)+y_^2,y_,x_,1):= (C_1* sqrt(exp(x))*besselJ(1,2 sqrt(exp(x))) + 2 sqrt(exp(x))* besselY(1,2 sqrt(exp(x))))/(C_1*besselJ(0,2 sqrt(exp(x))) + 2besselY(0,2 sqrt(exp(x))));
dsolve(exp(x_)+exp(y_),y_,x_,1):= exp(x) - log(C_1 - Ei(exp(x)));
dsolve(exp(x_*(1+y_))*1/x_-1/x_*y_,y_,x_,1):= -log(C_1 - e^x)/x;

#dsolve(c_+b_*x_+y_^2, y_,x_,1) := if(isconstant(b,c), cbrt(b)* AiPrime((-b* x - c)/abs(b)^(2/3))/Ai((-b *x - c)/abs(b)^(2/3)) and cbrt(b) *BiPrime((-b *x - c)/abs(b)^(2/3))/Bi((-b *x - c)/abs(b)^(2/3)) );
#dsolve(c_+x_+y_^2, y_,x_,1) := if(isconstant(c), AiPrime(-x-c)/Ai(-x-c) and BiPrime(-x-c)/Bi(-x-c) );
dsolve(y_^2-x_^2+c_, y_,x_,1) := if(isconstant(c), exp(c*x^2)/(C_1 - 1/2 sqrt(pi)* erfi(sqrt(c)*x)/sqrt(c)) + c*x );
dsolve(a_+y_+c_*y_^2, y_,x_,1) := if(has(a,x),block(u:=gsolution(1+d(c,x)/c,-a*c,y,x,2), -d(u,x)/u/c) );
#dsolve(a_+b_*y_+y_^2, y_,x_,1) := if(has(a,x) or has(b,x),if(b*b==4a, sqrt(b/2/x)*tan(C_1+sqrt(b/2/x)*x)-b/2, block(u:=gsolution(b,-a,y,x,2), -d(u,x)/u) ));
dsolve(a_+b_*y_+y_^2, y_,x_,p_) := if(b*b==4a, ode((y+b/2)^2,y,x,p), if(p==1, block(u:=gsolution(b,-a,y,x,2), -d(u,x)/u) ));
dsolve(a_+b_*y_+c_*y_^2, y_,x_,p_) := if(b*b==4a*c, ode(c*(y+b/2/c)^2,y,x,p), if(p==1, block(u:=gsolution(b+d(c,x)/c,-a*c,y,x,2), -d(u,x)/u/c) ));
#dsolve(a_+b_*y_+c_*y_^2, y_,x_,1) := if(has(a,x) or has(b,x),if(b*b==4a*c, sqrt(b/(2c*x))*tan(C_1+sqrt(b/(2c*x))*x)-b/2*sqrt(1/c), block(u:=gsolution(b+d(c,x)/c,-a*c,y,x,2), -d(u,x)/u/c) ));
dsolve(b_*x_+y_+y_^2, y_,x_,1) := if(isconstant(b), AiPrime((-b*x+1/4)/abs(b)^(2/3))/Ai((-b*x+1/4)/abs(b)^(2/3))*cbrt(b)-1/2 and BiPrime((-b*x+1/4)/abs(b)^(2/3))/Bi((-b*x+1/4)/abs(b)^(2/3))*cbrt(b)-1/2 );
dsolve(x_+y_+y_^2, y_,x_,1) := AiPrime(-x+1/4)/Ai(-x+1/4)-1/2 and BiPrime(-x+1/4)/Bi(-x+1/4)-1/2 ;

dsolve(c_+x_+y_+y_^2, y_,x_,1) := if(isconstant(c), AiPrime(-x+1/4-c)/Ai(-x+1/4-c)-1/2 and BiPrime(-x+1/4-c)/Bi(-x+1/4-c)-1/2 );
dsolve(c_+b_*x_+y_+y_^2, y_,x_,1) := if(isconstant(b,c), AiPrime((-b*x+1/4-c)/abs(b)^(2/3))/Ai((-b*x+1/4-c)/abs(b)^(2/3))*cbrt(b)-1/2 and BiPrime((-b*x+1/4-c)/abs(b)^(2/3))/Bi((-b*x+1/4-c)/abs(b)^(2/3))*cbrt(b)-1/2 );

dsolve(m_/n_,y_,x_,1) :=if(has(n,x) and has(n,y), if(hasnot(m,x),expand(dsolve(expand(replace(n,x,zz)/m),zz,y,1))=x, block(g:=(d(m,y)+d(n,x)), 
	if(g==0,block(p:=int(n,y),p-int(m+d(p,x),x)+C_1=0), 
	if(hasnot(simplify(-g/n),y),block(f:=exp(int(simplify(-g/n),x)),p:=int(f*n,y), p-int(f*m+d(p,x),x)+C_1=0), 
	if(hasnot(simplify(-g/m),x),block(f:=exp(int(simplify(-g/m),y)),p:=int(f*n,y), p-int(f*m+d(p,x),x)+C_1=0) 
	)) ))));
dsolve(n_/m_,y_,x_,1) :=if(has(m,x) and has(m,y), if(hasnot(n,x), expand(dsolve(expand(replace(m,x,zz)/n),zz,y,1))=x, block(g:=(d(n,y)+d(m,x)), 
	if(g==0,block(p:=int(m,y),p-int(n+d(p,x),x)+C_1=0), 
	if(hasnot(simplify(-g/m),y),block(f:=exp(int(simplify(-g/m),x)),p:=int(f*m,y), p-int(f*n+d(p,x),x)+C_1=0), 
	if(hasnot(simplify(-g/n),x),block(f:=exp(int(simplify(-g/n),y)),p:=int(f*m,y), p-int(f*n+d(p,x),x)+C_1=0) 
	)) ))));
dsolve(exp(y_)*xx_, y_,x_,1) := if(hasnot(xx,y), -log((C_1)-int(xx,x)));
dsolve(b_*y_^n_, y_,x_,1) := If(isfree(b,n,y) and isconstant(n), block(g:=(C_1-(n-1)*integrate(b,x))^(1/(1-n)), if(isodd(n), -g and g, g) ));
dsolve(b_*y_, y_,x_,1) := If(isfree(b,y),C_1*exp(integrate(b,x))  );

dsolve(y_^n_, y_,x_,1) := if(isconstant(n),if(isodd(n), -(C_1+(1-n)*x)^(1/(1-n)) and (C_1+(1-n)*x)^(1/(1-n)), (C_1+(1-n)*x)^(1/(1-n)) ));
dsolve((c_+y_)^2,y_,x_,1):= block(f:=d(c,x),if(f==0,1/(C_1-x)-c,if(hasnot(f,x), if(f<0,sqrt(-f)*tanh(C_1-sqrt(-f)*x)-c, sqrt(f)*tan(C_1+sqrt(f)*x)-c))));

dsolve((a_*x_+b_*y_)/(c_*x_+d_*y_),y_,x_,1) := if(isconstant(a,b,c,d), expand(replace(int(expand((c+d*yx)/(a+(b-c)*yx-d*yx^2)),yx),yx,y/x))+C_1+log(x)=0);
dsolve((x_-y_)/(x_+y_),y_,x_,1):= -sqrt(C_1+2x^2)-x and sqrt(C_1+2x^2)-x;
dsolve((x_+y_)/(x_-y_),y_,x_,1):= log(x^2+y^2)/2-atan(y/x)+C_1=0;
dsolve(x_/(x_+y_), y_,x_,1) := (-1/2-sqrt(5)/2)*x and (-1/2+sqrt(5)/2)*x;
dsolve(x_/(c_+x_+y_), y_,x_,1) := if(isconstant(c), (-1/2-sqrt(5)/2)*x-c and (-1/2+sqrt(5)/2)*x-c);
dsolve(x_/(a_*x_+b_*y_), y_,x_,1) := if(isconstant(a,b), -(a+sqrt(a*a+4b))/2/b*x and (-a+sqrt(a*a+4b))/2/b*x );
dsolve(x_/y_, y_,x_,1) := -sqrt(C_1+x^2) and sqrt(C_1+x^2);


#dsolve(c_+y(1,x_),y_,x_,2):=if(isfree(c,y),if(isfree(c,x), C_1+C_2*e^x-c*x, C_1+C_2*exp(x)-integrate(exp(x)*int(exp(-x)*c,x),x) ));
#dsolve(c_+b_*y(1,x_),y_,x_,2):=if(isfree(b,c,y),if(isfree(b,c,x), C_1+C_2*exp(b*x)-c*x/b, integrate(dsolve(c+b*y,y,x),x)+C_2 ));
#dsolve(c_+b_*y(1,x_),y_,x_,2):=if(isfree(b,c,y),if(isfree(b,c,x), C_1+C_2*exp(b*x)-c*x/b, C_1+C_2*exp(b*x)-integrate(exp(b*x)*int(exp(-b*x)*c,x),x) ));
#dsolve(c_+b_*y(1,x_),y_,x_,2):=if(isfree(b,c,y),if(isfree(b,c,x), C_1+C_2*exp(b*x)-c*x/b, C_1+C_2*int(exp(int(b,x)),x)+psolution(b,1,0,c,y,x,2) ));


dsolve(y(p_,x_)^n_+c_,y_,x_,q_):= if(hasnot(c,y),int(dsolve(y^n+c,y,x,q-p),x,p)+gsolution(0,y,x,p) );
dsolve(b_*y(p_,x_)^n_+c_,y_,x_,q_):= if(isfree(b,c,y),if(isfree(b,x),int(dsolve(b*y^n+c,y,x,q-p),x,p)+gsolution(0,y,x,p) )); 
dsolve(y(p_,x_)+c_,y_,x_,q_):= if(hasnot(c,y), int(dsolve(y+c,y,x,q-p),x,p)+gsolution(0,y,x,p));
dsolve(b_*y(p_,x_)+c_,y_,x_,q_):=  if(isconstant(b) and hasnot(c,y), int(dsolve(b*y+c,y,x,q-p),x,p)+gsolution(0,y,x,p) );
dsolve(b_*y(1,x_)+c_,y_,x_,q_):=if(isfree(b,c,y),if(isfree(b,c,x), C_1+C_2*exp(b*x)-c*x/b, gsolution(b,1,0,y,x,q)+psolution(b,1,0,c,y,x,q)),
	if(hasnot(b,x) and hasnot(c,y),  dsolve(int(b,y)+int(c,x),y,x,q-1) ));
dsolve(a_*y(n_,x_)*y(p_,x_)+c_,y_,x_,q_) := if(hasnot(a,c,y), block(g:=int(dsolve(a*y*y(p-n,x)+c,y,x,q-n),x,n),if(hasnot(g,y) and hasnot(g,C_1), gsolution(0,y,x,n)+g) ));
dsolve(b_*xy_*y(1,x_)+c_, xy_,x_,2):= if(hasnot(c,y), dsolve(b/2*y^2+int(c,x)+C_2,y,x,1) );
dsolve(y(1,x_)+c_,y_,x_,q_):=if(isfree(c,y),if(isfree(c,x), C_1+C_2*exp(x)-c*x, C_1+C_2*exp(x)+psolution(1,1,0,c,y,x,q) ));
dsolve(1/x_*y(p_,x_)+c_,y_,x_,q_):= if(hasnot(c,y) and q-p==1, int(int(c/x,x)*x,x,p)+C_1*x^q+C_2 );
#dsolve(xy_*y(1,x_)+c_, xy_,x_,2):= if(hasnot(c,y), dsolve(1/2y^2+int(c,x)+C_2,y,x,1) );
dsolve(xy_*y(1,x_)+c_, xy_,x_,2):= if(hasnot(c,y), block(u:=gsolution(-int(c,x)/2+C_2,y,x,2), -2d(u,x)/u) );
dsolve(y(n_,x_)*y(p_,x_)+c_,y_,x_,q_) := if(hasnot(c,y), block(g:=int(dsolve(y*y(p-n,x)+c,y,x,q-n),x,n),if(hasnot(g,y) and hasnot(g,C_1), gsolution(0,y,x,n)+g) ));

dsolve(-y(1,x_)^2/xy_+c_/xy_, xy_,x_,2):= -sqrt(2int(c,x,2)+C_1*x+C_2) and sqrt(2int(c,x,2)+C_1*x+C_2);
dsolve(y(1,x_)^2/xy_+c_/xy_, xy_,x_,2):= if(isconstant(c),-(sqrt(c)/c_1* tanh(c_2 + c_1* x))/sqrt(tanh( c_2 + c_1* x )^2 - 1) and (sqrt(c)/c_1* tanh(c_2 + c_1* x))/sqrt(tanh( c_2 + c_1* x)^2 - 1) );
dsolve(y(1,x_)^2/xy_+exp(x_)/xy_, xy_,x_,2):= -exp(x/2)*cosh(C_1*x+C_2)/C_1 and exp(x/2)*cosh(C_1*x+C_2)/C_1;
dsolve(y(1,x_)^2/xy_+exp(a_*x_)/xy_, xy_,x_,2):= if(hasnot(a,x),-exp(a*x/2)*cosh(C_1*x+C_2)/C_1 and exp(a*x/2)*cosh(C_1*x+C_2)/C_1);
dsolve(y(1,x_)^2+xy_, xy_,x_,2):= -sqrt(2)*int(1/(C_1*exp(2y)-2y-1),y)=C_2+x and sqrt(2)*int(1/(C_1*exp(2y)-2y-1),y)=C_2+x;

dsolve(a_*g_*y(p_,x_)+c_*g_, y_,x_,q_):= if(isconstant(a,c),-c/a*x^p/p!);
dsolve(y(p_,x_)+xy_, xy_,x_,q_) := gsolution(1,p,1,xy,x,q);
dsolve(y(p_,x_)+c_*xy_, xy_,x_,q_) := 	gsolution(1,p,c,xy,x,q);
dsolve(b_*y(p_,x_)+xy_, xy_,x_,q_) := 	if(hasnot(b,y),gsolution(b,p,1,xy,x,q));
dsolve(b_*y(p_,x_)+c_*xy_, xy_,x_,q_) := if(hasnot(b,y),gsolution(b,p,c,xy,x,q));
dsolve(b_*y(n_,x_)+c_*y(m_,x_),y_,x_,p_):= if(hasnot(b,c,y),gsolution(b,n,c,m,0,y,x,p) );
dsolve(a_*x_+b_*x_*y(p_,x_), y_,x_,q_) := If(isfree(a,b,x) and isfree(a,b,y), C_1*int(mittag(q-p,b*x^(1+q-p)/(1+q-p)!)-a/b,x,p));

dsolve(y(p_,x_)+xy_+d_, xy_,x_,q_) := if(isfree(d,xy),if(isfree(d,x), gsolution(1,p,1,xy,x,q)-d, gsolution(1,p,1,xy,x,q)+psolution(1,p,1,d,xy,x,q) ));
dsolve(y(p_,x_)+c_*xy_+d_, xy_,x_,q_) := if(isfree(c,d,xy),if(isfree(c,d,x), gsolution(1,p,c,xy,x,q)-d/c, gsolution(1,p,c,xy,x,q)+psolution(1,p,c,d,xy,x,q) ));
dsolve(b_*y(p_,x_)+xy_+d_, xy_,x_,q_) := if(isfree(b,d,xy),if(isfree(b,d,x), gsolution(b,p,1,xy,x,q)-d, gsolution(b,p,1,xy,x,q)+psolution(b,p,1,d,xy,x,q) ));
dsolve(b_*y(p_,x_)+c_*xy_+d_, xy_,x_,q_) := if(isfree(b,c,d,xy),if(isfree(b,c,d,x), gsolution(b,p,c,xy,x,q)-d/c, gsolution(b,p,c,xy,x,q)+psolution(b,p,c,d,xy,x,q) ));

dsolve(b_*y(p_,x_),y_,x_,q_):= gsolution(b,p,0,y,x,q);
dsolve(b_*y(p_,x_)^n_,y_,x_,q_):= if(isfree(b,y),if(isfree(b,x),
	if(n== -1, gsolution(0,y,x,p)+sqrt(b/(fallingfactorial((p+q)/2,p)*fallingfactorial((p+q)/2,q)))*(C_1+x)^((p+q)/2),
	b^(1/(1-n))*int(dsolve(y^n,y,x,q-p),x,p)+gsolution(0,y,x,p) ),
	if(q-p==1,block(g:=int(C_1-(n-1)*int(b,x)^(1/(1-n)),x,p), if(isodd(n), -g and g, g)+gsolution(0,y,x,p)) )),
	if(hasnot(b,x) and n<>2, int(((2-n)*int(b,y,q-p)+C_1)^(1/(n-q)),y,p)=(C_2+x)	
	));
#dsolve(b_*y(p_,x_)^n_,y_,x_,q_):= if(isfree(b,y),if(isfree(b,x),
	if(n== -1, gsolution(0,y,x,p)+sqrt(b/(fallingfactorial((p+q)/2,p)*fallingfactorial((p+q)/2,q)))*(C_1+x)^((p+q)/2),
	b^(1/(1-n))*int(dsolve(y^n,y,x,q-p),x,p)+gsolution(0,y,x,p) ),
	if(hasnot(b,x) and n<>2, int(((2-n)*int(b,y,q-p)+C_1)^(1/(n-q)),y,p)=(C_2+x)	
	)));
dsolve(b_/y(p_,x_),y_,x_,q_):= if(isconstant(b),C_2+sqrt(b/(fallingfactorial((p+q)/2,p)*fallingfactorial((p+q)/2,q)))*(C_1+x)^((p+q)/2),
	if(has(b,y) and hasnot(b,x), int(1/cbrt(C_1+int(b,y)),y)=cbrt(3)*x+C_1, 
	if(has(b,x) and hasnot(b,y), C_1-int(sqrt(C_1+2int(b,x)),x,p) and C_1+int(sqrt(C_1+2int(b,x)),x,p) 
	)));
dsolve(b_*xy_/y(p_,x_),xy_,x_,q_) := if(isconstant(b), b/(fallingfactorial(p+q,p)*fallingfactorial(p+q,q))*(C_1+x)^(p+q) );
dsolve(a_*y(n_,x_)*y(p_,x_),xy_,x_,q_) := if(hasnot(a,xy), int(dsolve(a*xy*y(p-n,x),xy,x,q-n),x,n)+gsolution(0,xy,x,n) );

dsolve(b_*y(1,x_)^2,y_,x_,2):=if(hasnot(b,y),if(hasnot(b,x), C_2 - log(C_1 - b*x)/b and C_2 - log(C_1 + b*x)/b, int(1/(C_1-int(b,x)),x)+C_2),
	int(exp(-int(b,y)),y)=C_2+C_1*x);
#dsolve(f_*y(1,x_)^2,y_,x_,2):= if(hasnot(f,y), int(1/(C_1-int(f,x)),x)+C_2, int(exp(-int(f,y)),y)=C_2+C_1*x );
dsolve(b_/y(1,x_),y_,x_,3):= if(isconstant(b), C_2-sqrt(-b*pi)/2erf(sqrt(2)*inverseerf(C_1+x*sqrt(2sgn(b))/sqrt(pi))) and C_2+sqrt(-b*pi)/2erf(sqrt(2)*inverseerf(C_1+x*sqrt(2sgn(b))/sqrt(pi))) );
dsolve(b_/xy_*y(1,x_)^2,xy_,x_,2) := if(isconstant(b), (c_2 + C_1 *x)^(1/(1 - b)) );

dsolve(exp(x_)*y(p_,x_)^n_,y_,x_,q_) := if(n== -1,sqrt(2^(p+q))*exp(x/2),(-1)^(q-n*p)*exp(-x))+gsolution(0,y,x,p);
#dsolve(xy_^m_*y(p_,x_)^n_,xy_,x_,q_) := if(m+n-1==0,C_1*exp(x), (fallingfactorial((n*p-q)/(n+m-1),p)^(-n)*fallingfactorial((n*p-q)/(n+m-1),q))^(1/(n+m-1))*(C_1+x)^((n*p-q)/(n+m-1)) );
dsolve(xy_^m_*y(p_,x_)^n_,xy_,x_,q_) := if(m+n-1==0,C_1*exp(x), (C_1+(fallingfactorial((n*p-q)/(n+m-1),p)^(-n)*fallingfactorial((n*p-q)/(n+m-1),q))^(1/(n*p-q))*x)^((n*p-q)/(n+m-1)) );
dsolve(xy_*y(1,x_)^2,xy_,x_,2):= sqrt(2)*inverseerf(C_1*x+C_2);
dsolve(xy_/y(p_,x_),xy_,x_,q_) := 1/(fallingfactorial(p+q,p)*fallingfactorial(p+q,q))*(C_1+x)^(p+q);
dsolve(exp(y_)/y(p_,x_),y_,x_,q_) :=log(d(-(p+q)*log(C_1+x),x,q)*d(-(p+q)*log(C_1+x),x,p) );
#dsolve(exp(x_)/y(p_,x_),y_,x_,q_) :=  sqrt(2^(p+q))*exp(x/2)+gsolution(0,y,x,p);
dsolve(y(n_,x_)^a_*y(p_,x_)^b_,xy_,x_,q_) := if(a+b-1==0, C_2*exp(C_1*x),int(dsolve(xy^a*y(p-n,x)^b,xy,x,q-n),x,n) )+gsolution(0,xy,x,n);
dsolve(y(n_,x_)*y(p_,x_),xy_,x_,q_) := int(dsolve(xy*y(p-n,x),xy,x,q-n),x,n)+gsolution(0,xy,x,n);
dsolve((c_+xy_)^m_/y(p_,x_),xy_,x_,q_) := if(isconstant(c),if(m==2,C_1*exp(x)-c, (fallingfactorial((p+q)/(2-m),p)*fallingfactorial((p+q)/(2-m),q))^(1/(m-2))*(C_1+x)^((p+q)/(2-m))-c ));

#dsolve(xy_^n_*y(1,x_)^m_,xy_,x_,2):= C_1 and x^((m-2)/(m+n-1)) ;
#dsolve(a_*y(p_,x_)^(-1),y_,x_,q_) := if(hasnot(a,y),-int(sqrt(C_1+2int(a,x,q-p)),x,p)+gsolution(0,y,x,p) and int(sqrt(C_1+2int(a,x,q-p)),x,p)+gsolution(0,y,x,p));
#dsolve(b_*xy_*y(1,x_)^2,y_,x_,2):= if(isconstant(b), (2/b)^0.5* inverseerf(x) );


dsolve(x_^n_*xy_^m_/y(p_,x_),xy_,x_,q_) := (fallingfactorial((n+p+q)/(2-m),p)*fallingfactorial((n+p+q)/(2-m),q))^(1/(m-2))*x^((n+p+q)/(2-m));
dsolve(x_^n_*xy_/y(p_,x_),xy_,x_,q_) := 1/(fallingfactorial(n+p+q,p)*fallingfactorial(n+p+q,q))*x^(n+p+q);
dsolve(x_*xy_^m_/y(p_,x_),xy_,x_,q_) := (fallingfactorial((1+p+q)/(2-m),p)*fallingfactorial((1+p+q)/(2-m),q))^(1/(m-2))*x^((1+p+q)/(2-m));
dsolve(x_*xy_/y(p_,x_),xy_,x_,q_) := 1/(fallingfactorial(1+p+q,p)*fallingfactorial(1+p+q,q))*x^(1+p+q);
dsolve(1/xy_*y(n_,x_)*y(p_,x_),xy_,x_,q_) := C_1*exp(C_2*x);
dsolve(exp(x_)*xy_/y(p_,x_),xy_,x_,q_) :=  exp(x);
dsolve(exp(x_)/y(p_,x_)*y(n_,x_),xy_,x_,q_) :=  exp(x)+gsolution(0,y,x,min(p,n));
dsolve(exp(x_)/y(m_,x_)*y(n_,x_),xy_,x_,q_) :=  exp(x)+gsolution(0,y,x,min(m,n));
dsolve(exp(x_)/y(m_,x_)/y(n_,x_),xy_,x_,q_) := cbrt(C_1+3^(m+n+q)*exp(x))+gsolution(0,y,x,min(m,n));
dsolve(exp(x_)*y(n_,x_)*y(p_,x_),xy_,x_,q_) :=  (-1)^(q-p-n)*exp(-x)+gsolution(0,y,x,n);

#dsolve(a_*b_*c_*y(p_,x_)^(-1),y_,x_,q_) := if(isconstant(a), a*dsolve(b*c/y(p,x),y,x,q) );
dsolve(a_*exp(x_)*y(n_,x_)*y(p_,x_),y_,x_,q_) := if(isconstant(a), (-1)^(q-p-n)/a*exp(-x)+gsolution(0,y,x,n) );
dsolve(a_*exp(x_)*xy_/y(p_,x_),xy_,x_,q_) := if(isconstant(a), a*exp(x));

dsolve((x_+n_)/x_*y(1,x_),y_,x_,2):= if(isconstant(n), C_1*Gamma(1+n,-x)+C_2);
dsolve(-(x_+n_)/x_*y(1,x_),y_,x_,2):= if(isconstant(n), C_1*Gamma(1-n,x)+C_2);


#dsolve(b_/x_*y(1,x_)+y(1,x_),y_,x_,2):= if(isconstant(b), C_1*Gamma(1+b,-x)+C_2);
#dsolve(a_*y(1,x_)+b_/x_*y(1,x_),y_,x_,2):= if(isconstant(a,b), C_1*(-a)^(-1-b)*Gamma(1+b,-a*x)+C_2);
#dsolve(xy_*y(1,x_)+x_, xy_,x_,p_):= (-1)^(p+1)*(2*(p-1))!/(p-1)!/(C_1+x)^(p-1)-i*x and (-1)^(p+1)*(2*(p-1))!/(p-1)!/(C_1+x)^(p-1)+i*x;
#dsolve(xy_*y(1,x_)+c_*x_, xy_,x_,p_):= if(c<0,-sqrt(-c)*x and sqrt(-c)*x,  2C_1*tan(C_1*x+C_2)-sqrt(-c)*x  and -2C_2/(C_1+C_2*x)  +sqrt(-c)*x );
#dsolve(xy_*y(1,x_)+c_*x_, xy_,x_,p_):= if(c<0,-sqrt(-c)*x and sqrt(-c)*x-2/x,  -2C_2/(C_1+C_2*x)-sqrt(-c)*x  and -2C_2/(C_1+C_2*x)  +sqrt(-c)*x );
#dsolve(b_*xy_*y(1,x_)+x_, xy_,x_,p_):= if(b<0,-sqrt(-1/b)*x and sqrt(-1/b)*x-2/b/x,  -2/b*C_2/(C_1+C_2*x)-sqrt(-1/b)*x  and -2/b*C_2/(C_1+C_2*x) +sqrt(-1/b)*x );
#dsolve(b_*xy_*y(1,x_)+c_*x_, xy_,x_,p_):= if(isconstant(b,c),if(sgn(b)==sgn(c), (-1)^(p-1)*(2*(p-1))!/(p-1)!*(C_2/(C_1+C_2*x))^(p-1)/(b)-sqrt(-c/b)*x  and (-1)^(p-1)*(2*(p-1))!/(p-1)!*(C_2/(C_1+C_2*x))^(p-1)/(b)+sqrt(-c/b)*x, 
	-sqrt(-c/b)*x and sqrt(-c/b)*x-(-1)^p*(2*(p-1))!/(p-1)!/b/x^(p-1) ));

dsolve(b_*y(1,x_)+xy_,xy_,x_,2):= if(hasnot(b,y),gsolution(b,1,1,xy,x,2) );
dsolve(b_*y(1,x_)+c_*xy_,xy_,x_,2):= if(hasnot(b,y),gsolution(b,1,c,xy,x,2) );
dsolve(b_*xy_*y(1,x_)+x_, xy_,x_,p_):= if(isconstant(b), if(p==2,sqrt(-1/b)*x and sqrt(-1/b)*x-2/b/x,  -sqrt(-1/b)*x and sqrt(-1/b)*x  ));
dsolve(b_*xy_*y(1,x_)+c_*x_, xy_,x_,p_):= if(isconstant(b,c), if(p==2,sqrt(-c/b)*x and sqrt(-c/b)*x-2/b/x,  -sqrt(-c/b)*x and sqrt(-c/b)*x ));
dsolve(xy_*y(1,x_)+c_*x_, xy_,x_,p_):= if(isconstant(c), if(p==2,sqrt(-c)*x and sqrt(-c)*x-2/x, -sqrt(-c)*x and sqrt(-c)*x ));
dsolve(xy_*y(1,x_)+x_, xy_,x_,p_):= if(p==2,i*x and i*x-2/x, -i*x and i*x);
dsolve(y(1,x_)+xy_,xy_,x_,2):= C_1*exp((1-sqrt(5))/2*x)+C_2*exp((1+sqrt(5))/2*x);
dsolve(y(1,x_)+c_*xy_,xy_,x_,2):= gsolution(1,1,c,xy,x,2) ;

dsolve(b_*y(1,x_)+xy_+d_,xy_,x_,2):= if(hasnot(d,xy), gsolution(b,1,1,xy,x,2) +psolution(b,1,1,d,xy,x,2));
dsolve(b_*y(1,x_)+c_*xy_+d_,xy_,x_,2):= if(hasnot(d,xy), gsolution(b,1,c,xy,x,2)+psolution(b,1,c,d,xy,x,2));
dsolve(y(1,x_)+xy_+d_,xy_,x_,2):= gsolution(1,1,1,xy,x,2)+psolution(1,1,1,d,xy,x,2);
dsolve(y(1,x_)+c_*xy_+d_,xy_,x_,2):= gsolution(1,1,c,xy,x,2)+psolution(1,1,c,d,xy,x,2);


#dsolve(a_*y(p_,x_)+b_*xy_+c_*xy_,xy_,x_,q_):= gsolution(a,p,b+c,xy,x,q);
#dsolve(y(p_,x_)+b_*xy_+c_*xy_,xy_,x_,q_):= gsolution(1,p,b+c,xy,x,q);

#dsolve(y(p_,x_)^n_,y_,x_,q_):= if(n== -1, C_2+1/sqrt(fallingfactorial((p+q)/2,p)*fallingfactorial((p+q)/2,q))*(C_1+x)^((p+q)/2),
	int(dsolve(y^n,y,x,q-p),x,p)+gsolution(0,y,x,p));
dsolve(y(p_,x_)^n_,y_,x_,q_):= 	int(dsolve(y^n,y,x,q-p),x,p)+gsolution(0,y,x,p);
#dsolve(y(p_,x_)^n_,y_,x_,q_):= if(n>0,gsolution(0,y,x,min(p,q)) and gsolution(0,y,x,min(p,q))+int(dsolve(y^(n^(sgn(q-p))),y,x,abs(q-p)),x,min(p,q)), gsolution(0,y,x,min(p,q))+int(dsolve(y^(n^(sgn(q-p))),y,x,abs(q-p)),x,min(p,q)) );
#dsolve(y(p_,x_)^n_,y_,x_,p_):= if(n>0,gsolution(0,y,x,p)+(C_2+x)^p/p!, gsolution(0,y,x,p)+(C_2+x)^p/p!);
#dsolve(y(1,x_)^n_,y_,x_,2) := C_2 - ((1-n)* (C_1 + x))^(1/(1 - n) + 1)/(n - 2);
#dsolve(y(1,x_)^2,y_,x_,2):= C_2 - log(C_1 + x) and C_2 - log(C_1 - x);
dsolve(y(1,x_)^p_,y_,x_,p_):= C_2 -(p-1)!^(1/(p-1))* log(C_1 + x) and C_2 -(p-1)!^(1/(p-1))* log(C_1 - x);
#dsolve(1/y(p_,x_),y_,x_,q_) := C_2+1/sqrt(fallingfactorial((p+q)/2,p)*fallingfactorial((p+q)/2,q))*(C_1+x)^((p+q)/2);

dsolve(y(p_,x_),y_,x_,q_):= C_1*exp(x)+gsolution(0,y,x,p);
dsolve((xy_-x_*y(1,x_))^n_,xy_,x_,1):= -C_1^n-C_1*x;
dsolve((-xy_+x_*y(1,x_))^n_,xy_,x_,1):= C_1*x-C_1^n;
dsolve((c_+xy_-x_*y(1,x_))^n_,xy_,x_,1):= if(isconstant(c),-C_1^n-C_1*x-c);
dsolve((c_-xy_+x_*y(1,x_))^n_,xy_,x_,1):= if(isconstant(c),C_1*x-C_1^n+c);
dsolve(y(1,x_),y_,x_,3):= C_1+C_2*exp(x)+C_3*exp(-x);
dsolve(y(2,x_),y_,x_,3):= C_1+C_2*exp(x)+C_3*x;
dsolve(1/y(1,x_),y_,x_,3):= C_2-i*sqrt(pi)/2erf(sqrt(2)*inverseerf(C_1+x*sqrt(-2)/sqrt(pi))) and C_2+i*sqrt(pi)/2erf(sqrt(2)*inverseerf(C_1+x*sqrt(-2)/sqrt(pi)));


dsolve(f_+b_*y(n_,x_)+c_*y(m_,x_),y_,x_,p_):=if(hasnot(f,x) and hasnot(b,c,f,y), gsolution(b,n,c,m,0,y,x,p)-if(hasnot(c,x),f/c*x^m/m!,0) );

dsolve(a_*y(p_,x_)+y(m_,x_)+c_*xy_, xy_,x_,q_) := gsolution(a,p,1,m,c,xy,x,q);
dsolve(a_*y(p_,x_)+b_*y(m_,x_)+xy_, xy_,x_,q_) := gsolution(a,p,b,m,1,xy,x,q);
dsolve(a_*y(p_,x_)+b_*y(m_,x_)+c_*xy_, xy_,x_,q_) := gsolution(a,p,b,m,c,xy,x,q);

dsolve(f_+a_*xy_+b_*y(n_,x_)+c_*y(m_,x_),xy_,x_,p_):=if(hasnot(b,c,f,xy), gsolution(b,n,c,m,a,xy,x,p)+psolution(b,n,c,m,a,f,xy,x,p) );


dsolve(a_=b_,y_):=dsolve(replace(toy(a-b),y(x),y),y);
dsolve(a_=b_,y(x_)):=dsolve(replace(toy(a-b),y(x),y),y);
#dsolve(b_=a_*y(n_,x_),y_):=dsolve(b/a,y,x,n);
#dsolve(a_*y(n_,x_)=b_,y_):=dsolve(b/a,y,x,n);
#dsolve(y(n_,x_)=b_,y_):=dsolve(b,y,x,n);
dsolve(x(n_,t_)=a_,x_) := dsolve(tox(a),x,t,n);


#dsolve(b_+x(n_,t_),d_+y(n_,t_)):= (x=pdsolve(-b,x,t,n),  y=pdsolve(-d,y,t,n));
#dsolve(b_+x(n_,t_),d_+x_+y(n_,t_)):= block(xx:=pdsolve(-b,x,t,n), x=xx,  y=pdsolve(-d-xx,y,t,n));
#dsolve(b_+x(n_,t_),d_+y(m_,t_)):= if(hasnot(b,y(t)) and hasnot(d,x(t)), (x=pdsolve(-b,x(t),t,n), y=pdsolve(-d,y(t),t,m)), 
	if(hasnot(b,y(t)),block(xx:=pdsolve(-b,x(t),t,n), (x=xx, y=pdsolve(-replace(d,x(t),xx),y(t),t,m))),
	if(hasnot(d,x(t)),block(xx:=pdsolve(-d,y(t),t,m), (y=xx, x=pdsolve(-replace(b,y(t),xx),x(t),t,n))) )));


dsolve(a_=b_,c_=d_):= dsolve(tox(a=b),toyt(c=d),x,y);
dsolve(ds(x_,t_)=f_,ds(y_,t_)=g_):= dsolve(x(1,t)=f,y(1,t)=g,x,y);

dsolve(x(m_,t_)=a_,y(n_,t_)=c_,x_,y_):= if(hasnot(a,y), block(xx:=solve(dsolve(a,x,t,m),x), (x=xx,  y=dsolve(replace(c,x,xx),y,t,n))), 
	if(hasnot(c,x), block(yy:=solve(dsolve(c,y,t,n),y), (x=dsolve(replace(a,y,yy),x,t,m),y=yy)),
	if(hasnot(a,c,t) and m==n, dsolve(c/a,y,x,n),
	if(hasnot(c,t), y=dsolve(c,y,x,n), 
	if(hasnot(a,t), swap(y=dsolve(swap(a,x,y),y,x,m),x,y) 
	)))));
dsolve(x(m_,t_)=c1_,y(n_,t_)=c2_+x_,x_,y_):= block(yy:=ode(d(c2,t,m)+replace(c1,x,y(n,t)-c2),y,t,m+n), (x=d(yy,t,n)-c2, y=yy) );
dsolve(x(m_,t_)=c1_,y(n_,t_)=c2_+x_+y_,x_,y_):= block(yy:=ode(y(m,t)+d(c2,t,m)+replace(c1,x,y(n,t)-y-c2),y,t,m+n), (x=d(yy,t,n)-c2-yy, y=yy) );
dsolve(x(m_,t_)=c1_,y(n_,t_)=c2_+x_+b_*y_,x_,y_):= block(yy:=ode(b*y(m,t)+d(c2,t,m)+replace(c1,x,y(n,t)-b*y-c2),y,t,m+n), (x=d(yy,t,n)-c2-b*yy, y=yy) );
dsolve(x(m_,t_)=c1_,y(n_,t_)=x_+b_*y_,x_,y_):= block(yy:=ode(b*y(m,t)+replace(c1,x,y(n,t)-b*y),y,t,m+n), (x=d(yy,t,n)-b*yy, y=yy) );
dsolve(x(m_,t_)=c1_+y_,y(n_,t_)=a_+y_,x_,y_):= block(xx:=ode(y(m,t)+toyt(d(replace(c1,x,y(t)),t,n))+replace(a-c1,x,y),y,t,m+n), (x=xx,y=expand(d(xx,t,m)-replace(c1,x,xx))) );
dsolve(x(m_,t_)=c1_+y_,y(n_,t_)=a_+b_*y_,x_,y_):= block(xx:=ode(b*y(m,t)+toyt(d(replace(c1,x,y(t)),t,n))+replace(a-b*c1,x,y),y,t,m+n), (x=xx,y=expand(d(xx,t,m)-replace(c1,x,xx))) );
dsolve(x(m_,t_)=c1_+f_*y_,y(n_,t_)=a_+y_,x_,y_):= block(xx:=ode(y(m,t)+toyt(d(replace(c1,x,y(t)),t,n))+replace(f*a-c1,x,y),y,t,m+n), (x=xx,y=expand(d(xx,t,m)/f-replace(c1,x,xx)/f)) );
dsolve(x(m_,t_)=c1_+f_*y_,y(n_,t_)=a_+b_*y_,x_,y_):= block(xx:=ode(b*y(m,t)+toyt(d(replace(c1,x,y(t)),t,n))+replace(f*a-b*c1,x,y),y,t,m+n), (x=xx,y=expand(d(xx,t,m)/f-replace(c1,x,xx)/f)) );

#dsolve(x(m_,t_)=a_*x_+b_*y_,y(m_,t_)=c_*x_+d_*y_,x_,y_):= if(isconstant(a,b,c,d) and (a+d)^2==4*(a*d-b*c),(x=C_1*mittag(m,(a+d)/2*t^m), y=C_1*mittag(m,(a+d)/2*t^m)), ode((c*x+d*y)/(a*x+b*y),y,x,m) );
dsolve(x(m_,t_)=a_*x_+y_,y(m_,t_)=x_+a_*y_,x_,y_):= if(isconstant(a),(x=C_1*mittag(m,(a+1)*t^m), y=C_1*mittag(m,(a+1)*t^m)) );
dsolve(x(m_,t_)=x_+b_*y_,y(m_,t_)=b_*x_+y_,x_,y_):= if(isconstant(b),(x=C_1*mittag(m,(1+b)*t^m), y=C_1*mittag(m,(1+b)*t^m)) );
dsolve(x(m_,t_)=a_*x_+b_*y_,y(m_,t_)=b_*x_+a_*y_,x_,y_):= if(isconstant(a,b),(x=C_1*mittag(m,(a+b)*t^m), y=C_1*mittag(m,(a+b)*t^m)) );
dsolve(x(m_,t_)=a_*x_+a_*y_,y(m_,t_)=b_*x_+b_*y_,x_,y_):= if(isconstant(a,b),(x=C_1*(a*mittag(m,(a+b)*t^m)+b), y=b*C_1*(mittag(m,(a+b)*t^m)-1)) );


dsolve(a_+b_*y(q_,x_),y_):= dsolve(-a/b,y,x,q);
dsolve(a_+y(q_,x_),y_):= dsolve(-a,y,x,q);

#dsolve(a_+b_*y(q_,x_),y_):= if(hasnot(a,y) and hasnot(b,x), d(b,y,-q)=gsolution(0,y,x,q)-d(a,x,-q), if(q<0, dsolve(y+expand(d(a/b,x,-q)),y), dsolve(-expand(a/b),y,x,q) ));
#dsolve(y(p_,x_)^n_+b_,y_):= dsolve(solve(replace(b,y(p,x),yy)+yy^n,yy),y,x,p);
#dsolve(a_*y(p_,x_)^n_+b_,y_):= dsolve(solve(expand(replace(b,y(p,x),yy)/a)+yy^n,yy),y,x,p);

#dsolve(integrate(f_,t_,a_,b_)=z_,y_) := dsolve((replace(f,t,b)-replace(f,t,a)-d(z,x)));
dsolve(c_+integrate(f_,t_,a_,b_),y_) := if(hasnot(c,y(x)),solve(replace(f,t,b)-replace(f,t,a)+d(c,x),y(x)), dsolve(replace(f,t,b)-replace(f,t,a)+d(c,x)));
dsolve(integrate(f_,t_,a_,b_)+z_,y_) := if(hasnot(z,y(x)),solve(replace(f,t,b)-replace(f,t,a)+d(z,x),y(x)), dsolve(replace(f,t,b)-replace(f,t,a)+d(z,x)));
dsolve(c_+integrates(f_,t_,a_,b_),y_) := if(hasnot(c,y(x)),solve(replace(f,t,b)-replace(f,t,a)+d(c,x),y), dsolve(replace(f,t,b)-replace(f,t,a)+d(c,x)));
dsolve(integrates(f_,t_,a_,b_)+z_,y_) := if(hasnot(z,y(x)),solve(replace(f,t,b)-replace(f,t,a)+d(z,x),y), dsolve(replace(f,t,b)-replace(f,t,a)+d(z,x)));
dsolve(c_+integrates(f_,t_,a_,b_)+z_,y_) := if(hasnot(c,z,y(x)),solve(replace(f,t,b)-replace(f,t,a)+d(c+z,x),y), dsolve(replace(f,t,b)-replace(f,t,a)+d(c+z,x)));
dsolve(c_+integrate(f_,t_,a_,b_)+z_,y_) := if(hasnot(c,z,y(x)),solve(replace(f,t,b)-replace(f,t,a)+d(c+z,x),y(x)), dsolve(replace(f,t,b)-replace(f,t,a)+d(c+z,x)));

#dsolve(y(p_,x_)^m_+b_,y_):= if(hasnot(b,y), block(f:=int((-b)^(1/m),x,p)+gsolution(0,y,x,p), if(iseven(m), -f and f,f)),
	dsolve((-b)^(1/m),y,x,p) );
#dsolve(a_*y(p_,x_)^m_+b_,y_):= if(hasnot(a,b,y), block(f:=int((-b/a)^(1/m),x,p)+gsolution(0,y,x,p), if(iseven(m), -f and f,f)),
	dsolve((-b/a)^(1/m),y,x,p) );
dsolve(a_*y(q_,x_)^m_+b_*y(p_,x_),xy_):=if(hasnot(b,xy),int(dsolve((-b/a*xy)^(1/m),xy,x,q-p),x,p), dsolve((-b/a)^(1/m)*y(p,x)^(1/m),xy,x,q) );
dsolve(y(q_,x_)^m_+b_*y(p_,x_),xy_):=   if(hasnot(b,xy),int(dsolve((-b*xy)^(1/m),xy,x,q-p),x,p), dsolve((-b)^(1/m)*y(p,x)^(1/m),xy,x,q) );
dsolve(y(p_,x_)^m_+b_*y(p_,x_),y_):=if(hasnot(b,y),if(m>0,gsolution(0,y,x,p) and gsolution(0,y,x,p)+int((-b)^(1/(m-1)),x,p), 
	gsolution(0,y,x,p)+int((-b)^(1/(m-1)),x,p)), dsolve((-b)^(1/(m-1)),y,x,p) );
dsolve(y(q_,x_)^2-xy_*y(p_,x_),xy_):= C_1*exp(x) and C_1;

dsolve(y(p_,x_)^m_+y(p_,x_)+c_,y_):= C_1+int(solve(y^m+y+c,y),x,p);
dsolve(y(p_,x_)^m_+b_*y(p_,x_)+c_,y_):= C_1+ode(solve(yy^m+b*yy+c,yy),y,x,p);
dsolve(y(p_,x_)^2+b_*y(p_,x_)+c_,y_):= if(hasnot(b,c,y),C_1-int(b/2+sqrt(b*b-4c)/2,x,p) and C_1+int(-b/2+sqrt(b*b-4c)/2,x,p), 
	expand(dsolve((-b-(b*b-4c)^0.5)/2,y,x,p)) );
dsolve(y(p_,x_)^2+b_*xy_*y(p_,x_)+c_*xy_^2,xy_):= if(b*b==4c, C_1*mittag(p,-p!/2*int(b,x,p)), 
	C_1*mittag(p,p!/2*int((-b+sqrt(b*b-4c)),x,p)) and C_2*mittag(p,-p!/2*int((b+(b*b-4c)^0.5),x,p)) );
dsolve(y(1,x_)^2+x_*y(1,x_)+c_*xy_,xy_):= if(isconstant(c),-(2+c)/4x^2);
dsolve(y(1,x_)^2+b_*x_*y(1,x_)+xy_,xy_):= if(isconstant(b),-(2b+1)/4x^2);
dsolve(y(1,x_)^2+b_*x_*y(1,x_)+c_*xy_,xy_):= if(isconstant(b,c),-(2b+c)/4x^2);
dsolve(y(1,x_)^2+b_*x_*y(1,x_)+c_*xy_+f_,xy_):= if(isconstant(b,c,f),-(2b+c)/4x^2-f/c);
dsolve(y(1,x_)^m_-x_*y(1,x_)+xy_,xy_):= -C_1^m+C_1*x;
dsolve(y(1,x_)^m_+x_*y(1,x_)-xy_,xy_):= C_1^m+C_1*x;
dsolve(y(1,x_)^2-x_*y(1,x_)+xy_,xy_):= -C_1^2+C_1*x and x^2/4;
dsolve(y(1,x_)^2+x_*y(1,x_)-xy_,xy_):= C_1^2+C_1*x and -x^2/4;

dsolve(y(p_,x_)^m_+b_*y(p_,x_)+c_*xy_+f_,xy_):=if(isconstant(f),dsolve(y(p_,x_)^m_+b_*y(p_,x_)+c_*xy_,xy_)-f/c);

dsolve(a_):=dsolve(replace(toy(a),y(x),y),y);