Scilab, cald function

function [d] = cald(T, y)
    exec('**/rich.sci', -1)

    tol = 1e-12

    n = size(y, 1)
    a = zeros(n, 1)
    b = zeros(n, 1)
    c = zeros(n, 1)
    d = zeros(n, 1)
    u = zeros(n, 1)

    h = zeros(n, 1)
    // vecteur l pour contenir la différence y(i+1)-y(i)
    l = zeros(n, 1)

    h = T(2:n) - T(1:n-1)
    l(1:n-1) = y(2:n) - y(1:n-1) // l(n) = 0 car de toute façon c(n) = 0

    for i = 1:n-1
        if (abs(h(i)) < tol)then
            error("h(i) trop petit")
        end
        // a(1) = 0
        a(i+1) = 1 / h(i)
        // c(n) = 0
        c(i) = 1 / h(i)
    end

    b = 2 * a(1:n) + 2 * c(1:n)
    u(1) = l(1) .* c(1) .* c(1)
    u(2:n) = l(1:n-1) .* (a(2:n).*a(2:n)) + l(2:n) .* (c(2:n).*c(2:n))

    d = rich(a, b, c, u)

endfunction