If you are using Maple input, and the code is all in one input block, then you can use the multi-line Maple comment, (* ... *)

> (*
>    These lines are commented out.
> *)

This won't be the easiest approach, but it may do the job.  I've arbitrarily assigned expressions for f, g1, and g2 in terms of x, y1, and y2.

alias(x=x(y1,y2)):

f := x^2 + y1/y2:
g1 := x*y1:
g2 := x+y2:

eq1 := f = 0:

d1 := solve(diff(eq1,y1), {diff(x,y1)});
                                   d          1
                           d1 := {--- x = - ------}
                                  dy1       2 x y2

d2 := solve(diff(eq1,y2), {diff(x,y2)});
                                    d        y1
                            d2 := {--- x = -------}
                                   dy2           2
                                           2 x y2

d3 := solve(diff(g1,y1)=0, {diff(x,y1)});
                                    d         x
                            d3 := {--- x = - ----}
                                   dy1        y1


d4 := diff(subs(d1,d3),y2);
                                 d                   d
                                --- x               --- x
                                dy2        1        dy2
                     d4 := {1/2 ----- + ------- = - -----}
                                 2            2      y1
                                x  y2   2 x y2

d5 := frontend(solve, [d4, diff(x,y2)], [{Not(specfunc(anything,diff))},{}]);
                              d              x y1
                      d5 := {--- x = - -----------------}
                             dy2              2
                                       y2 (2 x  y2 + y1)

d6 := {diff(g2,y2)=0};
                                   / d   \
                            d6 := {|--- x| + 1 = 0}
                                   \dy2  /


d7 := subs(d5, d6);
                                     x y1
                      d7 := {- ----------------- + 1 = 0}
                                      2
                               y2 (2 x  y2 + y1)

frontend(solve, [d7,y2], [{Not(identical(x))},{}]);
                      2      3    1/2                   2      3    1/2
             -y1 + (y1  + 8 x  y1)              y1 + (y1  + 8 x  y1)
       {y2 = ------------------------}, {y2 = - -----------------------}
                          2                                 2
                       4 x                               4 x

Alternatively, you can use declare the type of the parameter.  While that may not always be possible, it is generally more convenient to let Maple do the parameter processing.  For example

f := proc(x::integer) x^2; end proc:
f(1/2);
Error, invalid input: f expects its 1st argument, x, to be of type integer, but received 1/2

Here's a simple helper procedure that implements Robert's suggestion

p1 := plot(sin):
p2 := plot(cos):

recolor := proc(p, color)
    subsindets(p
               , 'specfunc(anything, {COLOR, COLOUR})'
               , clr -> COLOUR(RGB,op(ColorTools:-NameToRGB(convert(color,string))/255.0))
              )
end proc:

plots:-display(recolor~([p1,p2],[red,blue]));

Did you try the help page?  As an example,
 

