# psolution is particular solution of linear differential equation;

#14;
psolution(a_,n_,b_,m_,c_,q_,d_,f_,y_,t_,x_,u_,v_,p_):=if(isconstant(f),if(d==0, -f*t^p/p!, -f/d),
	if(hasnot(f,x) and hasnot(f,u) and hasnot(f,v) and has(f,t),psolution(-d,-f,y,t,p), 
	if(hasnot(f,t) and hasnot(f,u) and hasnot(f,v) and has(f,x),psolution(-d/a,-f/a,y,x,n),
	if(hasnot(f,t) and hasnot(f,x) and hasnot(f,v) and has(f,u),psolution(-d/b,-f/b,y,u,m),
	if(hasnot(f,t) and hasnot(f,x) and hasnot(f,u) and has(f,v),psolution(-d/c,-f/c,y,v,q)
	)))));
psolution(a_,n_,b_,m_,c_,q_,d_,c1_+c2_,y_,t_,x_,u_,v_,p_):=	psolution(a,n,b,m,c,q,d,c1,y,t,x,u,v,p)+psolution(a,n,b,m,c,q,d,c2,y,t,x,u,v,p);

#11;
psolution(a_,n_,b_,m_,c_,f_,y_,t_,x_,u_,p_):=if(isconstant(f),if(c==0, -f*t^p/p!, -f/c),
	if(hasnot(f,x) and hasnot(f,u) and has(f,t),psolution(-c,-f,y,t,p), 
	if(hasnot(f,t) and hasnot(f,u) and has(f,x),psolution(-c/a,-f/a,y,x,n),
	if(hasnot(f,t) and hasnot(f,x) and has(f,u),psolution(-c/b,-f/b,y,u,m)
	))));
psolution(a_,n_,b_,m_,c_,c1_+c2_,y_,t_,x_,u_,p_):= psolution(a,n,b,m,c,c1,y,t,x,u,p)+psolution(a,n,b,m,c,c2,y,t,x,u,p);


#10;
psolution1(a_,n_,b_,m_,c_,f_,y_,t_,x_,p_):=if(isconstant(f),if(c==0,-f*t^p/p!, -f/c),
	if(hasnot(f,x) and has(f,t),psolution(-c,-f,y,t,p), 
	if(hasnot(f,t) and has(f,x),psolution(-b/a,m,-c/a,-f/a,y,x,n) )));

psolution(a_,n_,b_,m_,c_,exp(tx_),y_,t_,x_,p_):= if(isconstant(a,b,c) and has(tx,t) and has(tx,x), -exp(tx)/(d(tx,t)^p+a*d(tx,x)^n+b*d(tx,x)^m+c) );
psolution(a_,n_,b_,m_,c_,c1_*c2_,y_,t_,x_,p_):=if(isconstant(c1),c1*psolution(a,n,b,m,c,c2,y,t,x,p),if(isconstant(c2),c2*psolution(a,n,b,m,c,c1,y,t,x,p) ));
#psolution(a_,n_,b_,m_,c_,t_*ax_,y_,t_,x_,p_):= if(p>1, psolution(b/a,m,c/a,ax/a,y,x,n)*t, if(c==0,psolution(a,n-1,b,t*ax,y,t,p) ));
#psolution(a_,n_,b_,m_,c_,tt_*x_,y_,t_,x_,p_):= if(p>1, psolution(-b/a,m,-c/a,-x/a,y,x,n)*t, psolution(c,tt,y,t,p)*x+if(m==1,c*b*int(tt,t,p)/2,0) );
#psolution(a_,n_,b_,m_,c_,tt_*x_,y_,t_,x_,p_):= if(m>=1, psolution(-b,m,c,tt*x,y,t,x,p) );
psolution(a_,n_,b_,m_,c_,t_*x_,y_,t_,x_,p_):=  if(p>1, psolution(b/a,m,c/a,x/a,y,x,n)*t, psolution(-b,m,-c,-t*x,y,t,x,p)  );
#psolution(a_,n_,b_,m_,c_,tt_*x_,y_,t_,x_,p_):=  if(n>1, if(p>1, psolution(-c/b,-x/b,y,x,m)*tt,  psolution(b,m,c,tt*x,y,t,x,p)  ));
psolution(a_,n_,b_,m_,c_,exp(x_)*tt_,y_,t_,x_,p_):= if(a-b-c==0 or isconstant(a,b,c), -psolution(-a-b-c,tt,y,t,p)*exp(x) );
psolution(a_,n_,b_,m_,c_,exp(x_)*at_,y_,t_,x_,p_):= if(a-b-c==0 or isconstant(a,b,c), -psolution(-a-b-c,at,y,t,p)*exp(x) );
psolution(a_,n_,b_,m_,c_,exp(t_)*xx_,y_,t_,x_,p_):= if(isconstant(a,b,c) and t<>x,psolution(-b/a,m,(-1-c)/a,-xx/a,y,x,n)*exp(t),0 );
psolution(a_,n_,b_,m_,c_,exp(t_)*ax_,y_,t_,x_,p_):= if(isconstant(a,b,c) and t<>x,psolution(-b/a,m,(-1-c)/a,-ax/a,y,x,n)*exp(t),0 );
#psolution(a_,n_,b_,m_,c_,exp(t)*x_,y_,t_,x_,p_):= if(isconstant(a,b,c), if(c==0 and m==1, (b+x)*exp(t), psolution(c-1,x,y,x,m)*exp(t) ));
#psolution(a_,n_,b_,m_,c_,exp(t)*x_,y_,t_,x_,p_):= if(isconstant(a,b,c), if(c==0 and m==1, (b+x)*exp(t), 
	if(c==1, psolution(b,m,0,x,y,x,n)*exp(t), if(m>=1, psolution(-b,m,c,exp(t)*x,y,t,p) ))));
psolution(a_,n_,b_,m_,c_,c1_*c2_*c3_,y_,t_,x_,p_):= if(isfree(c1,x) and hasnot(c1,t), c1*psolution(a,n,b,m,c,c2*c3,y,t,x,p), 
	if(isfree(c3,x) and hasnot(c3,t), c3*psolution(a,n,b,m,c,c1*c2,y,t,x,p) ));

psolution(a_,n_,b_,m_,c_,c1_+c2_,y_,t_,x_,p_):=	psolution(a,n,b,m,c,c1,y,t,x,p)+psolution(a,n,b,m,c,c2,y,t,x,p);


#9;
psolution(a_,n_,b_,m_,c_,d_,y_,x_,p_):= if(isconstant(p),if(c<>0 and isconstant(d),-d/c, if(has(d,t), psolution(b,m,c,d,y,t,p), if(has(d,x), psolution(b/a,m,c/a,d/a,y,x,n) ))),
	if(hasnot(a,c,d,x), psolution(-c/a,-d/a,y,t,n), if(hasnot(b,c,d,t), psolution(-c/b,-d/b,y,x,m) )));
psolution(a_,n_,b_,m_,c_,exp(x_),y_,x_,p_):= if(isconstant(a,b,c), -exp(x)/(a+b+c-1) );

psolution(a_,n_,b_,m_,c_,c1_*c2_,y_,t_,p_):= if(isconstant(c1), c1*psolution(a,n,b,m,c,c2,y,t,p), if(isconstant(c2), c2*psolution(a,n,b,m,c,c1,y,t,p) ));
psolution(a_,n_,b_,m_,c_,c1_*c2_*c3_,y_,t_,p_):= if(isconstant(c1), c1*psolution(a,n,b,m,c,c2*c3,y,t,p), if(isconstant(c2), c2*psolution(a,n,b,m,c,c1*c3,y,t,p) ));

psolution(a_,n_,b_,m_,c_,c1_+c2_,y_,x_,p_):=psolution(a,n,b,m,c,c1,y,x,p)+psolution(a,n,b,m,c,c2,y,x,p);

psolution(3,2,-3,1,1,exp(x_)*xx_,y_,x_,3):=exp(x)*int(xx,x,3);
psolution(3,2,-3,1,1,exp(x_),y_,x_,3):=exp(x)*x^3/6;


#8;
psolution2(b_,m_,c_,d_,y_,t_,x_,p_):= if(isconstant(d),if(c==0,-d*t^p/p!, -d/c),
	if(hasnot(d,x),psolution(-c,-d,y,t,p),
	if(has(d,x), if(hasnot(b,d,t),psolution(-c/b,-d/b,y,x,m), psolution(-c,-d,y,x,m)
	))));
psolution(b_,m_,0,d_,y_,t_,x_,m_):= 
	if(has(d,t) and has(d,x), if(isfunction(d), -d(argument(d,1),t)^m*int(d,t,m)/(d(argument(d,1),t)^m+b*d(argument(d,1),x)^m) ),
	if(hasnot(d,x), -int(d,t,m),
	if(has(d,x), if(hasnot(b,d,t), -int(d/b,x,m), -int(d,x,m)),
	-int(replace(replace(d,t,z),x,x-b*t+b*z),z,0,t) )));

psolution(b_,m_,c_,exp(tx_),y_,t_,x_,p_):= if(d(tx,t,2)==0 and d(tx,x,2)==0 and isconstant(b,c), if(d(tx,t)^p+b*d(tx,x)^m+c==0, -exp(tx)*t/p, -exp(tx)/(d(tx,t)^p+b*d(tx,x)^m+c) ));

psolution(b_,m_,c_,c1_*c2_,y_,t_,x_,p_):=if(isconstant(c1), c1*psolution(b,m,c,c2,y,t,x,p), if(isconstant(c2), c2*psolution(b,m,c,c1,y,t,x,p) ));

