@jud I am using version 2023, which has it, but 2022 doesn't.

In response to your email request, here is my version.

Download congruent.mw

The code in the worksheet should work in older versions of Maple. In 2023 it can be done a little shorter as:

congruent := proc(n::posint)
	local nn, An, Bn, c1, c2, c3, c4, x, y, p, q, mx, my, mz;
	description "Determine if n is a congruent number using Tunnell's theorem";
	# divide out largest square and use squarefree residue
	nn := mul(map(x -> x[1]^(x[2] mod 2), ifactors(n)[2]));
	if nn::odd then
  		c1, c2, c3, c4 := 2, 1, 32, 8
	else
		c1, c2, c3, c4 := 8, 2, 64, 16
	end if;
	An := 0;
	Bn := 0;
	for x from 0 while (p := nn - c1*x^2) >= 0 do
		for y from 0 while (q := p - c2*y^2) >= 0 do
			if issqr(q/c3) then
				# count all solutions e.g., for x nonzero there two possibilities of opposite sign
				mx := signum(x) + 1;
				my := signum(y) + 1;
				mz := signum(q) + 1;
				An := An + mx*my*mz;
			end if;
			if issqr(q/c4) then
				mx := signum(x) + 1;
				my := signum(y) + 1;
				mz := signum(q) + 1;
				Bn := Bn + mx*my*mz;
   			end if;
		end do
	end do:
	evalb(2*An = Bn);
end proc;

@mary120 

1- How can solve these two ODE(s) and  plot phi vs x in one phi-x coordinate?

Just solve them separately and generate two plots, in variables p1 and p2. Then display(p1,p2) will put them on the same plot.

2- How can I plot phi-x in a symmetric interval of x (for example, from x=-100..100).

The solver can go forwards and backwards from the initial point. But you didn't explain the significance of x=0 so it is not clear how the x axis is defined.

3- How can I export data of plot of these ODEs in the form of two distinct ascii table? (I couldn’t export data of the last figure in your last uploaded worksheet in the form of a ascii table)

I don't know why the export didn't work. If you upload a worksheet containing the error I can take a look. 

@mary120 if you assign the odeplot to a variable p, then M:=plottools:-getdata(p)[3]; gives a Matrix with the data. The matrix can then be exported with ExportMatrix to varrous  file types, e.g., target = csv will be comma separated. 

@mary120 To plot phi vs x just use

odeplot(ans1,[x(phi),phi],1e-2..0.5);

This goes through the origin because of the initial condition x(0)=0. But your plot doesn't go through the origin so I don't understand what you really want. Note the slope is wrong, so you probably need diff(x(phi),phi)=-integrand

Integral-New.mw

Here is something, but it is unclear how the x axis is constructed - what does x=0 mean on your plot?

something.mw

@Carl Love My original example was a graph with some randomly chosen edges. Later @Anthrazit wanted it to be a complete graph, and so I gave a version that does the complete graph for the same points, which does give the same path and distance as you find, see the worksheet TravelingSalesmanCompleteGraph.mw above.

I was also assuming that the start and end vertices were prescribed, and were to be given as the first and last in the list (since @Anthrazit mentioned something about nearest and outermost bolts), rather than just the shortest of all Hamiltonian paths. 

@Saha   You have gamma(Theta - Theta[a]) instead of gamma*(Theta - Theta[a])

@sand15 Yes! That syntax error can be avoided by diff(..., [eta$(m-2)]) , which returns the undifferentiated function. But probably not what the OP intended.

If I comment out with # the line setting values to various variables, then look at the output of the main loop I notice:

rps1.mw

1. Theta' is being interpreted as a derivative with respect to x. Probably you want to use declare(prime=eta)

2. You have multiplication before square brackets [] in two places, which Maple does not know how to interpret - if these are just for grouping you should use ().

3. The derivative with respect to eta m-2 times has not worked - the "d" should appear as upright, not italics. Use the calculus palette to form this, or use diff( ... , eta $ (m-2))

4. You have both Theta and Theta[a] and Theta with other subscripts. If it has subscripts you should not use Theta without subscripts as well.

5. If you are intending Theta without subscripts to be a function of eta you need to write Theta(eta)  or use declare so that Maple does not think it is a constant when it is differentiating.

6. use add instead of sum for just adding a finite number of terms.

Perhaps after fixing these, it will be clearer what you want to do, but at the moment I do not understand.

In Windows, if you click on the plot, the controls appear.

@Anthrazit Here's the CompleteGraph version (usually in traveling salesman problems there are not roads between every combination of cities). If you want shortest paths between two vertices that don't have to visit all vertices, DijkstrasAlgorithm or BellmanFordAlgorithm can be used - both can find shortest paths from vertex 1 to all other vertices, and then you can use an undirected graph - a path will not revisit an edge). If that is what you want I can set it up.

TravelingSalesmanCompleteGraph.mw

@C_R Yes, I forgot the abs; it works OK with it. I only used it without abs for a single "period" where I knew the sign, when I was playing around with it. I agree that the integral to infinity without abs is not relevant to @Mariusz Iwaniuk's integral, but I'm not sure whether or not the answer is wrong.

@Rouben Rostamian  PDEtools:-dchange gives the simpler result that you found, but comparing this with the IntegrationTools:-Change result shows that they are equivalent, so this is not a bug.

I1:=Int(sin(x^4)/(sqrt(x)+x^2),x=0..infinity);
PDEtools:-dchange(x=u^(1/4),I1);

dchange.mw

warnlevel = 0 works for me in Maple 2023.2.1

warnlevel.mw

@jalal To draw an arc representing the angle A in the triangle ABS, you can first find the cross product of the vector AB and AS to use as a rotation axis. Then find the rotation matrix tor rotation about this axis by arbitrary angle t. Apply this to a shortened version of vector AB, and its components (offset by A) are a parametric curve for the arc with parameter t, where t runs from zero to the angle A.

triangle(ABS, [A, B, S]);
angleSAB := FindAngle(ABS, A);
vAB := Vector(coordinates(B) - coordinates(A));
vAS := Vector(coordinates(S) - coordinates(A));
normalA := LinearAlgebra:-CrossProduct(vAB, vAS);
arcBAS := convert(0.2*((Student:-LinearAlgebra:-RotationMatrix(t, normalA)) . vAB), list);
plotarcBAS := spacecurve(coordinates(A) + arcBAS, t = 0 .. angleSAB, color = red, thickness = 2);
display(draw(ABS), plotarcBAS);

The code is at the end of the document.

rotation.mw