.gitignore

@@ -43,5 +43,4 @@
 src/*
 !src/Makefile
 !src/*.cpp
-!src/*.hpp
 !src/*.py

README.md

@@ -3,36 +3,3 @@
 ## Practical information
 
 - [Project](https://anderkve.github.io/FYS3150/book/projects/project1.html)
-
-## How to compile C++ code
-
-Make sure that you are inside the **src** directory before compiling the code.
-
-Now you can execute the command shown under to compile:
-
-```
-    make
-```
-
-This will create object files and link them together into 2 executable files.
-These files are called **main** and **analyticPlot**. 
-To run them, you can simply use the commands below:
-
-```
-    ./main
-```
-
-```
-    ./analyticPlot
-```
-
-## How to generate plots
-
-For generating the plots, there are 4 Python scripts. 
-You can run each one of them by using this command:
-
-```
-    python <PythonFile>
-```
-
-The plots will be saved inside **latex/images**.

latex/assignment_1.tex

@@ -80,7 +80,7 @@
 \begin{document}
 
 \title{Project 1}      % self-explanatory
-\author{Cory Alexander Balaton \& Janita Ovidie Sandtrøen Willumsen}          % self-explanatory
+\author{Cory Balaton \& Janita Willumsen}          % self-explanatory
 \date{\today}                             % self-explanatory
 \noaffiliation                            % ignore this, but keep it.
 
@@ -107,6 +107,4 @@
 
 \input{problems/problem9}
 
-\input{problems/problem10}
-
 \end{document}

latex/problems/problem1.tex

@@ -40,5 +40,5 @@ Using the values that we found for $c_1$ and $c_2$, we get
 
 \begin{align*}
     u(x) &= -e^{-10x} + (e^{-10} - 1) x + 1 \\
-         &= 1 - (1 - e^{-10})x - e^{-10x}
+         &= 1 - (1 - e^{-10}) - e^{-10x}
 \end{align*}

latex/problems/problem10.tex

@@ -1,8 +1,3 @@
 \section*{Problem 10}
 
 % Time and show result, and link to relevant files
-\begin{figure}[H]
-    \centering 
-    \includegraphics[width=0.7\linewidth]{images/problem10.pdf}
-    \caption{Timing of general algorithm vs. special for step sizes $n_{steps}$}
-\end{figure}

latex/problems/problem2.tex

@@ -6,8 +6,6 @@ Point generator code: https://github.uio.no/FYS3150-G2-203/Project-1/blob/main/s
 
 Plotting code: https://github.uio.no/FYS3150-G2-2023/Project-1/blob/main/src/analyticPlot.py
 
-\begin{figure}[H]
+Here is the plot of the analytical solution for $u(x)$.
-    \centering
+
-    \includegraphics[width=0.8\linewidth]{images/analytical_solution.pdf}
+\includegraphics[scale=.5]{analytical_solution.pdf}
-    \caption{Plot of the analytical solution $u(x)$.}
-\end{figure}

latex/problems/problem5.tex

@@ -1,6 +1,6 @@
 
 \section*{Problem 5}
 
-\subsection*{a \& b)}
+\subsection*{a)}
 
-$n = m - 2$ since when solving for $\vec{v}$, we are finding the solutions for all the points that are in between the boundaries and not the boundaries themselves. $\vec{v}^*$ on the other hand includes the boundary points.
+\subsection*{b)}

latex/problems/problem6.tex

@@ -17,7 +17,6 @@ Following Thomas algorithm for gaussian elimination, we first perform a forward
         \State $n \leftarrow$ length of $\vec{b}$
         \State $\vec{\hat{b}}$, $\vec{\hat{g}} \leftarrow$ vectors of length $n$.
         \State $\hat{b}_{1} \leftarrow b_{1}$ \Comment{Handle first element in main diagonal outside loop}
-        \State $\hat{g}_{1} \leftarrow g_{1}$
         \For{$i = 2, 3, ..., n$}
         \State $d \leftarrow \frac{a_{i}}{\hat{b}_{i-1}}$ \Comment{Calculating common expression}
         \State $\hat{b}_{i} \leftarrow b_{i} - d \cdot c_{i-1}$
@@ -44,4 +43,4 @@ Following Thomas algorithm for gaussian elimination, we first perform a forward
 Counting the number of FLOPs for the general algorithm by looking at one procedure at a time. 
 For every iteration of i in forward sweep we have 1 division, 2 multiplications, and 2 subtractions, resulting in $5(n-1)$ FLOPs. 
 For backward sweep we have 1 division, and for every iteration of i we have 1 subtraction, 1 multiplication, and 1 division, resulting in $3(n-1)+1$ FLOPs. 
-Total FLOPs for the general algorithm is $8(n-1)+1$.
+Total FLOPs for the general algorithm is $8(n-1)+1$.

latex/problems/problem7.tex

@@ -1,14 +1,3 @@
 \section*{Problem 7}
 
-\subsection*{a)}
 % Link to relevant files on gh and possibly add some comments
-The code can be found at https://github.uio.no/FYS3150-G2-2023/Project-1/blob/main/src/generalAlgorithm.cpp
-
-\subsection*{b)}
-Increasing the number of steps results in an approximation close to the exact solution.
-\begin{figure}[H]
-    \centering 
-    \includegraphics[width=0.8\linewidth]{images/problem7.pdf}
-    \caption{Plot showing the numeric solution of $u_{approx}$ for $n_{steps}$ and the exact solution $u_{exact}$.}
-\end{figure}
-

latex/problems/problem8.tex

@@ -1,27 +1,3 @@
 \section*{Problem 8}
 
 %link to relvant files and show plots
-\subsection*{a)}
-Increasing number of steps result in a decrease of absolute error.
-\begin{figure}[H]
-    \centering 
-    \includegraphics[width=0.8\linewidth]{images/problem8_a.pdf}
-    \caption{Absolute error for different step sizes $n_{steps}$.}
-\end{figure}
-
-\subsection*{b)}
-Increasing number of steps also result in a decrease of absolute error.
-\begin{figure}[H]
-    \centering 
-    \includegraphics[width=0.8\linewidth]{images/problem8_b.pdf}
-    \caption{Relative error for different step sizes $n_{steps}$.}
-\end{figure}
-
-\subsection*{c)}
-Increasing number of steps result in a decrease of maximum relative error, up to a certain number of steps. At $n_{steps} \approx 10^{5}$ the maximumrelative error increase. 
-This can be related to loss of numerical precicion when step size is small.
-\begin{figure}[H]
-    \centering 
-    \includegraphics[width=0.7\linewidth]{images/problem8_c.pdf}
-    \caption{Maximum relative error for each step sizes $n_{steps}$.}
-\end{figure}

latex/problems/problem9.tex

@@ -1,55 +1,3 @@
 \section*{Problem 9}
 
-\subsection*{a)}
+% Show the algorithm, then calculate FLOPs, then link to relevant files
-% Specialize algorithm
-The special algorithm does not require the values of all $a_{i}$, $b_{i}$, $c_{i}$. 
-We find the values of $\hat{b}_{i}$ from simplifying the general case
-\begin{align*}
-    \hat{b}_{i} &= b_{i} - \frac{a_{i} \cdot c_{i-1}}{\hat{b}_{i-1}} \\
-    \hat{b}_{i} &= 2 - \frac{1}{\hat{b}_{i-1}}
-\end{align*} 
-Calculating the first values to see a pattern 
-\begin{align*}
-    \hat{b}_{1} &= 2 \\
-    \hat{b}_{2} &= 2 - \frac{1}{2} = \frac{3}{2} \\
-    \hat{b}_{3} &= 2 - \frac{1}{\frac{3}{2}} = \frac{4}{3} \\
-    \hat{b}_{4} &= 2 - \frac{1}{\frac{4}{3}} = \frac{5}{4} \\
-    \vdots & \\
-    \hat{b}_{i} &= \frac{i+1}{i} && \text{for $i = 1, 2, ..., n$}
-\end{align*}
-
-
-\begin{algorithm}[H]
-    \caption{Special algorithm}\label{algo:special}
-    \begin{algorithmic}
-        \Procedure{Forward sweep}{$\vec{b}$}
-        \State $n \leftarrow$ length of $\vec{b}$
-        \State $\vec{\hat{b}}$, $\vec{\hat{g}} \leftarrow$ vectors of length $n$.
-        \State $\hat{b}_{1} \leftarrow 2$ \Comment{Handle first element in main diagonal outside loop}
-        \State $\hat{g}_{1} \leftarrow g_{1}$
-        \For{$i = 2, 3, ..., n$}
-        \State $\hat{b}_{i} \leftarrow \frac{i+1}{i}$
-        \State $\hat{g}_{i} \leftarrow g_{i} + \frac{\hat{g}_{i-1}}{\hat{b}_{i-1}}$
-        \EndFor
-        \Return $\vec{\hat{b}}$, $\vec{\hat{g}}$
-        \EndProcedure
-
-        \Procedure{Backward sweep}{$\vec{\hat{b}}$, $\vec{\hat{g}}$}
-        \State $n \leftarrow$ length of $\vec{\hat{b}}$
-        \State $\vec{v} \leftarrow$ vector of length $n$.
-        \State $v_{n} \leftarrow \frac{\hat{g}_{n}}{\hat{b}_{n}}$
-        \For{$i = n-1, n-2, ..., 1$}
-        \State $v_{i} \leftarrow \frac{\hat{g}_{i} + v_{i+1}}{\hat{b}_{i}}$ 
-        \EndFor
-        \Return $\vec{v}$
-        \EndProcedure
-    \end{algorithmic}
-\end{algorithm}
-
-
-\subsection*{b)}
-% Find FLOPs
-For every iteration of i in forward sweep we have 2 divisions, and 2 additions, resulting in $4(n-1)$ FLOPs. 
-For backward sweep we have 1 division, and for every iteration of i we have 1 addition, and 1 division, resulting in $2(n-1)+1$ FLOPs. 
-Total FLOPs for the special algorithm is $6(n-1)+1$.
-

src/Makefile

@@ -1,30 +1,17 @@
 CC=g++
 
-CCFLAGS= -std=c++11
+.PHONY: clean
 
-OBJS=generalAlgorithm.o specialAlgorithm.o funcs.o
+all: simpleFile analyticPlot
 
-EXEC=main analyticPlot
+simpleFile: simpleFile.o
-
+	$(CC) -o $@ $^
-.PHONY: clean create_dirs
-
-all: create_dirs $(EXEC)
-
-main: main.o $(OBJS)
-	$(CC) $(CCFLAGS) -o $@ $^
 
 analyticPlot: analyticPlot.o
-	$(CC) $(CCFLAGS) -o $@ $^
+	$(CC) -o $@ $^
 
 %.o: %.cpp
-	$(CC) $(CCFLAGS) -c -o $@ $^
+	$(CC) -c $< -o $@
 
 clean:
 	rm *.o
-	rm $(EXEC)
-	rm -r output
-
-create_dirs:
-	mkdir -p output/general
-	mkdir -p output/special
-	mkdir -p output/error

src/analyticPlot.cpp

@@ -7,7 +7,7 @@
 #include <iomanip>
 
 #define RANGE 1000
-#define FILENAME "output/analytical_solution.txt"
+#define FILENAME "analytical_solution.txt"
 
 double u(double x);
 void generate_range(std::vector<double> &vec, double start, double stop, int n);

src/analyticPlot.py

@@ -2,11 +2,11 @@ import numpy as np
 import matplotlib.pyplot as plt
 
 def main():
-    FILENAME = "../latex/images/analytical_solution.pdf"
+    FILENAME = "analytical_solution.pdf"
     x = []
     v = []
 
-    with open('output/analytical_solution.txt') as f:
+    with open('analytical_solution.txt') as f:
         for line in f:
             a, b = line.strip().split()
             x.append(float(a))

src/funcs.cpp

@@ -1,38 +0,0 @@
-#include "funcs.hpp"
-
-double f(double x) {
-    return 100*std::exp(-10*x);
-}
-
-double u(double x) {
-    return 1. - (1. - std::exp(-10.))*x - std::exp(-10.*x);
-}
-
-void build_g_vec(int n_steps, arma::vec& g_vec) {
-    g_vec.resize(n_steps-1);
-
-    double step_size = 1./ (double) n_steps;
-    for (int i=0; i < n_steps-1; i++) {
-        g_vec(i) = step_size*step_size*f((i+1)*step_size);
-    }
-}
-
-void build_arrays(
-    int n_steps, 
-    arma::vec& sub_diag, 
-    arma::vec& main_diag, 
-    arma::vec& sup_diag,
-    arma::vec& g_vec
-) 
-{
-    sub_diag.resize(n_steps-2);
-    main_diag.resize(n_steps-1);
-    sup_diag.resize(n_steps-2);
-
-    sub_diag.fill(-1);
-    main_diag.fill(2);
-    sup_diag.fill(-1);
-
-    build_g_vec(n_steps, g_vec);
-}
-

src/funcs.hpp

@@ -1,24 +0,0 @@
-#ifndef __FUNCS__
-#define __FUNCS__
-
-#include <armadillo>
-#include <cmath>
-
-#define PRECISION 8
-
-#define N_STEPS_EXP 7
-
-double f(double x);
-
-double u(double x);
-
-void build_g_vec(int n_steps, arma::vec& g_vec);
-
-void build_arrays(
-    int n_steps, 
-    arma::vec& sub_diag, 
-    arma::vec& main_diag, 
-    arma::vec& sup_diag,
-    arma::vec& g_vec
-);
-#endif

src/generalAlgorithm.cpp

@@ -1,84 +0,0 @@
-#include "funcs.hpp"
-#include "generalAlgorithm.hpp"
-#include <cmath>
-
-arma::vec& general_algorithm(
-    arma::vec& sub_diag,
-    arma::vec& main_diag,
-    arma::vec& sup_diag,
-    arma::vec& g_vec
-)
-{
-    int n = g_vec.n_elem;
-    double d;
-
-    for (int i = 1; i < n; i++) {
-        d = sub_diag(i-1) / main_diag(i-1);
-        main_diag(i) -= d*sup_diag(i-1);
-        g_vec(i) -= d*g_vec(i-1);
-    }
-
-    g_vec(n-1) /= main_diag(n-1);
-
-    for (int i = n-2; i >= 0; i--) {
-        g_vec(i) = (g_vec(i) - sup_diag(i) * g_vec(i+1)) / main_diag(i);
-    }
-    return g_vec;
-}
-
-void general_algorithm_main() 
-{
-    arma::vec sub_diag, main_diag, sup_diag, g_vec, v_vec;
-    std::ofstream ofile;
-    int steps;
-    double step_size;
-
-    for (int i = 0; i < N_STEPS_EXP; i++) {
-        steps = std::pow(10, i+1);
-        step_size = 1./(double) steps;
-
-        build_arrays(steps, sub_diag, main_diag, sup_diag, g_vec);
-
-        v_vec = general_algorithm(sub_diag, main_diag, sup_diag, g_vec);
-
-        ofile.open("output/general/out_" + std::to_string(steps) + ".txt");
-
-        for (int j=0; j < v_vec.n_elem; j++) {
-            ofile << std::setprecision(PRECISION) << std::scientific << step_size*(j+1) << ","
-                << std::setprecision(PRECISION) << std::scientific << v_vec(j) << std::endl;
-        }
-        ofile.close();
-    }
-}
-
-void general_algorithm_error() 
-{
-    arma::vec sub_diag, main_diag, sup_diag, g_vec, v_vec;
-    std::ofstream ofile;
-    int steps;
-    double step_size, abs_err, rel_err, u_i, v_i;
-
-    for (int i=0; i < N_STEPS_EXP; i++) {
-        steps = std::pow(10, i+1);
-        step_size = 1./(double) steps;
-
-        build_arrays(steps, sub_diag, main_diag, sup_diag, g_vec);
-
-        v_vec = general_algorithm(sub_diag, main_diag, sup_diag, g_vec);
-
-        ofile.open("output/error/out_" + std::to_string(steps) + ".txt");
-
-        for (int j=0; j < v_vec.n_elem; j++) {
-            u_i = u(step_size*(j+1));
-            v_i = v_vec(j);
-            abs_err = u_i - v_i;
-            ofile << std::setprecision(PRECISION) << std::scientific 
-                << step_size*(j+1) << ","
-                << std::setprecision(PRECISION) << std::scientific 
-                << std::log10(std::abs(abs_err)) << ","
-                << std::setprecision(PRECISION) << std::scientific 
-                << std::log10(std::abs(abs_err/u_i)) << std::endl;
-        }
-        ofile.close();
-    }
-}

src/generalAlgorithm.hpp

@@ -1,18 +0,0 @@
-#ifndef __GENERAL_ALG__
-#define __GENERAL_ALG__
-
-#include <armadillo>
-#include <iomanip>
-
-arma::vec& general_algorithm(
-    arma::vec& sub_diag,
-    arma::vec& main_diag,
-    arma::vec& sup_diag,
-    arma::vec& g_vec
-);
-
-void general_algorithm_main();
-
-void general_algorithm_error();
-
-#endif

src/main.cpp

@@ -1,68 +1,24 @@
+#include "GeneralAlgorithm.hpp"
 #include <armadillo>
-#include <cmath>
+#include <iostream>
-#include <ctime>
-#include <fstream>
-#include <iomanip>
-#include <ios>
-#include <string>
 
-#include "funcs.hpp"
+double f(double x) {
-#include "generalAlgorithm.hpp"
+    return 100. * std::exp(-10.*x);
-#include "specialAlgorithm.hpp"
-
-#define TIMING_ITERATIONS 5
-
-void timing() {
-    arma::vec sub_diag, main_diag, sup_diag, g_vec;
-    int n_steps;
-
-    std::ofstream ofile;
-    ofile.open("output/timing.txt");
-
-    // Timing
-    for (int i=1; i < N_STEPS_EXP; i++) {
-        n_steps = std::pow(10, i);
-        clock_t g_1, g_2, s_1, s_2;
-        double g_res = 0, s_res = 0;
-      
-        // Repeat a number of times to take an average
-        for (int j=0; j < TIMING_ITERATIONS; j++) {
-            
-            build_arrays(n_steps, sub_diag, main_diag, sup_diag, g_vec);
-
-            g_1 = clock();
-
-            general_algorithm(sub_diag, main_diag, sup_diag, g_vec);
-
-            g_2 = clock();
-
-            g_res += (double) (g_2 - g_1) / CLOCKS_PER_SEC;
-            // Rebuild g_vec for the special alg
-            build_g_vec(n_steps, g_vec);
-
-            s_1 = clock();
-
-            special_algorithm(-1., 2., -1., g_vec);
-
-            s_2 = clock();
-
-            s_res += (double) (s_2 - s_1) / CLOCKS_PER_SEC;
-
-        }
-        // Write the average time to file
-        ofile 
-            << n_steps << ","
-            << g_res / (double) TIMING_ITERATIONS << ","
-            << s_res / (double) TIMING_ITERATIONS << std::endl;
-    }
-
-    ofile.close();
 }
 
-int main() 
+double a_sol(double x) {
-{
+    return 1. - (1. - std::exp(-10)) * x - std::exp(-10*x);
-    timing();
+}
-    general_algorithm_main();
+
-    general_algorithm_error();
+int main() {
-    special_algorithm_main();
+    arma::mat A = arma::eye(3,3);
+
+    GeneralAlgorithm ga(3, &A, f, a_sol, 0., 1.);
+
+    ga.solve();
+    std::cout << "Time: " << ga.time(5) << std::endl;
+    ga.error();
+
+    return 0;
+
 }

src/plot_general_alg.py

@@ -1,30 +0,0 @@
-import matplotlib.pyplot as plt
-import numpy as np
-
-analytical_func = lambda x: 1 - (1 - np.exp(-10))*x - np.exp(-10*x)
-
-def main():
-    for i in range(7):
-        x = []
-        y = []
-        x.append(0.)
-        y.append(0.)
-        with open(f"output/general/out_{10**(i+1)}.txt", "r") as f:
-            lines = f.readlines()
-            for line in lines:
-                x_i, y_i = line.strip().split(",")
-                x.append(float(x_i))
-                y.append(float(y_i))
-
-        x.append(1.)
-        y.append(0.)
-
-        plt.plot(x, y, label=f"n$_{{steps}} = 10^{i+1}$")
-
-    x = np.linspace(0, 1, 1001)
-    plt.plot(x, analytical_func(x), label=f"u$_{{exact}}$")
-    plt.legend()
-    plt.savefig("../latex/images/problem7.pdf")
-
-if __name__ == "__main__":
-    main()

src/plot_general_alg_error.py

@@ -1,40 +0,0 @@
-import matplotlib.pyplot as plt
-
-# plt.rc('text', usetex=True)
-# plt.rc('font', family='serif')
-
-def main():
-    for i in range(7):
-        x = []
-        abs_err = []
-        rel_err = []
-        with open(f"output/error/out_{10**(i+1)}.txt", "r") as f:
-            lines = f.readlines()
-            for line in lines:
-                x_i, abs_err_i, rel_err_i = line.strip().split(",")
-                x.append(float(x_i))
-                abs_err.append(float(abs_err_i))
-                rel_err.append(float(rel_err_i))
-
-        plt.figure(1)
-        plt.plot(x, abs_err, label=f"n$_{{steps}} = 10^{i+1}$")
-        plt.figure(2)
-        plt.plot(x, rel_err, label=f"n$_{{steps}} = 10^{i+1}$")
-
-        plt.figure(3)
-        plt.plot(i+1, max(rel_err), marker="o", markersize=10)
-    
-    plt.figure(1)
-    plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)
-    plt.figure(2)
-    plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)
-    
-    plt.figure(1)
-    plt.savefig("../latex/images/problem8_a.pdf", bbox_inches="tight")
-    plt.figure(2)
-    plt.savefig("../latex/images/problem8_b.pdf", bbox_inches="tight")
-    plt.figure(3)
-    plt.savefig("../latex/images/problem8_c.pdf", bbox_inches="tight")
-
-if __name__ == "__main__":
-    main()

src/simpleFile.cpp

@@ -0,0 +1,162 @@
+#include <armadillo>
+#include <cmath>
+#include <ctime>
+#include <fstream>
+#include <iomanip>
+#include <ios>
+#include <string>
+
+#define TIMING_ITERATIONS 5
+
+arma::vec* general_algorithm(
+    arma::vec* sub_diag,
+    arma::vec* main_diag,
+    arma::vec* sup_diag,
+    arma::vec* g_vec
+)
+{
+    int n = g_vec->n_elem;
+    double d;
+
+    for (int i = 1; i < n; i++) {
+        d = (*sub_diag)(i-1) / (*main_diag)(i-1);
+        (*main_diag)(i) -= d*(*sup_diag)(i-1);
+        (*g_vec)(i) -= d*(*g_vec)(i-1);
+    }
+
+    (*g_vec)(n-1) /= (*main_diag)(n-1);
+
+    for (int i = n-2; i >= 0; i--) {
+        (*g_vec)(i) = ((*g_vec)(i) - (*sup_diag)(i) * (*g_vec)(i+1)) / (*main_diag)(i);
+    }
+    return g_vec;
+}
+
+
+arma::vec* special_algorithm(
+    double sub_sig,
+    double main_sig,
+    double sup_sig,
+    arma::vec* g_vec
+)
+{
+    int n = g_vec->n_elem;
+    arma::vec diag = arma::vec(n);
+
+    for (int i = 1; i < n; i++) {
+        // Calculate values for main diagonal based on indices
+        diag(i-1) = (double)(i+1) / i;
+        (*g_vec)(i) += (*g_vec)(i-1) / diag(i-1);
+    }
+    // The last element in main diagonal has value (i+1)/i = (n+1)/n
+    (*g_vec)(n-1) /= (double)(n+1) / (n);
+
+    for (int i = n-2; i >= 0; i--) {
+        (*g_vec)(i) = ((*g_vec)(i) + (*g_vec)(i+1))/ diag(i);
+    }
+
+    return g_vec;
+}
+
+void error(
+    std::string filename,
+    arma::vec* x_vec,
+    arma::vec* v_vec,
+    arma::vec* a_vec
+)
+{
+    std::ofstream ofile;
+    ofile.open(filename);
+
+    if (!ofile.is_open()) {
+        exit(1);
+    }
+
+    for (int i=0; i < a_vec->n_elem; i++) {
+        double sub = (*a_vec)(i) - (*v_vec)(i);
+        ofile << std::setprecision(8) << std::scientific << (*x_vec)(i)
+            << std::setprecision(8) << std::scientific << std::log10(std::abs(sub))
+            << std::setprecision(8) << std::scientific << std::log10(std::abs(sub/(*a_vec)(i)))
+            << std::endl;
+    }
+    
+    ofile.close();
+
+}
+
+double f(double x) {
+    return 100*std::exp(-10*x);
+}
+
+void build_array(
+    int n_steps, 
+    arma::vec* sub_diag, 
+    arma::vec* main_diag, 
+    arma::vec* sup_diag,
+    arma::vec* g_vec
+) 
+{
+    sub_diag->resize(n_steps-2);
+    main_diag->resize(n_steps-1);
+    sup_diag->resize(n_steps-2);
+
+    sub_diag->fill(-1);
+    main_diag->fill(2);
+    sup_diag->fill(-1);
+
+    g_vec->resize(n_steps-1);
+
+    double step_size = 1./ (double) n_steps;
+    for (int i=0; i < n_steps-1; i++) {
+        (*g_vec)(i) = f((i+1)*step_size);
+    }
+
+}
+
+void timing() {
+    arma::vec sub_diag, main_diag, sup_diag, g_vec;
+    int n_steps;
+
+    std::ofstream ofile;
+    ofile.open("timing.txt");
+
+    // Timing
+    for (int i=1; i <= 8; i++) {
+        n_steps = std::pow(10, i);
+        clock_t g_1, g_2, s_1, s_2;
+        double g_res = 0, s_res = 0;
+      
+        for (int j=0; j < TIMING_ITERATIONS; j++) {
+            build_array(n_steps, &sub_diag, &main_diag, &sup_diag, &g_vec);
+
+            g_1 = clock();
+
+            general_algorithm(&sub_diag, &main_diag, &sup_diag, &g_vec);
+
+            g_2 = clock();
+
+            g_res += (double) (g_2 - g_1) / CLOCKS_PER_SEC;
+            build_array(n_steps, &sub_diag, &main_diag, &sup_diag, &g_vec);
+
+            s_1 = clock();
+
+            special_algorithm(-1., 2., -1., &g_vec);
+
+            s_2 = clock();
+
+            s_res += (double) (s_2 - s_1) / CLOCKS_PER_SEC;
+
+        }
+        ofile 
+            << n_steps << ","
+            << g_res / (double) TIMING_ITERATIONS << ","
+            << s_res / (double) TIMING_ITERATIONS << std::endl;
+    }
+
+    ofile.close();
+}
+
+int main() 
+{
+    timing();
+}

src/specialAlgorithm.cpp

@@ -1,52 +0,0 @@
-#include "funcs.hpp"
-#include "specialAlgorithm.hpp"
-
-arma::vec& special_algorithm(
-    double sub_sig,
-    double main_sig,
-    double sup_sig,
-    arma::vec& g_vec
-)
-{
-    int n = g_vec.n_elem;
-    arma::vec diag = arma::vec(n);
-
-    for (int i = 1; i < n; i++) {
-        // Calculate values for main diagonal based on indices
-        diag(i-1) = (double)(i+1) / i;
-        g_vec(i) += g_vec(i-1) / diag(i-1);
-    }
-    // The last element in main diagonal has value (i+1)/i = (n+1)/n
-    g_vec(n-1) /= (double)(n+1) / (n);
-
-    for (int i = n-2; i >= 0; i--) {
-        g_vec(i) = (g_vec(i) + g_vec(i+1))/ diag(i);
-    }
-
-    return g_vec;
-}
-
-void special_algorithm_main() 
-{
-    arma::vec g_vec, v_vec;
-    std::ofstream ofile;
-    int steps;
-    double sub_sig, main_sig, sup_sig, step_size;
-
-    for (int i = 0; i < N_STEPS_EXP; i++) {
-        steps = std::pow(10, i+1);
-        step_size = 1./(double) steps;
-        build_g_vec(steps, g_vec);
-
-        v_vec = special_algorithm(sub_sig, main_sig, sup_sig, g_vec);
-
-        ofile.open("output/special/out_" + std::to_string(steps) + ".txt");
-
-        for (int j=0; j < v_vec.n_elem; j++) {
-            ofile << std::setprecision(PRECISION) << std::scientific << step_size*(j+1) << ","
-                << std::setprecision(PRECISION) << std::scientific << v_vec(j) << std::endl;
-        }
-        ofile.close();
-    }
-}
-

src/specialAlgorithm.hpp

@@ -1,16 +0,0 @@
-#ifndef __SPECIAL_ALG__
-#define __SPECIAL_ALG__
-
-#include <armadillo>
-#include <iomanip>
-
-arma::vec& special_algorithm(
-    double sub_sig,
-    double main_sig,
-    double sup_sig,
-    arma::vec& g_vec
-);
-
-void special_algorithm_main();
-
-#endif

src/timing.py

@@ -1,22 +0,0 @@
-import matplotlib.pyplot as plt
-
-def main():
-    x = []
-    gen_alg = []
-    spec_alg = []
-    with open(f"output/timing.txt", "r") as f:
-        lines = f.readlines()
-        for line in lines:
-            x_i, gen_i, spec_i = line.strip().split(",")
-            x.append(float(x_i))
-            gen_alg.append(float(gen_i))
-            spec_alg.append(float(spec_i))
-
-        plt.plot(x, gen_alg, label=f"General algorithm")
-        plt.plot(x, spec_alg, label=f"Special algorithm")
-
-    plt.legend()
-    plt.savefig("../latex/images/problem10.pdf")
-
-if __name__ == "__main__":
-    main()