Compare commits
No commits in common. "main" and "9-solve-problem-9" have entirely different histories.
main
...
9-solve-pr
1
.gitignore
vendored
1
.gitignore
vendored
@ -43,5 +43,4 @@
|
||||
src/*
|
||||
!src/Makefile
|
||||
!src/*.cpp
|
||||
!src/*.hpp
|
||||
!src/*.py
|
||||
|
||||
33
README.md
33
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**.
|
||||
|
||||
Binary file not shown.
@ -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}
|
||||
|
||||
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
@ -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*}
|
||||
|
||||
@ -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}
|
||||
@ -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]
|
||||
\centering
|
||||
\includegraphics[width=0.8\linewidth]{images/analytical_solution.pdf}
|
||||
\caption{Plot of the analytical solution $u(x)$.}
|
||||
\end{figure}
|
||||
Here is the plot of the analytical solution for $u(x)$.
|
||||
|
||||
\includegraphics[scale=.5]{analytical_solution.pdf}
|
||||
|
||||
@ -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)}
|
||||
|
||||
@ -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}
|
||||
|
||||
|
||||
@ -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}
|
||||
25
src/Makefile
25
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
|
||||
|
||||
.PHONY: clean create_dirs
|
||||
|
||||
all: create_dirs $(EXEC)
|
||||
|
||||
main: main.o $(OBJS)
|
||||
$(CC) $(CCFLAGS) -o $@ $^
|
||||
simpleFile: simpleFile.o
|
||||
$(CC) -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
|
||||
|
||||
@ -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);
|
||||
|
||||
@ -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))
|
||||
@ -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);
|
||||
}
|
||||
|
||||
@ -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
|
||||
@ -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();
|
||||
}
|
||||
}
|
||||
@ -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
|
||||
78
src/main.cpp
78
src/main.cpp
@ -1,68 +1,24 @@
|
||||
#include "GeneralAlgorithm.hpp"
|
||||
#include <armadillo>
|
||||
#include <cmath>
|
||||
#include <ctime>
|
||||
#include <fstream>
|
||||
#include <iomanip>
|
||||
#include <ios>
|
||||
#include <string>
|
||||
#include <iostream>
|
||||
|
||||
#include "funcs.hpp"
|
||||
#include "generalAlgorithm.hpp"
|
||||
#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;
|
||||
double f(double x) {
|
||||
return 100. * std::exp(-10.*x);
|
||||
}
|
||||
|
||||
ofile.close();
|
||||
double a_sol(double x) {
|
||||
return 1. - (1. - std::exp(-10)) * x - std::exp(-10*x);
|
||||
}
|
||||
|
||||
int main()
|
||||
{
|
||||
timing();
|
||||
general_algorithm_main();
|
||||
general_algorithm_error();
|
||||
special_algorithm_main();
|
||||
int 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;
|
||||
|
||||
}
|
||||
|
||||
@ -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()
|
||||
@ -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()
|
||||
0
src/problem2.cpp
Normal file
0
src/problem2.cpp
Normal file
162
src/simpleFile.cpp
Normal file
162
src/simpleFile.cpp
Normal file
@ -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();
|
||||
}
|
||||
@ -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();
|
||||
}
|
||||
}
|
||||
|
||||
@ -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
|
||||
@ -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()
|
||||
Loading…
Reference in New Issue
Block a user