diff --git a/Project.toml b/Project.toml index 7535448..4ef661c 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "MaterialPointVisualizer" uuid = "9ce2fbfb-c269-402f-8683-a675189e795c" authors = ["ZenanH "] -version = "0.1.11" +version = "0.2.0" [deps] ColorSchemes = "35d6a980-a343-548e-a6ea-1d62b119f2f4" @@ -12,18 +12,20 @@ HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f" Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7" +ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" WGLMakie = "276b4fcb-3e11-5398-bf8b-a0c2d153d008" WriteVTK = "64499a7a-5c06-52f2-abe2-ccb03c286192" [compat] -ColorSchemes = "3.30" +ColorSchemes = "3.31" DelimitedFiles = "1" -FastPointQuery = "0.1" +FastPointQuery = "0.2" HDF5 = "0.17" PrecompileTools = "1.2" +ProgressMeter = "1.10" StatsBase = "0.34" -WGLMakie = "0.11" +WGLMakie = "0.13" WriteVTK = "1.21" julia = "1.11" diff --git a/README.md b/README.md index 727281f..af9cf6d 100644 --- a/README.md +++ b/README.md @@ -3,7 +3,7 @@ [![CI](https://github.com/LandslideSIM/MaterialPointVisualizer.jl/actions/workflows/ci.yml/badge.svg)](https://github.com/LandslideSIM/MaterialPointVisualizer.jl/actions/workflows/ci.yml) [![Stable](https://img.shields.io/badge/docs-stable-blue.svg?logo=quicklook)](https://LandslideSIM.github.io/MaterialPointVisualizer.jl/stable) [![Dev](https://img.shields.io/badge/docs-dev-red.svg?logo=quicklook)](https://LandslideSIM.github.io/MaterialPointVisualizer.jl/dev) -[![Version](https://img.shields.io/badge/version-v0.1.11-pink)]() +[![Version](https://img.shields.io/badge/version-v0.2.0-pink)]() With this package, we can convert the MPM simulation results (HDF5 files from ***[MaterialPointSolver.jl](https://github.com/LandslideSIM/MaterialPointSolver.jl)*** ) into `.vtp` files or create ParaView-compatible animations. Additionally, it includes some post-processing functionalities. diff --git a/src/MaterialPointVisualizer.jl b/src/MaterialPointVisualizer.jl index d9e71b6..796a3fe 100644 --- a/src/MaterialPointVisualizer.jl +++ b/src/MaterialPointVisualizer.jl @@ -11,7 +11,8 @@ module MaterialPointVisualizer -using ColorSchemes, Dates, DelimitedFiles, HDF5, Logging, PrecompileTools, Printf, WGLMakie, WriteVTK +using ColorSchemes, Dates, DelimitedFiles, HDF5, Logging, PrecompileTools, Printf, + ProgressMeter, #=WGLMakie,=# WriteVTK import StatsBase: sample import FastPointQuery: trimesh, splashsurf, np, meshio @@ -26,10 +27,10 @@ import FastPointQuery: trimesh, splashsurf, np, meshio @sprintf("%02d days: %s", days, Dates.format(time, "HH:MM:SS")) end -include(joinpath(@__DIR__, "hdf2pvd.jl" )) -include(joinpath(@__DIR__, "pts2vtp.jl" )) -include(joinpath(@__DIR__, "pts2surf.jl" )) -include(joinpath(@__DIR__, "plot/display.jl")) +include(joinpath(@__DIR__, "hdf2pvd.jl")) +include(joinpath(@__DIR__, "pts2vtp.jl")) +include(joinpath(@__DIR__, "pts2surf.jl")) +# include(joinpath(@__DIR__, "plot/display.jl")) quiet(f) = redirect_stdout(devnull) do redirect_stderr(devnull) do @@ -39,6 +40,6 @@ quiet(f) = redirect_stdout(devnull) do end end -include(joinpath(@__DIR__, "precompile.jl")) +# include(joinpath(@__DIR__, "precompile.jl")) end \ No newline at end of file diff --git a/src/hdf2pvd.jl b/src/hdf2pvd.jl index af257eb..bd01f68 100644 --- a/src/hdf2pvd.jl +++ b/src/hdf2pvd.jl @@ -39,7 +39,7 @@ If your HDF5 file is generated by other means, please refer to [WriteVTK.jl](htt h5_str = joinpath(conf.prjdst, "$(conf.prjname).h5") h5_file = isfile(h5_str) ? h5_str : throw(ArgumentError("Expected a h5 file path at $h5_str")) anim_path = joinpath(conf.prjdst, "animation"); mkpath(anim_path) - pvd = paraview_collection(joinpath(conf.prjdst, "$(conf.prjname).pvd")) + #pvd = paraview_collection(joinpath(conf.prjdst, "$(conf.prjname).pvd")) fid = h5open(h5_file, "r") required = ("time", "ξ") groups = Tuple{String,Float64}[] @@ -51,34 +51,32 @@ If your HDF5 file is generated by other means, please refer to [WriteVTK.jl](htt end; sort!(groups, by = last) isempty(groups) && @warn "No valid time-step groups found in $h5_str" - p_total, p_iters = length(groups), 0 - t1, t_start = time(), time() - @inbounds for (step, (gname, t)) in enumerate(groups) - grp = fid[gname] - coords = read(grp["ξ"]) - vtk_path = joinpath(anim_path, @sprintf("%08d", step)) + p = Progress(length(1:1:length(groups)) - 1; dt=3.0, + desc = "\e[1;36m[ Info:\e[0m writing", + barlen = 20, + barglyphs = BarGlyphs(" ■■ ") + ) + @inbounds Threads.@threads for igroup in eachindex(groups) + grp = fid[groups[igroup][1]] + coords = permutedims(read(grp["ξ"])) + vtk_path = joinpath(anim_path, @sprintf("%08d", igroup)) vtp_cls = [MeshCell(PolyData.Verts(), [i]) for i in 1:size(coords, 2)] vtk_grid(vtk_path, coords, vtp_cls; compress=true, append=false, ascii=false) do vtk for name in keys(grp) name in ("ξ", "time") && continue - vtk[name] = read(grp[name]) + data = read(grp[name]) + data = ndims(data) == 2 ? permutedims(data) : data + vtk[name] = data end - pvd[t] = vtk # current simulation time - end - - t2 = time(); if t2 - t1 > 3.0 - percen = "animation (.pvd): " * @sprintf("%6.2f%%", (p_iters / p_total) * 100) - eta = "eta: " * format_seconds((p_total - p_iters) / (p_iters / (t2 - t_start))) - invo_str = " \e[1;32m⇌\e[0m " - info_con = "\e[1;32m[ Info: \e[0m" - print(stdout, "\r\e[2K"); print(stdout, info_con * percen * invo_str * eta) - t1 = t2 end - p_iters += 1; p_iters + 1 ≥ p_total && println() + next!(p) end - close(pvd) # close Paraview collection _grid2vtr(fid, joinpath(conf.prjdst, "grid.vtr")) close(fid) # close HDF5 文件 + times = map(last, groups) # 取出排序后的真实时间 + pvd_path = joinpath(conf.prjdst, "$(conf.prjname).pvd") + _write_pvd(pvd_path, times; basedir="animation", ext=".vtp", part=0) + return nothing end @views function _grid2vtr(fid, out_file::String) @@ -97,4 +95,25 @@ end end vtr["h"] = h close(vtr) +end + +function _write_pvd(outpath::AbstractString, times::AbstractVector{<:Real}; + basedir::AbstractString = "animation", + ext::AbstractString = ".vtp", + part::Integer = 0) + open(outpath, "w") do io + println(io, """""") + println(io, """""") + println(io, " ") + # 用 %.17g 写 timestep，既短又精确（避免 0.14749999999999996 这种串太长） + for (i, t) in enumerate(times) + step = @sprintf("%08d", i) + tstr = @sprintf("%.17g", float(t)) + fileattr = string(basedir, "/", step, ext) # 强制使用正斜杠，跨平台稳妥 + println(io, """ """) + end + println(io, " ") + println(io, "") + end + return outpath end \ No newline at end of file diff --git a/src/precompile.jl b/src/precompile.jl index 889bfe1..fe9abd9 100644 --- a/src/precompile.jl +++ b/src/precompile.jl @@ -1,8 +1,8 @@ @setup_workload begin @compile_workload begin quiet() do - scatter(rand(10, 3)) - pts2vtp(rand(2, 10); vtp_file=joinpath(tempdir(), "output.vtp"), data=(x=rand(10), y=rand(10))) + #scatter(rand(10, 3)) + # pts2vtp(rand(2, 10); vtp_file=joinpath(tempdir(), "output.vtp"), data=(x=rand(10), y=rand(10))) end end end \ No newline at end of file diff --git a/src/pts2surf.jl b/src/pts2surf.jl index 8f7dce2..60a7b83 100644 --- a/src/pts2surf.jl +++ b/src/pts2surf.jl @@ -62,8 +62,8 @@ function pts2surf(coords::AbstractArray, output_file::String; subdomain_num_cubes_per_dim::Int=64, output_mesh_smoothing_weights::Bool=true ) - @assert size(coords, 1) ∈ (2, 3) "coords should be a 2D or 3D array" - particles = np.array(coords', dtype=np.float64) + @assert size(coords, 2) ∈ (2, 3) "coords should be a 2D or 3D array" + particles = np.array(coords, dtype=np.float64) @assert isfile(output_file) "output_file should be a valid file path" file_name = endswith(lowercase(output_file), ".obj") ? output_file : output_file * ".obj" mesh_with_data, reconstruction = splashsurf.reconstruction_pipeline(particles, @@ -177,7 +177,7 @@ function hdf2surf(conf::NamedTuple; t1, t_start = time(), time() @inbounds for (step, (gname, t)) in enumerate(groups) grp = fid[gname] - coords = read(grp["ξ"])' # 读取粒子坐标 + coords = read(grp["ξ"]) # 读取粒子坐标 obj_filename = joinpath(obj_path, "$(@sprintf("%08d", step)).obj") pts2surf(coords, obj_filename; radius=radius, diff --git a/src/pts2vtp.jl b/src/pts2vtp.jl index 27a7b70..e0f4304 100644 --- a/src/pts2vtp.jl +++ b/src/pts2vtp.jl @@ -19,12 +19,14 @@ Description: Generates a `.vtp` file by passing custom fields. """ function pts2vtp(coords; vtp_file="output.vtp", data::NamedTuple=NamedTuple()) - pts_num = size(coords, 2) - size(coords, 1) in [2, 3] || throw(ArgumentError("The input coordinates should be 2D or 3D")) + pts_num = size(coords, 1) + size(coords, 2) in [2, 3] || throw(ArgumentError("The input coordinates should be 2D or 3D")) vtp_cls = [MeshCell(PolyData.Verts(), [i]) for i in 1:pts_num] - vtk_grid(vtp_file, coords, vtp_cls; compress=true, append=false, ascii=false) do vtk + vtk_grid(vtp_file, coords', vtp_cls; compress=true, append=false, ascii=false) do vtk keys(data) ≠ () && for vtp_key in keys(data) - vtk[string(vtp_key)] = getfield(data, vtp_key) + vdata = getfield(data, vtp_key) + vdata = ndims(vdata) == 2 ? permutedims(vdata) : vdata + vtk[string(vtp_key)] = vdata end end end \ No newline at end of file