@Markiyan Hirnyk No offense, but the timing difference is largely irrelevant.  That is, if this operation were going to be executed many times, it would most likely be inside a procedure---particularly if speed were a concern.  But the statement (seq(p, x=sol1)) doesn't work inside a procedure, not unless p is declared as global, which is unlikely and poor practice. I'm biased because of the way I use Maple. The method you propose is handy at the top-level because it's terse, but not useful for Maple coding.

@Markiyan Hirnyk An alternative to Acer's use of eval is to request that the solve procedure return equations rather than values.  Here that is done by enclosing the unknown, x, in a set.  Thus

proc()
local p,sols,s,x;
    p:=12+10*x^4+7*x^3-89*x^2+65*x:
    sols := [fsolve(p,{x})]:
    eval~(p,sols);
end proc();
                          -6         -8        -7        -6
                   [0.5 10  , -0.7 10  , 0.3 10  , 0.1 10  ]

The anonymous procedure was used to show it works with local variables.

Equations can be easily constructed with

 x =~ [ 1, 2, 3];
                                    [x = 1, x = 2, x = 3]

@Markiyan Hirnyk An alternative to Acer's use of eval is to request that the solve procedure return equations rather than values.  Here that is done by enclosing the unknown, x, in a set.  Thus

proc()
local p,sols,s,x;
    p:=12+10*x^4+7*x^3-89*x^2+65*x:
    sols := [fsolve(p,{x})]:
    eval~(p,sols);
end proc();
                          -6         -8        -7        -6
                   [0.5 10  , -0.7 10  , 0.3 10  , 0.1 10  ]

The anonymous procedure was used to show it works with local variables.

Equations can be easily constructed with

 x =~ [ 1, 2, 3];
                                    [x = 1, x = 2, x = 3]

The problem is not with the integration.  It is with your function U. You are doing (essentially)

U := proc(d)
   if d = 0 then (* do something *)
   else  ... U(d-1) ....
   fi
end proc:

That cannot work if d is not a nonnegative integer because the recursion will never terminate. That is, if you called U with 1/2, then you will get recursive calls U(-1/2), U(-3/2), ... If called with a symbolic variable, k, you get U(k), U(k-1), U(k-2), ... 

The problem is not with the integration.  It is with your function U. You are doing (essentially)

U := proc(d)
   if d = 0 then (* do something *)
   else  ... U(d-1) ....
   fi
end proc:

That cannot work if d is not a nonnegative integer because the recursion will never terminate. That is, if you called U with 1/2, then you will get recursive calls U(-1/2), U(-3/2), ... If called with a symbolic variable, k, you get U(k), U(k-1), U(k-2), ... 

@Markiyan Hirnyk Yes, that works, but note that the assignment to y has no effect.  You can get your original to work in a procedure by declaring y as global. Thus

 proc() global y; local x; y := x^2; seq(y, x=[1,2,3]) end ();
                                           1, 4, 9

That works because globals get fully evaluated, even in seq, etc. Normally you don't want stuff to fully evaluate in a procedure because it can be expensive. Also, the assignment to y is a side-effect, since y is global.

@Markiyan Hirnyk Yes, that works, but note that the assignment to y has no effect.  You can get your original to work in a procedure by declaring y as global. Thus

 proc() global y; local x; y := x^2; seq(y, x=[1,2,3]) end ();
                                           1, 4, 9

That works because globals get fully evaluated, even in seq, etc. Normally you don't want stuff to fully evaluate in a procedure because it can be expensive. Also, the assignment to y is a side-effect, since y is global.

A problem with this approach is that it typically doesn't work inside a procedure. More specifically, it doesn't work if y is a local variable. For example

proc() local y; y := x^2; seq(y, x=[1,2,3]) end proc();
                                  2   2   2
                                  x , x , x

y := x^2: seq(y, x=[1,2,3]);
                                    1, 4, 9

The reason it fails in the procedure is that ?seq has special evalution rules, so does not evaluate its first argument. To work around that, delay the evaluation of seq

proc() local y; y := x^2; eval('seq'(y, x=[1,2,3])); end proc();  
                                    1, 4, 9

A problem with this approach is that it typically doesn't work inside a procedure. More specifically, it doesn't work if y is a local variable. For example

proc() local y; y := x^2; seq(y, x=[1,2,3]) end proc();
                                  2   2   2
                                  x , x , x

y := x^2: seq(y, x=[1,2,3]);
                                    1, 4, 9

The reason it fails in the procedure is that ?seq has special evalution rules, so does not evaluate its first argument. To work around that, delay the evaluation of seq

proc() local y; y := x^2; eval('seq'(y, x=[1,2,3])); end proc();  
                                    1, 4, 9

Alas, that still isn't sufficient because there is no relation defining J(n).

@davidp That doesn't make sense.  If J(n) = 1/2 (i.e. is a constant function equal to 1/2 for all n), then  we have S(n+1) = 1/4 (i.e. is a constant function for all n). But S(0) is given as 0.45. 

The description of the recursive relations is incomplete. Missing is an equation for A(n).

@pagan Here's another method, though a bit slower than I expected

(nops(L)-add(sign(x),x=L))/2

I'm surprised sign, a builtin, is significantly slower than `if`. With a slight modification we get the speed back:

(nops(L)-add(CopySign(1,x),x=L)/2

@pagan Here's another method, though a bit slower than I expected

(nops(L)-add(sign(x),x=L))/2

I'm surprised sign, a builtin, is significantly slower than `if`. With a slight modification we get the speed back:

(nops(L)-add(CopySign(1,x),x=L)/2

That works fine for me.  Possibly you have assigned one of the variables in the input?

While it shouldn't matter in this particular usage, note that gamma is a symbolic constant in Maple, akin to Pi.  Thus

evalf(gamma);
                       0.5772156649