123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273 |
- """Utility functions for classifying and solving
- ordinary and partial differential equations.
- Contains
- ========
- _preprocess
- ode_order
- _desolve
- """
- from sympy.core import Pow
- from sympy.core.function import Derivative, AppliedUndef
- from sympy.core.relational import Equality
- from sympy.core.symbol import Wild
- def _preprocess(expr, func=None, hint='_Integral'):
- """Prepare expr for solving by making sure that differentiation
- is done so that only func remains in unevaluated derivatives and
- (if hint doesn't end with _Integral) that doit is applied to all
- other derivatives. If hint is None, don't do any differentiation.
- (Currently this may cause some simple differential equations to
- fail.)
- In case func is None, an attempt will be made to autodetect the
- function to be solved for.
- >>> from sympy.solvers.deutils import _preprocess
- >>> from sympy import Derivative, Function
- >>> from sympy.abc import x, y, z
- >>> f, g = map(Function, 'fg')
- If f(x)**p == 0 and p>0 then we can solve for f(x)=0
- >>> _preprocess((f(x).diff(x)-4)**5, f(x))
- (Derivative(f(x), x) - 4, f(x))
- Apply doit to derivatives that contain more than the function
- of interest:
- >>> _preprocess(Derivative(f(x) + x, x))
- (Derivative(f(x), x) + 1, f(x))
- Do others if the differentiation variable(s) intersect with those
- of the function of interest or contain the function of interest:
- >>> _preprocess(Derivative(g(x), y, z), f(y))
- (0, f(y))
- >>> _preprocess(Derivative(f(y), z), f(y))
- (0, f(y))
- Do others if the hint doesn't end in '_Integral' (the default
- assumes that it does):
- >>> _preprocess(Derivative(g(x), y), f(x))
- (Derivative(g(x), y), f(x))
- >>> _preprocess(Derivative(f(x), y), f(x), hint='')
- (0, f(x))
- Don't do any derivatives if hint is None:
- >>> eq = Derivative(f(x) + 1, x) + Derivative(f(x), y)
- >>> _preprocess(eq, f(x), hint=None)
- (Derivative(f(x) + 1, x) + Derivative(f(x), y), f(x))
- If it's not clear what the function of interest is, it must be given:
- >>> eq = Derivative(f(x) + g(x), x)
- >>> _preprocess(eq, g(x))
- (Derivative(f(x), x) + Derivative(g(x), x), g(x))
- >>> try: _preprocess(eq)
- ... except ValueError: print("A ValueError was raised.")
- A ValueError was raised.
- """
- if isinstance(expr, Pow):
- # if f(x)**p=0 then f(x)=0 (p>0)
- if (expr.exp).is_positive:
- expr = expr.base
- derivs = expr.atoms(Derivative)
- if not func:
- funcs = set().union(*[d.atoms(AppliedUndef) for d in derivs])
- if len(funcs) != 1:
- raise ValueError('The function cannot be '
- 'automatically detected for %s.' % expr)
- func = funcs.pop()
- fvars = set(func.args)
- if hint is None:
- return expr, func
- reps = [(d, d.doit()) for d in derivs if not hint.endswith('_Integral') or
- d.has(func) or set(d.variables) & fvars]
- eq = expr.subs(reps)
- return eq, func
- def ode_order(expr, func):
- """
- Returns the order of a given differential
- equation with respect to func.
- This function is implemented recursively.
- Examples
- ========
- >>> from sympy import Function
- >>> from sympy.solvers.deutils import ode_order
- >>> from sympy.abc import x
- >>> f, g = map(Function, ['f', 'g'])
- >>> ode_order(f(x).diff(x, 2) + f(x).diff(x)**2 +
- ... f(x).diff(x), f(x))
- 2
- >>> ode_order(f(x).diff(x, 2) + g(x).diff(x, 3), f(x))
- 2
- >>> ode_order(f(x).diff(x, 2) + g(x).diff(x, 3), g(x))
- 3
- """
- a = Wild('a', exclude=[func])
- if expr.match(a):
- return 0
- if isinstance(expr, Derivative):
- if expr.args[0] == func:
- return len(expr.variables)
- else:
- order = 0
- for arg in expr.args[0].args:
- order = max(order, ode_order(arg, func) + len(expr.variables))
- return order
- else:
- order = 0
- for arg in expr.args:
- order = max(order, ode_order(arg, func))
- return order
- def _desolve(eq, func=None, hint="default", ics=None, simplify=True, *, prep=True, **kwargs):
- """This is a helper function to dsolve and pdsolve in the ode
- and pde modules.
- If the hint provided to the function is "default", then a dict with
- the following keys are returned
- 'func' - It provides the function for which the differential equation
- has to be solved. This is useful when the expression has
- more than one function in it.
- 'default' - The default key as returned by classifier functions in ode
- and pde.py
- 'hint' - The hint given by the user for which the differential equation
- is to be solved. If the hint given by the user is 'default',
- then the value of 'hint' and 'default' is the same.
- 'order' - The order of the function as returned by ode_order
- 'match' - It returns the match as given by the classifier functions, for
- the default hint.
- If the hint provided to the function is not "default" and is not in
- ('all', 'all_Integral', 'best'), then a dict with the above mentioned keys
- is returned along with the keys which are returned when dict in
- classify_ode or classify_pde is set True
- If the hint given is in ('all', 'all_Integral', 'best'), then this function
- returns a nested dict, with the keys, being the set of classified hints
- returned by classifier functions, and the values being the dict of form
- as mentioned above.
- Key 'eq' is a common key to all the above mentioned hints which returns an
- expression if eq given by user is an Equality.
- See Also
- ========
- classify_ode(ode.py)
- classify_pde(pde.py)
- """
- if isinstance(eq, Equality):
- eq = eq.lhs - eq.rhs
- # preprocess the equation and find func if not given
- if prep or func is None:
- eq, func = _preprocess(eq, func)
- prep = False
- # type is an argument passed by the solve functions in ode and pde.py
- # that identifies whether the function caller is an ordinary
- # or partial differential equation. Accordingly corresponding
- # changes are made in the function.
- type = kwargs.get('type', None)
- xi = kwargs.get('xi')
- eta = kwargs.get('eta')
- x0 = kwargs.get('x0', 0)
- terms = kwargs.get('n')
- if type == 'ode':
- from sympy.solvers.ode import classify_ode, allhints
- classifier = classify_ode
- string = 'ODE '
- dummy = ''
- elif type == 'pde':
- from sympy.solvers.pde import classify_pde, allhints
- classifier = classify_pde
- string = 'PDE '
- dummy = 'p'
- # Magic that should only be used internally. Prevents classify_ode from
- # being called more than it needs to be by passing its results through
- # recursive calls.
- if kwargs.get('classify', True):
- hints = classifier(eq, func, dict=True, ics=ics, xi=xi, eta=eta,
- n=terms, x0=x0, hint=hint, prep=prep)
- else:
- # Here is what all this means:
- #
- # hint: The hint method given to _desolve() by the user.
- # hints: The dictionary of hints that match the DE, along with other
- # information (including the internal pass-through magic).
- # default: The default hint to return, the first hint from allhints
- # that matches the hint; obtained from classify_ode().
- # match: Dictionary containing the match dictionary for each hint
- # (the parts of the DE for solving). When going through the
- # hints in "all", this holds the match string for the current
- # hint.
- # order: The order of the DE, as determined by ode_order().
- hints = kwargs.get('hint',
- {'default': hint,
- hint: kwargs['match'],
- 'order': kwargs['order']})
- if not hints['default']:
- # classify_ode will set hints['default'] to None if no hints match.
- if hint not in allhints and hint != 'default':
- raise ValueError("Hint not recognized: " + hint)
- elif hint not in hints['ordered_hints'] and hint != 'default':
- raise ValueError(string + str(eq) + " does not match hint " + hint)
- # If dsolve can't solve the purely algebraic equation then dsolve will raise
- # ValueError
- elif hints['order'] == 0:
- raise ValueError(
- str(eq) + " is not a solvable differential equation in " + str(func))
- else:
- raise NotImplementedError(dummy + "solve" + ": Cannot solve " + str(eq))
- if hint == 'default':
- return _desolve(eq, func, ics=ics, hint=hints['default'], simplify=simplify,
- prep=prep, x0=x0, classify=False, order=hints['order'],
- match=hints[hints['default']], xi=xi, eta=eta, n=terms, type=type)
- elif hint in ('all', 'all_Integral', 'best'):
- retdict = {}
- gethints = set(hints) - {'order', 'default', 'ordered_hints'}
- if hint == 'all_Integral':
- for i in hints:
- if i.endswith('_Integral'):
- gethints.remove(i[:-len('_Integral')])
- # special cases
- for k in ["1st_homogeneous_coeff_best", "1st_power_series",
- "lie_group", "2nd_power_series_ordinary", "2nd_power_series_regular"]:
- if k in gethints:
- gethints.remove(k)
- for i in gethints:
- sol = _desolve(eq, func, ics=ics, hint=i, x0=x0, simplify=simplify, prep=prep,
- classify=False, n=terms, order=hints['order'], match=hints[i], type=type)
- retdict[i] = sol
- retdict['all'] = True
- retdict['eq'] = eq
- return retdict
- elif hint not in allhints: # and hint not in ('default', 'ordered_hints'):
- raise ValueError("Hint not recognized: " + hint)
- elif hint not in hints:
- raise ValueError(string + str(eq) + " does not match hint " + hint)
- else:
- # Key added to identify the hint needed to solve the equation
- hints['hint'] = hint
- hints.update({'func': func, 'eq': eq})
- return hints
|