#psolution1(b_,1,0,d_,y_,t_,x_,1):= if(has(d,t) and has(d,x), -int(replace(replace(d,t,z),x,x-b*t+b*z),z,0,t) );
psolution1(a_,1,0,tt_*xx_,y_,t_,x_,1):= if(isconstant(tt), tt*psolution(a,1,0,xx,y,t,x,1), 
	if(hasnot(tt/a,t), -int(tt*xx/a,x),
	if(has(tt,t) and has(xx,x) and hasnot(tt,x) and hasnot(xx,t), -int(replace(xx,x,x-a*t+a*z)*replace(tt,t,z),z,0,t) )));

#psolution(b_,m_,c_,ax_*t,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, int(ax,x,m)*t/b-if(p==1,1/b/b*int(ax,x,2m),0),-1/2/c*x*exp(t)-1/4/c*exp(t)-exp(t)/2 ));
#psolution(b_,m_,c_,ax_*t,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, int(ax,x,m)*t/b-if(p==1,1/b/b*int(ax,x,2m),0), if(p>1, psolution(c/b,ax/b,y,x,m)*t, psolution(c/b,ax/b,y,x,m)*t+c*int(ax/b,x,m)/2 )));
#psolution(b_,m_,c_,tt_*x,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, int(tt,t,p)*x-if(m==1,b*int(tt,t,2p),0), 1/b*exp(c*t) * integrate(exp(-c*t)* tt,t)*x+integrate(tt,t)/b/c/2 ));
#psolution(b_,m_,c_,tt_*x,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, int(tt,t,p)*x-if(m==1,b*int(tt,t,2p),0), if(m>1, psolution(c,tt,y,t,p)*x, -b*exp(c*t) * integrate(exp(-c*t)*tt,t)*x )));
#psolution(b_,m_,c_,tt_*x_,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, -int(tt,t,p)*x+if(m>1,gsolution(0,y,x,m),b*int(tt,t,2p)), if(m>1, psolution(-c,-tt,y,t,p)*x, psolution(-c,-tt,y,t,p)*x+psolution(-c,psolution(-c,-tt,y,t,p),y,t,p) )));

#psolution(b_,m_,c_,tt_*x_^n_,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, -int((x-b*t+b*z)^n*replace(tt,t,z),z,0,t) ));
psolution(b_,m_,c_,t_*x_^n_,y_,t_,x_,p_):= if(isconstant(b,c), if(c==0, -x^(m+n)*t/b/(n+m)+if(p>1,0,1/b/b*x^(n+2m)/(n+2m)/(n+m) )));
psolution(b_,m_,c_,t_*x_,y_,t_,x_,p_):= if(isconstant(b,c), if(m>1, if(p>1, -1/c*t*x, psolution(-c,-t,y,t,p)*x), if(p>1, psolution(c,x,y,x,m)*t/b, (1-c*t)*x/c^2+(-2+c*t)*b/c^3 )),0);
psolution(b_,m_,c_,exp(t_)*ax_,y_,t_,x_,p_):= if(isconstant(b,c) and has(ax,x), -psolution(-1/b-c/b,ax/b,y,x,m)*exp(t) );
psolution(b_,m_,c_,exp(t_)*xx_,y_,t_,x_,p_):= if(isconstant(b,c) and has(xx,x), -psolution(-1/b-c/b,xx/b,y,x,m)*exp(t) );
psolution(b_,m_,c_,exp(x_)*at_,y_,t_,x_,p_):= if(isconstant(b,c) or c+b==0 and has(at,t), -psolution(-b-c,at,y,t,p)*exp(x) );
psolution(b_,m_,c_,exp(x_)*tt_,y_,t_,x_,p_):= if(isconstant(b,c) or c+b==0 and has(tt,t), -psolution(-b-c,tt,y,t,p)*exp(x) );
#psolution(b_,m_,c_,exp(x_)*tt_,y_,t_,x_,p_):= if(b+c==0, -exp(x)*int(tt,t,p) );
#psolution(b_,m_,-1,exp(t_)*xx_,y_,t_,x_,p_):= if(isconstant(b), -exp(t)*int(xx,x,m)/b );

psolution(b_,m_,0,tt_*x_^n_,y_,t_,x_,1):= if(isconstant(b), -int((x-b*t+b*z)^n*replace(tt,t,z),z,0,t) );
#psolution(a_,1,0,c_*tt_*xx_,y_,t_,x_,1):= if(isconstant(a,c) and has(tt,t) and has(xx,x) and hasnot(tt,x) and hasnot(xx,t), expand(-c/a*int(replace(tt,t,t-x/a+z/a)*replace(xx,x,z),z,0,x)) );
#psolution(b_,m_,0,ax_*t_,y_,t_,x_,p_):= if(isconstant(b) and has(ax,x), if(p==1,-int(ax,x,m)*t/b+int(ax,x,2m)/b^2,-int(ax,x,m)*t/b ));
#psolution(b_,m_,0,tt_*x_,y_,t_,x_,p_):= if(isconstant(b) and has(tt,t), if(m==1,-int(tt,t,p)*x+b*int(tt,t,2p),-int(tt,t,p)*x ));
psolution(b_,m_,0,t_*x_,y_,t_,x_,p_):= if(isconstant(b), if(p>1, if(m>1, -1/2*t^(1+p)/(1+p)!*x+1/2*x^(1+m)/(1+m)!*t/b, -t/b*x^(1+m)/(1+m)!+if(p>1,0,int(x,x,2m)/b^2 )), -t^(1+p)/(1+p)!*x+if(m>1,0,b*int(t,t,2p)) ));
#psolution(b_,m_,0,t*x,y_,t_,x_,p_):= if(isconstant(b), if(p>1, if(m>1, 1/2*int(t,t,p)*x+1/2*int(x,x,m)*t/b,  psolution(0,x/b,y,x,m)*t), if(p<1,int(t,t,p)*x-b*int(t,t,2p),  if(m<1,int(x,x,m)*t/b-int(x,x,2m)/b^2 ))));
#psolution(b_,m_,0,t_*x_,y_,t_,x_,1):= if(isconstant(b), if(m>1, b*t^2*x/2,int(x,x,m)*t/b-int(x,x,2m)/b^2 ));
#psolution(b_,1,0,t_*x_,y_,t_,x_,p_):= if(isconstant(b), if(p>1, x^2*t/2/b,int(t,t,p)*x-b*int(t,t,2p) ));
#psolution(b_,p_,0,t_*x_,y_,t_,x_,p_):= if(isconstant(b), int(t,t,p)*x/2-int(x,x,p)*t/2/b+if(p>1,0, -int(t,t,p)/2+int(x,x,p)/2/b ));
psolution(b_,1,0,1/t_*x_,y_,t_,x_,1):= if(isconstant(b),(b*t-x)*log(b*t)-x);
psolution(b_,p_,0,tx_^n_,y_,t_,x_,p_):= if(d(tx,t,2)==0 and d(tx,x,2)==0 and isconstant(b), if(d(tx,t)^p+b*d(tx,x)^p==0, -t*(tx)^(n+p-1)/(n+p-1)!, -(tx)^(n+p)/fallingfactorial(n+p,p)/(d(tx,t)^p+b*d(tx,x)^p) ));
psolution(b_,1,0,t_*1/x_,y_,t_,x_,1):= if(isconstant(b),(x-b*t)*log(x)/b^2-x/b^2);
psolution(b_*t_,1,0,x_,y_,t_,x_,1):= if(isconstant(b), b*t^3/3-t*x );
psolution(b_*t_,1,0,x_^2,y_,t_,x_,1) := if(isconstant(b), -4/30b^2*t^5+2/3b*t^3*x-x^2*t );
psolution(b_*t_^n_,1,0,x_,y_,t_,x_,1):= if(isconstant(b), b*t^(n+2)/(n+2)-t*x );
psolution(x1_*tt_,1,0,xx_*tt_,y_,t_,x_,p_):= if(hasnot(xx/x1,t), int(xx/x1,x), if(hasnot(xx*tt,x), int(xx*tt,t,p) ));
psolution(t_,1,0,x_,y_,t_,x_,1):= t^3/3-t*x;
psolution(t_,1,0,x_^2,y_,t_,x_,1):= -4/30t^5+2/3t^3*x-x^2*t;
psolution(t_^n_,1,0,x_,y_,t_,x_,1):= t^(n+2)/(n+2)-t*x;

psolution(b_,m_,c_,c1_*c2_*c3_,y_,t_,x_,p_):=if(isconstant(c1), c1*psolution(b,m,c,c2*c3,y,t,x,p), if(isconstant(c2), c2*psolution(b,m,c,c1*c3,y,t,x,p),0 ));

psolution(xx_,q_,0,xx_,y_,t_,x_,p_):= -x^q/q!;
psolution(tt_,q_,0,tt_*xx_,y_,t_,x_,p_):= if(hasnot(xx,t), -int(xx,x,q) );
psolution(tt_,q_,0,tt_*cc_,y_,t_,x_,p_):= if(hasnot(c,t), -int(cc,x,q) );
psolution(x_,1,0,t_*x_,y_,t_,x_,1):= x-t*x;

psolution(b_,m_,c_,c1_+c2_,y_,t_,x_,n_):= psolution(b,m,c,c1,y,t,x,n)+psolution(b,m,c,c2,y,t,x,n);


