restart;
i2 := (x,y) -> -(1/2)*I*(exp(I*x)*(sin(x)/x)-exp(I*y)*(sin(y)/y))/(x-y);
i3_r := -(1/2)*I*(i2(y,z)-i2(y,x))/(z-x);
normal(eval(series(eval(i3_r,[x=t*x,y=t*y,z=t*z]),t=0),t=1));

On my rather old computer, for any of your graphs (having 21 vertices and 76 edges), VertexConnectivity needs about 0.5 sec.
So, it is easy to estimate the total time for your almost 15000 graphs.

It should be easy to write a compiled version for VertexConnectivity. It will be much faster (and a skilled user can do it), but probably the designers did not anticipate that someone will need it for so many graphs.

`simplify/size` is responsible for this. It calls `simplify/size/size` to compute the "math complexity" (see ?simplify/size) and chooses the simpler one. 

ex:=[t^2*x^2+t^2*y^2, t^2*x^2-t^2*y^2]:
`simplify/size/size`~(%);
#         [56, 59]

simplify~(ex);
#         [t^2*(x^2+y^2), t^2*x^2-t^2*y^2]

`simplify/size/size`~(%);
#        [38, 59]

Note that simplify(..., size=false) keeps both expressions unchanged.

If you want the rational exponent of x^r, where x is a variable and r is rational (including r=1)
then r=1 must be treated separately.

Exponent := proc(x::{name, name ^ rational}) 
  if x::name then 1 else op(2, x) fi 
end proc:
Exponent(x^(-3/2));   # -3/2
Exponent(Y[1,3]);     # 1

Note that r=0 is not possible because x^0 is automatically simplified to 1.

remove removes operands (1st level) of an expression.
x has one operand (x itself)
diff(y(x), x) has two operands (y(x) and x)  [so, you can remove only y(x) or x or both].
2*diff(y(x), x) has two operands (2 and diff(y(x), x))

Wind:=proc(L::listlist, P::list, isclosed:=true)
  local p,u,k, x0:=P[1],y0:=P[2];
  if member(P,L) then return undefined fi;
  u:=seq(arctan(p[2]-y0,p[1]-x0),p=L),  `if`(isclosed=true,arctan(L[1][2]-y0,L[1][1]-x0), NULL);
  (u[-1]-u[1])/(2*Pi) + add(piecewise( u[k+1]-u[k]<-Pi,+1, u[k+1]-u[k]>Pi,-1, abs(u[k+1]-u[k])=Pi, undefined, 0), k=1..nops([u])-1)
end:

It returns 0 if the point is outside, and undefined if it is on the boundary.

L:=[[0,0],[200,0],[200,300],[250,400], [500,300],[500,500],[0,500]]:
Wind(L, [300,300]), Wind(L, [100,200]), Wind(L, [0,50]);
                        0, 1, undefined

N depends on the order in which the numbers are eliminated.
The largest value of N is N_max = 2^(n-1) * (3*n - 4) + 2 
(you have n = 2012).
Try your hand for the smallest value!

It's obvious that the smartview option cannot guess all user's intentions.
If you don't know the relevant y-domain (used by view=...) you may try a slope criterion. E.g.

restart;
f:=exp(-3^(1/2)*(cos(x)-1)/sin(x)):
f1:= diff(f, x):
maxslope:=100:
ff:=piecewise(abs(f1) < maxslope, f, undefined):
plot(ff, x=0..2*Pi);

restart;
Digits := 15: 
nu_0 := sqrt(k^2 - k_0^2): nu_m := sqrt(k^2 - k_m^2): 
pref2 := 2*rho_me*exp(-nu_0*(zs + z))*BesselJ(0, k*sqrt((xs - x)^2 + (ys - y)^2))*k/(rho_me*nu_0 + rho_0*nu_m*tanh(nu_m*d)):

f:= (eval(0.5*(((eval(pref2, [z = 0, k_0 = 0.018371886863098, xs = 0, ys = 0, zs = 10, rho_0 = 1.2, d = 0.04, rho_me = 1.528516816439260 - 1235.297048886680*I, k_m = 0.490806242885258 - 0.490314205803914*I])))), [z = 0, k_0 = 0.018371886863098, xs = 0, ys = 0, zs = 10, rho_0 = 1.2, d = 0.04, rho_me = 1.528516816439260 - 1235.297048886680*I, k_m = 0.490806242885258 - 0.490314205803914*I]) ):

dom := [k=0..100, x = -2/2 .. 2/2, y = -3/2 .. 3/2]:  # 100 is enough
evalf(Int(f, dom, digits=10));

                               0.5867621705 - 0.1097638086*I

The answer is correct, a being a name and non-constant.
eval(diff(y(x,t), x$3), x=Pi);  and  eval(diff(y(x,t), x$3), x=a+1);

will look as you expect.

You may declare a as a constant. In this case the result will appear as for Pi:
constants := constants, a;
eval(diff(y(x,t), x$3), x=a);
            

restart;
N:=100:  dx:=evalf(2*Pi/N):
X:=<seq(k*dx, k=-N/2..N/2)>:
Y:= sin~(X):
f:=unapply(CurveFitting:-Spline(X, Y, x, degree=1),x):
Z:= [seq( evalf(Int(f, X[1]..x, epsilon=1e-5)), x=X)]:
plots:-display(plot(X,Y),  plot(X,Z), color=[red,blue]); 

Use the definition:

P := (x-6)^2 + (y+3)^2 = (x+2*y-1)^2 / 5;

with(plots): display(implicitplot([P, x+2*y-1], x=-10..10,y=-10..10), pointplot([6,-3], color=red));

The singularities of F1 are on the Jordan curve

plots:-implicitplot(1/F1, x=0..k - varepsilon,y=0..2*Pi, view=0..6.5, color=red);

Note that the integrand is even in x, so, one may integrate over x = 0 ... 1 - eps.

Unfortunately Maple cannot compute numerically a CPV integral, but it is easy to transform it into a regular one. If the singular point is c, then 

Int(f(x), x =c-r .. c+r, CPV) = Int(f(c+t) + f(c-t), t = 0 .. r)

Strange implementation (bug?) for numboccur.

restart;
L:=[1., HFloat(1.), 1., 0., 0.];
#                   L := [1., 1., 1., 0., 0.]

numboccur(L, 1.);
#                               2

numboccur(L, 1.0);
#                               0

numboccur(L, HFloat(1.0));
#                               1

evalb(L[1]=L[2]);
#                              true

select(`=`, L, 1.); select(`=`, L, 1.0);
#                          [1., 1., 1.]

#                          [1., 1., 1.]

dismantle(%);

#LIST(2)
#   EXPSEQ(4)
#      FLOAT(3): 1.
#         INTPOS(2): 1
#         INTPOS(2): 0
#      HFLOAT(2): 1.
#      FLOAT(3): 1.
#         INTPOS(2): 1
#         INTPOS(2): 0