%% Introduction to symbolic expressions in Matlab (Symbolic Math Toolbox)
% You should always start with this:

clearvars;          % clear all variables in the work space

% The Symbolic Math Toolbox allows to define symbolic expressions in terms of symbolic variables such as t,x,y. All results are exact symbolic results.
% If you type

format long g                      % show 15 digits of numerical results
[ sin(pi), sin(pi/3), sin(pi/7) ]

% Matlab uses a numerical approximation with 16 digits accuracy for \pi, and then computes numerical results.
% We can define

Pi = sym(pi)       % gives symbolic version of pi

% now

[ sin(Pi), sin(Pi/3), sin(Pi/7) ]

% gives symbolic results. We can then evaluate this using vpa (uses by default 32 digits):

vpa(ans)

% Note that we have to use sym(...) for integers to get symbolic results:

[ 1/3, exp(1) ]
[ 1/sym(3) , exp(sym(1)) ]

%% Defining a symbolic expression g depending on x
% Example: Define v=x\cdot\sin x

syms x              % define x as symbolic variable
v = x*sin(x);

% Note that v is not a function, but just denotes the expression x\cdot \sin x.
% Find the derivative   \displaystyle
% \frac{d}{dx}v

diff(v,x)

% Find the fifth derivative \displaystyle \left(\frac{d}{dx}\right)^5 v

diff(v,x,5)

% Find an antiderivative: \displaystyle  \int v\,dx

int(v,x)

% Use int(v,x,a,b) to find the definite integral \int_a^b v\,dx

int(v,x,0,Pi/3)    % we defined Pi = sym(pi) 

% Find   \displaystyle y =
% \left.
% v
% \right|_{\displaystyle x=3}, i.e., in v plug in x=3

y = subs(v,x,3)

% Find numerical value using vpa
% vpa(y) finds a numerical value. By default this shows 32 digits, you can change this with digits(n)

vpa(y)

%% Plot the symbolic expression v using fplot
% Plot the function g for 0\le x\le 20:

fplot(v,[0,20])
xlabel('x')                            % label x-axis
title('function g(x) = x\cdot sin x')  % title on top of graph
grid on                                % show grid

%% Find a solution x  of the nonlinear equation v=0 or v=a
% In our example there are infinitely many solutions.

solx = solve(v,x)

% This returns 0 which is definitely one solution.
% We can try to find more solutions:

[solx,param,conditions] = solve(v,x,'ReturnConditions',true)

% This shows that 0 and \pi k are solutions, where k\in\mathbb{Z} can be arbitrary. (Showing 0 separately is redundant as it is included in \pi k with k=0)
% This only works in simple equations.
% Let us try to solve x\cdot\sin x=5:

sol = solve(v==5,x)

% This is indeed an approximation of a solution. But we may be interested in the solution near x=7. With vpasolve we can specify an initial guess:

vpasolve(v==5,x,7)

%
%% Using an expression v of two variables x,y
% We can also use functions of more than one variable.
% Example: Define the function v=x^4+y^4-4(x^2+y^2)+4

syms x y                   % declare symbolic variables x,y
v = x^4+y^4-4*(x^2+y^2)+4;

% We can find partial derivatives:
% Find \partial_x v

diff(v,x)

% Find \displaystyle \partial_{xxy}v
% = (\frac\partial{\partial x})^2\frac\partial{\partial y} v

diff(v,x,x,y)

% Evaluate \displaystyle
% \left.
% v
% \right|_{x=1,y=2}, i.e. in v plug in x=1,y=2

subs(v,[x,y],[1,2])

% Plot the graph of v for -2\le x\le2, -2\le y\le 2 as a surface over the xy-plane:

fsurf(v,[-2 2 -2 2])
xlabel('x'); ylabel('y'); title('function x^4+y^4-4(x^2+y^2)+4')

% Make a contour plot of v for -2\le x\le2, -2\le y\le 2

fcontour(v,[-2 2 -2 2])
colorbar     % show colorbar on the right which explains colors of contours
axis equal   % use the same unit length for x and y
xlabel('x'); ylabel('y'); title('contour plot of x^4+y^4-4(x^2+y^2)+4')

% Here Matlab chooses the contour levels automatically. If you want more or fewer contours you can specify the spacing of the contours by using fcontour(v,[x1,x2,y1,y2],'LevelStep',dz):

fcontour(v,[-2 2 -2 2],'LevelStep',0.5)
colorbar; axis equal

% If you only want contours for a few specific z-values you can use

fcontour(v,[-2.1 2.1 -2.1 2.1],'LevelList',[0 1])
axis equal; 
title('points with v=0 (blue) and G(x,y)=1 (yellow)')

% The blue points satisfy v=0, the yellow points satisfy v=1
% If you only want to see the points with v=0 you can use fcontour(v,[x1,x2,y1,y2],'LevelList',[0]) or fimplicit(v,[x1,x2,y1,y2]):

fimplicit(v,[-2 2 -2 2]);
axis equal; xlabel('x'); ylabel('y'); title('points satisfying x^4+y^4-4(x^2+y^2)+4=0')

% NOTE: The plotting commands fplot, fsurf, fcontour, fimplicit work for both Matlab @-functions and symbolic expressions.
%
%