(f, symbol, gens)
return list(ordered(soln))
Return solution to ``f`` if it is a Lambert-type expression
    else raise NotImplementedError.

    The equality, ``f(x, a..f) = a*log(b*X + c) + d*X - f = 0`` has the
    solution,  `X = -c/b + (a/d)*W(d/(a*b)*exp(c*d/a/b)*exp(f/a))`. There
    are a variety of forms for `f(X, a..f)` as enumerated below:

    1a1)
      if B**B = R for R not [0, 1] then
      log(B) + log(log(B)) = log(log(R))
      X = log(B), a = 1, b = 1, c = 0, d = 1, f = log(log(R))
    1a2)
      if B*(b*log(B) + c)**a = R then
      log(B) + a*log(b*log(B) + c) = log(R)
      X = log(B); d=1, f=log(R)
    1b)
      if a*log(b*B + c) + d*B = R then
      X = B, f = R
    2a)
      if (b*B + c)*exp(d*B + g) = R then
      log(b*B + c) + d*B + g = log(R)
      a = 1, f = log(R) - g, X = B
    2b)
      if -b*B + g*exp(d*B + h) = c then
      log(g) + d*B + h - log(b*B + c) = 0
      a = -1, f = -h - log(g), X = B
    3)
      if d*p**(a*B + g) - b*B = c then
      log(d) + (a*B + g)*log(p) - log(c + b*B) = 0
      a = -1, d = a*log(p), f = -log(d) - g*log(p)
Return solution to ``f`` if it is a Lambert-type expression
    else raise NotImplementedError.
Return solution to f if it is a Lambert - type expression else raise NotImplementedError .
def _solve_lambert(f, symbol, gens):
    """Return solution to ``f`` if it is a Lambert-type expression
    else raise NotImplementedError.

    The equality, ``f(x, a..f) = a*log(b*X + c) + d*X - f = 0`` has the
    solution,  `X = -c/b + (a/d)*W(d/(a*b)*exp(c*d/a/b)*exp(f/a))`. There
    are a variety of forms for `f(X, a..f)` as enumerated below:

    1a1)
      if B**B = R for R not [0, 1] then
      log(B) + log(log(B)) = log(log(R))
      X = log(B), a = 1, b = 1, c = 0, d = 1, f = log(log(R))
    1a2)
      if B*(b*log(B) + c)**a = R then
      log(B) + a*log(b*log(B) + c) = log(R)
      X = log(B); d=1, f=log(R)
    1b)
      if a*log(b*B + c) + d*B = R then
      X = B, f = R
    2a)
      if (b*B + c)*exp(d*B + g) = R then
      log(b*B + c) + d*B + g = log(R)
      a = 1, f = log(R) - g, X = B
    2b)
      if -b*B + g*exp(d*B + h) = c then
      log(g) + d*B + h - log(b*B + c) = 0
      a = -1, f = -h - log(g), X = B
    3)
      if d*p**(a*B + g) - b*B = c then
      log(d) + (a*B + g)*log(p) - log(c + b*B) = 0
      a = -1, d = a*log(p), f = -log(d) - g*log(p)
    """

    nrhs, lhs = f.as_independent(symbol, as_Add=True)
    rhs = -nrhs

    lamcheck = [tmp for tmp in gens
                if (tmp.func in [exp, log] or
                (tmp.is_Pow and symbol in tmp.exp.free_symbols))]
    if not lamcheck:
        raise NotImplementedError()

    if lhs.is_Mul:
        lhs = expand_log(log(lhs))
        rhs = log(rhs)

    lhs = factor(lhs, deep=True)
    # make sure we are inverted as completely as possible
    r = Dummy()
    i, lhs = _invert(lhs - r, symbol)
    rhs = i.xreplace({r: rhs})

    # For the first ones:
    # 1a1) B**B = R != 0 (when 0, there is only a solution if the base is 0,
    #                     but if it is, the exp is 0 and 0**0=1
    #                     comes back as B*log(B) = log(R)
    # 1a2) B*(a + b*log(B))**p = R or with monomial expanded or with whole
    #                              thing expanded comes back unchanged
    #     log(B) + p*log(a + b*log(B)) = log(R)
    #     lhs is Mul:
    #         expand log of both sides to give:
    #         log(B) + log(log(B)) = log(log(R))
    # 1b) d*log(a*B + b) + c*B = R
    #     lhs is Add:
    #         isolate c*B and expand log of both sides:
    #         log(c) + log(B) = log(R - d*log(a*B + b))

    soln = []
    if not soln:
        mainlog = _mostfunc(lhs, log, symbol)
        if mainlog:
            if lhs.is_Mul and rhs != 0:
                soln = _lambert(log(lhs) - log(rhs), symbol)
            elif lhs.is_Add:
                other = lhs.subs(mainlog, 0)
                if other and not other.is_Add and [
                        tmp for tmp in other.atoms(Pow)
                        if symbol in tmp.free_symbols]:
                    if not rhs:
                        diff = log(other) - log(other - lhs)
                    else:
                        diff = log(lhs - other) - log(rhs - other)
                    soln = _lambert(expand_log(diff), symbol)
                else:
                    #it's ready to go
                    soln = _lambert(lhs - rhs, symbol)

    # For the next two,
    #     collect on main exp
    #     2a) (b*B + c)*exp(d*B + g) = R
    #         lhs is mul:
    #             log to give
    #             log(b*B + c) + d*B = log(R) - g
    #     2b) -b*B + g*exp(d*B + h) = R
    #         lhs is add:
    #             add b*B
    #             log and rearrange
    #             log(R + b*B) - d*B = log(g) + h

    if not soln:
        mainexp = _mostfunc(lhs, exp, symbol)
        if mainexp:
            lhs = collect(lhs, mainexp)
            if lhs.is_Mul and rhs != 0:
                soln = _lambert(expand_log(log(lhs) - log(rhs)), symbol)
            elif lhs.is_Add:
                # move all but mainexp-containing term to rhs
                other = lhs.subs(mainexp, 0)
                mainterm = lhs - other
                rhs=rhs - other
                if (mainterm.could_extract_minus_sign() and
                    rhs.could_extract_minus_sign()):
                    mainterm *= -1
                    rhs *= -1
                diff = log(mainterm) - log(rhs)
                soln = _lambert(expand_log(diff), symbol)

    # 3) d*p**(a*B + b) + c*B = R
    #     collect on main pow
    #     log(R - c*B) - a*B*log(p) = log(d) + b*log(p)

    if not soln:
        mainpow = _mostfunc(lhs, Pow, symbol)
        if mainpow and symbol in mainpow.exp.free_symbols:
            lhs = collect(lhs, mainpow)
            if lhs.is_Mul and rhs != 0:
                soln = _lambert(expand_log(log(lhs) - log(rhs)), symbol)
            elif lhs.is_Add:
                # move all but mainpow-containing term to rhs
                other = lhs.subs(mainpow, 0)
                mainterm = lhs - other
                rhs = rhs - other
                diff = log(mainterm) - log(rhs)
                soln = _lambert(expand_log(diff), symbol)

    if not soln:
        raise NotImplementedError('%s does not appear to have a solution in '
            'terms of LambertW' % f)

    return list(ordered(soln))
