Numerical Algorithms Test Suite

These are symbolically-solved problems with analytical solutions (and in particular, initial-value problems) aimed at testing Project Elara’s numerical computing libraries. By comparing hand-solved exact solutions with the numerical results, we can determine the accuracy and robustness of the numerical algorithms we have implemented.

ODE Problem 1

IVP to solve:

$$\begin{array}{c}\frac{dy}{dt}=\frac{2y(t)-1}{t+3},y(0)=-1\end{array}$$

Analytical solution:

$$\begin{array}{rl}\frac{dy}{dt}=\frac{2y(t)-1}{t+3}& \Rightarrow \frac{dy}{2y(t)-1}=\frac{dt}{t+3}\\ & \Rightarrow \int \frac{dy}{2y(t)-1}=\int \frac{dt}{t+3}\\ & \Rightarrow \frac{1}{2}\ln |2y(t)-1|=\ln |t+3|+C\end{array}$$ $$\begin{array}{rl}\frac{1}{2}\ln |2y(t)-1|=\ln |t+3|+C& \Rightarrow \frac{1}{2}\ln |2(-1)-1|=\ln |0+3|+C\\ & \Rightarrow \frac{1}{2}\ln |-3|=\ln |3|+C\\ & \Rightarrow \frac{1}{2}\ln (3)=\ln (3)+C\\ & \Rightarrow -\frac{1}{2}\ln (3)=C\end{array}$$ $$\begin{array}{rl}\frac{1}{2}\ln |2y(t)-1|=\ln |t+3|-\frac{1}{2}\ln (3)& \\ & \Rightarrow \ln |2y(t)-1|=2\ln |t+3|-\ln (3)\\ & \Rightarrow \ln |2y(t)-1|=\ln |(t+3{)}^2|-\ln (3)\\ & \Rightarrow \ln |2y(t)-1|=\ln |\frac{(t+3{)}^2}{3}|\\ & \Rightarrow 2y(t)-1=\frac{(t+3{)}^2}{3}\\ & \Rightarrow y(t)=\frac{(t+3{)}^2}{6}-\frac{1}{2}\end{array}$$

ODE Problem 2

IVP to solve:

$$\begin{array}{c}(1+t^2)y^{\prime }(t)-2ty(t)=t,y(0)=0\end{array}$$

Analytical solution (with steps):

$$\begin{array}{rl}(1+t^2)y^{\prime }(t)-2ty(t)=t& \Rightarrow y^{\prime }(t)-\frac{2t}{1+t^2}y(t)=\frac{t}{1+t^2}\end{array}$$ $$\begin{array}{c}p(t)=-\frac{2t}{1+t^2},q(t)=\frac{t}{1+t^2}\end{array}$$ $$\begin{array}{rl}\mu & =\exp (\int p(t)dt)\\ & =\exp (\int -\frac{2t}{1+t^2}dt)\\ & =\exp (-\ln (1+t^2))\\ & =\exp (\ln (\frac{1}{1+t^2}))\\ & =\frac{1}{1+t^2}\end{array}$$ $$\begin{array}{rl}y(t)& =\frac{1}{\mu }[\int \mu q(t)dt+C]\\ & =(1+t^2)[\int \frac{1}{(1+t^2{)}^2}dt+C]\\ & =(1+t^2)[\frac{1}{2(1+t^2)}\ln |(1+t^2{)}^2|+C]\\ & =\frac{1}{2}\ln |(1+t^2)|+C(1+t^2)\end{array}$$ $$\begin{array}{c}y(t)=\frac{1}{2}\ln |(1+t^2)|\end{array}$$

ODE Problem 3

$$\begin{array}{c}{y}^{\prime }-\frac{2y}{x}=-x^2{y}^2,y(a)=b\end{array}$$ $$\begin{array}{rl}{y}^{\prime }-\frac{2y}{x}=-x^2{y}^2& \Rightarrow \frac{y^{\prime }}{y^2}-\frac{2}{xy}=-x^2\end{array}$$

Using $u=\frac{1}{y},u^{\prime }=\frac{-y^{\prime }}{y^2}$

$$\begin{array}{rl}\frac{y^{\prime }}{y^2}-\frac{2}{xy}=-x^2& \Rightarrow -u^{\prime }-\frac{2}{x}u=-x^2\\ & \Rightarrow u^{\prime }+\frac{2}{x}u=x^2\end{array}$$ $$\begin{array}{rl}\mu & =\exp (2\int \frac{1}{x}dx)\\ & =\exp (2\ln |x|dx)\\ & =x^2\end{array}$$ $$\begin{array}{rl}u& =\frac{1}{x^2}(\int x^2\ast x^{2}dx+C)\\ & =\frac{1}{x^2}(\int x^{4}dx+C)\\ & =\frac{1}{x^2}(\frac{x^5}{5}+C)\\ & =\frac{x^3}{5}+\frac{C}{x^2}\end{array}$$ $$\begin{array}{rl}u=\frac{x^3}{5}+\frac{C}{x^2}& \Rightarrow \frac{1}{y}=\frac{x^3}{5}+\frac{C}{x^2}\\ & \Rightarrow y=\frac{1}{\frac{x^3}{5}+\frac{C}{x^2}}\\ & \Rightarrow y=\frac{5x^2}{x^5+C}\end{array}$$ $$\begin{array}{rl}y=\frac{5x^2}{x^5+C}& \Rightarrow y(a)=\frac{5a^2}{a^5+C}=b\\ & \Rightarrow 5a^2=b(a^5+C)=b{a}^5+Cb\\ & \Rightarrow 5a^2-b{a}^5=Cb\\ & \Rightarrow \frac{5a^5}{b}-a^5=C\end{array}$$ $$\begin{array}{c}y=\frac{5x^2}{x^5+\frac{5a^5}{b}-a^5}\end{array}$$

ODE Problem 4

This solves the generalized Bernoulli differential equation, which is given by:

$$\begin{array}{c}{y}^{\prime }+P(x)y=Q(x)y^n\end{array}$$

The solution is as follows:

$$\begin{array}{rl}{y}^{\prime }+P(x)y=Q(x)y^n& \Rightarrow \frac{y^{\prime }}{y^n}+P(x)\frac{1}{y^{n-1}}=Q(x)\end{array}$$

Using $u=\frac{1}{y^{n-1}},u^{\prime }=-\frac{y^{\prime }}{y^n}$

$$\begin{array}{rl}\frac{y^{\prime }}{y^n}+P(x)\frac{1}{y^{n-1}}=Q(x)& \Rightarrow -u^{\prime }+P(x)u=Q(x)\\ & \Rightarrow u^{\prime }-P(x)u=-Q(x)\end{array}$$

We define the integrating factor $\mu =\mu (x)$ as follows:

$$\begin{array}{c}\mu =\exp (-\int P(x)dx)\end{array}$$

Therefore:

$$\begin{array}{rl}u=\exp (\int P(x)dx)(\int P(x)Q(x)dx+C)& \Rightarrow \frac{1}{y^{n-1}}=\exp (\int P(x)dx)(\int P(x)Q(x)dx+C)\\ & \Rightarrow y^{n-1}=\frac{1}{\exp (\int P(x)dx)(\int P(x)Q(x)dx+C)}\\ & \Rightarrow y^{n-1}=\frac{\exp (-\int P(x)dx)}{(\int P(x)Q(x)dx+C)}\\ & \Rightarrow y=\sqrt[n-1]{\frac{\exp (-\int P(x)dx)}{(\int P(x)Q(x)dx+C)}}\end{array}$$