Update to current state

src/Peridynamics.jl

@@ -9,7 +9,7 @@ end
 export BBMaterial, CKIMaterial, NOSBMaterial, OSBMaterial, Body, point_set!,
        failure_permit!, material!, velocity_bc!, velocity_ic!, forcedensity_bc!, precrack!,
        VelocityVerlet, MultibodySetup, contact!, Job, read_vtk, uniform_box, submit,
-       process_each_export, @eff
+       process_each_export
 
 const MPI_INITIALIZED = Ref(false)
 const MPI_RANK = Ref(-1)
@@ -34,7 +34,7 @@ const ELASTIC_KWARGS = (:E, :nu)
 const FRAC_KWARGS = (:Gc, :epsilon_c)
 const DEFAULT_POINT_KWARGS = (:horizon, :rho, ELASTIC_KWARGS..., FRAC_KWARGS...)
 const CONTACT_KWARGS = (:radius, :sc)
-const EXPORT_KWARGS = (:path, :freq, :eff)
+const EXPORT_KWARGS = (:path, :freq, :fields)
 const DEFAULT_EXPORT_FIELDS = (:displacement, :damage)
 const JOB_KWARGS = (EXPORT_KWARGS...,)
 const SUBMIT_KWARGS = (:quiet,)
@@ -43,19 +43,19 @@ const PROCESS_EACH_EXPORT_KWARGS = (:serial,)
 const DimensionSpec = Union{Integer,Symbol}
 
 function __init__()
-    if MPI_INITIALIZED[]
-        return nothing
+    if !MPI_INITIALIZED[]
+        MPI.Init(finalize_atexit=true)
+        MPI_RANK[] = MPI.Comm_rank(MPI.COMM_WORLD)
+        MPI_SIZE[] = MPI.Comm_size(MPI.COMM_WORLD)
+        MPI_INITIALIZED[] = true
     end
-    MPI.Init(finalize_atexit=true)
-    MPI_RANK[] = MPI.Comm_rank(MPI.COMM_WORLD)
-    MPI_SIZE[] = MPI.Comm_size(MPI.COMM_WORLD)
-    MPI_INITIALIZED[] = true
     if haskey(ENV, "MPI_LOCALRANKID")
         MPI_SIM[] = true
     end
     @static if Sys.islinux() && !MPI_SIM[]
         pinthreads(:cores; force=false)
     end
+    BLAS.set_num_threads(1)
     return nothing
 end
 
@@ -127,7 +127,7 @@ include("auxiliary/process_each_export.jl")
 try
     include("auxiliary/precompile_workload.jl")
 catch err
-    @error "precompilation errored - can be safely ignored!\n" exception=err
+    @error "precompilation errored\n" exception=err
 end
 
 end

src/VtkReader/VtkReader.jl

@@ -244,9 +244,17 @@ function get_vtu_files(file::AbstractString)
     pieces = get_elements_by_tagname(unstruct_grid, "Piece")
     isnothing(pieces) && error("UnstructuredGrid does not contain `Piece`!\n")
     vtu_files::Vector{String} = attribute.(pieces, "Source"; required=true)
+    add_pwd!(vtu_files, dirname(file))
     return vtu_files
 end
 
+function add_pwd!(vtu_files::Vector{String}, dirname_file::String)
+    for i in eachindex(vtu_files)
+        vtu_files[i] = joinpath(dirname_file, vtu_files[i])
+    end
+    return nothing
+end
+
 function _read_vtu(vtu_file::AbstractString)
     raw_xml_str, raw_data = get_raw_xml_and_data(vtu_file)
     position_da, point_das, field_das = parse_data_arrays(raw_xml_str)

src/auxiliary/io.jl

@@ -1,37 +1,31 @@
-const ExportField = Tuple{Symbol,DataType}
-
-struct ExportOptions{F<:Function}
+struct ExportOptions
     exportflag::Bool
     root::String
     vtk::String
     logfile::String
     freq::Int
-    eff::F
+    fields::Vector{Symbol}
+end
 
-    function ExportOptions(root::String, freq::Int, eff::F) where {F<:Function}
-        if isempty(root)
-            _eff = ()->()
-            return new{typeof(_eff)}(false, "", "", "", 0, _eff)
-        end
-        vtk = joinpath(root, "vtk")
-        logfile = joinpath(root, "logfile.log")
-        return new{F}(true, root, vtk, logfile, freq, eff)
-    end
+function ExportOptions(root::String, freq::Int, fields::Vector{Symbol})
+    vtk = joinpath(root, "vtk")
+    logfile = joinpath(root, "logfile.log")
+    return ExportOptions(true, root, vtk, logfile, freq, fields)
 end
 
-function (eo::ExportOptions{F})(s::AbstractStorage) where {F<:Function}
-    return eo.eff(s)
+function ExportOptions()
+    return ExportOptions(false, "", "", "", 0, Vector{Symbol}())
 end
 
 function get_export_options(::Type{S}, o::Dict{Symbol,Any}) where {S<:AbstractStorage}
     local root::String
     local freq::Int
 
     if haskey(o, :path) && haskey(o, :freq)
