If you don't understand what's going on here, have a look at
  <http://para.inria.fr/~doligez/primrec/>

We're working with natural numbers, so any number different
from 0 must be strictly greater than 0.

Each line has an explanation on the right.
<<{foo}>> means <<by the definition of foo>>
<<by "foo">> means <<by lemma "foo">> or <<by theorem "foo">>
<<by recursive hypothesis>> means just that

Theorem "zero"
  zero() = 0

Proof
  zero() = Z				{zero}
         = 0				{Z}
QED

Theorem "one"
  one() = 1

Proof
  one() = S*(Z)()			{one}
        = S(Z())			{*}
        = S(0)				{Z}
        = 1				{S}
QED

Theorem "two"
  two() = 2

Proof
  two() = S*(S*(Z))()			{two}
        = S(S*(Z)())			{*}
        = S(S(Z()))			{*}
        = S(S(0))			{Z}
        = S(1)				{S}
        = 2				{S}
QED

Theorem "add"
  add(x,y) = x+y

Proof
  Recursion on y:
  Base case:
    add(x,0) = f(I,S[3]3)(x,0)		{add}
             = I(x)			{f}
             = x			{I}
             = x + 0
  Recursive case:
    add(x,y+1) = f(I,S[3]3)(x,y+1)	{add}
               = S[3]3(x,y,add(x,y))	{f}
               = S(add(x,y))		{[]}
               = add(x,y)+1		{S}
               = x+y+1			by recursive hypothesis
QED

Theorem "mul"
  mul(x,y) = xy

Proof
  Recursion on y:
  Base case:
    mul(x,0) = f(C,add[1,3]3)(x,0)	{mul}
             = C(x)			{f}
             = 0			{C}
  Recursive case:
    mul(x,y+1) = f(C,add[1,3]3)(x,y+1)		{mul}
               = add[1,3]3(x,y,mul(x,y))	{f}
               = add(x,mul(x,y))		{[]}
               = x+mul(x,y)			by "add"
               = x(y+1)				by recursive hypothesis
QED

Theorem "pred0"
  pred(0) = 0

Proof
  pred(0) = f(Z,I[1]2)(0)		{pred}
          = Z()				{f}
          = 0				{Z}
QED

Theorem "pred1"
  pred(x+1) = x

Proof
  pred(x+1) = f(Z,I[1]2)(x+1)		{pred}
            = I[1]2(x,pred(x))		{f}
            = x				{[]}
QED

Theorem "not"
  if x = 0, then not(x) = 1
  if x != 0, then not(x) = 0

Proof
  not(0) = f(S*(Z),Z[]2)(0)		{not}
         = S*(Z)()			{f}
         = S(Z())			{*}
         = S(0)				{Z}
         = 1				{S}
  not(x+1) = f(S*(Z),Z[]2)(0)		{not}
           = Z[]2(x,not(x))		{f}
           = Z()			{[]}
           = 0				{Z}
QED

Theorem "neither"
  neither(x,y) != 0 if and only if x = y = 0

Proof
  neither(0,0) = f(not,Z[]3)(0,0)		{neither}
               = not(0)				{f}
               = 1				by "not"
  neither(x+1,0) = f(not,Z[]3)(x+1,0)		{neither}
                 = not(x+1)			{f}
                 = 0				by "not"
  neither(x,y+1) = f(not,Z[]3)(x,y+1)		{neither}
                 = Z[]3(x,y,neither(x,y))	{f}
                 = Z()				{[]}
                 = 0				{Z}
QED

Lemma "sub0"
  sub(x,0) = x-0

Proof
  sub(x,0) = f(I,pred[3]3)(x,0)		{sub}
           = I(x)			{f}
           = x				{I}
           = x-0
QED

Lemma "sub1"
  sub(x,y+1) = pred(sub(x,y))

Proof
  sub(x,y+1) = f(I,pred[3]3)(x,y+1)	{sub}
             = pred[3]3(x,y,sub(x,y))	{f}
             = pred(sub(x,y))		{[]}
QED

Theorem "sub"
  if x >= y, then sub(x,y) = x-y
  if x <= y, then sub(x,y) = 0

