use _lambertw_pvlib etc.

docs/sphinx/source/whatsnew/v0.15.1.rst

@@ -36,6 +36,8 @@ Enhancements
 * Include `ross` and `faiman_rad` in the allowed models within 
   :py:meth:`pvlib.pvsystem.PVSystem.get_cell_temperature` (:issue:`2625`, :pull:`2631`)  
 * Accelerate the internals of :py:func:`~pvlib.solarpostion.ephemeris`. (:pull:`2626`)
+* Accelerate the intervals of :py:func:`~pvlib.pvsystem.singlediode` when
+  `method='lambertw'`. (:pull:`2732`, :pull:`2723`)
 
 Documentation
 ~~~~~~~~~~~~~

pvlib/ivtools/sdm/_fit_desoto_pvsyst_sandia.py

@@ -5,10 +5,9 @@
 import numpy as np
 
 from scipy import optimize
-from scipy.special import lambertw
 
 from pvlib.pvsystem import singlediode, v_from_i
-from pvlib.ivtools.utils import rectify_iv_curve, _numdiff
+from pvlib.ivtools.utils import rectify_iv_curve, _numdiff, _lambertw_pvlib
 from pvlib.pvsystem import _pvsyst_Rsh
 
 
@@ -535,7 +534,7 @@ def _calc_theta_phi_exact(vmp, imp, iph, io, rs, rsh, nnsvth):
             nnsvth == 0,
             np.nan,
             rsh * io / nnsvth * np.exp(rsh * (iph + io - imp) / nnsvth))
-        phi = np.where(argw > 0, lambertw(argw).real, np.nan)
+        phi = np.where(argw > 0, _lambertw_pvlib(argw), np.nan)
 
     # NaN where argw overflows. Switch to log space to evaluate
     u = np.isinf(argw)
@@ -561,7 +560,7 @@ def _calc_theta_phi_exact(vmp, imp, iph, io, rs, rsh, nnsvth):
             np.nan,
             rsh / (rsh + rs) * rs * io / nnsvth * np.exp(
                 rsh / (rsh + rs) * (rs * (iph + io) + vmp) / nnsvth))
-        theta = np.where(argw > 0, lambertw(argw).real, np.nan)
+        theta = np.where(argw > 0, _lambertw_pvlib(argw), np.nan)
 
     # NaN where argw overflows. Switch to log space to evaluate
     u = np.isinf(argw)

pvlib/ivtools/sdm/desoto.py

@@ -2,9 +2,8 @@
 
 from scipy import constants
 from scipy import optimize
-from scipy.special import lambertw
 
-from pvlib.ivtools.utils import rectify_iv_curve
+from pvlib.ivtools.utils import rectify_iv_curve, _lambertw_pvlib
 from pvlib.ivtools.sde import _fit_sandia_cocontent
 
 from pvlib.ivtools.sdm._fit_desoto_pvsyst_sandia import (
@@ -454,7 +453,7 @@ def fit_desoto_batzelis(v_mp, i_mp, v_oc, i_sc, alpha_sc, beta_voc):
     # Equation numbers refer to [1]
     t0 = 298.15  # K
     del0 = (1 - beta_voc * t0) / (50.1 - alpha_sc * t0)  # Eq 9
-    w0 = np.real(lambertw(np.exp(1/del0 + 1)))
+    w0 = _lambertw_pvlib(np.exp(1/del0 + 1))
 
     # Eqs 11-15
     a0 = del0 * v_oc

pvlib/singlediode.py

@@ -5,9 +5,10 @@
 import numpy as np
 import pandas as pd
 from pvlib.tools import _golden_sect_DataFrame
+from pvlib.ivtools.utils import _lambertw_pvlib, _log_lambertw
 
 from scipy.optimize import brentq, newton
-from scipy.special import lambertw
+
 
 # newton method default parameters for this module
 NEWTON_DEFAULT_PARAMS = {
@@ -801,12 +802,10 @@ def _lambertw_v_from_i(current, photocurrent, saturation_current,
         with np.errstate(over='ignore'):
             argW = I0 / (Gsh * a) * np.exp((-I + IL + I0) / (Gsh * a))
 
-        # lambertw typically returns complex value with zero imaginary part
-        # may overflow to np.inf
-        lambertwterm = lambertw(argW).real
+        # Record indices where lambertw input overflowed
+        idx_inf = np.isinf(argW)
 
-        # Record indices where lambertw input overflowed output
-        idx_inf = np.isinf(lambertwterm)
+        lambertwterm = _lambertw_pvlib(argW)
 
         # Only re-compute LambertW if it overflowed
         if np.any(idx_inf):
@@ -815,15 +814,7 @@ def _lambertw_v_from_i(current, photocurrent, saturation_current,
                        np.log(a[idx_inf]) +
                        (-I[idx_inf] + IL[idx_inf] + I0[idx_inf]) /
                        (Gsh[idx_inf] * a[idx_inf]))
-
-            # Three iterations of Newton-Raphson method to solve
-            #  w+log(w)=logargW. The initial guess is w=logargW. Where direct
-            #  evaluation (above) results in NaN from overflow, 3 iterations
-            #  of Newton's method gives approximately 8 digits of precision.
-            w = logargW
-            for _ in range(0, 3):
-                w = w * (1. - np.log(w) + logargW) / (1. + w)
-            lambertwterm[idx_inf] = w
+            lambertwterm[idx_inf] = _log_lambertw(logargW)
 
         # Eqn. 3 in Jain and Kapoor, 2004
         #  V = -I*(Rs + Rsh) + IL*Rsh - a*lambertwterm + I0*Rsh
@@ -879,7 +870,7 @@ def _lambertw_i_from_v(voltage, photocurrent, saturation_current,
 
         # lambertw typically returns complex value with zero imaginary part
         # may overflow to np.inf
-        lambertwterm = lambertw(argW).real
+        lambertwterm = _lambertw_pvlib(argW)
 
         # Eqn. 2 in Jain and Kapoor, 2004
         #  I = -V/(Rs + Rsh) - (a/Rs)*lambertwterm + Rsh*(IL + I0)/(Rs + Rsh)
@@ -1031,7 +1022,7 @@ def batzelis(photocurrent, saturation_current, resistance_series,
         voc = a * np.log(Iph / Is)
 
     # Eqs 5-8
-    w = np.real(lambertw(np.e * Iph / Is))
+    w = _lambertw_pvlib(np.e * Iph / Is)
     # vmp = (1 + Rs/Rsh) * a * (w - 1) - Rs * Iph * (1 - 1/w)  # not needed
     with np.errstate(divide='ignore', invalid='ignore'):  # zero Iph -> zero w
         imp = Iph * (1 - 1/w) - a * (w - 1) / Rsh