ODE symbolic result plotting with fplot()

Hello,

So I am trying to solve an ODE using the symbolic toolbox.

syms phi(t) g l d m
dphi = diff(phi,t); % Derivative of phi(t)

phi(t) = dsolve(diff(phi,t,t)== -g/l*phi-d/(m*l^2)*diff(phi,t),... % ODE (governing the movement of a pendulum)
                phi(0) ==1,... % Initial condition 1
                dphi(0)==0) % Initial condition 2


phi(t) = subs(phi(t),{l,m,g,d},{10,5,9.81,50}) % Replacing with known values
 
fplot(phi,[1 100]) % symbolic plotting between 1 and 100 s

I get the expected result but when I plot it, I get numerous vertical lines that add noise to my plot and it is almost unreadable. I am using Matlab on a Mac. How do i get rid of those lines ?

Answer 1

[Star Strider]

I am not certain what the problem is with the original function.

If you simplify it first:

phi = vpa(simplify(phi, 'Steps',500),5)

the vertical lines (indicating infinite results) disappear. (I use vpa here to shorten the output so I can see all of it. It is not necessary for the code.)

[M.Many]

Thanks for the answer, it works indeed but the vpa function just evaluates the terms (fractions,square roots etc...) in the symbolic phi function, so it is not anymore the exact solution.

I think I figured out the problem, phi is actually a complex function that's why it doesn't plot so well. I have to take the real part of phi and then plot it to get the expected solution and now it works. Thanks !

fplot(real(phi),[1 100]) % symbolic plotting between 1 and 100 s

[Star Strider]

As always, my pleasure!

The vpa function affects only the display of the expression. MATLAB maintains full internal symbolic precision, unless you tell it otherswise (using round, fix, ceil, floor, for example, that may not apply to symbolic expressions, or convert it to an anonymous function with matlabFunction, that then uses full floating-point precision).

When I plotted both the real and imaginary parts of ‘phi’, the imaginary part was identically 0. Also, if it had not been so, fplot (and other ploltting functions) would have thrown a warning about ignoring the imaginary parts.

[M.Many]

Yes you are right ! The simplify() command somehow doesn't want to further simplify the complex expression of the real function phi, since, as you said, the imaginary part is zero.

[Star Strider]

Thank you!

[M.Many]

However, the documentation on the vpa function clearly tells it approximates the symbolic expression

"vpa(x) uses variable-precision floating-point arithmetic (VPA) to evaluate each element of the symbolic input x to at least d significant digits, where d is the value of the digits function. The default value of digits is 32."

[Star Strider]

Everything in a computer has a finite precision (consider evaluating any expression involving π otherwise), and beyond a certain number of significant digits, the results are usually not important. Note that the symbolic numbers use 32 decimal digits (that can be increased if necessary), while floating-point precision is limited to 16 decimal digits. Those are generally enough for most applications. The vpa function simply limits the display of an expression. Full precision, however you define it, is maintained intermally, so nothing is lost.

[M.Many]

Thanks for your answer, I do understand this.

But I do not agree with that sentence : 'Full precision, however you define it, is maintained intermally, so nothing is lost.'

Please have a look at the following code, that I use with my previous phi function evaluated at 1:

>> phi(1)-vpa(phi(1))
 
ans =
 
cos(978500^(1/2)/1000)*exp(-1/20) + (978500^(1/2)*sin(978500^(1/2)/1000)*exp(-1/20))/19570 - 0.562748412343890484348604524131142562316
 

Which gives the exact error we make by using vpa (because it is an expression where the cosines and square roots are not evaluated) and it is not 0. We can evaluate the error by applying the vpa function again :

>> vpa(phi(1)-vpa(phi(1)),100)
 
ans =
 
0.000000000000000000000000000000000000000000000001010157116631602241471408195383622793040240093593829057563725786616790742733299511920032906991801494

So full precision is not maintained internally

Again if we type the following :

>> phi(1)-phi(1)
 
ans =
 
0

we get an exact symbolic 0 because it is the result of comparing expressions and not evaluations of expressions.

Am I correct or am I missing something here ?

[Star Strider]

The vpa function converts fractions to their decimal-fraction equivalents, so there can be approximation errors in the conversion of finite fractions to finite decimal fractions. Nevertheless, the fractions retain their full internal precision, and the results of vpa retains its full internal precision. Because of the conversion, the original fractions and the decimal equivalents are not absolutely the same.

So will always be the same, however (however many finite terms it uses) will never be exactly equal to it.