From 18ab5292cdd8cab7ea018fbfd342afe96fb3030a Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Sun, 22 Mar 2026 17:02:46 -0700
Subject: [PATCH 01/12] add lambertw_pvlib

---
 pvlib/tools.py      | 56 +++++++++++++++++++++++++++++++++++++++++++++
 tests/test_tools.py | 25 ++++++++++++++++++++
 2 files changed, 81 insertions(+)

diff --git a/pvlib/tools.py b/pvlib/tools.py
index 6cb631f852..da90a98259 100644
--- a/pvlib/tools.py
+++ b/pvlib/tools.py
@@ -586,3 +586,59 @@ def _file_context_manager(filename_or_object, mode='r', encoding=None):
         # otherwise, assume a filename or path
         context = open(str(filename_or_object), mode=mode, encoding=encoding)
     return context
+
+
+def _log_lambertw(logx):
+    r'''Computes W(x) starting from log(x).
+
+    Parameters
+    ----------
+    logx : numeric
+        Log(x) of
+
+    Returns
+    -------
+    numeric
+        Lambert's W(x)
+
+    '''
+    # handles overflow cases, but results in nan for x <= 1
+    w = logx - np.log(logx)  # initial guess, w = log(x) - log(log(x))
+
+    for _ in range(0, 3):
+        # Newton's. Halley's is not substantially faster or more accurate
+        # because f''(w) = -1 / (w**2) is small for large w
+        w = w * (1. - np.log(w) + logx) / (1. + w)
+    return w
+
+
+def lambertw_pvlib(x):
+    r'''Lambert's W function.
+
+    Parameters
+    ----------
+    x : numeric
+        Must be real numbers.
+
+    Returns
+    -------
+    numeric
+        Lambert's W(x). Principal branch only.
+
+    '''
+    w = np.full_like(x, np.nan)
+    small = x <= 10
+    # for large x, solve 0 = f(w) = w + log(w) - log(x) using Newton's
+    w[~small] = _log_lambertw(np.log(x[~small]))
+
+    # w will contain nan for these numbers due to log(w) = log(log(x))
+    # for small x, solve 0 = g(w) = w * exp(w) - x using Halley's method
+    if any(small):
+        z = x[small]
+        g = np.log(x[small] + 1) - np.log(np.log(x[small] + 1) + 1)
+        for _ in range(0, 3):
+            expg = np.exp(g)
+            g = g - (g*expg - z) * (g + 1) / (expg * (g + 1)**2 - 0.5*(g + 2)*(expg*g - z))
+        w[small] = g
+
+    return w
diff --git a/tests/test_tools.py b/tests/test_tools.py
index 4b733ad711..d053b24c08 100644
--- a/tests/test_tools.py
+++ b/tests/test_tools.py
@@ -283,3 +283,28 @@ def test__file_context_manager():
     buffer = StringIO("test content")
     with tools._file_context_manager(buffer) as obj:
         assert obj.read() == "test content"
+
+
+def test_lambertw_pvlib():
+    test_exp = np.arange(-10., 300, step=10)
+    test_x = 10.**test_exp
+    # known solution from scipy.special.lambertw
+    expected = np.array([
+        9.9999999989999997e-11, 5.6714329040978384e-01,
+        2.0028685413304952e+01, 4.2306755091738395e+01,
+        6.4904633770046118e+01, 8.7630277151947183e+01,
+        1.1042491882731335e+02, 1.3326278259180333e+02,
+        1.5613026581351718e+02, 1.7901931374150624e+02,
+        2.0192476320084489e+02, 2.2484310644511851e+02,
+        2.4777185185877809e+02, 2.7070916610249782e+02,
+        2.9365366103997610e+02, 3.1660426041503479e+02,
+        3.3956011295458728e+02, 3.6252053376149752e+02,
+        3.8548496362161768e+02, 4.0845294003314166e+02,
+        4.3142407612718210e+02, 4.5439804503371403e+02,
+        4.7737456808796901e+02, 5.0035340579834485e+02,
+        5.2333435083468805e+02, 5.4631722251791496e+02,
+        5.6930186244110166e+02, 5.9228813095427859e+02,
+        6.1527590431628334e+02, 6.3826507236734335e+02,
+        6.6125553661218726e+02])
+    result = tools.lambertw_pvlib(test_x)
+    assert np.allclose(result, expected, rtol=1e-14)

