标签: decomposition

I want to implement my own LU decomposition P,L,U = my_lu(A), so that given a matrix A, computes the LU decomposition with partial pivoting. But I only know how to do it without pivoting. Can anyone help to do the partial pivoting?

def lu(A):

    import numpy as np

    # Return an error if matrix is not square
    if not A.shape[0]==A.shape[1]:
        raise ValueError("Input matrix must be square")

    n = A.shape[0] 

    L = np.zeros((n,n),dtype='float64') 
    U = np.zeros((n,n),dtype='float64') 
    U[:] = A 
    np.fill_diagonal(L,1) # …

他们中的任何一个都暗示另一个吗？

我的逻辑是，如果保留所有依赖项，则不会丢失信息，同样，如果分解是无损的，则一定不会违反功能依赖项。

所以本质上，依赖保留是一种确保您的分解无损的方法。

我很难接受/拒绝它。那么这两者是相互保证的，还是有一种情况可以在没有另一个的情况下实现？

我最近使用 JVisualVM 分析了一些代码，发现一种特定方法占用了大量执行时间，原因是频繁调用和执行时间缓慢。该方法由一大块 if 语句组成，如下所示：（在实际方法中大约有 30 个）

    EcState c = candidate;

    if (waypoints.size() > 0)
    {
        EcState state = defaultDestination();
        for (EcState s : waypoints)
        {
            state.union(s);
        }
        state.union(this);
        return state.isSatisfied(candidate);
    }

    if (c.var1 < var1)
        return false;
    if (c.var2 < var2)
        return false;
    if (c.var3 < var3)
        return false;
    if (c.var4 < var4)
        return false;
    if ((!c.var5) & var5)
        return false;
    if ((!c.var6) & var6)
        return false;
    if ((!c.var7) & var7)
        return false;
    if ((!c.var8) & var8)
        return false; …

真的很困惑为什么使用RcppArmadillo的QR输出与R的QR输出不同; 犰狳文献也没有给出明确的答案.基本上当我给R一个n*q的矩阵Y(比如1000 X 20)时,我得到的Q是1000 X 20和R 20 X 1000.这就是我需要的.但是当我在Armadillo中使用QR求解器时,它会让我回到Q 1000 X 1000和R 1000 X 20.我可以调用R的qr函数吗？我需要Q有维度nxq,而不是qx q.下面的代码是我正在使用的(它是更大功能的一部分).

如果有人可以建议如何在RcppEigen中做到这一点,那也会有所帮助.

library(inline)
library(RcppArmadillo)

src <- '
    Rcpp::NumericMatrix Xr(Xs);
    int q = Rcpp::as<int>(ys);

    int n = Xr.nrow(), k = Xr.ncol();
    arma::mat X(Xr.begin(), n, k, false);

    arma::mat G, Y, B;

    G = arma::randn(n,q);

    Y = X*G;

    arma::mat Q, R;
    arma::qr(Q,R,Y);

    return Rcpp::List::create(Rcpp::Named("Q")=Q,Rcpp::Named("R")=R,Rcpp::Named("Y")=Y);'


rsvd <- cxxfunction(signature(Xs="numeric", ys="integer"), body=src, plugin="RcppArmadillo")

我正在解决一个编程练习,并遇到了一个问题,我无法满意地找到解决方案.问题如下:

Print all unique integer partitions given an integer as input.
Integer partition is a way of writing n as a sum of positive integers.

例如:输入= 4然后输出应为输出=

  1 1 1 1
  1 1 2
  2 2
  1 3
  4

我该如何考虑解决这个问题？我想知道使用递归.任何人都可以为我提供这个问题的算法吗？或暗示解决方案.对这类问题的任何解释都是受欢迎的.(我是编程世界的初学者)谢谢!!

