46 lines
815 B
CoffeeScript
46 lines
815 B
CoffeeScript
#-----------------------------------------------------------------------------
|
|
#
|
|
# Create a Hilbert matrix
|
|
#
|
|
# Input: Dimension on stack
|
|
#
|
|
# Output: Hilbert matrix on stack
|
|
#
|
|
# Example:
|
|
#
|
|
# > hilbert(5)
|
|
# ((1,1/2,1/3,1/4),(1/2,1/3,1/4,1/5),(1/3,1/4,1/5,1/6),(1/4,1/5,1/6,1/7))
|
|
#
|
|
#-----------------------------------------------------------------------------
|
|
|
|
|
|
|
|
#define A p1
|
|
#define N p2
|
|
|
|
#define AELEM(i, j) A->u.tensor->elem[i * n + j]
|
|
|
|
hilbert = ->
|
|
i = 0
|
|
j = 0
|
|
n = 0
|
|
save()
|
|
p2 = pop()
|
|
push(p2)
|
|
n = pop_integer()
|
|
if (n < 2)
|
|
push_symbol(HILBERT)
|
|
push(p2)
|
|
list(2)
|
|
restore()
|
|
return
|
|
push_zero_matrix(n, n)
|
|
p1 = pop()
|
|
for i in [0...n]
|
|
for j in [0...n]
|
|
push_integer(i + j + 1)
|
|
inverse()
|
|
p1.tensor.elem[i * n + j] = pop()
|
|
push(p1)
|
|
restore()
|