# Tag Info

52

Python (as of 2.6 and 3.0) now searches in the ~/.local directory for local installs, which do not require administrative privileges to install, so you just need to point your installer to that directory. If you have already downloaded the package foo and would like to install it manually, type: cd path/to/foo python setup.py install --user If you are ...

40

I'm going to break up my answer into three parts. Profiling, speeding up the python code via c, and speeding up python via python. It is my view that Python has some of the best tools for looking at what your code's performance is then drilling down to the actual bottle necks. Speeding up code without profiling is about like trying to kill a deer with an ...

32

First, if your undergraduates are like ours and had no prior introduction to computers, expect to spend some time teaching them how to use basic stuff like using a proper editor (i.e., not MS Word), the command line, etc. I think the answer somewhat depends on where you set the focus of your course (or what you are required to teach). For example: How ...

31

This is indeed called catastrophic cancellation. In fact, this particular case is very easy: rewrite the function using the equivalent, numerically stable expression $$\frac{t}{1+\sqrt{1-t^2}}.$$ Since you probably need a reference, this is discussed in most numerical methods textbooks in relation to the formula for solving quadratic equations (that ...

29

Ease of learning Python and Fortran are both relatively easy-to-learn languages. It's probably easier to find good Python learning materials than good Fortran learning materials because Python is used more widely, and Fortran is currently considered a "specialty" language for numerical computing. I believe the transition from Python to Fortran would be ...

24

You have to subclass the rv_continuous class in scipy.stats import scipy.stats as st class my_pdf(st.rv_continuous): def _pdf(self,x): return 3*x**2 # Normalized over its range, in this case [0,1] my_cv = my_pdf(a=0, b=1, name='my_pdf') now my_cv is a continuous random variable with the given PDF and range [0,1] Note that in this example ...

23

I am not super familiar with f2py internals, but I am very familiar with wrapping Fortran. F2py just automates some or all of the things below. You first need to export to C using the iso_c_binding module, as described for example here: http://fortran90.org/src/best-practices.html#interfacing-with-c Disclaimer: I am the main author of the fortran90.org ...

22

In 2014, I would've said Python. In 2017, I wholeheartedly believe that the language to teach undergraduates is Julia. Teaching is always about a tradeoff. On one hand, you want to choose something that is simple enough that it is easy to grasp. But secondly, you want to teach something that has staying power, i.e. something that can grow with you. The ...

22

Joblib does what you want. The basic usage pattern is: from joblib import Parallel, delayed def myfun(arg): do_stuff return result results = Parallel(n_jobs=-1, verbose=verbosity_level, backend="threading")( map(delayed(myfun), arg_instances)) where arg_instances is list of values for which myfun is computed in parallel. The main ...

20

Going from MATLAB to Python does introduce quite a bit of syntax overhead. One way to quantify it is the nice QuantEcon cheatsheet which showcases how there's a lot of extra "stuff" going on when trying to write simple linear algebra commands in Python. The verbose NumPy syntax is really just a symptom of how it was not developed as a technical computing ...

19

You can use the Python builtin ctypes module as described on fortran90.org. It is pretty straight forward and doesn't require any external dependencies. Also, the ndpointer arg type helper is very handy.

19

I will try to answer your question considering that you are asking for Python specifically. I will describe my own method of tackling a simulation problem. Strategies for faster simulations are given in this description. First, I prototype new simulations in Python. Of course, I try to make use of NumPy and SciPy as much as I can. Whereas NumPy provides a ...

19

Let me try and break down your requirements: Maintainability Reading/writing text data Strong interfaces/capability for LU factorizations Sparse linear solvers Performance and scalability to large data From this list, I would consider the following languages: C, C++, Fortran, Python, MATLAB, Java Julia is a promising new language, but the community is ...

18

There are two issues that you are likely to be encountering. Ill-conditioning First, the problem is ill-conditioned, but if you only provide a residual, Newton-Krylov is throwing away half your significant digits by finite differencing the residual to get the action of the Jacobian: $$J[x] y \approx \frac{F(x+\epsilon y) - F(x)}{\epsilon}$$ If you ...

17

Here is the Numba solution. On my machine the Numba version is >1000x faster than the python version without the decorator (for a 200x200 matrix, 'k' and 200-length vector 'a'). You can also use the @autojit decorator which adds about 10 microseconds per call so that the same code will work with multiple types. from numba import jit, autojit @jit('f8[:]...

15

You should check out sympy.stats. It provides an interface to deal with random variables. The following example provides a random variable X defined on the unit interval with density 2x In [1]: from sympy.stats import * In [2]: x = Symbol('x') In [3]: X = ContinuousRV(x, 2*x, Interval(0, 1)) In [4]: P(X>.5) Out[4]: 0.750000000000000 In [5]: Var(X) # ...

15

There are a number of issues in your question. Do not use Gaussian Elimination (LU factorization) to calculate the numerical rank of a matrix. LU factorization is unreliable for this purpose in floating-point arithmetic. Instead, use a rank-revealing QR decomposition (such as xGEQPX or xGEPQY in LAPACK, where x is C, D, S, or Z, though those routines are ...

15

A difficulty with any of these types of questions is that the answer is highly community-dependent. To answer some of your questions in haphazard order: MATLAB is used a lot both in academia and in industry. One of the reasons it's used quite a bit in industry is because it is taught in academia. I know for a fact that MATLAB is used at Lincoln Laboratory ...

15

The question has two very different subquestions. I will address the first one only. Matlab's version runs on average 24 times faster than my python equivalent! The second one is subjective. I would say that letting know the user that there is some problem with the integral is a good thing and this SciPy behavior outperforms the Matlab`s one to keep it ...

14

For the first part of my question, I found this very useful comparison for performance of different linear interpolation methods using python libraries: http://nbviewer.ipython.org/github/pierre-haessig/stodynprog/blob/master/stodynprog/linear_interp_benchmark.ipynb Below is list of methods collected so far. Standart interpolation, structured grid: http:/...

14

Here is R1, as computed in MATLAB: 1.0e+07 * -7.382605957465515 -9.599867106092937 -2.830412177259742 -0.000000000002830 -0.000000000002830 -1.230434326244253 -1.599977851015490 -0.471735362876624 -0.000000000000472 -0.000000000000472 3.691302978732758 4.799933553046468 1.415206088629871 0.000000000001415 0.000000000001415 -5....

14

First, see Mark L. Stone's answers, which is completely correct. Second, realize that this is the reason why people told you to use relative errors in your numerical analysis class. :) Third, the real question here is why the results do not coincide exactly, since both languages call some BLAS library functions for their computations. There are several very ...

13

Algorithms for rank-1 updates of the SVD (also called incremental SVD) do exist, but I haven't been able to find a LAPACK-like implementation anywhere. The one I've seen mentioned repeatedly is that of Brand (2003). Judging from this website, it seems as though Brand's algorithm is relatively simple to implement using existing LAPACK and BLAS routines as ...

13

I will address only the comparison of C to C++. While it is true that anything written in C can be ported to C++ with a few syntactic touch-ups, the communities have different values. The C library community, more than almost any other, values binary stability. Binary stability is critical for low-level libraries to avoid inflicting constant pain on those ...

13

What you're looking for is Numba, which can auto parallelize a for loop. From their documentation from numba import jit, prange @jit def parallel_sum(A): sum = 0.0 for i in prange(A.shape[0]): sum += A[i] return sum

13

To the best of my knowledge, Numpy does not support independent streams. Indeed, getting independent streams from the Mersenne Twister (Pythons RNG) is notoriously difficult although it can be done. Consider using the RandomGen package. It is fully compatible with Numpy, and provides you with the PCG64 generator, supporting up to $2^{63}$ independent ...

13

The problem you are encountering is likely not a consequence of your choice of algorithm, but in fact a consequence of the resulting dynamical system after applying time reversal. Per the definition of an attractor, all points in some neighborhood of the attractor will converge to the attractor under the flow of the dynamical system as $t\to\infty$. However, ...

13

There are libraries that you can use in Python that will give you all (or at least nearly all) of the functionality of MATLAB. For example, scipy.integrate.solve_ivp() supports a number of methods for ODE integration that are comparable to what you can get with the various odexxx() functions in MATLAB. So no, you wouldn't have to write your own ODE ...

12

For logging that allows full reproducibility, I highly recommend the Sumatra python package. It nicely links the version control commit number, machine state, and output files to each program run and has a django web interface to interact with the database of run info. The python API makes it very easy to include logging in my scripts.

12

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 eigenvalues and eigenvectors of Hermitian matrices; it's what scipy.eigsh() uses, through a call to the ARPACK library. In exact arithmetic, the Lanczos algorithm ...

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