Here is a little Maple. Note that using "sum" on the inside causes it to find an algebraic expression for the sum over $\ell$ and using "Sum" on the outside tells it to not try to sum that algebraically over $k$.

I'll request 120 digits then round to 100 for confidence but I didn't need to. Maple increases the working precision automatically, and this is confirmed by the fact that exactly the same answers are obtained if 100 digits are requested directly.

The documentation says that Levin's u-transform is used, with these references:

Fessler, T.; Ford, W.F.; and Smith, D.A. "HURRY: An acceleration algorithm for scalar sequences and series." ACM Trans. Math. Software, Vol. 9, (1983): 346-354.

Levin, D. "Development of non-linear transformations for improving convergence of sequences". Internat. J. Comput. Math, Vol. B3, (1973): 371-388.

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^2
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^2,l=0..infinity),
       k=0..infinity)):
Z22 := evalf[100](%);
Z22 := 4.7059056441271748413982167408282763463686149628226287913\
         07611234496885020388225363260883986399848630 + 0. I

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^1
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^3,l=0..infinity),
       k=0..infinity)):
Z13 := evalf[100](%);
Z13 := 1.6470471779364125793008916295858618268156925185886856043\
         46194382751557489805887318369558006800075685 + 0. I

evalf[100](2*Z22 + 4*Z13);
16.0000000000000000000000000000000000000000000000000000000000000\
   0000000000000000000000000000000000000 + 0. I

Here is a little Maple. Note that using "sum" on the inside causes it to find an algebraic expression for the sum over $\ell$ and using "Sum" on the outside tells it to not try to sum that algebraically over $k$.

I'll request 120 digits then round to 100 for confidence but I didn't need to. Maple increases the working precision automatically, and this is confirmed by the fact that exactly the same answers are obtained if 100 digits are requested directly.

The documentation says that Levin's u-transform is used, with these references:

Fessler, T.; Ford, W.F.; and Smith, D.A. "HURRY: An acceleration algorithm for scalar sequences and series." ACM Trans. Math. Software, Vol. 9, (1983): 346-354.

Levin, D. "Development of non-linear transformations for improving convergence of sequences". Internat. J. Comput. Math, Vol. B3, (1973): 371-388.

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^2
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^2,l=0..infinity),
       k=0..infinity)):
Z22 := evalf[100](%);
Z22 := 4.7059056441271748413982167408282763463686149628226287913\
         07611234496885020388225363260883986399848630 + 0. I

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^1
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^3,l=0..infinity),
       k=0..infinity)):
Z13 := evalf[100](%);
Z13 := 1.6470471779364125793008916295858618268156925185886856043\
         46194382751557489805887318369558006800075685 + 0. I

evalf[100](2*Z22 + 4*Z13);
16.0000000000000000000000000000000000000000000000000000000000000\
   0000000000000000000000000000000000000 + 0. I

Here is a little Maple. Note that using "sum" on the inside causes it to find an algebraic expression for the sum over $\ell$ and using "Sum" on the outside tells it to not try to sum that algebraically over $k$.

I'll request 120 digits then round to 100 for confidence but I didn't need to. Maple increases the working precision automatically, and this is confirmed by the fact that exactly the same values are obtained if 100 digits are requested directly.

I don't know what method Maple uses but the fact that it gives the result as a complex number might be a clue.

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^2
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^2,l=0..infinity),
       k=0..infinity)):
Z22 := evalf[100](%);
Z22 := 4.7059056441271748413982167408282763463686149628226287913\
         07611234496885020388225363260883986399848630 + 0. I

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^1
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^3,l=0..infinity),
       k=0..infinity)):
Z13 := evalf[100](%);
Z13 := 1.6470471779364125793008916295858618268156925185886856043\
         46194382751557489805887318369558006800075685 + 0. I

evalf[100](2*Z22 + 4*Z13);
16.0000000000000000000000000000000000000000000000000000000000000\
   0000000000000000000000000000000000000 + 0. I

Here is a little Maple. Note that using "sum" on the inside causes it to find an algebraic expression for the sum over $\ell$ and using "Sum" on the outside tells it to not try to sum that algebraically over $k$.

I'll request 120 digits then round to 100 for confidence but I didn't need to. Maple increases the working precision automatically, and this is confirmed by the fact that exactly the same answers are obtained if 100 digits are requested directly.

The documentation says that Levin's u-transform is used, with these references:

Fessler, T.; Ford, W.F.; and Smith, D.A. "HURRY: An acceleration algorithm for scalar sequences and series." ACM Trans. Math. Software, Vol. 9, (1983): 346-354.

Levin, D. "Development of non-linear transformations for improving convergence of sequences". Internat. J. Comput. Math, Vol. B3, (1973): 371-388.

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^2
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^2,l=0..infinity),
       k=0..infinity)):
Z22 := evalf[100](%);
Z22 := 4.7059056441271748413982167408282763463686149628226287913\
         07611234496885020388225363260883986399848630 + 0. I

evalf[120](Sum((2*k+3)/((k+2)*(k+1)/2)^1
    * sum( (2*l+3)/((k+2)*(k+1)/2 + (l+2)*(l+1)/2)^3,l=0..infinity),
       k=0..infinity)):
Z13 := evalf[100](%);
Z13 := 1.6470471779364125793008916295858618268156925185886856043\
         46194382751557489805887318369558006800075685 + 0. I

evalf[100](2*Z22 + 4*Z13);
16.0000000000000000000000000000000000000000000000000000000000000\
   0000000000000000000000000000000000000 + 0. I