As you no doubt know, the left side of your identity is Zeta(s). Since Maple also knows this, I'll use Zeta(s) for it for the rest of this article.

Constraints can be placed on variables by using an assuming clause, as in

is(Zeta(s) = product(ithprime(j)^s/(ithprime(j)^s - 1), j= 1..infinity)) assuming s > 1;

However, Maple doesn't have enough symbolic knowledge of ithprime to be able to solve this.

You can see all the symbolic knowledge (such as identities) that Maple knows about Zeta by doing

FunctionAdvisor(Zeta(s));

As far as I know, Maple has no symbolic knowledge about ithprime.

We can explore your identity through numeric computation like this (I've taken the logarithm of both sides of the identity):

nP:= 999: #number of primes to test
Ps:= hfarray([seq(ithprime(j), j= 1..nP)]):
F:= s-> evalhf(add(s*ln(Ps[k])-ln(Ps[k]^s-1), k= 1..nP)):
plot(ln@Zeta - F, 2..9);

The plot shows that the difference of the two sides of the identity is close to 0 when the infinite product is truncated.

A wrote a module to do this. It may be possible to improve the speed on this by using compiled code or by using evalhf differently. (That's not to say that this code is particularly slow!)

Pi_BaileyBorweinPlouffe:= module()
#Prints hexadecimal digits of Pi starting at an arbitrary position
local
   ModuleApply:= proc(d::posint)
   #d is the starting digit position. 9 digits will be printed because
   #that's the limit of hardware-float arithmetic according to the paper.
   option hfloat;
   local x, s, k, dm1:= d-1;
      Digits:= trunc(evalhf(Digits));
      x:= Frac(4*SSj(dm1,1) - 2*SSj(dm1,4) - SSj(dm1,5) - SSj(dm1,6));
      Digits:= 36;
      s:= sprintf("%38.36f", convert(x, binary));
      cat(seq(Hex(s[4*k-1..4*k+2]), k= 1..9))
   end proc,
      
   SSj:= proc(d::nonnegint, j::identical(1,4,5,6))
   option hfloat;
   local k;
      Frac(
         # *1
         Frac(`if`(d=0, 1./j, add((16 &^ (d-k) mod (8*k+j))/evalf(8*k+j), k= 0..d))) 
         + 
         evalhf(add(16^(d-k)/(8*k+j), k= d+1..d+14))
      )
   end proc,

   Hex:= proc(x::string)
   #Converts a string of exactly 4 binary digits to a hexadecimal digit
   option remember;
      convert(convert(x, decimal, 2), hexadecimal)
   end proc,

   Frac:= proc(x::numeric)  # *2
   option hfloat;
   local f:= x - trunc(x);
      `if`(f >= 0, f, 1+f)
   end proc
;
end module:

(* Footnotes:

*1: d=0 needs to be treated as a special case due to a bug in `mod/&^` such that
    16 &^ 0 mod m 
 returns an error. 

*2: Maple's endogenous `frac` won't work here because it may return negative. 
*)

Example usage:

Digits:= 15:
Pi_BaileyBorweinPlouffe(1);
                           243F6A888
evalf(convert(%, decimal, 16)/16^9);
                       0.141592653584667
frac(evalf(Pi));
                        0.14159265358979
Pi_BaileyBorweinPlouffe(1536);
                           362FB1341

The system and ssystem commands do the same thing as !, but they give you more flexibility for situations like this.

system(
   "curl -m 3 -o /Users/John/Documents/P\\ Admin/a.dat\"
   " --insecure https://www.google.com/?gws_rd=ssl"
);

Note carefully how I split the string over lines.

You had a semicolon at the end. I don't know whether this is required by MacOS or if you included it because you thought that it was required by Maple. If MacOS requires it, then put it at the end of the string above.

If you need the result of this command returned to Maple (as a string), then use ssystem instead.

Here's an example of generating random integers in a do loop and accumulating their lowest common multiple:

LCM:= 1:
to 9 do
   LCM:= ilcm(LCM, rand())
end do:
LCM;

The lowest common multiple is an associative operation. The technique above will work for any associative operation, and is commonly used for sums. The initial value (in this case, the 1) must be the identity element for the operation.

Having a variable, a name, that is assumed to be integer is not the same as having a value that is actually an integer (i.e., has type integer). The command numtheory:-bigomega wants an actual integer.

Now let's learn about premature evaluation, which is the essential difference between

sum(numtheory[bigomega](k), k= 2..9)

and

add(numtheory[bigomega](k), k= 2..9)

The command sum is like the vast, vast majority of other Maple commands: It evaluates its arguments before they are passed. That means that the above sum command attempts to evaluate numtheory[bigomega](k) before k has been assigned any numeric value. That causes an error.

On the other hand, the commands seq, add, and mul have special evaluation rules: Their index variable is assigned values from the prescribed range before the first argument is evaluated.

To some extent, premature evaluation can be controlled with unevaluation quotes; however, those are often difficult or impossible to use. In the above case, sum('numtheory[bigomega](k)', k= 2..9) works correctly, but it's better to just use add.

When you get an error message that says that a certain argument is expected to be a certain type, that always means that that argument must actually be that type; that error can never be corrected by making assumptions.

Your definition of N is correct. It's obvious that the midrange is approximately 10 and that the supposed answer is way off.

To get the min and max of a 2D plot, first assign the plot to a variable:

P:= plot(N, 0..20);

then

(min,max)(op([1,1], P)[.., 2]);

      9.49146814577784781, 11.3050481854892997

You can simply do

1 +~ ListTools:-PartialSums([1, 3, 4, 50, 10]);

You can't use the same variable, i, as both the index for the for loop and the index for the sum command. Simply use a different variable.

There are several other improvements that you could make, which I'll detail in a Reply.

Do

Matrix((2,2), (i,j)-> eval(p[j], solution[i]));

a) If you have actual data points, you should use Statistics:-PowerFit. It'll do the log transformation and find the closest exponent for you.

b) If R = a*P^(2/3), then R/P = a*P^(-1/3). From there, the two limits are trivial.

