Merge pull request #236 from kaipartmann/custom_export

Custom export fields

Author: kaipartmann

src/Peridynamics.jl

@@ -79,7 +79,7 @@ abstract type AbstractBondSystemMaterial{Correction} <: AbstractMaterial end
 abstract type AbstractBondBasedMaterial{CM} <: AbstractBondSystemMaterial{CM} end
 abstract type AbstractCorrespondenceMaterial{CM,ZEM} <: AbstractBondSystemMaterial{ZEM} end
 abstract type AbstractRKCMaterial{CM,C} <: AbstractBondSystemMaterial{C} end
-abstract type AbstractBondAssociatedSystemMaterial <: AbstractMaterial end
+abstract type AbstractBondAssociatedSystemMaterial <: AbstractBondSystemMaterial{Nothing} end
 abstract type AbstractConstitutiveModel end
 abstract type AbstractStressIntegration end
 abstract type AbstractZEMStabilization <: AbstractCorrection end

src/auxiliary/io.jl

@@ -102,9 +102,11 @@ end
 function check_export_fields(::Type{S}, fields::Vector{Symbol}) where {S}
     allowed_fields = point_data_fields(S)
     for f in fields
-        if !in(f, allowed_fields)
+        if !in(f, allowed_fields) && !custom_field(f)
             msg = "unknown point data field `:$(f)` specified for export!\n"
-            msg *= "All point data fields of $S:\n"
+            msg *= "If you intend to export a custom field, please define a method:"
+            msg *= "\n   `custom_field(::Val{:$(f)}) = true`\n"
+            msg *= "Otherwise, see here all available point data fields of $S:\n"
             for allowed_name in allowed_fields
                 msg *= "  - $allowed_name\n"
             end
@@ -114,6 +116,11 @@ function check_export_fields(::Type{S}, fields::Vector{Symbol}) where {S}
     return nothing
 end
 
+custom_field(field::Symbol) = custom_field(Val(field))
+# this function can be specialized to indicate custom fields
+# if this is not done, the export_field function will error out!
+custom_field(::Val{field}) where {field} = false
+
 function get_vtk_filebase(body::AbstractBody, root::AbstractString)
     body_name = replace(string(get_name(body)), " " => "_")
     filebase = isempty(body_name) ? "timestep" : body_name * "_timestep"
