9

This is a linear integral equation (because the right hand side is linear in $h$) and presumably the problem is ill-posed (for one, because of the multiplication by $\sin(\theta)$ but also because the "integral kernel" is smooth). Methods to solve this (approximately) are general Galerkin methods or collocation methods. For Galerkin methods you take a ...


6

Seems that you have a duplicate eigenvalue. Thus, you have two eigenpairs $(\lambda_1, x_1)$ and $(\lambda_2, x_2)$ where $\lambda_1 = \lambda_2$. Denote $\lambda = \lambda_1 = \lambda_2$. Let $\alpha$ and $\beta$ be arbitrary complex numbers. Then $$A (\alpha x_1 + \beta x_2) = \alpha A x_1 + \beta A x_2 = \alpha \lambda x_1 + \beta \lambda x_2 = \lambda (\...


6

Preamble: "integers" -> Integers in your problem statement, Mathematica uses CamelType. The method depends on your problem, naturally. In your case it appears to be the Contejean–Devie method. A [buried] excerpt from the Mathematica documentation: For polynomial systems, Reduce uses cylindrical algebraic decomposition for real domains and Gröbner basis ...


6

Actuallly you can scale the variable of the problem into a different interval. Transform $ x \in [0,a]$ where ($a=10^6$) to $\xi \in [0,1]$ by the formula $ \xi = \frac{x}{a} $. Thus, using the chain rule we have $\frac{d y }{ d x } = \frac{d y}{ d \xi } \frac{d \xi }{ d x } = \frac{ 1}{ a } \frac{d y}{ d \xi } $. Hence, replacing $\frac{d y }{ d x }$ by $\...


5

MATLAB always uses the LAPACK libraries to calculate eigenvectors which works on double precision, floating point numbers - MATLAB's default data type. Mathematica's method depends on its input type. For example, when you do TestMatrix = {{1, 2, 3}, {3, 1, 2}, {2, 3, 1}} Eigenvectors[TestMatrix] You'll get an exact answer involving Sqrt[3] and so on. ...


5

I want to mention another idea in case it helps. The truncation error of replacing $\int_{-\infty}^{\infty}$ with $\int_{-L}^L$ is on the order of $$\begin{aligned} \int_L^\infty \cos(tx)g(x)\,\mathrm{d}x &= \frac{1}{t}\sin(tx)g(x)\big|_{L}^{\infty} - \int_L^\infty\frac{\sin tx}{t}g'(x)\,\mathrm{d}x \\&= -\frac1L g(L)\sin(tL) + \text{asymptotically ...


5

Note: I'm somewhat worried at this point that the integral values Mathematica gives me are bogus. I thought it was working because it gave a sensible-looking result in a short time, but it might be the case that the method it tries to use is buggy or that I did something wrong. So it might be that the code below isn't working at all, I don't know, sorry. ...


5

Julia's macros and string macro are made for doing this kind of thing. You can define a domain-specific language (DSL) within a macro scope and use that as a way to define problems. One example of this is the DifferentialEquations.jl @ode_def macro. You may have to use a string macro (due to parsability into ASTs) but this sort of tool is what allows you to ...


4

@tqviet already says this in a some way, but let me phrase it differently: $10^6$ is not a large end time by itself. It seems like a large number because it has 6 zeros at the end, but it all depends on (i) which units you compute it, (ii) how fast the effects are that you are trying to model. To give examples: If I told you that you had to compute to $10^...


4

First of all, it is impossible to completely rewrite NIntegrate because usually, it evaluates (parts of) the integrand symbolically to check for certain properties. This means, you probably have to re-implement Mathematica completely. Nevertheless, if you are tackling one specific integral, I assume chances are very good that you can get similar numerical ...


3

Are you asking for a plot of $z(d)$ and $w(d)$? Wolfram Alpha might not like the long query, but Mathematica simplifies your number/letter soup just fine: $$w(d)=-\sqrt{3}\frac{ d^2 (x+2 y-1)-2 d \left(x (2 y-1)+y^2\right)+x (4 y-3)-2 y^2+6 y-3}{d^2-2 d y+4 y^2-6 y+3} = \frac{1}{2} \left(\sqrt{3} d-3 \sqrt{(d-1)^2}+\sqrt{3}\right)$$ $$z(d)=\frac{x \left(d^...


3

First, you're solving the Schrödinger equation not in a box, but rather you apply periodic boundary conditions. That is, it's more like solving the TDSE on a ring. But this just as a comment and shouldn't matter here. What's happening here is that you evolve a Gaussian wave packet in time. In a suitable hamonic oscillator potential, it would be stable and ...


3

Since you already know Matlab and Mathematica, I would go with SciPy Lecture notes. It is a course on using Python for Scientific purposes. Gael Varoquaux, Valentin Haenel, Pierre de Buyl, Gert-Ludwig Ingold, Emmanuelle Gouillart, Michael Hartmann, … João Felipe Santos. (2017, October 4). scipy-lectures/scipy-lecture-notes: Release 2017.1 (Version 2017.1). ...


3

The short answer is "It depends". Certainly, you can try to find the Helmholtz decomposition on your sampled data, and find your irrotational and solenoidal components. However, there are certain requirements on your original vector field you started with. In general, it is that the vector field you are trying to decompose has to be sufficiently smooth and ...


2

Seeing what are your boundary conditions, you should use an spectral(Galerkin) method. I am trying to solve sometrhing similar, but in 9 dimentions, And the way to treat the field fading in infinity as a boundary condition is with galerkin methods. Probably using a basis like $F_k(x)=e^{-\alpha x^2}\sin(k x)$. Though you should look more into the probable ...


2

Here is something I threw together, it is definitely incorrect, as the program is throwing errors. In particular, a SolutionVariableNumberError is raised. Likely this is because the equations are written incorrectly. The likely suspect, I'm guessing is in the SourceTerm's. This is a start though for how one might go about solving these equations in FiPy. ...


2

There are two possible issues going on here. Firstly as noted in the answer by knl any linear combination of eigenvectors corresponding to degenerate eigenvalues is also an eigenvector, as shown in that answer. However even for the non-degenerate case eigenvectors may vary between different invocations of a diagonaliser even if the input matrix is the same. ...


2

One the most useful things I found myself when converting from matlab to numpy/scipy was this page NumPy for MATLAB users. I even look at this every once in awhile to remind myself even to this day. I also highly recommend taking a look at the matplotlib gallery as what you probably want to do visually in either MATLAB or Mathematica is hopefully ...


1

Unfortunately I am not a Mathematica whiz. and in general do not encourage the use of Mathematica for extensive numerical simulations like Ising model or the likes. So I will post some answers pertaining to C++ and the general equations. First I am going to assume, that your $E_t$ is actually the spin-exchange energy $-JS_iS_j$ where $S_i$ are the spins and ...


1

It is impossible that finite element or finite difference violate their boundary condition that you imposed. It is not related to the physics of your problem and it's just an issue with your implementation. No matter your ODE is consistent with your boundary conditions or something else. The only thing matters is that you enforce your boundary conditions ...


1

One problem with converting Tex syntax to equations is that it can be ambiguous - the main goal with TeX is typesetting the look, not the actual meaning. If you have many formulas already in TeX, then looking to parse them into a symbolic form is okay. But if starting from scratch, a symbolic mathematics package would be a better choice. Another issue ...


1

You got good advices in comments, but as it is only 1D problem, you can manage it using fully explicit method with very small time steps. I checked your code in Mathematica and there are several misunderstandings there. After correcting them, the code works for me. Similarly, your description in the question contains, I suppose, some mistakes. As I guess ...


1

The classic book in B-Splines is de Boor, C. (1978). A Practical Guide to Splines. New York: Springer-Verlag, but I spent many weeks on it and was not able to make his algorithms work. The algorithm in Phillips, G. M., & Taylor, P. J. (1996). Splines and Other Approximations Chapter 6 of Theory and Applications of Numerical Analysis (2nd ed., pp. 131-...


1

If you want to diagonalize significantly larger matrices on large HPC clusters, I would look at Elemental, which is a state-of-the-art distributed dense linear algebra library. It's a better alternative than other libraries out there (e.g., ScaLAPACK, PLAPACK), and builds upon BLAS and LAPACK implementations.


1

How about the C++ (with python bindings) Gamma library for NMR simulations: The Gamma Library. In the GammaDetailedDescription, the authors mention the various components of the library, which in turn would solve many aspects of your question.


1

Given the apparent dearth of quasi Monte Carlo integrators for C++ (or C), I wrote my own implementation to be used with GSL. It's not especially well-tested, nor does it implement Mathematica's algorithm, but it should be better than nothing.


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