diff --git a/.gitignore b/.gitignore index 666f4a0..6aa0d5b 100644 --- a/.gitignore +++ b/.gitignore @@ -1,2 +1,2 @@ -**.DS_Store +**/.DS_Store **__pycache__ diff --git a/_site/plots/.DS_Store b/_site/plots/.DS_Store deleted file mode 100644 index d86e87f..0000000 Binary files a/_site/plots/.DS_Store and /dev/null differ diff --git a/main.py b/main.py index fb2502c..3e80cb5 100644 --- a/main.py +++ b/main.py @@ -50,7 +50,12 @@ def run(args): tf = max(0, args.Tf) cfl = min(1, args.cfl) - prefix = f"{args.pde}-{tf}-{args.ic}-{args.filter}-g{args.ggb}" + prefix = f"{args.pde}-{tf}-{args.ic}-{args.filter}-g" + if args.ggb: + prefix += str(args.ggblambda) + else: + prefix += "False" + print("PDE :", args.pde) print("TF :", tf) print("CFL :", cfl) @@ -75,7 +80,6 @@ def run(args): times = [0] us = [u_hat_init] - print(tf) if (tf > 0): output = solve_ivp(fun = rhs, t_span = [0, tf], @@ -95,8 +99,8 @@ def run(args): # Conveniently for us, the shock stays put at pi left = np.where(x < pi, True, False) right = np.where(x > pi, True, False) - ggbleft = ggb.recon_fft(x[left], uf, args.L) - ggbright = ggb.recon_fft(x[right], uf, args.L) + ggbleft = ggb.recon_fft(x[left], uf, args.ggblambda) + ggbright = ggb.recon_fft(x[right], uf, args.ggblambda) v[left] = ggbleft v[right] = ggbright sols.append((times[-1], us[0], uf, v)) @@ -118,26 +122,30 @@ def run(args): parser.add_argument('--Tf', type=float, default=1.0, help = "Final time") parser.add_argument('--pde', choices=('linadv', 'burgers'), default='burgers', help = "PDE to solve") # parser.add_argument('--add_visc', type=bool, default=False) - parser.add_argument('--filter', choices=('no_filter', 'exponential', 'cesaro', 'raisedcos', 'lanczos'), - default='no_filter', help = "Which filter, if any") + parser.add_argument('--filter', choices=('no_filter', 'exponential', 'cesaro', 'raisedcos', 'lanczos'), default='no_filter', help = "Which filter, if any") # parser.add_argument('--filterp', type=int, default=1, "p value for the exponential") - parser.add_argument('--ggb', type=bool, default=False, action='store_true' + parser.add_argument('--ggb', action='store_true', default=False, help = "Whether to reconstruct the analytic part of the solution") - parser.add_argument('-L', type=int, default=3, help = "Number of elements in the GGB basis") + parser.add_argument('---ggblambda', type=int, default=3, help = "Number of elements in the GGB basis") # parser.add_argument('--max_lgN', type=int, default=7, help = "Largest power of 2 to calculate until") #parser.add_argument('--integrator', choices=('solve_ivp', 'elrk4'), default='elrk4') - parser.add_argument('--show_markers', type=bool, default=False, action='store_true' + parser.add_argument('--show_markers',default=False, action='store_true', help = "Show the original sample values in red crosses") parser.add_argument('--exact', type=str, default=None, help = "Exact solution file") args = parser.parse_args() if args.N == None: - args.N = [16, 64] + args.N = [16, 64, 128] sols, prefix, initial_condition = run(args) plotname = prefix+'.svg' # PLOT! nn = 2048 fig, ax = plt.subplots(nrows=1, ncols=2, figsize = (14, 8), width_ratios=[3,1]) + # Exact + if args.exact != None: + exact = np.loadtxt(args.exact).transpose() + exact[0] *= 2 * pi + ax[0].plot(exact[0], exact[1], 'ko', markersize=0.8, label="Exact") # Initial x = np.linspace(xmin, xmax, 1000) ax[0].plot(x, initial_condition(x), linewidth=1, color='k', label="Init") @@ -145,8 +153,7 @@ def run(args): for (t, uhinit, uf, v) in sols: n = len(uf) label = str(n) + f", t={np.round(t, 3)}" - plots.smoothplot(uf, ax[0], label=label, linewidth=1) - plots.plot_resolution(uf, ax[1], linewidth=0.5, markersize=0.5) + # plots.smoothplot(uf, ax[0], label=label) if args.ggb: label += ",ggb" if args.filter != 'no_filter': @@ -154,17 +161,15 @@ def run(args): ax[0].plot(cgrid(n), v, label=label) elif args.filter != 'no_filter': label += "," + args.filter - plots.smoothplot(uf, ax[0]) + plots.smoothplot(uf, ax[0], label=label) + else: + # So no filter, no ggb + plots.smoothplot(uf, ax[0], label=label) if args.show_markers: ax[0].plot(cgrid(n), ifft(uf), "+", color= "red", markersize=5) - # Exact - if args.exact != None: - exact = np.loadtxt(args.exact).transpose() - ax[0].plot(exact[0] * 2 * np.pi, exact[1], 'ko', markersize=0.8, label="Exact") ax[0].grid(visible=True) ax[0].legend() fig.tight_layout() fig.savefig(plotname) - plt.close() print("Saved plot to "+plotname) diff --git a/plots.py b/plots.py index 1a0939d..48bc0c4 100644 --- a/plots.py +++ b/plots.py @@ -34,6 +34,19 @@ def plot_resolution(c, ax, **kwargs): ax.semilogy(k, np.abs(fftshift(c)), **kwargs) return +def plot_error_fft(uh, exact, ax, **kwargs): + x = exact[0] + if len(uh) < len(x): + u = ifft_at(x, uh) + e = np.abs(u - exact[1]) + ax[1].plot(x, e, **kwargs) + else: + print("WARNING: Interpolating exact not supported yet.\nSupply exact solution on finer grid.") + return + +def plot_error(c, exact, ax, **kwargs): + print("print_error not implemented yet.") + # Fourier interpolate the IFFT def smoothplot(v, ax,nn=2048, **plotargs): n = len(v) diff --git a/presentation.mplstyle b/presentation.mplstyle index 8ea6b28..1ed934e 100644 --- a/presentation.mplstyle +++ b/presentation.mplstyle @@ -1,5 +1,6 @@ axes.titlesize : 24 axes.labelsize : 20 +axes.linewidth : 1.5 lines.linewidth : 3 lines.markersize : 10 xtick.labelsize : 16