From 171505a4e078c7271e838c35f98f3df7f4cf77de Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Sat, 16 Dec 2023 02:05:34 +0100 Subject: [PATCH 1/6] implement --- src/callbacks/solution_saving.jl | 9 +++++++-- src/visualization/write2vtk.jl | 26 ++++++++++++++++++++------ 2 files changed, 27 insertions(+), 8 deletions(-) diff --git a/src/callbacks/solution_saving.jl b/src/callbacks/solution_saving.jl index 54eea689a..cbd6b558f 100644 --- a/src/callbacks/solution_saving.jl +++ b/src/callbacks/solution_saving.jl @@ -28,6 +28,7 @@ To ignore a custom quantity for a specific system, return `nothing`. - `custom_quantities...`: Additional user-defined quantities. - `write_meta_data`: Write meta data. - `verbose=false`: Print to standard IO when a file is written. +- `max_coordinates=2^15` Particles with absolute coordinates larger than this value will cause an error. # Examples ```julia @@ -57,6 +58,7 @@ struct SolutionSavingCallback{I, CQ} verbose :: Bool output_directory :: String prefix :: String + max_coordinates :: Float64 custom_quantities :: CQ latest_saved_iter :: Vector{Int} end @@ -66,6 +68,7 @@ function SolutionSavingCallback(; interval::Integer=0, dt=0.0, save_final_solution=true, output_directory="out", append_timestamp=false, prefix="", verbose=false, write_meta_data=true, + max_coordinates=32768.0, custom_quantities...) if dt > 0 && interval > 0 throw(ArgumentError("Setting both interval and dt is not supported!")) @@ -82,7 +85,8 @@ function SolutionSavingCallback(; interval::Integer=0, dt=0.0, solution_callback = SolutionSavingCallback(interval, save_initial_solution, save_final_solution, write_meta_data, verbose, output_directory, - prefix, custom_quantities, [-1]) + prefix, max_coordinates, custom_quantities, + [-1]) if dt > 0 # Add a `tstop` every `dt`, and save the final solution. @@ -137,7 +141,7 @@ end # affect! function (solution_callback::SolutionSavingCallback)(integrator) (; interval, output_directory, custom_quantities, write_meta_data, - verbose, prefix, latest_saved_iter) = solution_callback + verbose, prefix, latest_saved_iter, max_coordinates) = solution_callback vu_ode = integrator.u semi = integrator.p @@ -162,6 +166,7 @@ function (solution_callback::SolutionSavingCallback)(integrator) output_directory=output_directory, prefix=prefix, write_meta_data=write_meta_data, + max_coordinates=max_coordinates, custom_quantities...) # Tell OrdinaryDiffEq that u has not been modified diff --git a/src/visualization/write2vtk.jl b/src/visualization/write2vtk.jl index c6eaace56..a9f5d404c 100644 --- a/src/visualization/write2vtk.jl +++ b/src/visualization/write2vtk.jl @@ -15,6 +15,7 @@ Convert Trixi simulation data to VTK format. - `prefix`: Prefix for output files. Defaults to an empty string. - `write_meta_data`: Write meta data. - `custom_quantities...`: Additional custom quantities to include in the VTK output. TODO. +- `max_coordinates=Inf` Particles with absolute coordinates larger than this value will cause an error. # Example @@ -24,7 +25,7 @@ trixi2vtk(sol.u[end], semi, 0.0, iter=1, output_directory="output", prefix="solu TODO: example for custom_quantities """ function trixi2vtk(vu_ode, semi, t; iter=nothing, output_directory="out", prefix="", - write_meta_data=true, custom_quantities...) + write_meta_data=true, max_coordinates=Inf, custom_quantities...) (; systems) = semi v_ode, u_ode = vu_ode.x @@ -45,13 +46,15 @@ function trixi2vtk(vu_ode, semi, t; iter=nothing, output_directory="out", prefix trixi2vtk(v, u, t, system, periodic_box; output_directory=output_directory, system_name=filenames[system_index], iter=iter, prefix=prefix, - write_meta_data=write_meta_data, custom_quantities...) + write_meta_data=write_meta_data, max_coordinates=max_coordinates, + custom_quantities...) end end # Convert data for a single TrixiParticle system to VTK format function trixi2vtk(v, u, t, system, periodic_box; output_directory="out", prefix="", iter=nothing, system_name=vtkname(system), write_meta_data=true, + max_coordinates=Inf, custom_quantities...) mkpath(output_directory) @@ -71,6 +74,15 @@ function trixi2vtk(v, u, t, system, periodic_box; output_directory="out", prefix points = periodic_coords(current_coordinates(u, system), periodic_box) cells = [MeshCell(VTKCellTypes.VTK_VERTEX, (i,)) for i in axes(points, 2)] + if abs(maximum(points)) > max_coordinates || abs(minimum(points)) > max_coordinates + throw(DomainError(max_coordinates, + "At least one particle's absolute coordinates exceed `max_coordinates`")) + end + + println("max ", abs(maximum(points))) + println("min ", abs(minimum(points))) + println("max_cord", max_coordinates) + vtk_grid(file, points, cells) do vtk write2vtk!(vtk, v, u, t, system, write_meta_data=write_meta_data) @@ -107,10 +119,12 @@ end Convert coordinate data to VTK format. # Arguments -- `coordinates`: Coordinates to be saved. -- `output_directory` (optional): Output directory path. Defaults to `"out"`. -- `prefix` (optional): Prefix for the output file. Defaults to an empty string. -- `filename` (optional): Name of the output file. Defaults to `"coordinates"`. +- `coordinates`: Coordinates to be saved. + +# Keywords +- `output_directory='out'`: Output directory path. +- `prefix=''`: Prefix for the output file. +- `filename='coordinates'`: Name of the output file. # Returns - `file::AbstractString`: Path to the generated VTK file. From d676a6ec09449c089adfbc604cc8aa2420e24c16 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Sat, 16 Dec 2023 02:08:23 +0100 Subject: [PATCH 2/6] cleanup --- src/visualization/write2vtk.jl | 4 ---- 1 file changed, 4 deletions(-) diff --git a/src/visualization/write2vtk.jl b/src/visualization/write2vtk.jl index a9f5d404c..c1bc1bd8a 100644 --- a/src/visualization/write2vtk.jl +++ b/src/visualization/write2vtk.jl @@ -79,10 +79,6 @@ function trixi2vtk(v, u, t, system, periodic_box; output_directory="out", prefix "At least one particle's absolute coordinates exceed `max_coordinates`")) end - println("max ", abs(maximum(points))) - println("min ", abs(minimum(points))) - println("max_cord", max_coordinates) - vtk_grid(file, points, cells) do vtk write2vtk!(vtk, v, u, t, system, write_meta_data=write_meta_data) From c02ee97ce11ee66a1faf624693fa5eaecd6b78f3 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Mon, 18 Dec 2023 14:20:10 +0100 Subject: [PATCH 3/6] fix --- examples/n_body/n_body_solar_system.jl | 2 +- src/callbacks/solution_saving.jl | 5 ++--- 2 files changed, 3 insertions(+), 4 deletions(-) diff --git a/examples/n_body/n_body_solar_system.jl b/examples/n_body/n_body_solar_system.jl index f063d0599..d4444d5d6 100644 --- a/examples/n_body/n_body_solar_system.jl +++ b/examples/n_body/n_body_solar_system.jl @@ -39,7 +39,7 @@ tspan = (0.0, 10year) ode = semidiscretize(semi, tspan) info_callback = InfoCallback(interval=100000) -saving_callback = SolutionSavingCallback(dt=10day) +saving_callback = SolutionSavingCallback(dt=10day, max_coordinates=Inf) callbacks = CallbackSet(info_callback, saving_callback) diff --git a/src/callbacks/solution_saving.jl b/src/callbacks/solution_saving.jl index cbd6b558f..a874dde30 100644 --- a/src/callbacks/solution_saving.jl +++ b/src/callbacks/solution_saving.jl @@ -2,7 +2,7 @@ SolutionSavingCallback(; interval::Integer=0, dt=0.0, save_initial_solution=true, save_final_solution=true, - output_directory="out", append_timestamp=false, + output_directory="out", append_timestamp=false, max_coordinates=2^15, custom_quantities...) Callback to save the current numerical solution in VTK format in regular intervals. @@ -68,8 +68,7 @@ function SolutionSavingCallback(; interval::Integer=0, dt=0.0, save_final_solution=true, output_directory="out", append_timestamp=false, prefix="", verbose=false, write_meta_data=true, - max_coordinates=32768.0, - custom_quantities...) + max_coordinates=32768.0, custom_quantities...) if dt > 0 && interval > 0 throw(ArgumentError("Setting both interval and dt is not supported!")) end From 483af4f4a13f02f8f42e48c8406d3647c80c63ff Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Thu, 11 Jan 2024 11:52:36 +0100 Subject: [PATCH 4/6] review --- src/callbacks/solution_saving.jl | 2 +- src/visualization/write2vtk.jl | 12 +++++++----- 2 files changed, 8 insertions(+), 6 deletions(-) diff --git a/src/callbacks/solution_saving.jl b/src/callbacks/solution_saving.jl index a874dde30..449ddf37c 100644 --- a/src/callbacks/solution_saving.jl +++ b/src/callbacks/solution_saving.jl @@ -68,7 +68,7 @@ function SolutionSavingCallback(; interval::Integer=0, dt=0.0, save_final_solution=true, output_directory="out", append_timestamp=false, prefix="", verbose=false, write_meta_data=true, - max_coordinates=32768.0, custom_quantities...) + max_coordinates=Float64(2^15), custom_quantities...) if dt > 0 && interval > 0 throw(ArgumentError("Setting both interval and dt is not supported!")) end diff --git a/src/visualization/write2vtk.jl b/src/visualization/write2vtk.jl index c1bc1bd8a..737208447 100644 --- a/src/visualization/write2vtk.jl +++ b/src/visualization/write2vtk.jl @@ -75,8 +75,10 @@ function trixi2vtk(v, u, t, system, periodic_box; output_directory="out", prefix cells = [MeshCell(VTKCellTypes.VTK_VERTEX, (i,)) for i in axes(points, 2)] if abs(maximum(points)) > max_coordinates || abs(minimum(points)) > max_coordinates - throw(DomainError(max_coordinates, - "At least one particle's absolute coordinates exceed `max_coordinates`")) + println("Warning: At least one particle's absolute coordinates exceed `max_coordinates`") + for i in eachindex(points) + points[i] = clamp(points[i], -max_coordinates, max_coordinates) + end end vtk_grid(file, points, cells) do vtk @@ -118,9 +120,9 @@ Convert coordinate data to VTK format. - `coordinates`: Coordinates to be saved. # Keywords -- `output_directory='out'`: Output directory path. -- `prefix=''`: Prefix for the output file. -- `filename='coordinates'`: Name of the output file. +- `output_directory="out"`: Output directory path. +- `prefix=""`: Prefix for the output file. +- `filename="coordinates"`: Name of the output file. # Returns - `file::AbstractString`: Path to the generated VTK file. From e2e59bdcba07234af6ce9cbe61ffb2c02b94e350 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Thu, 11 Jan 2024 15:32:58 +0100 Subject: [PATCH 5/6] update --- src/callbacks/solution_saving.jl | 2 +- src/visualization/write2vtk.jl | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/callbacks/solution_saving.jl b/src/callbacks/solution_saving.jl index 449ddf37c..683636b7f 100644 --- a/src/callbacks/solution_saving.jl +++ b/src/callbacks/solution_saving.jl @@ -28,7 +28,7 @@ To ignore a custom quantity for a specific system, return `nothing`. - `custom_quantities...`: Additional user-defined quantities. - `write_meta_data`: Write meta data. - `verbose=false`: Print to standard IO when a file is written. -- `max_coordinates=2^15` Particles with absolute coordinates larger than this value will cause an error. +- `max_coordinates=2^15` Particles with absolute coordinates larger than this will be clipped. # Examples ```julia diff --git a/src/visualization/write2vtk.jl b/src/visualization/write2vtk.jl index 737208447..750ade0a6 100644 --- a/src/visualization/write2vtk.jl +++ b/src/visualization/write2vtk.jl @@ -15,7 +15,7 @@ Convert Trixi simulation data to VTK format. - `prefix`: Prefix for output files. Defaults to an empty string. - `write_meta_data`: Write meta data. - `custom_quantities...`: Additional custom quantities to include in the VTK output. TODO. -- `max_coordinates=Inf` Particles with absolute coordinates larger than this value will cause an error. +- `max_coordinates=Inf` Particles with absolute coordinates larger than this will be clipped. # Example @@ -75,7 +75,7 @@ function trixi2vtk(v, u, t, system, periodic_box; output_directory="out", prefix cells = [MeshCell(VTKCellTypes.VTK_VERTEX, (i,)) for i in axes(points, 2)] if abs(maximum(points)) > max_coordinates || abs(minimum(points)) > max_coordinates - println("Warning: At least one particle's absolute coordinates exceed `max_coordinates`") + println("Warning: At least one particle's absolute coordinates exceed `max_coordinates` and have been clipped") for i in eachindex(points) points[i] = clamp(points[i], -max_coordinates, max_coordinates) end From de3e9bc79ef0bc3818443162a8e8a8351efa60f9 Mon Sep 17 00:00:00 2001 From: Sven Berger Date: Thu, 11 Jan 2024 15:52:28 +0100 Subject: [PATCH 6/6] fix --- src/callbacks/solution_saving.jl | 2 +- src/visualization/write2vtk.jl | 6 ++++-- 2 files changed, 5 insertions(+), 3 deletions(-) diff --git a/src/callbacks/solution_saving.jl b/src/callbacks/solution_saving.jl index 683636b7f..884deb482 100644 --- a/src/callbacks/solution_saving.jl +++ b/src/callbacks/solution_saving.jl @@ -28,7 +28,7 @@ To ignore a custom quantity for a specific system, return `nothing`. - `custom_quantities...`: Additional user-defined quantities. - `write_meta_data`: Write meta data. - `verbose=false`: Print to standard IO when a file is written. -- `max_coordinates=2^15` Particles with absolute coordinates larger than this will be clipped. +- `max_coordinates=2^15` The coordinates of particles will be clipped if their absolute values exceed this threshold. # Examples ```julia diff --git a/src/visualization/write2vtk.jl b/src/visualization/write2vtk.jl index 750ade0a6..8f568cebc 100644 --- a/src/visualization/write2vtk.jl +++ b/src/visualization/write2vtk.jl @@ -15,7 +15,7 @@ Convert Trixi simulation data to VTK format. - `prefix`: Prefix for output files. Defaults to an empty string. - `write_meta_data`: Write meta data. - `custom_quantities...`: Additional custom quantities to include in the VTK output. TODO. -- `max_coordinates=Inf` Particles with absolute coordinates larger than this will be clipped. +- `max_coordinates=Inf` The coordinates of particles will be clipped if their absolute values exceed this threshold. # Example @@ -75,7 +75,9 @@ function trixi2vtk(v, u, t, system, periodic_box; output_directory="out", prefix cells = [MeshCell(VTKCellTypes.VTK_VERTEX, (i,)) for i in axes(points, 2)] if abs(maximum(points)) > max_coordinates || abs(minimum(points)) > max_coordinates - println("Warning: At least one particle's absolute coordinates exceed `max_coordinates` and have been clipped") + println("Warning: At least one particle's absolute coordinates exceed `max_coordinates`" + * + " and have been clipped") for i in eachindex(points) points[i] = clamp(points[i], -max_coordinates, max_coordinates) end