#7;
psolution1(b_,m_,c_,d_,y_,x_,n_):= if(hasnot(d,y),if(c==0,-int(d,x,m),
	if(hasnot(d/c,x),  if(m==1,-d/c*x, -d/c),
	if(hasnot(b/c,x) and hasnot(d/c/x,x),if(m==1, b/c^2*d/x-d/c, -d/c),
	if(b+c*x<>0 and hasnot(d/(b+c*x),x),-d*x/(b+c*x) 
	)))));
psolution1(b_,1,c_,d_,y_,x_,2):= if(c==0,-int(d,x),
	if(hasnot(d/c,x),  -d/c*x, 
	if(hasnot(b/c,x) and hasnot(d/c/x,x),  b/c^2*d/x-d/c, 
	if(b+c*x<>0 and hasnot(d/(b+c*x),x), -d*x/(b+c*x), 
	if(hasnot(b,c,y) and hasnot(b,c,x), 
	if(b*b+4c==0,  (exp(b/2*x)*x*int(exp(-b/2*x)*d,x)-exp(b/2*x)*int(d*exp(-b/2*x)*x,x)),
	((-exp((b-sqrt(b^2+4c))/2*x)*int(d*exp((-b+sqrt(b^2+4c))/2*x),x)+exp((b+sqrt(b^2+4c))/2*x)*int(d*exp((-b-sqrt(b^2+4c))/2*x),x))/sqrt(b*b+4c)) 
	))))));

#psolution(b_,m_,c_,c1_*c2_,y_,t_,n_):=if(hasnot(c1,t) and hasnot(c1,x), c1*psolution(b,m,c,c2,y,t,n), if(hasnot(c2,t) and hasnot(c2,x), c2*psolution(b,m,c,c1,y,t,n),
	if(m==1 and n==1 and has(c1,t) and hasnot(c1,x), if(c==0, int(replace(c1,t,t-x/b+z/b)*replace(c2,x,z),z,1,x)/b,
	int(exp(-c*z/b)*replace(c1,t,t-x/b+z/b)*replace(c2,x,z),z,1,x)*exp(c*x/b)/b ))));
psolution1(f_,p_,0,g_,y_,x_,q_):= int(psolution(f,g,y,x,q-p),x,p);


psolution(b_,1,c_/x_,f_,y_,x_,p_):= if(isconstant(b,c,f), if(b+c==0, f/b^2-f*log(x)*x/b, -f/(b+c)*x ));
psolution(b_,1,1/x_,f_,y_,x_,p_):= if(isconstant(b,f), if(b==-1, f+f*log(x)*x, -f/(b+1)*x ));
psolution(b_,1,c_*x_,f_,y_,x_,2):= if(isconstant(b,c), Bi((c*x+b*b/4)*c^(-2/3))*exp(b*x/2)*int(Ai((c*x+b*b/4)*c^(-2/3))*exp(-b*x/2)*f,x)*pi/cbrt(c)
	-Ai((c*x+b*b/4)*c^(-2/3))*exp(b*x/2)*int(Bi((c*x+b*b/4)*c^(-2/3))*exp(-b*x/2)*f,x)*pi/cbrt(c) );
psolution(b_,1,x_,f_,y_,x_,2):= if(isconstant(b),Bi(x+b*b/4)*exp(b*x/2)*int(Ai(x+b*b/4)*exp(-b*x/2)*f,x)*pi
	-Ai(x+b*b/4)*exp(b*x/2)*int(Bi(x+b*b/4)*exp(-b*x/2)*f,x)*pi );

psolution(b_/x_,1,c_*x_^(-2),d_,y_,x_,2):=if(isconstant(b,c,d),if(2-2b-c==0,d*log(x)*x^2, d/(2-2b-c)*x^2 ));
psolution(b_/x_,1,x_^(-2),d_,y_,x_,2):=if(isconstant(b,d),d/(1-2b)*x^2 );
psolution(1/x_,1,c_*x_^(-2),d_,y_,x_,2):=if(isconstant(c,d),-d/c*x^2 );
psolution(1/x_,1,x_^(-2),d_,y_,x_,2):=if(isconstant(d),-d*x^2 );
psolution(1/x_,1,x_^2,d_,y_,x_,2):=if(isconstant(d),1/2d*chi(x^2/2)*sinh(x^2/2)-1/2d*shi(x^2/2)*cosh(x^2/2) );
psolution(1/x_,1,x_,f_,y_,x_,2):= -AiPrime(x)*int(BiPrime(x)*f/x,x)*pi+BiPrime(x)*int(AiPrime(x)*f/x,x)*pi;

psolution(b_,m_,c_,c1_*c2_,y_,t_,n_):=if(hasnot(c1,t) and hasnot(c1,x), c1*psolution(b,m,c,c2,y,t,n), 
	if(hasnot(c2,t) and hasnot(c2,x), c2*psolution(b,m,c,c1,y,t,n),
	if(has(c1,t) and hasnot(c1,x), int(exp(-c*z/b)*replace(c1,t,t-x/b+z/b)*replace(c2,x,z),z,m,1,x)*exp(c*x/b)/b )));
psolution(b_,m_,c_,c1_*c2_*c3_,y_,x_,n_):=if(isconstant(c1), c1*psolution(b,m,c,c2*c3,y,x,n), if(isconstant(c2), c2*psolution(b,m,c,c1*c3,y,x,n) ));
psolution1(b_,m_,c_,c1_*c2_*c3_,y_,x_,n_):=if(isconstant(c1), c1*psolution(b,m,c,c2*c3,y,x,n), if(isconstant(c2), c2*psolution(b,m,c,c1*c3,y,x,n) ));

psolution(b_,m_,c_,c1_+c2_,y_,x_,n_):= psolution(b,m,c,c1,y,x,n)+psolution(b,m,c,c2,y,x,n);
psolution(b_,1,c_*x_^m_,f_+g_*x_^n_,y_,x_,p_):=if(isconstant(b,c,f,g) and n-m==1 and b*g/c==f, -g/c*x );
psolution(b_,1,x_^m_,f_+g_*x_^n_,y_,x_,p_):=if(isconstant(b,f,g) and n-m==1 and b*g==f, -g*x );
psolution(b_,1,c_*x_,f_+g_*x_^2,y_,x_,p_):=if(isconstant(b,c,f,g) and b*g/c==f, -g/c*x );
psolution(b_,1,x_,f_+g_*x_^2,y_,x_,p_):=if(isconstant(b,f,g) and b*g==f, -g*x );


psolution(b_,m_,c_,sinh(x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), psolution(b,m,c,exp(x)/2-exp(-x)/2,y,x,n));
psolution(b_,m_,c_,cosh(x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), psolution(b,m,c,exp(x)/2+exp(-x)/2,y,x,n));
psolution(b_,m_,c_,sin(x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), expand(totrig(psolution(b,m,c,i*exp(-i*x)/2-i*exp(i*x)/2,y,x,n))) );
psolution(b_,m_,c_,cos(x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), expand(totrig(psolution(b,m,c,exp(i*x)/2+exp(-i*x)/2,y,x,n))) );
psolution(b_,m_,c_,sinh(k_*x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,k,x), psolution(b,m,c,exp(k*x)/2-exp(-k*x)/2,y,x,n));
psolution(b_,m_,c_,cosh(k_*x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,k,x), psolution(b,m,c,exp(k*x)/2+exp(-k*x)/2,y,x,n));
psolution(b_,m_,c_,sin(k_*x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,k,x), expand(totrig(psolution(b,m,c,i*exp(-i*k*x)/2-i*exp(i*k*x)/2,y,x,n))) );
psolution(b_,m_,c_,cos(k_*x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,k,x), expand(totrig(psolution(b,m,c,exp(i*k*x)/2+exp(-i*k*x)/2,y,x,n))) );
psolution(b_,m_,c_,sin(f_+x_),y_,x_,n_):=if(isfree(b,c,f,y) and isfree(b,c,f,x), expand(totrig(psolution(b,m,c,i*exp(f-i*x)/2-i*exp(i*x+f)/2,y,x,n))) );
psolution(b_,m_,c_,cos(f_+x_),y_,x_,n_):=if(isfree(b,c,f,y) and isfree(b,c,f,x), expand(totrig(psolution(b,m,c,exp(i*x+f)/2+exp(f-i*x)/2,y,x,n))) );
psolution(b_,m_,c_,sin(f_+k_*x_),y_,x_,n_):=if(isfree(b,c,k,f,y) and isfree(b,c,k,f,x), expand(totrig(psolution(b,m,c,i*exp(f-i*k*x)/2-i*exp(i*k*x+f)/2,y,x,n))) );
psolution(b_,m_,c_,cos(f_+k_*x_),y_,x_,n_):=if(isfree(b,c,k,f,y) and isfree(b,c,k,f,x), expand(totrig(psolution(b,m,c,exp(i*k*x+f)/2+exp(f-i*k*x)/2,y,x,n))) );
#psolution(b_,m_,c_,exp(x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), if(n==b*m, 1/n!*exp(x)*x^n, if(1-b-c==0, exp(x)*x/(n-b*m), exp(x)/(1-b-c) )));
#psolution(b_,m_,c_,exp(a_*x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), if(a^n-b*a^m==c, exp(a*x)*x/(n*a^(n-1)-b*m*a^(m-1)), exp(a*x)/(a^n-b*a^m-c) ));
psolution(b_,m_,c_,exp(a_*x_),y_,x_,n_):=if(isfree(a,b,c,y),if(isfree(a,b,c,x), if(a^n==b*a^m+c, if(b*b== -4c, a/n/m*exp(a*x)*x^2, exp(a*x)*x/(n*a^(n-1)-b*m*a^(m-1)) ), exp(a*x)/(a^n-b*a^m-c) ), 
	if(b==c and a==-1 and isodd(m) and iseven(n),exp(-x) )));