Proof
  Assume x >= y
    Recursion on y:
    Base case: sub(x,0) = x-0		by "sub0"
    Recursive case: y>0, let y' = y-1, so x > y'
      sub(x,y) = sub(x,y'+1)		{y'}
               = pred(sub(x,y'))	by "sub1"
               = pred(x-y')		by recursion hypothesis
               = x-y'-1			by "pred1"
               = x-y
  Assume x <= y
    Recursion on y-x:
    Base case: sub(x,x) = 0		by "sub0"
    Recursive case: y>x, let y' = y-1, so y' >= x, and y'-x = y-x-1
      sub(x,y) = sub(x,y'+1)		{y'}
               = pred(sub(x,y'))	by "sub1"
               = pred(0)		by recursion hypothesis
               = 0			by "pred0"
QED

Theorem "greater"
  greater(x,y) != 0 if and only if x > y

Proof
  greater(x,y) = sub(x,y)		{greater}
  Assume x > y
    greater(x,y) = x-y			by "sub"
    greater(x,y) != 0
  Assume x <= y
    greater(x,y) = 0			by "sub"
QED

Theorem "equal"
  if x = y, then equal(x,y) = 1
  if x != y, then equal(x,y) = 0

Proof
  equal(x,y) = neither*(greater,greater)[1,2,2,1]2(x,y)		{equal}
             = neither*(greater,greater)(x,y,y,x)		{[]}
             = neither(greater(x,y),greater(y,x))		{*}
  Assume x = y
    equal(x,x) = neither(greater(x,x),greater(x,x))
               = neither(0,0)					by "greater"
               = 1						by "neither"
  Assume x > y
    equal(x,y) = neither(t,0)	with t != 0			by "greater"
               = 0						by "neither"
  Assume x < y
    equal(x,y) = neither(0,t)	with t != 0 			by "greater"
               = 0						by "neither"
QED

Theorem "if"
  if(0,y,z) = y
  if(x+1,y,z) = z

Proof
  if(0,y,z) = f(I[2]2,I[1]4)[2,3,1]3(0,y,z)	{if}
            = f(I[2]2,I[1]4)(y,z,0)		{[]}
            = I[2]2(y,z)			{f}
            = I(z)				{[]}
            = z					{I}
  if(x+1,y,z) = f(I[2]2,I[1]4)[2,3,1]3(x+1,y,z)	{if}
              = f(I[2]2,I[1]4)(y,z,x+1)		{[]}
              = I[1]4(y,z,x,if(x,y,z))		{f}
              = I(y)				{[]}
              = y				{I}
QED

Theorem "divides"
  if there exists t < z such that x = yt, then divides(x,y,z) = 1
  if not, then divides (x,y,z) = 0

Proof
  By recursion on z:
  Base case:
    divides(x,y,0) = f(Z[]2,if*(equal*(I,mul),one,I))(x,y,0)	{divides}
                   = Z[]2(x,y)					{f}
                   = Z()					{[]}
                   = 0						{Z}
  Recursive case:
    divides(x,y,z+1)
        = f(Z[]2,if*(equal*(I,mul),one,I))(x,y,z+1)		{divides}
        = if*(equal*(I,mul),one,I)(x,y,z,divides(x,y,z))	{f}
        = if(equal*(I,mul)(x,y,z),one,I(divides(x,y,z)))	{*}
        = if(equal(I(x),mul(y,z)),one,I(divides(x,y,z)))	{*}
        = if(equal(x,mul(y,z)),one,divides(x,y,z))		{I}
        = if(equal(x,yz),one,divides(x,y,z))			by "mul"
        = if(equal(x,yz),1,divides(x,y,z))			by "one"
    Assume x = yz
      divides(x,y,z+1)
          = if(1,1,divides(x,y,z))				by "equal"
          = 1							by "if"
    Assume x != yz
      divides(x,y,z+1)
          = if(0,1,divides(x,y,z))				by "equal"
          = divides(x,y,z)					by "if"
      Assume there exists t < z+1 such that x = yt
        t != z (because we assumed x != yz),
        so there exists t < z such that x = yt,
        so divides(x,y,z) = 1 (by recursion hypothesis)
        thus divides(x,y,z+1) = 1
      Assume there is no t < z+1 such that x = yt
        there is no t < z such that x = yt
        so divides(x,y,z) = 0 (by recursion hypothesis)
        thus divides(x,y,z+1) = 0
QED

Theorem "multiple"
  if x > 0 and y > 1 and x is a multiple of y, then multiple(x,y) = 1
  if x > 0 and (y <= 1 or x is not a multiple of y), then multiple(x,y) = 0

Proof
  multiple(x,y) = divides[1,2,1]2(x,y)		{multiple}
                = divides(x,y,x)		{[]}
  Assume x > 0 and y > 1
    Assume x is a multiple of y
      there exists t < x such that x = yt
      multiple(x,y) = 1				by "divides"
    Assume x is not a multiple of y
      there is no t < x such that x = yt
      multiple(x,y) = 0				by "divides"
  Assume x > 0 and y = 1
    there is no t < x such that x = t
    multiple(x,y) = 0				by "divides"
  Assume x > 0 and y = 0
    there is no t < x such that x = 0
    multiple(x,y) = 0				by "divides"
QED

Theorem "primaux"
  if x > 1 and there is a divisor of x in [2,y+1], then primaux(x,y) = 0
  if x > 1 and there is no divisor of x in [2,y+1], then primaux(x,y) = 1

Proof
  Assume x > 1
    By recursion on y:
    Base case:
      primaux(x,0)
          = f(S*(C),if*(multiple*(I,S*(S)),zero,I))(x,0)	{primaux}
          = S*(C)(x)						{f}
          = S(C(x))						{*}
          = S(0)						{C}
          = 1							{S}
    Recursive case:
      primaux(x,y+1)
          = f(S*(C),if*(multiple*(I,S*(S)),zero,I))(x,y+1)	{primaux}
          = if*(multiple*(I,S*(S)),zero,I)(x,y,primaux(x,y))	{f}
          = if(multiple*(I,S*(S))(x,y),zero,I(primaux(x,y)))	{*}
          = if(multiple(I(x),S*(S)(y)),zero,I(primaux(x,y)))	{*}
          = if(multiple(x,S*(S)(y)),zero,I(primaux(x,y)))	{I}
          = if(multiple(x,S(S(y))),zero,I(primaux(x,y)))	{*}
          = if(multiple(x,S(y+1)),zero,I(primaux(x,y)))		{S}
          = if(multiple(x,y+2),zero,I(primaux(x,y)))		{S}
          = if(multiple(x,y+2),0,I(primaux(x,y)))		by "zero"
          = if(multiple(x,y+2),0,primaux(x,y))			{I}
      Assume x is a multiple of y+2
        there exists n in [2,y+2] that divides x (n is y+2)
        primaux(x,y+1) = if(1,0,primaux(x,y))			by "multiple"
                       = 0					by "if"
      Assume x is not a multiple of y+2
        primaux(x,y+1) = if(0,0,primaux(x,y))			by "multiple"
                       = primaux(x,y)				by "if"
        Assume there exists n in [2,y+1] that divides x
          there exists n in [2,y+2] that divides x
          primaux(x,y+1) = 0			by recursion hypothesis
        Assume there is no n in [2,y+1] that divides x
          there is no n in [2,y+1] that divides x
          primaux(x,y+1) = 1			by recursion hypothesis
QED

Theorem "isprime"
  if x > 1 and x is a prime number, then isprime(x) = 1
  if x > 1 and x is not a prime number, then isprime(x) = 0

Proof
  Assume x >= 2
    isprime(x) = primaux*(I,pred*(pred))[1,1]1(x)		{isprime}
               = primaux*(I,pred*(pred))(x,x)			{[]}
               = primaux(I(x),pred*(pred)(x))			{[]}
               = primaux(x,pred*(pred)(x))			{I}
               = primaux(x,pred(pred(x)))			{*}
               = primaux(x,pred(x-1))				by "pred1"
               = primaux(x,x-2)					by "pred1"
    Assume x is not a prime number
      there is a divisor of x in [2,x-1]
      primaux(x,x-2) = 0					by "primaux"
    Assume x is a prime number
      there is no divisor of x in [2,x-1]
      primaux(x,x-2) = 1					by "primaux"
QED

Theorem "case"
  if c1 != 0, then case(e,t3,c3,t2,c2,t1,c1) = t1
  if c1 = 0 and c2 != 0 then case(e,t3,c3,t2,c2,t1,c1) = t2
  if c1 = c2 = 0 and c3 != 0 then case(e,t3,c3,t2,c2,t1,c1) = t3
  if c1 = c2 = c3 = 0 then case(e,t3,c3,t2,c2,t1,c1) = e

Proof
  Assume c1 = x+1
    case(e,t3,c3,t2,c2,t1,x+1)
        = f(f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)[1,2,3,4,5]6,I[6]8)
                                        (e,t3,c3,t2,c2,t1,x+1)	{case}
        = I[6]8(e,t3,c3,t2,c2,t1,x,case(e,t3,c3,t2,c2,t1,x))	{t}
        = I(t1)							{[]}
        = t1							{I}
  Assume c1 = 0 and c2 = x+1
    case(e,t3,c3,t2,x+1,t1,0)
        = f(f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)[1,2,3,4,5]6,I[6]8)
                                         (e,t3,c3,t2,x+1,t1,0)	{case}
        = f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)[1,2,3,4,5]6(e,t3,c3,t2,x+1,t1)
								{f}
        = f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)(e,t3,c3,t2,x+1)	{[]}
        = I[4]6(e,t3,c3,t2,x,f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)(e,t3,c3,t2,x)
								{f}
        = I(t2)							{[]}
        = t2							{I}
  Assume c1 = c2 = 0 and c3 = x+1
    case(e,t3,x+1,t2,0,t1,0)
        = f(f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)[1,2,3,4,5]6,I[6]8)
                                          (e,t3,x+1,t2,0,t1,0)	{case}
        = f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)[1,2,3,4,5]6(e,t3,x+1,t2,0,t1)
								{f}
        = f(f(I[1]2,I[2]4)[1,2,3]4,I[4]6)(e,t3,x+1,t2,0)	{[]}
        = f(I[1]2,I[2]4)[1,2,3]4(e,t3,x+1,t2)			{f}
        = f(I[1]2,I[2]4)(e,t3,x+1)				{[]}
        = I[2]4(e,t3,x,f(I[1]2,I[2]4)(e,t3,x))			{f}
        = I(t3)							{[]}
        = t3							{I}
