12
votes
Accepted
How to document math formulations in scientific computing codes?
I prefer using doxygen that supports C++ and LaTeX comments, both inline and as separate equations.
This way, you will keep your comments, including, say, the rigorous mathematical formulation of the ...
- 8,602
11
votes
Accepted
Memory Management - Why certain initializations order are faster?
Such an effect happens because of how the data of the int** a is stored in memory (as per C/C++). This question on StackOverflow has answers with some more details (...
- 8,602
10
votes
Accepted
Numerical evaluation of the Exponential Integral Ei by rational Chebyshev approximations fails
There is a prime missing both in your code and in the expression in your question. In the original paper, the expression is:
$$ \log(x/x_0) + (x-x_0) \frac{\sum^\prime_{0\leq j\leq n}p_j T_j^*(x/6)}{\...
- 11.4k
9
votes
Accepted
C standard for computational science
In theory, as the original authors, you're free to pick and name a standard, then expect others to follow it. In practise, if you're supporting an HPC system, then your choice is likely to be ...
- 2,229
8
votes
Accepted
What PRNG function is this?
A little playing with the sequence of numbers generated by the C code shows that the sequence is
$z_{i+1}=5z_{i}+273 \mod 2^{16}$
This is a linear congruential generator (LCG). It's easy to show ...
- 18.5k
8
votes
Accepted
Suming up the series $1+1/2+...+1/n$ in C
I see several problems with your code:
type error in scanf (n is an integer, not a float)
i=0 and not i==0 in the for loop
...
- 196
7
votes
Numerical integration giving incorrect sign
This might be an accuracy problem in computing the second term, because of those large exponentials when $x \gg 1$.
I would first work on that term: gather $e^x$ out from numerator and denominator and ...
- 11.1k
6
votes
Looking for Runge-Kutta 8th order in C/C++
summarizing some points:
If it's a long-term integration of a non-dissapative model, a symplectic integrator is what you're looking for.
Otherwise, since it's an equation of motion, Runge-Kutta ...
- 12.3k
6
votes
Do BLAS routines compute their respective operations with minimum error?
BLAS routines do not typically use stable summation algorithms.
In the case of gsl, you can look up its source code online - the source of gsl's sdot is contained ...
- 11.4k
6
votes
Accepted
Using GSL for basic operations
Of course it makes sense to use the GSL (or another library for that matter) for several reasons:
Don't reinvent the wheel. The work has been done, you can spend your time on more useful things.
If ...
- 6,169
6
votes
I've developed a derivative-free optimization method, looking for comments
I would like to hear comments from users that have some practical models (e.g. black-box hyperparameter optimization) which are still needed to be solved acceptably - whether this method works or not ...
- 2,176
6
votes
C standard for computational science
You should definitely jump to C99, or newer(!). The C99 standard introduced
the restrict keyword. Loosely speaking,
with this keyword you can inform the compiler ...
- 561
6
votes
Accepted
What language are best for creating data compression algorithms except for C/C++?
In my general experience, the following are required in order to get "reasonable" performance from a programming language:
The ability to specify static types. It's OK if the language also ...
- 2,391
5
votes
Kummer's confluent hypergeometric for complex arguments in C/C++?
One possibility is to wrap the arbitrary precision interval implementation in Arb. This will not be as fast as a dedicated double precision implementation, but it might still be fast enough.
Note 1: ...
5
votes
Accepted
How to deal with big numbers in intermediate calculations?
I think there's a simple way to do this. You have a rational function of identical cosh/sinh terms, where every expression is a homogeneous polynomial in cosh/sinh, and the only problem is that these ...
- 11.4k
5
votes
Accepted
Recommendations for ODE solvers for stiff equations
This is a huge open ended question, but I'll copy the recommendation section from the current release of DifferentialEquations.jl:
Stiff Problems
For stiff problems at high tolerances (>1e-2?) it ...
- 12.3k
5
votes
float in C: Size of the exponent
$1.2\times 10^{-38}$ is a decimal value. You need to convert it to binary exponential notation to see how the bits are set.
Another small mistake: $2^6+...2^1+2^0$ does not equal $2^{27}$, it equals $...
- 11.1k
4
votes
solver linear system equation
You haven't told us whether your matrix is positive definite or not- this is a very important factor in selecting a solver. You also haven't told us whether the matrix is well conditioned or poorly ...
- 18.5k
4
votes
Accepted
What do C, C++ and Java have that Fortran 2003 don't?
One thing that C++ includes that Fortran doesn't have is extensive support for generic programming and compile-time code ...
- 534
4
votes
Jacobi iteration for finite difference: when to stop?
For solving systems you shouldn't be comparing the results between each iteration but rather computing the residual.
If you consider the matrix representation to be in the standard form:
$ A\cdot y ...
- 407
4
votes
Accepted
C or fortran library to solve linear 2D/3D elliptic PDE
Most of the widely used finite element libraries are written in C++. If all you really care for -- and if all you will ever care for -- is solving an elliptic PDE on a rectangle, then it's probably ...
- 54.2k
4
votes
Accepted
Number of decimal of float and double in C
Prof. Bangerth's comment is completely correct. To add more detail, we can refer to IEEE754 standard which defines the floats as
...
- 2,411
3
votes
Accepted
Linear solve using CHLOMOD in C
Here is how I call CHOLMOD from a C program. It uses as input my own C++ structure for a sparse matrix (in standard Compressed Row Storage format).
There are several gotchas: understanding the ...
- 2,305
3
votes
Accepted
Implementing velocity verlet for harmonic oscillator in C gives error wrt conservation of energy
Your second acceleration ao in the else statement does not make sense. Since it assumes that the calculated position and ...
- 450
3
votes
Programming language for HPC
It is probably a good idea to stick with Python if you have expertise with that. I am assuming that this means you are using numpy/scipy. Depending on what you exactly want to do, well-written numpy/...
- 131
3
votes
Accepted
Sum over very small exponentials: Underflow
You may use Kahan summation algorithm [1]
The idea is to reschedule the sum operations in such a way precision loss is limited. The code is very simple (reproduced from [1] below).
If this does not ...
- 2,305
3
votes
Looking for Runge-Kutta 8th order in C/C++
I would like to add that while what Geoff Oxberry suggests for long-term integration (using symplectic integrators) is true, in some cases it won't work.
More specifically, if you have dissipative ...
- 45
3
votes
Decrease execution time using openMP
The algorithm you are using spawns extremely cheap tasks, especially at small $n$. As you don't have a luxury of changing the algorithm, the only reasonable (relatively) way to avoid the overhead of ...
- 8,602
3
votes
Accepted
How to deal with very low numerical values in C?
The best thing to do would be to re-scale your problem so that way you don't have such small values. If you can, do this and call it a day.
If every number is small, then you can still be fine ...
- 12.3k
3
votes
Accepted
Simple open-source Quantum chemistry or DFT code in C/C++
Quantum Chemistry codes can get very complicated very fast. Even if you limit yourself to DFT, there are many functionals to support. There's going to be a trade off here for you. You can get a simple,...
- 181
Only top scored, non community-wiki answers of a minimum length are eligible