New answers tagged

0

I believe breaking your monstrous matrix into smaller blocks won't help you that much cause as I mentioned in my comment, using inverse of a huge sparse matrix is a really bad idea cause generally inverse of a huge sparse matrix could be a dense matrix, which probably you can't afford to even store it in your memory. But, in order to show you how you can ...


4

I think you already know the answer- find the inverse of each of the small diagonal blocks. The details of how you might break the matrix up into the block in Python are something out of scope for this stackexchange group. As others have suggested you most likely don't really need the inverse of this matrix and could probably do everything you need with a ...


2

By definition, the Householder transformation is a reflection about a plane or hyperplane. The plane is described by its unit normal vector u and the transformation is then H = I - 2uu'. If u is a unit vector, 2uu' cannot be zero and thus H cannot be the identity matrix. You can find plenty of non-Householder transformations H that, for certain arguments, ...


3

Short answer: When he says that when $x=0$, any $u$ will do, he means any $u$ that satisfies $\|u\|_2 = \sqrt{2}$. Otherwise, it wouldn't be a Householder reflector. So yes, you have found a transformation that "introduces" zeros in $x$, but you haven't proven that that transformation is a Householder reflector (spoiler: it isn't ;) ).


2

As I mentioned in my comment, due to that you are searching for a method based on Python and ideally available in NumPy or SciPy, my suggestion is to use scipy.linalg.expm. It's really easy to use basically you have the $X$ matrix and you just pass it to the scipy.linalg.expm class and it would give you the exponential of $X$ (i.e. $e^{X}$). Due to that you ...


2

If changes in pivoting are an issue, then yes, to my knowledge there is no obvious $O(n^3)$ solution. However, you could consider switching to the QR factorization. That factorization costs twice as much as LU for the initial computation, but it can be updated in a stable way in $O(n^2)$. See Golub and Van Loan Matrix Computations for algorithms, or https://...


19

Your code is limited by memory bandwidth. For trivial math, it's often better to count memory accesses rather than flops. You'll get the following table: operation memory reads/writes matrix + matrix 3n² matrix * vector 2n²+n (if vector is not cached) matrix * vector n²+2n (if vector is only read once) vector + vector ...


60

Except for code which does a significant number of floating-point operations on data that are held in cache, most floating-point intensive code is performance limited by memory bandwidth and cache capacity rather than by flops. $v$ and the products $Av$ and $Bv$ are all vectors of length 2000 (16K bytes in double precision), which will easily fit into a ...


2

The Newton-Raphson method can be used to solve non-linear systems of equations. The first step is to write your system as a root finding problem: $$ f(T_n) = \left( \frac{C}{\Delta t} + K \right) T_n - Q_n - \frac{C}{\Delta t} T_{n-1} = 0 $$ Taylor expand this equation about an initial guess $T^0_n$, keeping only the linear term: $$ f(T_n) \approx f(T^0_n)...


2

First of all, I find the name "Smoothed Aggregation" a bit misleading, because the method - as I understood it - consists of both smoothing a tentative prolongation operator and implicitly considering the (near-)null-space for constructing the tentative prolongation operator. In standard algebraic multigrid with linear elasticity, corrections from coarser ...


4

You might want to use the fact that: $$ ||A||_2=\sigma_\max(A) $$ where $\sigma_\max$ is the largest singular value. If you are interested in details, this Math SO question should be interesting. Thus, $$ ||A^{-1}||_2=\frac{1}{\sigma_\min(A)} $$ where $\sigma_\min$ is the smallest singular value. You certainly want to avoid the actual calculation of the ...


Top 50 recent answers are included