# Julia vs NumPy broadcasting

**URL:** <https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908>\
**Category:** New to Julia\
**Created:** [July 21, 2025, 2:09pm UTC](https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908 "2025-07-21T14:09:05Z")\
**Posts on this page:** 4\
**Page:** 2

<div class="post-metadata">

**Author:** ![giordano](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/giordano/32/2166_2.png) [@giordano](https://discourse.julialang.org/u/giordano)\
**Post date:** [July 22, 2025, 11:07am UTC](https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908/21 "2025-07-22T11:07:31Z")

</div>

It works on macOS, too. Only major operating system currently missing is Windows.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [July 22, 2025, 1:16pm UTC](https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908/22 "2025-07-22T13:16:15Z")

</div>

> [@Donald](#):
>
> Is this approach of repeatedly multiplying by e^A recommended when computing e^{At} for integer t? I imagine you would get a build-up of rounding errors.

The roundoff errors should grow as O(\sqrt{n}) on average. In double precision, with n of a few thousand, it’s probably not worth worrying about?

(It’s the same error growth as computing a sum by repeated addition, and most people don’t worry about that either, though in a library function like `sum` we try to do better.)

That being said, if it’s a problem you could always recompute e^{kA} explicitly via `exp(k*A)` every N steps for some value of N, or there are other schemes to trade off computation for accuracy.

> [@Donald](#):
>
> However, as part of a larger optimisation problem, it’s possible that the code could be run with unfortunate inputs such that A is nearly singular. I was using some regularization on eigenvalues of A for this reason (also to keep them from getting too close to each other).

One option in that case is @mikmoore’s augmented-matrix approach, which handles arbitrary singular (and defective) A but is more expensive.

However, if you are willing to regularize the eigenvalues (i.e. change the solution when A is nearly singular, since I guess this is only an intermediate step in the optimization), maybe you could apply a regularization directly to the matrix (or your parameterization thereof) or to your system of equations. For example, you could replace `A \ B` in my code above with `[A; α*one(A)] \ [B; 0*B]` for some small value of `α`, e.g. `α = 1e-8` — this is a [Tikhonov regularization](https://en.wikipedia.org/wiki/Ridge_regression), and is equivalent to `A \ B` for \alpha \ll \Vert A \Vert.

(For optimization, even derivative-free optimization, you really want to keep things smooth as a function of the parameters. For example, a Tikhonov regularization is smooth as a function of A, whereas arbitrarily changing the eigenvalues when they get too small may not be, depending on how you do it.)

---

<div class="post-metadata">

**Author:** ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)\
**Post date:** [July 22, 2025, 2:22pm UTC](https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908/23 "2025-07-22T14:22:47Z")

</div>

> [@stevengj](#):
>
> The roundoff errors should grow as O(\sqrt{n}) on average. In double precision, with n of a few thousand, it’s probably not worth worrying about?

If A is a stable matrix (`real(eigvals(A)) < 0`), any numerical error or difference in initial condition etc. should be forgotten exponentially with a rate given by the real part of the eigenvalues of A?

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [July 22, 2025, 2:51pm UTC](https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908/24 "2025-07-22T14:51:54Z")

</div>

> [@baggepinnen](#):
>
> any numerical error or difference in initial condition etc. should be forgotten exponentially

Only in an absolute-error sense — the relative error still accumulates if you compute e^{Ak} by repeated multiplication.

[Previous page](https://discourse.julialang.org/t/julia-vs-numpy-broadcasting/130908.md?page=1)
