Sum of a series

Compute the nth partial sum of a series.

For this task, use S(x) = 1/x², from 1 to 1000. (This approximates the Riemann zeta function. The Basel problem solved this: ζ(2) = p²/6.)

ActionScript

<lang ActionScript>function partialSum(n:uint):Number { var sum:Number = 0; for(var i:uint = 1; i <= n; i++) sum += 1/(i*i); return sum; } trace(partialSum(1000));</lang>

ALGOL 68

Works with: ALGOL 68 version Standard - no extensions to language used

Works with: ALGOL 68G version Any - tested with release mk15-0.8b.fc9.i386

Works with: ELLA ALGOL 68 version Any (with appropriate job cards) - tested with release 1.8.8d.fc9.i386

<lang algol68>MODE RANGE = STRUCT(INT lwb, upb);

PROC sum = (PROC (INT)LONG REAL f, RANGE range)LONG REAL:(

 LONG REAL sum := LENG 0.0;
 FOR i FROM lwb OF range TO upb OF range DO
    sum := sum + f(i)
 OD;
 sum

);

RANGE range = (1,100); PROC f = (INT x)LONG REAL: LENG REAL(1) / LENG REAL(x)**2;

print(("Sum of f(x) from", lwb OF range, " to ",upb OF range," is ", SHORTEN sum(f,range),".", new line))</lang> Output: <lang algol68>Sum of f(x) from +1 to +100 is +1.63498390018489e +0.</lang>

AutoHotkey allows the precision of floating point numbers generated by math operations to be adjusted via the SetFormat command. The default is 6 decimal places. <lang autohotkey>SetFormat, FloatFast, 0.15 While A_Index <= 1000

sum += 1/A_Index**2

MsgBox,% sum ;1.643934566681554</lang>

AWK

<lang awk>$ awk 'BEGIN{for(i=1;i<=1000;i++)s+=1/(i*i);print s}' 1.64393</lang>

BASIC

Works with: QuickBasic version 4.5

<lang qbasic>function s(x%)

  s = 1 / x ^ 2

end function

function sum(low%, high%)

  ret = 0
  for i = low to high
     ret = ret + s(i)
  next i
  sum = ret

end function print sum(1, 1000)</lang>

C

<lang c>#include <stdio.h>

double Invsqr(double n) { return 1 / (n*n); }

int main (int argc, char *argv[]) { int i, start = 1, end = 1000; double sum = 0.0;

for( i = start; i <= end; i++) sum += Invsqr((double)i);

printf("%16.14f\n", sum);

return 0; }</lang>

C++

<lang cpp>#include <iostream>

double f(double x);

int main() {

   unsigned int start = 1;
   unsigned int end = 1000;
   double sum = 0;

   for( unsigned int x = start; x <= end; ++x )
   {
       sum += f(x);
   }
   std::cout << "Sum of f(x) from " << start << " to " << end << " is " << sum << std::endl;
   return 0;

}

double f(double x) {

   return ( 1.0 / ( x * x ) );

}</lang>

C#

<lang csharp>class Program {

   static void Main(string[] args)
   {
       // Create and fill a list of number 1 to 1000

       List<double> myList = new List<double>();
       for (double i = 1; i < 1001; i++)
       {
           myList.Add(i);
       }
       // Calculate the sum of 1/x^2

       var sum = myList.Sum(x => 1/(x*x));

       Console.WriteLine(sum);
       Console.ReadLine();
   }

}</lang>

An alternative approach using Enumerable.Range() to generate the numbers.

<lang csharp>class Program {

   static void Main(string[] args)
   {
       double sum = Enumerable.Range(1, 1000).Sum(x => 1.0 / (x * x));

       Console.WriteLine(sum);
       Console.ReadLine();
   }

}</lang>

Clojure