QED

Lemma "Chebychev"
  if x >= 2, then there is a prime number in [x+1,2x-1]

Proof
  left as an exercise to the reader (very hard !)
  see "Bertrand's conjecture" in Hardy and Wright,
  Introduction to Number Theory
QED

Definition
  findprime = f(C,case*(Z,I,isprime,I,I,Z,not*(greater))[2,2,3,3,2,1]3)

Lemma "findprime"
  if x >= 2 and there is a prime number in [x+1,y-1]
  then findprime(x,y) is the least such prime
  if x >= 2 and there is no prime in [x+1,y-1], then findprime(x,y) = 0

Proof
  Assume x >= 2
    By recursion on y:
    Base case:
      findprime(x,0)
          = f(C,case*(Z,I,isprime,I,I,Z,not*(greater))[2,2,3,3,2,1]3)(x,0)
								{findprime}
          = C(x)						{f}
          = 0							{C}
    Recursive case:
      findprime(x,y+1)
          = f(C,case*(Z,I,isprime,I,I,Z,not*(greater))[2,2,3,3,2,1]3)(x,y+1)
								    {findprime}
          = case*(Z,I,isprime,I,I,Z,not*(greater))[2,2,3,3,2,1]3
                                            (x,y,findprime(x,y))    {f}
          = case*(Z,I,isprime,I,I,Z,not*(greater))
                         (y,y,findprime(x,y),findprime(x,y),y,x)    {[]}
          = case(Z(),I(y),isprime(y),I(findprime(x,y)),
                 I(findprime(x,y)),Z(),not*(greater)(y,x))	    {*}
          = case(0,y,isprime(y),findprime(x,y),findprime(x,y),
                 0,not(greater(y,x)))				    {Z,I,*}
      Assume x >= y
        [x+1,y] is empty, so there is no prime in it
        not(greater(y,x)) = 1				by "greater", "not"
        findprime(x,y+1) = 0				by "case"
      Assume x < y
        not(greater(y,x)) = 0				by "greater", "not"
        Assume there is a prime in [x+1,y-1]
          findprime(x,y) is the least prime in [x+1,y-1]
							by recursion hypothesis
          findprime(x,y+1) = findprime(x,y)		by "case"
          findprime(x,y+1) is the least prime in [x+1,y]
        Assume there is no prime in [x+1,y-1]
          Assume y is prime
            isprime(y) = 1				by "isprime"
            findprime(x,y+1) = y			by "case"
            findprime(x,y+1) is the least prime in [x+1,y]
          Assume y is not prime
            isprime(y) = 0				by "isprime"
            findprime(x,y+1) = 0
