Skip to content

BUG: amos.h: fix bugs in buni, uni1, uoik#145

Merged
steppi merged 1 commit into
scipy:mainfrom
JaRoSchm:amos_buni_uni1_uoik_fix
May 14, 2026
Merged

BUG: amos.h: fix bugs in buni, uni1, uoik#145
steppi merged 1 commit into
scipy:mainfrom
JaRoSchm:amos_buni_uni1_uoik_fix

Conversation

@JaRoSchm

Copy link
Copy Markdown
Contributor

Reference issue

What does this implement/fix?

Fixes translation bugs in the functions buni, uni1, and uoik in amos.h. These led to wrong results for besj with z = -50.0 - 50.0j, fnu = 0.0, kode = 1, n = 1000. I have also added a test inspired by #92 (comment), which fails on main but passes on this branch. I replaced the test from #138 with this approach here and added the new test case, as this is much simpler going forward compared to the old approach.

Additional information

AI Generation Disclosure

Comment thread include/xsf/amos/amos.h
xx = std::real(z);
yy = std::imag(z);
ax = std::fabs(xx) + std::sqrt(3.);
ax = std::fabs(xx) * std::sqrt(3.);

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Comment thread include/xsf/amos/amos.h
/* 80 */
if (nd > 2) {
rz = 1.0 / z;
rz = 2.0 / z;

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Comment thread include/xsf/amos/amos.h
zn = -zr * std::complex<double>(0, 1);
if (yy <= 0.) {
zn = std::conj(zn);
zn.real(-zn.real());

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Comment thread include/xsf/amos/amos.h
if (rcz > -elim) {
ascle = 1e3 * d1mach[0] / tol;
cz = std::log(phi);
cz += std::log(phi);

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Comment thread include/xsf/amos/amos.h
ax = std::exp(rcz) / tol;
ay = std::imag(cz);
cz = ax * (std::cos(ay) + std::sin(ay * std::complex<double>(0, 1)));
cz = ax * (std::cos(ay) + std::sin(ay) * std::complex<double>(0, 1));

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@steppi steppi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks @JaRoSchm! The links to the original Fortran source were very helpful.

@steppi steppi merged commit 409ce44 into scipy:main May 14, 2026
4 of 6 checks passed
@JaRoSchm

Copy link
Copy Markdown
Contributor Author

Thank you! If time allows would you please have a look at #92. After this is merged, I could also add a test to #93.

@JaRoSchm JaRoSchm deleted the amos_buni_uni1_uoik_fix branch May 15, 2026 12:40
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants