Skip to main content
The 2024 Developer Survey results are live! See the results
21 votes

What programming language should I choose and why?

You should definitely check out Julia. Julia is a programming language which is similar to Python or MATLAB but utilizes a strong type-inference algorithm + JIT in order to optimize code. If types can ...
Chris Rackauckas's user avatar
13 votes
Accepted

Why is the central difference method dispersing my solution?

I'll write the equation short as $$\ddot x(t)+c\dot x(t)=a(t,x(t))$$ to separate the "easy" linear parts from the non-linear and forcing terms. On the first method The claimed order of the ...
Lutz Lehmann's user avatar
  • 6,129
12 votes
Accepted

When and when not to use automatic differentiation

Given code that computes a function $f(x)$, automatic differentiation tools produce a code that can compute $f(x)$ and its derivatives at the same time. Solving a differential equation is an entirely ...
Brian Borchers's user avatar
12 votes

Finite-difference software for solving custom equations

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 ...
Daniel Shapero's user avatar
11 votes

What programming language should I choose and why?

What is it you want to achieve? If you want to develop methods/algorithms you might prefer a language that is flexible, and that you are familiar with. As you stated in your question, the Fortran ...
Dohn Joe's user avatar
  • 559
11 votes

What programming language should I choose and why?

Start simple. Learn Python. I have been paid to write programms for over forty years and I have used all the languages mentioned in other answers (except Julia - I had never heard of it before now). ...
Paul Smith's user avatar
11 votes

Is using iterative methods to solve a linear system always superior to inversing the matrix?

First off, there are basically no scenarios where one would ever actually compute and store $A^{-1}$ in memory, even for small problems. An LU factorization offers both superior efficiency and ...
whpowell96's user avatar
  • 2,731
10 votes

Numerical derivative and finite difference coefficients: any update of the Fornberg method?

Overview Good question. There is a paper entitled "Improving the accuracy of the matrix differentiation method for arbitrary collocation points" by R. Baltensperger. It's no big deal in my opinion, ...
davidhigh's user avatar
  • 3,187
10 votes

Non-hermitian discretizations in quantum mechanics

If your discrete equations represent your problem properly you should have the right eigenvalues, although numerically you can get a small imaginary part when solving the eigenvalue problem. As an ...
nicoguaro's user avatar
  • 8,582
10 votes

Finite difference for 1D wave equation: why the spike initial data results in a noisy output?

I think you're probably seeing artifacts that are due to numerical dispersion. In brief, in the discrete case different (spatial) frequencies of a wave function will propagate at different phase/group ...
rchilton1980's user avatar
  • 4,946
10 votes

Why FVM can handle unstructured meshes while FDM cannot?

The difference between finite volumes and finite differences is really more about the form of the equations solved. In typical FV methods, the conservative form is ...
Darren Engwirda's user avatar
9 votes
Accepted

Correct eigenfunctions of Laplace operator by Finite Differences

You should specify the eigenvalues you want with which="SM", for example. Check the following snippet. I also changed the solver, since your system is symmetric. <...
nicoguaro's user avatar
  • 8,582
8 votes

ODE $x''(t)+\eta x'(t)+x(t)=0$ with the $\eta$ extremely small

The problem you have is what is called "stiff": it has two time scales, namely one for the oscillation (which is ${\cal O}(1)$) and one for the damping (which is ${\cal O}(\eta^{-1})$). If $\eta$ is ...
Wolfgang Bangerth's user avatar
8 votes
Accepted

Solving Ax = b with sparse A and sparse b

When looking at the solution of your system, you will find that almost all entries of $x$ are nonzero although the right-hand side is "sparse". Hence, whatever algorithm you use, it'll have to visit ...
Nico Schlömer's user avatar
8 votes

Finite Difference and Finite Volume as special cases of Finite Element

