diff --git a/projects/project1/main.pdf b/projects/project1/main.pdf index 58de6bb..1267e38 100644 Binary files a/projects/project1/main.pdf and b/projects/project1/main.pdf differ diff --git a/projects/project1/main.tex b/projects/project1/main.tex index e56a89e..f22e9ad 100644 --- a/projects/project1/main.tex +++ b/projects/project1/main.tex @@ -243,7 +243,7 @@ Plotting the absolute error for different values of $N$ reveals the said converg \end{figure} \subsection*{Subproblem (b)} -Plotting the relative error for different values of $N$ reveals the increase of relative error as $u$ approaches 0 in the vicinity of the boundaries. This behavior is prevalent in all solutions as the limit of small divisors is approached. This behavior is quite logical, but it also shows, that not all uncertainties can be correctly minimized with an increase in resolution. The results are shown in \autoref{fig:rel_error_plot}. To overcome the division by zero issue for the relative error, only the values of $u > \num{1e-10}$ are considered. +Plotting the relative error for different values of $N$ reveals a constant value for the relative error accross the entire domain of $x$. The relative error decreases by two orders of magnitude for every order of magnitude increase in $N$. The results are summarized in \autoref{fig:rel_error_plot}. To avoiud issues with divisions by zero, the points where $u(x) \leq \num{e-10}$ were omitted from the plot. \begin{figure} \centering @@ -253,7 +253,7 @@ Plotting the relative error for different values of $N$ reveals the increase of \end{figure} \subsection*{Subproblem (c)} -The maximum relative errors are calculated directly using C++ as the output size scales linearly with $N$ and gets overwhelmingly large for bigger $N$. The results are output using the \texttt{std::cout} command and retrieved using GNU/Linux command line utilities. The code for calculating the maximum relative error $\max_i \frac{\delta v_i}{u_i}$ (see \eqref{eq:relerror}) is included in the C++ implementation. The results are summarized in \autoref{tab:rel_errors}. We see a steady decrease of the relative error up to $N = \num{100000}$, after which the relative error increases again. This will most likely be due to the artifacts seen in \autoref{fig:rel_error_plot} starting at $N = \num{10000}$. The source of these artifacts is not entirely clear. +The maximum relative errors are calculated directly using C++ as the output size scales linearly with $N$ and gets overwhelmingly large for bigger $N$. The results are output using the \texttt{std::cout} command and retrieved using GNU/Linux command line utilities. The code for calculating the maximum relative error $\max_i \frac{\delta v_i}{u_i}$ (see \eqref{eq:relerror}) is included in the C++ implementation. The results are summarized in \autoref{tab:rel_errors}. We see a steady decrease of the relative error up to $N = \num{100000}$, after which the relative error increases again. This will most likely be due to the limited numerical precision of floating point numbers in C++ and associated numerical errors. \begin{table}[H] \centering @@ -332,3 +332,4 @@ To evaluate the performance of the different algorithms, the algorithms were exe \bibliography{include/citations} \end{document} + \ No newline at end of file diff --git a/src/project1/poisson_solver.cpp b/src/project1/poisson_solver.cpp index 27c23e8..a752bc4 100644 --- a/src/project1/poisson_solver.cpp +++ b/src/project1/poisson_solver.cpp @@ -43,11 +43,11 @@ int main(int argc, char* argv[]){ // Output the solution string filename = "numerical_solution_" + to_string(N) + ".csv"; ofstream ofs(filename); - ofs << "x,u(x)" << endl; + ofs << "x,v(x),u(x)" << endl; for (int i = 0; i <= N; i++) { - ofs << fixed << setprecision(decimal_places) << scientific << x[i] << "," << solution[i] << endl; - } - ofs.close(); + ofs << fixed << setprecision(decimal_places) << scientific << x[i] << "," << solution[i] << "," << u[i] << endl; + } + ofs.close(); } diff --git a/src/project1/python/poisson_plotter.py b/src/project1/python/poisson_plotter.py index 287dd9d..2660edc 100644 --- a/src/project1/python/poisson_plotter.py +++ b/src/project1/python/poisson_plotter.py @@ -13,12 +13,43 @@ args = parser.parse_args() FILENAME = args.filename REFERENCE = args.reference +def get_rc_params(): + colors = ["FF220C", "70D6FF", "8AAA79", "666370", "1C1F33"] + rcParams = plt.rcParams + # Use LaTeX for rendering + # Setup fonts + rcParams["text.usetex"] = True + rcParams["font.family"] = "serif" + rcParams["font.size"] = 12 + rcParams["axes.labelsize"] = 12 + rcParams["axes.titlesize"] = 12 + rcParams["legend.fontsize"] = 10 + rcParams["xtick.labelsize"] = 10 + rcParams["ytick.labelsize"] = 10 + # Figure size and resolution + rcParams["figure.figsize"] = (6, 4) + rcParams["figure.dpi"] = 300 + # Use colors from the palette + rcParams["axes.prop_cycle"] = plt.cycler(color=[f"#{color}" for color in colors]) + # Grid + rcParams["axes.grid"] = True + rcParams["grid.alpha"] = 0.5 + rcParams["grid.linestyle"] = "--" + # Point ticks to the inside of the axes + rcParams["xtick.direction"] = "in" + rcParams["ytick.direction"] = "in" + rcParams["xtick.top"] = True + rcParams["ytick.right"] = True + return rcParams + + def singular_plot(): + plt.rcParams.update(get_rc_params()) df = pd.read_csv(FILENAME) fig, ax = plt.subplots(figsize=(6,3)) - ax.plot(df["x"], df["u(x)"]) + ax.plot(df["x"], df[df.columns[1]]) ax.set_xlabel(r"$x$") - ax.set_ylabel(r"$u(x)$") + ax.set_ylabel(rf"${df.columns[1]}$") fig.tight_layout() @@ -27,22 +58,21 @@ def singular_plot(): fig.savefig(os.path.join(output_directory, output_file)) def reference_plot(): + plt.rcParams.update(get_rc_params()) df = pd.read_csv(FILENAME).sort_values("x") - df_ref = pd.read_csv(REFERENCE).sort_values("x") fig, (ax1, ax2) = plt.subplots(nrows=2, sharex=True, height_ratios=[3, 1], figsize=(4, 8)) - ax1.plot(df["x"], df["u(x)"], label="Numerical Solution") - ax1.plot(df_ref["x"], df_ref["u(x)"], label="Reference Solution", linestyle="--") + ax1.plot(df["x"], df["v(x)"], label="Numerical Solution") + ax1.plot(df["x"], df["u(x)"], label="Reference Solution", linestyle="--") ax1.set_xlabel(r"$x$") ax1.set_ylabel(r"$u(x)$") ax1.legend() - combined_df = pd.merge_asof(df, df_ref, on="x", suffixes=("_num", "_ref")) - ratio = (combined_df["u(x)_num"] - combined_df["u(x)_ref"]) / combined_df["u(x)_ref"] + ratio = (df["v(x)"] - df["u(x)"]) / df["u(x)"] - ax2.plot(combined_df["x"], ratio.abs()) + ax2.plot(df["x"], ratio.abs()) ax2.set_xlabel(r"$x$") - ax2.set_ylabel(r"$|\frac{u(x)_{num} - u(x)_{ref}}{u(x)_{ref}}|$") + ax2.set_ylabel(r"$|\frac{v(x) - u(x)}{u(x)}|$") ax2.axhline(y=1, color="gray", linestyle="--") ax2.set_yscale("log") @@ -53,30 +83,29 @@ def reference_plot(): fig.savefig(os.path.join(output_directory, output_file)) def multiple_reference_plot(): + plt.rcParams.update(get_rc_params()) # Filename is wildcard -> Find files import glob files = glob.glob(FILENAME) files.sort() fig, (ax1, ax2) = plt.subplots(nrows=2, sharex=True, height_ratios=[3, 1], figsize=(6,5)) - df_ref = pd.read_csv(REFERENCE).sort_values("x") for file in files: df = pd.read_csv(file).sort_values("x") label = f"Num. (N = {os.path.basename(file).replace(".csv", "").split("_")[-1]})" - ax1.plot(df["x"], df["u(x)"], label=label) + ax1.plot(df["x"], df["v(x)"], label=label) - combined_df = pd.merge_asof(df, df_ref, on="x", suffixes=("_num", "_ref")) - ratio = (combined_df["u(x)_num"] - combined_df["u(x)_ref"]) / combined_df["u(x)_ref"] - ax2.plot(combined_df["x"], ratio.abs(), label=label) - - ax1.plot(df_ref["x"], df_ref["u(x)"], label="Reference Solution", linestyle="--", color="gray") + ratio = (df["v(x)"] - df["u(x)"]) / df["u(x)"] + ax2.plot(df["x"], ratio.abs(), label=label) + + ax1.plot(df["x"], df["u(x)"], label="Reference Solution", linestyle="--", color="gray") ax1.set_xlabel(r"$x$") ax1.set_ylabel(r"$u(x)$") ax1.legend() ax2.axhline(y=1, color="gray", linestyle="--") ax2.set_xlabel(r"$x$") - ax2.set_ylabel(r"$|\frac{u(x)_{num} - u(x)_{ref}}{u(x)_{ref}}|$") + ax2.set_ylabel(r"$|\frac{v(x) - u(x)}{u(x)}|$") ax2.set_yscale("log") fig.tight_layout() @@ -92,6 +121,7 @@ def error_plots(): def abs_error_plot(): + plt.rcParams.update(get_rc_params()) import glob files = glob.glob(FILENAME) files.sort() @@ -101,12 +131,12 @@ def abs_error_plot(): fig, ax = plt.subplots(figsize=(6,3)) for file in files: df = pd.read_csv(file).sort_values("x") - combined_df = pd.merge_asof(df, df_ref, on="x", suffixes=("_num", "_ref")) - log_abs_err = np.log10((combined_df["u(x)_num"] - combined_df["u(x)_ref"]).abs()) - ax.plot(combined_df["x"], log_abs_err, label=f"N={os.path.basename(file).replace('.csv','').split('_')[-1]}") + # combined_df = pd.merge_asof(df, df_ref, on="x", suffixes=("_num", "_ref")) + log_abs_err = np.log10((df["u(x)"] - df["v(x)"]).abs()) + ax.plot(df["x"], log_abs_err, label=f"N={os.path.basename(file).replace('.csv','').split('_')[-1]}") ax.set_xlabel(r"$x$") - ax.set_ylabel(r"$\log_{10}(|u(x)_{num} - u(x)_{ref}|)$") + ax.set_ylabel(r"$\log_{10}(|u(x) - v(x)|)$") ax.legend() fig.tight_layout() @@ -117,6 +147,7 @@ def abs_error_plot(): def rel_error_plot(): + plt.rcParams.update(get_rc_params()) import glob files = glob.glob(FILENAME) files.sort() @@ -126,14 +157,14 @@ def rel_error_plot(): fig, ax = plt.subplots(figsize=(6,3)) for file in files: df = pd.read_csv(file).sort_values("x") - combined_df = pd.merge_asof(df, df_ref, on="x", suffixes=("_num", "_ref")) + # combined_df = pd.merge_asof(df, df_ref, on="x", suffixes=("_num", "_ref")) # Drop rows where u(x)_ref is zero - combined_df = combined_df[combined_df["u(x)_ref"].abs() > 1e-10] - log_rel_err = np.log10((combined_df["u(x)_num"] - combined_df["u(x)_ref"]).abs() / combined_df["u(x)_ref"].abs()) - ax.plot(combined_df["x"], log_rel_err, label=f"N={os.path.basename(file).replace('.csv','').split('_')[-1]}") + df = df[df["u(x)"] > 1e-10] + log_rel_err = np.log10((df["u(x)"] - df["v(x)"]).abs() / df["u(x)"].abs()) + ax.plot(df["x"], log_rel_err, label=f"N={os.path.basename(file).replace('.csv','').split('_')[-1]}") ax.set_xlabel(r"$x$") - ax.set_ylabel(r"$\log_{10}(|\frac{u(x)_{num} - u(x)_{ref}}{u(x)_{ref}}|)$") + ax.set_ylabel(r"$\log_{10}(|\frac{u(x) - v(x)}{u(x)}|)$") ax.legend() fig.tight_layout() diff --git a/src/project1/python/timing_plotter.py b/src/project1/python/timing_plotter.py index f7dcf35..6082f46 100644 --- a/src/project1/python/timing_plotter.py +++ b/src/project1/python/timing_plotter.py @@ -4,10 +4,40 @@ import matplotlib.pyplot as plt import argparse import os +def get_rc_params(): + colors = ["FF220C", "70D6FF", "8AAA79", "666370", "1C1F33"] + rcParams = plt.rcParams + # Use LaTeX for rendering + # Setup fonts + rcParams["text.usetex"] = True + rcParams["font.family"] = "serif" + rcParams["font.size"] = 12 + rcParams["axes.labelsize"] = 12 + rcParams["axes.titlesize"] = 12 + rcParams["legend.fontsize"] = 10 + rcParams["xtick.labelsize"] = 10 + rcParams["ytick.labelsize"] = 10 + # Figure size and resolution + rcParams["figure.figsize"] = (6, 4) + rcParams["figure.dpi"] = 300 + # Use colors from the palette + rcParams["axes.prop_cycle"] = plt.cycler(color=[f"#{color}" for color in colors]) + # Grid + rcParams["axes.grid"] = True + rcParams["grid.alpha"] = 0.5 + rcParams["grid.linestyle"] = "--" + # Point ticks to the inside of the axes + rcParams["xtick.direction"] = "in" + rcParams["ytick.direction"] = "in" + rcParams["xtick.top"] = True + rcParams["ytick.right"] = True + return rcParams + parser = argparse.ArgumentParser(description="Plot timing results from CSV file.") parser.add_argument("csv_file", help="Path to the CSV file containing timing results.") args = parser.parse_args() +plt.rcParams.update(get_rc_params()) df = pd.read_csv(args.csv_file) for column in df.columns[1:]: df[column] = df[column] / df["N"] # Normalize by N