#psolution(b_,m_,c_,exp(x_+f_),y_,x_,n_):=if(isfree(b,c,f,y) and isfree(b,c,f,x), if(b+c==1, if(b*b== -4c, 1/n/m*exp(x+f)*x^2, exp(x+f)*x/(n-b*m) ), exp(x+f)/(1-b-c) ));
psolution(b_,m_,c_,exp(f_+x_),y_,x_,n_):=if(isfree(b,c,f,y),if(isfree(b,c,x), if(b+c==1, if(b*b== -4c and 2m==n, 1/m/n*exp(x+f)*x^2, exp(x+f)*x/(n-b*m) ), exp(x+f)/(1-b-c)), if(b+c==0,exp(f+x) )));
psolution(b_,m_,c_,exp(f_+a_*x_),y_,x_,n_):=if(isfree(a,b,c,f,y),if(isfree(a,b,c,f,x), if(a^n==b*a^m+c, if(b*b== -4c, a/n/m*exp(a*x+f)*x^2, exp(a*x+f)*x/(n*a^(n-1)-b*m*a^(m-1)) ), exp(a*x+f)/(a^n-b*a^m-c) ), 
	if(b==c and a==-1 and isodd(m) and iseven(n),exp(f-x) )));
psolution(b_,m_,c_,exp(x_),y_,x_,n_):=if(isfree(b,c,y),if(b+c==0,exp(x),if(isfree(b,c,x), if(b+c==1, if(b*b== -4c and 2m==n, 1/m/n*exp(x)*x^2, exp(x)*x/(n-b*m) ), exp(x)/(1-b-c) ))));

#psolution(b_,m_,c_,x_^n_,y_,x_,p_):=if(c<>0 and isfree(b,c,y) and isfree(c,x), if(m==n, -x^n/c+b/c^2, if(m>n, -x^n/c, psolution(c,x^n,y,x,p)-psolution(c/b,x^n,y,x,m)-b*x^n/c^2 )));
#psolution(b_,m_,c_,x_^n_,y_,x_,p_):=if(isfree(b,c,y) and isfree(c,x),if(m>n and isinteger(m-n), -x^n/c, if(isinteger(n) or isinteger(p-n) and m==n, -x^n/c+b/c^2*n! )));
psolution(b_,m_,c_,x_^n_,y_,x_,p_):=if(n>0 and isfree(b,c,y) and isfree(c,x),if(m>n and isinteger(m-n),-x^n/c, if(p>=n,-psolution(-c/b,x^n/b,y,x,m)-if(p==n,n!/c^2,0) )));
psolution(b_,m_,c_,x_,y_,x_,p_):=if(isfree(b,c,y) and isfree(c,x),if(m>1,-x/c, if(m==1,-x/c+b/c^2, b*x^m/m!/c^2-x/c-b^2/c^3 )));

psolution(b_,1,c_,sin(x_),y_,x_,2):=if(isfree(b,c,y) and isfree(b,c,x), - (c+1)*sin(x)/(b^2 + (c + 1)^2) + b*cos(x)/(b^2 + (c + 1)^2) );
psolution(b_,1,c_,cos(x_),y_,x_,2):=if(isfree(b,c,y) and isfree(b,c,x), - b*sin(x)/(b^2 + (c + 1)^2) - (c+1)*cos(x)/(b^2 + (c + 1)^2) );
psolution(b_,1,c_,sin(k_*x_),y_,x_,2):=if(b==0 and c== -k*k, sgn(k)*cos(k*x)*x/c, if(isfree(b,c,k,y) and isfree(b,c,k,x), b*k*cos(k*x)/(b^2*k^2 + c^2 + 2c*k^2+k^4) - (c+k^2)*sin(k*x)/(b^2*k^2 + c^2 + 2c*k^2+k^4)  ));
psolution(b_,1,c_,cos(k_*x_),y_,x_,2):=if(b==0 and c== -k*k, -sin(k*x)*x/c, if(isfree(b,c,k,y) and isfree(b,c,k,x), -b*k*sin(k*x)/(b^2*k^2 + c^2 + 2c*k^2+k^4) - (c+k^2)*cos(k*x)/(b^2*k^2 + c^2 + 2c*k^2+k^4)  ));
psolution(b_,1,c_,sin(k_+x_),y_,x_,2):=if(isfree(b,c,k,y) and isfree(b,c,k,x), - (c+1)*sin(x+k)/(b^2 + (c + 1)^2) + b*cos(x+k)/(b^2 + (c + 1)^2) );
psolution(b_,1,c_,cos(k_+x_),y_,x_,2):=if(isfree(b,c,k,y) and isfree(b,c,k,x), - b*sin(x+k)/(b^2 + (c + 1)^2) - (c+1)*cos(x+k)/(b^2 + (c + 1)^2) );
psolution(b_,1,c_,sin(f_+k_*x_),y_,x_,2):=if(b==0 and c== -k*k, sgn(k)*cos(f+k*x)*x/c, if(isfree(b,c,k,f,y) and isfree(b,c,k,f,x),  b*k*cos(k*x+f)/(b^2*k^2 + c^2 + 2c*k^2+k^4) - (c+k^2)*sin(k*x+f)/(b^2*k^2 + c^2 + 2c*k^2+k^4)  ));
psolution(b_,1,c_,cos(f_+k_*x_),y_,x_,2):=if(b==0 and c== -k*k, -sin(f+k*x)*x/c, if(isfree(b,c,k,f,y) and isfree(b,c,k,f,x), -b*k*sin(k*x+f)/(b^2*k^2 + c^2 + 2c*k^2+k^4) - (c+k^2)*cos(k*x+f)/(b^2*k^2 + c^2 + 2c*k^2+k^4)  ));

psolution(b_,k_,x_^m_,x_^n_,y_,x_,p_):= if(k>n-m and n-m>0,-x^(n-m), if(n-m==1,-psolution(b,1,x^m,b,y,x,p)-x, if(p>n-m,-psolution(-x^m/b,x^n/b,y,x,k) )));
psolution(b_,k_,x_^m_,x_,y_,x_,p_):= if(k>1-m and 1-m>0,-x^(1-m), if(p>1-m,-psolution(-x^m/b,x/b,y,x,k) ));
psolution(b_,k_,x_,x_^n_,y_,x_,p_):= if(k>n-1 and n>1,-x^(n-1), if(p>n-1,-psolution(-x/b,x^n/b,y,x,k) ));
psolution(b_,1,x_,x_^2,y_,x_,p_):=  if(b==x^2, -x/2, -psolution(b,1,x,b,y,x,p)-x );
psolution(b_,1,x_,x_^3,y_,x_,2):= if(isconstant(b),psolution(b,1,x,2,y,x,2)+2b-x^2 );
psolution(b_,1,x_,x_^4,y_,x_,2):= if(isconstant(b),psolution(b,1,x,3b*b,y,x,2)-6+3b*x-x^3 );


psolution(x_,1,c_,x_,y_,x_,p_) := if(isconstant(c) and c<> -1, -x/(c+1) );
psolution(x_,1,c_,x_^2,y_,x_,2):= if(isconstant(c) and c<> -2, -2/((c+2)*c)-x^2/(c+2) );
psolution(x_,1,-2,x_^2,y_,x_,2):= 1-x^2/2;
psolution(x_,1,1,x_^3,y_,x_,2):= -(3x+x^3)/4;
psolution(x_,1,1,x_^4,y_,x_,2):= -(8+4x^2+x^4)/5;
psolution(x_,1,1,x_^5,y_,x_,2):= -(15x+5x^3+x^5)/6;
#psolution(x_,1,1,x_^n_,y_,x_,2):= -psolution(1/x,x^(n-1),y,x,2)/(1+n);

psolution(b_,m_,c_,bb_*c_,y_,x_,n_):=if(isfree(bb,c,y) and isfree(bb,x),  -bb);
psolution(b_,m_,a_*c_,c_,y_,x_,n_):=if(isfree(a,c,y) and isfree(a,x),  -1/a);
psolution(b_,m_,a_*c_,bb_*c_,y_,x_,n_):=if(isfree(a,bb,c,y) and isfree(a,bb,x), -bb/a);
psolution(b_,mm_,c_,c_*x_^m_,y_,x_,p_):=if(isinteger(m) or isinteger(p-m) and p>m and m>0,if(mm>=m,-x^m+if(mm==m,m!, 0), -sum(sgn(-b/c)^k*m!/(m-k)!*x^(m-k),k,0,m,mm) ));
psolution(b_,m_,c_,c_*x_,y_,x_,n_):= if(m>1,-x, if(m==1,-x+b/c ));
psolution(b_,m_,c_,c_,y_,x_,n_):=  if(m==1,-x,if(m>0, -1 ));
psolution(b_,m_,c_,-c_,y_,x_,n_):= if(m==1,x, if(m>0, 1 ));
psolution(b_,m_,-c_,c_,y_,x_,n_):= if(m==1,x, if(m>0, 1 ));


psolution(b_,1,x_,x_^3,y_,x_,3):= if(isconstant(b),2b-x^2);

