When solving differential equations with Heaviside functions, a meijerg function appears in the outcome, while a (simpler) expression is possible as well. For example a beam problem:
import sympy as sp
w = sp.symbols('w', cls=sp.Function)
x = sp.symbols('x')
C1, C2, C3, C4 = sp.symbols('C1 C2 C3 C4')
Q, a, b, EI = sp.symbols('Q, a, b, EI')
q = Q*(1 - sp.Heaviside(x-a))+F*sp.DiracDelta(x-(a+b))
diffeq = sp.Eq(EI*sp.diff(w(x),x,4),q)
w = sp.dsolve(diffeq)
w = w.rhs
print(w)
gives:
C1 + x**3*(C4 - Q*Piecewise((Piecewise((0, Abs(x/a) < 1), (a*meijerg(((2, 1), ()), ((), (1, 0)), x/a), True)), a > 0), (Integral(Heaviside(-a + x), x), True))/(6*EI) + 35*Heaviside(-a - b + x)/(6*EI)) + x**2*(C3 + Q*Piecewise((Piecewise((0, Abs(x/a) < 1), (a**2*meijerg(((3, 1), ()), ((), (2, 0)), x/a), True)), a > 0), (Integral(x*Heaviside(-a + x), x), True))/(2*EI) - 35*a*Heaviside(-a - b + x)/(2*EI) - 35*b*Heaviside(-a - b + x)/(2*EI)) + x*(C2 - Q*Piecewise((Piecewise((0, Abs(x/a) < 1), (a**3*meijerg(((4, 1), ()), ((), (3, 0)), x/a), True)), a > 0), (Integral(x**2*Heaviside(-a + x), x), True))/(2*EI) + 35*(a + b)**2*Heaviside(-a - b + x)/(2*EI)) + Q*x**4/(24*EI) + Q*Piecewise((Piecewise((0, Abs(x/a) < 1), (a**4*meijerg(((5, 1), ()), ((), (4, 0)), x/a), True)), a > 0), (Integral(x**3*Heaviside(-a + x), x), True))/(6*EI) - 35*(a + b)**3*Heaviside(-a - b + x)/(6*EI)
Thus including meijerg functions.
If the same result is found by integrating multiple times with meijerg=False:
import sympy as sp
w = sp.symbols('w', cls=sp.Function)
x = sp.symbols('x')
C1, C2, C3, C4 = sp.symbols('C1 C2 C3 C4')
Q, a, b, EI = sp.symbols('Q, a, b, EI')
q = Q*(1 - sp.Heaviside(x-a))+F*sp.DiracDelta(x-(a+b))
V = sp.integrate(-q,x,meijerg=False)+C1
M = sp.integrate(V,x,meijerg=False)+C2
kappa = M / EI
phi = sp.integrate(kappa,x,meijerg=False)+C3
w = sp.integrate(-phi,x,meijerg=False)+C4
print(w)
the following result is obtained:
-C3*x + C4 - (C1*x**3/6 + C2*x**2/2 - Q*(x**4/24 - (a**4/24 - a**3*x/6 + a**2*x**2/4 - a*x**3/6 + x**4/24)*Heaviside(-a + x)) - 35*(-a*x**2/2 + a*x*(a + b) - a*(a + b)**2/2 - b*x**2/2 + b*x*(a + b) - b*(a + b)**2/2 + x**3/6 - x*(a + b)**2/2 + (a + b)**3/3)*Heaviside(-a - b + x))/EI
Would it be possible to implement a meijerg=false option for dsolve as well?
Kind regards,
Tom van Woudenberg