taylor(x^3+1, x=2);
                                              2          3
                    9 + 12 (x - 2) + 6 (x - 2)  + (x - 2

Followup:  taylor returns a Maple type series.   Use convert/polynom to convert back to a polynomial.

convert(%,polynom);
                                            2          3
                      -15 + 12 x + 6 (x - 2)  + (x - 2)

Note that Maple automatically expands the linear term.

expand(%);
                                     3
                                    x  + 1

Jakubi's approach is simpler, but this may be of interest.  The basic idea is to use simplify/polar; alas, it does not currently handle assumptions so we will surgically repair it (temporarily).  Here's the limitation

simplify(polar(r,theta)^k) assuming k::real;
                                           polar(r,theta)^k

Here's the minor surgery

`simplify/polar` := subs('{numeric, realcons}' = 'satisfies'(rcurry(is,real))
                         , eval(`simplify/polar`)
                        ):

Now

simplify(polar(r,theta)^k) assuming k::real;
                                                  polar(r^k,theta*k)

With that fixed, we can now convert to polar form and use the new simplification.

zr := a+b*I;
zp := polar(zr);
zc := conjugate(zp);
y1 := zp^k + zc^k;
simplify(y1,polar) assuming real;
evalc(%);
                 2*((a^2+b^2)^(1/2))^k*cos(arctan(b,a)*k)

With Maple 13 you can append the tilde to a symbol and it will map over elements.  Thus

fg~(L);

would also work.  See ?operators/elementwise.

What do you expect to get?  Does p actually represent a general expression in a(t) and b(t)?

Consider

p := (a,b) -> a^2 + sin(b):
y := p(a(t),b(t)):

The frontend procedure can be used here to compute the derivative with respect to a(t):
 

frontend(diff, [y,a(t)]);
                                                    2 a(t)

Actually, that only worked because the sin didn't contain a(t). Consider
 

y := a(t) + b(t)*sin(a(t)+b(t)):

To use frontend effectively here you need to use its optional third argument that specifies what not to freeze.

frontend(diff, [y, a(t)], [{`*`,`+`, And(function,Not(identical(a(t),b(t))))},{}]);
                                          1 + b(t) cos(a(t) + b(t))

Suppose that your expression contains an inert function of a(t) and b(t).

frontend(diff, [f(a(t),b(t)), a(t)], [{`*`,`+`, And(function,Not(identical(a(t),b(t))))},{}]);
                                          diff(f(a(t), b(t)), a(t))

%;
Error, invalid input: diff received a(t), which is not valid for its 2nd argument

That raises an error because the diff procedure cannot deal with a function as the second argument. The proper way to express the derivative is with the D operator.  That can be done by using  convert/D inside the call to frontend,

frontend(`convert/D`@diff, [f(a(t),b(t)), a(t)], [{`*`,`+`, And(function,Not(identical(a(t),b(t))))},{}]);
                 D[1](f)(a(t), b(t))

b[1] := b1:
b1_func := unapply(b[1]+t,t):
b1_func(x);
                                    b1 + x

b[1] := b2:
b1_func(x);
                                    b1 + x

Another approach, less efficient than nested seq's but possibly easier to comprehend, is to use nested loops and store the intermediate result in a table that you later convert to a set. To reduce the size I'll only save ordered pairs.

pairs := table():
cnt := 0:
for i from 100 to 110 do
    for j from i+1 to 110 do
        if igcd(i,j) = 1 then
            cnt := cnt+1;
            pairs[cnt] := [i,j];
        end if;
    end do;
end do:

convert(pairs, set);

I assume regexprep is your own procedure.  Why are you passing true/false to ArrayTools:-AllNonZero?

Followup: I'm slow today. I hadn't realize that was the output of FromMFile.  The problem there is that the Matlab local variable fi is a Maple reserved word. Change it, say to ffi.  That won't fix the other issues with this translation...

Maybe your expectation of what Eval does is incorrect?  Eval is an inert version of eval. Contrast

Eval(x,x=1);
                                                      x|
                                                       |x = 1

eval(x,x=1);
                                                         1

Oh, now I understand. You want to show the unevaluated version of your sin. First off, I don't believe that your technique will work, at least not if done at top-level.  The reason is that, at top-level, :-sin is the same as sin, so your assignment leads to an infinite recursion.  It will work if done inside a procedure or module.  Here I'll use the name Sin which avoid the issue.

Sin := proc(x)
    if x :: numeric then
        :-sin(x*Pi/180);
    else
        'procname'(x)
    end if;
end proc:

Eval(y=Sin(x), x=30);
                              (y = Sin(x))|
                                          |x = 30

value(%);
                                    y = 1/2

A disadvantage with this method is that it doesn't allow, say, Sin(44) to exist unevaluated.  An alternative is to leave Sin inert but extend value so that it know how to evaluate it.

restart;
`value/Sin` := x -> sin(x*Pi/180):
Eval(y = Sin(x), x=30);
                              (y = Sin(x))|
                                          |x = 30

value(%);
                                    y = 1/2

Sin(30);
                                    Sin(30)

value(%);
                                      1/2

Another possibility is to use %Sin as the inert version:

restart;
Sin := proc(x)
    if x :: realcons then
        :-sin(x*Pi/180);
    else
        'procname'(x)
    end if;
end proc:

Eval(y = %Sin(x), x=30);
                             (y = %Sin(x))|
                                           |x = 30

value(%);
                                    y = 1/2

That works because value converts %XXX to XXX

That method will not work.  The reason is that the ideal op-amp is modeled as one would by hand, that it, it impresses a virtual-short between its two input terminals.  The only possible solution for this circuit with that constraint (V+ = V-) is to have V+ = V-  = 0, which is not oscillatilng.

If you substitute the ideal op-amp with limited gain and set the voltages sources to something reasonable, and set an initial condition on the capacitor, you can see the output start to move, but the model will hit a singularity and terminate. Proving that that must be the case would make a nice exercise.

What you really want to use is a comparator rather than an op-amp.  I modeled one by using a VCVS driving a controlled commutating switch, with one pole connected to ground and the other to a positive supply.  The common terminal is the output.  That does oscillate.

Do you want the symbols (or names) from a particular expression, or all expressions in the worksheet?

The command for extracting symbols or names (or other types) from an expression is indets.  Check out its help page, it takes an optional second argument that specifies the type. 

Collecting names from the worksheet as a whole is a bit trickier. The commands of use are anames and unames.

Use the normal integrator and divide the result by the time. That will give an average power that is continuously updated.  You may need to avoid dividing by zero at t=0, that can be done in a few ways.

Another possibility is to make the time period a parameter and then use it to compute the average as well as set the period of the driving sources (assuming that there is a driving source that sets the period).