# Numerics

**URL:** https://discourse.julialang.org/c/domain/numerics/20.md?page=2

[Latest](https://discourse.julialang.org/latest.md) · [Categories](https://discourse.julialang.org/categories.md) · [Tags](https://discourse.julialang.org/tags.md)

**Page:** 3

---

## [Why doesn't 1//0 throw an error?](https://discourse.julialang.org/t/why-doesnt-1-0-throw-an-error/131354)

<div class="topic-metadata">

**Author:** [@HMegh](https://discourse.julialang.org/u/HMegh)\
**Replies:** 5\
**Last updated:** [August 5, 2025, 4:21am UTC](https://discourse.julialang.org/t/why-doesnt-1-0-throw-an-error/131354 "2025-08-05T04:21:30Z")

</div>

Hi, I have noticed that 1//0 does not throw any errors, which is a bit unexpected. Even more, inv(0//1) is 1//0. On the other hand, 0//0 throws an argument error (understandably). Is this the intended behavior? My cur…

---

## [Reverse rule in Enzyme for an implicit function](https://discourse.julialang.org/t/reverse-rule-in-enzyme-for-an-implicit-function/131327)

<div class="topic-metadata">

**Author:** [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Replies:** 2\
**Last updated:** [August 4, 2025, 1:19pm UTC](https://discourse.julialang.org/t/reverse-rule-in-enzyme-for-an-implicit-function/131327 "2025-08-04T13:19:08Z")

</div>

Suppose H is the cdf of a distribution, A, B are real numbers, and p is defined by p = H(p A + (1 - p) B) This can be implemented simply in Julia as using Roots, Distributions struct Problem{T,S} H::T A::S …

---

## [Finalizer triggers "task switch not allowed" when used under multi-threaded GC](https://discourse.julialang.org/t/finalizer-triggers-task-switch-not-allowed-when-used-under-multi-threaded-gc/131232)

<div class="topic-metadata">

**Author:** [@Cyunkitlau](https://discourse.julialang.org/u/Cyunkitlau)\
**Replies:** 7\
**Last updated:** [August 4, 2025, 4:39am UTC](https://discourse.julialang.org/t/finalizer-triggers-task-switch-not-allowed-when-used-under-multi-threaded-gc/131232 "2025-08-04T04:39:03Z")

</div>

Hi there: I’m working on parameter identification using ModelingToolkit-based ODE models. I construct a symbolic system using ModelingToolkit.jl and build a custom loss function, which is optimized under multithreaded…

---

## [Solving for transcendental equation](https://discourse.julialang.org/t/solving-for-transcendental-equation/131229)

<div class="topic-metadata">

**Author:** [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Replies:** 15\
**Last updated:** [July 31, 2025, 4:14pm UTC](https://discourse.julialang.org/t/solving-for-transcendental-equation/131229 "2025-07-31T16:14:16Z")

</div>

Hi, I have the function f\_n(x) = e^x T\_n(x), where T\_n are Touchard’s polynomials. The function is bijective from R to R and everything is Ok, just that the inverse is unknown… It has the nice taylor exapnsion : f\_n(x)…

---

## [DifferentiationInterface do-block best practices](https://discourse.julialang.org/t/differentiationinterface-do-block-best-practices/131125)

<div class="topic-metadata">

**Author:** [@afleming](https://discourse.julialang.org/u/afleming)\
**Replies:** 1\
**Last updated:** [July 28, 2025, 3:15pm UTC](https://discourse.julialang.org/t/differentiationinterface-do-block-best-practices/131125 "2025-07-28T15:15:18Z")

</div>

Is it better / recommended to use do-blocks that capture variables or should I use the Context mechanism provided by DifferentiationInterface? Even though this pack = (y0 = 1.0, y1 = 2.0) x1 = \[2.0, 4.0\] jac = Differen…

---

## [Support for adjoint sensitivities with \`SciMLStructures\` in ODE problems](https://discourse.julialang.org/t/support-for-adjoint-sensitivities-with-scimlstructures-in-ode-problems/131049)

<div class="topic-metadata">

**Author:** [@albangossard](https://discourse.julialang.org/u/albangossard)\
**Replies:** 2\
**Last updated:** [July 25, 2025, 6:15pm UTC](https://discourse.julialang.org/t/support-for-adjoint-sensitivities-with-scimlstructures-in-ode-problems/131049 "2025-07-25T18:15:02Z")

</div>

Hi all, In complex geoscientific models, it’s common to have many simulation parameters that need to be passed into the ODE function. At the same time, we often want to compute gradients with respect to some of these pa…

---

## [Numerical artifact in solution of system of ODEs](https://discourse.julialang.org/t/numerical-artifact-in-solution-of-system-of-odes/130941)

<div class="topic-metadata">

**Author:** [@Gako](https://discourse.julialang.org/u/Gako)\
**Replies:** 12\
**Last updated:** [July 24, 2025, 10:09am UTC](https://discourse.julialang.org/t/numerical-artifact-in-solution-of-system-of-odes/130941 "2025-07-24T10:09:46Z")

</div>

I am currently solving the following system of ODEs: \\begin{cases} \\ddot\\phi=-3H\\dot{\\phi}-\\phi\\\\ \\dot{a}=Ha\\\\ \\ddot\\sigma=-H\\dot{\\sigma}-\\frac{m\_{Eff}^2}{a^2}\\sigma\\\\ \\ddot{\\rho}\_0=\\frac{4\\rho\_2}{a^2}-H\\dot{\\rho}\_{0}-…

---

## [Solving Large Sparse SPD Linear Systems with Hybrid Direct-Iterative Methods in Julia](https://discourse.julialang.org/t/solving-large-sparse-spd-linear-systems-with-hybrid-direct-iterative-methods-in-julia/130650)

<div class="topic-metadata">

**Author:** [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Replies:** 1\
**Last updated:** [July 11, 2025, 8:50pm UTC](https://discourse.julialang.org/t/solving-large-sparse-spd-linear-systems-with-hybrid-direct-iterative-methods-in-julia/130650 "2025-07-11T20:50:39Z")

</div>

I’m exploring the idea of solving large, sparse symmetric positive definite (SPD) linear systems using hybrid approaches that combine direct methods (e.g., Cholesky) with iterative methods (e.g., CG with incomplete Chole…

---

## [\`\`\`julia -t\`\`\` Does not respect the number of threads](https://discourse.julialang.org/t/julia-t-does-not-respect-the-number-of-threads/130533)

<div class="topic-metadata">

**Author:** [@Dario-Rosa85](https://discourse.julialang.org/u/Dario-Rosa85)\
**Replies:** 13\
**Last updated:** [July 7, 2025, 3:49pm UTC](https://discourse.julialang.org/t/julia-t-does-not-respect-the-number-of-threads/130533 "2025-07-07T15:49:54Z")

</div>

Hi all, I am running some large-scale multithreading numerical code on a couple of different servers. In one of them (and only in one) I am seeing that the machine does not obey the restriction on the number of threads …

---

## [How to specify the nonlinear solver algorithm for Shooting and MultipleShooting BVP solvers](https://discourse.julialang.org/t/how-to-specify-the-nonlinear-solver-algorithm-for-shooting-and-multipleshooting-bvp-solvers/130392)

<div class="topic-metadata">

**Author:** [@ThermionicsGuy](https://discourse.julialang.org/u/ThermionicsGuy)\
**Replies:** 5\
**Last updated:** [July 3, 2025, 4:40pm UTC](https://discourse.julialang.org/t/how-to-specify-the-nonlinear-solver-algorithm-for-shooting-and-multipleshooting-bvp-solvers/130392 "2025-07-03T16:40:05Z")

</div>

SciML’s documentation suggests that the boundary value problem solvers Shooting() and MultipleShooting() should be able to use any of the NonlinearSolve.jl algorithms, but I haven’t been able to find a way to specify the…

---

## [Stochastic differential equations: how to improve the exactness of the solution?](https://discourse.julialang.org/t/stochastic-differential-equations-how-to-improve-the-exactness-of-the-solution/130379)

<div class="topic-metadata">

**Author:** [@Marco\_Masi](https://discourse.julialang.org/u/Marco_Masi)\
**Replies:** 6\
**Last updated:** [July 3, 2025, 12:27pm UTC](https://discourse.julialang.org/t/stochastic-differential-equations-how-to-improve-the-exactness-of-the-solution/130379 "2025-07-03T12:27:09Z")

</div>

I’m trying to find out how to integrate stochastic differential equations. I would like to see how an ODE that I have already solved behaves when subjected to a noisy signal. I wrote the code below. I don’t know if it is…

---

## [Solve linear system under Hermitian constraint](https://discourse.julialang.org/t/solve-linear-system-under-hermitian-constraint/130158)

<div class="topic-metadata">

**Author:** [@araujoms](https://discourse.julialang.org/u/araujoms)\
**Replies:** 25\
**Last updated:** [June 30, 2025, 1:26pm UTC](https://discourse.julialang.org/t/solve-linear-system-under-hermitian-constraint/130158 "2025-06-30T13:26:05Z")

</div>

Is there a nice way to solve the linear system XA = B under the constraint X = X^\\dagger? Here A and B are complex matrices that are not necessarily square. What I came up with is to solve YA = B and A^\\dagger Y = B^\\da…

---

## [Solving a large sparse ODE using Exponential methods](https://discourse.julialang.org/t/solving-a-large-sparse-ode-using-exponential-methods/125159)

<div class="topic-metadata">

**Author:** [@albertomercurio](https://discourse.julialang.org/u/albertomercurio)\
**Replies:** 8\
**Last updated:** [June 25, 2025, 7:49pm UTC](https://discourse.julialang.org/t/solving-a-large-sparse-ode-using-exponential-methods/125159 "2025-06-25T19:49:29Z")

</div>

Hello, I need to integrate a differential equation of the form \\dot{u} = L(t) u \\, , where u is a vector, L(t) = L\_0 + c(t) \* L\_1 \\, , L\_0 and L\_1 large sparse matrices, and c(t) a time dependent coefficient. I k…

---

## [Solvers for Block-Tridiagonal](https://discourse.julialang.org/t/solvers-for-block-tridiagonal/130013)

<div class="topic-metadata">

**Author:** [@rsenne](https://discourse.julialang.org/u/rsenne)\
**Replies:** 14\
**Last updated:** [June 24, 2025, 6:41pm UTC](https://discourse.julialang.org/t/solvers-for-block-tridiagonal/130013 "2025-06-24T18:41:55Z")

</div>

Block Tridiagonal Linear System Performance: Custom vs Package Implementations I recently submitted an issue to BlockBandedMatrices.jl about this, but thought it might be interesting to the broader community. Problem I …

---

## [Matrix multiplication and vector multiplication give different results!](https://discourse.julialang.org/t/matrix-multiplication-and-vector-multiplication-give-different-results/130177)

<div class="topic-metadata">

**Author:** [@tomtom](https://discourse.julialang.org/u/tomtom)\
**Replies:** 4\
**Last updated:** [June 24, 2025, 5:04pm UTC](https://discourse.julialang.org/t/matrix-multiplication-and-vector-multiplication-give-different-results/130177 "2025-06-24T17:04:11Z")

</div>

I found that (slightly) different results coming from matrix-matrix vs. matrix-vector multiplication! julia\> m1 = randn(200, 300); julia\> m2 = randn(300, 400); julia\> m12 = m1 \* m2; # matrix-matrix multiplicati…

---

## [ExponentialUtilities.jl expv!() allocations](https://discourse.julialang.org/t/exponentialutilities-jl-expv-allocations/129917)

<div class="topic-metadata">

**Author:** [@random\_guy](https://discourse.julialang.org/u/random_guy)\
**Replies:** 11\
**Last updated:** [June 23, 2025, 4:56pm UTC](https://discourse.julialang.org/t/exponentialutilities-jl-expv-allocations/129917 "2025-06-23T16:56:05Z")

</div>

Hi everyone, I’m trying to use Julia’s ExponentialUtilities.jl package to perform matrix-exponential-matrix multiplication. I have a sparse NxN Hamiltonian H, a small time step dt, and an NxN matrix A. The goal is to co…

---

## [Issue when using previous solution for initial guess with MIRK Boundary Value Problem solver](https://discourse.julialang.org/t/issue-when-using-previous-solution-for-initial-guess-with-mirk-boundary-value-problem-solver/129948)

<div class="topic-metadata">

**Author:** [@ThermionicsGuy](https://discourse.julialang.org/u/ThermionicsGuy)\
**Replies:** 3\
**Last updated:** [June 21, 2025, 7:26pm UTC](https://discourse.julialang.org/t/issue-when-using-previous-solution-for-initial-guess-with-mirk-boundary-value-problem-solver/129948 "2025-06-21T19:26:57Z")

</div>

I’ve found an odd issue while using a previous solution as the initial guess in solving a boundary value problem using the SciML MIRK solvers and TwoPointBVProblem (part of the BoundaryValueDiffEq.jl package. According t…

---

## [Repeated high (6–7) dimensional interpolation](https://discourse.julialang.org/t/repeated-high-6-7-dimensional-interpolation/130042)

<div class="topic-metadata">

**Author:** [@JADekker](https://discourse.julialang.org/u/JADekker)\
**Replies:** 3\
**Last updated:** [June 19, 2025, 9:21pm UTC](https://discourse.julialang.org/t/repeated-high-6-7-dimensional-interpolation/130042 "2025-06-19T21:21:51Z")

</div>

Thank you, I’ll check it out! I can sketch the problem: We are interested in approximating a (costly) function f(x\_1, ..., x\_p) by simulating a grid of points (x\_1^{(n)},...,x\_p^{(n)}),~n = 1,...,N, evaluating the funct…

---

## [Discrepancy of ODE Sensitivity Analysis paper results with benchmarks](https://discourse.julialang.org/t/discrepancy-of-ode-sensitivity-analysis-paper-results-with-benchmarks/129853)

<div class="topic-metadata">

**Author:** [@mra](https://discourse.julialang.org/u/mra)\
**Replies:** 31\
**Last updated:** [June 18, 2025, 2:01pm UTC](https://discourse.julialang.org/t/discrepancy-of-ode-sensitivity-analysis-paper-results-with-benchmarks/129853 "2025-06-18T14:01:19Z")

</div>

I have tried to reproduce results from the paper “A Comparison of Automatic Differentiation and Continuous Sensitivity Analysis for Derivatives of Differential Equation Solutions” (https://arxiv.org/pdf/1812.01892, Fig. …

---

## [How to get upper triangular Cholesky factor of sparse matrix?](https://discourse.julialang.org/t/how-to-get-upper-triangular-cholesky-factor-of-sparse-matrix/129979)

<div class="topic-metadata">

**Author:** [@slwu89](https://discourse.julialang.org/u/slwu89)\
**Replies:** 3\
**Last updated:** [June 18, 2025, 1:22pm UTC](https://discourse.julialang.org/t/how-to-get-upper-triangular-cholesky-factor-of-sparse-matrix/129979 "2025-06-18T13:22:34Z")

</div>

Hi all. I read the documentation here Linear Algebra · The Julia Language for what is returned when you call cholesky on a SparseMatrixCSC object but I still can’t figure out how to get back the upper triangular factor o…

---

## [Rank-deficient jacobian error in Newton-Raphson method for unstable limit cycle refinement](https://discourse.julialang.org/t/rank-deficient-jacobian-error-in-newton-raphson-method-for-unstable-limit-cycle-refinement/129975)

<div class="topic-metadata">

**Author:** [@R.Kanakov](https://discourse.julialang.org/u/R.Kanakov)\
**Replies:** 3\
**Last updated:** [June 17, 2025, 10:20pm UTC](https://discourse.julialang.org/t/rank-deficient-jacobian-error-in-newton-raphson-method-for-unstable-limit-cycle-refinement/129975 "2025-06-17T22:20:57Z")

</div>

I am trying to refine the position of an unstable limit cycle within a chaotic attractor using the Newton-Raphson method, but I am running into a confusing error related to a rank-deficient matrix. To avoid degeneracy o…

---

## [Sparse GPU linear solve from documentation fails](https://discourse.julialang.org/t/sparse-gpu-linear-solve-from-documentation-fails/129834)

<div class="topic-metadata">

**Author:** [@Jonas208](https://discourse.julialang.org/u/Jonas208)\
**Replies:** 2\
**Last updated:** [June 17, 2025, 10:57am UTC](https://discourse.julialang.org/t/sparse-gpu-linear-solve-from-documentation-fails/129834 "2025-06-17T10:57:22Z")

</div>

To solve the Poisson equation on the GPU, I would like to use the interface from LinearSolve.jl to do a sparse linear solve. However, the following minimal reproducible example (also provided in the documentation) gives …

---

## [Block matrix linear solve performance](https://discourse.julialang.org/t/block-matrix-linear-solve-performance/129353)

<div class="topic-metadata">

**Author:** [@Karajan](https://discourse.julialang.org/u/Karajan)\
**Replies:** 6\
**Last updated:** [June 15, 2025, 5:41pm UTC](https://discourse.julialang.org/t/block-matrix-linear-solve-performance/129353 "2025-06-15T17:41:51Z")

</div>

Hi there! I have located a bottleneck in my simulations I could use some help with. I think this is mainly a linear algebra problem rather than pure Julia performance, but I know very little about that field so I’m not …

---

## [Integrators with analytic derivatives w/r/t initial conditions (especially the independent variable)](https://discourse.julialang.org/t/integrators-with-analytic-derivatives-w-r-t-initial-conditions-especially-the-independent-variable/129832)

<div class="topic-metadata">

**Author:** [@stevenhughes](https://discourse.julialang.org/u/stevenhughes)\
**Replies:** 2\
**Last updated:** [June 13, 2025, 12:43pm UTC](https://discourse.julialang.org/t/integrators-with-analytic-derivatives-w-r-t-initial-conditions-especially-the-independent-variable/129832 "2025-06-13T12:43:37Z")

</div>

Hello, I am working on an application that requires derivatives odes to compute partials for solvers. The approach I am taking is to prefer analytic partials (for performance, which is key), and fallback to AD or finit…

---

## [Eigen produces non-normalized eigenvectors](https://discourse.julialang.org/t/eigen-produces-non-normalized-eigenvectors/129829)

<div class="topic-metadata">

**Author:** [@nvenkov1](https://discourse.julialang.org/u/nvenkov1)\
**Replies:** 5\
**Last updated:** [June 12, 2025, 11:58am UTC](https://discourse.julialang.org/t/eigen-produces-non-normalized-eigenvectors/129829 "2025-06-12T11:58:48Z")

</div>

I realized that, when using LinearAlgebra.jl’s eigen function to solve a general eigenvalue problem for a matrix pencil (A, B), the eigenvectors returned do not seem to be normalized in any way. The standard in the numer…

---

## [Speeding up combination of Funs](https://discourse.julialang.org/t/speeding-up-combination-of-funs/129400)

<div class="topic-metadata">

**Author:** [@nilsbecker](https://discourse.julialang.org/u/nilsbecker)\
**Replies:** 15\
**Last updated:** [June 11, 2025, 1:28pm UTC](https://discourse.julialang.org/t/speeding-up-combination-of-funs/129400 "2025-06-11T13:28:48Z")

</div>

I’m building up an ApproxFun.Fun by successively multiplying a lot of individual factors along a graph. The individual factors are themselves also Funs. The situation is similar to this fake code factors = \[\] for n in …

---

## [Code fails using BoundaryValueDiffEq Version 5.18.0](https://discourse.julialang.org/t/code-fails-using-boundaryvaluediffeq-version-5-18-0/129793)

<div class="topic-metadata">

**Author:** [@ThermionicsGuy](https://discourse.julialang.org/u/ThermionicsGuy)\
**Replies:** 5\
**Last updated:** [June 11, 2025, 3:45am UTC](https://discourse.julialang.org/t/code-fails-using-boundaryvaluediffeq-version-5-18-0/129793 "2025-06-11T03:45:57Z")

</div>

I had a code working using the MIRK4 boundary value problem solver from SciML. Yesterday I updated BoundaryValueDiffEq to Version 5.18.0 and now the same code errors out. I’m using Julia 1.11.5. Here is a code fragment, …

---

## [Quick & dirty online interpolation in Rn](https://discourse.julialang.org/t/quick-dirty-online-interpolation-in-rn/129698)

<div class="topic-metadata">

**Author:** [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Replies:** 2\
**Last updated:** [June 6, 2025, 12:11pm UTC](https://discourse.julialang.org/t/quick-dirty-online-interpolation-in-rn/129698 "2025-06-06T12:11:55Z")

</div>

I am estimating the model which requires mapping a set of parameters x to an implicit solution y = f(x), defined by some F(x, y) = 0. I am solving this problem about 40000 times in a typical run. The cost of calculating …

---

## [Heat/Burgers using rectangular collocation via ClassicalOrthogonalPolynomials.jl and OrdinaryDiffEq.jl](https://discourse.julialang.org/t/heat-burgers-using-rectangular-collocation-via-classicalorthogonalpolynomials-jl-and-ordinarydiffeq-jl/129520)

<div class="topic-metadata">

**Author:** [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Replies:** 7\
**Last updated:** [June 4, 2025, 7:24pm UTC](https://discourse.julialang.org/t/heat-burgers-using-rectangular-collocation-via-classicalorthogonalpolynomials-jl-and-ordinarydiffeq-jl/129520 "2025-06-04T19:24:50Z")

</div>

I recently created an example of solving time-evolution PDEs with rectangular collocation and thought this might be helpful for others. This uses ClassicalOrthogonalPolynomials.jl, which has similar underlying mathematic…

---

## [Computing the adjugate/cofactors of a matrix](https://discourse.julialang.org/t/computing-the-adjugate-cofactors-of-a-matrix/125184)

<div class="topic-metadata">

**Author:** [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Replies:** 1\
**Last updated:** [June 3, 2025, 12:04am UTC](https://discourse.julialang.org/t/computing-the-adjugate-cofactors-of-a-matrix/125184 "2025-06-03T00:04:19Z")

</div>

For pedagogy (in our matrix calculus course at MIT), I wanted to have a routine to compute the adjugate of a matrix, since this arises in the gradient of the determinant. This is defined by: \\mathrm{adj}(A) = A^{-1} …

[Previous page](https://discourse.julialang.org/c/domain/numerics/20.md?page=1)

[Next page](https://discourse.julialang.org/c/domain/numerics/20.md?page=3)