我正在使用choleski分解来计算半正定矩阵的逆矩阵.然而,当我的矩阵变得非常大并且其中有零时,我得知我的矩阵不再(从计算机的角度来看数字)正定.因此,为了解决这个问题,我使用pivot = TRUE了choleski命令中的选项R.但是,(如下所示)两者返回相同的输出,但行和列或矩阵重新排列.我试图找出是否有一种方法(或转换)使它们相同.这是我的代码:

X = matrix(rnorm(9),nrow=3)
A = X%*%t(X)

inv1 = function(A){
    Q = chol(A)
    L = t(Q)
    inverse = solve(Q)%*%solve(L)
    return(inverse)
}


inv2 = function(A){
    Q = chol(A,pivot=TRUE)
    L = t(Q)
    inverse = solve(Q)%*%solve(L)
    return(inverse)
}

运行时导致:

 > inv1(A)
              [,1]      [,2]      [,3]
    [1,]  9.956119 -8.187262 -4.320911
    [2,] -8.187262  7.469862  3.756087
    [3,] -4.320911  3.756087  3.813175
    > 
    > inv2(A)
              [,1]      [,2]      [,3]
    [1,]  7.469862  3.756087 -8.187262
    [2,]  3.756087  3.813175 -4.320911
    [3,] -8.187262 -4.320911  9.956119

有没有办法让两个答案匹配？我想inv2()从中回答答案 …

我最近读到了如何使用Choleski分解计算QR分解的R矩阵.关系是:

R = Choleski分解(A ^ TA)

例:

> A=matrix(c(1,2,3,2,3,5,1,3,2), nrow=3)
> A
     [,1] [,2] [,3]
[1,]    1    2    1
[2,]    2    3    3
[3,]    3    5    2

> AtA = t(A)%*%A
> AtA
     [,1] [,2] [,3]
[1,]   14   23   13
[2,]   23   38   21
[3,]   13   21   14

现在计算QR和Choleski分解:

> chol(AtA)
         [,1]     [,2]       [,3]
[1,] 3.741657 6.147009  3.4743961
[2,] 0.000000 0.462910 -0.7715167
[3,] 0.000000 0.000000  1.1547005

> qr_A = qr(A)
> qr.R(qr_A)
          [,1]      [,2]       [,3]
[1,] -3.741657 -6.147009 -3.4743961
[2,] …

我正在尝试使用 PCA 分解我的列。

我在如何使用 scikit learn in python 选择函数 PCA 的 n_components 方面遇到了一些困难。我做了这个

sc = StandardScaler()
Z = sc.fit_transform(X)
pca = PCA(n_components = 5')

你能解释一下吗？

一个简单的矩阵基准测试表明,Revolution Analytics R 2.13.2的LU分解比矩阵乘法慢近5倍.理论和多年的实践表明,LU应该是1/3到2/3的时间A*A.

Revo R和Matlab正在使用英特尔的数学核心进行此测试. R 2.14.1没有使用内核.一切都是64位.

异常现象见下表2.这是表1的标准化A*A.还有其他(明显的)异常,但LU是最明显的异常.

                        Table 1 (secs)

                    A*A     LU     A\b    Det   Inv
----------------------------------------------------
R 2.14.1           0.757   0.43   0.45   0.20  1.11
Revo R 2.13.2      0.063   0.35   0.11   0.03  0.14
Matlab 2011b       0.062   0.08   0.10   0.07  0.16
----------------------------------------------------
Averaged over 20 runs on a 1000x1000 random matrix


                       Table 2 (normalized)

                    A*A     LU     A\b    Det   Inv
----------------------------------------------------
R 2.14.1             1     0.57   0.19   0.26  1.47
Revol R 2.13.2       1     4.67* …

WF 在这里给出了一个现在删除的答案我的问题使用了这一行:

auto [x, y] = div_t{ 1, 0 };

从答案的代码看起来这像是tie对div_t结构.我希望有人能解释这里发生的事情.完整的功能代码如下:

constexpr bool first_quot() {
    auto [x, y] = std::div_t{1, 0};
    (void)y;
    return x;
}