From 33be0e3721dc25d85c9491eccf10a11776dd0c97 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Mon, 23 Mar 2026 08:22:51 -0700
Subject: [PATCH 02/12] whatsnew

---
 docs/sphinx/source/whatsnew/v0.15.1.rst | 7 +++++--
 1 file changed, 5 insertions(+), 2 deletions(-)

diff --git a/docs/sphinx/source/whatsnew/v0.15.1.rst b/docs/sphinx/source/whatsnew/v0.15.1.rst
index f86fed987b..504efb588b 100644
--- a/docs/sphinx/source/whatsnew/v0.15.1.rst
+++ b/docs/sphinx/source/whatsnew/v0.15.1.rst
@@ -23,9 +23,12 @@ Bug fixes
 
 Enhancements
 ~~~~~~~~~~~~
-* Use ``k`` and ``cap_adjustment`` from :py:func:`pvlib.pvsystem.Array.module_parameters` in :py:func:`pvlib.pvsystem.PVSystem.pvwatts_dc`
+* Use ``k`` and ``cap_adjustment`` from :py:func:`pvlib.pvsystem.Array.module_parameters`
+  in :py:func:`pvlib.pvsystem.PVSystem.pvwatts_dc`
   (:issue:`2714`, :pull:`2715`)
-  
+* Add :py:func:`pvlib.tools.lambertw_pvlib` to speed up calculations when using
+  LambertW single diode equation methods.
+  (:pull:`2723`)  
 
 Documentation
 ~~~~~~~~~~~~~

From e2fa011874a47ce26aa70fa00abf6fc409aeef16 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Mon, 23 Mar 2026 08:25:55 -0700
Subject: [PATCH 03/12] formatting

---
 pvlib/tools.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/pvlib/tools.py b/pvlib/tools.py
index da90a98259..c0aecc6cca 100644
--- a/pvlib/tools.py
+++ b/pvlib/tools.py
@@ -638,7 +638,8 @@ def lambertw_pvlib(x):
         g = np.log(x[small] + 1) - np.log(np.log(x[small] + 1) + 1)
         for _ in range(0, 3):
             expg = np.exp(g)
-            g = g - (g*expg - z) * (g + 1) / (expg * (g + 1)**2 - 0.5*(g + 2)*(expg*g - z))
+            g = g - (g*expg - z) * (g + 1) / \
+                (expg * (g + 1)**2 - 0.5*(g + 2)*(expg*g - z))
         w[small] = g
 
     return w

From e0bc631ff6e1a0bc6659bab717ad27c0933e3943 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Mon, 23 Mar 2026 17:58:57 -0700
Subject: [PATCH 04/12] Apply suggestions from code review

Co-authored-by: Echedey Luis <80125792+echedey-ls@users.noreply.github.com>
---
 docs/sphinx/source/whatsnew/v0.15.1.rst | 2 +-
 pvlib/tools.py                          | 7 +++----
 2 files changed, 4 insertions(+), 5 deletions(-)

diff --git a/docs/sphinx/source/whatsnew/v0.15.1.rst b/docs/sphinx/source/whatsnew/v0.15.1.rst
index 504efb588b..a4969597dd 100644
--- a/docs/sphinx/source/whatsnew/v0.15.1.rst
+++ b/docs/sphinx/source/whatsnew/v0.15.1.rst
@@ -28,7 +28,7 @@ Enhancements
   (:issue:`2714`, :pull:`2715`)
 * Add :py:func:`pvlib.tools.lambertw_pvlib` to speed up calculations when using
   LambertW single diode equation methods.
-  (:pull:`2723`)  
+  (:discussion:`2720`, :pull:`2723`)  
 
 Documentation
 ~~~~~~~~~~~~~
diff --git a/pvlib/tools.py b/pvlib/tools.py
index c0aecc6cca..baee899f64 100644
--- a/pvlib/tools.py
+++ b/pvlib/tools.py
@@ -613,17 +613,16 @@ def _log_lambertw(logx):
 
 
 def lambertw_pvlib(x):