-        root = string(o[:path])
+        root = abspath(string(o[:path]))
         freq = Int(o[:freq])
     elseif haskey(o, :path) && !haskey(o, :freq)
-        root = string(o[:path])
+        root = abspath(string(o[:path]))
         freq = 10
     elseif !haskey(o, :path) && haskey(o, :freq)
         msg = "if `freq` is spedified, the keyword `path` is also needed!\n"
@@ -42,42 +36,46 @@ function get_export_options(::Type{S}, o::Dict{Symbol,Any}) where {S<:AbstractSt
     end
     freq < 0 && throw(ArgumentError("`freq` should be larger than zero!\n"))
 
-    eff = get_exported_fields_function(S, o)
+    fields = get_export_fields(S, o)
 
-    return ExportOptions(root, freq, eff)
-end
+    if isempty(root)
+        eo = ExportOptions()
+    else
+        eo = ExportOptions(root, freq, fields)
+    end
 
-function export_disp_and_dmg(s::AbstractStorage)
-    return (("displacement", s.displacement), ("damage", s.damage))
+    return eo
 end
 
-macro eff(vars...)
-    local args = Expr[]
-    for var in vars
-        var isa QuoteNode || error("only symbols are allowed args!\n")
-        local name = string(var)[begin+1:end]
-        push!(args, :($(name), Core.getfield(x, $(esc(var)))))
+function get_export_fields(::Type{S}, o::Dict{Symbol,Any}) where {S}
+    local _fields::NTuple{N,Symbol} where {N}
+    if haskey(o, :fields)
+        _fields = o[:fields]
+    else
+        _fields = DEFAULT_EXPORT_FIELDS
     end
-    return Expr(:->, :x, Expr(:tuple, args...))
+
+    fields = [field for field in _fields]
+    check_storage_fields(S, fields)
+
+    return fields
 end
 
-function get_exported_fields_function(::Type{S}, o::Dict{Symbol,Any}) where {S}
-    if haskey(o, :eff)
-        eff = o[:eff]
-    else
-        eff = export_disp_and_dmg
+function check_storage_fields(::Type{S}, fields::Vector{Symbol}) where {S}
+    allowed_fields = fieldnames(S)
+    for f in fields
+        if !in(f, allowed_fields)
+            msg = "unknown field $(name) specified for export!\n"
+            msg *= "Allowed fields for $S:\n"
+            for allowed_name in fieldnames(S)
+                msg *= "  - $allowed_name\n"
+            end
+            throw(ArgumentError(msg))
+        end
     end
-    return eff
+    return nothing
 end
 
-# function unknown_fieldname_error(::Type{S}, name::Symbol) where {S<:AbstractStorage}
-#     msg = "unknown field $(name) specified for export!\n"
-#     msg *= "Allowed fields for $S:\n"
-#     for allowed_name in fieldnames(S)
-#         msg *= "  - $allowed_name\n"
-#     end
-#     throw(ArgumentError(msg))
-# end
 
 function export_results(dh::AbstractDataHandler, options::ExportOptions, chunk_id::Int,
                         timestep::Int, time::Float64)
@@ -98,22 +96,17 @@ end
 
 function _export_results(b::AbstractBodyChunk, chunk_id::Int, n_chunks::Int,
                          options::ExportOptions, n::Int, t::Float64)
-    filename = joinpath(options.vtk, @sprintf("timestep_%05d", n))
+    filename = @sprintf("timestep_%05d", n)
     position = get_loc_position(b)
     pvtk_grid(filename, position, b.cells; part=chunk_id, nparts=n_chunks) do vtk
-        for (name, field) in options(b.store)
-            vtk[name] = field
+        for field in options.fields
+            vtk[string(field), VTKPointData()] = get_storage_field(b.store, field)
         end
         vtk["time", VTKFieldData()] = t
     end
     return nothing
 end
 
-# @inline function get_export_field(s::AbstractStorage, name::Symbol, ::Type{V}) where {V}
-#     export_field::V = getfield(s, name)
-#     return export_field
-# end
-
 @inline function get_loc_position(b::AbstractBodyChunk)
     return @views b.store.position[:, 1:b.ch.n_loc_points]
 end

src/conditions/boundary_conditions.jl

@@ -2,27 +2,18 @@
 """
 
 """
-struct SingleDimBC{F<:Function,G<:Function} <: AbstractCondition
+struct SingleDimBC{F<:Function} <: AbstractCondition
     fun::F
-    cff::G
     field::Symbol
     point_set::Symbol
     dim::UInt8
 end
 
-@inline function (b::SingleDimBC{F,G})(t::Float64) where {F<:Function,G<:Function}
+@inline function (b::SingleDimBC{F})(t::Float64) where {F}
     value::Float64 = b.fun(t)
     return value
 end
 
-@inline function (b::SingleDimBC{F,G})(s::AbstractStorage) where {F<:Function,G<:Function}
-    return b.cff(s)
-end
-
-@inline get_cff_velocity_half(s::AbstractStorage) = s.velocity_half
-@inline get_cff_velocity(s::AbstractStorage) = s.velocity
-@inline get_cff_b_ext(s::AbstractStorage) = s.b_ext
-
 function override_eachother(a::SingleDimBC, b::SingleDimBC)
     same_field = a.field === b.field
     same_point_set = a.point_set === b.point_set
@@ -34,14 +25,17 @@ function apply_bcs!(b::AbstractBodyChunk, time::Float64)
     for bc in b.sdbcs
         apply_bc!(b.store, b.psets, bc, time)
     end
+    for bc in b.pdsdbcs
+        apply_bc!(b.store, b.psets, bc, b.dscr.position, time)
+    end
     return nothing
 end
 
 function apply_bc!(s::AbstractStorage, psets::Dict{Symbol,Vector{Int}},
-                   bc::SingleDimBC{F,G}, time::Float64) where {F,G}
+                   bc::SingleDimBC{F}, time::Float64) where {F}
     value = bc(time)
     isnan(value) && return nothing