psolution(b_,1,c_/x_,x_,y_,x_,2):= if(isconstant(b,c), if(b+c==0,-x^2/b+2x/b^2+2/b^3-2x*log(x)/b^2, -1/(2b+c)*x^2-2*x/(6*(b+c)) ));
psolution(b_,1,b_/x_,x_,y_,x_,2):= if(isconstant(b), -x^2/(3b)-x/(3b*b) );
psolution(b_,1,1/x_,x_^n_,y_,x_,p_):=if(p>n+1 and isconstant(b), -x^(n+1)/(b*n+b+1) );
psolution(b_,1,1/x_,x_,y_,x_,p_):= if(isconstant(b) and p>2,-x^2/(2b+1) );
psolution(b_,1,1/x_,x_,y_,x_,2):= if(isconstant(b),-x^2/(2b+1)-b/(3+6b)*x );

psolution(b_/x_,1,x_^(-2),log(x_)*x_^(-2),y_,x_,2):=if(isconstant(b),(b+1)-log(x) );
psolution(b_/x_,1,c_*x_^(-2),log(x_)*x_^(-2),y_,x_,2):=if(isconstant(b,c),(b+1)/c^2-1/c*log(x) );
psolution(b_/x_,1,c_*x_^(-2),x_^n_,y_,x_,2):=if(isconstant(b,c), x^(2+n)/((n+2)*(-b+(n+2)-1)-c) );
psolution(b_/x_,1,c_*x_^(-2),x_,y_,x_,2):=if(isconstant(b,c), x^3/(3*(-b+2)-c) );
psolution(b_/x_,1,c_*x_^(-2),x_^(-3),y_,x_,2):=if(isconstant(b,c),log(x)/x/(c-b-1) );

psolution(2,1,c_,exp(x_)*x_,y_,x_,2):= if(c+1==0, exp(x)*x^3/6,-exp(x)*x/(1+c) );
psolution(-2,1,c_,exp(-x_)*x_,y_,x_,2):= if(c+1==0, exp(-x)*x^3/6,-exp(-x)*x/(1+c) );
psolution(1,1,c_/x_,x_,y_,x_,2):= if(isconstant(c), -1/(2+c)*x^2-c*x/(6*(1+c)) );
psolution(2,1,1/x_,x_,y_,x_,2):= -1/5*x^2-2/15*x ;
psolution(1,1,-1/x_,x_,y_,x_,2):= -x^2+2x+2-2x*log(x) ;
psolution(1,1,1/x_,x_,y_,x_,2):= -1/3*x^2-1/3*x ;
psolution(1,1,2/x_,x_,y_,x_,2):= -1/4*x^2-1/6*x ;
psolution(1,1,3/x_,x_,y_,x_,2):= -1/5*x^2-1/10*x ;
psolution(1,1,4/x_,x_,y_,x_,2):= -1/6*x^2-1/15*x ;
psolution(1,1,5/x_,x_,y_,x_,2):= -1/7*x^2-1/21*x ;
psolution(2,1,-1,exp(x_)*xx_,y_,x_,2):= exp(x)*int(xx,x,2);
psolution(-2,1,-1,exp(-x_)*xx_,y_,x_,2):= exp(-x)*int(xx,x,2);

psolution(x_^m_,1,x_^m_,x_^n_,y_,x_,p_):= if(n-m==1,1-x, if(p>n, -psolution(-1,x^(n-m),y,x,1) ));
#psolution(x_,1,x_,x_^n_,y_,x_,p_):= if(n==2,1-x, if(p>=n,0, -psolution(x,1,x,(-1)^(n+p)*(n-1)!,y,x,p))-psolution(-1,x^(n-1),y,x,1) );
psolution(x_,1,x_,x_^n_,y_,x_,p_):= if(n==2,1-x, if(p>=n,0, -psolution(x,1,x,(-1)^(n+p)*(n-1)!,y,x,p))-sum((-1)^k*(n-1)!/(n-k-1)!*x^(n-k-1),k,0,n-1) );

psolution(1/x_,1,1,x_^2,y_,x_,2):= -x^2;
psolution(1/x_,1,x_^(-2),x_^n_,y_,x_,2):=x^(2+n)/((n+2)^2-2*(n+2)-1);
psolution(1/x_,1,x_^(-2),x_,y_,x_,2):=x^3/2;
psolution(1/x_,1,1/x_,x_^n_,y_,x_,2):= if(n>0,-sum((n+1)!/(n+1-k)!*(n-1)!/(n-1-k)!*x^(n+1-k),k,0,n-1) );
psolution(1/x_,1,1/x_,x_,y_,x_,2):= -x^2;
psolution(1/x_,1,x_,x_^2,y_,x_,2):= -psolution(1/x,1,x,1/x,y,x,2)-x;
psolution(1/x_,1,x_,x_^3,y_,x_,2):= -x^2;
psolution(1/x_,1,x_,x_^4,y_,x_,2):= -3-x^3;
psolution(1/x_,1,x_,x_^6,y_,x_,2):= -15x^2-x^5;
psolution(1/x_,1,x_,x_^7,y_,x_,2):= -x^6-24x^3-72;
psolution(1/x_,1,x_,x_^9,y_,x_,2):= -x^8-48x^5-720x^2;

psolution(x_,1,x_,x_^5,y_,x_,3):= -psolution(x,1,x,24,y,x,3)-48+24*x-12*x^2+4*x^3-x^4;
psolution(x_,1,x_,x_^4,y_,x_,3):= -psolution(x,1,x,-6,y,x,3)+6-6*x+3*x^2-x^3;
psolution(x_/a_,1,1/a_,x_/a_,y_,x_,2):= -x/2;
psolution(x_,1,x_,x_^5,y_,x_,2):= -psolution(x,1,x,-24,y,x,2)+12+12x-12*x^2+4*x^3-x^4;
psolution(x_,1,x_,x_^4,y_,x_,2):= -psolution(x,1,x,6,y,x,2)-6*x+3*x^2-x^3;
psolution(x_,1,x_,x_^3,y_,x_,2):= -psolution(x,1,x,-2,y,x,2)-1+2*x-x^2;
psolution(x_,1,x_,x_^2,y_,x_,p_):= 1-x;

psolution(b_,1,0,x_^2,y_,x_,n_):= if(isconstant(b) and n>2, -x/b-x^3/b/3 );
psolution(b_,1,0,x_,y_,x_,p_):= if(isconstant(b) and p>2, -1/b/2*x^2);
psolution(b_,1,0,x_^2,y_,x_,2):= if(isconstant(b), -2x/b-x^2/b-x^3/b/3 );
psolution(b_,1,0,x_,y_,x_,2):= if(isconstant(b), -1/b/2*x^2-x/b^2 );
psolution(b_,1,0,exp(b_*x_),y_,x_,2):= exp(b*x)*x/b;
psolution(b_/x_,1,0,x_^(-2),y_,x_,2):= if(b==-1,log(x)^2/2,log(x)/(-1-b) );
psolution(b_,1,0,b_,y_,x_,p_):= -x;
psolution(b_,1,0,-b_,y_,x_,p_):= x;
psolution(-b_,1,0,b_,y_,x_,p_):= x;
psolution(1/x_,1,0,x_^(-2),y_,x_,2):= -log(x)/2;
psolution(1,1,0,sin(x_),y_,x_,2):= cos(x)/2-sin(x)/2;
psolution(1,1,0,cos(x_),y_,x_,2):= -cos(x)/2-sin(x)/2;
psolution(1,1,0,exp(x_),y_,x_,2):= exp(x)*x;


#5;
#psolution1(c_,d_,y_,x_,n_):=if(isfree(c,d,y), if(isfree(c,d,x), -d/c, if(n==2, if(isfree(c,x), 
	if(c<0, expand((cos(sqrt(-c)*x)*int(d*sin(sqrt(-c)*x),x)-sin(sqrt(-c)*x)*int(d*cos(sqrt(-c)*x),x))/2/sqrt(-c)),
		expand((exp(sqrt(c)*x)*int(d*exp(-sqrt(c)*x),x)-exp(-sqrt(c)*x)*int(d*exp(sqrt(c)*x),x))/2/sqrt(c)) )),
	if(isinteger(n), block(dsolve:=integrate(c,x), expand(exp(dsolve)*integrate(d*exp(-dsolve),x))),
	 int((x-z)^(n-1)*replace(d,x,z)*mittag(n,n,replace(int(c,x,n),x,x-z)*n!), z,0,x) ))));
#psolution1(c_,d_,y_,x_,n_):=if(isfree(c,d,y) and n<1, block(gs:=gsolution(c,y,x,n)/C_1, int(simplify(d/gs),x,n)*gsolution(c,y,x,n)/C_1) );
#psolution1(c_,d_,y_,x_,n_):=if(isfree(c,d,y) and n<=1, block(gs:=gsolution(c,y,x,n)/C_1, int(simplify(d/gs),x,n)*gs) );
#psolution(c_,d_,y_,x_,2):=if(isfree(c,d,y) and isfree(c,x),if(c<0, cos(sqrt(-c)*x)*int(d*sin(sqrt(-c)*x),x)/2/sqrt(-c)-sin(sqrt(-c)*x)*int(d*cos(sqrt(-c)*x),x)/2/sqrt(-c),
	expand(exp(sqrt(c)*x)*int(d*exp(-sqrt(c)*x),x)/2/sqrt(c)-exp(-sqrt(c)*x)*int(d*exp(sqrt(c)*x),x)/2/sqrt(c)) ));

