GSoC 2016 Application Kshitij Saraogi: Solvers: Extending Solveset - sympy/sympy GitHub Wiki

About Me

Basic Information

Name: Kshitij Saraogi
University: Indian Institute of Technology, Kharagpur
Major: Mathematics and Computing
Email: [email protected]
Github/IRC : kshitij10496
Blog : #TODO
Time Zone : IST (UTC + 5:30)

Personal Background

I am Kshitij Saraogi, a second year undergraduate student at IIT Kharagpur, India. I'm pursuing a degree in Mathematics and Computing.I have taken courses on Discrete Mathematics, Set Theory, Higher Algebra, Linear Algebra, Numerical Methods, ODE and PDE. I am also familiar with real and complex analysis. In addition to this, I have a basic idea of abstract algebra as well.

Programming Details

I work on Ubuntu 15.10 with Vim as my primary editor. I like Vim because of the simple reason that it is a modal editor. The different modes make navigating the files with the keyboard a lot faster. Also the fact that each key on the keyboard has a different meaning in each mode intrigues me.

I was introduced to Object Oriented Programming at an early age and I have been coding for nearly 6 years now. I started with Java and then followed it up with C. I switched to Python last summer due to it’s intuitive syntax.

Favourite features of Python:

High emphasis on code readability (PEP 008).
Given my past experiences with Java and C, I find Python’s dynamic type system quite amazing.
Python can be used for a broad range of programming tasks, from little shell scripts to building a web applications to scientific uses. It may not be as good at any of those as a purpose-built programming language but it can do all of them, and do them well. This is why I appreciate the power of this interpreted, dynamic language. I started to find the fun in programming because of Python.

SymPy has so many wonderful features. Out of all the cool features of SymPy, what amazes me the most are the ConditionSet and ImageSet classes.

The ConditionSet object represents a set of elements, in a given super set, which satisfy a given condition. Apart from the wide range of application this can be used for, I was intrigued by its most basic application : containing an variable in an Interval. For example, the complete solution of Abs(x) - n = 0 we need something like {- n, n} given n >= 0 This can be achieved by ConditionSet with something like: {x | x ∊ {-n, n} ∧ (n ∈ [0, ∞))}
An ImageSet object represents the image of a Set under a mathematical function. This can be used for such wide variety of operations on a whole set of numbers. For example, if we want an object to represent all the rational numbers, it is as simple as ImageSet(Lambda((x, y), x/y), S.Integers * S.Naturals). Without an ImageSet object, doing this would have been a lot complex.

Also, I am familiar with using Git version control system and Github. I have been using them for the past 6 months and learnt them through practice.
I have set up my development environment and am familiar with SymPy’s workflow.
If I get stuck, I search it up, read about it at Stack Overflow, and come up with a solution.

Contributions to SymPy

Pull Requests: Bugs Reported:

The Project

Brief Overview

Solving equations is a quintessential feature of any Computer Algebra System. And being able to solve a variety of equations, with accuracy, adds new dimensions to the system. In order to understand what needs to be done, we need to know exactly what has already been done. There have been 2 Google Summer of Code Projects on Solvers in the last couple of years.

In 2014, Harsh Gupta during his Google Summer of Code project improved Solvers and wrote a new submodule named solveset. The aim of this new module was to replace solve when it is mature enough. For this, Harsh rewrote the univariate solver, and built the basic Set infrastructure to support the solutions returned for solveset.

In 2015, Amit Kumar extended the solveset module to support solving System of Linear Equations by adding the linsolve solver. Amit also built a ComplexPlane class, the Complex Set infrastructure, in order to represent the regions of Complex domain. He started replacing solve with solveset in the codebase, with the idea of making a smooth transition from solve to solveset. However, Harsh figured out that the solveset is vunerable to API changes, so it's not advisable to replace solve with solveset and the work was halted.

The major issues with merging solveset with solve are that the lack of solvers for univariate Transcendental Equations and Multivariate equations. I, on my part, intend to implement these and a few more solvers to extend the solveset module and make it more robust.

The Plan

We have two different APIs for solving equations - solve and solveset. The solve module is capaable of solving wide range of equation but its API is obsolete. Currently, the new solveset module can be used only for solving univarate equations(single equation for a single variable) and a system of linear equations (with N variables and M equations). In order to extend solveset to atleast replace the old solve module, we need to implement the transcendental equation solver (equations solvable by LambertW function) and a non-linear multivariate system of equations solver. Also, for solving ineqaulities, we need to incorporate the solve_univariate_inequality according to the new solveset API.

Basic Idea

A robust framework for solving equations is a salient feature of any Computer Algebra System. Two important ideas in this regard are Rewriting and Decomposition. Basically, by using these techniques, we try to reduce the given equation to a set of equations which can be easily solved. An example here would be the best way to establish the significance of the above ideas.

If we have to solve for a function f(x) = 0, given by:

$f(x) = 4^{x} - 3 (2^{x}) + 2$

By the method of Decomposition, we can decompose f(x) as :

$f(x) = g(x)^{2} - 3*g(x) + 2$

where g(x) is given by:

$g(x) = 2^{x}$

