@mehdibaghaee 

The idea is that if you are solving, say, dy(x)/dx = -k*y(x) with y(0)=1 (y(x) is "unknown") the solution is y(x)=exp(-k*x), and if you change variables x=2*X, then the solution is y(X)=exp(-2*k*X) in which y(X) is a different function than y(x). But if you have a known function like sin(x) then you want sin(x) to change to sin(2*X)
 

Download inttt3.mw

@mehdibaghaee  add known=u0 to get this.

@mehdibaghaee Looks like a bug

Using dchange gives the expected double integral
 

Download inttt2.mw

aside from the main issue raised by @rlopez your loop should be "from 0 to 99" not "by 0 to 99"

@umar khan  

fsolve(y, x = 0 .. 2e6);

works for me.

@leslieglasser Still not clear what the possible inputs are. If it is a monic polynomial and the other coefficents are the variables (as in your example), then the following works:

I don't have any suggestions for more complicated cases.
 

Download Companion3.mw

@leslieglasser  Well you can make a 3x3 one by bordering the 2x2 one with zeroes except for one in the 1,1 entry. But I'm guessing that isn't what you want. Can you be more specifc about your exact requirements? - I don't know what separated variables mean or nomograph in this context.

int(Heaviside(1-x^2-y^2), [x=0..1,y=0..1]); #should be Pi/4

gives -Pi/4

(I thought using even symmetry might help, but it went negative.)

@dharr If you set x=Pi/2 rather than solving for it, then you get the a, b, c values that you want.

ab should be a*b, but if I've interpreted your equation correctly, then odetest suggests that your solution is not correct.

does not return zero, meaning it is not a solution. Using dsolve gives something with Heun functions, but you probably need to check you have the right equation first.

@sunit I'm not entirely sure what you mean. You could certainly write a procedure that tested whether there were minus signs to its argument and adjusted things before passing to sin. But covering all the unexpected cases (the above one was certainly unexpected for me; I just mistakenly typed Phi instead of phi) is probably going to take more effort than it is worth. As I noted, Maple isn't human, and often outputs things in forms that are unexpected (but mathematically correct)

Download expand2.mw

@sunit 

I misread the correct answer at line (13) followed by the comment "The results are worse than single iteration using jacobian" to mean that you didn't think the default solver was working. I see you have now clarified that with an edit.

As for correctly coding lsode I don't have much comment; I have tried to move away low level coding.

FunctionAdvisor(hypergeom([a1,a2],[b1,b1],z)); returns "The system is unable to compute the "aymptotic expansion" for hypergeom([a1,a2],[b1,b2],z)" 

The  NIST handbook gives a formula but it's pretty messy:

http://dlmf.nist.gov/16.11#E7

I programmed mgear in FORTRAN long ago, and then used it in Maple, then Maple's lsode, then Maple's rosenbrock, and I haven't noticed any difficulties.

I can't tell exactly what you are doing, but it seems as if you know what the solution should be? What is that solution and why do you think the 0.3496  answer is wrong? There doesn't seem to be anything pathological about the DE. At least locally, the answer is right, as this shows:

Download NonlinearDE.mw