Skip to content

Commit

Permalink
Improved log1p implementation for complex inputs (#107100)
Browse files Browse the repository at this point in the history
This PR improves the implementation of `torch.log1p` for complex inputs as mentioned in issue #107022. The new implementation is based on the insights provided in numpy/numpy#22611 (comment).

Pull Request resolved: #107100
Approved by: https://github.com/lezcano
  • Loading branch information
tringwald authored and pytorchmergebot committed Aug 16, 2023
1 parent 35cca79 commit 4de5e17
Showing 1 changed file with 18 additions and 0 deletions.
18 changes: 18 additions & 0 deletions c10/util/complex_math.h
Original file line number Diff line number Diff line change
Expand Up @@ -291,6 +291,11 @@ C10_HOST_DEVICE inline c10::complex<T> atanh(const c10::complex<T>& x) {

template <typename T>
C10_HOST_DEVICE inline c10::complex<T> log1p(const c10::complex<T>& z) {
#if defined(__APPLE__) || defined(__MACOSX)
// For Mac, the new implementation yielded a high relative error. Falling back
// to the old version for now.
// See https://github.com/numpy/numpy/pull/22611#issuecomment-1667945354

// log1p(z) = log(1 + z)
// Let's define 1 + z = r * e ^ (i * a), then we have
// log(r * e ^ (i * a)) = log(r) + i * a
Expand All @@ -316,6 +321,19 @@ C10_HOST_DEVICE inline c10::complex<T> log1p(const c10::complex<T>& z) {
T z0 = std::hypot(x + 1, y);
return {std::log(z0), theta};
}
#else
// Based on https://github.com/numpy/numpy/pull/22611#issuecomment-1667945354
c10::complex<T> u = z + T(1);
if (u == T(1)) {
return z;
} else {
auto log_u = log(u);
if (u - T(1) == z) {
return log_u;
}
return log_u * (z / (u - T(1)));
}
#endif
}

template <typename T>
Expand Down

0 comments on commit 4de5e17

Please sign in to comment.