It is difficult, and often impossible, to use an elementwise operator when one of the operands is a container (set, list, table, or rtable) that you want to be treated as a solitary operand.

Here's one more way, which I present as a curiosity rather than as a practical solution:

eval([{1,2}, {2,3}] minus~ S, S= {2})

My favorite way has already been presented by Joe:

map(`minus`, [{1,2}, {2,3}], {2})

Not only is this a bit briefer than

map(u-> u minus {2}, [{1,2}, {2,3}])

it's also more efficient. Keep this is mind whenever you want to map a multi-argument function or binary operator.

You need to define f in your code. From your related Question, "Decrypt Message", I know that f can be defined via

f:= a-> x-> a*x mod 256;

Like this:

for i to 3 do f[i]:= subs(_i= i, x-> P(_i,1,1,x)) od;

Here's a procedure to estimate the order of convergence. The first argument is the converging sequence as a list. The points in the list could be real numbers, complex numbers, vectors, or even more-complicated structures.

OrderOfConvergence:= proc(S::list, {Norm::appliable:= abs}, {diff::appliable:= `-`})
local E:= remove(`=`, map(Norm, map(diff, S[..-2], S[-1])), 0.);
   Statistics:-PowerFit(E[..-2], E[2..])[2]
end proc:

Here's an example of using the procedure to estimate the order of convergence of using Newton's method to approximate sqrt(2):

f:= x-> x - (x^2-2)/(2*x):
Digits:= 200:
S[1]:= 1.:
for k to 9 do S[k+1]:= f(S[k]) od:
ord:= OrderOfConvergence(convert(S,list)):
evalf[10](ord);
                          2.000367876

If the sequence is a sequence of vectors, then the second argument should be LinearAlgebra:-Norm:

OrderOfConvergence(S, LinearAlgebra:-Norm);

0 = ln(1) = ln((-1)*(-1)) = ln(-1) + ln(-1) = 2*I*Pi ?  So, clearly the rule ln(x) + ln(y) = ln(x*y) is not true for all x and all y. However, if the rule is incorrectly applied, the most that you'll be wrong by is an integer multiple of 2*I*Pi. This is mentioned in the Wikipedia article that you linked to, if you read further.

You should show the command that you used for "what I get".

This is pretty close to your desired fieldplot:

plots:-fieldplot([-(x+1)*(x-3)*y, x*(y+2)*(y-1)], x= -1..3, y= -2..1);

Unfortunately, implicitplot3d is one (and I believe the only) plot command from which one cannot extract the data. But your equation can be easily solved for f0 or sigma, each yielding a two-branched solution. That's the way to go.

Any map R^2 -> R^2 can be viewed as a map C -> C and thus plotted with plots:-conformal. It is irrelevant that the complex map is neither conformal nor analytic; the plot command doesn't care.

restart:
#Convert any map R^2->R^2 to the equivalent C->C:
ToComplex:= Z-> unapply(inner([Z((Re,Im)(_z))],[1,I]), _z):
Z:= ToComplex((r,t)-> (r*exp(t), r*exp(-t)));
plots:-conformal(Z, 0..1+I, labels= [r,t], grid= [13,13], numxy= [51,51]);

Regarding linspace: It is now not needed because its purpose is handled by the grid option. However, it is a generally useful command from Matlab. It can be implemented in Maple by

linspace:= proc(a,b,n::And(posint, Not(1)))
local k, h:= (b-a)/(n-1);
   seq(a+k*h, k= 0..n-1)
end proc:

A more Mapleish implementation uses ranges:

linspace:= proc(ab::range(algebraic), n::And(posint, Not(1)))
local k, a:= op(1,ab), h:= (op(2,ab)-a)/(n-1);
   seq(a+k*h, k= 0..n-1)
end proc:

So it can be called via

linspace(0..1, 11);

It is possible to solve your system symbolically and obtain 7 explicit solutions for p, v, and w in terms of t that are not terribly complicated:

restart:
eq1:= 3*v^2-v*t*(v^2+3)-3*w^2+w*t*(w^2+3)=0:
eq2:= v*(p*v^2+3)/(3*v^2+1)-w*(p*w^2+3)/(3*w^2+1) = 0:
eq3:= -v*(p*v^2+3)/(3*v^2+1)-w*(p*w^2+3)/(3*w^2+1) = 0:
_EnvExplicit:= true:
Sols:= [eliminate({eq1,eq2,eq3},{p,v,w})]:
simplify(Sols);

Use a parametric plot3d for the plane:

plot3d([u,u,v], u= 0..0.01, v= 0..100, transparency= .8);

Here's another example. As Preben noted, the key point is to include the output= Array(...) option to dsolve, where the array contains the independent-variable values. In this example, I FFT the famous Lorenz Attractor.

restart:
N:= 12: #Use 2^N t-values.
(Tmin,Tmax):= (0,20):
xyz0:= [1,1,1]:
sys:= {
   diff(x(t),t) = sigma*(y(t) - x(t)),
   diff(y(t),t) = x(t)*(rho - z(t)) - y(t),
   diff(z(t),t) = x(t)*y(t) - beta*z(t),
   ([x,y,z](Tmin) =~ xyz0)[]
}:
params:= [sigma= 10, rho= 28, beta= 8/3]:
h:= (Tmax-Tmin)/(2^N-1):
T:= Array(1..2^N, k-> Tmin + h*(k-1)):
Sol:= dsolve(eval(sys, params), numeric, maxfun= 0, output= T);
plots:-odeplot(Sol, [x,y,z](t));
#When you modify this example for your system, do NOT change ANY of the numbers 
#in the next line.
XYZ:= DiscreteTransforms:-FourierTransform(Sol[2,1][.., 2..], 1); 
for P in [Re,Im] do
   print(plots:-pointplot3d(P~(XYZ), style= line))
od;

Yes, to solve systems of polynomial equations/inequalities over the reals, you should use RegularChains, not Groebner. The RegularChains package is massive, perhaps the largest in Maple, but you probably don't need to access it directly. There are several user-friendly access points from outside the package. I think that the most appropriate in your case is SolveTools:-SemiAlgebraic.

Please post your system of equations/inequalities using plaintext Maple code, not MathJax. I would've tried solving your system, but I don't fully understand it.

To differentiate with respect to a function (such as another derivative), use Physics:-diff:

y:= diff(tau(t), t)^2:
Physics:-diff(y, diff(tau(t),t));

FirstHit:= proc(L::{Array,list}(float), C::realcons)
local c:= evalf(C), k;
   for k to numelems(L) do
      if L[k] > c then return k fi
   od;
   FAIL
end proc:

What you have written is a fine answer for part b of the problem, not part a. To answer part a, simply change plot to line in your Tangent command.

The animation has two parts: the part that changes with each frame (the tangent line) and the part that doesn't (the function's plot). The part that doesn't change is in the background option below:

f:= x-> 2/(1+exp(-x)):
plots:-animate(
   plot,
   [
       Student:-Calculus1:-Tangent(f(x), x= x0, output= line),
       x= -6..6,
       color= "DarkGreen", thickness= 3
   ],
   x0= -5..5,
   background= plot(f, -6..6, thickness= 3)
); 

To play the animation, you need to click on it and start it with the toolbar controls or the context menu.