The files may not be written to disk on the OS until you either restart the session, or fclose them, or close the worksheet.

So, you can add this line inside the end of the loop.
   fclose(pout);

I am presuming that you will also uncomment the line that executes the plot. You may also want it to be,
  print(plot(...));
if your looping is further nested, or if all this is done within a procedure.

For example,

for j from 1 to 3 do
pout := cat("C:/Users/Eli/Documents/Animation/", "file", j, ".bmp"):
print(pout):
plotsetup(bmp, plotoutput = pout):
print(plot(...));
fclose(pout);
od:

You don't need differential equations for this. You just need to know how to differentiate and substitute.

For brevity, there is this one-liner:

eval( D(2/f)(-1), [D(f)(-1)=2, f(-1)=4] );

                -1
                --
                 4

Or, slightly longer,

eval( D(h=2/f)(-1), [D(f)(-1)=2, f(-1)=4] );

            D(h)(-1) = -1/4

In 2D Input that could be done tersely as follows,

Download calc1.mw

Download calc2.mw

Here are a few more ways, done using 2D Input. (It's a little awkward, since the 2D Input of  f'(-1) is hooked into the D form rather than the eval(diff(...)...) form.

Download calcexample.mw

And, using context-menus,

Download calcexample2.mw

Even in old Maple 17 the Help pages for radsimp states, "The original radsimp has been removed and replaced by a wrapper to simplify(...,radical,symbolic)."

If you don't want to operate under assumptions about the signum of the squared term under the radical then don't use the symbolic option (or radsimp).

Download radsimp_example.mw

You are using Maple 18 (and not Maple 2018 as your Question's header was incorrectly marked by you). That is why Tom's earlier suggestion did not work for you. I have corrected the Question's header.

Please don't submit duplicates queries for this problem. Post your followups here.

(You should also figure out how to search Mapleprimes -- since I suspect that previous students of your course/supervisor have been tasked with very similar examples and asked here before.)

You can change the name of the color palette (assigned at the top of this attachment). There are a few other choices. (Perhaps see here too.)

Download lines_ac.mw

You can augment the ODE system with a new differential equation involving g(t) such that g(1) equals the integral of Y0 from 0 to 1. For example, diff(g(t),t)=y0(t) and initial condition g(0)=0.

For example,

restart;
ICs := {y0(0)=3}:
ODESys := {diff(y0(t),t)=y0(t)}:
F := dsolve(ODESys union ICs union {diff(g(t),t)=y0(t), g(0)=0},
            {y0(t),g(t)}, numeric, output=listprocedure):

G := eval(g(t),F):

G(1);
            5.15484531216027

Or you can numerically integrate a Y0 procedure using evalf(Int(....)). For example,

restart;
ICs := {y0(0)=3}:
ODESys := {diff(y0(t),t)=y0(t)}:
F := dsolve(ODESys union ICs,
            {y0(t)}, numeric, output=listprocedure):

Y0 := eval(y0(t),F):

evalf(Int(Y0(t),t=0..1));
            5.15484514547602

I believe that using the augmented ODE system is more often somewhat better (numerically speaking) than numerically integrating the result procedure.

You don't have to use the listprocedure option, however (though I prefer it). There are other valid alternatives for picking off y0(t).

You could also do it with just a slight modification of your original code, using operator-form for the integrand and a floating-point end-point for the range of integration. (I consider this way of picking off the Y0 procedure to be fragile and inferior, and I also prefer my earlier syntax and methodology for the numeric integration.) For example,

restart;
ICs := {y0(0)=3}:
ODESys := {diff(y0(t),t)=y0(t)}:
F := dsolve(ODESys union ICs, {y0(t)}, numeric):

Y0 := t -> rhs(op(2, F(t))):

int(Y0,0..1.0);
            5.15484514547602

Download sqr_ex.mw

Naturally, the robustness of any technique can depend on the examples provided up front. I don't know what other examples you might have.

I am supposing that you really do want the factor of 1/2 distributed throughout the sum in brackets in the result, otherwise you could take combine(subexp) which returns with the 1/2 factored out front.

Download combine_trig_ex.mw

note. The order of the additive terms in the bracket depend only on the order in which they appear when subexp2 is first formed.

I could have used `*` instead of op(0,subexp) (since I can see that it is so) but I wanted it more generally programmatically valid and less ad hoc. Ie,

`*`(combine~([selectremove(type,subexp,trig)])[]);

Alternatives -- getting more precise -- might include,

subsindets(subexp,`*`,
           u->`*`(combine~([selectremove(type,subexp,trig)],trig)[]));

subsindets(subexp,And(`*`,satisfies(u->nops(select(type,u,trig))>1)),
           u->`*`(combine~([selectremove(type,subexp,trig)],trig)[]));

Download fit_solmodule.mw

There are various "Paragraph Styles" which you can modify.

There is a "Paragraph Style" names "Maple Output", whose alignment property can be changed from "Center" to "Left".

Here's what you do:

 - On the main menubar, Format -> Styles...

 - In the popup, left combo-box, select "P Maple Output", then click the "Modify" button

 - Change the "Alignment" dropmenu from "Center" to "Left"

Then execute and compute in your Document, and output should appear left-aligned.

You can even save an empty Document with this change as a so-called style sheet, and make it your default for New Documents.

See in the Help system, here and here. 

You can use the right-click Unit Formatting facility, or the Units:-UseUnit command.

In the actual Maple GUI the units appear in upright roman without the double-bracing. They appear that way when the worksheet is inlined here on Mapleprimes because that's the "old" way that they used to appear in Maple.

Download unitpref.mw

Starting with a RealRange, you can convert to relations (inequalities) and use the solve command. Union and intersection can be expressed pretty well using inert And and Or.

You can use Open ends, or not, and get strict or nonstrict inequailities accordingly.

Of course, there is nothing stopping you from beginning with inequalities. You don't have to begin from RealRange structures.

Depending on how you express the variables parameter to solve, the result can be in RealRange or relation form.

And so on.

Download solveRealRange.mw

You have used both y(t) and y, and it complains.

Perhaps you intended this:

eq1 := diff(y(t), t) = 1/2*cos(t - 2*y(t));
with(DEtools);
DEplot(eq1, y(t), t = -2 .. 2, y(t) = -2 .. 2);

If you want to substitute for diff(x(t), t, t) you don't have to assign to it.

(I find you intention slightly confusing, since it appears that the substitution is merely the equation from isolating `Eq__λ1`=0 for diff(x(t), t, t), hence the result is zero. If you intended something else then please explain.)

For example,

Download subst_example.mw

Outside the loop (right after calling dsolve) you can make this assignment,

   Y:=eval(y(t),lambdaE):

Then, replace all your calls lambdaE(tCompute[i]) by Y(tCompute[i]) in the later code.

The above works because you used the output=listprocedure option when you called the dsolve command.

If you had omitted the output=listprocedure option then you could get avoid the assignment to Y above, and instead replace all your calls lambdaE(tCompute[i]) by eval(y(t),lambdaE(tCompute[i])) in the subsequent code.

I prefer the way with listprocedure, but that's a stylistic coding preference.

[edit] I could add that it can be needlessly inefficient to call lambdaE(tCompute[i])  (or whatever you replace it with) multiple times for the same numeric value of tCompute[i]. Instead, inside your loop, you could simply call it once and assign the result to a temporary variable.
    temp := lambdaE(tCompute[i]);  # or whatever
And then use temp in the multiple locations your formulas.

Here are two ways.

The first has Nu[a] in italics (because I used single left-quotes to make it a name).

The second has Nu[a] in upright roman.

I specified the axes font, only so that the labels renders more clearly on this site. Adjust as you wish.
 

Download indexedname_label.mw

[edit] A couple of simpler ways (square brackets matching the bold, which you are free remove)
   labels=[eta,Typesetting:-mtext("Nu[a]")(eta)]
and,
   labels=[eta,`#mtext("Nu[a]");`(eta)]