# differential equation solver; # solve differential equation and fractional differential equation for function 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_+b_, y_,x_,q_) := If(isfree(a,b,y), d(a+b,x,-q)+C_1+if(q>1, C_2*x,0), If(isfree(a,b,x), solve(d(1/(a+b),y,-q)=x^q/q!,y)+if(q>1, C_2*x^2+C_1*x,C_1+x), If(isfree(a,y) and isfree(d(b,y),y), block(dsolve:=d(d(b,y),x,-q),simplify(mittag(q,dsolve)*d(a*mittag(q,-dsolve),x,-q))+mittag(q,dsolve)), If(isfree(b,y) and isfree(d(a,y),y), block(dsolve:=d(d(a,y),x,-q),simplify(mittag(q,dsolve)*d(b*mittag(q,-dsolve),x,-q))+mittag(q,dsolve)) )))); dsolve(a_+y, y_,x_,q_) := If(isfree(a,y), If(isfree(a,x), If(q>0,C_1*exp(x)-a, C_1*mittag(-q,x^(-q)/(-q)!) ), simplify(mittag(q,x^q/q!)*d(a*mittag(q,-x^q/q!),x,-q))+mittag(q,x^q/q!) )); dsolve(a_+b_*y, y_,x_,q_) := If(isfree(a,b,y), if(isfree(a,b,x), If(q<0, C_1*mittag(-q,1/b*x^(-q)/(-q)!), If(b>0, C_1*exp(b^(1/q)*x)-a/b, C_1*mittag(q,b*x^q/q!)-a/b )), If(a==b, if(q>0,C_1*mittag(q,d(a,x,-q))-1, C_1*mittag(-q,d(a,x,-q))/d(a,x,-q) ) ))); dsolve(a_+a_*y, y_,x_,q_) := If(isfree(a,y), If(isfree(a,x), If(q>0, C_1*exp(a^(1/q)*x)-1, C_1*mittag(-q,1/a*x^(-q)/(-q)!) ), If(q>0,C_1*mittag(q,d(a,x,-q))-1, C_1*mittag(-q,d(a,x,-q))*d(a,x,-q) ) )); #dsolve(x_^q_+y_, y_,x_,q_) := -x^q-q!; #dsolve(x_^n_+b_*y_, y_,x_,q_) := If(isfree(b,x) and isfree(b,y), If(n==q, -x^q/b-q!/b^2, If(q>0, n!*x^(n+q)*mittag(q,1+n+q,b*x^q), gamma(n+1)*x^n*mittag(-q,n+1,1/b*x^(-q)) ))); #dsolve(x_^n_+b_*y_, y_,x_,q_) := If(isfree(b,x) and isfree(b,y), If(n==q, -x^n/b-q!/b^2, n!*x^(n+q)*mittag(q,1+n+q,b*x^q) )); #dsolve(a_*x_^n_+b_*y_, y_,x_,q_) := If(isfree(a,b,x) and isfree(a,b,y), If(n==q, -a*x^n/b-a*q!/b^2, a*n!*x^(n+q)*mittag(q,1+n+q,b*x^q) )); dsolve(x_^n_+y, y_,x_,q_) := If(n==q, -x^n-n!, n!*x^(n+q)*mittag(q,1+n+q,x^q)); dsolve(a_*x_^n_+y, y_,x_,q_) := If(isfree(a,x) and isfree(a,y), If(n==q, -a*x^n-a*n!, a*n!*x^(n+q)*mittag(q,1+n+q,x^q) )); dsolve(x_^n_+y+c_, y_,x_,q_) := If(isfree(c,x) and isfree(c,y), If(n==q, -x^n-n!-c, n!*x^(n+q)*mittag(q,1+n+q,x^q)-c)); dsolve(a_*x_^n_+y+c_, y_,x_,q_) := If(isfree(a,c,x) and isfree(a,c,y), If(n==q, -a*x^n-a*n!-c, a*n!*x^(n+q)*mittag(q,1+n+q,x^q)-c )); dsolve(x_^n_+b_*y_+c_, y_,x_,q_) := If(isfree(b,c,x) and isfree(b,c,y), If(n==q, -x^n/b-q!/b^2-c/b, n!*x^(n+q)*mittag(q,1+n+q,b*x^q)-c/b )); dsolve(a_*x_^n_+b_*y_+c_, y_,x_,q_) := If(isfree(a,b,c,x) and isfree(a,b,c,y), If(n==q, -a*x^n/b-a*q!/b^2-c/b, a*n!*x^(n+q)*mittag(q,1+n+q,b*x^q)-c/b )); dsolve(a_*x_+b_*y_, y_,x_,q_) := If(isfree(a,b,x) and isfree(a,b,y), a*x^(1+q)*mittag(q,2+q,b*x^q)); dsolve(a_*x_+y_, y_,x_,q_) := If(isfree(a,x) and isfree(a,y), a*x^(1+q)*mittag(q,2+q,x^q)); dsolve(x_+b_*y_, y_,x_,q_) := If(isfree(b,x) and isfree(b,y), x^(1+q)*mittag(q,2+q,b*x^q) ); dsolve(x_+y, y_,x_,q_) := x^(1+q)*mittag(q,2+q,x^q); dsolve(a_*x_+y_+c_, y_,x_,q_) := If(isfree(a,c,x) and isfree(a,c,y), If(q>0,a*x^(1+q)*mittag(q,2+q,x^q)-c,a*x^(1+q)*mittag(q,2+q,x^q) )); dsolve(x_+b_*y_+c_, y_,x_,q_) := If(isfree(b,c,x) and isfree(c,b,y), If(q>0,x^(1+q)*mittag(q,2+q,b*x^q)-c/b,x^(1+q)*mittag(q,2+q,b*x^q) )); dsolve(x_+y_+c_, y_,x_,q_) := If(isfree(c,x) and isfree(c,y), If(q>0,x^(1+q)*mittag(q,2+q,x^q)-c,x^(1+q)*mittag(q,2+q,x^q) )); dsolve(a_*x_+b_*y_+c_, y_,x_,q_) := If(isfree(a,b,c,x) and isfree(a,b,c,y), If(q>0, a*x^(1+q)*mittag(q,2+q,b*x^q)-c/b,a*x^(1+q)*mittag(q,2+q,b*x^q) )); dsolve(a_*x_+x_*y_, y_,x_,q_) := If(isfree(a,x) and isfree(a,y), mittag(q,x^(1+q)/(1+q)!)-a); dsolve(a_*x_+b_*x_*y_, y_,x_,q_) := If(isfree(a,b,x) and isfree(a,b,y), mittag(q,b*x^(1+q)/(1+q)!)/b-a/b); dsolve(x_+x_*y_, y_,x_,q_) := mittag(abs(q),x^(1+abs(q))/(1+abs(q))!)+if(q>0, -1,0) ; dsolve(b_*y^n_, y_,x_,q_) := If(isfree(b,y), If(isfree(b,x), b^(1-n)*fallingfactorial(q/(1-n),q)^(1/(n-1))*(x+C_1)^(q/(1-n)), (1/fallingfactorial(-n,q)*d(b,x,-q))^(1/(q-n)) )); dsolve(y_^n_, y_,x_,q_) := fallingfactorial(q/(1-n),q)^(1/(n-1))*(x+C_1)^(q/(1-n)) ; dsolve(b_*y_, y_,x_,q_) := If(b>0, C_1*exp(b^(1/q)*x)+if(iseven(q),C_2*exp(-b^(1/q)*x),0), If(b<0, C_1*mittag(q,b*x^q/q!), C_1*mittag(q,d(b,x,-q)) )); dsolve(b_*x_*y_, y_,x_,q_) := If(b>0, C_1*exp(b^(1/q)*q/(1+q)*x^(1/q+1)), If(b<0, If(q>0, C_1*mittag(q,b*x^(q+1)/(1+q)!), C_1*mittag(-q,1/b*x^(1-q)) ))); dsolve(x_*y_, y_,x_,q_) := If(q>0, C_1*exp(q/(1+q)*x^(1/q+1)), C_1*mittag(-q,x^(1-q)/(1-q)!)); dsolve(y_, y_,x_,q_) := if(iseven(q),C_1*exp(-x)+C_2*exp(x), C_1*exp(x)); dsolve(f_, y_,x_,1) :=dsolve(f, y,x); dsolve(y(p_,x_),y_,x_,q_):= C_1*e^x+if(p>0,C_2,0); dsolve(b_*y(p_,x_),y_,x_,q_):=if(isfree(b,y), if(b>0, C_1*exp(b^(1/(q-p))*x)+if(p>0,C_2,0), C_1*mittag(q-p,d(b,x,p-q))+if(p>0,C_2,0) )); dsolve(y(p_,x_)+c_,y_,x_,q_):= d(dsolve(c+y,y,x,q-p),x,-p)+if(p>0,C_2,0); dsolve(b_*y(p_,x_)+c_,y_,x_,q_):=if(isfree(b,y), d(dsolve(c+b*y,y,x,q-p),x,-p))+if(p>0,C_2,0); dsolve(y(p_,x_)+c_*y, y_,x_,q_) := if(isfree(c,y) and isfree(c,x), if(q==2p, C_1*e^(((1-sqrt(1+4c))/2)^(1/p)*x)+C_2*e^(((1+sqrt(1+4c))/2)^(1/p)*x) )); dsolve(y(p_,x_)+y, y_,x_,q_) := if(q==2p, C_1*e^(((1-sqrt(5))/2)^(1/p)*x)+C_2*e^(((1+sqrt(5))/2)^(1/p)*x) ); dsolve(b_*y(p_,x_)+c_*y, y_,x_,q_) := if(isfree(b,c,y) and isfree(b,c,x), if(q==2p, C_1*e^(((b-sqrt(b^2+4c))/2)^(1/p)*x)+C_2*e^(((b+sqrt(b^2+4c))/2)^(1/p)*x), if(b+c==1, C_1*e^(x) ))); dsolve(b_*y(p_,x_)+y, y_,x_,q_) := if(isfree(b,y) and isfree(b,x), if(q==2p, C_1*e^(((b-sqrt(b^2+4))/2)^(1/p)*x)+C_2*e^(((b+sqrt(b^2+4))/2)^(1/p)*x) )); dsolve(b_*y(p_,x_)+c_*y+d_, y_,x_,q_) := if(isfree(b,c,d,y) and isfree(b,c,d,x), dsolve(b*y(p,x)+c*y, y,x,q)-d/c); dsolve(b_*y(p_,x_)+y+d_, y_,x_,q_) := if(isfree(b,d,y) and isfree(b,d,x), dsolve(b*y(p,x)+y, y,x,q)-d); dsolve(y(p_,x_)+c_*y+d_, y_,x_,q_) := if(isfree(c,d,y) and isfree(c,d,x), dsolve(y(p,x)+c*y, y,x,q)-d/c); dsolve(y(p_,x_)+y+d_, y_,x_,q_) := if(isfree(d,y) and isfree(d,x), dsolve(y(p,x)+y, y,x,q)-d); dsolve(y(1,x_),y_,x_,3):=C_1+C_2*e^(x)+C_3*e^(-x); dsolve(y(2,x_),y_,x_,3):=C_1+C_2*e^(x)+C_3*x; dsolve(a_+b_*y_, y_,x_,2) := if(isfree(a,b,y), if(isfree(a,b,x), C_1*exp(sqrt(b)*x)+C_2*exp(-sqrt(b)*x)-a/b, if(a==b,log(1+y)*(1+y)-y=int(int(a,x),x)+C_1*x+C_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)-1/2*exp(-x)*int(exp(x)*a,x)+1/2*exp(x)*int(exp(-x)*a,x) )); dsolve(x_+y_, y_,x_,2) := C_1*exp(x)+C_2*exp(-x)-x; dsolve(x_+b_*y_, y_,x_,2) := if(isfree(b,x), C_1*exp(sqrt(b)*x)+C_2*exp(-sqrt(b)*x)-x/b); dsolve(a_*x_+b_*y_, y_,x_,2) := if(isfree(a,b,x), C_1*exp(sqrt(b)*x)+C_2*exp(-sqrt(b)*x)-a*x/b); dsolve(a_*x_+b_*y_+c_, y_,x_,2) := if(isfree(a,b,c,x), C_1*exp(sqrt(b)*x)+C_2*exp(-sqrt(b)*x)-a*x/b-c/b); dsolve(b_*y_, y_,x_,2) := if(isfree(b,y),if(isfree(b,x), if(b<0, C_1*cos(sqrt(-b)*x)+C_2*sin(sqrt(-b)*x), C_1*exp(sqrt(b)*x)+C_2*exp(-sqrt(b)*x) ), log(y)*y-y=int(int(b,x),x)+C_1*x+C_2 )); dsolve(x_*y_, y_,x_,2) := C_1*AiryAi(x)+C_2*AiryBi(x); dsolve(a_*x_*y_, y_,x_,2) := if(isfree(a,x) and isfree(a,y), C_1*AiryAi(cbrt(a)*x)+C_2*AiryBi(cbrt(a)*x)); dsolve(b_*y(1,x_),y_,x_,2):=if(isfree(b,y),if(isfree(b,x), C_1+C_2*e^(b*x), C_1*integrate(exp(int(b,x)),x)+C_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*e^(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*e^(b*x)/b-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*e^(b*x)/b-c*x/b, C_1+C_2*e^(b*x)/b+integrate(exp(b*x)*int(exp(-b*x)*c,x),x) )); dsolve(y+y(1,x_),y_,x_,2):= C_1*e^((1-sqrt(5))/2*x)+C_2*e^((1+sqrt(5))/2*x); dsolve(b_*y(1,x_)+y,g_,x_,2):=if(isfree(b,g) and isfree(b,x), C_1*e^((b-sqrt(b^2+4))/2*x)+C_2*e^((b+sqrt(b^2+4))/2*x)); dsolve(y(1,x_)+c_*y,g_,x_,2):=if(isfree(c,g) and isfree(c,x), C_1*e^((1-sqrt(1+4c))/2*x)+C_2*e^((1+sqrt(1+4c))/2*x)); dsolve(b_*y(1,x_)+c_*y,g_,x_,2):=if(isfree(b,c,g) and isfree(b,c,x), if(b*b+4c==0, C_1*exp(b*x/2)+C_2*exp(b*x/2)*x, if(b*b+4c<0, C_1*e^(b*x/2)*cos(sqrt(-b^2-4c)/2*x)+C_2*e^(b*x/2)*sin(sqrt(-b^2-4c)/2*x),C_1*e^((b-sqrt(b^2+4c))/2*x)+C_2*e^((b+sqrt(b^2+4c))/2*x) ))); dsolve(b_*y(1,x_)+y+d_,g_,x_,2):=if(isfree(b,d,g) and isfree(b,d,x),dsolve(b*y(1,x)+y,g,x,2)-d*x-d); dsolve(y(1,x_)+c_*y+d_,g_,x_,2):=if(isfree(c,d,g) and isfree(c,d,x),dsolve(y(1,x)+c*y,g,x,2)-d*x/c-d/c^2); dsolve(b_*y(1,x_)+c_*y+d_,g_,x_,2):=if(isfree(b,c,d,g) and isfree(b,c,d,x),dsolve(b*y(1,x)+c*g,g,x,2)-d/c, if(b*b+4c==0, C_1*exp(b/2*x)+C_2*exp(b/2*x)*x+expand(exp(b/2*x)*x*int(exp(-b/2*x)*d,x)-exp(b/2*x)*int(d*exp(-b/2*x)*x,x)), C_1*e^((b-sqrt(b^2+4c))/2*x)+C_2*e^((b+sqrt(b^2+4c))/2*x)+ expand(-e^((b-sqrt(b^2+4c))/2*x)*int(d*e^((-b+sqrt(b^2+4c))/2*x),x)+e^((b+sqrt(b^2+4c))/2*x)*int(d*e^((-b-sqrt(b^2+4c))/2*x),x)) )); dsolve(y_, y_,x_,2) := C_1*exp(-x)+C_2*exp(x); #dsolve(a_+b_, y_,x_) := If(isfree(a,b,y), integrate(a+b,x)+C_1, If(has(a,y) and has(b,y) and isfree(d(a+b,y),y), exp(integrate(d(a+b,y),x))*C_1, If(isfree(b,y) and isfree(d(a,y),y), block(dsolve:=int(d(a,y),x),exp(dsolve)*int(b*exp(-dsolve),x)+exp(dsolve)*C_1), If(isfree(a,y) and isfree(d(b,y),y), block(dsolve:=int(d(b,y),x),exp(dsolve)*int(a*exp(-dsolve),x)+exp(dsolve)*C_1) )))); dsolve(a_+b_*y_,y_,x_):=If(isfree(a,b,y), If(isfree(a,b,x), C_1*exp(b*x)-a/b, block(dsolve:=integrate(b,x), expand(exp(dsolve)*integrate(a*exp(-dsolve),x))+exp(dsolve)*C_1 ) )); dsolve(a_+a_*y_,y_,x_):=If(isfree(a,y),C_1*exp(integrate(a,x))-1); dsolve(a_+y_,y_,x_):=If(isfree(a,y), If(isfree(a,x), C_1*exp(x)-a, expand(exp(x)*integrate(a*exp(-x),x))+exp(x)*C_1 )); dsolve(a_*y_+b_*y_^n_,y_,x_):=If(isfree(a,b,y) , block(dsolve:=integrate(a,x), exp(dsolve)*(1+(1-n)*integrate(b*exp((n-1)*dsolve),x))^(1/(1-n)) ) ); dsolve(a_*y_+y_^n_,y_,x_):=If(isfree(a,n,y) , block(dsolve:=integrate(a,x), exp(dsolve)*(1+(1-n)*integrate(exp((n-1)*dsolve),x))^(1/(1-n)) ) ); dsolve(y_+y_^n_,y_,x_):= exp(x)*(C_1-exp((n-1)*x))^(1/(1-n)); dsolve(x_+y_, y_,x_) := C_1*exp(x)-x-1; #dsolve(a_*x_+y_, y_,x_) := If(isfree(a, x) and isfree(a, y), -a*x-a); #dsolve(x_+b_*y_,y_,x_) := If(isfree(b, x) and isfree(b, y), -1/b*x-1/b^2); #dsolve(a_*x_+b_*y_,y_,x_):=If(isfree(a,b, x) and isfree(a,b, y), -a/b*x-a/b^2); #dsolve(x_+y_+c_, y_,x_) := If(isfree(c,x) and isfree(c,y), -x-1-c); #dsolve(x_+b_*y_+c_, y_,x_) := If(isfree(b,c,x) and isfree(b,c,y), -1/b*x-1/b^2-c/b); #dsolve(a_*x_+y_+c_, y_,x_) := If(isfree(a,c,x) and isfree(a,c,y), -a*x-a-c); #dsolve(a_*x_+b_*y_+c_, y_,x_) := If(isfree(a,b,c,x) and isfree(a,b,c,y), -a/b*x-a/b^2-c/b); dsolve(g_*y_, y_,x_) := If(isfree(g,y), C_1*exp(integrate(g,x))); dsolve(g_*y_^n_, y_,x_) := If(isfree(g,y), ((1-n)*integrate(g,x)+C_1)^(1/(1-n)) ); dsolve(y_,y_,x_) := C_1*exp(x); #dsolve(y(n_,x_),y_,x_) := C_1*exp(x); #dsolve(f_,y(x_)) := dsolve(replace(d(y(x),x)-f,y(x),y),y); #dsolve(a_=b_,y_):= if(has(a,ds(y,x,2)),dsolve(ds(y,x,2)-a+b,y,x,2), if(has(a,ds(y)),dsolve(ds(y)-a+b,y,x), if(has(a,dy),dsolve(dy/dx-a+b,y,x), dsolve(a-b,y) ))); #dsolve(y'=f_,y) := dsolve(f,y,x); #dsolve(y''=f_,y) := dsolve(f,y,x,2); #dsolve(dy/dx=f_,y) := dsolve(f,y,x); #dsolve(d(y(x_),x_)=f_,y(x_)) := dsolve(replace(f,y(x),y),y,x); #dsolve(ds(y_)=a_,y_) := dsolve(a,y,x); #dsolve(ds(y_,x_)=a_,y_) := dsolve(a,y,x); #dsolve(ds(y_,x_,q_)=a_,y_):= dsolve(a,y,x,q); #dsolve(a_+ds(y_,x_,q_),y_):= dsolve(-a,y,x,q); #dsolve(ds(y_,x_,q_)+z_,y_):= dsolve(-z,y,x,q); #dsolve(ds(y_,x_,q_)+f_+z_,y_):= dsolve(-f-z,y,x,q); #dsolve(a_+ds(y_,x_,q_)+z_,y_):= dsolve(-a-z,y,x,q); dsolve(a_+y(q_,x_),y_):= dsolve(-a,y,x,q); dsolve(a_+b_*y(q_,x_),y_):= dsolve(-expand(a/b),y,x,q); dsolve(a_+y(q_,x_),y(x_)):= dsolve(-replace(a,y(x),y),y,x,q); dsolve(a_=b_,y_):=dsolve(toyn(a-b),y); dsolve(y(n_,x_)=a_,y_) := dsolve(toyn(a),y,x,n); dsolve(a_):=dsolve(toyn(a),y); dsolve(d(y(x_),x_)=f_) := dsolve(replace(f,y(x),y),y,x); y':=y(1,x); y'':=y(2,x);