text: #
  # This is the module for ODE solver classes for single ODEs.
  #
  import typing
  
  if typing.TYPE_CHECKING:
      from typing import ClassVar
  from typing import Dict, Iterable, List, Optional, Type
  
  from sympy.core import S
  from sympy.core.exprtools import factor_terms
  from sympy.core.expr import Expr
  from sympy.core.function import AppliedUndef, Derivative, Function, expand, Subs
  from sympy.core.relational import Equality, Eq
  from sympy.core.symbol import Symbol, Dummy, Wild
  from sympy.functions import exp, sqrt, tan, log
  from sympy.integrals import Integral
  from sympy.polys.polytools import cancel, factor, factor_list
  from sympy.simplify import simplify
  from sympy.simplify.radsimp import fraction
  from sympy.utilities import numbered_symbols
  
  from sympy.solvers.solvers import solve
  from sympy.solvers.deutils import ode_order, _preprocess
  
  
  class ODEMatchError(NotImplementedError):
      """Raised if a SingleODESolver is asked to solve an ODE it does not match"""
      pass
  
  
  def cached_property(func):
      '''Decorator to cache property method'''
      attrname = '_' + func.__name__
      def propfunc(self):
          val = getattr(self, attrname, None)
          if val is None:
              val = func(self)
              setattr(self, attrname, val)
          return val
      return property(propfunc)
  
  
  class SingleODEProblem:
      """Represents an ordinary differential equation (ODE)
  
      This class is used internally in the by dsolve and related
      functions/classes so that properties of an ODE can be computed
      efficiently.
  
      Examples
      ========
  
      This class is used internally by dsolve. To instantiate an instance
      directly first define an ODE problem:
  
      >>> from sympy import Eq, Function, Symbol
      >>> x = Symbol('x')
      >>> f = Function('f')
      >>> eq = f(x).diff(x, 2)
  
      Now you can create a SingleODEProblem instance and query its properties:
  
      >>> from sympy.solvers.ode.single import SingleODEProblem
      >>> problem = SingleODEProblem(f(x).diff(x), f(x), x)
      >>> problem.eq
      Derivative(f(x), x)
      >>> problem.func
      f(x)
      >>> problem.sym
      x
      """
  
      # Instance attributes:
      eq = None  # type: Expr
      func = None  # type: AppliedUndef
      sym = None  # type: Symbol
      _order = None  # type: int
      _eq_expanded = None  # type: Expr
      _eq_preprocessed = None  # type: Expr
  
      def __init__(self, eq, func, sym, prep=True):
          assert isinstance(eq, Expr)
          assert isinstance(func, AppliedUndef)
          assert isinstance(sym, Symbol)
          assert isinstance(prep, bool)
          self.eq = eq
          self.func = func
          self.sym = sym
          self.prep = prep
  
      @cached_property
      def order(self) -> int:
          return ode_order(self.eq, self.func)
  
      @cached_property
      def eq_preprocessed(self) -> Expr:
          return self._get_eq_preprocessed()
  
      @cached_property
      def eq_expanded(self) -> Expr:
          return expand(self.eq_preprocessed)
  
      def _get_eq_preprocessed(self) -> Expr:
          if self.prep:
              process_eq, process_func = _preprocess(self.eq, self.func)
              if process_func != self.func:
                  raise ValueError
          else:
              process_eq = self.eq
          return process_eq
  
      def get_numbered_constants(self, num=1, start=1, prefix='C') -> List[Symbol]:
          """
          Returns a list of constants that do not occur
          in eq already.
          """
          ncs = self.iter_numbered_constants(start, prefix)
          Cs = [next(ncs) for i in range(num)]
          return Cs
  
      def iter_numbered_constants(self, start=1, prefix='C') -> Iterable[Symbol]:
          """
          Returns an iterator of constants that do not occur
          in eq already.
          """
          atom_set = self.eq.free_symbols
          func_set = self.eq.atoms(Function)
          if func_set:
              atom_set |= {Symbol(str(f.func)) for f in func_set}
          return numbered_symbols(start=start, prefix=prefix, exclude=atom_set)
  
      # TODO: Add methods that can be used by many ODE solvers:
      # order
      # is_linear()
      # get_linear_coefficients()
      # eq_prepared (the ODE in prepared form)
  
  
  class SingleODESolver:
      """
      Base class for Single ODE solvers.
  
      Subclasses should implement the _matches and _get_general_solution
      methods. This class is not intended to be instantiated directly but its
      subclasses are as part of dsolve.
  
      Examples
      ========
  
      You can use a subclass of SingleODEProblem to solve a particular type of
      ODE. We first define a particular ODE problem:
  
      >>> from sympy import Eq, Function, Symbol
      >>> x = Symbol('x')
      >>> f = Function('f')
      >>> eq = f(x).diff(x, 2)
  
      Now we solve this problem using the NthAlgebraic solver which is a
      subclass of SingleODESolver:
  
      >>> from sympy.solvers.ode.single import NthAlgebraic, SingleODEProblem
      >>> problem = SingleODEProblem(eq, f(x), x)
      >>> solver = NthAlgebraic(problem)
      >>> solver.get_general_solution()
      [Eq(f(x), _C*x + _C)]
  
      The normal way to solve an ODE is to use dsolve (which would use
      NthAlgebraic and other solvers internally). When using dsolve a number of
      other things are done such as evaluating integrals, simplifying the
      solution and renumbering the constants:
  
      >>> from sympy import dsolve
      >>> dsolve(eq, hint='nth_algebraic')
      Eq(f(x), C1 + C2*x)
      """
  
      # Subclasses should store the hint name (the argument to dsolve) in this
      # attribute
      hint = None  # type: ClassVar[str]
  
      # Subclasses should define this to indicate if they support an _Integral
      # hint.
      has_integral = None  # type: ClassVar[bool]
  
      # The ODE to be solved
      ode_problem = None  # type: SingleODEProblem
  
      # Cache whether or not the equation has matched the method
      _matched = None  # type: Optional[bool]
  
      # Subclasses should store in this attribute the list of order(s) of ODE
      # that subclass can solve or leave it to None if not specific to any order
      order = None  # type: None or list
  
      def __init__(self, ode_problem):
          self.ode_problem = ode_problem
  
      def matches(self) -> bool:
          if self.order is not None and self.ode_problem.order not in self.order:
              self._matched = False
              return self._matched
  
          if self._matched is None:
              self._matched = self._matches()
          return self._matched
  
      def get_general_solution(self, *, simplify: bool = True) -> List[Equality]:
          if not self.matches():
              msg = "%s solver can not solve:\n%s"
              raise ODEMatchError(msg % (self.hint, self.ode_problem.eq))
          return self._get_general_solution()
  
      def _matches(self) -> bool:
          msg = "Subclasses of SingleODESolver should implement matches."
          raise NotImplementedError(msg)
  
      def _get_general_solution(self, *, simplify: bool = True) -> List[Equality]:
          msg = "Subclasses of SingleODESolver should implement get_general_solution."
          raise NotImplementedError(msg)
  
  
  class SinglePatternODESolver(SingleODESolver):
      '''Superclass for ODE solvers based on pattern matching'''
  
      def wilds(self):
          prob = self.ode_problem
          f = prob.func.func
          x = prob.sym
          order = prob.order
          return self._wilds(f, x, order)
  
      def wilds_match(self):
          match = self._wilds_match
          return [match.get(w, S.Zero) for w in self.wilds()]
  
      def _matches(self):
          eq = self.ode_problem.eq_expanded
          f = self.ode_problem.func.func
          x = self.ode_problem.sym
          order = self.ode_problem.order
          df = f(x).diff(x)
  
          if order != 1:
              return False
  
          pattern = self._equation(f(x), x, 1)
  
          if not pattern.coeff(df).has(Wild):
              eq = expand(eq / eq.coeff(df))
          eq = eq.collect([f(x).diff(x)], func = cancel)
  
          self._wilds_match = match = eq.match(pattern)
          if match is not None:
              return self._verify(f(x))
          return False
  
      def _verify(self, fx) -> bool:
          return True
  
      def _wilds(self, f, x, order):
          msg = "Subclasses of SingleODESolver should implement _wilds"
          raise NotImplementedError(msg)
  
      def _equation(self, fx, x, order):
          msg = "Subclasses of SingleODESolver should implement _equation"
          raise NotImplementedError(msg)
  
  
  class NthAlgebraic(SingleODESolver):
      r"""
      Solves an `n`\th order ordinary differential equation using algebra and
      integrals.
  
      There is no general form for the kind of equation that this can solve. The
      the equation is solved algebraically treating differentiation as an
      invertible algebraic function.
  
      Examples
      ========
  
      >>> from sympy import Function, dsolve, Eq
      >>> from sympy.abc import x
      >>> f = Function('f')
      >>> eq = Eq(f(x) * (f(x).diff(x)**2 - 1), 0)
      >>> dsolve(eq, f(x), hint='nth_algebraic')
      [Eq(f(x), 0), Eq(f(x), C1 - x), Eq(f(x), C1 + x)]
  
      Note that this solver can return algebraic solutions that do not have any
      integration constants (f(x) = 0 in the above example).
      """
  
      hint = 'nth_algebraic'
      has_integral = True  # nth_algebraic_Integral hint
  
      def _matches(self):
          r"""
          Matches any differential equation that nth_algebraic can solve. Uses
          `sympy.solve` but teaches it how to integrate derivatives.
  
          This involves calling `sympy.solve` and does most of the work of finding a
          solution (apart from evaluating the integrals).
          """
          eq = self.ode_problem.eq
          func = self.ode_problem.func
          var = self.ode_problem.sym
  
          # Derivative that solve can handle:
          diffx = self._get_diffx(var)
  
          # Replace derivatives wrt the independent variable with diffx
          def replace(eq, var):
              def expand_diffx(*args):
                  differand, diffs = args[0], args[1:]
                  toreplace = differand
                  for v, n in diffs:
                      for _ in range(n):
                          if v == var:
                              toreplace = diffx(toreplace)
                          else:
                              toreplace = Derivative(toreplace, v)
                  return toreplace
              return eq.replace(Derivative, expand_diffx)
  
          # Restore derivatives in solution afterwards
          def unreplace(eq, var):
              return eq.replace(diffx, lambda e: Derivative(e, var))
  
          subs_eqn = replace(eq, var)
          try:
              # turn off simplification to protect Integrals that have
              # _t instead of fx in them and would otherwise factor
              # as t_*Integral(1, x)
              solns = solve(subs_eqn, func, simplify=False)
          except NotImplementedError:
              solns = []
  
          solns = [simplify(unreplace(soln, var)) for soln in solns]
          solns = [Equality(func, soln) for soln in solns]
  
          self.solutions = solns
          return len(solns) != 0
  
      def _get_general_solution(self, *, simplify: bool = True):
          return self.solutions
  
      # This needs to produce an invertible function but the inverse depends
      # which variable we are integrating with respect to. Since the class can
      # be stored in cached results we need to ensure that we always get the
      # same class back for each particular integration variable so we store these
      # classes in a global dict:
      _diffx_stored = {}  # type: Dict[Symbol, Type[Function]]
  
      @staticmethod
      def _get_diffx(var):
          diffcls = NthAlgebraic._diffx_stored.get(var, None)
  
          if diffcls is None:
              # A class that behaves like Derivative wrt var but is "invertible".
              class diffx(Function):
                  def inverse(self):
                      # don't use integrate here because fx has been replaced by _t
                      # in the equation; integrals will not be correct while solve
                      # is at work.
                      return lambda expr: Integral(expr, var) + Dummy('C')
  
              diffcls = NthAlgebraic._diffx_stored.setdefault(var, diffx)
  
          return diffcls
  
  
  class FirstExact(SinglePatternODESolver):
      r"""
      Solves 1st order exact ordinary differential equations.
  
      A 1st order differential equation is called exact if it is the total
      differential of a function. That is, the differential equation
  
      .. math:: P(x, y) \,\partial{}x + Q(x, y) \,\partial{}y = 0
  
      is exact if there is some function `F(x, y)` such that `P(x, y) =
      \partial{}F/\partial{}x` and `Q(x, y) = \partial{}F/\partial{}y`.  It can
      be shown that a necessary and sufficient condition for a first order ODE
      to be exact is that `\partial{}P/\partial{}y = \partial{}Q/\partial{}x`.
      Then, the solution will be as given below::
  
          >>> from sympy import Function, Eq, Integral, symbols, pprint
          >>> x, y, t, x0, y0, C1= symbols('x,y,t,x0,y0,C1')
          >>> P, Q, F= map(Function, ['P', 'Q', 'F'])
          >>> pprint(Eq(Eq(F(x, y), Integral(P(t, y), (t, x0, x)) +
          ... Integral(Q(x0, t), (t, y0, y))), C1))
                      x                y
                      /                /
                     |                |
          F(x, y) =  |  P(t, y) dt +  |  Q(x0, t) dt = C1
                     |                |
                    /                /
                    x0               y0
  
      Where the first partials of `P` and `Q` exist and are continuous in a
      simply connected region.
  
      A note: SymPy currently has no way to represent inert substitution on an
      expression, so the hint ``1st_exact_Integral`` will return an integral
      with `dy`.  This is supposed to represent the function that you are
      solving for.
  
      Examples
      ========
      >>> from sympy import Symbol,Function
      >>> from sympy import E,dsolve,Eq,log,sin,cos
      >>> x=Symbol('x')
      >>> f=Function('f')
      >>> eq=Eq(0,(sin(x)*sin(f(x))-x*E**f(x))*f(x).diff(x)-(E**f(x)+cos(x)*cos(f(x))))
      >>> dsolve(eq,hint='1st_exact')
      Eq(x*exp(f(x)) + sin(x)*cos(f(x)), C1)
      >>> dsolve(cos(f(x)) - (x*sin(f(x)) - f(x)**2)*f(x).diff(x),f(x), hint='1st_exact')
      Eq(x*cos(f(x)) + f(x)**3/3, C1)
  
      References
      ==========
  
      - "Differential Equations With Applications And Historical Notes",
         Second Edition, George.F.Simmons
      - https://en.wikipedia.org/wiki/Exact_differential_equation
  
      # indirect doctest
  
      """
      hint = "1st_exact"
      has_integral = True
      order = [1]
  
      def _wilds(self, f, x, order):
          P = Wild('P', exclude=[f(x).diff(x)])
          Q = Wild('Q', exclude=[f(x).diff(x)])
          return P, Q
  
      def _equation(self, fx, x, order):
          P, Q = self.wilds()
          return P + Q*fx.diff(x)
  
      def _verify(self, fx) -> bool:
          P, Q = self.wilds()
          x = self.ode_problem.sym
          y = Dummy('y')
  
          m, n = self.wilds_match()
  
          try:
              if simplify(m) != 0:
                  m = m.subs(fx,y)
                  n = n.subs(fx,y)
                  numerator= simplify(m.diff(y)-n.diff(x))
                  # The following few conditions try to convert a non-exact
                  # differential equation into an exact one.
                  # References :
                  #1. Differential equations with applications
                  # and historical notes - George E. Simmons
                  #2. https://math.okstate.edu/people/binegar/2233-S99/2233-l12.pdf
                  if numerator:
                      # If (dP/dy - dQ/dx) / Q = f(x)
                      # then exp(integral(f(x))*equation becomes exact
                      factor = simplify(numerator/n)
                      variables = factor.free_symbols
                      if len(variables) == 1 and x == variables.pop():
                          factor = exp(Integral(factor).doit())
                          m *= factor
                          n *= factor
                          self._wilds_match[P] = m.subs(y, fx)
                          self._wilds_match[Q] = n.subs(y, fx)
                          return True
                      else:
                          # If (dP/dy - dQ/dx) / -P = f(y)
                          # then exp(integral(f(y))*equation becomes exact
                          factor = simplify(-numerator / m)
                          variables = factor.free_symbols
                          if len(variables) == 1 and y == variables.pop():
                              factor = exp(Integral(factor).doit())
                              m *= factor
                              n *= factor
                              self._wilds_match[P] = m.subs(y, fx)
                              self._wilds_match[Q] = n.subs(y, fx)
                              return True
                  else:
                      return True
          except NotImplementedError:
              # Differentiating the coefficients might fail because of things
              # like f(2*x).diff(x).  See issue 4624 and issue 4719.
              pass
  
          return False
  
      def _get_general_solution(self, *, simplify: bool = True):
          m, n = self.wilds_match()
          fx = self.ode_problem.func
          x = self.ode_problem.sym
          (C1,) = self.ode_problem.get_numbered_constants(num=1)
          y = Dummy('y')
  
          m = m.subs(fx, y)
          n = n.subs(fx, y)
  
          gen_sol = Eq(Subs(Integral(m, x)
                            + Integral(n - Integral(m, x).diff(y), y), y, fx), C1)
          return gen_sol
  
  
  class FirstLinear(SinglePatternODESolver):
      r"""
      Solves 1st order linear differential equations.
  
      These are differential equations of the form
  
      .. math:: dy/dx + P(x) y = Q(x)\text{.}
  
      These kinds of differential equations can be solved in a general way.  The
      integrating factor `e^{\int P(x) \,dx}` will turn the equation into a
      separable equation.  The general solution is::
  
          >>> from sympy import Function, dsolve, Eq, pprint, diff, sin
          >>> from sympy.abc import x
          >>> f, P, Q = map(Function, ['f', 'P', 'Q'])
          >>> genform = Eq(f(x).diff(x) + P(x)*f(x), Q(x))
          >>> pprint(genform)
                      d
          P(x)*f(x) + --(f(x)) = Q(x)
                      dx
          >>> pprint(dsolve(genform, f(x), hint='1st_linear_Integral'))
                  /       /                   \
                  |      |                    |
                  |      |         /          |     /
                  |      |        |           |    |
                  |      |        | P(x) dx   |  - | P(x) dx
                  |      |        |           |    |
                  |      |       /            |   /
          f(x) = |C1 +  | Q(x)*e           dx|*e
                  |      |                    |
                  \     /                     /
  
  
      Examples
      ========
  
      >>> f = Function('f')
      >>> pprint(dsolve(Eq(x*diff(f(x), x) - f(x), x**2*sin(x)),
      ... f(x), '1st_linear'))
      f(x) = x*(C1 - cos(x))
  
      References
      ==========
  
      - https://en.wikipedia.org/wiki/Linear_differential_equation#First_order_equation
      - M. Tenenbaum & H. Pollard, "Ordinary Differential Equations",
        Dover 1963, pp. 92
  
      # indirect doctest
  
      """
      hint = '1st_linear'
      has_integral = True
      order = [1]
  
      def _wilds(self, f, x, order):
          P = Wild('P', exclude=[f(x)])
          Q = Wild('Q', exclude=[f(x), f(x).diff(x)])
          return P, Q
  
      def _equation(self, fx, x, order):
          P, Q = self.wilds()
          return fx.diff(x) + P*fx - Q
  
      def _get_general_solution(self, *, simplify: bool = True):
          P, Q = self.wilds_match()
          fx = self.ode_problem.func
          x = self.ode_problem.sym
          (C1,)  = self.ode_problem.get_numbered_constants(num=1)
          gensol = Eq(fx, (((C1 + Integral(Q*exp(Integral(P, x)),x))
              * exp(-Integral(P, x)))))
          return [gensol]
  
  
  class AlmostLinear(SinglePatternODESolver):
      r"""
      Solves an almost-linear differential equation.
  
      The general form of an almost linear differential equation is
  
      .. math:: a(x) g'(f(x)) f'(x) + b(x) g(f(x)) + c(x)
  
      Here `f(x)` is the function to be solved for (the dependent variable).
      The substitution `g(f(x)) = u(x)` leads to a linear differential equation
      for `u(x)` of the form `a(x) u' + b(x) u + c(x) = 0`. This can be solved
      for `u(x)` by the `first_linear` hint and then `f(x)` is found by solving
      `g(f(x)) = u(x)`.
  
      See Also
      ========
      :meth:`sympy.solvers.ode.single.FirstLinear`
  
      Examples
      ========
  
      >>> from sympy import Function, Derivative, pprint, sin, cos
      >>> from sympy.solvers.ode import dsolve, classify_ode
      >>> from sympy.abc import x
      >>> f = Function('f')
      >>> d = f(x).diff(x)
      >>> eq = x*d + x*f(x) + 1
      >>> dsolve(eq, f(x), hint='almost_linear')
      Eq(f(x), (C1 - Ei(x))*exp(-x))
      >>> pprint(dsolve(eq, f(x), hint='almost_linear'))
                          -x
      f(x) = (C1 - Ei(x))*e
      >>> example = cos(f(x))*f(x).diff(x) + sin(f(x)) + 1
      >>> pprint(example)
                          d
      sin(f(x)) + cos(f(x))*--(f(x)) + 1
                          dx
      >>> pprint(dsolve(example, f(x), hint='almost_linear'))
                      /    -x    \             /    -x    \
      [f(x) = pi - asin\C1*e   - 1/, f(x) = asin\C1*e   - 1/]
  
  
      References
      ==========
  
      - Joel Moses, "Symbolic Integration - The Stormy Decade", Communications
        of the ACM, Volume 14, Number 8, August 1971, pp. 558
      """
      hint = "almost_linear"
      has_integral = True
      order = [1]
  
      def _wilds(self, f, x, order):
          P = Wild('P', exclude=[f(x).diff(x)])
          Q = Wild('Q', exclude=[f(x).diff(x)])
          return P, Q
  
      def _equation(self, fx, x, order):
          P, Q = self.wilds()
          return P*fx.diff(x) + Q
  
      def _verify(self, fx):
          a, b = self.wilds_match()
          c, b = b.as_independent(fx) if b.is_Add else (S.Zero, b)
          # a, b and c are the function a(x), b(x) and c(x) respectively.
          # c(x) is obtained by separating out b as terms with and without fx i.e, l(y)
          # The following conditions checks if the given equation is an almost-linear differential equation using the fact that
          # a(x)*(l(y))' / l(y)' is independent of l(y)
  
          if b.diff(fx) != 0 and not simplify(b.diff(fx)/a).has(fx):
              self.ly = factor_terms(b).as_independent(fx, as_Add=False)[1] # Gives the term containing fx i.e., l(y)
              self.ax = a / self.ly.diff(fx)
              self.cx = -c  # cx is taken as -c(x) to simplify expression in the solution integral
              self.bx = factor_terms(b) / self.ly
              return True
  
          return False
  
      def _get_general_solution(self, *, simplify: bool = True):
          x = self.ode_problem.sym
          (C1,)  = self.ode_problem.get_numbered_constants(num=1)
          gensol = Eq(self.ly, (((C1 + Integral((self.cx/self.ax)*exp(Integral(self.bx/self.ax, x)),x))
                  * exp(-Integral(self.bx/self.ax, x)))))
  
          return [gensol]
  
  
  class Bernoulli(SinglePatternODESolver):
      r"""
      Solves Bernoulli differential equations.
  
      These are equations of the form
  
      .. math:: dy/dx + P(x) y = Q(x) y^n\text{, }n \ne 1`\text{.}
  
      The substitution `w = 1/y^{1-n}` will transform an equation of this form
      into one that is linear (see the docstring of
      :py:meth:`~sympy.solvers.ode.single.FirstLinear`).  The general solution is::
  
          >>> from sympy import Function, dsolve, Eq, pprint
          >>> from sympy.abc import x, n
          >>> f, P, Q = map(Function, ['f', 'P', 'Q'])
          >>> genform = Eq(f(x).diff(x) + P(x)*f(x), Q(x)*f(x)**n)
          >>> pprint(genform)
                      d                n
          P(x)*f(x) + --(f(x)) = Q(x)*f (x)
                      dx
          >>> pprint(dsolve(genform, f(x), hint='Bernoulli_Integral'), num_columns=110)
                                                                                                                -1
                                                                                                               -----
                                                                                                               n - 1
                 //         /                                /                           \                    \
                 ||        |                                |                            |                    |
                 ||        |                 /              |                 /          |            /       |
                 ||        |                |               |                |           |           |        |
                 ||        |       (1 - n)* | P(x) dx       |       (1 - n)* | P(x) dx   |  (n - 1)* | P(x) dx|
                 ||        |                |               |                |           |           |        |
                 ||        |               /                |               /            |          /         |
          f(x) = ||C1 - n* | Q(x)*e                   dx +  | Q(x)*e                   dx|*e                  |
                 ||        |                                |                            |                    |
                 \\       /                                /                             /                    /
  
  
      Note that the equation is separable when `n = 1` (see the docstring of
      :py:meth:`~sympy.solvers.ode.ode.ode_separable`).
  
      >>> pprint(dsolve(Eq(f(x).diff(x) + P(x)*f(x), Q(x)*f(x)), f(x),
      ... hint='separable_Integral'))
      f(x)
          /
      |                /
      |  1            |
      |  - dy = C1 +  | (-P(x) + Q(x)) dx
      |  y            |
      |              /
      /
  
  
      Examples
      ========
  
      >>> from sympy import Function, dsolve, Eq, pprint, log
      >>> from sympy.abc import x
      >>> f = Function('f')
  
      >>> pprint(dsolve(Eq(x*f(x).diff(x) + f(x), log(x)*f(x)**2),
      ... f(x), hint='Bernoulli'))
                      1
      f(x) =  -----------------
              C1*x + log(x) + 1
  
      References
      ==========
  
      - https://en.wikipedia.org/wiki/Bernoulli_differential_equation
  
      - M. Tenenbaum & H. Pollard, "Ordinary Differential Equations",
        Dover 1963, pp. 95
  
      # indirect doctest
  
      """
      hint = "Bernoulli"
      has_integral = True
      order = [1]
  
      def _wilds(self, f, x, order):
          P = Wild('P', exclude=[f(x)])
          Q = Wild('Q', exclude=[f(x)])
          n = Wild('n', exclude=[x, f(x), f(x).diff(x)])
          return P, Q, n
  
      def _equation(self, fx, x, order):
          P, Q, n = self.wilds()
          return fx.diff(x) + P*fx - Q*fx**n
  
      def _get_general_solution(self, *, simplify: bool = True):
          P, Q, n = self.wilds_match()
          fx = self.ode_problem.func
          x = self.ode_problem.sym
          (C1,) = self.ode_problem.get_numbered_constants(num=1)
          if n==1:
              gensol = Eq(log(fx), (
              (C1 + Integral((-P + Q),x)
          )))
          else:
              gensol = Eq(fx**(1-n), (
                  (C1 - (n - 1) * Integral(Q*exp(-n*Integral(P, x))
                              * exp(Integral(P, x)), x)
                  ) * exp(-(1 - n)*Integral(P, x)))
              )
          return [gensol]
  
  
  class Factorable(SingleODESolver):
      r"""
          Solves equations having a solvable factor.
  
          This function is used to solve the equation having factors. Factors may be of type algebraic or ode. It
          will try to solve each factor independently. Factors will be solved by calling dsolve. We will return the
          list of solutions.
  
          Examples
          ========
  
          >>> from sympy import Function, dsolve, Eq, pprint, Derivative
          >>> from sympy.abc import x
          >>> f = Function('f')
          >>> eq = (f(x)**2-4)*(f(x).diff(x)+f(x))
          >>> pprint(dsolve(eq, f(x)))
                                          -x
          [f(x) = 2, f(x) = -2, f(x) = C1*e  ]
  
  
          """
      hint = "factorable"
      has_integral = False
  
      def _matches(self):
          eq = self.ode_problem.eq
          f = self.ode_problem.func.func
          x = self.ode_problem.sym
          order =self.ode_problem.order
          df = f(x).diff(x)
          self.eqs = []
          eq = eq.collect(f(x), func = cancel)
          eq = fraction(factor(eq))[0]
          roots = factor_list(eq)[1]
          if len(roots)>1 or roots[0][1]>1:
              for base,expo in roots:
                  if base.has(f(x)):
                      self.eqs.append(base)
              if len(self.eqs)>0:
                  return True
          roots = solve(eq, df)
          if len(roots)>0:
              self.eqs = [(df - root) for root in roots]
              if len(self.eqs)==1:
                  if order>1:
                      return False
                  return fraction(factor(self.eqs[0]))[0]-eq!=0
              return True
          return False
  
  
      def _get_general_solution(self, *, simplify: bool = True):
          func = self.ode_problem.func.func
          x = self.ode_problem.sym
          eqns = self.eqs
          sols = []
          for eq in eqns:
              try:
                  sol = dsolve(eq, func(x))
              except NotImplementedError:
                  continue
              else:
                  if isinstance(sol, list):
                      sols.extend(sol)
                  else:
                      sols.append(sol)
  
          if sols == []:
              raise NotImplementedError("The given ODE " + str(eq) + " cannot be solved by"
                  + " the factorable group method")
          return sols
  
  
  class RiccatiSpecial(SinglePatternODESolver):
      r"""
      The general Riccati equation has the form
  
      .. math:: dy/dx = f(x) y^2 + g(x) y + h(x)\text{.}
  
      While it does not have a general solution [1], the "special" form, `dy/dx
      = a y^2 - b x^c`, does have solutions in many cases [2].  This routine
      returns a solution for `a(dy/dx) = b y^2 + c y/x + d/x^2` that is obtained
      by using a suitable change of variables to reduce it to the special form
      and is valid when neither `a` nor `b` are zero and either `c` or `d` is
      zero.
  
      >>> from sympy.abc import x, y, a, b, c, d
      >>> from sympy.solvers.ode import dsolve, checkodesol
      >>> from sympy import pprint, Function
      >>> f = Function('f')
      >>> y = f(x)
      >>> genform = a*y.diff(x) - (b*y**2 + c*y/x + d/x**2)
      >>> sol = dsolve(genform, y)
      >>> pprint(sol, wrap_line=False)
              /                                 /        __________________       \\
              |           __________________    |       /                2        ||
              |          /                2     |     \/  4*b*d - (a + c)  *log(x)||
             -|a + c - \/  4*b*d - (a + c)  *tan|C1 + ----------------------------||
              \                                 \                 2*a             //
      f(x) = ------------------------------------------------------------------------
                                              2*b*x
  
      >>> checkodesol(genform, sol, order=1)[0]
      True
  
      References
      ==========
  
      1. http://www.maplesoft.com/support/help/Maple/view.aspx?path=odeadvisor/Riccati
      2. http://eqworld.ipmnet.ru/en/solutions/ode/ode0106.pdf -
         http://eqworld.ipmnet.ru/en/solutions/ode/ode0123.pdf
      """
      hint = "Riccati_special_minus2"
      has_integral = False
      order = [1]
  
      def _wilds(self, f, x, order):
          a = Wild('a', exclude=[x, f(x), f(x).diff(x), 0])
          b = Wild('b', exclude=[x, f(x), f(x).diff(x), 0])
          c = Wild('c', exclude=[x, f(x), f(x).diff(x)])
          d = Wild('d', exclude=[x, f(x), f(x).diff(x)])
          return a, b, c, d
  
      def _equation(self, fx, x, order):
          a, b, c, d = self.wilds()
          return a*fx.diff(x) + b*fx**2 + c*fx/x + d/x**2
  
      def _get_general_solution(self, *, simplify: bool = True):
          a, b, c, d = self.wilds_match()
          fx = self.ode_problem.func
          x = self.ode_problem.sym
          (C1,) = self.ode_problem.get_numbered_constants(num=1)
          mu = sqrt(4*d*b - (a - c)**2)
  
          gensol = Eq(fx, (a - c - mu*tan(mu/(2*a)*log(x) + C1))/(2*b*x))
          return [gensol]
  
  
  # Avoid circular import:
  from .ode import dsolve
  