psolution1(c_,d_,y_,x_,n_):=if(isfree(d/c,x), -d/c, 
	if(n<=1,block(dsolve:=n!*integrate(c,x,n), expand(integrate(d/mittag(n,dsolve),x,n)*mittag(n,dsolve))) ));

#psolution1(b_,x_^n_, y_,x_,q_) := If(isfree(b,x) and isfree(b,y) and n>0, if(q==n,-x^n/b-n!, 
	if(q-n>0,-x^n/b, 
	if(isinteger(n) or isinteger(n-q),-sum(b^(-1-k/q)*n!/(n-k)!*x^(n-k),k,0,n,q),
	-gamma(n+1,b*x)*exp(b*x) ))));


psolution1(c_,g_,y_,x_,-1):=if(hasnot(g,c,x),-g/c*erf(x)*exp(x^2), if(hasnot(c,g,E(0.5,x)) and hasnot(c,x), block(f:= -psolution(1/c,d(g,x)/c,y,x,1),
	if(hasnot(f,int(y(x),x)),f) )));
psolution1(c_,g_,y_,x_,-0.5):=if(hasnot(g,c,x),-g/c*erf(x^0.5)*exp(x), if(hasnot(g,f(x)) and hasnot(c,g,E(0.5,x)) and hasnot(c,x), 
	block(f:= -psolution(1/c^2,d(g,x)/c+d(g,x,0.5)/c^2,y,x,1),
	if(hasnot(f,int(y(x),x)),f) )));
psolution(c_,g_,y_,x_,0.5):=if(hasnot(g/c,x),-g/c, if(hasnot(g,f(x)) and hasnot(c,g,E(0.5,x)) and hasnot(c,x), 
	block(f:=psolution(c*c,d(g,x,0.5)+c*g,y,x,1),
	if(hasnot(f,int(y(x),x)),f)) ));
#psolution1(c_,g_,y_,x_,1/2):=if(hasnot(g/c,x),-g/c, if(hasnot(c,g,E(0.5,x)), psolution(c*c+d(c,x,0.5),d(g,x,0.5)+c*g,y,x,1) ));
psolution1(c_,d_,y_,x_,4):=if(isfree(c,x), gsolution1(c,y,x,4)*int(d*gsolution2(c,y,x,4),x)/(4*(c)^(3/4))
	-gsolution2(c,y,x,4)*int(d*gsolution1(c,y,x,4),x)/(4*(c)^(3/4))+gsolution3(c,y,x,4)*int(d*gsolution4(c,y,x,4),x)/(2*(c)^(3/4))
	-gsolution4(c,y,x,4)*int(d*gsolution3(c,y,x,4),x)/(2*(c)^(3/4))  );
psolution1(c_,d_,y_,x_,3):=if(isfree(c,x), (gsolution1(c,y,x,3)*int(d*gsolution2(c,y,x,3),x)+gsolution2(c,y,x,3)*int(d*gsolution3(c,y,x,3),x)
	+gsolution3(c,y,x,3)*int(d*gsolution1(c,y,x,3),x))/(3cbrt(3)*c^(2/3)) );
psolution1(c_,d_,y_,x_,2):=if(isfree(c,x),if(isfree(d,x), -d/c,  
	if(c<0, expand((cos(sqrt(-c)*x)*int(d*sin(sqrt(-c)*x),x)-sin(sqrt(-c)*x)*int(d*cos(sqrt(-c)*x),x))/sqrt(-c)),
	expand((exp(sqrt(c)*x)*int(d*exp(-sqrt(c)*x),x)-exp(-sqrt(c)*x)*int(d*exp(sqrt(c)*x),x))/2/sqrt(c)) )));

psolution1(c_,g_,y_,x_,1):=if(hasnot(g/c,x),-g/c, if(isfree(c,g,y), 
	block(dsolve:=integrate(c,x), exp(dsolve)*(integrate(g*exp(-dsolve),x))) ));
#psolution(b_*x_,c_,y_,x_,1) := if(isconstant(c),if(b>0, c*Dawsonminus(sqrt(b/2)*x)/sqrt(b/2), if(b<0, c*Dawson(sqrt(-b/2)*x)/sqrt(-b/2) )));
psolution(b_*x_,c_,y_,x_,1) := if(isconstant(c) and b<0, c*Dawson(sqrt(-b/2)*x)/sqrt(-b/2) );

psolution1(exp(x_),d_,y_,x_,2):= 2 besselI(0,2 sqrt(exp(x)))* int(besselK(0,2 sqrt(exp(x)))* d,x) -2 besselK(0,2 sqrt(exp(x))) *int(besselI(0,2 sqrt(exp(x)))* d,x);
psolution1(1/x_,d_,y_,x_,2):= 2sqrt(x)*besselI(1,2sqrt(x))*int(besselK(1,2sqrt(x))*sqrt(x)*d,x)-2sqrt(x)*besselK(1,2sqrt(x))*int(besselI(1,2sqrt(x))*sqrt(x)*d,x);
psolution1(x_^(-2),d_,y_,x_,2):= x^(1/2 + sqrt(5)/2)* int(x^(1/2 - sqrt(5)/2) *d,x)/sqrt(5) - x^(1/2 - sqrt(5)/2)* int(x^(1/2 + sqrt(5)/2)*d,x)/sqrt(5);
psolution1(b_*x_^(-2),d_,y_,x_,2):=  if(isconstant(b),-i*x^(1/2 + 1/2 i *sqrt(-4 b - 1))* int(( x^(1/2 - 1/2 i* sqrt(-4 b - 1))* f(x)),x)/sqrt(-4 b - 1) 
	+ i*x^(1/2 - 1/2 i* sqrt(-4 b - 1))* int( (x^(1/2 i* sqrt(-4 b - 1) + 1/2)* f(x)),x)/sqrt(-4 b - 1) );
psolution1(c_*x_,f_,y_,x_,2):= if(isconstant(c),pi*Bi(cbrt(c)*x)*int(Ai(cbrt(c)*x)*f,x)/cbrt(c)-pi*Ai(cbrt(c)*x)*int(Bi(cbrt(c)*x)*f,x)/cbrt(c));
psolution1(x_,f_,y_,x_,2):= pi*Bi(x)*int(Ai(x)*f,x)-pi*Ai(x)*int(Bi(x)*f,x);


psolution(b_*x_^m_,c_,y_,x_,q_) := if(isconstant(b,c),if(m+q==0,
	if(b==q!,c*log(x)*x^q/(fallingfactorial(4,q-1)-1),c/x^m/(q!-b)), 
	if(m<0 and re(m+q)>0,-c/b/x^m )));
psolution(b_/x_,c_,y_,x_,n_):= if(re(n)>1 and isconstant(b,c),-c/b*x );
psolution(x_^m_,c_,y_,x_,q_) := if(isconstant(c),if(m+q==0,c/x^m/(q!-1), if(m<0 and re(m+q)>0,-c/x^m )));
psolution(1/x_,c_,y_,x_,n_):= if(re(n)>1 and isconstant(c),-c*x );
psolution(f_/(x_+b_),c_,y_,x_,n_):= if(re(n)>1 and isconstant(b,c,f),-c*x/f-b*c/f );
psolution(f_/(a_*x_+b_),c_,y_,x_,n_):= if(re(n)>1 and isconstant(a,b,c,f),-a*c*x/f-b*c/f );
psolution(1/(x_+b_),c_,y_,x_,n_):= if(re(n)>1 and isconstant(b,c),-c*x-b*c );
psolution(1/(a_*x_+b_),c_,y_,x_,n_):= if(re(n)>1 and isconstant(a,b,c),-a*c*x-b*c );

psolution(2*x_^(-2),c_,y_,x_,2) := if(isconstant(c),c*log(x)*x^2/3 );
psolution(6*x_^(-3),c_,y_,x_,3) := if(isconstant(c),c*log(x)*x^3/11 );
psolution(24*x_^(-4),c_,y_,x_,4):= if(isconstant(c),c*log(x)*x^4/50 );

psolution(c_,c1_*c2_,y_,x_,n_):=if(isfree(c1,x),c1*psolution(c,c2,y,x,n), 
	if(isfree(c2,x),c2*psolution(c,c1,y,x,n) ));
psolution(c_,c1_*c2_*c3_,y_,x_,n_):=if(isfree(c1,x),c1*psolution(c,c2*c3,y,x,n), 
	if(isfree(c2,x),c2*psolution(c,c1*c3,y,x,n), 
	if(isfree(c3,x),c3*psolution(c,c1*c2,y,x,n) 
	)));


psolution(c_,c1_+c2_,y_,x_,p_):=psolution(c,c1,y,x,p)+psolution(c,c2,y,x,p);

psolution(a_*c_,c_,y_,x_,p_):=if(isfree(a,x),if(p>0,  -1/a, 1/(1-a) ));
psolution(c_,b_*c_,y_,x_,p_):=if(isfree(b,x),if(p>0,  -b, 1-b ));
psolution(a_*c_,b_*c_,y_,x_,p_):=if(isfree(a,b,x), if(p>0,-b/a, b/(1-a) ));
psolution(a_*z_,b_*x_*z_,y_,x_,p_):=if(isfree(a,b,x), if(p>0,-b/a*x, b/(1-a) ));

