Make some changes
This commit is contained in:
parent
85e469f101
commit
54f5d4ab5d
@ -1,17 +1,19 @@
|
|||||||
CC=g++
|
CC=g++
|
||||||
|
|
||||||
|
OBJS= generalAlgorithm.o specialAlgorithm.o
|
||||||
|
|
||||||
.PHONY: clean
|
.PHONY: clean
|
||||||
|
|
||||||
all: simpleFile analyticPlot
|
all: main analyticPlot
|
||||||
|
|
||||||
simpleFile: simpleFile.o
|
main: main.o $(OBJS)
|
||||||
$(CC) -o $@ $^
|
$(CC) -o $@ $^
|
||||||
|
|
||||||
analyticPlot: analyticPlot.o
|
analyticPlot: analyticPlot.o
|
||||||
$(CC) -o $@ $^
|
$(CC) -o $@ $^
|
||||||
|
|
||||||
%.o: %.cpp
|
%.o: %.cpp
|
||||||
$(CC) -c $< -o $@
|
$(CC) -c -o $@ $^
|
||||||
|
|
||||||
clean:
|
clean:
|
||||||
rm *.o
|
rm *.o
|
||||||
|
|||||||
61
src/generalAlgorithm.cpp
Normal file
61
src/generalAlgorithm.cpp
Normal file
@ -0,0 +1,61 @@
|
|||||||
|
#include "generalAlgorithm.hpp"
|
||||||
|
#include <cmath>
|
||||||
|
#include <fstream>
|
||||||
|
#include <iomanip>
|
||||||
|
#include <ios>
|
||||||
|
#include <string>
|
||||||
|
|
||||||
|
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 < 6; i++) {
|
||||||
|
steps = std::pow(10, i+1);
|
||||||
|
step_size = 1./(double) steps;
|
||||||
|
main_diag.resize(steps - 1);
|
||||||
|
g_vec.resize(steps - 1);
|
||||||
|
sub_diag.resize(steps - 2);
|
||||||
|
sup_diag.resize(steps - 2);
|
||||||
|
|
||||||
|
sub_diag.fill(-1.);
|
||||||
|
sup_diag.fill(-1.);
|
||||||
|
main_diag.fill(2.);
|
||||||
|
|
||||||
|
v_vec = general_algorithm(&sub_diag, &main_diag, &sup_diag, &g_vec);
|
||||||
|
|
||||||
|
ofile.open("general_algorithm_" + std::to_string(steps) + ".txt");
|
||||||
|
|
||||||
|
for (int j=0; j < v_vec->n_elem; j++) {
|
||||||
|
ofile << std::setprecision(4) << std::scientific << step_size*(i+1) << ","
|
||||||
|
<< std::setprecision(4) << std::scientific << (*v_vec)(i) << std::endl;
|
||||||
|
}
|
||||||
|
ofile.close();
|
||||||
|
}
|
||||||
|
}
|
||||||
15
src/generalAlgorithm.hpp
Normal file
15
src/generalAlgorithm.hpp
Normal file
@ -0,0 +1,15 @@
|
|||||||
|
#ifndef __GENERAL_ALG__
|
||||||
|
#define __GENERAL_ALG__
|
||||||
|
|
||||||
|
#include <armadillo>
|
||||||
|
|
||||||
|
arma::vec* general_algorithm(
|
||||||
|
arma::vec* sub_diag,
|
||||||
|
arma::vec* main_diag,
|
||||||
|
arma::vec* sup_diag,
|
||||||
|
arma::vec* g_vec
|
||||||
|
);
|
||||||
|
|
||||||
|
void general_algorithm_main();
|
||||||
|
|
||||||
|
#endif
|
||||||
124
src/main.cpp
124
src/main.cpp
@ -1,24 +1,124 @@
|
|||||||
#include "GeneralAlgorithm.hpp"
|
|
||||||
#include <armadillo>
|
#include <armadillo>
|
||||||
#include <iostream>
|
#include <cmath>
|
||||||
|
#include <ctime>
|
||||||
|
#include <fstream>
|
||||||
|
#include <iomanip>
|
||||||
|
#include <ios>
|
||||||
|
#include <string>
|
||||||
|
|
||||||
|
#include "generalAlgorithm.hpp"
|
||||||
|
#include "specialAlgorithm.hpp"
|
||||||
|
|
||||||
|
#define TIMING_ITERATIONS 5
|
||||||
|
|
||||||
|
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) {
|
double f(double x) {
|
||||||
return 100. * std::exp(-10.*x);
|
return 100*std::exp(-10*x);
|
||||||
}
|
}
|
||||||
|
|
||||||
double a_sol(double x) {
|
void build_g_vec(int n_steps, arma::vec* g_vec) {
|
||||||
return 1. - (1. - std::exp(-10)) * x - std::exp(-10*x);
|
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);
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
int main() {
|
void build_array(
|
||||||
arma::mat A = arma::eye(3,3);
|
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);
|
||||||
|
|
||||||
GeneralAlgorithm ga(3, &A, f, a_sol, 0., 1.);
|
sub_diag->fill(-1);
|
||||||
|
main_diag->fill(2);
|
||||||
|
sup_diag->fill(-1);
|
||||||
|
|
||||||
ga.solve();
|
build_g_vec(n_steps, g_vec);
|
||||||
std::cout << "Time: " << ga.time(5) << std::endl;
|
}
|
||||||
ga.error();
|
|
||||||
|
|
||||||
return 0;
|
|
||||||
|
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_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;
|
||||||
|
|
||||||
|
}
|
||||||
|
ofile
|
||||||
|
<< n_steps << ","
|
||||||
|
<< g_res / (double) TIMING_ITERATIONS << ","
|
||||||
|
<< s_res / (double) TIMING_ITERATIONS << std::endl;
|
||||||
|
}
|
||||||
|
|
||||||
|
ofile.close();
|
||||||
|
}
|
||||||
|
|
||||||
|
void error_file() {
|
||||||
|
|
||||||
}
|
}
|
||||||
|
|
||||||
|
int main()
|
||||||
|
{
|
||||||
|
timing();
|
||||||
|
general_algorithm_main();
|
||||||
|
special_algorithm_main();
|
||||||
|
}
|
||||||
|
|||||||
@ -1,162 +0,0 @@
|
|||||||
#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();
|
|
||||||
}
|
|
||||||
52
src/specialAlgorithm.cpp
Normal file
52
src/specialAlgorithm.cpp
Normal file
@ -0,0 +1,52 @@
|
|||||||
|
#include "specialAlgorithm.hpp"
|
||||||
|
#include <iomanip>
|
||||||
|
|
||||||
|
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 < 6; i++) {
|
||||||
|
steps = std::pow(10, i+1);
|
||||||
|
step_size = 1./(double) steps;
|
||||||
|
g_vec.resize(steps - 1);
|
||||||
|
|
||||||
|
v_vec = special_algorithm(sub_sig, main_sig, sup_sig, &g_vec);
|
||||||
|
|
||||||
|
ofile.open("special_algorithm_" + std::to_string(steps) + ".txt");
|
||||||
|
|
||||||
|
for (int j=0; j < v_vec->n_elem; j++) {
|
||||||
|
ofile << std::setprecision(4) << std::scientific << step_size*(i+1) << ","
|
||||||
|
<< std::setprecision(4) << std::scientific << (*v_vec)(i) << std::endl;
|
||||||
|
}
|
||||||
|
ofile.close();
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
15
src/specialAlgorithm.hpp
Normal file
15
src/specialAlgorithm.hpp
Normal file
@ -0,0 +1,15 @@
|
|||||||
|
#ifndef __SPECIAL_ALG__
|
||||||
|
#define __SPECIAL_ALG__
|
||||||
|
|
||||||
|
#include <armadillo>
|
||||||
|
|
||||||
|
arma::vec* special_algorithm(
|
||||||
|
double sub_sig,
|
||||||
|
double main_sig,
|
||||||
|
double sup_sig,
|
||||||
|
arma::vec* g_vec
|
||||||
|
);
|
||||||
|
|
||||||
|
void special_algorithm_main();
|
||||||
|
|
||||||
|
#endif
|
||||||
Loading…
Reference in New Issue
Block a user