原创

Jacobi并行拆解

版权声明:本文为博主原创文章,遵循 CC 4.0 BY 版权协议,转载请附上原文出处链接和本声明。
本文链接:https://wqqeqeqwe.blog.csdn.net/article/details/80479036

作者:桂。

时间:2018-04-23  21:12:02

链接:http://www.cnblogs.com/xingshansi/p/8921815.html 


前言

本文主要是复数矩阵分解的拆解思路(矩阵求逆/特征值分解)一文的补充。

一、并行拆解思路

   回顾前文,对于8X8的实数矩阵:

仿真:

clc;clear all;
X = rand(8);
R = X+X';
Iteration = 20;
[D,U] = Jac_sweep(R,Iteration);
[u,s,v] = svd(R);
[sort(diag(s)),sort(abs(diag(D)))]

  Jac_sweep并行代码:

function [D,U] = Jac_sweep(A,Iteration)
iter = 0;
n = size(A,1);
U = eye(n);
A_c = [1,2;3,4;5,6;7,8;...
    2,3;4,5;6,7;1,8;...
    1,3;2,4;5,7;6,8;...
    3,5;4,6;1,7;2,8;...
    1,4;3,6;5,8;2,7;...
    2,5;4,7;1,6;3,8;...
    1,5;2,6;3,7;4,8];

while iter <Iteration
    iter = iter+1;
    for flag = 1:size(A_c,1)/4
        T = eye(n);
        for t = 1:4
            p = A_c((flag-1)*4+t,1);q = A_c((flag-1)*4+t,2);
            y = 2*A(p,q);
            x = A(p,p)-A(q,q);
            phi = atan2(y,x)/2;

            T(p,p) = cos(phi);
            T(q,q) = cos(phi);
            T(p,q) = -sin(phi);
            T(q,p) = sin(phi);
        end
        D = T'*A*T;
        U = U*T;
        A = D;
    end
end

  

打印结果与SVD分解的结果一致,并行思路可行。其他维度依次类推。对于维度N的矩阵,

  • N为偶数,可并行N/2路;
  • N为奇数,可并行[N-1]/2路;

二、改进思路

  每一次sweep,需要1次Cordic:phi = atan2(y,x)/2,和两次Cordic(两次可并行):cos(phi) / sin(phi),二者串行。对于atan2操作,可借鉴复数相位近似估计一文的思路,即对于atan的计算,考虑到CORDIC耗时较长,内存资源充足的情况下,1)直接查表;若内存相对紧张,2)多项式逼近。二者较CORDIC均减少运算时间。

  另外,关于定点仿真可调用fi工具包,其中CORDIC对应指令:cordicatan2、cordiccos......

文章最后发布于: 2018-04-23 21:43:00
展开阅读全文
评论将由博主筛选后显示,对所有人可见 | 还能输入1000个字符

没有更多推荐了,返回首页

分享到微信朋友圈

×

扫一扫,手机浏览