Skip to content

Instantly share code, notes, and snippets.

@aplund
Created March 31, 2022 00:26
Show Gist options
  • Save aplund/888da79a96fb295e7d614f43780119f9 to your computer and use it in GitHub Desktop.
Save aplund/888da79a96fb295e7d614f43780119f9 to your computer and use it in GitHub Desktop.
[ins] In [1]: import sympy
[ins] In [2]: x=sympy.Function('x')
...: t=sympy.Function('t')
...: s=sympy.Symbol('s',real=True)
...: de_system=(sympy.Eq(x(s).diff(s),t(s).diff(s)), sympy.Eq(t(s).diff(s),1))
...: sympy.dsolve(de_system,funcs=(x,t),t=s,ics={x(0):1,t(0):0})
---------------------------------------------------------------------------
NotAlgebraic Traceback (most recent call last)
Input In [2], in <cell line: 5>()
3 s=sympy.Symbol('s',real=True)
4 de_system=(sympy.Eq(x(s).diff(s),t(s).diff(s)), sympy.Eq(t(s).diff(s),1))
----> 5 sympy.dsolve(de_system,funcs=(x,t),t=s,ics={x(0):1,t(0):0})
File /usr/lib/python3.10/site-packages/sympy/solvers/ode/ode.py:559, in dsolve(eq, func, hint, simplify, ics, xi, eta, x0, n, **kwargs)
553 # This may have to be changed in future
554 # when we have weakly and strongly
555 # connected components. This have to
556 # changed to show the systems that haven't
557 # been solved.
558 try:
--> 559 sol = dsolve_system(eq, funcs=func, ics=ics, doit=True)
560 return sol[0] if len(sol) == 1 else sol
561 except NotImplementedError:
File /usr/lib/python3.10/site-packages/sympy/solvers/ode/systems.py:2131, in dsolve_system(eqs, funcs, t, ics, doit, simplify)
2129 if ics:
2130 constants = Tuple(*sol).free_symbols - variables
-> 2131 solved_constants = solve_ics(sol, funcs, constants, ics)
2132 sol = [s.subs(solved_constants) for s in sol]
2134 if simplify:
File /usr/lib/python3.10/site-packages/sympy/solvers/ode/ode.py:793, in solve_ics(sols, funcs, constants, ics)
791 # TODO: Use solveset here
792 try:
--> 793 solved_constants = solve(subs_sols, constants, dict=True)
794 except NotImplementedError:
795 solved_constants = []
File /usr/lib/python3.10/site-packages/sympy/solvers/solvers.py:1108, in solve(f, *symbols, **flags)
1106 solution = _solve(f[0], *symbols, **flags)
1107 else:
-> 1108 solution = _solve_system(f, symbols, **flags)
1110 #
1111 # postprocessing
1112 ###########################################################################
1113 # Restore masked-off objects
1114 if non_inverts:
File /usr/lib/python3.10/site-packages/sympy/solvers/solvers.py:1768, in _solve_system(exprs, symbols, **flags)
1766 subsyms = list(sorted(subsyms, key = lambda x: sym_indices[x]))
1767 flags['_split'] = False # skip split step
-> 1768 subsol = _solve_system(subexpr, subsyms, **flags)
1769 if not isinstance(subsol, list):
1770 subsol = [subsol]
File /usr/lib/python3.10/site-packages/sympy/solvers/solvers.py:1799, in _solve_system(exprs, symbols, **flags)
1796 failed.append(g)
1797 continue
-> 1799 poly = g.as_poly(*symbols, extension=True)
1801 if poly is not None:
1802 polys.append(poly)
File /usr/lib/python3.10/site-packages/sympy/core/expr.py:1116, in Expr.as_poly(self, *gens, **args)
1113 from sympy.polys.polytools import Poly
1115 try:
-> 1116 poly = Poly(self, *gens, **args)
1118 if not poly.is_Poly:
1119 return None
File /usr/lib/python3.10/site-packages/sympy/polys/polytools.py:181, in Poly.__new__(cls, rep, *gens, **args)
179 return cls._from_poly(rep, opt)
180 else:
--> 181 return cls._from_expr(rep, opt)
File /usr/lib/python3.10/site-packages/sympy/polys/polytools.py:311, in Poly._from_expr(cls, rep, opt)
309 """Construct a polynomial from an expression. """
310 rep, opt = _dict_from_expr(rep, opt)
--> 311 return cls._from_dict(rep, opt)
File /usr/lib/python3.10/site-packages/sympy/polys/polytools.py:255, in Poly._from_dict(cls, rep, opt)
252 domain = opt.domain
254 if domain is None:
--> 255 domain, rep = construct_domain(rep, opt=opt)
256 else:
257 for monom, coeff in rep.items():
File /usr/lib/python3.10/site-packages/sympy/polys/constructor.py:367, in construct_domain(obj, **args)
364 coeffs = [obj]
366 coeffs = list(map(sympify, coeffs))
--> 367 result = _construct_simple(coeffs, opt)
369 if result is not None:
370 if result is not False:
File /usr/lib/python3.10/site-packages/sympy/polys/constructor.py:65, in _construct_simple(coeffs, opt)
62 max_prec = max(c._prec for c in float_numbers) if float_numbers else 53
64 if algebraics:
---> 65 domain, result = _construct_algebraic(coeffs, opt)
66 else:
67 if floats and complexes:
File /usr/lib/python3.10/site-packages/sympy/polys/constructor.py:105, in _construct_algebraic(coeffs, opt)
102 trees = build_trees(coeffs)
103 exts = list(ordered(exts))
--> 105 g, span, H = primitive_element(exts, ex=True, polys=True)
106 root = sum([ s*ext for s, ext in zip(span, exts) ])
108 domain, g = QQ.algebraic_field((g, root)), g.rep.rep
File /usr/lib/python3.10/site-packages/sympy/polys/numberfields/subfield.py:365, in primitive_element(extension, x, ex, polys)
362 return cls(g), coeffs
364 gen, coeffs = extension[0], [1]
--> 365 f = minimal_polynomial(gen, x, polys=True)
366 K = QQ.algebraic_field((f, gen)) # incrementally constructed field
367 reps = [K.unit] # representations of extension elements in K
File /usr/lib/python3.10/site-packages/sympy/polys/numberfields/minpoly.py:701, in minimal_polynomial(ex, x, compose, polys, domain)
697 raise GeneratorsError("the variable %s is an element of the ground "
698 "domain %s" % (x, domain))
700 if compose:
--> 701 result = _minpoly_compose(ex, x, domain)
702 result = result.primitive()[1]
703 c = result.coeff(x**degree(result, x))
File /usr/lib/python3.10/site-packages/sympy/polys/numberfields/minpoly.py:618, in _minpoly_compose(ex, x, dom)
616 res = _minpoly_rootof(ex, x)
617 else:
--> 618 raise NotAlgebraic("%s doesn't seem to be an algebraic element" % ex)
619 return res
NotAlgebraic: Integral(1, (s, 0)) doesn't seem to be an algebraic element
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment