add copy(::Network) and improve metadata dealiasing

NEWS.md

@@ -4,7 +4,8 @@
 
 - **ModelingToolkit v11 compatibility** ([#344](https://github.com/JuliaDynamics/NetworkDynamics.jl/pull/344)):
   - Raises compat from `ModelingToolkit` v10 to v11. MTK v11 split into `ModelingToolkitBase` (MIT-licensed) and `ModelingToolkit` (AGPL-licensed); NetworkDynamics.jl now only requires the MIT-licensed `ModelingToolkitBase` as a direct dependency, dropping the AGPL dependency.
-  - A custom structural simplification pipeline is added as the new default (`mtkcompile=false`). With `mtkcompile=true`, MTK's own `mtkcompile` is called: full structural simplification when `ModelingToolkit` is loaded, or the stub from `ModelingToolkitBase` (no structural simplification) when it is not.
+  - A custom structural-simplification pipeline replaces MTK's `mtkcompile` as the new default (`mtkcompile=false`). It performs alias/linear-state elimination, breaks algebraic and nonlinear loops, and handles simple DAE index reduction. Set the global default via `NetworkDynamics.set_mtkcompile!`; use `mtkcompile=:compare` to run both backends and print a comparison for debugging. With `mtkcompile=true`, MTK's own `mtkcompile` is called: full structural simplification when `ModelingToolkit` is loaded, or the stub from `ModelingToolkitBase` (no structural simplification) when it is not.
+  - MTK models now warn when they contain discrete variables (unsupported).
   - New extension dependencies: `Moshi`, `Hungarian`.
   - **Dependency version bumps**: `SymbolicUtils` ≥ 4, `Symbolics` ≥ 7.
   - `implicit_output` macro updated to use `Symbolics.@register_symbolic` instead of `ModelingToolkit.@register_symbolic`.
@@ -14,8 +15,10 @@
   - `find_fixpoint(nw, x0::NWState, p::NWParameter)` is deprecated; use `find_fixpoint(nw, x0::NWState)`.
   - The default initial state for fixpoint search now applies guesses and formulas automatically (`guess=true, apply_formulas=true`).
 - **`NWState` constructor** gains `guess`, `apply_formulas`, and `verbose` keyword arguments. With `default=true` (the default), values are filled in order: defaults/inits → `InitFormula`s → guesses (if `guess=true`) → `GuessFormula`s. `apply_formulas=true` by default; `guess=false` by default.
+- `initialize_component` / `initialize_componentwise` gain a `warn` keyword (default `true`) to silence initialization warnings.
 - `doctor` check: added smoketest for the observable function (`obsf`) of each component.
-- `SciMLBase@3` and `OrdinaryDfifEq@7` compat. Most notably, upstream `DiffEq` changed its default `initializealg` to *check* rather than *reinit* inconsistent initial conditions. To preserve the previous behavior, the `ODEProblem(nw::Network, ...)` constructor now defaults to `initializealg=BrownFullBasicInit()` (deliberately differing from the OrdinaryDiffEq default). Pass `initializealg=...` to the constructor to override. Also, the interface of the `VectorContinuousComponentCallback` changed to fit that of `VectorContinuousCallback`: the affect now receives a buffer with per-element up or down crossing information rather than a single event idx.
+- `SciMLBase@3` and `OrdinaryDiffEq@7` compat. Most notably, upstream `DiffEq` changed its default `initializealg` to *check* rather than *reinit* inconsistent initial conditions. To preserve the previous behavior, the `ODEProblem(nw::Network, ...)` constructor now defaults to `initializealg=BrownFullBasicInit()` (deliberately differing from the OrdinaryDiffEq default). Pass `initializealg=...` to the constructor to override. Also, the interface of the `VectorContinuousComponentCallback` changed to fit that of `VectorContinuousCallback`: the affect now receives a buffer with per-element up or down crossing information rather than a single event idx.
+- new `copy(::Network)` method (dealiases all components and buffers)
 
 ## v0.10.17 Changelog
 - **Open-loop linearization** ([#341](https://github.com/JuliaDynamics/NetworkDynamics.jl/pull/341)):

ext/NetworkDynamicsMTKExt.jl

@@ -30,6 +30,12 @@ include("MTKExt_simplification.jl")
 import NetworkDynamics: implicit_output, RHSDifferentialsError
 Symbolics.@register_symbolic implicit_output(x)
 
+# `System`s and `Equation`s are large symbolic compile artifacts that are treated as
+# immutable once stored in component metadata (e.g. `:odesystem`, `:equations`). Opt
+# them into sharing so copying/reconstructing a component does not deep-copy them (which
+# is orders of magnitude slower). See `NetworkDynamics.dealias_metadata`.
+NetworkDynamics.dealias_metadata(x::Union{System,Equation}) = x
+
 const MTKCOMPILE_DEFAULT = Ref{Union{Symbol,Bool}}(false)
 function NetworkDynamics.set_mtkcompile!(value)
     MTKCOMPILE_DEFAULT[] = value

src/component_functions.jl

@@ -613,6 +613,19 @@ Base.@nospecializeinfer function _construct_comp(::Type{T}, @nospecialize(kwargs
     return c
 end
 
+"""
+    dealias_metadata(x)
+
+Recursively copy `x` so the result can be mutated without affecting the original,
+while avoiding costly deep copies of large immutable objects. Used to dealias a
+component's `metadata` and `symmetadata` when it is copied or reconstructed (see
+`copy(::ComponentModel)` and `_reconstruct_comp`).
+"""
+dealias_metadata(@nospecialize(x)) = deepcopy(x)
+dealias_metadata(x::ComponentModel) = copy(x)
+dealias_metadata(x::AbstractDict{Symbol}) = typeof(x)(k => dealias_metadata(v) for (k, v) in x)
+dealias_metadata(x::Union{AbstractVector,Tuple,NamedTuple}) = map(dealias_metadata, x)
+
 Base.@nospecializeinfer function _reconstruct_comp(::Type{T}, cf::ComponentModel, @nospecialize(kwargs)) where {T}
     fields = fieldnames(T)
     # fields to copy
@@ -623,9 +636,8 @@ Base.@nospecializeinfer function _reconstruct_comp(::Type{T}, cf::ComponentModel
     for f in nfields
         dict[f] = getproperty(cf, f)
     end
-    for f in cfields
-        dict[f] = deepcopy(getproperty(cf, f))
-    end
+    dict[:metadata]    = dealias_metadata(getproperty(cf, :metadata))
+    dict[:symmetadata] = dealias_metadata(getproperty(cf, :symmetadata))
 
     # if there is no change in symbols and it had a clash before, keep the setting
     if !haskey(kwargs, :allow_output_sym_clash) && !haskey(kwargs, :sym) && !haskey(kwargs, :outsym)
@@ -1114,8 +1126,10 @@ end
 """
     copy(c::NetworkDynamics.ComponentModel)
 
-Shallow copy of the component model. Creates a deepcopy of `metadata` and `symmetadata`
-but references the same objects everywhere else.
+Shallow copy of the component model. `metadata` and `symmetadata` are dealiased via
+`dealias_metadata` (fresh container spine, values deep-copied by default,
+immutable artifacts like MTK `System`s shared) but all other fields reference the same
+objects.
 """
 @generated function Base.copy(c::ComponentModel)
     fields = fieldnames(c)
@@ -1125,7 +1139,8 @@ but references the same objects everywhere else.
     nfields = setdiff(fields, cfields)
     assign = Expr(:block,
         (:($(field) = c.$field) for field in nfields)...,
-        (:($(field) = deepcopy(c.$field)) for field in cfields)...)
+        :(metadata = dealias_metadata(c.metadata)),
+        :(symmetadata = dealias_metadata(c.symmetadata)))
     construct = Expr(:call, c, [:($field) for field in fields]...)
 
     quote

src/construction.jl

@@ -523,3 +523,14 @@ function Network(nw::Network;
 
     return new_nw
 end
+
+"""
+    copy(nw::Network)
+
+Copy a `Network`. The vertex and edge component models are copied (see
+[`copy(::ComponentModel)`](@ref) for how their metadata is dealiased) while the
+underlying graph is shared. Much cheaper than `deepcopy(nw)`, which would additionally
+clone the (immutable) graph, caches and buffers as well as the components' symbolic
+compile artifacts (e.g. the MTK `System`s stored in metadata).
+"""
+Base.copy(nw::Network) = Network(nw)

test/construction_test.jl

@@ -348,6 +348,109 @@ end
     @test v1 != v2
 end
 
+# Two metadata-value stand-ins: `DealiasImmutable` mimics an arbitrary third-party
+# immutable-wrapper-of-mutable-state (no dealias override -> must be defensively
+# deep-copied). `DealiasShared` mimics an MTK `System`/`Equation`: an immutable compile
+# artifact that opts into sharing via a `dealias_metadata` method, exactly as the MTK
+# extension does. See NetworkDynamics.dealias_metadata.
+struct DealiasImmutable
+    payload::Vector{Int}
+end
+struct DealiasShared
+    payload::Vector{Int}
+end
+NetworkDynamics.dealias_metadata(x::DealiasShared) = x
+
+@testset "component metadata dealiasing" begin
+    # copy(::ComponentModel)/reconstruction dealiases metadata & symmetadata via
+    # `dealias_metadata`: rebuild the container spine, deepcopy leaves by default
+    # (defensive), copy nested ComponentModels, and share only types that opt in.
+    dealias = NetworkDynamics.dealias_metadata
+
+    # defensive default: an unknown immutable-wrapper is deep-copied, NOT shared
+    imm = DealiasImmutable([1,2,3])
+    @test dealias(imm) !== imm
+    @test dealias(imm).payload !== imm.payload         # inner mutable state independent
+    # opt-in sharing (the System/Equation pattern): a type with a method is shared as-is
+    shared = DealiasShared([1,2,3])
+    @test dealias(shared) === shared
+
+    # mutable leaf -> deepcopy; ComponentModel -> copy (independent but equal)
+    let v = [1,2,3]
+        @test dealias(v) !== v && dealias(v) == v
+    end
+    inner = VertexModel(; g=(o,i,p,t)->nothing, outdim=1, name=:inner, check=false)
+    @test dealias(inner) !== inner && dealias(inner) == inner
+
+    # Symbol-keyed containers: fresh spine of the same type, elements dealiased
+    d = Dict{Symbol,Any}(:a => shared, :b => [1,2,3])
+    d2 = dealias(d)
+    @test d2 isa Dict{Symbol,Any} && d2 !== d
+    @test d2[:a] === d[:a]                              # opted-in element shared
+    @test d2[:b] !== d[:b]                              # mutable element deep-copied
+    nt = (x = shared, y = [1,2,3])
+    nt2 = dealias(nt)
+    @test nt2 isa NamedTuple{(:x,:y)} && nt2.x === nt.x && nt2.y !== nt.y
+    @test dealias([shared, shared])[1] === shared       # vector: fresh, shared elements
+
+    # the dict rule is restricted to AbstractDict{Symbol}: a non-Symbol dict is deep-
+    # copied wholesale (defensive), so even an opted-in element inside is NOT shared
+    @test dealias(Dict{String,Any}("a" => shared))["a"] !== shared
+
+    # through copy(::ComponentModel)
+    f = (dv, v, ein, p, t) -> nothing
+    g = (out, in, p, t) -> nothing
+    md = Dict{Symbol,Any}(:artifact=>shared, :user_data=>[10,20,30], :nested=>inner)
+    v = VertexModel(; f, g, sym=[:x=>(;default=1)], outdim=2, metadata=md,
+                    name=:vp, check=false)
+    v2 = copy(v)
+    m, m2 = NetworkDynamics.metadata(v), NetworkDynamics.metadata(v2)
+
+    @test m2 !== m                          # fresh dict
+    @test m2[:artifact] === m[:artifact]    # opted-in artifact shared (no deepcopy!)
+    @test m2[:user_data] !== m[:user_data]  # mutable user data deep-copied, independent
+    m2[:user_data][1] = -1
+    @test m[:user_data][1] == 10            # original untouched
+    @test m2[:nested] !== m[:nested]        # nested ComponentModel: copied, not aliased
+    @test m2[:nested] == m[:nested]
+
+    # symmetadata: fresh inner dicts, mutation of the copy is isolated
+    sm, sm2 = NetworkDynamics.symmetadata(v), NetworkDynamics.symmetadata(v2)
+    @test sm2 !== sm
+    @test sm2[:x] !== sm[:x]
+    set_default!(v2, :x, 99)
+    @test get_default(v, :x) == 1
+    @test get_default(v2, :x) == 99
+end
+
+@testset "copy(::Network)" begin
+    using NetworkDynamics: getcomp
+    vertexf = (dv, u, esum, p, t) -> (dv[1] = esum[1]; nothing)
+    edgeg   = (out, src, dst, p, t) -> (out[1] = src[1] - dst[1]; nothing)
+    v1 = VertexModel(; f=vertexf, g=[1], sym=[:x=>(;default=1.0)], insym=[:i], vidx=1, name=:v1, check=false)
+    v2 = VertexModel(; f=vertexf, g=[1], sym=[:x=>(;default=2.0)], insym=[:i], vidx=2, name=:v2, check=false)
+    v3 = VertexModel(; f=vertexf, g=[1], sym=[:x=>(;default=3.0)], insym=[:i], vidx=3, name=:v3, check=false)
+    e1 = EdgeModel(; g=AntiSymmetric(edgeg), outsym=[:f], src=1, dst=2, check=false)
+    e2 = EdgeModel(; g=AntiSymmetric(edgeg), outsym=[:f], src=1, dst=3, check=false)
+    e3 = EdgeModel(; g=AntiSymmetric(edgeg), outsym=[:f], src=2, dst=3, check=false)
+    nw = Network([v1,v2,v3], [e1,e2,e3])
+
+    nw2 = copy(nw)
+    @test nw2 isa Network
+    @test nv(nw2) == nv(nw)
+    @test ne(nw2) == ne(nw)
+    # components are copied (independent objects) but structurally equal
+    @test getcomp(nw2, VIndex(1)) !== getcomp(nw, VIndex(1))
+    @test getcomp(nw2, VIndex(1)) == getcomp(nw, VIndex(1))
+    @test getcomp(nw2, EIndex(1)) !== getcomp(nw, EIndex(1))
+    # graph is shared (immutable), equivalent to but cheaper than deepcopy
+    @test nw2.im.g === nw.im.g
+    # mutating a component of the copy must not leak into the original
+    set_default!(nw2, VIndex(1, :x), 42.0)
+    @test get_default(getcomp(nw2, VIndex(1)), :x) == 42.0
+    @test get_default(getcomp(nw, VIndex(1)), :x) == 1
+end
+
 @testset "test network-remake constructor" begin
     g = (out, in, p, t) -> nothing
     ge = (out, src, dst, p, t) -> nothing