-    r'''Lambert's W function.
+    r'''Lambert's W function principal branch, :math:`W_0(x)`, for :math:`x` real valued.
 
     Parameters
     ----------
-    x : numeric
+    x : np.array
         Must be real numbers.
 
     Returns
     -------
-    numeric
-        Lambert's W(x). Principal branch only.
+    np.array
 
     '''
     w = np.full_like(x, np.nan)

From 935ff3544f070f1e0f3b807c1d617fc970f75f5e Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Mon, 23 Mar 2026 18:02:16 -0700
Subject: [PATCH 05/12] from review

---
 pvlib/tools.py      | 3 ++-
 tests/test_tools.py | 1 +
 2 files changed, 3 insertions(+), 1 deletion(-)

diff --git a/pvlib/tools.py b/pvlib/tools.py
index baee899f64..e0f7cc43ae 100644
--- a/pvlib/tools.py
+++ b/pvlib/tools.py
@@ -613,7 +613,8 @@ def _log_lambertw(logx):
 
 
 def lambertw_pvlib(x):
-    r'''Lambert's W function principal branch, :math:`W_0(x)`, for :math:`x` real valued.
+    r'''Lambert's W function principal branch, :math:`W_0(x)`, for :math:`x`
+    real valued.
 
     Parameters
     ----------
diff --git a/tests/test_tools.py b/tests/test_tools.py
index d053b24c08..f316cb4605 100644
--- a/tests/test_tools.py
+++ b/tests/test_tools.py
@@ -289,6 +289,7 @@ def test_lambertw_pvlib():
     test_exp = np.arange(-10., 300, step=10)
     test_x = 10.**test_exp
     # known solution from scipy.special.lambertw
+    # scipy 1.7.1, python 3.13.1, numpy 2.3.5
     expected = np.array([
         9.9999999989999997e-11, 5.6714329040978384e-01,
         2.0028685413304952e+01, 4.2306755091738395e+01,

From c22c7636c36db9203a980c98d44f6042430e19ed Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Mon, 23 Mar 2026 18:05:26 -0700
Subject: [PATCH 06/12] move to ivtools, make private

---
 pvlib/ivtools/utils.py      | 57 +++++++++++++++++++++++++++++++++++++
 pvlib/tools.py              | 57 -------------------------------------
 tests/ivtools/test_utils.py | 27 ++++++++++++++++++
 tests/test_tools.py         | 26 -----------------
 4 files changed, 84 insertions(+), 83 deletions(-)

diff --git a/pvlib/ivtools/utils.py b/pvlib/ivtools/utils.py
index cde50655dc..25ed4f2c84 100644
--- a/pvlib/ivtools/utils.py
+++ b/pvlib/ivtools/utils.py
@@ -544,3 +544,60 @@ def astm_e1036(v, i, imax_limits=(0.75, 1.15), vmax_limits=(0.75, 1.15),
     result['mp_fit'] = mp_fit
 
     return result
+
+
+def _log_lambertw(logx):
+    r'''Computes W(x) starting from log(x).
+
+    Parameters
+    ----------
+    logx : numeric
+        Log(x) of
+
+    Returns
+    -------
+    numeric
+        Lambert's W(x)
+
+    '''
+    # handles overflow cases, but results in nan for x <= 1
+    w = logx - np.log(logx)  # initial guess, w = log(x) - log(log(x))
+
+    for _ in range(0, 3):
+        # Newton's. Halley's is not substantially faster or more accurate
+        # because f''(w) = -1 / (w**2) is small for large w
+        w = w * (1. - np.log(w) + logx) / (1. + w)
+    return w
+
+
+def _lambertw_pvlib(x):
+    r'''Lambert's W function principal branch, :math:`W_0(x)`, for :math:`x`
+    real valued.
+
+    Parameters
+    ----------
+    x : np.array
+        Must be real numbers.
+
+    Returns
+    -------
+    np.array
+
+    '''
+    w = np.full_like(x, np.nan)
+    small = x <= 10
+    # for large x, solve 0 = f(w) = w + log(w) - log(x) using Newton's
+    w[~small] = _log_lambertw(np.log(x[~small]))
+
+    # w will contain nan for these numbers due to log(w) = log(log(x))
+    # for small x, solve 0 = g(w) = w * exp(w) - x using Halley's method
+    if any(small):
+        z = x[small]
+        g = np.log(x[small] + 1) - np.log(np.log(x[small] + 1) + 1)
+        for _ in range(0, 3):
+            expg = np.exp(g)
+            g = g - (g*expg - z) * (g + 1) / \
+                (expg * (g + 1)**2 - 0.5*(g + 2)*(expg*g - z))
+        w[small] = g
+
+    return w
diff --git a/pvlib/tools.py b/pvlib/tools.py
index e0f7cc43ae..6cb631f852 100644
--- a/pvlib/tools.py
+++ b/pvlib/tools.py
@@ -586,60 +586,3 @@ def _file_context_manager(filename_or_object, mode='r', encoding=None):
         # otherwise, assume a filename or path
         context = open(str(filename_or_object), mode=mode, encoding=encoding)
     return context
-
-
-def _log_lambertw(logx):
-    r'''Computes W(x) starting from log(x).
-
-    Parameters
-    ----------
-    logx : numeric
-        Log(x) of
-
-    Returns
-    -------
-    numeric
-        Lambert's W(x)
-
-    '''
-    # handles overflow cases, but results in nan for x <= 1
-    w = logx - np.log(logx)  # initial guess, w = log(x) - log(log(x))
-
-    for _ in range(0, 3):
-        # Newton's. Halley's is not substantially faster or more accurate
-        # because f''(w) = -1 / (w**2) is small for large w
-        w = w * (1. - np.log(w) + logx) / (1. + w)
-    return w
-
-
-def lambertw_pvlib(x):
-    r'''Lambert's W function principal branch, :math:`W_0(x)`, for :math:`x`
-    real valued.
-
-    Parameters
-    ----------
-    x : np.array
-        Must be real numbers.
-
-    Returns
-    -------
-    np.array
-
-    '''
-    w = np.full_like(x, np.nan)
-    small = x <= 10
-    # for large x, solve 0 = f(w) = w + log(w) - log(x) using Newton's
-    w[~small] = _log_lambertw(np.log(x[~small]))
-
-    # w will contain nan for these numbers due to log(w) = log(log(x))
-    # for small x, solve 0 = g(w) = w * exp(w) - x using Halley's method
-    if any(small):
-        z = x[small]
-        g = np.log(x[small] + 1) - np.log(np.log(x[small] + 1) + 1)
-        for _ in range(0, 3):
-            expg = np.exp(g)
-            g = g - (g*expg - z) * (g + 1) / \
-                (expg * (g + 1)**2 - 0.5*(g + 2)*(expg*g - z))
-        w[small] = g
-
-    return w
diff --git a/tests/ivtools/test_utils.py b/tests/ivtools/test_utils.py
index 24c5cd81ad..0e799faaee 100644
--- a/tests/ivtools/test_utils.py
+++ b/tests/ivtools/test_utils.py
@@ -3,6 +3,7 @@
 import pytest
 from pvlib.ivtools.utils import _numdiff, rectify_iv_curve, astm_e1036
 from pvlib.ivtools.utils import _schumaker_qspline
+from pvlib.ivtools.utils import _lambertw_pvlib
 
 from tests.conftest import TESTS_DATA_DIR
 
@@ -171,3 +172,29 @@ def test_astm_e1036_fit_points(v_array, i_array):
                 'ff': 0.7520255886236707}
     result.pop('mp_fit')
     assert result == pytest.approx(expected)
+
+
+def test_lambertw_pvlib():
+    test_exp = np.arange(-10., 300, step=10)
+    test_x = 10.**test_exp
+    # known solution from scipy.special.lambertw
+    # scipy 1.7.1, python 3.13.1, numpy 2.3.5
+    expected = np.array([
+        9.9999999989999997e-11, 5.6714329040978384e-01,
+        2.0028685413304952e+01, 4.2306755091738395e+01,
+        6.4904633770046118e+01, 8.7630277151947183e+01,
+        1.1042491882731335e+02, 1.3326278259180333e+02,
+        1.5613026581351718e+02, 1.7901931374150624e+02,
+        2.0192476320084489e+02, 2.2484310644511851e+02,
+        2.4777185185877809e+02, 2.7070916610249782e+02,
+        2.9365366103997610e+02, 3.1660426041503479e+02,
+        3.3956011295458728e+02, 3.6252053376149752e+02,
+        3.8548496362161768e+02, 4.0845294003314166e+02,
+        4.3142407612718210e+02, 4.5439804503371403e+02,
+        4.7737456808796901e+02, 5.0035340579834485e+02,
+        5.2333435083468805e+02, 5.4631722251791496e+02,
+        5.6930186244110166e+02, 5.9228813095427859e+02,
+        6.1527590431628334e+02, 6.3826507236734335e+02,
+        6.6125553661218726e+02])
+    result = _lambertw_pvlib(test_x)
+    assert np.allclose(result, expected, rtol=1e-14)
diff --git a/tests/test_tools.py b/tests/test_tools.py
index f316cb4605..4b733ad711 100644
--- a/tests/test_tools.py
+++ b/tests/test_tools.py
@@ -283,29 +283,3 @@ def test__file_context_manager():
     buffer = StringIO("test content")
     with tools._file_context_manager(buffer) as obj:
         assert obj.read() == "test content"
-
-
-def test_lambertw_pvlib():
-    test_exp = np.arange(-10., 300, step=10)
-    test_x = 10.**test_exp
-    # known solution from scipy.special.lambertw
-    # scipy 1.7.1, python 3.13.1, numpy 2.3.5
-    expected = np.array([
-        9.9999999989999997e-11, 5.6714329040978384e-01,
-        2.0028685413304952e+01, 4.2306755091738395e+01,
-        6.4904633770046118e+01, 8.7630277151947183e+01,
-        1.1042491882731335e+02, 1.3326278259180333e+02,
-        1.5613026581351718e+02, 1.7901931374150624e+02,
-        2.0192476320084489e+02, 2.2484310644511851e+02,
-        2.4777185185877809e+02, 2.7070916610249782e+02,
-        2.9365366103997610e+02, 3.1660426041503479e+02,
-        3.3956011295458728e+02, 3.6252053376149752e+02,
-        3.8548496362161768e+02, 4.0845294003314166e+02,
-        4.3142407612718210e+02, 4.5439804503371403e+02,
-        4.7737456808796901e+02, 5.0035340579834485e+02,
-        5.2333435083468805e+02, 5.4631722251791496e+02,
-        5.6930186244110166e+02, 5.9228813095427859e+02,
-        6.1527590431628334e+02, 6.3826507236734335e+02,
-        6.6125553661218726e+02])
-    result = tools.lambertw_pvlib(test_x)
-    assert np.allclose(result, expected, rtol=1e-14)

From f0897ba0d054667ac041a5228d6fddc97b3fc722 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Tue, 24 Mar 2026 09:08:34 -0700
Subject: [PATCH 07/12] Apply suggestions from code review

Co-authored-by: Kevin Anderson <kevin.anderso@gmail.com>
---
 docs/sphinx/source/whatsnew/v0.15.1.rst | 4 ++--
 pvlib/ivtools/utils.py                  | 9 ++++++---
 2 files changed, 8 insertions(+), 5 deletions(-)

diff --git a/docs/sphinx/source/whatsnew/v0.15.1.rst b/docs/sphinx/source/whatsnew/v0.15.1.rst
index a4969597dd..0d2a7bac55 100644
--- a/docs/sphinx/source/whatsnew/v0.15.1.rst
+++ b/docs/sphinx/source/whatsnew/v0.15.1.rst
@@ -26,9 +26,9 @@ Enhancements
 * Use ``k`` and ``cap_adjustment`` from :py:func:`pvlib.pvsystem.Array.module_parameters`
   in :py:func:`pvlib.pvsystem.PVSystem.pvwatts_dc`
   (:issue:`2714`, :pull:`2715`)
-* Add :py:func:`pvlib.tools.lambertw_pvlib` to speed up calculations when using
+* Add :py:func:`pvlib.tools._lambertw_pvlib` to speed up calculations when using
   LambertW single diode equation methods.
-  (:discussion:`2720`, :pull:`2723`)  
+  (:discuss:`2720`, :pull:`2723`)  
 
 Documentation
 ~~~~~~~~~~~~~
diff --git a/pvlib/ivtools/utils.py b/pvlib/ivtools/utils.py
index 25ed4f2c84..bda59899b4 100644
--- a/pvlib/ivtools/utils.py
+++ b/pvlib/ivtools/utils.py
@@ -593,11 +593,14 @@ def _lambertw_pvlib(x):
     # for small x, solve 0 = g(w) = w * exp(w) - x using Halley's method
     if any(small):
         z = x[small]
-        g = np.log(x[small] + 1) - np.log(np.log(x[small] + 1) + 1)
+        temp = np.log(x[small] + 1)
+        g = temp - np.log(temp + 1)
         for _ in range(0, 3):
             expg = np.exp(g)
-            g = g - (g*expg - z) * (g + 1) / \
-                (expg * (g + 1)**2 - 0.5*(g + 2)*(expg*g - z))
+            g_expg_z = g*expg - z
+            g_p1 = g + 1
+            g = g - g_expg_z * g_p1 / \
+                (expg * g_p1**2 - 0.5*(g + 2)*g_expg_z)
         w[small] = g
 
     return w

From 610bbb205b1ebd06eb519229bdbb8efb10ee0120 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Tue, 24 Mar 2026 09:39:40 -0700
Subject: [PATCH 08/12] accept float, restrict to x>=0, simplify test

---
 docs/sphinx/source/whatsnew/v0.15.1.rst |  3 ---
 pvlib/ivtools/utils.py                  | 26 ++++++++++-----------
 tests/ivtools/test_utils.py             | 30 ++++++++++---------------
 3 files changed, 25 insertions(+), 34 deletions(-)

diff --git a/docs/sphinx/source/whatsnew/v0.15.1.rst b/docs/sphinx/source/whatsnew/v0.15.1.rst
index 0d2a7bac55..44539d4a5c 100644
--- a/docs/sphinx/source/whatsnew/v0.15.1.rst
+++ b/docs/sphinx/source/whatsnew/v0.15.1.rst
@@ -26,9 +26,6 @@ Enhancements
 * Use ``k`` and ``cap_adjustment`` from :py:func:`pvlib.pvsystem.Array.module_parameters`
   in :py:func:`pvlib.pvsystem.PVSystem.pvwatts_dc`
   (:issue:`2714`, :pull:`2715`)
-* Add :py:func:`pvlib.tools._lambertw_pvlib` to speed up calculations when using
-  LambertW single diode equation methods.
-  (:discuss:`2720`, :pull:`2723`)  
 
 Documentation
 ~~~~~~~~~~~~~
diff --git a/pvlib/ivtools/utils.py b/pvlib/ivtools/utils.py
index bda59899b4..e845ab3e6a 100644
--- a/pvlib/ivtools/utils.py
+++ b/pvlib/ivtools/utils.py
@@ -552,7 +552,6 @@ def _log_lambertw(logx):
     Parameters
     ----------
     logx : numeric
-        Log(x) of
 
     Returns
     -------
@@ -571,29 +570,30 @@ def _log_lambertw(logx):
 
 
 def _lambertw_pvlib(x):
-    r'''Lambert's W function principal branch, :math:`W_0(x)`, for :math:`x`
-    real valued.
+    r'''Lambert's W function principal branch, :math:`W_0(x)`, for
+    :math:`x>=0`.
 
     Parameters
     ----------
-    x : np.array
+    x : float or np.array
         Must be real numbers.
 
     Returns
     -------
-    np.array
+    float or np.array
 
     '''