psolution(c_,E(n_,x_^n_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x), if(c==1, 2/sqrt(pi)*E(n,x^n)*x^n, E(n,x^n)/(1-c) ));
psolution(c_,E(n_,b_*x_^n_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), if(c==b, 2/sqrt(pi)*E(n,b*x^n)*x^n,E(n,b*x^n)/(b-c) ));
psolution(c_,sinh(x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x),psolution(c,exp(x)/2-exp(-x)/2,y,x,n));
psolution(c_,cosh(x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x),psolution(c,exp(x)/2+exp(-x)/2,y,x,n));
#psolution(c_,sin(x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,i*exp(-i*x)/2-i*exp(i*x)/2,y,x,n))) );
#psolution(c_,cos(x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,exp(i*x)/2+exp(-i*x)/2,y,x,n))) );
psolution(c_,sinh(k_*x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x),psolution(c,exp(k*x)/2-exp(-k*x)/2,y,x,n));
psolution(c_,cosh(k_*x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x),psolution(c,exp(k*x)/2+exp(-k*x)/2,y,x,n));
#psolution(c_,sin(k_*x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,i*exp(-i*k*x)/2-i*exp(i*k*x)/2,y,x,n))) );
#psolution(c_,cos(k_*x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,exp(i*k*x)/2+exp(-i*k*x)/2,y,x,n))) );
#psolution(c_,sin(f_+x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,i*exp(-i*x-f*i)/2-i*exp(i*x+f*i)/2,y,x,n))) );
#psolution(c_,cos(f_+x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,exp(i*x+f*i)/2+exp(-i*x-f*i)/2,y,x,n))) );
#psolution(c_,sin(f_+k_*x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,i*exp(-i*k*x-f*i)/2-i*exp(i*k*x+f*i)/2,y,x,n))) );
#psolution(c_,cos(f_+k_*x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x) and n<>1,expand(totrig(psolution(c,exp(i*k*x+f*i)/2+exp(-i*k*x-f*i)/2,y,x,n))) );
psolution(c_,exp(a_*x_),y_,x_,n_):=if(isfree(a,c,y) and isfree(c,x), if(a^n==c,exp(a*x)*x/n/a^(n-1), exp(a*x)/(a^n-c) ));
psolution(c_,exp(a_+x_),y_,x_,n_):=if(isfree(a,c,y) and isfree(a,c,x), if(c==1,exp(a+x)*x/n, exp(a+x)/(1-c) ));
psolution(c_,exp(b_+a_*x_),y_,x_,n_):=if(isfree(b,c,y) and isfree(b,c,x), if(a^n==c,exp(a*x+b)*x/n/a^(n-1), exp(a*x+b)/(a^n-c) ));
psolution(c_,exp(x_)*x_,y_,x_,n_):=if(isfree(c,y) and isfree(c,x), if(c==1,exp(x)*x^2/2/n-(n-1)*exp(x)*x/2/n, -n*exp(x)/(c-1)^2-exp(x)*x/(c-1) ));
psolution(c_,exp(a_*x_)*x_,y_,x_,n_):=if(isfree(a,c,y) and isfree(a,c,x), if(a^n==c,exp(a*x)*x^2/(2n*a^(n-1))-(n-1)*exp(a*x)*x/(2n*a^n), -a^(n-1)*n*exp(a*x)/(c-a^n)^2-exp(a*x)*x/(c-a^n) ));
psolution(c_,exp(d_+x_)*x_,y_,x_,n_):=if(isfree(c,d,y) and isfree(c,d,x), if(c==1,exp(d+x)*x^2/2/n-(n-1)*exp(d+x)*x/2/n, -n*exp(d+x)/(c-1)^2-exp(d+x)*x/(c-1) ));
psolution(c_,exp(b_+a_*x_)*x_,y_,x_,n_):=if(isfree(a,b,c,y) and isfree(a,b,c,x), if(a^n==c,exp(a*x+b)*x^2/(2n*a^(n-1))-(n-1)*exp(a*x+b)*x/(2n*a^n), -a^(n-1)*n*exp(a*x+b)/(c-a^n)^2-exp(a*x+b)*x/(c-a^n) ));
psolution(c_,exp(x_^2),y_,x_,n_):=if(isfree(c,y) and isfree(c,x),1/2*erfi(1-n,-c/2+x)*exp(-c^2/4+c*x)*pi^(1/2) );
psolution(c_,exp(x_)/x_,y_,x_,n_):=if(isfree(c,y) and isfree(c,x), if(c==1,exp(x)*ln(1-n,x),exp(c/n*x)*Ei(1-n,(1-c)*x) ));
psolution(c_,exp(x_),y_,x_,n_):=if(isfree(c,y) and isfree(c,x), if(c==1,exp(x)*x/n, exp(x)/(1-c) ));

psolution(c_,sin(x_),y_,x_,2):= if(c== -1,-cos(x)*x/2, if(isfree(c,x), -sin(x)/(c + 1) ));
psolution(c_,cos(x_),y_,x_,2):= if(c== -1,sin(x)*x/2, if(isfree(c,x), -cos(x)/(c + 1) ));
psolution(c_,tan(x_),y_,x_,2):= if(c== -1,-cos(x)*log(sec(x)+tan(x)) );
psolution(c_,cot(x_),y_,x_,2):= if(c== -1,sin(x)*log(sec(x)+tan(x)) );
psolution(c_,sin(d_*x_),y_,x_,2):=if(c== -d*d, sgn(d)*cos(d*x)*x/c, if(isfree(c,d,x), -sin(d*x)/(c+d^2) ));
psolution(c_,cos(d_*x_),y_,x_,2):=if(c== -d*d,-sin(d*x)*x/c, if(isfree(c,d,x), -cos(d*x)/(c+d^2) ));
psolution(c_,sin(b_+x_),y_,x_,2):=if(isfree(b,c,x), if(c== -1,-cos(x+b)*x/2, -sin(x+b)/(c + 1) ));
psolution(c_,cos(b_+x_),y_,x_,2):=if(isfree(b,c,x), if(c== -1,sin(x+b)*x/2, -cos(x+b)/(c + 1) ));
psolution(c_,sin(f_+d_*x_),y_,x_,2):=if(c== -d*d,sgn(d)*cos(f+d*x)*x/c, if(isfree(c,d,f,x), -sin(d*x+f)/(c+d^2) ));
psolution(c_,cos(f_+d_*x_),y_,x_,2):=if(c== -d*d,-sin(f+d*x)*x/c, if(isfree(c,d,f,x), -cos(d*x+f)/(c+d^2) ));
psolution(c_,sin(x_)*x_,y_,x_,2):= if(c== -1,-cos(x)*x^2/4+sin(x)*x/4, if(isfree(c,x), c/2*cos(x)/(c+1)^2-sin(x)*x/(c + 1) ));
psolution(c_,cos(x_)*x_,y_,x_,2):= if(c== -1,sin(x)*x^2/4+cos(x)*x/4, if(isfree(c,x), -c/2*sin(x)/(c+1)^2-cos(x)*x/(c + 1) ));

psolution(c_,x_^n_, y_,x_,q_) := If(isfree(c,n,x) and isfree(c,y) and re(n)>0,if(re(q-n)>0,-x^n/c, if(q==n,-x^n/c-n!/c^2, if(isinteger(2*(n-q)),-sum((c)^round(-1-k/q)*n!/(n-k)!*x^(n-k),k,0,n,q) ))));
psolution(c_,x_,y_,x_,n_):= if(isfree(c,y) and isfree(c,x),  if(re(n)>1,-x/c, if(n==1,-x/c-1/c^2, if(n==-0.5,1+x^(1+n)/(1+n)!, if(n==0.5,-x/c-x^n/n!/c^2-1/c^3 )))));
psolution(c_,c_*x_^m_,y_,x_,p_):=if(isinteger(m) or isinteger(p-m) and p>m and m>0,-x^m, if(p==m and isconstant(c),-x^m-p!/c ));
psolution(c_,c_*x_,y_,x_,p_):=if(re(p)>1,-x, if(re(p)==1,-x-1/c ));
psolution(a_*c_,c_*x_^m_,y_,x_,p_):=if(isinteger(m) or isinteger(p-m) and p>m and m>0 isconstant(a),-x^m/a);

psolution(b_*x_,c_*x_^n_,y_,x_,p_):=if(isconstant(c/b),if(re(p)-n+1>0 and isinteger(n), -c/b*x^(n-1),if(p==n-1,-c/b*x^(n-1)+psolution(b*x,c/b*p!,y,x,p) )));
psolution(b_*x_^m_,x_^n_,y_,x_,p_) := if(isconstant(b),if(re(p)-n+m>0 and isinteger(n-m), -x^(n-m)/b, if(p==n-m,-1/b*x^(n-m)+psolution(b*x^m,1/b*p!,y,x,p) )));
#psolution(x_^m_,c_*x_^n_,y_,x_,p_):=if(isconstant(c),if(re(p)-n+m>0 and isinteger(n-m), -c*x^(n-m), if(p==n-m,-c*x^(n-m)+psolution(x^m,c*p!,y,x,p) )));
#psolution(b_*x_^m_,c_*x_^n_,y_,x_,p_):=if(isconstant(c/b),if(re(p)-n+m>0 and isinteger(n-m), -c/b*x^(n-m), if(p==n-m,-c/b*x^(n-m)+psolution(b*x^m,c/b*p!,y,x,p) )));

psolution(b_*x_^m_,c_*x_^n_,y_,x_,q_) := if(isconstant(c/b),
	if(m+q==0,c*x^(n-m)/(fallingfactorial(n-m,q)-b), 
	if(q-n+m>0 and n-m>0,-x^(n-m)*c/b, 
	if(q==n-2m,-x^(n-m)*c/b-c*fallingfactorial(n-m,q)/b^q, 
	if(q==n-m,-x^(n-m)*c/b+psolution(b*x^m,q!*c/b,y,x,q) )))));
