Optimize and fix thermal_dm for nbar=0

[yuyichao]

Checklist

• Please read Contributing to Quantum Toolbox in Julia.
• Any code changes were done in a way that does not break public API.
• Appropriate tests were added and tested locally by running: make test.
• Any code changes should be julia formatted by running: make format.
• All documents (in docs/ folder) related to code changes were updated and able to build locally by running: make docs.
• (If necessary) the CHANGELOG.md should be updated (regarding to the code changes) and built by running: make changelog.

Description

Fix thermal_dm for nbar=0, which previous returns NaN results.

julia> thermal_dm(3, 0.0)

Quantum Object:   type=Operator()   dims=[3]   size=(3, 3)   ishermitian=false
3×3 Matrix{ComplexF64}:
 NaN+NaN*im  0.0+0.0im   0.0+0.0im
 0.0+0.0im   NaN+NaN*im  0.0+0.0im
 0.0+0.0im   0.0+0.0im   NaN+NaN*im

The main fix is to avoid evaluating exp(-β * 0) when filling the diagonal data array. Other changes that improves performances includes

• pre-computing the normalization factor to avoid looping through the array twice
• use precomputed value/specialized expressions for small j (this mainly matters for BigFloat results and doesn't hurt for machine size numbers)
• Use a loop directly since the syntax of using branches in comprehension, while allowed and totally works, is not as elegant anymore. It also reduces reliance on inlining a little bit.

Here's the measured performance difference compared to master. For dense matrix at larger size the performance is going to be dominated by filling the matrix but for most other cases the performance difference is ~10-30%.

--- master      2025-12-15 17:55:40.759069248 -0500
+++ pr  2025-12-15 17:55:50.030156180 -0500
@@ -1,35 +1,35 @@
 julia> @btime thermal_dm($(3), $(0.2));
-  97.823 ns (6 allocations: 448 bytes)
+  62.755 ns (4 allocations: 336 bytes)

 julia> @btime thermal_dm($(30), $(0.2));
-  824.015 ns (7 allocations: 15.27 KiB)
+  738.646 ns (5 allocations: 14.71 KiB)

 julia> @btime thermal_dm($(300), $(0.2));
-  33.833 μs (9 allocations: 1.38 MiB)
+  32.459 μs (6 allocations: 1.38 MiB)

 julia> @btime thermal_dm($(3), $(big(0.2)));
-  8.222 μs (73 allocations: 6.34 KiB)
+  7.458 μs (37 allocations: 3.04 KiB)

 julia> @btime thermal_dm($(30), $(big(0.2)));
-  73.375 μs (641 allocations: 69.90 KiB)
+  65.209 μs (281 allocations: 35.77 KiB)

 julia> @btime thermal_dm($(300), $(big(0.2)));
-  780.091 μs (6361 allocations: 1.91 MiB)
+  696.049 μs (2761 allocations: 1.58 MiB)

 julia> @btime thermal_dm($(3), $(0.2); sparse=$(Val(true)));
-  312.675 ns (25 allocations: 1.17 KiB)
+  222.500 ns (23 allocations: 1.06 KiB)

 julia> @btime thermal_dm($(30), $(0.2); sparse=$(Val(true)));
-  825.236 ns (25 allocations: 4.94 KiB)
+  564.070 ns (23 allocations: 4.38 KiB)

 julia> @btime thermal_dm($(300), $(0.2); sparse=$(Val(true)));
-  6.017 μs (37 allocations: 40.95 KiB)
+  4.375 μs (34 allocations: 36.20 KiB)

 julia> @btime thermal_dm($(3), $(big(0.2)); sparse=$(Val(true)));
-  8.250 μs (82 allocations: 6.13 KiB)
+  7.458 μs (46 allocations: 2.84 KiB)

 julia> @btime thermal_dm($(30), $(big(0.2)); sparse=$(Val(true)));
-  72.292 μs (595 allocations: 53.56 KiB)
+  63.876 μs (235 allocations: 19.44 KiB)

 julia> @btime thermal_dm($(300), $(big(0.2)); sparse=$(Val(true)));
-  723.216 μs (5785 allocations: 526.22 KiB)
+  640.298 μs (2185 allocations: 184.15 KiB)

[codecov]

Codecov Report

❌ Patch coverage is 0% with 14 lines in your changes missing coverage. Please review.
✅ Project coverage is 19.33%. Comparing base (0c77d51) to head (f0cb438).
⚠️ Report is 4 commits behind head on main.

Files with missing lines	Patch %	Lines
src/qobj/states.jl	0.00%	14 Missing ⚠️

❗ There is a different number of reports uploaded between BASE (0c77d51) and HEAD (f0cb438). Click for more details.

HEAD has 8 uploads less than BASE

Flag	BASE (0c77d51)	HEAD (f0cb438)
	10	2

Additional details and impacted files

@@             Coverage Diff             @@
##             main     #613       +/-   ##
===========================================
- Coverage   92.28%   19.33%   -72.95%     
===========================================
  Files          50       50               
  Lines        3615     3589       -26     
===========================================
- Hits         3336      694     -2642     
- Misses        279     2895     +2616     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[ytdHuang]

Thank you for reporting this issue.

For n = 0, we can just return fock_dm(N, 0). However, your code doesn work for n = Inf.

Thus, I have made an extended fix in PR #614 .