All Questions
11,162
questions
0
votes
0
answers
18
views
Constructing Mesh In FreeFem++
I am new to building 3D meshes in FreeFem++. Is it possible to build the mesh below? The domain is 1+1/8 in length, 1 height, and 1 width. There is a 1/4 depth and 1/8 wide channel that cuts though it ...
1
vote
0
answers
12
views
Numerical calculation of Lyapunov exponents using SciPy's built-in solve_ivp
I have previously successfully implemented the QR decomposition method in MATLAB to calculate Lyapunov exponents for Lorenz equations. See here.
This method integrates the stacked system, i.e. the ...
3
votes
1
answer
357
views
In linear programming, is there a way to constrain two variables to not have opposite sign
Say I have two sets of variables $x$ and $y$ of equal size. $x$'s have a lower bound $x_{min}<0$, and $y$'s have a lower bound $0$.
Is there a linear way to constrain that $x_i\geq0$ if the ...
0
votes
0
answers
17
views
Discretization of 2D advection equation with non-constant speed
Suppose I have a 2D advection equation $$\frac{\partial \rho}{\partial t}=-\nabla\cdot(\vec{w}\rho)$$ with $\vec{w}=(u,v)$ non-constant and having zero divergence. I want to numerically solve this ...
0
votes
0
answers
21
views
How to use a custom OdeSolver in Scipy's solve_ivp
In Scipy's solve_ivp documentation, we see the method argument can be either a string or a user-defiend OdeSolver inherited from ...
1
vote
1
answer
74
views
accuracy problem for a null space calculation on a sparse rectangular matrix
I have been using the QR-based approach on this link to find the null space of rectangular matrices, and possibly sparse matrices, that emerge as a result of some coupling conditions of different ...
-1
votes
1
answer
39
views
Gmsh problems in Google Colab, visualize mesh
I'm trying to implement mesh in Google colab from gmesh tutorials.
I have an error: Exception: Fltk not available
My code is:
...
1
vote
1
answer
63
views
Numerically computing envelope of Gibbs oscillation
If I numerically compute the envelope of $\sin(\pi t)$ using a Hilbert transform, I obtain exactly what I expect:
If I do the same for $\mathrm{sinc}(t)$, still I obtain an envelope which agrees with ...
2
votes
1
answer
146
views
Where am I making a mistake in solving the heat equation using the spectral method (Chebyshev's differentiation matrix)?
I would like to numerically solve the following heat equation problem:
$$ u_t = \Bigg(2{a \over l}\Bigg)^2 u_{xx} \tag 1$$
$$ x \in [ -1, 1 ] \tag 2$$
$$ u(x, 0) = 0 \tag 3$$
$$ u(1, t) = A \sin \Bigg(...
0
votes
0
answers
49
views
ENO-Runge-Kutta discretization
One beginner's question about discretization of a Hamilton-Jacobi equation(non-linear)
$$
u_t = H(u_x)
$$
$u_x$ is discreated with 2nd order ENO-FD
1st order: $D_1^{\pm}u = \pm [u_{x\pm1} - u_x ] / \...
0
votes
0
answers
24
views
Help with inferring Network topology from Spectral templates
I am trying to use matlab and YALMIP to solve a graph learning problem of recovering eigenvalues from the eigenvectors of the covariance of sampled graph signal data. This is to implement the ...
0
votes
0
answers
24
views
Recommendations for some new books about computational contact mechanics in solid mechanics
I want to simulate some frictional contact problems, but I'm not familiar with this field, could you please recommend some new books as introductions?thank you
1
vote
0
answers
50
views
Overlap matrix and its inverse matrix
Now, we consider a non-orthonormal basis:
$$\mathcal{S}_N=\{|\alpha\rangle,a^\dagger|\alpha\rangle,a^{\dagger 2}|\alpha\rangle,\ldots,a^{\dagger N}|\alpha\rangle\},$$
where $|\alpha\rangle$ is the ...
0
votes
0
answers
43
views
3D Quadrature schemes with points on boundary
In one dimension there are two types of quadrature schemes.
asymmetric rules like Newton-Cotes like formulas (Trapezodi, Simpson), and Clenshaw-Curtis place sampling points on boundary of the ...
0
votes
0
answers
35
views
Methods for delaying the "break" in non-linear least squares optimisation when the step size gets too small?
I am using the Levenberg-Marquardt method for calibration purposes. Typically, the RMSE of my calibration looks like:
I want to break the algorithm when the algorithm step-updates start to slow down, ...
0
votes
1
answer
75
views
Are there any established direct eigensolvers for sparse hermitian matrices?
I have experience with LAPACK (direct solvers) and ARPACK (sparse iterative solvers), but are there any sparse direct solvers? I am concerned more with preserving space than with fast solutions. ...
0
votes
0
answers
38
views
Which numerical method can I use to solve this system of hyperbolic PDEs?
Backround
The mathematical model I am trying to numerically solve models wave propagation inside a cylinder with specific material properties suited for dynamic loading. The cylinder's upper base is ...
1
vote
1
answer
45
views
Finding the Vector $v$ for a Given Householder Matrix Transformation of Non-Collinear Vectors $a$ and $b$
Consider a vector $v$ in $\mathbb{R}^{n\times1}$. The Householder matrix is defined as follows:
$$H(v)=I-\dfrac{2vv^T}{v^Tv}.$$
It can be demonstrated that $H(v)$ is symmetric and orthonormal. The ...
2
votes
1
answer
76
views
Solving for $X$ in $\sum_{a,b} b a^T b^T X a = Y$
Suppose I have $k$ pairs of $(a,b)$ where $a$ and $b$ are vectors in $\mathbb{R}^d$, $Y$ is $d\times d$ and I need least squares solution for $X$ in the following
$$\sum_{(a,b)}^k b a^T (b^T X a) = Y$...
0
votes
0
answers
49
views
eigenvalues of inhomogeneous Helmholtz equation violate superposition using FEM
I am trying to solve the in-homogeneous Helmholtz equation with damping and forcing using finite element method (FEM) with FEniCSx. The equation is;
$c^2\nabla^2p - c^2\omega i d \nabla^2p + \omega^2 ...
1
vote
1
answer
139
views
Solving PDEs using FEM using cubic Hermite polynomials
everyone.
I am a beginner in Numerical mathematics, I have some idea of how to use Galerkin method to solve PDEs numerically, but so far I had no luck finding an example of how to solve a simple PDE ...
0
votes
1
answer
96
views
Step size constraint in Euler backward
I am dealing with an assignment in MATLAB. It has to do with 'self-driving' cars which are driving in-front/behind eachother. Assuming M cars on a single-lane road, each car adjusts its speed based on ...
2
votes
1
answer
128
views
Large set of nonlinear equations in Sympy
I have a set of 6 nonlinear equations, and using Sympy I find the values of the 6 unknowns. This works perfectly and it directly gives the exact solution, using sympy.solve to be specific. Now I ...
5
votes
2
answers
438
views
Number of function calls and jacobian calls in scipy.root
Just as an exercise, I am numerically solving the following system of equations:
$$
\begin{equation}
\begin{cases}
x^2 + y^2 = 32 \\
3x + 7y = 15
\end{cases}
\end{equation}
$$
...
2
votes
1
answer
89
views
references for optimization in the context of parameter identification with finite elements
i am performing parameter identification for a non-linear partial differential equation (elasticity) that I solve with finite elements.
My optimization problem is a non-linear least squares data-...
3
votes
1
answer
244
views
Finite difference problem
I have a problem to resolve with the Finite Difference method in $[a,b]$:
$$-\frac{d}{dx}(\alpha(x)\frac{du}{dx})= g(x),$$
with $\alpha(x) \in L^{\infty}$ continuous in $]a,c[$ and $]c,b[$ and ...
3
votes
0
answers
70
views
How is the Alternating Schwarz Method used as a Preconditioner to a Krylov Method?
I am reading "Domain Decomposition: Parallel Multilevel Methods for Elliptic Partial Differential Equations" (Smith 1996), and I am confused as to how the below Alternating Schwarz algorithm ...
1
vote
1
answer
86
views
Symmetry axis boundary condition
I was wondering about the symmetry axis boundary condition in commercial CFD solvers such as ANSYS Fluent.
If the problem is the flow through a round pipe or out of a round nozzle, it is natural to ...
2
votes
1
answer
156
views
Solving non-linear partial differential equation numerically: $u_{xx}+u_{yy}=\mathrm{e}^{u}$
To start with, I need to solve this partial equation numerically, but I do not know how to do that. If I try a finite difference method, I face a problem that $u_{i,j}$ is also located in exponential, ...
2
votes
1
answer
113
views
How to formulate a convex expression to minimize the difference between Frobenius norm of a positive semidefinite matrix and a positive value
So what I am trying to do is to minimize the distance between the Frobenius norm of a PSD matrix and a real positive value, which can be formulated as
$$\min \left|\|\textbf{P}\|_F - J\right|^2$$
...
2
votes
0
answers
42
views
How to use a preconditioner estimated from a subset of data?
Suppose I'm solving $Ax=b$ using row-action method like Kaczmarz for $m\times n$ matrix A with $m\approx \infty$ and have $H_k=\frac{1}{k}A_k^T A_k$ which is an estimate of the Hessian obtained from ...
3
votes
1
answer
185
views
Solving underdetermined Lyapunov equation?
I'm solving the following for $X$ with $A,B$ singular positive semidefinite matrices.
$$AX + XA = B$$
Because $A$, $B$ are singular, standard Lyapunov solver fails
However, if I heuristically skip ...
0
votes
1
answer
44
views
draw a log-log plot of MSD (mean square displacement) versus `t` of a movement of the polymer chain
Cross-posted on MMSE (Matter Modeling Stack Exchange).
The following are the movements of the center of mass of a polymer chain over time in a monte carlo simulation.
...
1
vote
2
answers
96
views
Implementing matrix term version of Gauss-seidel
I am trying to implement the below description from Ch. 11 of Heath's "Scientific Computing An Introductory Survey" of the Gauss-Seidel iterative method for solving a system of linear ...
0
votes
1
answer
67
views
How to plot the power spectrum
I have an array of data whose columns are solution vectors to a system of ODEs at a specific time. I want to plot the power spectrum of a solution at a specific time, but when I attempt this I get ...
5
votes
0
answers
114
views
Single precision vs double precision conjugate gradients
I tested my conjugate gradients implementation with float and double precision and contrary to my guess the double code was twice faster than the single precision code. The reason is that I need many ...
0
votes
0
answers
12
views
Simulating a dataset from model output when model includes multiple binary deviation-coded variables
I am trying to simulate data using parameters from a glmer() model output. The model, which comes from a published paper, is as follows: DV ~ 1 + group* sex *verb type + trial number + (1 |participant)...
2
votes
0
answers
57
views
Numerical integration in Fourier space over 3D grid
I am attempting to implement a model outlined in this paper:
General magnetostatic shape–shape interactions
Background
This model allows the calculation of magnetostatic interaction energies between ...
4
votes
0
answers
120
views
FVM for non-regular domain with triangular mesh
Setup
The 1D convection-diffusion equation is given by:
\begin{equation}\tag{1}
\frac{\partial u}{\partial t} + v \frac{\partial u}{\partial x} - \mu \frac{\partial^2 u}{\partial x^2} = 0,
\end{...
0
votes
0
answers
81
views
A confusion about the bubble function in lumped mass FEM
I am studying knowledge related to lumped mass finite elements. As is well known, lumped mass finite element methods higher than 2nd order on simplex mesh require the construction of new function ...
0
votes
0
answers
46
views
Is a sort of "z-drift" the result of numerical precision errors in FDM?
Upon solving the 2D wave equation with Neumann boundary conditions $u_x = u_y = 0$ on a rectangular $10 \times 10 \times 10$ grid, I noticed something odd - $u$ seemed to shift upwards with time. This ...
1
vote
1
answer
83
views
Adding a diffusion term to the MUSCL - Kurganov and Tadmor central scheme
Im currently using a MUSCL scheme with a rusanov flux and Van Leer limiter to simulate the 2d euler equations:
$$ \frac{\partial \rho}{\partial t} + \frac{\partial \rho v_x}{\partial x} + \frac{\...
1
vote
0
answers
42
views
Solution to the Liouville-Gibbs equation
What would be the approach to numerically solve for $\rho(x,t)$ the following equation with some initial conditions
$$\frac{\partial\rho}{\partial t}
+\sum_{i=1}^n\left(\frac{\partial(\rho g_i)}{\...
0
votes
1
answer
53
views
Local truncation error of given implicit 1-step scheme
I'm given the 1-step implicit scheme $$y_{n+1} = y_n + \frac{h}{6}[4f(t_n, y_n) + 2f(t_{n+1}, y_{n+1}) + hf'(t_n, y_n)],$$
where $y'(t) = f(t, y)$, and I'm seeking the scheme's local truncation error. ...
1
vote
0
answers
35
views
Implementation of operator splitting method for Wigner equation
I am dealing with the integro-differential equation for Wigner function,
$$\frac{\partial f}{\partial t}+p\frac{\partial f}{\partial x}+\\+\frac{1}{\chi}\left\{\int_{-\pi}^{+\pi}dy\,\int_{-\infty}^{+\...
3
votes
1
answer
175
views
Stability of Euler forward method
I am trying to solve a linear system of ODEs of the form:
$$ \frac{du}{dt} = A u, \quad u(0)=k$$
where $A$ is a 2x2 matrix and $u(t)$ is a 2x1 column vector. I want to solve this numerically, using ...
0
votes
0
answers
26
views
Doubt of wavelets about the return of comand plot() of package Wavethresh in R language
I take a plot of father wavelet with this code
library(wavethresh)
y <- c(1,1,7,9,2,8,8,6)
ywd <- wd(y, filter.number=1, family='DaubExPhase')
plot(ywd)
But ...
1
vote
1
answer
116
views
Calculating Hypergeometric1F1 for large arguments
Cross posted on StackOverflow
I am trying to use the gsl library for calculating 1F1. I have some C code.
The following works and matches Mathematica's results for ...
0
votes
0
answers
16
views
How do you determine the Mott-Insulator to Superfluid transition in the Bose Hubbard System
I am doing some simulations on various systems expressed in 2nd quantization and one of the points of interest of mine was Phase transitions in the Bose-Hubbard model
$$
H = \sum_{k} \{ t_k(b^\dagger_{...
1
vote
1
answer
82
views
derivative matrix and the Dirac delta distribution
For a project I'm working on, I was working with the following equation
$$
w(x) = \int k(x,y)v(y)dy
$$
I noticed that if I choose
$$
k(x,y) = -\delta'(x-y)
$$
Then we probably get (I haven't touched ...