Substituting g(x) with a dummy variable t, we get a polynomial in f(t) given as:

$f(t) = (t - 1)* (t - 2)$

Resubstituting the value of t as g(x) and further composing f(x) from g(x), we get:

$f(x) = (2^{x} - 1)*(2^{x} - 2)$

We can now represent f(x) as the product of two simpler functions by the idea of Rewriting.

$f(x) = k(x) * l(x)$

Calculating the solutions of these functions is comparatively easier with respect to the original equation. We solve for these functions separately and get individual solutions as:

$2^{x} - 1 = 0$ or $2^{x} - 2 = 0$

$\Rightarrow x = 0$ or
$\Rightarrow x = 1$

Thus, we get two possible solutions of the given transcendental equation through the methods of Rewriting and Decomposition.

1. Equations solvable by LambertW function( Transcendental Equation Solver)

A transcendental equation is an equation which cannot be expressed in terms of a finite sequence of the algebraic operations of addition, multiplication, and root extraction, in the variable being solved for. Generally, such equations donot have a closed-form solution, i.e, the solution cannot be evaluated in a finite number of operations. In some cases, special functions (like Lambert W function) can be used to write the solutions to transcendental equations in closed form.

Brief Overview of the LambertW Function:

The Lambert W function, also called the omega function or product logarithm, is a set of functions, namely the branches of the inverse relation of the function
$f(z) = ze^{z}$ , where $e^{z}$ is the exponential function and $z$ is any complex number.
For a real number x, we call write the above as:

$x = LambertW(xe^{x})$

In general, we can write the LambertW function in the form:

$LambertW( foo \cdot e^{foo}) = foo$

where foo can be any function of the real variable x.

We can denote the LambertW(x) by W(x). And on further substituting W(x) with t, we get:

$t = W(t) \cdot e^{W(t)}$

Diving both sides by W(t), we get:

$e^{W(t)} = \frac{t}{W(t)}$

This result can also be written in a subtle general form:

$e^{W(foo)} = \frac{foo}{W(foo)}$

The above results, in their general form, are useful identities for solving exponential equations. We will try to manipulate the given equation to either form for solving it with the help of LambertW function.

Branches of LambertW function

The equation $z = W(z)\cdot e^{W(z)}$ has infinitely many solutions on the complex plane. These solutions are represented by $W_{k}(x)$ with the branch index k ranging over the integers.
If x is real, then for -1/e <= x < 0, there are two possible values of W(x) or LambertW(x) (Refer Graph) :

The branch satisfying -1 <= W(x), known as the Principal Branch, is denoted by $W_{0}(x)$ .
The branch satisfying W(x) <= -1 is denoted by $W_{-1}(x)$

Summary

The general equation of the following form:

$p^{ax+b} = cx + d$

can be solved by the use of LambertW function through suitable substitution and decomposition.

The final solution of the above general equation is

$x = -\frac{W(-\frac{a \log_e{p}}{c} p ^{b - \frac{ad}{c}})}{a \log_e{p}}-\frac{d}{c}$

Thus, the LambertW Function can be used to write the solutions of an extensive class of exponential equations in closed form. Currently, only the principal branch of LambertW function( $W_{0}(x)$ ) is implemented in old solve, which is one of the reasons for loss of solutions[] and the solveset doesn't support this as of yet.

The Execution

Equations solvable by LambertW function

The general approach for solving a transcendental equation using LambertW function (principal branch) :

Manipulate all occurrences of the unknown variable x into an expression of the form $f(x)\cdot e^{f(x)}$ .
Taking LambertW function of both the LHS and the RHS.
Using the following identity for simplification:

$LambertW( foo \cdot e^{foo}) = foo$

Solving for x in terms of the LambertW function.

A basic implementation of the above approach for solving the equation $5t - 2^{t} = 0$ can be:

This can also be written as : $5t = 2^{t}$

Step 01 : Manipulating the LHS until the expression resembles the desired form :

$5t \cdot 2^{-t} = 1$

$5t \cdot e^{-t\cdot \log_{e}2} = 1$

$- \log_{e}2 \cdot t \cdot e^{-t\cdot \log_{e}2} = - \frac{\log_{e}2}{5}$

Step 02 : We take LambertW functions of both the sides of the equation as :

$LambertW(- \log_{e}2 \cdot t \cdot e^{-t\cdot \log_{e}2}) = LambertW(- \frac{\log_{e}2}{5})$

Step 03 : Using the identity, we get:

$- \log_{e}2 \cdot t = LambertW(- \frac{\log_{e}2}{5})$

Step 04 : Solving for x and simplifying, the result is given as:

$t = - \frac{LambertW(- \frac{\log_{e}2}{5})}{\log_{e}2}$

We already have a _tsolve function in the solve module for solving transcendental equations. I intend to extend it to new solveset module. After completing the task, the final API should look like :

In []: solveset(x + log(x), x)
Out[]: {LambertW(1)}

In []: solveset(x + exp(x**2), x)
Out[]: {I*sqrt(LambertW(-2)/2)}

In []: solveset(2**x - 5*x, x)
Out[]:-LambertW(-ln(2)/5)/ln(2))