c) You have a lot of questions. Are you trying to do your homework for the entire semester in one night? You should indicate/show that you understand the answers that I already gave before proceeding.

Like this

restart:
K:= 100:  N__0:= 20:
NP:= R-> rsolve({N(t+1) = R*N(t)/(1+(R-1/K)*N(t)), N(0)= N__0}, N(t), makeproc):
plot(
   [seq(['[k, NP(R)(k)]' $ k= 1..10], R= [2, 5, 10])],
   legend= (R=~ [2, 5, 10]),
   labels= [t,N]
);

If it appears as a straight line on a log-log plot, then the equation must be as I've given as Meq, for some yet-to-be-determined values of k and C. The two given data points are enough to solve for these.

restart:
Meq:= t = exp(k*ln(T)+C):
S:= solve({eval(Meq, [t=3, T=20]), eval(Meq, [t=20, T=5])});  
plots:-loglogplot(eval(t, eval(Meq, S)), T= 1..30);

Do

restart:
R:= k*(a-x)*(b-2*x)^2:
a:= 5:  b:= 6:  k:= 0.3:
plot(R, x= 0..min(a,b/2), labels= [x,'R']);

You can also plot x as a function of time like this:

Sol:= dsolve({diff(x(t),t) = subs(x= x(t), R), x(0)=0}, numeric):
plots:-odeplot(Sol, [t,x(t)], t= 0..min(a,b/2));

In Maple's plaintext input (aka 1D input or Maple Input), the juxtaposition of parenthesized expressions does not imply their multiplication. In other words (a-1)(a-2) doesn't mean (a-1)*(a-2). Unfortunately, it does mean something, so it's not a syntax error. What it means is rather obscure and not worth discussing here.

Do you really need a Groebner basis? Or do you just need a solution? The equations can be easily solved by solve. Here, I just use equations 17-44, because I didn't feel like correcting the first 16, and it's trivial to see that the given solution satisfies the first 16.

solve({e||(17..44)});
      {a = 2, b = 3, c = 4, d = 1, e = 1, f = 4, g = 3, h = 2, i = 3, j = 2, k = 1, l = 4, m = 4, n = 1, o = 2, p = 3}

Here's a pair of procedures for it. They'll work in any base, although the results will always be reported in base 10. The second parameter, n, of KaprekarSeq is the maximum number of terms of the sequence to report. This is needed because it's not guaranteed that the sequence will converge to a fixed point. It's likely that it will, but not guaranteed.

KaprekarSeq:= proc(x::posint, n::posint, b::posint:= 10)
local j, r, kx;
   r[0]:= x;
   for j to n-1 do 
      kx:= Kaprekar(r[j-1], b);
      if kx = r[j-1] then break end if;
      r[j]:= kx
   end do;
   convert(r, list)
end proc:

Kaprekar:= proc(x::posint, b::posint:= 10)
local X:= convert(x, base, b), N:= b^nops(X), S:= sort(X);
   convert(S, base, b, N)[] - convert(ListTools:-Reverse(S), base, b, N)[]
end proc:

Example usage:

KaprekarSeq(1234, 9);

       [1234, 3087, 8352, 6174]