All Questions
9,348
questions
0
votes
0answers
29 views
“This DAE appears to be of index greater than 1” daeic12 (line76) error code
Hi I am trying to solve a set of pde converted into ODE and DAE using central finite difference method. I have used the MATLAB 'solve' command to determine the coefficients of fictitious nodes for ...
3
votes
1answer
55 views
Understand the need for Welford's online algorithm
I am puzzled by the Wikipedia entry discussing many online algorithms for computing the sample variance, including the Welford's online algorithm.
In particular, the sample variance $s_n^2$ can be ...
-1
votes
0answers
32 views
Calling scipy.integrate.odeint, solve_ivp, or any possible python solver in matlab
Do you know any way to call scipy.integrate.odeint, solve_ivp, or any python solver directly within a Matlab code? I have tried ...
7
votes
4answers
18k views
Plot integral function with scipy and matplotlib
I want to plot a numerical integral function of some function $f$ using scipy and matplotlib. How can I do this?
I tried the ...
3
votes
1answer
649 views
Interpolation using compactly supported radial basis function
I have been struggling for two days with the following problem.
I would like to do a $d$-dimensional interpolation over some data. I tried first
to use polyharmonic splines, but when the size of data ...
4
votes
1answer
55 views
Backwards Difference Implicit Method for Nonlinear Parabolic PDE in Python
Original Stack Overflow Question: https://stackoverflow.com/questions/65683788/indexerror-index-31-is-out-of-bounds-for-axis-1-with-size-31?noredirect=1#comment116218335_65683788
PDE: u_t = u_xx + u(...
0
votes
1answer
61 views
Technique to find the CFL condition using the Galerkin method in space and finite-difference in time?
I am using the Galerkin method (Discontinuous to be precise) to discretize in space the scalar linear wave equation and the explicit second order centered finite difference scheme to discretize in ...
14
votes
3answers
467 views
What guidelines should I follow for simulation software projects?
I am not sure whether this question belongs here, but I would like to give it a try and benefit from the experience of the people at scicomp.SE.
From my experience, the software quality in ...
7
votes
1answer
561 views
full rank update to cholesky decomposition for multivariate normal distribution
This question is a specialization of full rank update to cholesky decomposition, to which I hope to get a more positive answer.
When calculating the minus log of the multivariate normal distribution, ...
1
vote
0answers
36 views
Integrating a wavelike equation with absorbing boundary conditions
I am trying to numerically solve the following equation:
$\frac{\partial^{2} \phi}{\partial t^{2}}-\frac{\partial^{2} \phi}{\partial x^{2}}+V(x) \phi(x, t)=0$ On some domain, with:
$\phi(x, 0) = I(x)$
...
1
vote
1answer
61 views
Gauss Integration of $\sqrt(x)$
I want to construct a gauss integration for the weight function $w(x) = x^{1/2}$ for
$$\int_{0}^{1}x^{1/2}f(x)dx = a_{1}f(x_{1})+a_{2}f(x_{2})$$
Solving
\begin{align*}
a_{1}+a_{2} =& \int_{0}^{...
1
vote
0answers
50 views
3D Matrix (Tensor) Operating On a 1D Vector
Say you have a tensor $T$ and the components are represented by a 3 by 3 by 3 matrix. And you want to use that tensor to map a vector $u$ into a new vector $s$, both of which are 3 by 1 column vectors....
5
votes
1answer
106 views
Solving PDEs in parallel
I have read different approaches on how to solve pdes in parallel which are discretized using finite element method. For example:
Non-overlapping domain decomposition approach as mentioned in https://...
2
votes
0answers
127 views
Parametric nonlinear programming
I believe, I have a parametric nonlinear optimization problem.
The non-convex constraints depend on some parameters, and I seek a solution that satisfies these constraints for all parameters in a ...
7
votes
1answer
135 views
Which finite difference better approximates $uu'$?
I want to approximate $uu'$ with a finite difference. On the one hand, it seems to be
$$(uu')_i=u_i\frac{u_{i+1}-u_{i-1}}{2\Delta t}=\frac{u_iu_{i+1}-u_iu_{i-1}}{2\Delta t}$$
On the other hand,
$$(uu')...
-1
votes
0answers
21 views
Help in agglomerating a mesh with Metis?
I have got a mesh of triangles (for a finite element code) that I'd like to partition into polygons, much better if convex. I tried metis and it is very complicated but powerful if one knows about ...
0
votes
0answers
25 views
Split classroom sizes in half linear integer programming
I work for a school who is looking to split students into two groups for COVID related reasons. Ultimately we need to split the students by class section and try to keep families to be in the same ...
0
votes
1answer
85 views
R function or package for carrying out maximum likelihood techniques in random effect models
I am applying optim() function in R to obtain maximum likelihood estimates of the fixed effects and random effects in a model with bivariate random effects. The ...
1
vote
1answer
76 views
Is there any way to have a better guess for initial condition of an ODE coupled to CFD as a boundary condition?
I'm doing CFD simulations for blood flow in unstructured grids. My boundary condition at the outlets is called three-element Windkessel which basically calculates the pressure by solving this ODE:
$$...
2
votes
1answer
83 views
Efficiently compute a projection matrix from Householders reflectors
Let $A \in \mathbb{R}^{m \times n}$ where $m \geq n$.
Let $B$ and $\tau$ be the result of applying LAPACK's dgeqrfp method (R on the upper right triangle, and the ...
5
votes
3answers
888 views
implicit vs. explicit domain decomposition methods
I've been working on a finite element code on unstructured methods, which I've parallelized using the Schur complement method. Here's a summary of how I did it:
Assign each triangle of the mesh to a ...
3
votes
1answer
285 views
Givens rotation vs 2x2 Householder reflection
The usual story of Givens rotations vs Householder reflections is that Householder reflections are better if you want to map a long vector to $e_1$, while Givens is better if you want to map a 2-...
0
votes
1answer
45 views
Error too large in leapfrog method for solving the wave equation of a vibrating string
I have been trying to figure out what I did wrong for the last two days. I do not know if I actually did something wrong or if the error is supposed to be this large in usual leapfrog problems. I ...
1
vote
2answers
75 views
Efficient ODE steppers with query of $f$ and $\nabla f$ is efficient
Assume we have an IVP $y'(t) = f(t,y)$, and that $\partial_t f$ and $\nabla f$ are cheap to compute.
Assume further that more derivatives are not cheap to compute, or inaccessible for some reason, ...
0
votes
0answers
52 views
Interpolation of a Concave Mesh
I would like to have an algorithm that interpolates the values attached to nodes in a concave mesh mesh.
To be me precise, assume we have a point cloud P (e.g. in 3 dimensions) and a list of edges E ...
6
votes
2answers
1k views
How to properly calculate CPU and GPU FLOPS performance?
Problem
I'm trying to calculate CPU / GPU FLOPS performance but I'm not sure if I'm doing it correctly.
Let's say we have:
A Kaby Lake CPU (clock: 2.8 GHz, cores: 4, threads: 8)
A Pascal GPU (clock: ...
17
votes
2answers
1k views
Best practice for storing hierarchical simulation data
TL,DR
What is the accepted best practice in scientific computing circles for storing large quantities of hierarchically structured data? For example, SQL does not play nicely with large sparse ...
1
vote
1answer
37 views
Finite element analysis software for acoustic and electrostatic
I need to do a simulation for my thesis project involving some piezoelectric nanoparticles in a fluid beamed with ultrasounds. I'm looking for a software for such simulation and for now it seems me ...
0
votes
0answers
13 views
How to define "Characteristic length in Ansys CFX Power Syntax script?
I am trying to evaluate Reynolds number in Ansys CFX 2020 R2 using CCL/Perl script on different domains in a radial turbine. But for that, I need to define the Characteristic length. So I define the ...
0
votes
1answer
45 views
How can I compute the gradient of the noiseless image given by the compass operator?
I need to generate an image with 11x11 pixels having in the center of the image a square of 5x5 pixels, with the gray level of the background 0 and the gray level of the square is 50. I need to ...
0
votes
0answers
9 views
Different values of Gaussian Integration in Pracma library in R
I want to construct a gauss integration for the weight function:
$$w(x) = x^{1/2}$$ $\in(0,1)$ of the form:
$$\int_{0}^{1}x^{1/2}f(x)dx = w_{1}f(x_{1})+w_{2}f(x_{2})$$
\begin{align*}
w_{1}+w_{2} =...
0
votes
1answer
356 views
Imposing periodic boundary condition for linear advection equation - Node problem
I've spent the whole day trying to figure out what is the correct way to impose (and implement) periodic boundary conditions $u(0,t)=u(1,t)$ for all $t>0$ for the simple advection equation $u_t + v ...
2
votes
1answer
65 views
How to include negative number in the log-sum-exp?
I want to know summation of some small numbers, such as {e^-1000, -e^1001, e^1002...}
If all numbers are positive, I can use log-sum-exp algorithm. But unfortunately, negative numbers are also ...
0
votes
1answer
66 views
Smoothed Particle Hydrodynamics: Weird clustering of particles. Is that normal?
I implemented a rather simple SPH simulation using a cubic-spline-kernel and a simple non-iterative pressure solver as described in this PDF in equation 9. I followed algorithm 1 of that paper (...
0
votes
1answer
58 views
How to compute gradient of each node in finite element method?
I have a problem with computing gradient of each node in finite element method. I can get the value of each node. But how can I get the gradient?
I know $u = \sum u_i \phi_i$ where $\phi_i$ are the ...
1
vote
0answers
34 views
N-body correct scaling
I realized an usual way to scale an N-body problem for an N-body simulation is by choosing units such that gravitational constant $G = 1$, but I'm probably doing it the wrong way. Suppose I simply ...
0
votes
0answers
26 views
Largest triangle that contains a point
Given the location of $n$ points on a 2D plane ($P_1, P_2, \ldots, P_n$); and the location of a special point $X$.
Find three points $P_i,P_j,P_k$ ($i \neq j \neq k$) such that point $X$ is inside the ...
-1
votes
0answers
31 views
Gmsh undesired behaviour 3D mesh
I would like to generate a 3D hexahedral mesh starting from a STL file using Gmsh 4.8.0.
I followed these steps:
Opened a STL file: file->open
Added a volume: Modules -> Geometry -> ...
4
votes
1answer
110 views
Floating point and global error in Euler Method
Inspired by this answer, I tried to understand when floating point errors come into visibility and to check it also comparing the plot of the numerical solution with Explicit Euler with the analytical ...
-1
votes
0answers
26 views
Defacto Barnes Hut / N-Body MPI C implementation
I'm looking for a relatively straightforward but robust MPI-based parallel n-body application implementation in C/C++ -- ideally one that's off the shelf and has been used in other papers before. A ...
-1
votes
0answers
26 views
Is it ok to get the parameter of a random variable, from a random variable itself?
I want to generate the parameter, for a Poisson distribution, from a lognormal distribution. Has there been literature on if this is proper in probability & statistics, or are there nuances to ...
0
votes
1answer
90 views
Numerical solution of ill-conditioned differential equation
I want to solve the following Cauchy problem
\begin{equation}
y' = y^2 + \frac{t^4 - 6t^3 + 12t^2 - 14t + 9}{(1+t)^2}
\end{equation}
with initial condition: $y(0) = 2$ for $t \in [0,1.6]$ using a 3 ...
-1
votes
0answers
21 views
How could I find out how many FLOPs my calculator is capable of?
Looking at previous S.E answers, we have a formula for this:
The theoretical peak FLOP/s is given by
$$ \text{Number of Cores} * \text{Average frequency} * \text{Operations per cycle} $$
After ...
0
votes
0answers
27 views
Parameter sampling from a 3D isosurface in R
Before we start, a small disclaimer: I am not a computer scientist, and the field of isosurfaces is new to me, so hopefully, the question is phrased clearly:) Otherwise, let me know, and I will (try ...
0
votes
2answers
69 views
Problem of half-planes intersection
Consider the half-planes $\{x \leqslant 2\}$ and $\{x+y \leqslant 3\}$. These two half-planes are coded with the R package 'rcdd' as follows:
...
1
vote
0answers
48 views
Open-source Hexahedral meshing tools
I'm looking for an open-source automatic meshing tool or library able to generate full or dominant hexahedral mesh starting from a provided arbitrary model i.e. a surface triangular mesh stored in STL ...
3
votes
1answer
56 views
How does the diffusion of a finite volume method with a WENO scheme compare with that of spectral methods?
I know that, in general, finite volume (FV) methods are more (numerically) diffusive than spectral methods. However, I can't find any information on how the advection scheme changes that.
For example, ...
0
votes
0answers
21 views
How to limit the output value of an Linear Congruential Generator, but still maintain an extremely long period?
Let us suppose that I have a full cycle LCG, following the standard equation
Hence, as it has a full cycle, it will have a period of m, and will output values ranging from 0 to m. However, I want to ...
3
votes
2answers
132 views
Using the BDF and RK4 methods to solve this coupled system of ODEs in C++
I'm trying to solve a system of ODEs using the BDF order 4 method. I find the first 3 points using RK4, then for the implicit part of the BDF, I use Newton-Raphson iteration. Unfortunately my solution ...
3
votes
1answer
66 views
Is there a source/cookbook of equations that approximate geometric shapes?
I'm numerically modelling flows around various geometric 2D shapes. Is there a good source/cookbook of equations that approximate these? Some examples are
Rectangle:
$(x-a)^n+(y-b)^n < r^n$ where ...