-    apply_sdbc!(bc(s), value, bc.dim, psets[bc.point_set])
+    apply_sdbc!(get_storage_field(s, bc.field), value, bc.dim, psets[bc.point_set])
     return nothing
 end
 
@@ -52,3 +46,41 @@ end
     end
     return nothing
 end
+
+struct PosDepSingleDimBC{F<:Function} <: AbstractCondition
+    fun::F
+    field::Symbol
+    point_set::Symbol
+    dim::UInt8
+end
+
+@inline function (b::PosDepSingleDimBC{F})(p::AbstractVector, t::Float64) where {F}
+    value::Float64 = b.fun(p, t)
+    return value
+end
+
+function override_eachother(a::PosDepSingleDimBC, b::PosDepSingleDimBC)
+    same_field = a.field === b.field
+    same_point_set = a.point_set === b.point_set
+    same_dim = a.dim == b.dim
+    return same_field && same_point_set && same_dim
+end
+
+function apply_bc!(s::AbstractStorage, psets::Dict{Symbol,Vector{Int}},
+                   bc::PosDepSingleDimBC{F}, position::Matrix{Float64},
+                   time::Float64) where {F}
+    apply_pdsdbc!(get_storage_field(s, bc.field), position, bc, psets[bc.point_set], time)
+    return nothing
+end
+
+@inline function apply_pdsdbc!(field::Matrix{Float64}, pos::Matrix{Float64},
+                               bc::PosDepSingleDimBC{F}, point_ids::Vector{Int},
+                               t::Float64) where {F}
+    @simd for i in point_ids
+        value = bc(SVector{3}(pos[1,i], pos[2,i], pos[3,i]), t)
+        if !isnan(value)
+            field[bc.dim, i] = value
+        end
+    end
+    return nothing
+end

src/conditions/condition_checks.jl

@@ -8,8 +8,13 @@ end
 function check_condition_function(f::F) where {F<:Function}
     func_method = get_method_of_function(f)
     args = get_argument_names_of_function(func_method)
-    if length(args) != 1
-        error("too many arguments for condition! Only one (time) is allowed!\n")
+    if length(args) == 1
+        return :sdbc
+    elseif length(args) == 2 && args[1] === :p && args[2] === :t
+        return :pdsdbc
+    else
+        msg = "wrong arguments for position dependent condition function!\n"
+        throw(ArgumentError(msg))
     end
     return nothing
 end

src/conditions/initial_conditions.jl

@@ -1,18 +1,13 @@
 """
 
 """
-struct SingleDimIC{G<:Function} <: AbstractCondition
+struct SingleDimIC <: AbstractCondition
     value::Float64
-    cff::G
     field::Symbol
     point_set::Symbol
     dim::UInt8
 end
 
-@inline function (sdic::SingleDimIC{G})(s::AbstractStorage) where {G<:Function}
-    return sdic.cff(s)
-end
-
 function override_eachother(a::SingleDimIC, b::SingleDimIC)
     same_field = a.field === b.field
     same_point_set = a.point_set === b.point_set
@@ -22,7 +17,7 @@ end
 
 function apply_ic!(b::AbstractBodyChunk, ic::SingleDimIC)
     for point_id in b.psets[ic.point_set]
-        setindex!(ic(b.store), ic.value, ic.dim, point_id)
+        setindex!(get_storage_field(b.store, ic.field), ic.value, ic.dim, point_id)
     end
     return nothing
 end

src/core/submit.jl

@@ -35,6 +35,21 @@ function submit_threads(job::Job, n_chunks::Int)
     point_decomp = PointDecomposition(job.spatial_setup, n_chunks)
     tdh = ThreadsDataHandler(job.spatial_setup, job.time_solver, point_decomp)
     init_time_solver!(job.time_solver, tdh)
+    _pwd = init_export(job.options)
     solve!(tdh, job.time_solver, job.options)
+    finish_export(job.options, _pwd)
     return tdh
 end
+
+function init_export(options)
+    options.exportflag || return ""
+    _pwd = pwd()
+    mkpath(options.vtk)
+    cd(options.vtk)
+    return _pwd
+end
+
+function finish_export(options, _pwd)
+    options.exportflag && cd(_pwd)
+    return nothing
+end