How can I solve a differential equation including the Heaviside function in MATLAB R2014a?

I am trying to symbolically solve the following differential equation that contains a Heaviside function:

 

V'(t) = 10 - (1 + 3*heaviside(t-5))*V(t)^2

 

with the following initial conditions:

 

V(0) = 0

 

 

When I try to use "dsolve" and type the following in MATLAB:

 

syms V(t) t

dsolve(diff(V)==10-(1+3*heaviside(t-5))*V^2,V(0)==0)

 

 

I received the following:

 

Warning: Explicit solution could not be found. 

> In dsolve at 197 

 

ans =

 

[ empty sym ]

Answer 1

[MathWorks Support Team]

The “dsolve” function can be used to find solutions to a simple ODE with a piecewise constant right-hand side when the ODE can be reduced to an integration problem.

For instance, for an ODE of the form

 

V’(t) = a(t)*V(t)+b(t)

“dsolve” finds a solution even when a Heaviside function is present as it uses the function “int” under the hood. You can see for instance that the result of the following command:

 

dsolve(diff(V)==(-1+2*heaviside(t-1))*V,V(0)==1)

is

ans =

 

exp(2*heaviside(t - 1)*(t - 1) - t)

 

 

For more complex ODEs, “dsolve” cannot be used when the right-hand side contains a piecewise-constant function at this time.

There are two possible workarounds.

 

1. Split the ODE into its continuous parts

 

For instance, the ODE above can be solved for the two continuous parts of the Heaviside function:

syms V(t) t

% Solution that applies to the first constant piece 

sol1 = dsolve(diff(V)==10-V^2,V(0)==0); 

% Value of the solution at the discontinuity 

Vs = double(subs(sol1,t,5)); 

% Solution that applies after the discontinuity

sol2 = dsolve(diff(V)==10-4*V^2,V(5)==Vs); 

 

2. Use a constant to represent the piecewise function

For instance, use a constant named “HEAV” in place of the Heaviside function in the ODE. Solve using "dsolve" and substitute the constant for the actual Heaviside function in the solution:

 

syms V(t) t HEAV

sol3 = dsolve(diff(V)==10-(1+3*HEAV)*V^2,V(0)==0);

sol3 = subs(sol3, HEAV, heaviside(t-5));

 

This technique must be used with precaution since the “dsolve” function treats “HEAV” as a constant. As a consequence, the solution could potentially contain expressions of the form 1/HEAV + 1/(HEAV-1) that make the solution undefined. The result has to be checked for these scenarios.