Hot answers tagged

5

You could define a linear opearator and pass it to the function eigsh. Ideally, your matrices $L$ and $X$ are sparse so you can take advantage of the matrix-vector product. In your case, you would have something like the following. import numpy as np from scipy.sparse.linalg import LinearOperator, eigsh def mv(v): a = 2.3 return L@(L@v) + a*X@(X.T@v)...


5

There is no simple fix. For an ill-conditioned matrix $A$, the harm (loss of precision) is already done the moment you wrote those numbers in a numpy array, because that tiny $10^{-16}$ perturbation from the exact non-representable values is already harmful. You could increase your working precision; but at that point the question is if your matrix entries $...


4

Ideally you would have the matrix already in a sparse matrix data structure. But for this example we can do the conversion via In [31]: M = scipy.sparse.coo_matrix(np.array([[0,0,0,0],[1,0,0,0],[1,1,0,0],[1,1,1,0]])) Then you can do In [32]: Mcsr = M.tocsr() In [33]: np.split(Mcsr.indices, Mcsr.indptr) Out[33]: [array([], dtype=int32), array([...


4

If that is the case, what computational method can be used to confirm that the linked configuration is, in fact, a good spherical code? You can't, with the given data. With only 16 significant digits available, you will never know if that quantity is exactly 1 or not. Or even if you were given Float128s, or numbers with 2000 significant digits, for that ...


3

mats_row=np.array([[A1,A2,A3,...,A_n]]) #Create array of matrices with shape (1,M,N,N) mats_column=np.transpose(mats_row,(1,0,2,3)) #Make a copy with shape (M,1,N,N) block=.5*(mats_row+mats_column) #Add, broadcasting to a (M,M,N,N) array This works by utilizing Numpy's broadcasting capabilities. The basic idea is similar to if you added a matrix where each ...


2

Let's reconstruct this from first principles: Defining the ODE system In the method-of-lines discretization you solve an ODE system $\dot U=F(U)$, $U=(U_0,U_1,...,U_{M+1})$, $U_k(t)=u(x_k,t)$, and similarly $F=(F_0,F_1...,F_{M+1})$. Because of the boundary conditions $$ u(x, 0) = 40 · x^2 · (1 - x) / 3 \\ u(0, t) = u(1, t) = 0 $$ $U_0=U_{M+1}=0$ and ...


2

We consider two thin infinite straight parallel wires at distance $2 a$ apart and carrying equal currents $I$ in opposite directions. We need to find an approximate description of the system in terms of multipole expansion moments. Let's first make a few preliminary comments. The magnetic field is described by the magnetic vector ...


2

Well, you "overwrite" the table at every loop... As a result when the for loop is terminated the table has only the t,X,Y columns for theta=75. Try this: columns = [] for theta in range(30, 90, 15): tmax = v * np.sin(theta*np.pi/180)/g #Time to reach maximum height t = np.linspace(0, tmax,600) X = x(theta) #Maximum horizontal ...


2

This is the code I came up with. I still dont feel confident about the results, mostly because I didn't find the usual amount of errors yet. import numpy as np from scipy.special import kn, iv import matplotlib.pyplot as plt global size, mu_0, p, a, q, current, order mu_0 = 4 * np.pi * 1e-7 size = 10 order = 5 # order of bessel functions considered p = 1.8 ...


2

The program solves the problem Place n points on a sphere in d dimensions so as to maximize the minimal distance (or equivalently the minimal angle) between them. As an aside, let's take $n = 3$ and $d = 2$. Then the solution is three points on a circle 120° apart. But that means we have coordinates $(0,0)$,$(cos(2\pi/3),sin(2\pi/3))$ and $(cos(2\pi/3),-...


1

xprime = lambda x,k,alpha,rho,beta,t: is in the wrong argument order, first are the fixed arguments, then the parameters. xprime = lambda x,t,k,alpha,rho,beta: I'm not sure where the division-by-zero error originates, it could be a consequence of the solution diverging, doing strange things in the step-size controller.


1

This post is better suited for Stackoverflow I think. Anyway you can simply solve your problem by changing psi[:,ig]to psi[:,ig:ig+1]. Then the left-hand side is truly a nx1 matrix, not just a vector of size n. Or you can remove the line nodes = np.array([nodes]).T which is useless here, and causes Numpy to transform the array "nodes" (size n), ...


1

The exact solution is $$x(t)=\frac{3}{1-3t}.$$ This means that the domain of the solution is $(-\infty,\frac13)$. odeint with the integration interval $[0,5]$ stalls at this singularity, as it can, correctly, not move beyond this point.


1

In your Ifft function, you are taking the absolute value of the real part. Just take the real part without doing the absolute value.


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