-    w = np.full_like(x, np.nan)
-    small = x <= 10
+    localx = np.asarray(x, float)
+    w = np.full_like(localx, np.nan)
+    small = localx <= 10
     # for large x, solve 0 = f(w) = w + log(w) - log(x) using Newton's
-    w[~small] = _log_lambertw(np.log(x[~small]))
-
     # w will contain nan for these numbers due to log(w) = log(log(x))
+    w[~small] = _log_lambertw(np.log(localx[~small]))
+
     # for small x, solve 0 = g(w) = w * exp(w) - x using Halley's method
-    if any(small):
-        z = x[small]
-        temp = np.log(x[small] + 1)
+    if np.any(small):
+        z = localx[small]
+        temp = np.log(localx[small] + 1)
         g = temp - np.log(temp + 1)
         for _ in range(0, 3):
             expg = np.exp(g)
@@ -603,4 +603,4 @@ def _lambertw_pvlib(x):
                 (expg * g_p1**2 - 0.5*(g + 2)*g_expg_z)
         w[small] = g
 
-    return w
+    return w[0] if w.shape==1 else w
diff --git a/tests/ivtools/test_utils.py b/tests/ivtools/test_utils.py
index 0e799faaee..c4bb16ae5b 100644
--- a/tests/ivtools/test_utils.py
+++ b/tests/ivtools/test_utils.py
@@ -175,26 +175,20 @@ def test_astm_e1036_fit_points(v_array, i_array):
 
 
 def test_lambertw_pvlib():
