Ask Your Question

Integral equation solver

asked 2021-01-08 04:04:10 -0600

danielvolinski gravatar image

updated 2021-01-09 03:21:34 -0600

slelievre gravatar image

Can Sage solve integral equations (where the unknown is a function and the equation involves integrals), either analytically or numerically?


For example, can Sage solve a Volterra equation of the second kind, such as

y(x) == integrate(exp(x-t)*y(t), t, 0, x)

The solution is expected in the form $y(x) = ...$ where the right hand side is an explicit function of $x$ not involving $y$.

Notice in the example the function $y$ appears in two places, once in the left hand side and the other inside the integral sign.

I am looking for something like what you see in

i.e., a command that would give the solution once I give the equation.

edit retag flag offensive close merge delete


To make this page easier to read for anyone landing here, I edited the question to include all the bits and pieces spread around in comments, and I deleted the no longer relevant comments, hope that is okay for everyone.

slelievre gravatar imageslelievre ( 2021-01-09 03:26:38 -0600 )edit

1 answer

Sort by ยป oldest newest most voted

answered 2021-01-08 16:50:35 -0600

Emmanuel Charpentier gravatar image

updated 2021-01-09 09:40:16 -0600

Sagemath has currently no "canned" solution for this kind of integro-differential equations. However, such a problem can sometimes be transformed in (a system of) sagemath-solvable differential equations.

In the example given by Daniel Volinski in the second comment :

var("x, t")
y = function("y")
E1 = y(x) == integrate(e^(-t + x)*y(t), t, 0, x)
y(x) == integrate(e^(-t + x)*y(t), t, 0, x)

[ Note : This does not seem to fullfill the definition of a [Volterra equation](]

One can clean up the problem by noting that in the righ-hand size, $e^x$ is a constant not dependent of the integration variable $t$, hence expressed out of the integrand :

with assuming(x>0): E2=E1.subs(E1.rhs()==e^x*integrate(e^-t*y(t),t,0,x))
with assuming(x<0): E3=E1.subs(E1.rhs()==e^x*integrate(e^-t*y(t),t,0,x))

The identity of E2 and E3 is easily checked :

sage: E2
y(x) == e^x*integrate(e^(-t)*y(t), t, 0, x)
sage: E3
y(x) == -e^x*integrate(e^(-t)*y(t), t, x, 0)

Differentiating (the two sides of) this equation with respect to $x$ gives a differential equation solvable by Sagemath

with assuming(x>0): S2=desolve((E2/e^x).diff(x),y(x), ivar=x)
# Declare the integration constants :
[var(str(u)) for u in S2.variables()]
[_C, x]

However, this candidate solution satisfies the original equation only for the trivial case of the null function :

[_C == 0]

By the way :

sage: E1.substitute_function(y,S2.function(x)).unhold()
_C*e^(2*x) == _C*(e^(2*x) - e^x)
sage: E2.substitute_function(y,S2.function(x)).unhold()
_C*e^(2*x) == _C*(e^x - 1)*e^x

EDIT : As already told, there is, to the best of my knowledge", no pre-programmed function for this in Sagemath. Feel free to submit one...

For the hell of it, I verified that a similar procedure was able to solve the equation proposed in the Mathematica example :

var("x, t, lambda_")
y = function("y")
E1 = y(x) == x^3 + lambda_*integrate((t-x)*y(t), t, 0, x)`

One notes that Sagemath insists to reformat this equation "its way" :

$$ y\left(x\right) = x^{3} - \lambda {\left(\int_{0}^{x} -t y\left(t\right)\,{d t} + \int_{0}^{x} x y\left(t\right)\,{d t}\right)} $$

Isolating the integral terms and differentiating twice gives us a second-order differential equation :

-(6*x - diff(y(x), x, x))/lambda_ == -y(x)

Sagemath default ODE solver (i. eI Maxima's) depends on the sign of $\lambda$ :

with assuming(lambda_>0): Sp = desolve(((E1-x^3)/lambda_).diff(x, 2), y(x), ivar=x)
[var(str(u)) for u in Sp.variables()] # Declaring the integration constants...
[_K1, _K2, lambda_, x]
_K2*cos(sqrt(lambda_)*x) + _K1*sin(sqrt(lambda_)*x) + 6*x/lambda_

Determining the constants is done by substituting in the original equation, simplifying a bit and solving for them the resulting polynomial in $x$ :

Ep = E1.substitute_function(y,Sp.function(x)).unhold().expand()-_K1*sin(x*sqrt(lambda_))-_K2*cos(x*sqrt(lambda_))
Consp = solve((Ep.lhs()-Ep.rhs()).coefficients(x, sparse=False),[_K1, _K2], solution_dict=True)
Solp = [Sp.subs(u).function(x) for u in Consp]
sage: Ep
6*x/lambda_ == -_K1*sqrt(lambda_)*x - _K2
sage: Consp
[{_K1: -6/lambda_^(3/2), _K2: 0}]
sage: Solp
[x |--> 6*x/lambda_ - 6*sin(sqrt(lambda_)*x)/lambda_^(3/2)]

$$\left[x \ {\mapsto}\ \frac{6\,x}{\lambda} - \frac{6 \, \sin\left(\sqrt{\lambda} x\right)}{\lambda^{\frac{3}{2}}}\right]$$

This unique solution is identical to Mathematica's proposed solution.

Similarly, for $\lambda<0$, we get :

with assuming(lambda_<0): Sn = desolve(((E1-x^3)/lambda_).diff(x, 2), y(x), ivar=x)
En = E1.substitute_function(y,Sn.function(x)).unhold().expand()-_K1*e^(-sqrt(-lambda_)*x)-_K2*e^(sqrt(-lambda_)*x)
Consn = solve((En.lhs()-En.rhs()).coefficients(x, sparse=False),[_K1, _K2], solution_dict=True)
Soln = [Sn.subs(u).maxima_methods().demoivre().trig_reduce().function(x) for u in Consn]
sage: Soln
[x |--> 6*x/lambda_ - 6*sin(sqrt(lambda_)*x)/lambda_^(3/2)]

This solution is formally identical to the solution obtained for $\lambda>0$ ; however, its meaning is different, $\lambda$ being negative...


edit flag offensive delete link more


Thank you,

I was looking for something like what you see in Solve a Volterra Integral Equation , not something manual as you did but a command that would give the solution once I give the equation.

Daniel Volinski

danielvolinski gravatar imagedanielvolinski ( 2021-01-08 21:23:46 -0600 )edit

What I'd suggest now is to put this procedure and knowledge into a command so that there will be no need to do it manually. You will just need to write the Integral Equation and use this new command to obtain the required solution. Daniel Volinski

danielvolinski gravatar imagedanielvolinski ( 2021-01-10 04:07:26 -0600 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower


Asked: 2021-01-08 04:04:10 -0600

Seen: 102 times

Last updated: Jan 09