Sweet Snippet 之 Gram-Schmidt 正交化

悠悠 2021-10-23 13:22 285阅读 0赞

> Gram-Schmidt 正交化的简单实现

Gram-Schmidt(格拉姆-施密特) 正交化可以正交化一组给定的向量,使这些向量两两垂直,这里列出一份简单的实现(Lua):

-- 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

有兴趣的朋友可以试试这组向量的 Gram-Schmidt 正交化:

a = ( 1 , 0 , 0 , 1 ) b = ( 0 , 1 , 0 , 1 ) c = ( 0 , 0 , 1 , 1 ) d = ( 0 , 1 , 1 , 1 ) \\begin\{aligned\} & a = ( 1, 0, 0, 1 ) \\\\ & b = ( 0, 1, 0, 1 ) \\\\ & c = ( 0, 0, 1, 1 ) \\\\ & d = ( 0, 1, 1, 1 ) \\end\{aligned\} a=(1,0,0,1)b=(0,1,0,1)c=(0,0,1,1)d=(0,1,1,1)

#### 更多资料 ####

*  [wiki][]

[wiki]: https://zh.wikipedia.org/wiki/%E6%A0%BC%E6%8B%89%E5%A7%86-%E6%96%BD%E5%AF%86%E7%89%B9%E6%AD%A3%E4%BA%A4%E5%8C%96