QED

Lemma "nextprimedef"
  nextprime = findprime*(I,add)[1,1,1]1

Proof
  Expand the definitions of findprime and nextprime in the above equation.
QED

Theorem "nextprime"
  if x >= 2, then nextprime(x) is the least prime number greater than x

Proof
  nextprime(x) = findprime*(I,add)[1,1,1]1(x)		by "nextprimedef"
               = findprime*(I,add)(x,x,x)		{[]}
               = findprime(I(x),add(x,x))		{*}
               = findprime(x,add(x,x))			{I}
               = findprime(x,2x)			by "add"
  there is a prime number in [x+1,2x-1]			by "Chebychev"
  nextprime(x) is the least prime number in [x+1,2x-1]	by "findprime"
  nextprime(x) is the least prime number greater than x
QED

Theorem "nthprime"
  nthprime(n) is the prime number p such that there are n prime numbers
  less than p.

Proof
  By recursion on n:
  Base case:
    nthprime(0) = f(two,nextprime[2]2)(0)		{nthprime}
                = two					{f}
                = 2					by "two"
  Recursive case:
    nthprime(n+1) = f(two,nextprime[2]2)(n+1)		{nthprime}
                  = nextprime[2]2(n,nthprime(n))	{f}
                  = nextprime(nthprime(n))		{[]}
    let p = nthprime(n)
    p is prime and there are n primes less than p	by recursion hypothesis
    nthprime(n+1) = nextprime(p)
    nthprime(n+1) is the least prime number greater than p	by "nextprime"
    nthprime(n+1) is prime and there are n+1 primes less than it
QED

Expanding the above definitions in one big function is left as an
exercise for the reader's favorite text editor.