# Numerics

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

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

**Page:** 24

---

## [Unexpected behavior of QuadGK with extended precision floats](https://discourse.julialang.org/t/unexpected-behavior-of-quadgk-with-extended-precision-floats/72684)

<div class="topic-metadata">

**Author:** [@PeterSimon](https://discourse.julialang.org/u/PeterSimon)\
**Replies:** 4\
**Last updated:** [December 7, 2021, 9:42pm UTC](https://discourse.julialang.org/t/unexpected-behavior-of-quadgk-with-extended-precision-floats/72684 "2021-12-07T21:42:36Z")

</div>

This topic arose from a topic in the Fortran Discourse group. Calling the gauss function from the QuadGK package produces different results depending on how it was previously called. Here is an example: …

---

## [Approximation matching asymptotes](https://discourse.julialang.org/t/approximation-matching-asymptotes/72648)

<div class="topic-metadata">

**Author:** [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Replies:** 3\
**Last updated:** [December 6, 2021, 7:36pm UTC](https://discourse.julialang.org/t/approximation-matching-asymptotes/72648 "2021-12-06T19:36:46Z")

</div>

I have a function f: \\mathbb{R} \\to \\mathbb{R} that is expensive to calculate (implicitly defined, needs iterative methods). But otherwise f is well-behaved: continuous, differentiable, and strictly increasing on all the…

---

## [Leading principal submatrix of Cholesky factor](https://discourse.julialang.org/t/leading-principal-submatrix-of-cholesky-factor/72456)

<div class="topic-metadata">

**Author:** [@DanielVandH](https://discourse.julialang.org/u/DanielVandH)\
**Replies:** 3\
**Last updated:** [December 3, 2021, 1:03am UTC](https://discourse.julialang.org/t/leading-principal-submatrix-of-cholesky-factor/72456 "2021-12-03T01:03:09Z")

</div>

I have an algorithm that uses only the first k rows and k columns of a Cholesky factor L (i.e., the kth leading principal submatrix L = L\[1:k, 1:k\]). Currently I just extract it after the computation of the whole Cholesk…

---

## [NLsolve with StaticArrays](https://discourse.julialang.org/t/nlsolve-with-staticarrays/72465)

<div class="topic-metadata">

**Author:** [@wsshin](https://discourse.julialang.org/u/wsshin)\
**Replies:** 3\
**Last updated:** [December 2, 2021, 4:10pm UTC](https://discourse.julialang.org/t/nlsolve-with-staticarrays/72465 "2021-12-02T16:10:44Z")

</div>

I am wondering if there is a way to make NLsolve work with StaticArrays, because my function is written to take SVector. I thought if I use the keyword argument inplace=false as described in here, I would be able to use…

---

## [Element-wise operations in Tullio.jl on GPU](https://discourse.julialang.org/t/element-wise-operations-in-tullio-jl-on-gpu/72186)

<div class="topic-metadata">

**Author:** [@cuihantao](https://discourse.julialang.org/u/cuihantao)\
**Replies:** 2\
**Last updated:** [November 28, 2021, 3:34pm UTC](https://discourse.julialang.org/t/element-wise-operations-in-tullio-jl-on-gpu/72186 "2021-11-28T15:34:47Z")

</div>

Hello, I wonder what is the correct way to implement element-wise operations in Tullio.jl on GPU. The MWE is as follows: using CUDA, Tullio using BenchmarkTools mul(A, B, C) = @tullio C\[k\] = A\[k\] \* B\[k\] a = rand(8000…

---

## [Complex step differentiation method explained](https://discourse.julialang.org/t/complex-step-differentiation-method-explained/14647)

<div class="topic-metadata">

**Author:** [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)\
**Replies:** 8\
**Last updated:** [November 27, 2021, 11:33am UTC](https://discourse.julialang.org/t/complex-step-differentiation-method-explained/14647 "2021-11-27T11:33:20Z")

</div>

I wrote a brief post explaining the “complex step method” that was the first topic in Nick Higham’s JuliaCon 2018 talk. (Using Julia for demonstrations, of course !)

---

## [Options for Sparse Distributed Linear Algebra](https://discourse.julialang.org/t/options-for-sparse-distributed-linear-algebra/71744)

<div class="topic-metadata">

**Author:** [@lucasbanting](https://discourse.julialang.org/u/lucasbanting)\
**Replies:** 2\
**Last updated:** [November 19, 2021, 7:52am UTC](https://discourse.julialang.org/t/options-for-sparse-distributed-linear-algebra/71744 "2021-11-19T07:52:21Z")

</div>

What is the current state of sparse distributed linear algebra in Julia? I am interested in switching from C++ to Julia for finite element codes. In C PETSc seems to be the best option for sparse distributed linear alge…

---

## [Getting non-vectorized code to the speed of vectorized ones](https://discourse.julialang.org/t/getting-non-vectorized-code-to-the-speed-of-vectorized-ones/71525)

<div class="topic-metadata">

**Author:** [@raphaelchinchilla](https://discourse.julialang.org/u/raphaelchinchilla)\
**Replies:** 22\
**Last updated:** [November 18, 2021, 8:42pm UTC](https://discourse.julialang.org/t/getting-non-vectorized-code-to-the-speed-of-vectorized-ones/71525 "2021-11-18T20:42:48Z")

</div>

I am not new to Julia, but I am not really a regular user either. I am trying to make a code that is non vectorized to achieve the performance of a vectorized one, but I cannot do it. Here is the problem in general lines…

---

## [Eigvals pinv](https://discourse.julialang.org/t/eigvals-pinv/67373)

<div class="topic-metadata">

**Author:** [@Jorge\_Silva](https://discourse.julialang.org/u/Jorge_Silva)\
**Replies:** 3\
**Last updated:** [November 18, 2021, 6:18am UTC](https://discourse.julialang.org/t/eigvals-pinv/67373 "2021-11-18T06:18:32Z")

</div>

I´m trying to implement a relatively complex model that I found in literature. It involves to calculate the spectrum of a linearized system in order to evaluate the stability of the steady-state solution. Thus I´m using …

---

## [Zygote: Differentiate over Zygote operation?](https://discourse.julialang.org/t/zygote-differentiate-over-zygote-operation/70853)

<div class="topic-metadata">

**Author:** [@ThummeTo](https://discourse.julialang.org/u/ThummeTo)\
**Replies:** 4\
**Last updated:** [November 16, 2021, 9:03am UTC](https://discourse.julialang.org/t/zygote-differentiate-over-zygote-operation/70853 "2021-11-16T09:03:37Z")

</div>

Hello everyone, I wondered, if it is possible to differentiate with Zygote over another Zygote-Operation (I attached a minimal example). using Zygote function myFunction(a) Zygote.jacobian(sin, a)\[1\] end Zygote.j…

---

## [Possible bug in determinant calculation](https://discourse.julialang.org/t/possible-bug-in-determinant-calculation/20323)

<div class="topic-metadata">

**Author:** [@Daniel\_Scheinerman](https://discourse.julialang.org/u/Daniel_Scheinerman)\
**Replies:** 21\
**Last updated:** [November 15, 2021, 8:58pm UTC](https://discourse.julialang.org/t/possible-bug-in-determinant-calculation/20323 "2021-11-15T20:58:46Z")

</div>

I found a (not poorly conditioned) integer matrix that julia computes the determinant incorrectly for. using LinearAlgebra A=\[1 1 0 0 1 0 0 0; 1 0 1 0 0 1 0 0; 1 0 0 1 0 0 1 0; 0 1 1 1 0 0 0 0; 0 1 0 0 0 0 1 1; 0 0 1 0 …

---

## [Can I save ODESystem separately and call it every time I need?](https://discourse.julialang.org/t/can-i-save-odesystem-separately-and-call-it-every-time-i-need/70872)

<div class="topic-metadata">

**Author:** [@shaurya](https://discourse.julialang.org/u/shaurya)\
**Replies:** 6\
**Last updated:** [November 13, 2021, 10:08am UTC](https://discourse.julialang.org/t/can-i-save-odesystem-separately-and-call-it-every-time-i-need/70872 "2021-11-13T10:08:56Z")

</div>

Hi! I am using ModelingToolkit to build a symbolic ODESystem of 13 vector equations. However, the system is unusually large and therefore, takes a lot of time to compute every time I run the code, which is not required.…

---

## [Making the result of \`nullspace\` sparse](https://discourse.julialang.org/t/making-the-result-of-nullspace-sparse/71310)

<div class="topic-metadata">

**Author:** [@jchen975](https://discourse.julialang.org/u/jchen975)\
**Replies:** 1\
**Last updated:** [November 12, 2021, 12:22pm UTC](https://discourse.julialang.org/t/making-the-result-of-nullspace-sparse/71310 "2021-11-12T12:22:07Z")

</div>

Is it possible to make the result of LinearAlgebra.nullspace as sparse as possible? How about when using the package Nemo for sparse matrices over \\mathbb{R}? Regardless of whether or not this is possible, what would be…

---

## [How to find crossing point using Spline1D in julia?](https://discourse.julialang.org/t/how-to-find-crossing-point-using-spline1d-in-julia/71364)

<div class="topic-metadata">

**Author:** [@Luigi\_Marongiu](https://discourse.julialang.org/u/Luigi_Marongiu)\
**Replies:** 2\
**Last updated:** [November 12, 2021, 9:38am UTC](https://discourse.julialang.org/t/how-to-find-crossing-point-using-spline1d-in-julia/71364 "2021-11-12T09:38:08Z")

</div>

I have fitted a spline function using Spline1D from the package Dierckx to the data: x = \[0.0, 0.25167322597331027, 1.6518531853618088, 4.032506892599284, 6.871304792339401, 10.090447889235595, 13.686254426579751, 17.61…

---

## [Fast Legendre Transform?](https://discourse.julialang.org/t/fast-legendre-transform/68152)

<div class="topic-metadata">

**Author:** [@Veenty](https://discourse.julialang.org/u/Veenty)\
**Replies:** 3\
**Last updated:** [November 6, 2021, 7:13am UTC](https://discourse.julialang.org/t/fast-legendre-transform/68152 "2021-11-06T07:13:04Z")

</div>

Is anyone aware of a package that is able to do a Fast Legendre Transform similar to the one in chebfun (The fast discrete Legendre transform » Chebfun). The closest that I could find was something in FastTransform, whe…

---

## [Tridiagonal Toeplitz matrix complex eigenvalues instead of real](https://discourse.julialang.org/t/tridiagonal-toeplitz-matrix-complex-eigenvalues-instead-of-real/70892)

<div class="topic-metadata">

**Author:** [@jamblejoe](https://discourse.julialang.org/u/jamblejoe)\
**Replies:** 7\
**Last updated:** [November 4, 2021, 10:08am UTC](https://discourse.julialang.org/t/tridiagonal-toeplitz-matrix-complex-eigenvalues-instead-of-real/70892 "2021-11-04T10:08:17Z")

</div>

The spectrum of tridiagonal Toeplitz matrices is known analytically (Tridiagonal-Toeplitz-Matrix – Wikipedia) and is (for non-zero off-diagonals) simple and real. Julia’s eigen function returns for small non-symmetric To…

---

## [Algorithms for finding the self-intersection of a plane curve?](https://discourse.julialang.org/t/algorithms-for-finding-the-self-intersection-of-a-plane-curve/70192)

<div class="topic-metadata">

**Author:** [@liuyxpp](https://discourse.julialang.org/u/liuyxpp)\
**Replies:** 14\
**Last updated:** [October 29, 2021, 8:45pm UTC](https://discourse.julialang.org/t/algorithms-for-finding-the-self-intersection-of-a-plane-curve/70192 "2021-10-29T20:45:37Z")

</div>

Does anyone here happen to know about algorithms for finding the self-intersection points of a plane curve? For example, find the self-intersection (X, Y) = (0, 0) in the below plane curve. Edit: Add an actual exampl…

---

## [MPLAPACK (multiprecision blas+lapack)](https://discourse.julialang.org/t/mplapack-multiprecision-blas-lapack/70514)

<div class="topic-metadata">

**Author:** [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Replies:** 8\
**Last updated:** [October 28, 2021, 8:32pm UTC](https://discourse.julialang.org/t/mplapack-multiprecision-blas-lapack/70514 "2021-10-28T20:32:20Z")

</div>

MPLAPACK works with any of these multiprecision libs: GMP, MPFR+MPC, libquadmath, and David Bailey’s DD/QD library. Here is the user guide and the repository, All BLAS routines are available using multiple precision ar…

---

## [Implementation of Jacobi method not yielding acurate results](https://discourse.julialang.org/t/implementation-of-jacobi-method-not-yielding-acurate-results/70511)

<div class="topic-metadata">

**Author:** [@scaredFish](https://discourse.julialang.org/u/scaredFish)\
**Replies:** 2\
**Last updated:** [October 28, 2021, 9:17am UTC](https://discourse.julialang.org/t/implementation-of-jacobi-method-not-yielding-acurate-results/70511 "2021-10-28T09:17:37Z")

</div>

Hello, I am fairly new to the language and i have tried to implement the Jacobi method. However, this has not given me the expected results. The solution should be (1,2,-1,1) and the returned values at convergence are (…

---

## [Describing a particular second-order equation in DifferentialEquations.jl](https://discourse.julialang.org/t/describing-a-particular-second-order-equation-in-differentialequations-jl/53253)

<div class="topic-metadata">

**Author:** [@anC](https://discourse.julialang.org/u/anC)\
**Replies:** 8\
**Last updated:** [October 24, 2021, 7:49am UTC](https://discourse.julialang.org/t/describing-a-particular-second-order-equation-in-differentialequations-jl/53253 "2021-10-24T07:49:15Z")

</div>

I am using DifferentialEquations.jl . If I wanted to solve a simple second order equation of the form: a t^2 f’‘(t) + b t f’(t) +c = 0 (for a,b,c real constants) I think the way to describe the problem in dif…

---

## [Lyapunov Spectrum simple code problem](https://discourse.julialang.org/t/lyapunov-spectrum-simple-code-problem/70182)

<div class="topic-metadata">

**Author:** [@Gabriel\_Carita](https://discourse.julialang.org/u/Gabriel_Carita)\
**Replies:** 2\
**Last updated:** [October 22, 2021, 12:43am UTC](https://discourse.julialang.org/t/lyapunov-spectrum-simple-code-problem/70182 "2021-10-22T00:43:20Z")

</div>

Hello there, I’m trying to figure out how to solve my lyapunov spectrum handwrite code using the QR method for a simple example (lorenz system), there is something wrong with my code ? I’m trying to reproduce the values…

---

## [DifferentialEquations.jl and ReverseDiff.jl: MethodError: Cannot \`convert\` an object of type Vector{ReverseDiff.TrackedReal{](https://discourse.julialang.org/t/differentialequations-jl-and-reversediff-jl-methoderror-cannot-convert-an-object-of-type-vector-reversediff-trackedreal/70175)

<div class="topic-metadata">

**Author:** [@salazardetroya](https://discourse.julialang.org/u/salazardetroya)\
**Replies:** 1\
**Last updated:** [October 21, 2021, 7:43pm UTC](https://discourse.julialang.org/t/differentialequations-jl-and-reversediff-jl-methoderror-cannot-convert-an-object-of-type-vector-reversediff-trackedreal/70175 "2021-10-21T19:43:49Z")

</div>

I am trying to get gradients in a heat equation solved with DifferentialEquations.jl: using DiffEqOperators, OrdinaryDiffEq, Plots, Zygote, DiffEqSensitivity nknots = 10 h = 1.0 / (nknots + 1) knots = range(h, step=h, …

---

## [DiffEqOperators and gradients](https://discourse.julialang.org/t/diffeqoperators-and-gradients/70121)

<div class="topic-metadata">

**Author:** [@salazardetroya](https://discourse.julialang.org/u/salazardetroya)\
**Replies:** 2\
**Last updated:** [October 21, 2021, 4:00pm UTC](https://discourse.julialang.org/t/diffeqoperators-and-gradients/70121 "2021-10-21T16:00:05Z")

</div>

I am trying to solve the following problem \\begin{aligned} \\nabla \\cdot \\nabla \\phi\_{1} &= \\frac{\\partial\\left(\\phi\_{1}-\\phi\_{2}\\right)}{\\partial t}\\\\ \\nabla \\cdot \\nabla \\phi\_{2} &= -\\tau p\_1p\_2\\frac{\\partial\\left(\\…

---

## [Julia in fusion research](https://discourse.julialang.org/t/julia-in-fusion-research/26460)

<div class="topic-metadata">

**Author:** [@johnh](https://discourse.julialang.org/u/johnh)\
**Replies:** 10\
**Last updated:** [October 21, 2021, 1:18pm UTC](https://discourse.julialang.org/t/julia-in-fusion-research/26460 "2021-10-21T13:18:56Z")

</div>

Following on fromn the discussion on PIC codes. Is there much activity in using Julia for fusion and plasma simulations? Any pointers gratefully received.

---

## [Type instability in ForwardDiff.gradient](https://discourse.julialang.org/t/type-instability-in-forwarddiff-gradient/68973)

<div class="topic-metadata">

**Author:** [@ElOceanografo](https://discourse.julialang.org/u/ElOceanografo)\
**Replies:** 8\
**Last updated:** [October 20, 2021, 12:04am UTC](https://discourse.julialang.org/t/type-instability-in-forwarddiff-gradient/68973 "2021-10-20T00:04:14Z")

</div>

Working on a nested optimization problem, I noticed that a lot of time was being spent in gradient calls, and a little digging found that doing out-of-place gradients was introducing a type-instability. I followed the a…

---

## [Error using Differential Equations and Unitful](https://discourse.julialang.org/t/error-using-differential-equations-and-unitful/69909)

<div class="topic-metadata">

**Author:** [@MisterBiggs](https://discourse.julialang.org/u/MisterBiggs)\
**Replies:** 5\
**Last updated:** [October 17, 2021, 4:59pm UTC](https://discourse.julialang.org/t/error-using-differential-equations-and-unitful/69909 "2021-10-17T16:59:54Z")

</div>

I am using Differential Equations and Unitful together for the first time and I am running into the following error when trying to solve my OdeProblem: ArgumentError: zero(Unitful.Quantity{Float64}) not defined. My ODE …

---

## [\*Cache and \*ConstantCache versions of ODE solvers. What is the difference?](https://discourse.julialang.org/t/cache-and-constantcache-versions-of-ode-solvers-what-is-the-difference/69810)

<div class="topic-metadata">

**Author:** [@Andrey](https://discourse.julialang.org/u/Andrey)\
**Replies:** 1\
**Last updated:** [October 16, 2021, 1:30pm UTC](https://discourse.julialang.org/t/cache-and-constantcache-versions-of-ode-solvers-what-is-the-difference/69810 "2021-10-16T13:30:34Z")

</div>

I was going through Rodas4 solver source code and found that there are 2 different implementations of this solver - Rodas4Cache and Rodas4ConstantCache. Is there any principal difference between these 2 variants? It lo…

---

## [Interval Compiler for Sound Float Computations](https://discourse.julialang.org/t/interval-compiler-for-sound-float-computations/69857)

<div class="topic-metadata">

**Author:** [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Replies:** 0\
**Last updated:** [October 15, 2021, 9:38pm UTC](https://discourse.julialang.org/t/interval-compiler-for-sound-float-computations/69857 "2021-10-15T21:38:55Z")

</div>

For those interested in llvm mods for enhancing floats. IGen is a source-to-source compiler for interval arithmetic. IGen takes as input a C function doing floating-point computations, possibly using SIMD intrinsics, a…

---

## [Argument "tgrad" in ODEFunction appears to have no influence. Correct?](https://discourse.julialang.org/t/argument-tgrad-in-odefunction-appears-to-have-no-influence-correct/69236)

<div class="topic-metadata">

**Author:** [@Andrey](https://discourse.julialang.org/u/Andrey)\
**Replies:** 2\
**Last updated:** [October 5, 2021, 2:53pm UTC](https://discourse.julialang.org/t/argument-tgrad-in-odefunction-appears-to-have-no-influence-correct/69236 "2021-10-05T14:53:44Z")

</div>

\[Question replicated from DifferentialEquations Lobby at JuliaDiffEq/Lobby - Gitter \] Hi All. I have got a question to the “DifferentialEquations” and “ModelingToolbox” Pkg developers, probably the right addressees are…

---

## [Is there a package for fractional calculus?](https://discourse.julialang.org/t/is-there-a-package-for-fractional-calculus/19702)

<div class="topic-metadata">

**Author:** [@phelipe](https://discourse.julialang.org/u/phelipe)\
**Replies:** 10\
**Last updated:** [October 4, 2021, 5:46am UTC](https://discourse.julialang.org/t/is-there-a-package-for-fractional-calculus/19702 "2021-10-04T05:46:20Z")

</div>

Is there any package for fractional derivatives?

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

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