To compute A + B within float32 where both A and B is ~1e100, we can do:
log(A + B) = log(exp(a) + exp(b)) = a + softplus(b - a)where a = log(A) and b = log(B).
Created
April 27, 2022 09:31
-
-
Save xrq-phys/96b7e05f288e4f15624edc1ae5013463 to your computer and use it in GitHub Desktop.
Complex-valued Softplus in C++. Useful when summing two VERY LARGE numbers.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| /** | |
| * This Source Code Form is subject to the terms of the Mozilla Public | |
| * License, v. 2.0. If a copy of the MPL was not distributed with this | |
| * file, You can obtain one at https://mozilla.org/MPL/2.0/. | |
| */ | |
| #pragma once | |
| #include <complex> | |
| #include <algorithm> | |
| // Real softplus. | |
| template <typename T> | |
| T softplus(T x) | |
| { | |
| using namespace std; | |
| return log(1 + exp(-abs(x))) + max(x, T(0)); | |
| } | |
| // Complex softplus. | |
| template <typename T> | |
| std::complex<T> softplus(std::complex<T> z) | |
| { | |
| using namespace std; | |
| T r = real(z); | |
| T i = imag(z); | |
| T logval0 = softplus(2*r); | |
| T logval2; | |
| T ci = cos(i); | |
| if (ci > 0) { | |
| T logval1 = r + log(2*ci); | |
| logval2 = 0.5*(logval0 + softplus(logval1 - logval0)); | |
| } else { | |
| T logval1 = r + log(2*abs(ci)); | |
| logval2 = 0.5*(logval0 + log(1 - exp(logval1 - logval0))); | |
| } | |
| T arg = asin(exp(r - logval2) * sin(i)); | |
| if (ci < 0 && r > log(-1/ci)) { | |
| arg = M_PI - arg; | |
| } | |
| return std::complex<T>(logval2, arg); | |
| } | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment