@9009134 You are still focusing on the contents of the file 'u(y)-y.dat' .  That's the wrong place to look.  Before you put anything in that file, you should verify that there is something useful to put in it.  Did you watch the output of lprint(tt,u) on the screen? Is that what you expect? Probably not.

Fix that error first, and then we will talk about why 'u(y)-y.dat' is empty.

@9009134 Forget about the file 'u(y)-y.dat' for now. Did you do what I suggested?

Temporarily replace the fprintf line with lprint(tt,u) so that you can inspect the values of tt and u.  You will see that u does not evaluate to a number. Track down the error.

The Q[..,3]:=Q[..,3]/Q[1,3] calculation adjusts the third eigenvector.  In general to adjust all columns in a matrix Q so that their first components are 1 (if possible) we may do:
Q := Q . DiagonalMatrix(1/~Row(Q,1));

@leafgreen As an alternative to Kitonum's correction, you may consider defining f as an indexed variable, as in:

f := Array(antisymmetric, 1..3, 1..3,
    (i,j) -> -m[i]*m[j]*(x[i]-x[j])/((x[i]-x[j])^2+(y[i]-y[j])^2) );

Then f[1,2] will produce the x component of the gravitational force.  You need to define another function to produce the y component.

I should point out that your force formula assumes that the gravitational force is proportional to the inverse of the distance between the particles.  Do you really mean that?  Newton says that the gravitational force is proportional to the square of the inverse distance.

In either case, I am afraid that your worksheet takes you in a direction which is not suitable at all for implementing Verlet's algorithms.  What happens later when the particles have moved?  Their distances change, and the forces need to be recomputed.  The current setup is not quite suited to that.

Someone needs to grab you and turn you around by about 90 degrees to point you in the correct direction.  Your professor should be able to help with that.

@Franzs If you want the axes to rotate along with the object, define P1 as before, and then do:

p := Plot(P1, view = [-40 .. 40, -40 .. 40, -40 .. 40],
    color = "Orange", thickness = 2,
    viewpoint=["circleleft", frames=60], axes=normal);

I don't understand what you mean by "without the names of its spatial coordinates".

@rrbaldino  You keep saying that the second answer is the solution although you plainly see that it fails to satisfy your initial condition.

Think about this:  If I give you several candidates for solution, how would you pick the right one?  Do you know a way other than checking that the differential equation and the initial conditions are satisfied?  Is there some other magical way?

Let me repeat what I said before:  Your problem is ill-posed.  The way to fix it, is to shift the time of the application of the delta function to somewhere after the initial time.

What did you get for the critical point?  My calculation gives T = -45926.73299, v = 0.2514345875.

But this smells like a question coming from thermodynamics, and therefore that temperature looks a bit suspicious to me.

@rrbaldino  Did you plot the two "solutions" as I suggested?  Neither of them satisfies the initial condition dy/dt=0 at t=0 that you have asked for.  One or the other answer that you have obtained may be a solution to some problem, but certainly not to your problem.

If anything, Maple should warn that both answers are wrong.

@vv What you have here is a tour de force of recursive programming.  Awesome!

I think that you are asking how to convert ln(x)/ln(10) to log₁₀(x). So far as I can tell that's not possible—Maple automatically converts logarithms from any base to base e.

For displaying purposes you may temporarily prevent the conversion by putting the desired expression in quotations, as in 'log[10](x)'; which yields log₁₀(x). But if you refer to that result later, it will revert to natural logarithms.

@nikkarkian In mathematics, "unknown" is what is to be solved for.  In your equation, V is not "unknown" because you are not solving for it.  The unkown is q(t), so all you need is:

EulerLagrange(L, t, q(t));

There are two parts to your question.  The first part is to figure out the shape and measurements of the desired net.  The second part is to draw that in Maple.

I suggest that you do the first part.  That way you can honestly claim that you did a part of the work for your project.  Then post here a hand-drawing of the net.  Many participant will be happy to show you how to do that drawing in Maple.

@Kitonum Thank you very much for this answer.  It does exactly what I was looking for.

You are attempting to solve

The differebtial equation says that the graph of  g(y) is a straight line.  But a straight line cannot fit those boundary conditions.

There may be other issues with the rest of your equations but this is the most glaring one.