Yes, the condition number always matters in floating-point arithmetic, whether you choose to solve your system with an iterative or direct method. The relative accuracy of an approximate solution to $... View answer Accepted answer 19 votes To add to Lutz Lehmann's answer, you can look up the latency for the CPU instructions in this comprehensive table by Agner Fog. For example, on the Intel Ivy Bridge processors: FADD / FSUB (floating ... View answer Accepted answer 15 votes In general, there is no shortcut other than completely re-factoring the matrix. There have been a few similar questions on this SE that cover the topic in more depth than I can: Can diagonal plus ... View answer Accepted answer 14 votes Chapter 8 of Brenner and Scott's Mathematical Theory of Finite Element Methods is devoted to this subject. In particular, theorem 8.1.11 and the corollary give you that$ \|u - u_h\|_{W^1_\infty} \le ...

The conjugate gradient method is the provably fastest iterative solver, but only for symmetric, positive-definite systems. What would be awfully convenient is if there was an iterative method with ...

The degeneracy of some eigenvalues looks to me like the hallmark of the breakdown of the Lanczos algorithm. The Lanczos algorithm is one of the more commonly used methods to approximate the ...

I'm going to assume since you mention electrodynamics that you're interested in PDEs. You've already mentioned FEniCS in your comment. FEniCS offers a domain-specific language (DSL) called UFL. This ...

Suppose that the solution $u$ of the PDE lives in some function space $X$. We'll write the PDE as a bilinear form $A(u, v) = f(v)$ for all $v$ in $X$, where $f$ is some element of the dual space X^*... View answer 12 votes The conjugate gradient method is good for finding the minimum of a strictly convex functional. This is typical when you reformulate a nonlinear elliptic PDE as an optimization problem. If you want to ... View answer Accepted answer 11 votes If it's shock-capturing that you're interested in, I would suggest you use the finite volume method instead of the finite element method. When applied naively, FEM is actually notoriously bad at ... View answer 10 votes With your formulation of the constraints, you can only produce ellipses where the semi-major and semi-minor axes are aligned with the coordinate axes, but it's clear from the figure that you attached ... View answer Accepted answer 10 votes If the operation is as trivial as you say, and each node has all the information necessary to carry out the operation, then the communication will be substantially more expensive than recomputing ... View answer Accepted answer 10 votes You may find the expansion of a plane wave in spherical waves to be helpful here:$e^{i\mathbf{k}\cdot\mathbf{x}} = 4\pi\sum_{l=0}^\infty\sum_{m=-l}^l i^lj_l(kr)Y_{lm}(\theta,\phi)Y_{lm}^*(\... View answer 10 votes Generally speaking, you'll want to use an implicit method for parabolic equations (the diffusion part) -- explicit schemes for parabolic PDE need to have a very short timestep to be stable. Conversely,... View answer 9 votes This paper by Kirby and Mitchell describes the implementation of$C^1$elements in the Firedrake package*. One of the main use cases is biharmonic problems, which show up in the elastic deformation of ... View answer 8 votes That depends on how well you know the coordinates and velocities. If you have exact values, you can get a reasonable answer using Hermite interpolation. This will give you a degree-3 polynomial in ... View answer Accepted answer 8 votes There's a simple test to see if a point$(x, y)$is enclosed within a curve. Draw a ray from$(x, y)$to infinity, and count how many times it crosses the curve; if the count is odd, then$(x, y)$is ... View answer 8 votes Often$M^{-1}A$is not symmetric, even if$M$and$A$are. There are two common approaches to dealing with this: Find a Cholesky factorization of$M$into$LL^\top$, and instead solve for$L^{-1}AL^{-...

Using an Eigen matrix type where the number of rows and columns is encoded into the type at compile time gives you an edge over LAPACK, where the matrix size is known only at runtime. This extra ...

There is another approach called the adjoint method, which is commonly used in inverse problems for PDE and which is quite easy to generalize to other problems. This is going to be long. Your ...

Which one was faster when you tried using Eigen? You may also consider Trilinos if you're using C++. By "generic Laplacian matrix", do you mean a finite difference / finite element discretization of ...

FILTLAN is a C++ library for computing interior eigenvalues of sparse symmetric matrices. The fact that there is a whole package devoted to just this should tell you that it's a pretty hard problem. ...

Is the matrix sparse or dense? Is it symmetric? I'm assuming by "rearrange" you mean permute the entries, rather than apply some more general similarity transformation to the matrix. In that case, ...

If the only non-zero entries of $A_{ij}$ have $j$ in $\{i - 1, i, i + 1\}$, then $A$ is a banded matrix with bandwidth 1. More generally, you can talk about matrices of bandwidth $k$ where $k$ is any ...
As Brian Borchers pointed out, the question only makes sense up to factors of $2\pi$. I'll assume that you want an answer in $[0, 2\pi)$ then. If there are really big contrasts in the magnitudes $|... View answer 6 votes There are multiple ways of implementing Dirichlet boundary conditions when using the finite element method. For each unknown$i$of the system belonging to the Dirichlet boundary, you can zero out ... View answer Accepted answer 6 votes Short answer: yes, you have to use something different, e.g. a Neumann-Neumann method. A good reference is Widlund's book. Non-overlapping methods are based on the principle that, if$u$solves the ... View answer 6 votes Trefethen and Bau's book Numerical Linear Algebra has the Householder QR algorithm in chapter 10, and it's written considering general rectangular matrices. It's also in Matrix Computations by Golub ... View answer 6 votes The characteristic polynomial of a matrix$M$can be written as$\chi_M(z) = z^n + \textrm{trace(M)}z^{n-1}+\ldots + \det(M)$. Since you know ahead of time that all the entries of$M$are rational, ... View answer 6 votes Usually you would do one of two things: Compute$z = Ax$, then$y = y+Az$; If you anticipate needing to do this operation many times, you might compute the matrix$B = AA'$and then evaluate$y = y+...