Sweet Snippet 之 Gram-Schmidt 正交化

-- vector add
function add(a, b)
    if a and b and #a == #b then
        local ret = {}
        
        for i = 1, #a do
            table.insert(ret, a[i] + b[i])
        end
        
        return ret
    end
end

-- vector sub
function sub(a, b)
    if a and b and #a == #b then
        local ret = {}
        
        for i = 1, #a do
            table.insert(ret, a[i] - b[i])
        end
        
        return ret
    end
end

-- dot product
function dot(a, b)
    if a and b and #a == #b then
        local ret = 0
        
        for i = 1, #a do
            ret = ret + a[i] * b[i]
        end
        
        return ret
    end
end

-- magnitude
function mag(a)
    local val = dot(a, a)
    if val then
        return math.sqrt(val)
    end
end

-- normalize, do not change param
function norm(a)
    local magnitude = mag(a)
    if magnitude and magnitude ~= 0 then
        local normalize = {}
        
        for i = 1, #a do
            table.insert(normalize, a[i] / magnitude)
        end
        
        return normalize
    end
end

-- project a to b
function proj(a, b)
    if a and b and #a == #b then
        local norm_b = norm(b)
        local val = dot(a, norm_b)
        
        if val then
            local projection = {}
            
            for i = 1, #norm_b do
                table.insert(projection, norm_b[i] * val)
            end
            
            return projection
        end
    end
end

-- perpendicular a to b
function perp(a, b)
    local projection = proj(a, b)
    if projection then
        return sub(a - projection)
    end
end

-- gram schmidt
function gram_schmidt(...)
    local vecs = { ... }
    local ret = {}
    
    if #vecs > 0 then
        table.insert(ret, vecs[1])
    end
    
    for i = 2, #vecs do
        local base = vecs[i]
        for j = 1, i - 1 do
            base = sub(base, proj(vecs[i], vecs[j]))
        end
        table.insert(ret, base)
        vecs[i] = base
    end
    
    return table.unpack(ret)
end

-- use to check gram schmidt result
function check_perp(...)
    local vecs = { ... }
    for i = 1, #vecs - 1 do
        for j = i + 1, #vecs do
            local val = dot(vecs[i], vecs[j])
            if math.abs(val) > 0.001 then
                return false
            end
        end
    end
    
    return true
end