I guess, that Markiyan Hirnyk says: 

It is not a continous function in parameters in general.

@roman_pearce 

What does the "$" mean in the definition of the proc?

map(evalc, A):  # assumes, that variables are Reals
simplify(%): 
value(%): 
res:=simplify(%);

  res := abs(-1+a^2)*(signum(-1+a^2)^2+1)*csgn((a+1)/(a-1))*Pi/(-1+a^2)

That already ignores N and now proceed

convert(res, piecewise, a):  # silently assumes, that a is a Real 
simplify(%);

                       {  0           a = -1
                       {
                       {  0            a = 1
                       {
                       { 2 Pi        otherwise

map(evalc, A):  # assumes, that variables are Reals
simplify(%): 
value(%): 
res:=simplify(%);

  res := abs(-1+a^2)*(signum(-1+a^2)^2+1)*csgn((a+1)/(a-1))*Pi/(-1+a^2)

That already ignores N and now proceed

convert(res, piecewise, a):  # silently assumes, that a is a Real 
simplify(%);

                       {  0           a = -1
                       {
                       {  0            a = 1
                       {
                       { 2 Pi        otherwise

That sounds strange - I have no idea for that.

Please ask your local support, they may also check your installation
for user rights on the machine and configuration of the browser.

That sounds strange - I have no idea for that.

Please ask your local support, they may also check your installation
for user rights on the machine and configuration of the browser.

The following should give the result with correct relative
precision (for 'arbitrary' Digits):

S:=proc(a)
  local n, epsilon:= Float(1,-Digits), S0, Soo, S1, oldDigits;
  
  # Zeta(n) ~ 1 with relative error   
  n:= ceil( -ln(epsilon)/ln(2) ) + 3;
  
  # brute way to care for cancellation errors
  oldDigits:=Digits;
  Digits:=2*Digits;
   
  S0:=  evalf( Sum(Zeta(p-2)*a^p/GAMMA(1/2*p+1), p = 4 .. n) );
  
  # Sum(a^p/GAMMA(1/2*p+1), p = n+1 .. infinity) = Soo - S1
  S1:=  evalf( add(   1     *a^p/GAMMA(1/2*p+1), p = 0 .. n) );
  Soo:= evalf( (1+erf(a))*exp(a^2) );
  
  #print([S0, Soo, S1, Soo - S1]);  
  S0 + Soo - S1;
  return evalf[oldDigits](%);
end proc;

The following should give the result with correct relative
precision (for 'arbitrary' Digits):

S:=proc(a)
  local n, epsilon:= Float(1,-Digits), S0, Soo, S1, oldDigits;
  
  # Zeta(n) ~ 1 with relative error   
  n:= ceil( -ln(epsilon)/ln(2) ) + 3;
  
  # brute way to care for cancellation errors
  oldDigits:=Digits;
  Digits:=2*Digits;
   
  S0:=  evalf( Sum(Zeta(p-2)*a^p/GAMMA(1/2*p+1), p = 4 .. n) );
  
  # Sum(a^p/GAMMA(1/2*p+1), p = n+1 .. infinity) = Soo - S1
  S1:=  evalf( add(   1     *a^p/GAMMA(1/2*p+1), p = 0 .. n) );
  Soo:= evalf( (1+erf(a))*exp(a^2) );
  
  #print([S0, Soo, S1, Soo - S1]);  
  S0 + Soo - S1;
  return evalf[oldDigits](%);
end proc;

A:=map(convert, A, rational):
1/A: evalf[16](%): convert(%, list);

[4.587858207136615, -.9714476510691793, -2.036021455162406, 
  -.9561250876917093, -.5897670589620581, -.5183393016505515,
-.9714476510691793, 2.225705816897946, .9962413794921655e-2,
   -.6431880770231276, -.1283352063807487, -.1646368544191163,
-2.036021455162406, .9962413794921655e-2, 2.810538455018751,
   .8187140864632667e-2, -.3459189551452366, -.1888857297921770e-1,
-.9561250876917093, -.6431880770231276, .8187140864632667e-2,
   2.132711063728444, -.1288794640122821e-1, -.1285508973757944,
-.5897670589620581, -.1283352063807487, -.3459189551452366,
   -.1288794640122821e-1, 1.596704324983773, .3128566096228758e-1,
-.5183393016505515, -.1646368544191163, -.1888857297921770e-1,
   -.1285508973757944, .3128566096228758e-1, 1.392967817819466]

A:=map(convert, A, rational):
1/A: evalf[16](%): convert(%, list);

[4.587858207136615, -.9714476510691793, -2.036021455162406, 
  -.9561250876917093, -.5897670589620581, -.5183393016505515,
-.9714476510691793, 2.225705816897946, .9962413794921655e-2,
   -.6431880770231276, -.1283352063807487, -.1646368544191163,
-2.036021455162406, .9962413794921655e-2, 2.810538455018751,
   .8187140864632667e-2, -.3459189551452366, -.1888857297921770e-1,
-.9561250876917093, -.6431880770231276, .8187140864632667e-2,
   2.132711063728444, -.1288794640122821e-1, -.1285508973757944,
-.5897670589620581, -.1283352063807487, -.3459189551452366,
   -.1288794640122821e-1, 1.596704324983773, .3128566096228758e-1,
-.5183393016505515, -.1646368544191163, -.1888857297921770e-1,
   -.1285508973757944, .3128566096228758e-1, 1.392967817819466]

I think your c(k) writes as

H := k -> 2*3^k/GAMMA(k+1) * 
  Int(exp(-14/3*t)*(-GAMMA(k+1,7/3*t)+GAMMA(k+1,4/3*t)),t = 0 .. infinity)

I currently miss familiarity with the Gamma distribution to simplify more.

I think your c(k) writes as

H := k -> 2*3^k/GAMMA(k+1) * 
  Int(exp(-14/3*t)*(-GAMMA(k+1,7/3*t)+GAMMA(k+1,4/3*t)),t = 0 .. infinity)

I currently miss familiarity with the Gamma distribution to simplify more.

I get c(k) = 95/729, omitting x:=0.

Where I would use UpperCase Int in the definitions and finally would call 'value'
to get the result.

Sum(c(k), k = 0 .. 20); value(%);

  198750905161/73222472421

I get c(k) = 95/729, omitting x:=0.

Where I would use UpperCase Int in the definitions and finally would call 'value'
to get the result.

Sum(c(k), k = 0 .. 20); value(%);

  198750905161/73222472421

If you type in your answers then you see a green arrow in the menu bar.

It is for uploading (avoid fancy characters in file names). Just try, it works.