-    test_exp = np.arange(-10., 300, step=10)
-    test_x = 10.**test_exp
+    test_x = np.array([0., 1.e-10, 1., 10., 100., 1.e+10, 1.e+100, 1.e+300])
     # known solution from scipy.special.lambertw
     # scipy 1.7.1, python 3.13.1, numpy 2.3.5
     expected = np.array([
-        9.9999999989999997e-11, 5.6714329040978384e-01,
-        2.0028685413304952e+01, 4.2306755091738395e+01,
-        6.4904633770046118e+01, 8.7630277151947183e+01,
-        1.1042491882731335e+02, 1.3326278259180333e+02,
-        1.5613026581351718e+02, 1.7901931374150624e+02,
-        2.0192476320084489e+02, 2.2484310644511851e+02,
-        2.4777185185877809e+02, 2.7070916610249782e+02,
-        2.9365366103997610e+02, 3.1660426041503479e+02,
-        3.3956011295458728e+02, 3.6252053376149752e+02,
-        3.8548496362161768e+02, 4.0845294003314166e+02,
-        4.3142407612718210e+02, 4.5439804503371403e+02,
-        4.7737456808796901e+02, 5.0035340579834485e+02,
-        5.2333435083468805e+02, 5.4631722251791496e+02,
-        5.6930186244110166e+02, 5.9228813095427859e+02,
-        6.1527590431628334e+02, 6.3826507236734335e+02,
-        6.6125553661218726e+02])
+        0.0000000000000000e+00, 9.9999999989999997e-11, 5.6714329040978384e-01,
+        1.7455280027406994e+00, 3.3856301402900502e+00, 2.0028685413304952e+01,
+        2.2484310644511851e+02, 6.8424720862976085e+02])
     result = _lambertw_pvlib(test_x)
     assert np.allclose(result, expected, rtol=1e-14)
