From e45af966038bb2ae8cea4ef4e496028dae1d074d Mon Sep 17 00:00:00 2001 From: Matt Date: Thu, 5 Dec 2024 16:44:13 -0500 Subject: [PATCH] applying comments except for removing :diffusivity --- packages/algjulia-service/demos.txt | 7 - .../src/decapodes-service/DecapodesService.jl | 42 +- .../src/decapodes-service/geometry.jl | 6 +- .../decapodes-service/initial_conditions.jl | 6 +- .../src/decapodes-service/model.jl | 28 +- .../src/decapodes-service/plotting.jl | 28 -- .../src/decapodes-service/theory.jl | 19 +- packages/algjulia-service/src/decapodes.jl | 416 ------------------ packages/algjulia-service/test/runtests.jl | 102 ++--- 9 files changed, 90 insertions(+), 564 deletions(-) delete mode 100644 packages/algjulia-service/demos.txt delete mode 100644 packages/algjulia-service/src/decapodes-service/plotting.jl delete mode 100644 packages/algjulia-service/src/decapodes.jl diff --git a/packages/algjulia-service/demos.txt b/packages/algjulia-service/demos.txt deleted file mode 100644 index 952ceda6..00000000 --- a/packages/algjulia-service/demos.txt +++ /dev/null @@ -1,7 +0,0 @@ -Navier-Stokes -http://localhost:5173/analysis/01938869-4114-73c0-b412-42d86d941296 -https://tinyurl.com/catcolabns - -Diffusion -http://localhost:5173/analysis/0193937b-0735-72e1-96ee-0f4ae6d46e9f -https://tinyurl.com/catcolabdiffusion \ No newline at end of file diff --git a/packages/algjulia-service/src/decapodes-service/DecapodesService.jl b/packages/algjulia-service/src/decapodes-service/DecapodesService.jl index 0655d40e..311785f6 100644 --- a/packages/algjulia-service/src/decapodes-service/DecapodesService.jl +++ b/packages/algjulia-service/src/decapodes-service/DecapodesService.jl @@ -39,7 +39,7 @@ include("model.jl") ## model-building (is this actually a diagram of a model?) struct PodeSystem pode::SummationDecapode plotvar::Vector{Symbol} - scalars::Dict{Symbol, Any} # closures # TODO rename scalars => anons + scalars::Dict{Symbol, Any} # closures geometry::Geometry init::ComponentArray # TODO Is it always ComponentVector? generate::Any @@ -47,39 +47,42 @@ struct PodeSystem end export PodeSystem +function PodeSystem(json_string::String, args...) + json_object = JSON3.read(json_string) + PodeSystem(json_object, args...) +end + """ Construct a `PodeSystem` object from a JSON string. """ -function PodeSystem(json_string::String,hodge=GeometricHodge()) - json_object = JSON3.read(json_string); - +function PodeSystem(json_object::AbstractDict, hodge=GeometricHodge()) # converts the JSON of (the fragment of) the theory # into theory of the DEC, valued in Julia - theory = Theory(json_object[:model]); + theory = Theory(json_object[:model]) # this is a diagram in the model of the DEC. it wants to be a decapode! - diagram = json_object[:diagram]; + diagram = json_object[:diagram] # any scalars? - scalars = haskey(json_object, :scalars) ? json_object[:scalars] : []; + scalars = haskey(json_object, :scalars) ? json_object[:scalars] : [] # pode, anons, and vars (UUID => ACSetId) - decapode, anons, vars = Decapode(diagram, theory; scalars=scalars); + decapode, anons, vars = Decapode(diagram, theory; scalars=scalars) dot_rename!(decapode) - uuid2symb = uuid_to_symb(decapode, vars); + uuid2symb = uuid_to_symb(decapode, vars) # plotting variables - plotvars = [uuid2symb[uuid] for uuid in json_object[:plotVariables]]; + plotvars = [uuid2symb[uuid] for uuid in json_object[:plotVariables]] # extract the domain in order to create the mesh, dualmesh - geometry = Geometry(json_object); + geometry = Geometry(json_object) # initialize operators - ♭♯_m = ♭♯_mat(geometry.dualmesh); - wedge_dp10 = dec_wedge_product_dp(Tuple{1,0}, geometry.dualmesh); - dual_d1_m = dec_mat_dual_differential(1, geometry.dualmesh); + ♭♯_m = ♭♯_mat(geometry.dualmesh) + wedge_dp10 = dec_wedge_product_dp(Tuple{1,0}, geometry.dualmesh) + dual_d1_m = dec_mat_dual_differential(1, geometry.dualmesh) star0_inv_m = dec_mat_inverse_hodge(0, geometry.dualmesh, hodge) - Δ0 = Δ(0,geometry.dualmesh); + Δ0 = Δ(0,geometry.dualmesh) #fΔ0 = factorize(Δ0); function sys_generate(s, my_symbol, hodge=hodge) op = @match my_symbol begin @@ -97,7 +100,7 @@ function PodeSystem(json_string::String,hodge=GeometricHodge()) # end initialize # initial conditions - u0 = initial_conditions(json_object, geometry, uuid2symb); + u0 = initial_conditions(json_object, geometry, uuid2symb) # symbol => uuid. we need this to reassociate the var symb2uuid = Dict([v => k for (k,v) in pairs(uuid2symb)]) @@ -106,7 +109,7 @@ function PodeSystem(json_string::String,hodge=GeometricHodge()) end export PodeSystem -points(system::PodeSystem) = collect(values(system.geometry.mesh.subparts.point.m)); +points(system::PodeSystem) = collect(values(system.geometry.mesh.subparts.point.m)) indexing_bounds(system::PodeSystem) = indexing_bounds(system.geometry.domain) ## SIM HELPERS ## @@ -156,7 +159,7 @@ function state_at_time(soln::ODESolution, system::PodeSystem, plotvar::Symbol, t end function state_at_time(soln::ODESolution, domain::Rectangle, var::Symbol, t::Int) # last two args can be one Point2 - (x, y) = indexing_bounds(domain); + (x, y) = indexing_bounds(domain) [SVector(i, j, getproperty(soln.u[t], var)[x*(i-1) + j]) for i in 1:x+1, j in 1:y+1] end @@ -167,9 +170,8 @@ function grid(pt3::Point3, grid_size::Vector{Int}) end function state_at_time(soln::ODESolution, domain::Sphere, var::Symbol, t::Int, points) - l , _ = indexing_bounds(domain); # TODO this is hardcoded to return 100, 100 + l , _ = indexing_bounds(domain) # TODO this is hardcoded to return 100, 100 northern_indices = filter(i -> points[i][3] > 0, keys(points)) - # TODO we don't have access to points in this function yet. need to pass that in map(northern_indices) do n i, j = grid(points[n], [l, l]) # TODO SVector(i, j, getproperty(soln.u[t], var)[n]) diff --git a/packages/algjulia-service/src/decapodes-service/geometry.jl b/packages/algjulia-service/src/decapodes-service/geometry.jl index 98ec9921..0fef00b8 100644 --- a/packages/algjulia-service/src/decapodes-service/geometry.jl +++ b/packages/algjulia-service/src/decapodes-service/geometry.jl @@ -74,7 +74,7 @@ struct Geometry dualmesh::HasDeltaSet end -function Geometry(json_object::JSON3.Object) +function Geometry(json_object::AbstractDict) mesh_name = Symbol(json_object[:mesh]) domain = PREDEFINED_MESHES[mesh_name] Geometry(domain) @@ -95,14 +95,14 @@ end # function Geometry(r::Periodic, division::SimplexCenter=Circumcenter()) end function Geometry(m::Sphere, division::SimplexCenter=Circumcenter()) - s = loadmesh(Icosphere(m.dim, m.radius)); + s = loadmesh(Icosphere(m.dim, m.radius)) sd = EmbeddedDeltaDualComplex2D{Bool, Float64, Point3{Float64}}(s) subdivide_duals!(sd, division) Geometry(m, s, sd) end function Geometry(m::UV, division::SimplexCenter=Circumcenter()) - s, _, _ = makeSphere(m); + s, _, _ = makeSphere(m) sd = EmbeddedDeltaDualComplex2D{Bool, Float64, Point3{Float64}}(s) subdivide_duals!(sd, division) Geometry(m, s, sd) diff --git a/packages/algjulia-service/src/decapodes-service/initial_conditions.jl b/packages/algjulia-service/src/decapodes-service/initial_conditions.jl index ff7d37a4..dcd3251e 100644 --- a/packages/algjulia-service/src/decapodes-service/initial_conditions.jl +++ b/packages/algjulia-service/src/decapodes-service/initial_conditions.jl @@ -43,8 +43,8 @@ GaussianIC(r::Rectangle) = GaussianIC(r, GaussianData(r)) TaylorVortexIC(d::Sphere) = TaylorVortexIC(d, TaylorVortexData()) -function initial_conditions(json_object::JSON3.Object, geometry::Geometry, uuid2symb::Dict{String, Symbol}) - ic_specs = json_object[:initialConditions]; # this is "C" +function initial_conditions(json_object::AbstractDict, geometry::Geometry, uuid2symb::Dict{String, Symbol}) + ic_specs = json_object[:initialConditions] # this is "C" dict = Dict([uuid2symb[string(uuid)] => ic_specs[string(uuid)] for uuid ∈ keys(ic_specs)]...) initial_conditions(dict, geometry) # the resulting sim will only have (C,) as initial conditions end @@ -99,7 +99,7 @@ end function initial_conditions(ics::TaylorVortexIC, geometry::Geometry) # TODO prefer not to load `s0` here but che sara sara - s0 = dec_hodge_star(0, geometry.dualmesh, GeometricHodge()); + s0 = dec_hodge_star(0, geometry.dualmesh, GeometricHodge()) X = vort_ring(ics, geometry) du = s0 * X return du diff --git a/packages/algjulia-service/src/decapodes-service/model.jl b/packages/algjulia-service/src/decapodes-service/model.jl index eceee606..17110939 100644 --- a/packages/algjulia-service/src/decapodes-service/model.jl +++ b/packages/algjulia-service/src/decapodes-service/model.jl @@ -6,7 +6,7 @@ export add_to_pode! function add_to_pode!(d::SummationDecapode, vars::Dict{String, Int}, # mapping between UUID and ACSet ID theory::Theory, - content::JSON3.Object, + content::AbstractDict, nc::Vector{Int}, ::ObType) theory_elem = theory.data[content[:over][:content]] # indexes the theory by UUID @@ -24,7 +24,7 @@ function add_to_pode!(d::SummationDecapode, return d end -function Base.nameof(theory::Theory, content::JSON3.Object) +function Base.nameof(theory::Theory, content::AbstractDict) Symbol(theory.data[content[:over][:content]].name) end @@ -32,12 +32,13 @@ end function add_to_pode!(d::SummationDecapode, vars::Dict{String, Int}, # mapping between UUID and ACSet ID theory::Theory, - content::JSON3.Object, + content::AbstractDict, scalars::Any, anons::Dict{Symbol, Any}, ::HomType) - dom = content[:dom][:content]; cod = content[:cod][:content] + dom = content[:dom][:content] + cod = content[:cod][:content] # TODO we need a safe way to fail this if haskey(vars, dom) && haskey(vars, cod) # get the name of the Op1 and add it to the theory @@ -53,21 +54,23 @@ function add_to_pode!(d::SummationDecapode, scalar = scalars[Symbol(content[:over][:content])] push!(anons, op1 => x -> scalar * x) end - # TODO if scalars were typed correctly, we could probably do away with the !isempty check end d end -""" Decapode(jsondiagram::JSON3.Object, theory::Theory) => SummationDecapode +""" Decapode(diagram::AbstractVector{<:AbstractDict}, theory::Theory) => (::SummationDecapode, ::Dict{Symbol, Any}, ::Dict{String, Int}) -This returns a Decapode given a jsondiagram and a theory. +This returns + 1. a Decapode + 2. a dictionary of symbols mapped to anonymous functions + 3. a dictionary of JSON UUIDs mapped to symbols """ -function Decapode(diagram::AbstractVector{JSON3.Object}, theory::Theory; scalars=[]) +function Decapode(diagram::AbstractVector{<:AbstractDict}, theory::Theory; scalars=[]) # initiatize decapode and its mapping between UUIDs and ACSet IDs - pode = SummationDecapode(parse_decapode(quote end)); - vars = Dict{String, Int}(); # UUID => ACSetID - nc = [0]; # array is a mutable container - anons = Dict{Symbol, Any}(); + pode = SummationDecapode(parse_decapode(quote end)) + vars = Dict{String, Int}() # UUID => ACSetID + nc = [0] # array is a mutable container + anons = Dict{Symbol, Any}() # for each cell in the notebook, add it to the diagram foreach(diagram) do cell @match cell begin @@ -80,7 +83,6 @@ function Decapode(diagram::AbstractVector{JSON3.Object}, theory::Theory; scalars return pode, anons, vars end export Decapode -# the proper name for this constructor should be "SummationDecapode" function uuid_to_symb(decapode::SummationDecapode, vars::Dict{String, Int}) Dict([key => (subpart(decapode, vars[key], :name)) for key ∈ keys(vars)]) diff --git a/packages/algjulia-service/src/decapodes-service/plotting.jl b/packages/algjulia-service/src/decapodes-service/plotting.jl deleted file mode 100644 index 4fc14bf0..00000000 --- a/packages/algjulia-service/src/decapodes-service/plotting.jl +++ /dev/null @@ -1,28 +0,0 @@ -## PLOTTING CODE - -abstract type AbstractPlotType end - -struct Heatmap <: AbstractPlotType end - -# TODO make length a conditional value so we can pass it in if we want -function Base.reshape(::Heatmap, data) - l = floor(Int64, sqrt(length(data))) - reshape(data, l, l) -end - -# TODO size fixed -# function at_time(sr::SimResult, t::Int) -# [sr.state[t][i,j][3] for i in 1:51 for j in 1:51] -# end - -# function show_heatmap(sr::SimResult, t::Int) -# data = at_time(result, t) -# ℓ = floor(Int64, sqrt(length(data))); -# reshaped = reshape(data, ℓ, ℓ) -# Plots.heatmap(1:51, 1:51, reshaped, clims=(minimum(data), maximum(data)); palette=:redsblues) -# end - -# @gif for t ∈ 1:length(result.time) -# show_heatmap(result, t) -# end every 5 - diff --git a/packages/algjulia-service/src/decapodes-service/theory.jl b/packages/algjulia-service/src/decapodes-service/theory.jl index a13d1149..5e7dbb80 100644 --- a/packages/algjulia-service/src/decapodes-service/theory.jl +++ b/packages/algjulia-service/src/decapodes-service/theory.jl @@ -22,18 +22,12 @@ function to_theory(theory::ThDecapode, type::HomType, name::String) "Δ" => :Δ "Δ⁻¹" => :Δ⁻¹ "d*" || "d̃₁" => :dual_d₁ - # \star on LHS - "⋆" => :⋆₁ + "⋆" || "⋆₁" || "★₁" || "★1" => :⋆₁ "⋆⁻¹" || "⋆₀⁻¹" => :⋆₀⁻¹ - # => :⋆₀⁻¹ - # \bigstar on LHS "★" || "★⁻¹" => :⋆₁ - # => :⋆₀⁻¹ "diffusivity" => :diffusivity - # new "d" || "d₀" || "d01" => :d₀ "d12" => :d₁ - "⋆1" => :⋆₁ "⋆2" => :⋆₂ "♭♯" => :♭♯ "lamb" => :dpsw # dual-primal self-wedge @@ -99,19 +93,16 @@ struct Theory{T<:ThDecapode} end export Theory -# TODO engooden -Base.show(io::IO, theory::Theory) = show(io, theory.data) - Base.values(theory::Theory) = values(theory.data) function add_to_theory! end; export add_to_theory! -function add_to_theory!(theory::Theory{T}, content::Any, type::ObType) where T +function add_to_theory!(theory::Theory{T}, content::AbstractDict, type::ObType) where T push!(theory.data, content[:id] => TheoryElement(;name=to_theory(T(), type, content[:name]))) end -function add_to_theory!(theory::Theory{T}, content::Any, type::HomType) where T +function add_to_theory!(theory::Theory{T}, content::AbstractDict, type::HomType) where T push!(theory.data, content[:id] => TheoryElement(;name=to_theory(T(), type, content[:name]), val=HomData(dom=content[:dom][:content], @@ -122,8 +113,8 @@ end # ...an object, we convert its type to a symbol and add it to the theorydict # ...a morphism, we add it to the theorydict with a field for the ids of its # domain and codomain to its -function Theory(model::AbstractVector{JSON3.Object}) - theory = Theory(ThDecapode()); +function Theory(model::AbstractVector{<:AbstractDict}) + theory = Theory(ThDecapode()) foreach(model) do cell @match cell begin IsObject(content) => add_to_theory!(theory, content, ObType()) diff --git a/packages/algjulia-service/src/decapodes.jl b/packages/algjulia-service/src/decapodes.jl deleted file mode 100644 index fd7fef7a..00000000 --- a/packages/algjulia-service/src/decapodes.jl +++ /dev/null @@ -1,416 +0,0 @@ -# - parameterize by the theory. this is currently fixed to decapodes -# - switch to different meshes -# - use enum instead of val - -# NOTES: -# TODO "anonymous objects" are •n - -# algebraicjulia dependencies -using ACSets -using DiagrammaticEquations -using Decapodes -using Decapodes: dec_mat_dual_differential, dec_mat_inverse_hodge -using CombinatorialSpaces - -# dependencies -import JSON3 -using StaticArrays -using MLStyle -using LinearAlgebra -using ComponentArrays -using Distributions # for initial conditions - -# meshing -using CoordRefSystems -using GeometryBasics: Point2, Point3 -Point3D = Point3{Float64}; - -using OrdinaryDiffEq - -export infer_types!, evalsim, default_dec_generate, default_dec_matrix_generate, - DiagonalHodge, ComponentArray - -struct ImplError <: Exception - name::String -end - -Base.showerror(io::IO, e::ImplError) = print(io, "$(e.name) not implemented") - -# funcitons for geometry and initial conditions -include("geometry_and_init.jl") - -## THEORY BUILDING - -""" Functions to build a dictionary associating ids in the theory to elements in the model""" -function to_decapode_theory end -export to_decapode_theory - -""" Helper function to convert CatColab values (Obs) in Decapodes """ -function to_decapode_theory(::Val{:Ob}, name::String) - @match lowercase(name) begin - "0-form" => :Form0 - "1-form" => :Form1 - "2-form" => :Form2 - "primal 0-form" => :Form0 - "primal 1-form" => :Form1 - "primal 2-form" => :Form2 - "dual 0-form" => :DualForm0 - "dual 1-form" => :DualForm1 - "dual 2-form" => :DualForm2 - x => throw(ImplError(x)) - end -end - -""" Helper function to convert CatColab values (Homs) in Decapodes """ -function to_decapode_theory(::Val{:Hom}, name::String) - @match replace(name," " => "") begin - "∂t" => :∂ₜ - "∂ₜ" => :∂ₜ - "Δ" => :Δ - "Δ⁻¹" => :Δ⁻¹ - "d" => :d₀ - "d*" => :dual_d₁ - "d̃₁" => :dual_d₁ - # \star on LHS - "⋆" => :⋆₁ - "⋆⁻¹" => :⋆₀⁻¹ - "⋆₀⁻¹" => :⋆₀⁻¹ - # \bigstar on LHS - "★" => :⋆₁ - "★⁻¹" => :⋆₀⁻¹ - "diffusivity" => :diffusivity - # new - "d₀" => :d₀ - "d12" => :d₁ - "⋆₁" => :⋆₁ - "⋆2" => :⋆₂ - "♭♯" => :♭♯ - "lamb" => :dpsw # dual-primal self-wedge - "-" => :neg - x => throw(ImplError(x)) - end -end - -# Build the theory - -#= -@active patterns are MLStyle-implementations of F# active patterns that forces us to work in the Maybe/Option pattern. -Practically, yet while a matter of opinion, they make @match statements cleaner; a statement amounts to a helpful pattern -name and the variables we intend to capture. -=# -@active IsObject(x) begin - x[:tag] == "object" ? Some(x) : nothing -end - -@active IsMorphism(x) begin - x[:tag] == "morphism" ? Some(x) : nothing -end - -export IsObject, IsMorphism - -""" Obs, Homs """ -abstract type ElementData end - -""" Struct capturing the name of the object and its relevant information. ElementData may be objects or homs, each of which has different data. -""" -struct TheoryElement - name::Union{Symbol, Nothing} - val::Union{ElementData, Nothing} - function TheoryElement(;name::Symbol=nothing,val::Any=nothing) - new(name, val) - end -end -export TheoryElement - -Base.nameof(t::TheoryElement) = t.name - -struct ObData <: ElementData end -# TODO not being used right now but added for completeness. - -struct HomData <: ElementData - dom::Any - cod::Any - function HomData(;dom::Any,cod::Any) - new(dom,cod) - end -end -export HomData -# TODO type dom/cod - -"""Struct wrapping a dictionary""" -struct Theory - data::Dict{String, TheoryElement} - function Theory() - new(Dict{String, TheoryElement}()) - end -end -export Theory - -# TODO engooden -Base.show(io::IO, theory::Theory) = show(io, theory.data) - -Base.values(theory::Theory) = values(theory.data) - -function add_to_theory! end -export add_to_theory! - -function add_to_theory!(theory::Theory, content::Any, type::Val{:Ob}) - push!(theory.data, - content[:id] => TheoryElement(;name=to_decapode_theory(type, content[:name]))) -end - -function add_to_theory!(theory::Theory, content::Any, type::Val{:Hom}) - push!(theory.data, content[:id] => - TheoryElement(;name=to_decapode_theory(type, content[:name]), - val=HomData(dom=content[:dom][:content], - cod=content[:cod][:content]))) -end - -# for each cell, if it is... -# ...an object, we convert its type to a symbol and add it to the theorydict -# ...a morphism, we add it to the theorydict with a field for the ids of its -# domain and codomain to its -function Theory(model::AbstractVector{JSON3.Object}) - theory = Theory(); - foreach(model) do cell - @match cell begin - IsObject(content) => add_to_theory!(theory, content, Val(:Ob)) - IsMorphism(content) => add_to_theory!(theory, content, Val(:Hom)) - x => throw(ImplError(x)) - end - end - return theory -end -export Theory - -## MODEL BUILDING - -function add_to_pode! end -export add_to_pode! - -function add_to_pode!(d::SummationDecapode, - vars::Dict{String, Int}, # mapping between UUID and ACSet ID - theory::Theory, - content::JSON3.Object, - nc::Vector{Int}, - ::Val{:Ob}) - theory_elem = theory.data[content[:over][:content]] # indexes the theory by UUID - # checks if the cell is an anonymous (intermediate) variable. - # if so, we increment the intermediate variable counter and make an intermediate variable name. - # otherwise we use the existing name of the given content. - name = if isempty(content[:name]) - nc[1] += 1 - Symbol("•$(nc[1])") - else - Symbol(content[:name]) - end - id = add_part!(d, :Var, name=name, type=nameof(theory_elem)) - push!(vars, content[:id] => id) - return d -end - -function op1_name(theory::Theory, content::JSON3.Object) - Symbol(theory.data[content[:over][:content]].name) -end - -# TODO we are restricted to Op1 -function add_to_pode!(d::SummationDecapode, - vars::Dict{String, Int}, # mapping between UUID and ACSet ID - theory::Theory, - content::JSON3.Object, - scalars::Any, - anons::Dict{Symbol, Any}, - ::Val{:Hom}) - dom = content[:dom][:content]; cod = content[:cod][:content] - # TODO we need a safe way to fail this - if haskey(vars, dom) && haskey(vars, cod) - # get the name of the Op1 and add it to the theory - op1 = op1_name(theory, content) - add_part!(d, :Op1, src=vars[dom], tgt=vars[cod], op1=op1) - # we need to add an inclusion to the TVar table - if op1 == :∂ₜ - add_part!(d, :TVar, incl=vars[cod]) - end - # if the dom is anonymous, we treat it as a something which will receive x -> k * x. - # we store its value in another array - if !isempty(scalars) && haskey(scalars, Symbol(content[:over][:content])) - scalar = scalars[Symbol(content[:over][:content])] - push!(anons, op1 => x -> scalar * x) - end - # TODO if scalars were typed correctly, we could probably do away with the !isempty check - end - d -end - -""" Decapode(jsondiagram::JSON3.Object, theory::Theory) => SummationDecapode - -This returns a Decapode given a jsondiagram and a theory. -""" -function Decapode(diagram::AbstractVector{JSON3.Object}, theory::Theory; scalars=[]) - # initiatize decapode and its mapping between UUIDs and ACSet IDs - pode = SummationDecapode(parse_decapode(quote end)); - vars = Dict{String, Int}(); # UUID => ACSetID - nc = [0]; # array is a mutable container - anons = Dict{Symbol, Any}(); - # for each cell in the notebook, add it to the diagram - foreach(diagram) do cell - @match cell begin - # TODO merge nameless_count into vars - IsObject(content) => add_to_pode!(pode, vars, theory, content, nc, Val(:Ob)) - IsMorphism(content) => add_to_pode!(pode, vars, theory, content, scalars, anons, Val(:Hom)) - _ => throw(ImplError(cell[:content][:tag])) - end - end - return pode, anons, vars -end -export Decapode -# the proper name for this constructor should be "SummationDecapode" - -function init_conditions(vars::Vector{Symbol}, sd::HasDeltaSet) - c_dist = MvNormal([50, 50], LinearAlgebra.Diagonal(map(abs2, [7.5, 7.5]))) - c = [pdf(c_dist, [p[1], p[2]]) for p in sd[:point]] - u0 = ComponentArray(; Dict([var=>c for var in vars])...) - return u0 -end - -struct PodeSystem - pode::SummationDecapode - plotvar::Vector{Symbol} - scalars::Dict{Symbol, Any} # closures # TODO rename scalars => anons - mesh::HasDeltaSet - dualmesh::HasDeltaSet - init::ComponentArray # TODO Is it always ComponentVector? - generate::Any - uuiddict::Dict{Symbol, String} -end -export PodeSystem - -""" -Construct a `PodeSystem` object from a JSON string. -""" -function PodeSystem(json_string::String,hodge=GeometricHodge()) - json_object = JSON3.read(json_string); - - # converts the JSON of (the fragment of) the theory - # into theory of the DEC, valued in Julia - theory = Theory(json_object[:model]); - - # this is a diagram in the model of the DEC. it wants to be a decapode! - diagram = json_object[:diagram]; - - # any scalars? - scalars = haskey(json_object, :scalars) ? json_object[:scalars] : []; - - # pode, anons, and vars (UUID => ACSetId) - decapode, anons, vars = Decapode(diagram, theory; scalars=scalars); - dot_rename!(decapode) - uuid2symb = Dict( - [key => (subpart(decapode, vars[key], :name)) for key ∈ keys(vars)] - ); - # plotting variables - plotvars = [uuid2symb[uuid] for uuid in json_object[:plotVariables]]; - - # create the mesh - mesh_name = Symbol(json_object[:mesh]) - mesh_builder = predefined_meshes[mesh_name] - s, sd = create_mesh(mesh_builder) - - # initialize operators - ♭♯_m = ♭♯_mat(sd); - wedge_dp10 = dec_wedge_product_dp(Tuple{1,0}, sd); - dual_d1_m = dec_mat_dual_differential(1, sd); - star0_inv_m = dec_mat_inverse_hodge(0, sd, hodge) - Δ0 = Δ(0,sd); - #fΔ0 = factorize(Δ0); - # end initialize - - # initial conditions - ic_specs = json_object[:initialConditions]; - # Dict(:duu => "TaylorVortex") - dict = Dict([uuid2symb[string(uuid)] => ic_specs[string(uuid)] for uuid ∈ keys(ic_specs)]...) - u0 = initial_conditions(dict, mesh_builder, sd) - - function sys_generate(s, my_symbol,hodge=hodge) - op = @match my_symbol begin - sym && if sym ∈ keys(anons) end => anons[sym] - :♭♯ => x -> ♭♯_m * x # [1] - :dpsw => x -> wedge_dp10(x, star0_inv_m[1]*(dual_d1_m[1]*x)) - :Δ⁻¹ => x -> begin - y = Δ0 \ x - y .- minimum(y) - end - _ => default_dec_matrix_generate(s, my_symbol, hodge) - end - return (args...) -> op(args...) - end - - # symbol => uuid. we need this to reassociate the var - symb2uuid = Dict([v => k for (k,v) in pairs(uuid2symb)]) - - return PodeSystem(decapode, plotvars, anons, s, sd, u0, sys_generate, symb2uuid) -end -export PodeSystem - -function run_sim(fm, u0, t0, constparam) - prob = ODEProblem(fm, u0, (0, t0), constparam) - soln = solve(prob, Tsit5(), saveat=0.01) -end -export run_sim - -struct SimResult - time::Vector{Float64} - state::Dict{String, Vector{AbstractArray{SVector{3, Float64}}}} - x::Vector{Float64} # axis - y::Vector{Float64} -end -export SimResult - -function SimResult(sol::ODESolution, system::PodeSystem) - - points = collect(values(system.mesh.subparts.point.m)); - geometry = length(points[1]) == 2 ? :Rect : :Sphere - - function grid(pt3,grid_size) - pt2 = [(pt3[1]+1)/2,(pt3[2]+1)/2] - [round(Int,pt2[1]*grid_size[1]),round(Int,pt2[2]*grid_size[2])] - end - l = 0 - function at_time(sol::ODESolution, plotvar::Symbol, timeidx::Int) - if geometry == :Rect - l = 51 - return [SVector(i, j, getproperty(sol.u[timeidx], plotvar)[l*(i-1) + j]) for i in 1:l, j in 1:l] - else - northern_indices = filter(i -> points[i][3]>0,keys(points)) - l = 100 - return map(northern_indices) do n - i,j = grid(points[n], [l,l]) - SVector(i, j, getproperty(sol.u[timeidx], plotvar)[n]) - end - end - end - - function state_vals(plotvar::Symbol) - map(1:length(sol.t)) do i - at_time(sol, plotvar, i) - end - end - state_val_dict = Dict([(system.uuiddict[plotvar] => state_vals(plotvar)) for plotvar in system.plotvar]) - - # TODO engooden - # Dict("UUID1" => VectorMatrixSVectr...) - SimResult(sol.t, state_val_dict, 0:l-1, 0:l-1) -end -# TODO generalize to HasDeltaSet - -## PLOTTING CODE - -abstract type AbstractPlotType end - -struct Heatmap <: AbstractPlotType end - -# TODO make length a conditional value so we can pass it in if we want -function Base.reshape(::Heatmap, data) - l = floor(Int64, sqrt(length(data))) - reshape(data, l, l) -end -#∧ᵈᵖ₁₀(-,⋆d(-)) diff --git a/packages/algjulia-service/test/runtests.jl b/packages/algjulia-service/test/runtests.jl index ef72a559..c50bf13c 100644 --- a/packages/algjulia-service/test/runtests.jl +++ b/packages/algjulia-service/test/runtests.jl @@ -13,15 +13,12 @@ using StaticArrays using LinearAlgebra import OrdinaryDiffEq: ReturnCode -# visualization -#using Plots - const KEYS = Set([:mesh, :plotVariables, :initialConditions, :domain, :diagram, :model, :scalars]) # load data data = open(JSON3.read, joinpath(@__DIR__, "test_jsons", "diffusion_data.json"), "r") -diagram = data[:diagram]; -model = data[:model]; +diagram = data[:diagram] +model = data[:model] @testset "Text-to-Pode" begin @@ -55,31 +52,31 @@ end _ => false end - theory = Theory(ThDecapode()); + theory = Theory(ThDecapode()) @match model[1] begin IsObject(content) => add_to_theory!(theory, content, ObType()) _ => nothing end - _id = "019323fa-49cb-7373-8c5d-c395bae4006d"; + _id = "019323fa-49cb-7373-8c5d-c395bae4006d" @test theory.data[_id] == TheoryElement(;name=:Form0, val=nothing) end @testset "Making the Decapode" begin - theory = Theory(model); + theory = Theory(model) @test Set(nameof.(values(theory))) == Set([:Form0, :Form1, :Form2, :Δ, :∂ₜ]) - handcrafted_pode = SummationDecapode(parse_decapode(quote end)); - add_part!(handcrafted_pode, :Var, name=:C, type=:Form0); - add_part!(handcrafted_pode, :Var, name=Symbol("dC/dt"), type=:Form0); - add_part!(handcrafted_pode, :TVar, incl=2); - add_part!(handcrafted_pode, :Op1, src=1, tgt=2, op1=:∂ₜ); - add_part!(handcrafted_pode, :Op1, src=1, tgt=2, op1=:Δ); + handcrafted_pode = SummationDecapode(parse_decapode(quote end)) + add_part!(handcrafted_pode, :Var, name=:C, type=:Form0) + add_part!(handcrafted_pode, :Var, name=Symbol("dC/dt"), type=:Form0) + add_part!(handcrafted_pode, :TVar, incl=2) + add_part!(handcrafted_pode, :Op1, src=1, tgt=2, op1=:∂ₜ) + add_part!(handcrafted_pode, :Op1, src=1, tgt=2, op1=:Δ) # no scalars in second position - decapode, _, _ = Decapode(diagram, theory); + decapode, _, _ = Decapode(diagram, theory) @test decapode == handcrafted_pode @@ -87,13 +84,13 @@ end @testset "Simulation - Diffusion" begin - json_string = read(joinpath(@__DIR__, "test_jsons", "diffusion_data.json"), String); + json_string = read(joinpath(@__DIR__, "test_jsons", "diffusion_data.json"), String) @test Set(keys(JSON3.read(json_string))) == KEYS - system = PodeSystem(json_string); + system = PodeSystem(json_string) - simulator = evalsim(system.pode); - f = simulator(system.geometry.dualmesh, default_dec_generate, DiagonalHodge()); + simulator = evalsim(system.pode) + f = simulator(system.geometry.dualmesh, default_dec_generate, DiagonalHodge()) soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)); # returns ::ODESolution @@ -104,11 +101,11 @@ end @test soln.retcode == ReturnCode.Success - result = SimResult(soln, system); + result = SimResult(soln, system) @test typeof(result.state) == Dict{String, Vector{AbstractArray{SVector{3, Float64}}}} - jv = JsonValue(result); + jv = JsonValue(result) end @@ -126,11 +123,6 @@ we are trying to index `soln.u[t]` by `Ċ `, but only `C` is present. I think t # f = simulator(system.geometry.dualmesh, system.generate, DiagonalHodge()) # soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)); -# # returns ::ODESolution -# # - retcode -# # - interpolation -# # - t -# # - u::Vector{ComponentVector} # @test soln.retcode == ReturnCode.Success @@ -145,8 +137,8 @@ we are trying to index `soln.u[t]` by `Ċ `, but only `C` is present. I think t @testset "Parsing the Theory JSON Object - Diffusion Long Trip" begin data = open(JSON3.read, joinpath(@__DIR__, "test_jsons", "diffusion_long_trip.json"), "r") - diagram = data[:diagram]; - model = data[:model]; + diagram = data[:diagram] + model = data[:model] @test Set(keys(data)) == KEYS @test @match model[1] begin @@ -159,7 +151,7 @@ we are trying to index `soln.u[t]` by `Ċ `, but only `C` is present. I think t _ => false end - theory = Theory(ThDecapode()); + theory = Theory(ThDecapode()) @match model[1] begin IsObject(content) => add_to_theory!(theory, content, ObType()) _ => nothing @@ -172,16 +164,16 @@ end # # GOOD @testset "Simulation - Diffusion Long Trip" begin - json_string = read(joinpath(@__DIR__, "test_jsons", "diffusion_long_trip.json"), String); + json_string = read(joinpath(@__DIR__, "test_jsons", "diffusion_long_trip.json"), String) @test Set(keys(JSON3.read(json_string))) == KEYS - system = PodeSystem(json_string); + system = PodeSystem(json_string) simulator = evalsim(system.pode) - f = simulator(system.geometry.dualmesh, default_dec_matrix_generate, DiagonalHodge()); + f = simulator(system.geometry.dualmesh, default_dec_matrix_generate, DiagonalHodge()) # time - soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)); + soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)) # returns ::ODESolution # - retcode # - interpolation @@ -190,11 +182,11 @@ end @test soln.retcode == ReturnCode.Success - result = SimResult(soln, system); + result = SimResult(soln, system) @test typeof(result.state) == Dict{String, Vector{AbstractArray{SVector{3, Float64}}}} - jv = JsonValue(result); + jv = JsonValue(result) end @@ -202,9 +194,9 @@ end @testset "Parsing the Theory JSON Object - Diffusivity Constant" begin data = open(JSON3.read, joinpath(@__DIR__, "test_jsons", "diffusivity_constant.json"), "r") - diagram = data[:diagram]; - model = data[:model]; - scalars = data[:scalars]; + diagram = data[:diagram] + model = data[:model] + scalars = data[:scalars] @test Set(keys(data)) == KEYS @@ -218,7 +210,7 @@ end _ => false end - theory = Theory(ThDecapode()); + theory = Theory(ThDecapode()) @match model[1] begin IsObject(content) => add_to_theory!(theory, content, ObType()) _ => nothing @@ -230,47 +222,37 @@ end @testset "Simulation - Diffusivity Constant" begin - json_string = read(joinpath(@__DIR__, "test_jsons", "diffusivity_constant.json"), String); - system = PodeSystem(json_string); + json_string = read(joinpath(@__DIR__, "test_jsons", "diffusivity_constant.json"), String) + system = PodeSystem(json_string) simulator = evalsim(system.pode) - f = simulator(system.geometry.dualmesh, system.generate, DiagonalHodge()); + f = simulator(system.geometry.dualmesh, system.generate, DiagonalHodge()) - soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)); - # returns ::ODESolution - # - retcode - # - interpolation - # - t - # - u::Vector{ComponentVector} + soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)) @test soln.retcode == ReturnCode.Success - result = SimResult(soln, system); + result = SimResult(soln, system) @test typeof(result.state) == Dict{String, Vector{AbstractArray{SVector{3, Float64}}}} - jv = JsonValue(result); + jv = JsonValue(result) end @testset "Simulation - Navier-Stokes Vorticity" begin - json_string = read(joinpath(@__DIR__, "test_jsons", "ns_vorticity.json"), String); + json_string = read(joinpath(@__DIR__, "test_jsons", "ns_vorticity.json"), String) @test Set(keys(JSON3.read(json_string))) == KEYS - system = PodeSystem(json_string); + system = PodeSystem(json_string) simulator = evalsim(system.pode) - f = simulator(system.geometry.dualmesh, system.generate, DiagonalHodge()); + f = simulator(system.geometry.dualmesh, system.generate, DiagonalHodge()) - soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)); - # returns ::ODESolution - # - retcode - # - interpolation - # - t - # - u::Vector{ComponentVector} + soln = run_sim(f, system.init, 50.0, ComponentArray(k=0.5,)) @test soln.retcode == ReturnCode.Success @@ -278,6 +260,6 @@ end @test typeof(result.state) == Dict{String, Vector{AbstractArray{SVector{3, Float64}}}} - jv = JsonValue(result); + jv = JsonValue(result) end