All Questions
848
questions
7
votes
1answer
2k views
Solving for null space of a matrix with mkl LAPACK
I want to find a solution for $xA=0$, where $A$ is a square matrix. I know that most of the LAPACK routines solve for $Ax=b$. So I take $A^T$ as a, and set $b=0$. I have an additional restriction of $\...
6
votes
1answer
1k views
method of frozen coefficients and its relation to von Neumann stability analysis
I am considering two equations $$u_t=a(x)u_{xx}$$ and $$v_t=b(x)v_x$$ as classical representatives of the parabolic and hyperbolic family of equations. If $a(x)=a$ and $b(x)=b$ were constants, to show ...
6
votes
3answers
692 views
role of initial guess for iterative linear solver
Suppose we use a preconditioned iterative solver for a linear system. If the initial state for the solver can be chosen very close to the exact solution - does this reduce requirements for the ...
6
votes
0answers
465 views
How to do FEM in sector elements?
Suppose we have an element in the polar coordinate as $0\leq r\leq 1$, $0\leq\theta\leq \alpha$. What are the correct basis functions in this element? For this kind of mesh, there
are a lot of ...
6
votes
3answers
147 views
Accurate evaluation of the sign of a polynomial
Let $p$ be a polynomial with floating-point coefficients and let $a$ be a floating-point value.
Is there a method for accurately evaluating the sign of $p(a)$ in floating-point arithmetic?
I don't ...
5
votes
2answers
476 views
How to parallelize a banded direct solver?
I have a linear system whose matrix that is diagonally dominant, non-symmetric, but banded. Since the band-radius is 2 (producing only 5 variables per equation), a banded direct solver (gaussian ...
5
votes
1answer
2k views
PDE - Conservative form, conservative methods and discrete conservation
I cannot find a reference explaining clearly and rigorously the links between the notions of conservative form for a PDE, a conservative numerical method and discrete conservation.
I would be very ...
4
votes
1answer
300 views
adjoint method for reaction-diffusion problem
I'm trying to code a parameter estimation for a reaction-diffusion problem. Namely, knowing the distribution of tumor density $u$ at time $0$ and $T_f$ ($u^0$ and $u^f$), what are the best ...
4
votes
3answers
543 views
Quality Measures for Various Pseudo-Random Number Generators
According to this paper,
Ideally, a pseudorandom number generator would produce a stream of
numbers that:
are uniformly distributed,
are uncorrelated,
never repeats itself,
...
4
votes
2answers
7k views
Error in result of finite-difference approximation when refining
I have calculated the first derivative of following equation using Euler method (first order), Three point Finite Difference method (second order) and Four point Finite Difference method (third order)....
4
votes
1answer
336 views
Finite Elements Weak Formulation generalization
I am struggling with an equation that represents the Weak form of Galerkin method:
$ \phi^{T}F(\textbf{u})\sim \int_{\Omega}^{ } \phi.f_{0}(\mathit{u},\nabla \mathit{u}) + \nabla\phi:f_{1}(\mathit{u},...
4
votes
1answer
2k views
Is there a mesh generator that will generate zero thickness elements for interfaces?
I have a geometry and a finite element method where there are a couple internal interfaces across which I want to enforce some fairly simple jump conditions.
One recommendation I've gotten has been ...
3
votes
1answer
1k views
Implementation of gradient zero boundary conditon in advection-diffusion equation
My question is about Finite Element Method.
I want to know how to implement "gradient zero" conditions to advection-diffusion equations in conservative form like,
$\frac{\partial \rho}{\partial t} + ...
2
votes
1answer
3k views
Introduction to Lattice Boltzmann methods [closed]
I am trying to learn the Lattice-Boltzmann method and was looking for some good beginner resources explaining the method. I have been looking at some codes online, but have been having trouble ...
2
votes
1answer
223 views
Imposing boundary conditions for PDE quadratic eigenvalue problem
I have a quadratic eigenvalue problem of the form:
$$(A_2 s^2 + A_1 s + A_0)\hat{v} = 0$$
where $s$ is the eigenvalue. The matrices $A_i$ contain derivatives up to order six, and I have six boundary ...
2
votes
1answer
217 views
Physical interpretation of divergence theorem
In a diverging pipe section like the following,
the pipe of radius $r$ splits into two pipes of radius $r/2$.
Consider a solute transported by convection from node 1.
$$\frac{\partial C}{\partial t}...
2
votes
1answer
134 views
Calculation of the EFIE integral
I need help computing the following integral:
$$
\int_{}\frac{(1+jk|\vec{r}-\vec{r}^\prime|)e^{-jk|\vec{r}-\vec{r}^\prime|}}{|\vec{r}-\vec{r}^\prime|}d\vec{r}^\prime
$$
in this integral $\vec{r}$ ...
2
votes
0answers
188 views
Heisenberg Model python : Specific heat capacity for spin 2
I have the correct plot for specific heat capacity when I am using the formula which is $C_V$ = differentiation of entropy with respect to temperature. However, When I try to calculate $C_V$, by using ...
1
vote
1answer
429 views
FEM on tet10 element: negetive determinant at the Gauss point
I am trying to implement a fem code on tet10 elements. I follow the lecture notes for tet10 implementation given in
http://www.colorado.edu/engineering/CAS/courses.d/AFEM.d/AFEM.Ch10.d/AFEM.Ch10.pdf
...
1
vote
1answer
292 views
Calculate inverse of dense matrix with entries of very different magnitude
I need to calculate the inverse of a dense matrix, with some elements taking values as high as 1e9 and some around 1e2. What would be the best method to do it?
Note:
I am more concerned about the ...
0
votes
2answers
145 views
Adaptive numerical integration of a univariate vector integrand
Background & Problem formulation
I'm trying to write a simple program in C++ that performs adaptive numerical integration of vector valued integrands (in one variable), i.e.
$$\int_a^b \bar{f}(...
0
votes
0answers
2k views
Finite differences and Neumann boundary conditions
I am dealing with a highly nonlinear system of two PDEs. I already have a code to solve the system in case of Dirichlet boundary conditions. The explicit system is:
$$
\begin{eqnarray*}
\partial_{t}u ...
11
votes
1answer
347 views
For software submitted to ACM TOMS, how does the ACM software license agreement interact with other licenses?
The journal Association for Computing Machinery Transactions on Mathematical Software (ACM TOMS) publishes many articles on numerical algorithms that include software implementations. According to ...
8
votes
1answer
1k views
Shortley-Weller finite difference method
can you give me a link for a good and simple explanation of the Shortley-Weller finite-difference scheme? I tried to google it but all I get is (inaccessible) academic publications. I also tried ...
8
votes
2answers
225 views
How should I report profiling/timing information about my code?
I've seen a lot of publications in Computational Physics journals use different metrics for the performance of their code. Especially for GPGPU code, there seems to be a great variety of timing ...
7
votes
3answers
289 views
Is there a minimum angle requirement for cells in the finite volume method?
In his talk "What is a good linear finite element?", Shewchuk states that small dihedral angles in linear tetrahedra elements cause ill-conditioning of the stiffness matrix.
Do small dihedral angles ...
7
votes
3answers
541 views
C - OpenMP, MPI, Serial Program
I'm part of a Computational Science course and come from a non-programming background, so please forgive me my ignorance. I'm working on a set of code in C to numerically solve the Navier Stokes ...
7
votes
1answer
275 views
How to avoid overflow error in program that computes product of two numbers, such that when one is big enough to cause overflow, other is $0$?
Let us say that I have a function like so:
def f(x):
return g(x)*h(x)
Now, g(x) and ...
6
votes
3answers
4k views
How to find QR decomposition of a rectangular matrix in overdetermined linear system solution?
While trying to find cell-centered gradients in finite volume method computation of incompressible fluid flow I get over-determined linear system. This is a well known "cell based least-square" ...
6
votes
2answers
657 views
Perturbation of Cholesky decomposition for matrix inversion
I am looking for a computationally cheap way to compute $x$ such that $$(L L^T + \mu^2 I)x = y$$
where $L \in \mathbb{R}^{n \times n}$ is a lower triangular definite positive matrix (with some very ...
5
votes
2answers
326 views
Approximating and visualizing basins of attraction
I am working on estimating the position and orientation (pose) of a model (rigid object) from its silhouette in an image. For this, I have constructed an error measure between the model in its pose ...
5
votes
3answers
3k views
Discrete Poisson Equation with Pure Neumann Boundary Conditions
I'm trying to implement the Helmholtz-Hodge Decomposition in 2D, which states that a vector field is composed by a rotational free component, a divergence free component and a harmonic component.
...
5
votes
2answers
500 views
What does symmetrize mean? (imposing multifreedom constraints to stiffness matrix)
I have a small FEM implementation program. And I want to add imposing multifreedom constraints (MFC) feature to it. The theory of master-slave method is given here (page 10 for general case).
...
4
votes
1answer
2k views
Partial trace algorithm (original)
In general, is there a partial trace algorithm (ideally for systems of any size) that can be coded using basic matrix operations found in software like Mathematica or Maple? All of the methods I'm ...
4
votes
4answers
413 views
Role of weight function in Galerkin methods
I have difficulties in understanding the role of the weight function $w(x)$ that occurs in the solution of PDEs via the Galerkin approach. Consider a linear differential equation of the form
$$
\...
4
votes
1answer
598 views
LU Decomposition of PSD Matrix + Diagonal Matrix
If I have a psd, symmetric matrix $\mathbf{A}$ and I need to do LU decomps on $\mathbf{B_i}= \mathbf{A} + \mathbf{D_i}$ (where $\mathbf{D_i}$ is a diagonal psd matrix, where $\mathbf{D_i}$ changes ...
4
votes
1answer
93 views
Linear programming boundedness
Assume the optimal value of a primal problem is bounded. Is the following statement true?
If the primal problem is bounded, then its dual problem is bounded as well.
4
votes
1answer
372 views
Interpolating a mathematical function using a Hermite Cubic Finite Element Space
I have a Hermite Cubic Finite Element Space on a computer in the form of Matlab m-files. More specifically, I can evaluate four "shape functions" $N_1, N_2, N_3,$ and $N_4$, for which the following ...
4
votes
1answer
146 views
Fast and free server for computing
I have to calculate a huge differential equation. With my laptop, it's going to be computed for several days.
Is there a free (I need just for 3 days) fast server for scientific calculations?
My ...
4
votes
1answer
177 views
Compute hypergeometric function ratio: $\frac{_{2}F_{1}(a+1,b;c;x)}{_{2}F_{1}(a,b;c;x)}$?
I need a numerically stable way to compute the following ratio:
$$\frac{_{2}F_{1}(a+1,b;c;x)}{_{2}F_{1}(a,b;c;x)}$$
All the parameters are real numbers, with $a< 0$,$\ $ $b,c > 0$ and $0<x&...
3
votes
2answers
386 views
Get equation for a curve which intersects x at seemingly randomly distributed points?
Is there any type of function that when graphed would show a curve which intersects the x axis multiple times, with each point being an arbitrary distance from the last?
I mean, not like a trig ...
3
votes
1answer
213 views
Finding zeroes of an infinitely differentiable function of ~100 to ~1000 variables
I have a function that is not only infinitely differentiable, but it is also very cheap to calculate any of those derivatives. It looks like:
$f(\boldsymbol{C}, \boldsymbol{x})=\sum_{i} C_{i} \prod_{...
3
votes
1answer
5k views
Why a finite difference scheme would give second order of accuracy in norm L2 but 1.5 with L1 (while 1 with Linf)?
My finite difference scheme for the 2D Euler equations is second order accurate in theory, since all the terms are second order accurate, with the advective terms being third order. So I expect a rate ...
3
votes
1answer
258 views
non-smooth convex c++ solver
I happened to know that there are advanced established techniques for non-smooth convex optimization in research. For example, these two papers:
Nesterov, "Smooth minimization of non-smooth functions"...
3
votes
1answer
199 views
How many operations are needed for LAPACK's zgesv to solve a linear system?
I have a linear system of complex numbers. I am using LAPACK' zgesv (actually I am using intel MKL LAPACKE, but I am assuming the algorithm is the same). No assumption can be made about the system.
I ...
3
votes
2answers
543 views
DST using FFT routine
Please can you help me with my problem? On Wikipedia, in article Discrete sine transform, this is written (chapter Computation):
"Although the direct application of these formulas would require O(N2) ...
3
votes
2answers
735 views
NP-Completeness
Consider an instance of non-convexoptimization problem:
It seems that this problem is NP-complete. How can I find a suitable reduction for this?
3
votes
0answers
219 views
logsumexp with one very large term and many very small terms
I want to compute an expression of the form:
$$L = \ln\sum_i e^{x_i}$$
Suppose that there are many small terms, say $e^{x_i} \approx \epsilon$. If there are $N_\epsilon$ such terms, their ...
2
votes
2answers
169 views
Preconditioner for scalar laplacian system
Suppose that I have a large (on the order of 10^6 unknowns) 3D scalar Poisson system which I would like to precondition. The boundary conditions have been treated so that the system is SPD. I.e.,
$$\...
2
votes
2answers
197 views
Weak form of the Navier-Cauchy equation
I am trying to obtain the weak form of the Navier-Cauchy equation, which is
$$- \rho \omega ^2 \textbf{U} - \mu \nabla ^2 \textbf{U} - (\mu + \lambda) \nabla (\nabla \cdot \textbf{U}) = \textbf{F}$$
...