+    # with float input
+    for x, k in zip([1.e-10, 1.e+10], [1, 5]):
+        result = _lambertw_pvlib(x)
+        assert np.isclose(result, expected[k])
+    # with 1d array
+    for x, k in zip([1.e-10, 1.e+10], [1, 5]):
+        result = _lambertw_pvlib(np.array([x]))
+        assert np.isclose(result, expected[k])

From b8e1e3dda77824a74d7647ab86360f61cc33dbee Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Tue, 24 Mar 2026 09:41:13 -0700
Subject: [PATCH 09/12] formatting

---
 pvlib/ivtools/utils.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pvlib/ivtools/utils.py b/pvlib/ivtools/utils.py
index e845ab3e6a..969680f2bb 100644
--- a/pvlib/ivtools/utils.py
+++ b/pvlib/ivtools/utils.py
@@ -603,4 +603,4 @@ def _lambertw_pvlib(x):
                 (expg * g_p1**2 - 0.5*(g + 2)*g_expg_z)
         w[small] = g
 
-    return w[0] if w.shape==1 else w
+    return w[0] if w.shape == 1 else w

From 9a0cf8509536ae4357ee90bdaaca192b65b4a0ca Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Tue, 24 Mar 2026 09:42:03 -0700
Subject: [PATCH 10/12] small < 100