<lang lisp>(reduce + (map #(/ 1.0 % %) (range 1 1001)))</lang>

Common Lisp

<lang lisp>(loop for x from 1 to 1000 summing (/ (expt x 2)))</lang>

D

<lang d>import std.stdio, std.traits;

ReturnType!(TF) series(TF)(TF func, int end, int start=1) {

   ReturnType!(TF) sum = 0;
   for (int i = start; i <= end; i++)
       sum += func(i);
   return sum;

}

void main() {

   writefln("Sum: ", series((int n){return 1.0L / (n*n);}, 1_000));

}</lang>

Delphi

unit Form_SumOfASeries_Unit;

interface

uses

 Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
 Dialogs, StdCtrls;

type

 TFormSumOfASeries = class(TForm)
   M_Log: TMemo;
   B_Calc: TButton;
   procedure B_CalcClick(Sender: TObject);
 private
   { Private-Deklarationen }
 public
   { Public-Deklarationen }
 end;

var

 FormSumOfASeries: TFormSumOfASeries;

implementation

{$R *.dfm}

function Sum_Of_A_Series(_from,_to:int64):extended; begin

 result:=0;
 while _from<=_to do
 begin
   result:=result+1.0/(_from*_from);
   inc(_from);
 end;

end;

procedure TFormSumOfASeries.B_CalcClick(Sender: TObject); begin

 try
   M_Log.Lines.Add(FloatToStr(Sum_Of_A_Series(1,trunc(sqrt($7FFFFFFFFFFFFFFF-1)))));
 except
   M_Log.Lines.Add('Error');
 end;

end;

end.

</lang>

Pi^2/6 = 1.644934066_848226436472415166646

result = 1,644934066_51266

E

<lang e>pragma.enable("accumulator") accum 0 for x in 1..1000 { _ + 1 / x ** 2 }</lang>

Erlang

<lang erlang>lists:sum([1/math:pow(X,2) || X <- lists:seq(1,1000)]).</lang>

Factor

<lang factor>1000 [1,b] [ >float sq recip ] map-sum</lang>

Forth

<lang forth>: sum ( fn start count -- fsum )

 0e
 bounds do
   i s>d d>f dup execute f+
 loop drop ;

noname ( x -- 1/x^2 ) fdup f* 1/f ; ( xt )

1 1000 sum f. \ 1.64393456668156 pi pi f* 6e f/ f. \ 1.64493406684823</lang>

Fortran

In ISO Fortran 90 and later, use SUM intrinsic: <lang fortran>real, dimension(1000) :: a = (/ (1.0/(i*i), i=1, 1000) /) real :: result

result = sum(a);</lang>

F#

The following function will do the task specified. <lang fsharp>let rec f (x : float) =

   match x with
       | 0. -> x
       | x -> (1. / (x * x)) + f (x - 1.)</lang>

In the interactive F# console, using the above gives: <lang fsharp>> f 1000. ;; val it : float = 1.643934567</lang> However this recursive function will run out of stack space eventually (try 100000). A tail-recursive implementation will not consume stack space and can therefore handle much larger ranges. Here is such a version: <lang fsharp>#light let sum_series (max : float) =

   let rec f (a:float, x : float) = 
       match x with
           | 0. -> a
           | x -> f ((1. / (x * x) + a), x - 1.)
   f (0., max)

[<EntryPoint>] let main args =

   let (b, max) = System.Double.TryParse(args.[0])
   printfn "%A" (sum_series max)
   0</lang>

This block can be compiled using fsc --target exe filename.fs or used interactively without the main function.

Groovy

Start with smallest terms first to minimize rounding error: <lang groovy>println ((1000..1).collect { x -> 1/(x*x) }.sum())</lang>

Output:

1.6439345654

Haskell

With a list comprehension: <lang haskell>sum [1 / x ^ 2 | x <- [1..1000]]</lang> With higher-order functions: <lang haskell>sum $ map (\x -> 1 / x ^ 2) [1..1000]</lang> In point-free style: <lang haskell>(sum . map (1/) . map (^2)) [1..1000]</lang>

Icon

<lang icon>procedure main()

  local i, sum
  sum := 0 & i := 0
  every sum +:= 1.0/((| i +:= 1 ) ^ 2) \1000
  write(sum)

end</lang>

IDL

<lang idl>print,total( 1/(1+findgen(1000))^2)</lang>

J

<lang j> NB. sum of reciprocals of squares of first thousand positive integers

  +/ % *: >: i. 1000

1.64393

  (*:o.1)%6       NB. pi squared over six, for comparison

1.64493

  1r6p2           NB.  As a constant (J has a rich constant notation)

1.64493</lang>

Java

<lang java>public class Sum{

   public static double f(double x){
      return 1/(x*x);
   }

   public static void main(String[] args){
      double start = 1;
      double end = 1000;
      double sum = 0;

      for(double x = start;x <= end;x++) sum += f(x);

      System.out.println("Sum of f(x) from " + start + " to " + end +" is " + sum);
   }

}</lang>

JavaScript

<lang javascript>function sum(a,b,fn) {

  var s = 0;
  for ( ; a <= b; a++) s += fn(a);
  return s;

}

sum(1,1000, function(x) { return 1/(x*x) } )  // 1.64393456668156</lang>

Logo

<lang logo>to series :fn :a :b

 localmake "sigma 0
 for [i :a :b] [make "sigma :sigma + invoke :fn :i]
 output :sigma

end to zeta.2 :x

 output 1 / (:x * :x)

end print series "zeta.2 1 1000 make "pi (radarctan 0 1) * 2 print :pi * :pi / 6</lang>

Lua

<lang lua> sum = 0 for i = 1, 1000 do sum = sum + 1/i^2 end print(sum) </lang>

Lucid

<lang lucid>series = ssum asa n >= 1000

  where
        num = 1 fby num + 1;
        ssum = ssum + 1/(num * num)
  end;</lang>

Mathematica

This is the straightforward solution of the task: <lang mathematica>Sum[1/x^2, {x, 1, 1000}]</lang> However this returns a quotient of two huge integers (namely the exact sum); to get a floating point approximation, use N: <lang mathematica>N[Sum[1/x^2, {x, 1, 1000}]]</lang> or better: <lang mathematica>NSum[1/x^2, {x, 1, 1000}]</lang> Which gives a higher or equal accuracy/precision. Alternatively, get Mathematica to do the whole calculation in floating point by using a floating point value in the formula: <lang mathematica>Sum[1./x^2, {x, 1, 1000}]</lang> Other ways include (exact, approximate,exact,approximate): <lang mathematica>Total[Table[1/x^2, {x, 1, 1000}]] Total[Table[1./x^2, {x, 1, 1000}]] Plus@@Table[1/x^2, {x, 1, 1000}] Plus@@Table[1./x^2, {x, 1, 1000}]</lang>

Maxima

<lang maxima>(%i45) sum(1/x^2, x, 1, 1000);

      835459384831496894781878542648[806 digits]396236858699094240207812766449

(%o45) ------------------------------------------------------------------------

      508207201043258126178352922730[806 digits]886537101453118476390400000000

(%i46) sum(1/x^2, x, 1, 1000),numer; (%o46) 1.643934566681561</lang>

MAXScript

<lang maxscript>total = 0 for i in 1 to 1000 do (

   total += 1.0 / pow i 2

) print total</lang>

MMIX

<lang mmix>x IS $1 % flt calculations y IS $2 % id z IS $3 % z = sum series t IS $4 % temp var

LOC Data_Segment GREG @ BUF OCTA 0,0,0 % print buffer

LOC #1000 GREG @

// print floating point number in scientific format: 0.xxx...ey.. // most of this routine is adopted from: // http://www.pspu.ru/personal/eremin/emmi/rom_subs/printreal.html // float number in z GREG @ NaN BYTE "NaN..",0 NewLn BYTE #a,0 1H LDA x,NaN TRAP 0,Fputs,StdOut GO $127,$127,0

prtFlt FUN x,z,z % test if z == NaN BNZ x,1B CMP $73,z,0 % if necessary remember it's neg BNN $73,4F Sign BYTE '-' LDA $255,Sign TRAP 0,Fputs,StdOut ANDNH z,#8000 % make number pos // normalizing float number 4H SETH $74,#4024 % initialize mulfactor = 10.0 SETH $73,#0023 INCMH $73,#86f2 INCML $73,#6fc1 % FLOT $73,$73 % $73 = float 10^16 SET $75,16 % set # decimals to 16 8H FCMP $72,z,$73 % while z >= 10^16 do BN $72,9F % FDIV z,z,$74 % z = z / 10.0 ADD $75,$75,1 % incr exponent JMP 8B % wend 9H FDIV $73,$73,$74 % 10^16 / 10.0 5H FCMP $72,z,$73 % while z < 10^15 do BNN $72,6F FMUL z,z,$74 % z = z * 10.0 SUB $75,$75,1 % exp = exp - 1 JMP 5B NulPnt BYTE '0','.',#00 6H LDA $255,NulPnt % print '0.' to StdOut TRAP 0,Fputs,StdOut FIX z,0,z % convert float z to integer // print mantissa 0H GREG #3030303030303030 STO 0B,BUF STO 0B,BUF+8 % store print mask in buffer LDA $255,BUF+16 % points after LSD % repeat 2H SUB $255,$255,1 % move pointer down DIV z,z,10 % (q,r) = divmod z 10 GET t,rR % get remainder INCL t,'0' % convert to ascii digit STBU t,$255,0 % store digit in buffer BNZ z,2B % until q == 0 TRAP 0,Fputs,StdOut % print mantissa Exp BYTE 'e',#00 LDA $255,Exp % print 'exponent' indicator TRAP 0,Fputs,StdOut // print exponent 0H GREG #3030300000000000 STO 0B,BUF LDA $255,BUF+2 % store print mask in buffer CMP $73,$75,0 % if exp neg then place - in buffer BNN $73,2F ExpSign BYTE '-' LDA $255,ExpSign TRAP 0,Fputs,StdOut NEG $75,$75 % make exp positive 2H LDA $255,BUF+3 % points after LSD % repeat 3H SUB $255,$255,1 % move pointer down DIV $75,$75,10 % (q,r) = divmod exp 10 GET t,rR INCL t,'0' STBU t,$255,0 % store exp. digit in buffer BNZ $75,3B % until q == 0 TRAP 0,Fputs,StdOut % print exponent LDA $255,NewLn TRAP 0,Fputs,StdOut % do a NL GO $127,$127,0 % return

i IS $5 ;iu IS $6 Main SET iu,1000 SETH y,#3ff0 y = 1.0 SETH z,#0000 z = 0.0 SET i,1 for (i=1;i<=1000; i++ ) { 1H FLOT x,i x = int i FMUL x,x,x x = x^2 FDIV x,y,x x = 1 / x FADD z,z,x s = s + x ADD i,i,1 CMP t,i,iu PBNP t,1B } z = sum GO $127,prtFlt print sum --> StdOut TRAP 0,Halt,0</lang> Output:

~/MIX/MMIX/Rosetta> mmix sumseries
0.1643934566681562e1

Nial

<lang nial>|sum (1 / power (count 1000) 2) =1.64393</lang>

OCaml

<lang ocaml>let sum a b fn =

 let result = ref 0. in
 for i = a to b do
   result := !result +. fn i
 done;
 !result</lang>

# sum 1 1000 (fun x -> 1. /. (float x ** 2.))
- : float = 1.64393456668156124

or in a functional programming style: <lang ocaml>let sum a b fn =

 let rec aux i r =
   if i > b then r
   else aux (succ i) (r +. fn i)
 in
 aux a 0.

</lang>

Octave

Given a vector, the sum of all its elements is simply sum(vector); a range can be generated through the range notation: sum(1:1000) computes the sum of all numbers from 1 to 1000. To compute the requested series, we can simply write:

OpenEdge/Progress

Conventionally like elsewhere:

<lang openedge>def var dcResult as decimal no-undo. def var n as int no-undo.

do n = 1 to 1000 :

 dcResult = dcResult + 1 / (n * n)  .

end.

display dcResult .</lang>

or like this:

<lang openedge>def var n as int no-undo.

repeat n = 1 to 1000 :

 accumulate 1 / (n * n) (total).

end.

display ( accum total 1 / (n * n) ) .</lang>

Oz

With higher-order functions: <lang oz>declare

 fun {SumSeries S N}
    {FoldL {Map {List.number 1 N 1} S}
     Number.'+' 0.}
 end

 fun {S X}
    1. / {Int.toFloat X*X}
 end

in

 {Show {SumSeries S 1000}}</lang>

Iterative: <lang oz> fun {SumSeries S N}

    R = {NewCell 0.}
 in
    for I in 1..N do
       R := @R + {S I}
    end
    @R
 end</lang>

Perl

<lang perl>my $sum = 0; $sum += 1 / $_ ** 2 foreach 1..1000; print "$sum\n";</lang> or <lang perl>use List::Util qw(reduce); $sum = reduce { $a + 1 / $b ** 2 } 0, 1..1000; print "$sum\n";</lang>

Perl 6

Works with: Rakudo version #22 "Thousand Oaks"

In general, the $nth partial sum of a series whose terms are given by a unary function &f is

So what's needed in this case is

PL/I

<lang PL/I> /* sum the first 1000 terms of the series 1/n**2. */ s = 0;

do i = 1000 to 1 by -1;

  s = s + i/1000q0;

end;

put skip list (s); </lang>

Pop11

<lang pop11>lvars s = 0, j; for j from 1 to 1000 do

   s + 1.0/(j*j) -> s;

endfor;

s =></lang>

PowerShell

<lang powershell>$x = 1..1000 `

      | ForEach-Object { 1 / ($_ * $_) } `
      | Measure-Object -Sum

Write-Host Sum = $x.Sum</lang>

PureBasic

<lang PureBasic>Define i, sum.d

For i=1 To 1000

 sum+1.0/(i*i)

Next i

Debug sum</lang> Answer = 1.6439345666815615

Python

<lang python>print ( sum(1.0 / x ** 2 for x in range(1, 1001)) )</lang>

R

<lang r>print( sum( 1/seq(1000)^2 ) )</lang>

Ruby

<lang ruby>puts (1..1000).inject(0) {|sum, x| sum + 1.0 / x ** 2}</lang>

Scala

<lang scala>scala> 1 to 1000 map (x => 1.0 / (x * x)) sum res30: Double = 1.6439345666815615</lang>

Scheme

<lang scheme>(define (sum a b fn)

 (do ((i a (+ i 1))
      (result 0 (+ result (fn i))))
     ((> i b) result)))

(sum 1 1000 (lambda (x) (/ 1 (* x x)))) ; fraction (exact->inexact (sum 1 1000 (lambda (x) (/ 1 (* x x))))) ; decimal</lang>

Slate

Manually coerce it to a float, otherwise you will get an exact (and slow) answer:

<lang slate>((1 to: 1000) reduce: [|:x :y | x + (y squared reciprocal as: Float)]).</lang>

Smalltalk

<lang smalltalk>( (1 to: 1000) fold: [:sum :aNumber |

 sum + (aNumber squared reciprocal) ] ) asFloat displayNl.</lang>

Tcl

Works with: Tcl version 8.5

<lang tcl>package require Tcl 8.5

proc partial_sum {func - start - stop} {

   for {set x $start; set sum 0} {$x <= $stop} {incr x} {
       set sum [expr {$sum + [apply $func $x]}]
   }
   return $sum

}

set S {x {expr {1.0 / $x**2}}}

partial_sum $S from 1 to 1000 ;# => 1.6439345666815615</lang>

Or, this uses package struct::list from

Library: tcllib

<lang tcl>package require Tcl 8.5 package require struct::list

proc sum_of {lambda nums} {

   struct::list fold [struct::list map $nums [list apply $lambda]] 0 ::tcl::mathop::+

}

sum_of $S [range 1 1001] ;# ==> 1.6439345666815615</lang>

The helper range procedure is: <lang tcl># a range command akin to Python's proc range args {

   foreach {start stop step} [switch -exact -- [llength $args] {
       1 {concat 0 $args 1}
       2 {concat   $args 1}
       3 {concat   $args  }
       default {error {wrong # of args: should be "range ?start? stop ?step?"}}
   }] break
   if {$step == 0} {error "cannot create a range when step == 0"}
   set range [list]
   while {$step > 0 ? $start < $stop : $stop < $start} {
       lappend range $start
       incr start $step
   }
   return $range

}</lang>

TI-89 BASIC

UnixPipes

<lang bash>term() {

  b=$1;res=$2
  echo "scale=5;1/($res*$res)+$b" | bc

}

sum() {

 (read B; res=$1;
 test -n "$B" && (term $B $res) || (term 0 $res))

}

fold() {

 func=$1
 (while read a ; do
     fold $func | $func $a
 done)

}

(echo 3; echo 1; echo 4) | fold sum</lang>

Ursala

The expression plus:-0. represents a function returning the sum of any given list of floating point numbers, or zero if it's empty, using the built in reduction operator, :-, and the binary addition function, plus. The rest the expression constructs the series by inverting the square of each number in the list from 1 to 1000. <lang Ursala>#import flo

import nat

cast %e

total = plus:-0 div/*1. sqr* float*t iota 1001</lang> output:

1.643935e+00