This is actually a visual representation of a Riemann surface for Ln (= multi-valued log)

col:=[green, cyan, red, maroon, grey]:
plots:-display( seq( plot3d([ r*cos(u), r*sin(u), u],
u=-Pi+2*k*Pi .. Pi+2*k*Pi, r=0..1,  labels=[Rez,Imz, 'Im(Ln(z))'], grid=[120,30], color=col[k+3], lightmodel=none), k=-2..2));

A few considerations.

The equality 
arctan(z/b)/b = int(1/(b^2+u^2), u = 0 .. z)     (1)

is well known when z, b are real and b<>0 (note that I changed the dummy variable, it's not a good idea to use the same letter here).
Your substitution does not make much sense.

Now about the complex case in (1), we must be careful with complex integration. Maple is used for formal (symbolic) computations but it is supposed that the user knows the maths behind these computations.

A complex integral such as (1) is usually a "complex path integral". The integrand is analytic (holomorphic} in C \ {+-I*b}. If the path is missing (as here), it usually means that the integral is independent of the path (i.e. depends only on the endpoints 0 and z). This is true (for our analytic function) if it is considered in a simply connected domain containing 0 and z but not containing the poles +-I*b. So, we must restrict the computations to such domains; for example, C \ {+-I*b} is not simply connected. Ignoring this aspect leads to wrong results or nonsenses.
The conclusion is that you must have knowledge in complex analysis to address such problems.

restart;
Vol:=(u,v,w,U,V,W)->sqrt(-u^2*U^4 + ((-W^2 + u^2 + v^2)*V^2 + (u^2 + w^2)*W^2 - u^4 + (v^2 + w^2)*u^2 - v^2*w^2)*U^2 - v^2*V^4 + ((v^2 + w^2)*W^2 + (v^2 - w^2)*u^2 - v^4 + v^2*w^2)*V^2 - W^2*(W^2*w^2 + (v^2 - w^2)*u^2 - v^2*w^2 + w^4))/12:
m:=20; nr:=0: nrmax:=50;
for U to m do
for V from U to m do
for W from V to U+V-1  do
for u to m do
for v from abs(W-u)+1 to W+u-1  do
for w from max(abs(W-u),abs(U-v))+1 to min(W+u,U+v)-1  do
  vv:=Vol(u,v,w,U,V,W);
  if vv::posint then nr++; lprint(nr, [u,v,w],[U,V,W],vol=vv)  fi;
  if nr>=nrmax then break 6 fi
od:od:od:od:od:od:

# Your example appears at the position 441  (for m=20)

                            m := 20

                          nrmax := 50

1, [8, 8, 7], [2, 3, 4], vol = 6
2, [17, 16, 16], [2, 3, 4], vol = 15
3, [6, 7, 8], [2, 4, 4], vol = 6
4, [6, 8, 7], [2, 4, 4], vol = 6
5, [7, 5, 6], [2, 4, 4], vol = 6
6, [7, 6, 5], [2, 4, 4], vol = 6
7, [4, 8, 7], [2, 5, 6], vol = 6
8, [7, 4, 4], [2, 5, 6], vol = 6
9, [16, 14, 13], [2, 5, 6], vol = 12
10, [4, 6, 5], [2, 7, 8], vol = 6
11, [6, 4, 4], [2, 7, 8], vol = 6
12, [8, 4, 3], [2, 7, 8], vol = 6
13, [12, 14, 13], [2, 7, 8], vol = 24
14, [10, 11, 12], [2, 8, 8], vol = 21
15, [10, 12, 11], [2, 8, 8], vol = 21
16, [11, 9, 10], [2, 8, 8], vol = 21
17, [11, 10, 9], [2, 8, 8], vol = 21
18, [17, 11, 12], [2, 8, 8], vol = 21
19, [17, 12, 11], [2, 8, 8], vol = 21
20, [8, 12, 11], [2, 9, 10], vol = 21
21, [11, 8, 8], [2, 9, 10], vol = 21
22, [8, 10, 9], [2, 11, 12], vol = 21
23, [10, 8, 8], [2, 11, 12], vol = 21
24, [17, 8, 8], [2, 11, 12], vol = 21
25, [12, 8, 7], [2, 13, 14], vol = 24
26, [16, 6, 5], [2, 13, 14], vol = 12
27, [16, 22, 21], [2, 13, 14], vol = 60
28, [17, 3, 4], [2, 16, 16], vol = 15
29, [17, 4, 3], [2, 16, 16], vol = 15
30, [12, 19, 19], [2, 19, 19], vol = 72
31, [8, 24, 23], [2, 19, 20], vol = 42
32, [17, 16, 17], [3, 5, 7], vol = 30
33, [4, 7, 8], [3, 6, 6], vol = 9
34, [4, 8, 7], [3, 6, 6], vol = 9
35, [14, 16, 17], [3, 6, 6], vol = 36
36, [14, 17, 16], [3, 6, 6], vol = 36
37, [4, 6, 6], [3, 7, 8], vol = 9
38, [7, 11, 10], [3, 7, 8], vol = 24
39, [8, 4, 2], [3, 7, 8], vol = 6
40, [11, 7, 8], [3, 7, 8], vol = 24
41, [16, 16, 17], [3, 7, 8], vol = 48
42, [16, 8, 10], [3, 8, 10], vol = 21
43, [11, 15, 16], [3, 9, 9], vol = 42
44, [11, 16, 15], [3, 9, 9], vol = 42
45, [7, 8, 7], [3, 10, 11], vol = 24
46, [9, 15, 14], [3, 11, 12], vol = 48
47, [8, 14, 13], [3, 13, 14], vol = 48
48, [9, 12, 11], [3, 14, 15], vol = 48
49, [9, 17, 15], [3, 15, 16], vol = 42
50, [11, 9, 9], [3, 15, 16], vol = 42

Most "complex plot" commands are very simple, but are convenient. Basically they just take the real and imaginary parts and call a regular plot.
For example, for f:=sin(x+I):

complexplot(f,x=-Pi..Pi);

is equivqlent to:

plot([Re(f),Im(f),x=-Pi..Pi]);

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