psolution(b_/x_,c_*x_,y_,x_,q_):= if(isconstant(c/b),
	if(q==2, -2c/b^2*x-c/b*x^2,
	if(q>2, -c/b*x^2 )));

#psolution(b_*x_^m_,x_^n_,y_,x_,q_) := if(isconstant(b),if(m== -q,x^(n-m)/(fallingfactorial(n-m,q)-b), 
	-1/b*x^(n-m)-if(q>=n-m,0, fallingfactorial(n-m,q)/b^q ) ));
psolution(b_*x_^m_,x_^n_,y_,x_,q_) := if(isconstant(b),
	if(m+q==0,x^(n-m)/(fallingfactorial(n-m,q)-b), 
	if(q-n+m>0 and n-m>0,-x^(n-m)/b, 
	if(q==n-2m,-x^(n-m)/b-fallingfactorial(n-m,q)/b^q, 
	if(q==n-m,-x^(n-m)/b+psolution(b*x^m,q!/b,y,x,q) )))));
psolution(b_*x_^m_,x_,y_,x_,q_) := if(isconstant(b),
	if(m+q==0,x^(1+q)/((1+q)!-b), 
	if(q==2,-2/b^2*x-x^(1-m)/b, 
	if(q>1-m,-x^(1-m)/b  ))));
psolution(b_*x_,x_^n_,y_,x_,q_):= if(isconstant(b), 
	if(n>1 and n-q<1,-x^(n-1)/b,
	if(q>1 and n-q>1 and  n<2q+2,-fallingfactorial(n-1,q)*x^(n-2-q)/b^2-x^(n-1)/b )));
psolution(b_/x_,x_^n_,y_,x_,q_):= if(isconstant(b),
	if(n> -1 and q>1+n,-x^(1+n)/b, 
	if(q>2 and q-n>=0 and q-n<=1,-x^(1+n)/b-fallingfactorial(n+1,q)*x^(n+2-q)/b^2,  
	if(q>1 and n>0,-sum(n!*(n+1)!/(n+1-k)!/(n-k)!*b^round(-1-k)*x^(n+1-k),k,0,n,q-1) ))));
#psolution(c_*x_,x_^n_,y_,x_,q_):= if(isconstant(c), if(q>=n-1,-1/c*x^(n-1), if(n-q==2,-(n-1)!/c/c-x^3/c, 
	if(q==2,pi*Bi(cbrt(c)*x)*int(Ai(cbrt(c)*x)*x^n,x)/cbrt(c)-pi*Ai(cbrt(c)*x)*int(Bi(cbrt(c)*x)*x^n,x)/cbrt(c) ))));
#psolution(c_*x_,x_^n_,y_,x_,q_):= if(isconstant(c), if(q>=n-1,-1/c*x^(n-1), if(n-q==2,-(n-1)!/c^2-x^3/c )));

#psolution(1/x_,c_*x_,y_,x_,q_):= if(isconstant(c), 
	if(q==2,-2c*x-c*x^2, 
	if(q>2,-c*x^2 )));
#psolution(1/x_,c_*x_^n_,y_,x_,q_):= if(isconstant(c) and n <> -1,  if(q==1,c*x^(n+1)/n,
	if(n<q-1,-c*x^(1+n), 
	if(n>=2*q-3,-c*fallingfactorial(n+1,q)*fallingfactorial(n+2-q,q)*x^(n+3-2q)-c*fallingfactorial(n+1,q)*x^(n+2-q)-c*x^(1+n),
	if(n>=q-2,-c*fallingfactorial(n+1,q)*x^(n+2-q)-c*x^(1+n)	)))));

psolution(x_^m_,x_^n_,y_,x_,q_) := if(m+q==0,x^(n-m)/(fallingfactorial(n-m,q)-1), 
	if(q==n-2m,-x^(n-m)-fallingfactorial(n-m,q), 
	if(q==n-m,-x^(n-m)+psolution(x^m,q!,y,x,q),
	if(q>n-m and n-m>0,-x^(n-m) 
	))));
psolution(x_^m_,x_,y_,x_,q_) := if(m+q==0,x^(1-m)/((1-m)!-1), 
	if(q==2,-2x-x^(1-m), 
	if(q==1-2m,-x^(1-m)-fallingfactorial(1-m,q), 
	if(q==1-m,-x^(1-m)+psolution(x^m,q!,y,x,q),
	if(q>1-m and 1-m>0,-x^(1-m) 
	)))));
#psolution(x_,x_^n_,y_,x_,q_):= if(q-n== -2,-(n-1)!-x^(n-1), if(q-n== -1,-(n)!*x-x^(n-1), if(q-n>-1,-x^(n-1), if(q==2,pi*Bi(x)*int(Ai(x)*x^n,x)-pi*Ai(x)*int(Bi(x)*x^n,x) ))));
#psolution(x_,x_^n_,y_,x_,q_):= if(n>1 and n<1+q,-x^(n-1), if(q>1 and n-q>1 and n<2q+2,-fallingfactorial(n-1,q)*x^(n-2-q)-x^(n-1),  if(q>1 and n>1,-sum((n-1)!*(n)!/(n-1-k)!/(n-k)!*x^(n-1-k),k,0,n-1,q+1) )));
#psolution(1/x_,x_^n_,y_,x_,q_):= if(q==n,-(1+n)!*x^2-x^(1+n), if(q-n==1,-q!*x-x^q, if(n> -1 and q-n>1,-x^(1+n), if(q>1 and n>0,-sum((n+1)!/(n-2/3*k+1)!*x^(n-k+1),k,0,n+1,q-1) ))));
#psolution(1/x_,x_^n_,y_,x_,q_):= if(q>2 and q-n>=0 and q-n<=1,-fallingfactorial(n+1,q)*x^(n+2-q)-x^(1+n), if(n> -1 and q-n>1,-x^(1+n), if(q==2 and n>0,-sum((n)!/(n-k)!*(n+1)!/(n+1-k)!*x^(n+1-k),k,0,n,q-1) )));
#psolution(1/x_,x_^n_,y_,x_,q_):= if(q>2 and n/(q-1)<3 and q-n>= -2 and q-n<=1,-fallingfactorial(n+1,q)*fallingfactorial(n+2-q,q)*x^(n+3-2q)-fallingfactorial(n+1,q)*x^(n+2-q)-x^(1+n), ;
psolution(1/x_,x_^n_,y_,x_,q_):= if(n == -1,-1, if(q==1,x^(n+1)/n,
	if(q>n+1,-x^(1+n), 
	if(n>=2*q-3,-fallingfactorial(n+1,q)*fallingfactorial(n+2-q,q)*x^(n+3-2q)-fallingfactorial(n+1,q)*x^(n+2-q)-x^(n+1),
	if(q<=n+2,-fallingfactorial(n+1,q)*x^(n+2-q)-x^(1+n)	
	)))));
psolution(x_,x_^n_,y_,x_,q_):=  if(n>1 and isinteger(n),	
	if(q>n-1, -x^(n-1), 
	if(n>=2*q+3,-fallingfactorial(n-1,q)*fallingfactorial(n-2-q,q)*x^(n-3-2q)-fallingfactorial(n-1,q)*x^(n-2-q)-x^(n-1),
	if(q<=n-2,-fallingfactorial(n-1,q)*x^(n-2-q)-x^(n-1) 	))));


psolution(1/x_,x_^n_,y_,x_,2):=	if(n>0,-sum(n!/(n-k)!*(n+1)!/(n+1-k)!*x^(n+1-k),k,0,n) );
psolution(1/x_,x_^6,y_,x_,3):= (-75600.0)*x-12600*x^3-210*x^5-x^7;
psolution(x_,x_^3,y_,x_,2):= -2Ai(x)*int(Bi(x),x)*pi+2Bi(x)*int(Ai(x),x)*pi-x^2;
psolution(x_,x_^6,y_,x_,2):= -40Ai(x)*int(Bi(x),x)*pi+40Bi(x)*int(Ai(x),x)*pi-20x^2-x^5;
psolution(x_,x_^9,y_,x_,2):= -2240Ai(x)*int(Bi(x),x)*pi+2240Bi(x)*int(Ai(x),x)*pi-1120*x^2-56*x^5-x^8;

psolution(c_,c_,y_,x_,p_):= if(p>0,-1);
psolution(c_,-c_,y_,x_,p_):= if(p>0,1);
psolution(-c_,c_,y_,x_,p_):= if(p>0,1);
psolution(c_*xx_,xx_,y_,x_,p_):= if(p>0 and isconstant(c),-1/c);
psolution(xx_,c_*xx_,y_,x_,p_):= if(p>0 and isconstant(c),-c);
psolution(c_*xx_,d_*xx_,y_,x_,p_):= if(p>0 and isconstant(c,d),-d/c);

psolution(b_*x_,x_^n_,y_,x_,-1):= if(isconstant(b),x^(n-1)/(fallingfactorial(n-1,-1)-b) );
psolution(b_*x_,x_,y_,x_,-1):= if(hasnot(b,x),1/(1-b) );
psolution(x_,x_^n_,y_,x_,-1):= x^(n-1)/(fallingfactorial(n-1,-1)-1);

psolution(0,c_,y_,x_,n_):=if(isfree(c,y), int(c,x,n) );