A good example is if you use bilinear finite elements for the Laplace equation on a uniform mesh, and then approximate the integrals using the trapezoidal rule, then you get the usual 3-point stencil (...
Wolfgang Bangerth's user avatar
8 votes

Why is my simulation of a first-order wave equation not stable?

Numerical solution of the advection equation with centered differences in space and forward Euler in time is unconditionally unstable. So the behavior you are seeing is expected. Here is a nice ...
Bill Greene's user avatar
  • 6,229
8 votes
Accepted

Why do many people use FDM method to solve Stokes equations, i.e., saddle point matrix?

The preconditioner for the FDM method that corresponds to the one you outline for the FEM (i.e., the Sylvester-Wathen approach) will still contain the Schur complement of the FDM matrix. The Schur ...
Wolfgang Bangerth's user avatar
8 votes
Accepted

Computation of diffusion time

It's easy to derive that equation from Fick's law. You have this diffusion equation as: $$\frac{\partial C}{\partial t} = D \nabla^{2} C$$ The mean square displacement weighted by the concentration ...
Mithridates the Great's user avatar
8 votes

Finite-difference software for solving custom equations

You might want to check out DifferentialEquations.jl. It supports ODEs, PDEs, stochastic equations, delay equations, and basically everything else. It also has really good automatic sparsity detection ...
Oscar Smith's user avatar
8 votes
Accepted

3D laplacian operator

Yes, that finite difference is correct. You can obtain it using a finite difference in each direction for the Laplace operator in each coordinate. \begin{align} \nabla^2 u =& \frac{\partial^2 u}{\...
nicoguaro's user avatar
  • 8,582
7 votes

Derivatives Approximation on non uniform grid

So not quite sure what you're missing, but here's how you go about doing this sort of thing. So first, I am assuming this is 1D due to your description. Second, I'm assuming you know the relationship ...
spektr's user avatar
  • 4,278
7 votes

ODE $x''(t)+\eta x'(t)+x(t)=0$ with the $\eta$ extremely small

Let us write the equation in state-space form \begin{equation} \frac{d}{dt}\begin{bmatrix}x_{1}(t)\\ x_{2}(t) \end{bmatrix}+\begin{bmatrix}\eta/2 & -\omega\\ \omega & \eta/2 \end{bmatrix}\...
Richard Zhang's user avatar
7 votes

Advection equation with finite difference: importance of forward, backward or centered difference formula for the first derivative

Okay. Let's begin with the first situation. Your equation is: $$ \frac{\partial f}{\partial t}+v\frac{\partial f}{\partial x}=0\tag{1}$$ $\textbf{Previous comments}$ There are plenty of webs that ...
HBR's user avatar
  • 1,648
7 votes
Accepted

Finite difference method basic implementation on Octave

Summarizing and formulating the answer: (part of it was already given in the comments by user14717, Christian Clason, and Kirill) While performing numerical differentiation using finite differences, ...
Anton Menshov's user avatar
  • 8,702
7 votes
Accepted

What is the maximum attainable accuracy with a given set of $\alpha,\beta$?

If you are looking for general accuracy of your stencil, it could be extracted by using Taylor expansion. Basically, if I want to write down your stencil explicitly, it contains terms of $u(x_{i}-2h)$,...
Mithridates the Great's user avatar
7 votes

Time complexity of numerical finite differences

As pointed out in the comments, the cost of evaluating $f$ is critical, and in most practical cases will be the dominant cost. Lets suppose it takes $C$ operations to evaluate $f$. For nontrivial ...
Steven Roberts's user avatar
6 votes
Accepted

Thomas algorithm for 3D finite difference

The problem is that 2d and 3d discretizations are not block tridiagonal. The only tridiagonal decomposition would be a 2x2 decomposition that encompasses the entire matrix. Of course, this fact is a ...
Wolfgang Bangerth's user avatar
6 votes
Accepted

Data corruption when taking gradient of numerical data in python

The $x$-coordinates of your data points are not equally spaced (x[1:]-x[:-1] is not constant), so numpy.gradient is not ...
Kirill's user avatar
  • 11.4k

Only top scored, non community-wiki answers of a minimum length are eligible