@@ -134,7 +141,7 @@ function _export_results(options::AbstractJobOptions, chunk::AbstractBodyChunk,
     filename = get_filename(options, chunk.body_name, n)
     position = get_loc_position(chunk)
     pvtk_grid(filename, position, chunk.cells; part=chunk_id, nparts=n_chunks) do vtk
-        export_fields!(vtk, chunk, options.fields)
+        export_fields!(vtk, chunk, options.fields, t)
         vtk["time", VTKFieldData()] = t
     end
     return nothing
@@ -152,20 +159,41 @@ end
     return @sprintf("%s_%d", vtk_filebase[body_name], n)
 end
 
-function export_fields!(vtk, chunk, fields::Vector{Symbol})
+function export_fields!(vtk, chunk, fields::Vector{Symbol}, t)
+    (; mat, system, storage, paramsetup) = chunk
     for field in fields
-        point_data = get_loc_point_data(chunk.storage, chunk.system, field)
+        point_data = export_field(Val(field), mat, system, storage, paramsetup, t)
         vtk[string(field), VTKPointData()] = point_data
     end
     return nothing
 end
 
-function export_fields!(vtk, chunk, fields_spec::Dict{Symbol,Vector{Symbol}})
+function export_fields!(vtk, chunk, fields_spec::Dict{Symbol,Vector{Symbol}}, t)
     fields = fields_spec[chunk.body_name]
-    export_fields!(vtk, chunk, fields)
+    export_fields!(vtk, chunk, fields, t)
     return nothing
 end
 
+# this function can be specialized for each field, even custom export fields can be written!
+function export_field(::Val{field}, mat, system, storage, paramsetup, t) where {field}
+    return get_loc_point_data(storage, system, field)
+end
+
 @inline function get_loc_position(chunk::AbstractBodyChunk)
     return @views chunk.storage.position[:, 1:get_n_loc_points(chunk)]
 end
+
+function msg_export_fields(fields::Vector{Symbol}; indentation=2)
+    return msg_vec("exported fields", fields; continuation_label=" "^15, indentation)
+end
+
+function msg_export_fields(fields_spec::Dict{Symbol,Vector{Symbol}}; indentation=2)
+    msg = " "^indentation * "EXPORTED FIELDS PER BODY\n"
+    for (body_name, fields) in fields_spec
+        body_name_str = string(body_name)
+        n = length(body_name_str)
+        ind = indentation + 2
+        msg *= msg_vec(body_name_str, fields; continuation_label=" "^n, indentation=ind)
+    end
+    return msg
+end

src/auxiliary/logs.jl

@@ -214,6 +214,224 @@ function msg_qty(descr::Any, qty::Any; kwargs...)
     return msg_qty(string(descr), string(qty); kwargs...)
 end
 
+function msg_path(descr_raw::AbstractString, path_raw::AbstractString;
+                  linewidth::Union{Int,Nothing}=nothing,
+                  leftwidth::Union{Int,Nothing}=nothing, indentation::Int=2,
+                  filler::Char='.', separator::AbstractString="",
+                  delimiter::AbstractString=" ",
+                  continuation_label::AbstractString="(continued)")
+    descr, path = strip(descr_raw), strip(path_raw)
+
+    # Determine effective linewidth
+    effective_linewidth = something(linewidth, leftwidth !== nothing ? default_linewidth() : nothing, default_linewidth())
+
+    # Try to fit on one line
+    len_filling = get_len_filling(linewidth, leftwidth, indentation, descr, separator, path)
+    required_length = indentation + length(descr) + length(separator) + max(len_filling, 1) + length(path)
+
+    if required_length <= effective_linewidth
+        # Fits on one line - use msg_qty
+        return msg_qty(descr, path; linewidth=linewidth, leftwidth=leftwidth,
+                      indentation=indentation, filler=filler, separator=separator,
+                      delimiter=delimiter, newline=true)
+    end
+
+    # Path is too long - split into multiple lines
+    return _msg_path_multiline(descr, path, effective_linewidth, indentation, filler,
+                               separator, delimiter, continuation_label)
+end
+
+function _msg_path_multiline(descr, path, linewidth, indentation, filler, separator, delimiter, continuation_label)
+    min_filler = 2 * length(delimiter)
+
+    # First line: calculate budget for path
+    first_line_budget = linewidth - indentation - length(descr) - length(separator) - min_filler
+    break_idx = _find_path_break(path, first_line_budget)
+
+    path_first = break_idx > 0 ? path[1:break_idx] : ""
+    path_rest = break_idx > 0 ? path[break_idx+1:end] : path
+
+    # Build first line
+    first_line = if !isempty(path_first)
+        len_fill = linewidth - indentation - length(descr) - length(separator) - length(path_first)
+        _msg_qty(descr, path_first, len_fill, indentation, filler, separator, delimiter, true)
+    else
+        len_fill = linewidth - indentation - length(descr) - length(separator)
+        _msg_qty(descr, "", len_fill, indentation, filler, separator, delimiter, true)
+    end
+
+    # Build continuation lines
+    continuation_lines = _build_continuation_lines(path_rest, continuation_label, linewidth,
+                                                   indentation, filler, delimiter)
+
+    return first_line * continuation_lines
+end
+
+function _build_continuation_lines(remaining_path, label, linewidth, indentation, filler, delimiter)
+    isempty(remaining_path) && return ""
+
+    result = ""
+    min_filler = 2 * length(delimiter)
+    prefix_len = indentation + length(label)
+
+    while !isempty(remaining_path)
+        budget = linewidth - prefix_len - min_filler
+
+        if length(remaining_path) <= budget
+            # Last piece fits
+            len_fill = linewidth - prefix_len - length(remaining_path)
+            result *= _msg_qty(label, remaining_path, len_fill, indentation, filler, "", delimiter, true)
+            break
+        else
+            # Need another split
+            break_idx = _find_path_break(remaining_path, budget)
+            piece = break_idx > 0 ? remaining_path[1:break_idx] : remaining_path[1:min(budget, end)]
+            remaining_path = length(piece) < length(remaining_path) ? remaining_path[length(piece)+1:end] : ""
+
+            len_fill = linewidth - prefix_len - length(piece)
+            result *= _msg_qty(label, piece, len_fill, indentation, filler, "", delimiter, true)
+        end
+    end
+
+    return result
+end
+
+function _find_path_break(path::AbstractString, max_length::Int)
+    length(path) <= max_length && return length(path)
+    max_length <= 0 && return 0
+
+    # Find last directory separator within budget
+    last_sep = findlast(c -> c == '/' || c == '\\', @view path[1:min(max_length, end)])
+
+    return something(last_sep, min(max_length, length(path)))
+end
+
+function msg_vec(descr_raw::AbstractString, vec::AbstractVector;
+                 linewidth::Union{Int,Nothing}=nothing,
+                 leftwidth::Union{Int,Nothing}=nothing, indentation::Int=2,
+                 filler::Char='.', separator::AbstractString="",
+                 delimiter::AbstractString=" ",
+                 continuation_label::AbstractString="(continued)",
+                 vec_delimiter::AbstractString=", ",
+                 vec_brackets::Tuple{AbstractString,AbstractString}=("[", "]"))
+    descr = strip(descr_raw)
+
+    # Format vector as string
+    vec_str = _format_vector(vec, vec_delimiter, vec_brackets)
+
+    # Determine effective linewidth
+    effective_linewidth = something(linewidth, leftwidth !== nothing ? default_linewidth() : nothing, default_linewidth())
+
+    # Try to fit on one line
+    len_filling = get_len_filling(linewidth, leftwidth, indentation, descr, separator, vec_str)
+    required_length = indentation + length(descr) + length(separator) + max(len_filling, 1) + length(vec_str)
+
+    if required_length <= effective_linewidth
+        # Fits on one line - use msg_qty
+        return msg_qty(descr, vec_str; linewidth=linewidth, leftwidth=leftwidth,
+                      indentation=indentation, filler=filler, separator=separator,
+                      delimiter=delimiter, newline=true)
+    end
+
+    # Vector string is too long - split into multiple lines
+    return _msg_vec_multiline(descr, vec_str, effective_linewidth, indentation, filler,
+                              separator, delimiter, continuation_label)
+end
+
+function _format_vector(vec::AbstractVector, vec_delimiter::AbstractString,
+                       vec_brackets::Tuple{AbstractString,AbstractString})
+    isempty(vec) && return vec_brackets[1] * vec_brackets[2]
+    elements = _format_elements(vec)
+    vec_str = vec_brackets[1] * join(elements, vec_delimiter) * vec_brackets[2]
+    return vec_str
+end
+
+function _format_elements(vec::AbstractVector{<:Real})
+    return [@sprintf("%.7g", x) for x in vec]
+end
+function _format_elements(vec::AbstractVector{T}) where {T}
+    return [string(x) for x in vec]
+end
+
+function _msg_vec_multiline(descr, vec_str, linewidth, indentation, filler, separator,
+                            delimiter, continuation_label)
+    min_filler = 2 * length(delimiter)
+
+    # First line: calculate budget for vector string
+    first_line_budget = linewidth - indentation - length(descr) - length(separator) - min_filler
+    break_idx = _find_vec_break(vec_str, first_line_budget)
+
+    vec_first = break_idx > 0 ? vec_str[1:break_idx] : ""
+    vec_rest = break_idx > 0 ? vec_str[break_idx+1:end] : vec_str
+
+    # Build first line
+    first_line = if !isempty(vec_first)
+        len_fill = linewidth - indentation - length(descr) - length(separator) - length(vec_first)
+        _msg_qty(descr, vec_first, len_fill, indentation, filler, separator, delimiter, true)
+    else
+        len_fill = linewidth - indentation - length(descr) - length(separator)
+        _msg_qty(descr, "", len_fill, indentation, filler, separator, delimiter, true)
+    end
+
+    # Build continuation lines
+    continuation_lines = _build_vec_continuation_lines(vec_rest, continuation_label, linewidth,
+                                                       indentation, filler, delimiter)
+
+    return first_line * continuation_lines
+end
+
+function _build_vec_continuation_lines(remaining_vec, label, linewidth, indentation, filler, delimiter)
+    isempty(remaining_vec) && return ""
+
+    result = ""
+    min_filler = 2 * length(delimiter)
+    prefix_len = indentation + length(label)
+
+    while !isempty(remaining_vec)
+        budget = linewidth - prefix_len - min_filler
+
+        if length(remaining_vec) <= budget
+            # Last piece fits
+            len_fill = linewidth - prefix_len - length(remaining_vec)
+            result *= _msg_qty(label, remaining_vec, len_fill, indentation, filler, "", delimiter, true)
+            break
+        else
+            # Need another split
+            break_idx = _find_vec_break(remaining_vec, budget)
+            piece = break_idx > 0 ? remaining_vec[1:break_idx] : remaining_vec[1:min(budget, end)]
+            remaining_vec = length(piece) < length(remaining_vec) ? remaining_vec[length(piece)+1:end] : ""
+
+            len_fill = linewidth - prefix_len - length(piece)
+            result *= _msg_qty(label, piece, len_fill, indentation, filler, "", delimiter, true)
+        end
+    end
+
+    return result
+end
+
+function _find_vec_break(vec_str::AbstractString, max_length::Int)
+    length(vec_str) <= max_length && return length(vec_str)
+    max_length <= 0 && return 0
+
+    # Find last comma or space within budget (prefer comma)
+    last_comma = findlast(==(','), @view vec_str[1:min(max_length, end)])
+
+    if last_comma !== nothing
+        # Include the comma and any following spaces
+        break_point = last_comma
+        while break_point < min(max_length, length(vec_str)) &&
+              vec_str[break_point + 1] == ' '
+            break_point += 1
+        end
+        return break_point
+    end
+
+    # If no comma, look for space
+    last_space = findlast(==(' '), @view vec_str[1:min(max_length, end)])
+
+    return something(last_space, min(max_length, length(vec_str)))
+end
+
 function msg_fields_inline(obj::T) where {T}
     fields = fieldnames(T)
     values = Tuple(getfield(obj, field) for field in fields)

src/auxiliary/nans.jl

@@ -5,8 +5,8 @@ function containsnan(K::T) where {T<:AbstractArray}
     return false
 end
 
-function nancheck(chunk::AbstractBodyChunk, t, Δt)
-    if containsnan(chunk.storage.b_int)
+function nancheck(storage::AbstractStorage, t, Δt)
+    if containsnan(storage.b_int)
         n = (t > 0 && Δt > 0) ? Int(t ÷ Δt) : 0
         msg = "NaN's found in field `b_int` at simulation time $(t), step $(n)!\n"
         error(msg)

src/core/data_handler.jl

@@ -1,4 +1,5 @@
 function initialize!(dh, solver)
+    calc_force_density!(dh, 0.0, solver.Δt)
     return nothing
 end
 

src/core/job.jl

@@ -180,3 +180,16 @@ function get_root_and_freq(o::Dict{Symbol,Any})
 
     return root, freq
 end
+
+function log_job_options(options::JobOptions)
+    msg = "EXPORT OPTIONS:\n"
+    if options.export_allowed
+        msg *= msg_qty("export frequency (export every nth step)", options.freq)
+        msg *= msg_path("root path", options.root; continuation_label="         ")
+        msg *= msg_export_fields(options.fields)
+    else
+        msg *= msg_qty("export allowed", "false")
+    end
+    log_it(options, msg)
+    return nothing
+end

src/core/mpi_body_data_handler.jl

@@ -421,7 +421,8 @@ function export_reference_results(dh::MPIBodyDataHandler, options::AbstractJobOp
     return nothing
 end
 
-function initialize!(::AbstractMPIBodyDataHandler, ::AbstractTimeSolver)
+function initialize!(dh::AbstractMPIBodyDataHandler, solver::AbstractTimeSolver)
+    calc_force_density!(dh, 0.0, solver.Δt)
     return nothing
 end
 

src/core/submit.jl

@@ -48,6 +48,7 @@ function submit_mpi(job::Job)
             initialize!(dh, job.time_solver)
             log_create_data_handler_end()
             log_data_handler(job.options, dh)
+            log_job_options(job.options)
             log_timesolver(job.options, job.time_solver)
         end
         @timeit_debug TO "solve!" begin
@@ -72,6 +73,7 @@ function submit_threads(job::Job, n_chunks::Int)
         initialize!(dh, job.time_solver)
         log_create_data_handler_end()
         log_data_handler(job.options, dh)
+        log_job_options(job.options)
         log_timesolver(job.options, job.time_solver)
         solve!(dh, job.time_solver, job.options)
     end

src/core/threads_multibody_data_handler.jl

@@ -85,8 +85,9 @@ function update_caches!(dh::ThreadsMultibodyDataHandler)
     return nothing
 end
 
-function initialize!(dh::AbstractThreadsMultibodyDataHandler, ::AbstractTimeSolver)
+function initialize!(dh::AbstractThreadsMultibodyDataHandler, solver::AbstractTimeSolver)
     init_contact_nhs!(dh)
+    calc_force_density!(dh, 0.0, solver.Δt)
     return nothing
 end