def _solve_lambert ( f , symbol , gens ) : nrhs , lhs = f . as_independent ( symbol , as_Add = True ) rhs = - nrhs lamcheck = [ tmp for tmp in gens if ( tmp . func in [ exp , log ] or ( tmp . is_Pow and symbol in tmp . exp . free_symbols ) ) ] if not lamcheck : raise NotImplementedError ( ) if lhs . is_Mul : lhs = expand_log ( log ( lhs ) ) rhs = log ( rhs ) lhs = factor ( lhs , deep = True ) # make sure we are inverted as completely as possible r = Dummy ( ) i , lhs = _invert ( lhs - r , symbol ) rhs = i . xreplace ( { r : rhs } ) # For the first ones: # 1a1) B**B = R != 0 (when 0, there is only a solution if the base is 0, #                     but if it is, the exp is 0 and 0**0=1 #                     comes back as B*log(B) = log(R) # 1a2) B*(a + b*log(B))**p = R or with monomial expanded or with whole #                              thing expanded comes back unchanged #     log(B) + p*log(a + b*log(B)) = log(R) #     lhs is Mul: #         expand log of both sides to give: #         log(B) + log(log(B)) = log(log(R)) # 1b) d*log(a*B + b) + c*B = R #     lhs is Add: #         isolate c*B and expand log of both sides: #         log(c) + log(B) = log(R - d*log(a*B + b)) soln = [ ] if not soln : mainlog = _mostfunc ( lhs , log , symbol ) if mainlog : if lhs . is_Mul and rhs != 0 : soln = _lambert ( log ( lhs ) - log ( rhs ) , symbol ) elif lhs . is_Add : other = lhs . subs ( mainlog , 0 ) if other and not other . is_Add and [ tmp for tmp in other . atoms ( Pow ) if symbol in tmp . free_symbols ] : if not rhs : diff = log ( other ) - log ( other - lhs ) else : diff = log ( lhs - other ) - log ( rhs - other ) soln = _lambert ( expand_log ( diff ) , symbol ) else : #it's ready to go soln = _lambert ( lhs - rhs , symbol ) # For the next two, #     collect on main exp #     2a) (b*B + c)*exp(d*B + g) = R #         lhs is mul: #             log to give #             log(b*B + c) + d*B = log(R) - g #     2b) -b*B + g*exp(d*B + h) = R #         lhs is add: #             add b*B #             log and rearrange #             log(R + b*B) - d*B = log(g) + h if not soln : mainexp = _mostfunc ( lhs , exp , symbol ) if mainexp : lhs = collect ( lhs , mainexp ) if lhs . is_Mul and rhs != 0 : soln = _lambert ( expand_log ( log ( lhs ) - log ( rhs ) ) , symbol ) elif lhs . is_Add : # move all but mainexp-containing term to rhs other = lhs . subs ( mainexp , 0 ) mainterm = lhs - other rhs = rhs - other if ( mainterm . could_extract_minus_sign ( ) and rhs . could_extract_minus_sign ( ) ) : mainterm *= - 1 rhs *= - 1 diff = log ( mainterm ) - log ( rhs ) soln = _lambert ( expand_log ( diff ) , symbol ) # 3) d*p**(a*B + b) + c*B = R #     collect on main pow #     log(R - c*B) - a*B*log(p) = log(d) + b*log(p) if not soln : mainpow = _mostfunc ( lhs , Pow , symbol ) if mainpow and symbol in mainpow . exp . free_symbols : lhs = collect ( lhs , mainpow ) if lhs . is_Mul and rhs != 0 : soln = _lambert ( expand_log ( log ( lhs ) - log ( rhs ) ) , symbol ) elif lhs . is_Add : # move all but mainpow-containing term to rhs other = lhs . subs ( mainpow , 0 ) mainterm = lhs - other rhs = rhs - other diff = log ( mainterm ) - log ( rhs ) soln = _lambert ( expand_log ( diff ) , symbol ) if not soln : raise NotImplementedError ( '%s does not appear to have a solution in ' 'terms of LambertW' % f ) return list ( ordered ( soln ) )