---
 pvlib/ivtools/utils.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pvlib/ivtools/utils.py b/pvlib/ivtools/utils.py
index 969680f2bb..4b4bbf83b3 100644
--- a/pvlib/ivtools/utils.py
+++ b/pvlib/ivtools/utils.py
@@ -585,7 +585,7 @@ def _lambertw_pvlib(x):
     '''
     localx = np.asarray(x, float)
     w = np.full_like(localx, np.nan)
-    small = localx <= 10
+    small = localx <= 100
     # for large x, solve 0 = f(w) = w + log(w) - log(x) using Newton's
     # w will contain nan for these numbers due to log(w) = log(log(x))
     w[~small] = _log_lambertw(np.log(localx[~small]))

From eba17772749a0407e5c2b39a88b6398d48962f32 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Tue, 24 Mar 2026 16:29:55 -0700
Subject: [PATCH 11/12] Apply suggestions from code review

Co-authored-by: Echedey Luis <80125792+echedey-ls@users.noreply.github.com>
---
 pvlib/ivtools/utils.py      | 4 ++--
 tests/ivtools/test_utils.py | 4 ++--
 2 files changed, 4 insertions(+), 4 deletions(-)

diff --git a/pvlib/ivtools/utils.py b/pvlib/ivtools/utils.py
index 4b4bbf83b3..3662601f4d 100644
--- a/pvlib/ivtools/utils.py
+++ b/pvlib/ivtools/utils.py
@@ -593,8 +593,8 @@ def _lambertw_pvlib(x):
     # for small x, solve 0 = g(w) = w * exp(w) - x using Halley's method
     if np.any(small):
         z = localx[small]
-        temp = np.log(localx[small] + 1)
-        g = temp - np.log(temp + 1)
+        temp = np.log1p(localx[small])
+        g = temp - np.log1p(temp)
         for _ in range(0, 3):
             expg = np.exp(g)
             g_expg_z = g*expg - z
diff --git a/tests/ivtools/test_utils.py b/tests/ivtools/test_utils.py
index c4bb16ae5b..adbfcaaa51 100644
--- a/tests/ivtools/test_utils.py
+++ b/tests/ivtools/test_utils.py
@@ -185,10 +185,10 @@ def test_lambertw_pvlib():
     result = _lambertw_pvlib(test_x)
     assert np.allclose(result, expected, rtol=1e-14)
     # with float input
-    for x, k in zip([1.e-10, 1.e+10], [1, 5]):
+    for x, k in zip(test_x[[1, 5]], expected[[1, 5]]):
         result = _lambertw_pvlib(x)
         assert np.isclose(result, expected[k])
     # with 1d array
-    for x, k in zip([1.e-10, 1.e+10], [1, 5]):
+    for x, k in zip(test_x[[1, 5]], expected[[1, 5]]):
         result = _lambertw_pvlib(np.array([x]))
         assert np.isclose(result, expected[k])

From 189d8fbc97e1df3267daaca71c477f45d64adee7 Mon Sep 17 00:00:00 2001
From: Cliff Hansen <cwhanse@sandia.gov>
Date: Tue, 24 Mar 2026 16:40:31 -0700
Subject: [PATCH 12/12] fix test loop

---
 tests/ivtools/test_utils.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/tests/ivtools/test_utils.py b/tests/ivtools/test_utils.py
index adbfcaaa51..ae87d9c496 100644
--- a/tests/ivtools/test_utils.py
+++ b/tests/ivtools/test_utils.py
@@ -185,10 +185,10 @@ def test_lambertw_pvlib():
     result = _lambertw_pvlib(test_x)
     assert np.allclose(result, expected, rtol=1e-14)
     # with float input
-    for x, k in zip(test_x[[1, 5]], expected[[1, 5]]):
+    for x, known in zip(test_x[[1, 5]], expected[[1, 5]]):
         result = _lambertw_pvlib(x)
-        assert np.isclose(result, expected[k])
+        assert np.isclose(result, known)
     # with 1d array
-    for x, k in zip(test_x[[1, 5]], expected[[1, 5]]):
+    for x, known in zip(test_x[[1, 5]], expected[[1, 5]]):
         result = _lambertw_pvlib(np.array([x]))
-        assert np.isclose(result, expected[k])
+        assert np.isclose(result, known)