From aed87b90013c4c5d96567fe45d500dd950071b81 Mon Sep 17 00:00:00 2001 From: Braam van Dyk Date: Fri, 10 Nov 2023 23:09:45 +0200 Subject: [PATCH 1/2] Add data cleaning --- Manifest.toml | 193 +++++++++++++++-- Project.toml | 12 +- docs/Manifest.toml | 491 ++++++++++++++++++++++++++++++++++++++++++-- src/datacleaning.jl | 93 +++++++++ test/Manifest.toml | 2 +- 5 files changed, 758 insertions(+), 33 deletions(-) create mode 100644 src/datacleaning.jl diff --git a/Manifest.toml b/Manifest.toml index 8d187b9..9bb412c 100644 --- a/Manifest.toml +++ b/Manifest.toml @@ -1,21 +1,19 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.9.2" +julia_version = "1.9.3" manifest_format = "2.0" -project_hash = "193f8d3e099802ba629129d15a5ddcd7bb658227" +project_hash = "386e3531f9b1b93f33ab1ccde1ce009c6c3d47bb" [[deps.Adapt]] deps = ["LinearAlgebra", "Requires"] git-tree-sha1 = "76289dc51920fdc6e0013c872ba9551d54961c24" uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" version = "3.6.2" +weakdeps = ["StaticArrays"] [deps.Adapt.extensions] AdaptStaticArraysExt = "StaticArrays" - [deps.Adapt.weakdeps] - StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" - [[deps.ArgTools]] uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" version = "1.1.1" @@ -45,9 +43,31 @@ version = "7.4.11" [[deps.Artifacts]] uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" +[[deps.AxisAlgorithms]] +deps = ["LinearAlgebra", "Random", "SparseArrays", "WoodburyMatrices"] +git-tree-sha1 = "66771c8d21c8ff5e3a93379480a2307ac36863f7" +uuid = "13072b0f-2c55-5437-9ae7-d433b7a33950" +version = "1.0.1" + [[deps.Base64]] uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" +[[deps.Calculus]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "f641eb0a4f00c343bbc32346e1217b86f3ce9dad" +uuid = "49dc2e85-a5d0-5ad3-a950-438e2897f1b9" +version = "0.5.1" + +[[deps.ChainRulesCore]] +deps = ["Compat", "LinearAlgebra"] +git-tree-sha1 = "e0af648f0692ec1691b5d094b8724ba1346281cf" +uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" +version = "1.18.0" +weakdeps = ["SparseArrays"] + + [deps.ChainRulesCore.extensions] + ChainRulesCoreSparseArraysExt = "SparseArrays" + [[deps.CommonSubexpressions]] deps = ["MacroTools", "Test"] git-tree-sha1 = "7b8a93dba8af7e3b42fecabf646260105ac373f7" @@ -56,9 +76,9 @@ version = "0.3.0" [[deps.Compat]] deps = ["UUIDs"] -git-tree-sha1 = "5ce999a19f4ca23ea484e92a1774a61b8ca4cf8e" +git-tree-sha1 = "e460f044ca8b99be31d35fe54fc33a5c33dd8ed7" uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" -version = "4.8.0" +version = "4.9.0" weakdeps = ["Dates", "LinearAlgebra"] [deps.Compat.extensions] @@ -95,9 +115,9 @@ version = "1.15.0" [[deps.DataStructures]] deps = ["Compat", "InteractiveUtils", "OrderedCollections"] -git-tree-sha1 = "cf25ccb972fec4e4817764d01c82386ae94f77b4" +git-tree-sha1 = "3dbd312d370723b6bb43ba9d02fc36abade4518d" uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" -version = "0.18.14" +version = "0.18.15" [[deps.DataValueInterfaces]] git-tree-sha1 = "bfc1187b79289637fa0ef6d4436ebdfe6905cbd6" @@ -126,10 +146,37 @@ git-tree-sha1 = "23163d55f885173722d1e4cf0f6110cdbaf7e272" uuid = "b552c78f-8df3-52c6-915a-8e097449b14b" version = "1.15.1" +[[deps.Distances]] +deps = ["LinearAlgebra", "Statistics", "StatsAPI"] +git-tree-sha1 = "5225c965635d8c21168e32a12954675e7bea1151" +uuid = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7" +version = "0.10.10" +weakdeps = ["ChainRulesCore", "SparseArrays"] + + [deps.Distances.extensions] + DistancesChainRulesCoreExt = "ChainRulesCore" + DistancesSparseArraysExt = "SparseArrays" + [[deps.Distributed]] deps = ["Random", "Serialization", "Sockets"] uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" +[[deps.Distributions]] +deps = ["FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SpecialFunctions", "Statistics", "StatsAPI", "StatsBase", "StatsFuns"] +git-tree-sha1 = "a6c00f894f24460379cb7136633cef54ac9f6f4a" +uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" +version = "0.25.103" + + [deps.Distributions.extensions] + DistributionsChainRulesCoreExt = "ChainRulesCore" + DistributionsDensityInterfaceExt = "DensityInterface" + DistributionsTestExt = "Test" + + [deps.Distributions.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d" + Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" + [[deps.DocStringExtensions]] deps = ["LibGit2"] git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d" @@ -141,6 +188,12 @@ deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" version = "1.6.0" +[[deps.DualNumbers]] +deps = ["Calculus", "NaNMath", "SpecialFunctions"] +git-tree-sha1 = "5837a837389fccf076445fce071c8ddaea35a566" +uuid = "fa6b7ba4-c1ee-5f82-b5fc-ecf0adba8f74" +version = "0.6.8" + [[deps.FileWatching]] uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" @@ -168,24 +221,34 @@ version = "2.21.1" [[deps.ForwardDiff]] deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions"] -git-tree-sha1 = "00e252f4d706b3d55a8863432e742bf5717b498d" +git-tree-sha1 = "cf0fe81336da9fb90944683b8c41984b08793dad" uuid = "f6369f11-7733-5829-9624-2563aa707210" -version = "0.10.35" +version = "0.10.36" +weakdeps = ["StaticArrays"] [deps.ForwardDiff.extensions] ForwardDiffStaticArraysExt = "StaticArrays" - [deps.ForwardDiff.weakdeps] - StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" - [[deps.Future]] deps = ["Random"] uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" +[[deps.HypergeometricFunctions]] +deps = ["DualNumbers", "LinearAlgebra", "OpenLibm_jll", "SpecialFunctions"] +git-tree-sha1 = "f218fe3736ddf977e0e772bc9a586b2383da2685" +uuid = "34004b35-14d8-5ef3-9330-4cdb6864b03a" +version = "0.3.23" + [[deps.InteractiveUtils]] deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" +[[deps.Interpolations]] +deps = ["Adapt", "AxisAlgorithms", "ChainRulesCore", "LinearAlgebra", "OffsetArrays", "Random", "Ratios", "Requires", "SharedArrays", "SparseArrays", "StaticArrays", "WoodburyMatrices"] +git-tree-sha1 = "721ec2cf720536ad005cb38f50dbba7b02419a15" +uuid = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" +version = "0.14.7" + [[deps.IrrationalConstants]] git-tree-sha1 = "630b497eafcc20001bba38a4651b327dcfc491d2" uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" @@ -239,6 +302,12 @@ version = "7.2.0" deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +[[deps.Loess]] +deps = ["Distances", "LinearAlgebra", "Statistics", "StatsAPI"] +git-tree-sha1 = "a113a8be4c6d0c64e217b472fb6e61c760eb4022" +uuid = "4345ca2d-374a-55d4-8d30-97f9976e7612" +version = "0.6.3" + [[deps.LogExpFunctions]] deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"] git-tree-sha1 = "c3ce8e7420b3a6e071e0fe4745f5d4300e37b13f" @@ -302,6 +371,12 @@ version = "1.0.2" uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" version = "1.2.0" +[[deps.OffsetArrays]] +deps = ["Adapt"] +git-tree-sha1 = "2ac17d29c523ce1cd38e27785a7d23024853a4bb" +uuid = "6fe1bfb0-de20-5000-8ca7-80f57d26f881" +version = "1.12.10" + [[deps.OpenBLAS_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" @@ -329,6 +404,12 @@ git-tree-sha1 = "2e73fe17cac3c62ad1aebe70d44c963c3cfdc3e3" uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" version = "1.6.2" +[[deps.PDMats]] +deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"] +git-tree-sha1 = "f6f85a2edb9c356b829934ad3caed2ad0ebbfc99" +uuid = "90014a1f-27ba-587c-ab20-58faa44d9150" +version = "0.11.29" + [[deps.Parameters]] deps = ["OrderedCollections", "UnPack"] git-tree-sha1 = "34c0e9ad262e5f7fc75b10a9952ca7692cfc5fbe" @@ -346,6 +427,12 @@ git-tree-sha1 = "17275485f373e6673f7e7f97051f703ed5b15b20" uuid = "85a6dd25-e78a-55b7-8502-1745935b8125" version = "0.2.4" +[[deps.PrecompileTools]] +deps = ["Preferences"] +git-tree-sha1 = "9673d39decc5feece56ef3940e5dafba15ba0f81" +uuid = "aea7be01-6a6a-4083-8856-8a6e6704d82a" +version = "1.1.2" + [[deps.Preferences]] deps = ["TOML"] git-tree-sha1 = "7eb1686b4f04b82f96ed7a4ea5890a4f0c7a09f1" @@ -362,6 +449,12 @@ version = "2.2.7" deps = ["Unicode"] uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7" +[[deps.QuadGK]] +deps = ["DataStructures", "LinearAlgebra"] +git-tree-sha1 = "9ebcd48c498668c7fa0e97a9cae873fbee7bfee1" +uuid = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" +version = "2.9.1" + [[deps.REPL]] deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" @@ -370,6 +463,24 @@ uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" deps = ["SHA", "Serialization"] uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +[[deps.Ratios]] +deps = ["Requires"] +git-tree-sha1 = "1342a47bf3260ee108163042310d26f2be5ec90b" +uuid = "c84ed2f1-dad5-54f0-aa8e-dbefe2724439" +version = "0.4.5" + + [deps.Ratios.extensions] + RatiosFixedPointNumbersExt = "FixedPointNumbers" + + [deps.Ratios.weakdeps] + FixedPointNumbers = "53c48c17-4a7d-5ca2-90c5-79b7896eea93" + +[[deps.RecipesBase]] +deps = ["PrecompileTools"] +git-tree-sha1 = "5c3d09cc4f31f5fc6af001c250bf1278733100ff" +uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" +version = "1.3.4" + [[deps.Reexport]] git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b" uuid = "189a3867-3050-52da-a836-e630ba90ab69" @@ -381,6 +492,18 @@ git-tree-sha1 = "838a3a4188e2ded87a4f9f184b4b0d78a1e91cb7" uuid = "ae029012-a4dd-5104-9daa-d747884805df" version = "1.3.0" +[[deps.Rmath]] +deps = ["Random", "Rmath_jll"] +git-tree-sha1 = "f65dcb5fa46aee0cf9ed6274ccbd597adc49aa7b" +uuid = "79098fc4-a85e-5d69-aa6a-4863f24498fa" +version = "0.7.1" + +[[deps.Rmath_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "6ed52fdd3382cf21947b15e8870ac0ddbff736da" +uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f" +version = "0.4.0+0" + [[deps.SHA]] uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" version = "0.7.0" @@ -394,6 +517,10 @@ git-tree-sha1 = "e2cc6d8c88613c05e1defb55170bf5ff211fbeac" uuid = "efcf1570-3423-57d1-acb7-fd33fddbac46" version = "1.1.1" +[[deps.SharedArrays]] +deps = ["Distributed", "Mmap", "Random", "Serialization"] +uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383" + [[deps.Sockets]] uuid = "6462fe0b-24de-5631-8697-dd941f90decc" @@ -412,12 +539,20 @@ deps = ["IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_j git-tree-sha1 = "7beb031cf8145577fbccacd94b8a8f4ce78428d3" uuid = "276daf66-3868-5448-9aa4-cd146d93841b" version = "2.3.0" +weakdeps = ["ChainRulesCore"] [deps.SpecialFunctions.extensions] SpecialFunctionsChainRulesCoreExt = "ChainRulesCore" - [deps.SpecialFunctions.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" +[[deps.StaticArrays]] +deps = ["LinearAlgebra", "Random", "StaticArraysCore"] +git-tree-sha1 = "0adf069a2a490c47273727e029371b31d44b72b2" +uuid = "90137ffa-7385-5640-81b9-e52037218182" +version = "1.6.5" +weakdeps = ["Statistics"] + + [deps.StaticArrays.extensions] + StaticArraysStatisticsExt = "Statistics" [[deps.StaticArraysCore]] git-tree-sha1 = "36b3d696ce6366023a0ea192b4cd442268995a0d" @@ -441,6 +576,20 @@ git-tree-sha1 = "75ebe04c5bed70b91614d684259b661c9e6274a4" uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" version = "0.34.0" +[[deps.StatsFuns]] +deps = ["HypergeometricFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] +git-tree-sha1 = "f625d686d5a88bcd2b15cd81f18f98186fdc0c9a" +uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c" +version = "1.3.0" + + [deps.StatsFuns.extensions] + StatsFunsChainRulesCoreExt = "ChainRulesCore" + StatsFunsInverseFunctionsExt = "InverseFunctions" + + [deps.StatsFuns.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" + [[deps.StringManipulation]] git-tree-sha1 = "46da2434b41f41ac3594ee9816ce5541c6096123" uuid = "892a3eda-7b42-436c-8928-eab12a02cf0e" @@ -481,6 +630,12 @@ version = "1.10.0" deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +[[deps.TimeSeries]] +deps = ["Dates", "DelimitedFiles", "DocStringExtensions", "RecipesBase", "Reexport", "Statistics", "Tables"] +git-tree-sha1 = "8b9288d84da88ea44693ca8cf9c236da1778f274" +uuid = "9e3dc215-6440-5c97-bce1-76c03772f85e" +version = "0.23.2" + [[deps.UUIDs]] deps = ["Random", "SHA"] uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" @@ -493,6 +648,12 @@ version = "1.0.2" [[deps.Unicode]] uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" +[[deps.WoodburyMatrices]] +deps = ["LinearAlgebra", "SparseArrays"] +git-tree-sha1 = "de67fa59e33ad156a590055375a30b23c40299d3" +uuid = "efce3f68-66dc-5838-9240-27a6d6f5f9b6" +version = "0.5.5" + [[deps.Zlib_jll]] deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" diff --git a/Project.toml b/Project.toml index 40ea776..60868f7 100644 --- a/Project.toml +++ b/Project.toml @@ -7,21 +7,25 @@ version = "0.1.0" Base64 = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" Dates = "ade2ca70-3891-5945-98fb-dc099432e06a" DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab" +Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" Downloads = "f43a241f-c20a-4ad4-852c-f6b1247861c6" +Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" +Loess = "4345ca2d-374a-55d4-8d30-97f9976e7612" MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" +Missings = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28" Optim = "429524aa-4258-5aef-a3af-852621145aeb" PrettyTables = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +TimeSeries = "9e3dc215-6440-5c97-bce1-76c03772f85e" [compat] -julia = "1" DelimitedFiles = "1" Downloads = "1" -Statistics = "1" -Optim = "1.7" MacroTools = "0.5" +Optim = "1.7" PrettyTables = "2" - +Statistics = "1" +julia = "1" [extras] Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/docs/Manifest.toml b/docs/Manifest.toml index 06b4c13..5963728 100644 --- a/docs/Manifest.toml +++ b/docs/Manifest.toml @@ -1,75 +1,403 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.7.3" +julia_version = "1.9.2" manifest_format = "2.0" +project_hash = "2fc9abc4c3c6874720ed0957ebb2c422cd6263a5" [[deps.ANSIColoredPrinters]] git-tree-sha1 = "574baf8110975760d391c710b6341da1afa48d8c" uuid = "a4c015fc-c6ff-483c-b24f-f7ea428134e9" version = "0.0.1" +[[deps.Adapt]] +deps = ["LinearAlgebra", "Requires"] +git-tree-sha1 = "76289dc51920fdc6e0013c872ba9551d54961c24" +uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" +version = "3.6.2" + + [deps.Adapt.extensions] + AdaptStaticArraysExt = "StaticArrays" + + [deps.Adapt.weakdeps] + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + +[[deps.ArgTools]] +uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" +version = "1.1.1" + +[[deps.ArrayInterface]] +deps = ["Adapt", "LinearAlgebra", "Requires", "SparseArrays", "SuiteSparse"] +git-tree-sha1 = "f83ec24f76d4c8f525099b2ac475fc098138ec31" +uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" +version = "7.4.11" + + [deps.ArrayInterface.extensions] + ArrayInterfaceBandedMatricesExt = "BandedMatrices" + ArrayInterfaceBlockBandedMatricesExt = "BlockBandedMatrices" + ArrayInterfaceCUDAExt = "CUDA" + ArrayInterfaceGPUArraysCoreExt = "GPUArraysCore" + ArrayInterfaceStaticArraysCoreExt = "StaticArraysCore" + ArrayInterfaceTrackerExt = "Tracker" + + [deps.ArrayInterface.weakdeps] + BandedMatrices = "aae01518-5342-5314-be14-df237901396f" + BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" + CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" + GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" + StaticArraysCore = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" + Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" + +[[deps.Artifacts]] +uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" + [[deps.Base64]] uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" +[[deps.CommonSubexpressions]] +deps = ["MacroTools", "Test"] +git-tree-sha1 = "7b8a93dba8af7e3b42fecabf646260105ac373f7" +uuid = "bbf7d656-a473-5ed7-a52c-81e309532950" +version = "0.3.0" + +[[deps.Compat]] +deps = ["UUIDs"] +git-tree-sha1 = "e460f044ca8b99be31d35fe54fc33a5c33dd8ed7" +uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" +version = "4.9.0" +weakdeps = ["Dates", "LinearAlgebra"] + + [deps.Compat.extensions] + CompatLinearAlgebraExt = "LinearAlgebra" + +[[deps.CompilerSupportLibraries_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" +version = "1.0.5+0" + +[[deps.ConstructionBase]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "fe2838a593b5f776e1597e086dcd47560d94e816" +uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" +version = "1.5.3" + + [deps.ConstructionBase.extensions] + ConstructionBaseIntervalSetsExt = "IntervalSets" + ConstructionBaseStaticArraysExt = "StaticArrays" + + [deps.ConstructionBase.weakdeps] + IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953" + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + +[[deps.Crayons]] +git-tree-sha1 = "249fe38abf76d48563e2f4556bebd215aa317e15" +uuid = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f" +version = "4.1.1" + +[[deps.DataAPI]] +git-tree-sha1 = "8da84edb865b0b5b0100c0666a9bc9a0b71c553c" +uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a" +version = "1.15.0" + +[[deps.DataStructures]] +deps = ["Compat", "InteractiveUtils", "OrderedCollections"] +git-tree-sha1 = "3dbd312d370723b6bb43ba9d02fc36abade4518d" +uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" +version = "0.18.15" + +[[deps.DataValueInterfaces]] +git-tree-sha1 = "bfc1187b79289637fa0ef6d4436ebdfe6905cbd6" +uuid = "e2d170a0-9d28-54be-80f0-106bbe20a464" +version = "1.0.0" + [[deps.Dates]] deps = ["Printf"] uuid = "ade2ca70-3891-5945-98fb-dc099432e06a" +[[deps.DelimitedFiles]] +deps = ["Mmap"] +git-tree-sha1 = "9e2f36d3c96a820c678f2f1f1782582fcf685bae" +uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab" +version = "1.9.1" + +[[deps.DiffResults]] +deps = ["StaticArraysCore"] +git-tree-sha1 = "782dd5f4561f5d267313f23853baaaa4c52ea621" +uuid = "163ba53b-c6d8-5494-b064-1a9d43ac40c5" +version = "1.1.0" + +[[deps.DiffRules]] +deps = ["IrrationalConstants", "LogExpFunctions", "NaNMath", "Random", "SpecialFunctions"] +git-tree-sha1 = "23163d55f885173722d1e4cf0f6110cdbaf7e272" +uuid = "b552c78f-8df3-52c6-915a-8e097449b14b" +version = "1.15.1" + +[[deps.Distributed]] +deps = ["Random", "Serialization", "Sockets"] +uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" + [[deps.DocStringExtensions]] deps = ["LibGit2"] -git-tree-sha1 = "5158c2b41018c5f7eb1470d558127ac274eca0c9" +git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d" uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" -version = "0.9.1" +version = "0.9.3" [[deps.Documenter]] deps = ["ANSIColoredPrinters", "Base64", "Dates", "DocStringExtensions", "IOCapture", "InteractiveUtils", "JSON", "LibGit2", "Logging", "Markdown", "REPL", "Test", "Unicode"] -git-tree-sha1 = "ee945ed10767de73cf8b610951b98228cb65af80" +git-tree-sha1 = "39fd748a73dce4c05a9655475e437170d8fb1b67" uuid = "e30172f5-a6a5-5a46-863b-614d45cd2de4" -version = "0.27.22" +version = "0.27.25" + +[[deps.Downloads]] +deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] +uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" +version = "1.6.0" + +[[deps.FileWatching]] +uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" + +[[deps.FillArrays]] +deps = ["LinearAlgebra", "Random", "SparseArrays", "Statistics"] +git-tree-sha1 = "f372472e8672b1d993e93dada09e23139b509f9e" +uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" +version = "1.5.0" + +[[deps.FiniteDiff]] +deps = ["ArrayInterface", "LinearAlgebra", "Requires", "Setfield", "SparseArrays"] +git-tree-sha1 = "c6e4a1fbe73b31a3dea94b1da449503b8830c306" +uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" +version = "2.21.1" + + [deps.FiniteDiff.extensions] + FiniteDiffBandedMatricesExt = "BandedMatrices" + FiniteDiffBlockBandedMatricesExt = "BlockBandedMatrices" + FiniteDiffStaticArraysExt = "StaticArrays" + + [deps.FiniteDiff.weakdeps] + BandedMatrices = "aae01518-5342-5314-be14-df237901396f" + BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" [[deps.FlowsheetTools]] +deps = ["Base64", "Dates", "DelimitedFiles", "Downloads", "MacroTools", "Optim", "PrettyTables", "Statistics", "TimeSeries"] path = ".." uuid = "1024d294-f334-4c0a-aed8-72739789c799" version = "0.1.0" +[[deps.ForwardDiff]] +deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions"] +git-tree-sha1 = "cf0fe81336da9fb90944683b8c41984b08793dad" +uuid = "f6369f11-7733-5829-9624-2563aa707210" +version = "0.10.36" + + [deps.ForwardDiff.extensions] + ForwardDiffStaticArraysExt = "StaticArrays" + + [deps.ForwardDiff.weakdeps] + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + +[[deps.Future]] +deps = ["Random"] +uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" + [[deps.IOCapture]] deps = ["Logging", "Random"] -git-tree-sha1 = "f7be53659ab06ddc986428d3a9dcc95f6fa6705a" +git-tree-sha1 = "d75853a0bdbfb1ac815478bacd89cd27b550ace6" uuid = "b5f81e59-6552-4d32-b1f0-c071b021bf89" -version = "0.2.2" +version = "0.2.3" [[deps.InteractiveUtils]] deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" +[[deps.IrrationalConstants]] +git-tree-sha1 = "630b497eafcc20001bba38a4651b327dcfc491d2" +uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" +version = "0.2.2" + +[[deps.IteratorInterfaceExtensions]] +git-tree-sha1 = "a3f24677c21f5bbe9d2a714f95dcd58337fb2856" +uuid = "82899510-4779-5014-852e-03e436cf321d" +version = "1.0.0" + +[[deps.JLLWrappers]] +deps = ["Preferences"] +git-tree-sha1 = "abc9885a7ca2052a736a600f7fa66209f96506e1" +uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210" +version = "1.4.1" + [[deps.JSON]] deps = ["Dates", "Mmap", "Parsers", "Unicode"] -git-tree-sha1 = "3c837543ddb02250ef42f4738347454f95079d4e" +git-tree-sha1 = "31e996f0a15c7b280ba9f76636b3ff9e2ae58c9a" uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6" -version = "0.21.3" +version = "0.21.4" + +[[deps.LaTeXStrings]] +git-tree-sha1 = "f2355693d6778a178ade15952b7ac47a4ff97996" +uuid = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f" +version = "1.3.0" + +[[deps.LibCURL]] +deps = ["LibCURL_jll", "MozillaCACerts_jll"] +uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" +version = "0.6.3" + +[[deps.LibCURL_jll]] +deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"] +uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" +version = "7.84.0+0" [[deps.LibGit2]] deps = ["Base64", "NetworkOptions", "Printf", "SHA"] uuid = "76f85450-5226-5b5a-8eaa-529ad045b433" +[[deps.LibSSH2_jll]] +deps = ["Artifacts", "Libdl", "MbedTLS_jll"] +uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" +version = "1.10.2+0" + +[[deps.Libdl]] +uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" + +[[deps.LineSearches]] +deps = ["LinearAlgebra", "NLSolversBase", "NaNMath", "Parameters", "Printf"] +git-tree-sha1 = "7bbea35cec17305fc70a0e5b4641477dc0789d9d" +uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255" +version = "7.2.0" + +[[deps.LinearAlgebra]] +deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"] +uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" + +[[deps.LogExpFunctions]] +deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"] +git-tree-sha1 = "c3ce8e7420b3a6e071e0fe4745f5d4300e37b13f" +uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" +version = "0.3.24" + + [deps.LogExpFunctions.extensions] + LogExpFunctionsChainRulesCoreExt = "ChainRulesCore" + LogExpFunctionsChangesOfVariablesExt = "ChangesOfVariables" + LogExpFunctionsInverseFunctionsExt = "InverseFunctions" + + [deps.LogExpFunctions.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + ChangesOfVariables = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" + InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" + [[deps.Logging]] uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" +[[deps.MacroTools]] +deps = ["Markdown", "Random"] +git-tree-sha1 = "42324d08725e200c23d4dfb549e0d5d89dede2d2" +uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" +version = "0.5.10" + [[deps.Markdown]] deps = ["Base64"] uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" +[[deps.MbedTLS_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" +version = "2.28.2+0" + +[[deps.Missings]] +deps = ["DataAPI"] +git-tree-sha1 = "f66bdc5de519e8f8ae43bdc598782d35a25b1272" +uuid = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28" +version = "1.1.0" + [[deps.Mmap]] uuid = "a63ad114-7e13-5084-954f-fe012c677804" +[[deps.MozillaCACerts_jll]] +uuid = "14a3606d-f60d-562e-9121-12d972cd8159" +version = "2022.10.11" + +[[deps.NLSolversBase]] +deps = ["DiffResults", "Distributed", "FiniteDiff", "ForwardDiff"] +git-tree-sha1 = "a0b464d183da839699f4c79e7606d9d186ec172c" +uuid = "d41bc354-129a-5804-8e4c-c37616107c6c" +version = "7.8.3" + +[[deps.NaNMath]] +deps = ["OpenLibm_jll"] +git-tree-sha1 = "0877504529a3e5c3343c6f8b4c0381e57e4387e4" +uuid = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" +version = "1.0.2" + [[deps.NetworkOptions]] uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" +version = "1.2.0" + +[[deps.OpenBLAS_jll]] +deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] +uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" +version = "0.3.21+4" + +[[deps.OpenLibm_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "05823500-19ac-5b8b-9628-191a04bc5112" +version = "0.8.1+0" + +[[deps.OpenSpecFun_jll]] +deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1" +uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e" +version = "0.5.5+0" + +[[deps.Optim]] +deps = ["Compat", "FillArrays", "ForwardDiff", "LineSearches", "LinearAlgebra", "NLSolversBase", "NaNMath", "Parameters", "PositiveFactorizations", "Printf", "SparseArrays", "StatsBase"] +git-tree-sha1 = "e3a6546c1577bfd701771b477b794a52949e7594" +uuid = "429524aa-4258-5aef-a3af-852621145aeb" +version = "1.7.6" + +[[deps.OrderedCollections]] +git-tree-sha1 = "2e73fe17cac3c62ad1aebe70d44c963c3cfdc3e3" +uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" +version = "1.6.2" + +[[deps.Parameters]] +deps = ["OrderedCollections", "UnPack"] +git-tree-sha1 = "34c0e9ad262e5f7fc75b10a9952ca7692cfc5fbe" +uuid = "d96e819e-fc66-5662-9728-84c9c7592b0a" +version = "0.12.3" [[deps.Parsers]] -deps = ["Dates"] -git-tree-sha1 = "0044b23da09b5608b4ecacb4e5e6c6332f833a7e" +deps = ["Dates", "PrecompileTools", "UUIDs"] +git-tree-sha1 = "716e24b21538abc91f6205fd1d8363f39b442851" uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0" -version = "2.3.2" +version = "2.7.2" + +[[deps.Pkg]] +deps = ["Artifacts", "Dates", "Downloads", "FileWatching", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] +uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" +version = "1.9.2" + +[[deps.PositiveFactorizations]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "17275485f373e6673f7e7f97051f703ed5b15b20" +uuid = "85a6dd25-e78a-55b7-8502-1745935b8125" +version = "0.2.4" + +[[deps.PrecompileTools]] +deps = ["Preferences"] +git-tree-sha1 = "9673d39decc5feece56ef3940e5dafba15ba0f81" +uuid = "aea7be01-6a6a-4083-8856-8a6e6704d82a" +version = "1.1.2" + +[[deps.Preferences]] +deps = ["TOML"] +git-tree-sha1 = "7eb1686b4f04b82f96ed7a4ea5890a4f0c7a09f1" +uuid = "21216c6a-2e73-6563-6e65-726566657250" +version = "1.4.0" + +[[deps.PrettyTables]] +deps = ["Crayons", "LaTeXStrings", "Markdown", "Printf", "Reexport", "StringManipulation", "Tables"] +git-tree-sha1 = "ee094908d720185ddbdc58dbe0c1cbe35453ec7a" +uuid = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" +version = "2.2.7" [[deps.Printf]] deps = ["Unicode"] @@ -83,18 +411,157 @@ uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" deps = ["SHA", "Serialization"] uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +[[deps.RecipesBase]] +deps = ["PrecompileTools"] +git-tree-sha1 = "5c3d09cc4f31f5fc6af001c250bf1278733100ff" +uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" +version = "1.3.4" + +[[deps.Reexport]] +git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b" +uuid = "189a3867-3050-52da-a836-e630ba90ab69" +version = "1.2.2" + +[[deps.Requires]] +deps = ["UUIDs"] +git-tree-sha1 = "838a3a4188e2ded87a4f9f184b4b0d78a1e91cb7" +uuid = "ae029012-a4dd-5104-9daa-d747884805df" +version = "1.3.0" + [[deps.SHA]] uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" +version = "0.7.0" [[deps.Serialization]] uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b" +[[deps.Setfield]] +deps = ["ConstructionBase", "Future", "MacroTools", "StaticArraysCore"] +git-tree-sha1 = "e2cc6d8c88613c05e1defb55170bf5ff211fbeac" +uuid = "efcf1570-3423-57d1-acb7-fd33fddbac46" +version = "1.1.1" + [[deps.Sockets]] uuid = "6462fe0b-24de-5631-8697-dd941f90decc" +[[deps.SortingAlgorithms]] +deps = ["DataStructures"] +git-tree-sha1 = "c60ec5c62180f27efea3ba2908480f8055e17cee" +uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c" +version = "1.1.1" + +[[deps.SparseArrays]] +deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"] +uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" + +[[deps.SpecialFunctions]] +deps = ["IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] +git-tree-sha1 = "7beb031cf8145577fbccacd94b8a8f4ce78428d3" +uuid = "276daf66-3868-5448-9aa4-cd146d93841b" +version = "2.3.0" + + [deps.SpecialFunctions.extensions] + SpecialFunctionsChainRulesCoreExt = "ChainRulesCore" + + [deps.SpecialFunctions.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + +[[deps.StaticArraysCore]] +git-tree-sha1 = "36b3d696ce6366023a0ea192b4cd442268995a0d" +uuid = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" +version = "1.4.2" + +[[deps.Statistics]] +deps = ["LinearAlgebra", "SparseArrays"] +uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +version = "1.9.0" + +[[deps.StatsAPI]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "45a7769a04a3cf80da1c1c7c60caf932e6f4c9f7" +uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0" +version = "1.6.0" + +[[deps.StatsBase]] +deps = ["DataAPI", "DataStructures", "LinearAlgebra", "LogExpFunctions", "Missings", "Printf", "Random", "SortingAlgorithms", "SparseArrays", "Statistics", "StatsAPI"] +git-tree-sha1 = "75ebe04c5bed70b91614d684259b661c9e6274a4" +uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" +version = "0.34.0" + +[[deps.StringManipulation]] +git-tree-sha1 = "46da2434b41f41ac3594ee9816ce5541c6096123" +uuid = "892a3eda-7b42-436c-8928-eab12a02cf0e" +version = "0.3.0" + +[[deps.SuiteSparse]] +deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"] +uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" + +[[deps.SuiteSparse_jll]] +deps = ["Artifacts", "Libdl", "Pkg", "libblastrampoline_jll"] +uuid = "bea87d4a-7f5b-5778-9afe-8cc45184846c" +version = "5.10.1+6" + +[[deps.TOML]] +deps = ["Dates"] +uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" +version = "1.0.3" + +[[deps.TableTraits]] +deps = ["IteratorInterfaceExtensions"] +git-tree-sha1 = "c06b2f539df1c6efa794486abfb6ed2022561a39" +uuid = "3783bdb8-4a98-5b6b-af9a-565f29a5fe9c" +version = "1.0.1" + +[[deps.Tables]] +deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "OrderedCollections", "TableTraits", "Test"] +git-tree-sha1 = "1544b926975372da01227b382066ab70e574a3ec" +uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c" +version = "1.10.1" + +[[deps.Tar]] +deps = ["ArgTools", "SHA"] +uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" +version = "1.10.0" + [[deps.Test]] deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +[[deps.TimeSeries]] +deps = ["Dates", "DelimitedFiles", "DocStringExtensions", "RecipesBase", "Reexport", "Statistics", "Tables"] +git-tree-sha1 = "8b9288d84da88ea44693ca8cf9c236da1778f274" +uuid = "9e3dc215-6440-5c97-bce1-76c03772f85e" +version = "0.23.2" + +[[deps.UUIDs]] +deps = ["Random", "SHA"] +uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" + +[[deps.UnPack]] +git-tree-sha1 = "387c1f73762231e86e0c9c5443ce3b4a0a9a0c2b" +uuid = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" +version = "1.0.2" + [[deps.Unicode]] uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" + +[[deps.Zlib_jll]] +deps = ["Libdl"] +uuid = "83775a58-1f1d-513f-b197-d71354ab007a" +version = "1.2.13+0" + +[[deps.libblastrampoline_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" +version = "5.8.0+0" + +[[deps.nghttp2_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" +version = "1.48.0+0" + +[[deps.p7zip_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0" +version = "17.4.0+0" diff --git a/src/datacleaning.jl b/src/datacleaning.jl new file mode 100644 index 0000000..8abd7bc --- /dev/null +++ b/src/datacleaning.jl @@ -0,0 +1,93 @@ +using Dates, Loess, Interpolations, Missings, TimeSeries + +# These are used only during testing +using Plots, Distributions + +""" + calcHoL(timestamps) + +Calculates the hours on-line for a series of timestamps, assuming the first timestamp = 0 hours on-line. +""" +function calcHoL(timestamps) + starttime = timestamps[begin] + endtime = timestamps[end] + + return (timestamps .- starttime) ./ Hour(1) +end + + +# Generate dummy data for testing +function gendata(timestamps, period, fracfilled) + HoL = calcHoL(timestamps) + data = zeros(Union{Float64, Missing}, length(times)) + + for i in eachindex(HoL) + if rand() > fracfilled + data[i] = missing + else + data[i] = sin(π*period*(HoL[i]/(Hour(endtime - starttime)/Hour(1)))) + norm = Normal(0, 0.1*abs(data[i])) + data[i] += rand(norm) + end + end + + return data +end + +""" + smooth_and_fill(raw; loessspan=0.3, suggest_start=false, startval=0.0, suggest_end=false, endval=0.0) + +Smooth and fill a time series using loess, with suggested start and end values or linear extrapolations. +""" +function smooth_and_fill(raw; loessspan=0.3, suggest_start=false, startvals=Float64[], suggest_end=false, endvals=Float64[]) + HoL = calcHoL(timestamp(raw)) + data = similar(values(raw)) + + for (i, col) in enumerate(colnames(raw)) + # Drop missing data entries + _HoL, _data = collect.(skipmissings(HoL, values(raw[col]))) + + # Prepare a linear interpolator that will extrapolate to start and end values + extrap = linear_interpolation(_HoL, _data, extrapolation_bc = Line()) + + # Add start and end values, if missing + if !any(x -> x == HoL[begin], _HoL) + if suggest_start + pushfirst!(_HoL, HoL[begin]) + pushfirst!(_data, startvals[i]) + else + pushfirst!(_HoL, HoL[begin]) + pushfirst!(_data, extrap(HoL[begin])) + end + end + if !any(x -> x == HoL[end], _HoL) + if suggest_end + push!(_HoL, HoL[end]) + push!(_data, endvals[i]) + else + push!(_HoL, HoL[end]) + push!(_data, extrap(HoL[end])) + end + end + + # Fit a loess model to the data + model = loess(_HoL, _data; span=loessspan, degree=2) + data[:, i] .= predict(model, HoL) + end + + return TimeArray(timestamp(raw), data, colnames(raw)) +end + + +# Generate dummy data with missing values +starttime = DateTime(2023, 1, 1, 0, 0) +endtime = DateTime(2023, 1, 12, 24, 0) +times = starttime:Hour(6):endtime + +raw = TimeArray(times, [gendata(times, 2, 0.95) gendata(times, 2, 0.95)], [:raw1, :raw2]) +plot(raw, leg=:bottomleft) + +sf1 = smooth_and_fill(raw, suggest_start=true, suggest_end=true, startvals=[0.0, 0.0], endvals=[0.0, 0.0]) +sf2 = smooth_and_fill(raw) +plot!(sf1) +plot!(sf2) diff --git a/test/Manifest.toml b/test/Manifest.toml index ab6c466..974a54d 100644 --- a/test/Manifest.toml +++ b/test/Manifest.toml @@ -1,6 +1,6 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.8.2" +julia_version = "1.9.2" manifest_format = "2.0" project_hash = "71d91126b5a1fb1020e1098d9d492de2a4438fd2" From 20cc5741e8100e1e0e9bdf02bb2b5e54ab53b47c Mon Sep 17 00:00:00 2001 From: Braam van Dyk Date: Fri, 1 Dec 2023 20:20:28 +0200 Subject: [PATCH 2/2] Complete rebuild using TimeArrays for streams --- Manifest.toml | 181 ++++++--- Project.toml | 5 +- components/nitrogen.comp | 5 + demo2.jl | 46 +++ demos.jl | 297 +++++++------- src/FlowsheetTools.jl | 40 +- src/boundaries.jl | 257 +++--------- src/closure.jl | 162 +++----- src/components.jl | 39 +- src/datacleaning.jl | 58 +-- src/flowsheets.jl | 43 +- src/kpis.jl | 143 +++---- src/stepchange.jl | 2 + src/streams.jl | 594 ++++++++++------------------ src/unitops.jl | 274 ++++--------- streamhistories/C2.csv | 28 ++ streamhistories/FeedStreamShort.csv | 2 + streamhistories/Hydrogen.csv | 28 ++ streamhistories/Product.csv | 28 ++ streamhistories/Wrongstream.csv | 28 ++ 20 files changed, 1002 insertions(+), 1258 deletions(-) create mode 100644 components/nitrogen.comp create mode 100644 demo2.jl create mode 100644 streamhistories/C2.csv create mode 100644 streamhistories/FeedStreamShort.csv create mode 100644 streamhistories/Hydrogen.csv create mode 100644 streamhistories/Product.csv create mode 100644 streamhistories/Wrongstream.csv diff --git a/Manifest.toml b/Manifest.toml index 9bb412c..2629735 100644 --- a/Manifest.toml +++ b/Manifest.toml @@ -1,14 +1,14 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.9.3" +julia_version = "1.9.4" manifest_format = "2.0" -project_hash = "386e3531f9b1b93f33ab1ccde1ce009c6c3d47bb" +project_hash = "8960c8cb9a90acf90744e7273e593bbb51356eff" [[deps.Adapt]] deps = ["LinearAlgebra", "Requires"] -git-tree-sha1 = "76289dc51920fdc6e0013c872ba9551d54961c24" +git-tree-sha1 = "02f731463748db57cc2ebfbd9fbc9ce8280d3433" uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" -version = "3.6.2" +version = "3.7.1" weakdeps = ["StaticArrays"] [deps.Adapt.extensions] @@ -20,9 +20,9 @@ version = "1.1.1" [[deps.ArrayInterface]] deps = ["Adapt", "LinearAlgebra", "Requires", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "f83ec24f76d4c8f525099b2ac475fc098138ec31" +git-tree-sha1 = "247efbccf92448be332d154d6ca56b9fcdd93c31" uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" -version = "7.4.11" +version = "7.6.1" [deps.ArrayInterface.extensions] ArrayInterfaceBandedMatricesExt = "BandedMatrices" @@ -68,6 +68,22 @@ weakdeps = ["SparseArrays"] [deps.ChainRulesCore.extensions] ChainRulesCoreSparseArraysExt = "SparseArrays" +[[deps.ChangePointDetection]] +deps = ["LinearAlgebra", "Random"] +git-tree-sha1 = "b7e6b08f71aaf04ad3e2d97ca0e5674abc99d05f" +uuid = "68f453bf-a161-50e8-aafe-8ad75c8d6944" +version = "1.2.0" + +[[deps.Combinatorics]] +git-tree-sha1 = "08c8b6831dc00bfea825826be0bc8336fc369860" +uuid = "861a8166-3701-5b0c-9a16-15d98fcdc6aa" +version = "1.0.2" + +[[deps.CommonSolve]] +git-tree-sha1 = "0eee5eb66b1cf62cd6ad1b460238e60e4b09400c" +uuid = "38540f10-b2f7-11e9-35d8-d573e4eb0ff2" +version = "0.2.4" + [[deps.CommonSubexpressions]] deps = ["MacroTools", "Test"] git-tree-sha1 = "7b8a93dba8af7e3b42fecabf646260105ac373f7" @@ -76,9 +92,9 @@ version = "0.3.0" [[deps.Compat]] deps = ["UUIDs"] -git-tree-sha1 = "e460f044ca8b99be31d35fe54fc33a5c33dd8ed7" +git-tree-sha1 = "886826d76ea9e72b35fcd000e535588f7b60f21d" uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" -version = "4.9.0" +version = "4.10.1" weakdeps = ["Dates", "LinearAlgebra"] [deps.Compat.extensions] @@ -91,9 +107,9 @@ version = "1.0.5+0" [[deps.ConstructionBase]] deps = ["LinearAlgebra"] -git-tree-sha1 = "fe2838a593b5f776e1597e086dcd47560d94e816" +git-tree-sha1 = "c53fc348ca4d40d7b371e71fd52251839080cbc9" uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" -version = "1.5.3" +version = "1.5.4" [deps.ConstructionBase.extensions] ConstructionBaseIntervalSetsExt = "IntervalSets" @@ -148,9 +164,9 @@ version = "1.15.1" [[deps.Distances]] deps = ["LinearAlgebra", "Statistics", "StatsAPI"] -git-tree-sha1 = "5225c965635d8c21168e32a12954675e7bea1151" +git-tree-sha1 = "66c4c81f259586e8f002eacebc177e1fb06363b0" uuid = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7" -version = "0.10.10" +version = "0.10.11" weakdeps = ["ChainRulesCore", "SparseArrays"] [deps.Distances.extensions] @@ -198,10 +214,16 @@ version = "0.6.8" uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" [[deps.FillArrays]] -deps = ["LinearAlgebra", "Random", "SparseArrays", "Statistics"] -git-tree-sha1 = "f372472e8672b1d993e93dada09e23139b509f9e" +deps = ["LinearAlgebra", "Random"] +git-tree-sha1 = "28e4e9c4b7b162398ec8004bdabe9a90c78c122d" uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" -version = "1.5.0" +version = "1.8.0" +weakdeps = ["PDMats", "SparseArrays", "Statistics"] + + [deps.FillArrays.extensions] + FillArraysPDMatsExt = "PDMats" + FillArraysSparseArraysExt = "SparseArrays" + FillArraysStatisticsExt = "Statistics" [[deps.FiniteDiff]] deps = ["ArrayInterface", "LinearAlgebra", "Requires", "Setfield", "SparseArrays"] @@ -239,15 +261,27 @@ git-tree-sha1 = "f218fe3736ddf977e0e772bc9a586b2383da2685" uuid = "34004b35-14d8-5ef3-9330-4cdb6864b03a" version = "0.3.23" +[[deps.HypothesisTests]] +deps = ["Combinatorics", "Distributions", "LinearAlgebra", "Printf", "Random", "Rmath", "Roots", "Statistics", "StatsAPI", "StatsBase"] +git-tree-sha1 = "4b5d5ba51f5f473737ed9de6d8a7aa190ad8c72f" +uuid = "09f84164-cd44-5f33-b23f-e6b0d136a0d5" +version = "0.11.0" + [[deps.InteractiveUtils]] deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" [[deps.Interpolations]] deps = ["Adapt", "AxisAlgorithms", "ChainRulesCore", "LinearAlgebra", "OffsetArrays", "Random", "Ratios", "Requires", "SharedArrays", "SparseArrays", "StaticArrays", "WoodburyMatrices"] -git-tree-sha1 = "721ec2cf720536ad005cb38f50dbba7b02419a15" +git-tree-sha1 = "274ad8005db8b4ef6cc46d1392927083405813c2" uuid = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" -version = "0.14.7" +version = "0.15.0" + + [deps.Interpolations.extensions] + InterpolationsUnitfulExt = "Unitful" + + [deps.Interpolations.weakdeps] + Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" [[deps.IrrationalConstants]] git-tree-sha1 = "630b497eafcc20001bba38a4651b327dcfc491d2" @@ -260,25 +294,25 @@ uuid = "82899510-4779-5014-852e-03e436cf321d" version = "1.0.0" [[deps.JLLWrappers]] -deps = ["Preferences"] -git-tree-sha1 = "abc9885a7ca2052a736a600f7fa66209f96506e1" +deps = ["Artifacts", "Preferences"] +git-tree-sha1 = "7e5d6779a1e09a36db2a7b6cff50942a0a7d0fca" uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210" -version = "1.4.1" +version = "1.5.0" [[deps.LaTeXStrings]] -git-tree-sha1 = "f2355693d6778a178ade15952b7ac47a4ff97996" +git-tree-sha1 = "50901ebc375ed41dbf8058da26f9de442febbbec" uuid = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f" -version = "1.3.0" +version = "1.3.1" [[deps.LibCURL]] deps = ["LibCURL_jll", "MozillaCACerts_jll"] uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" -version = "0.6.3" +version = "0.6.4" [[deps.LibCURL_jll]] deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"] uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" -version = "7.84.0+0" +version = "8.4.0+0" [[deps.LibGit2]] deps = ["Base64", "NetworkOptions", "Printf", "SHA"] @@ -287,7 +321,7 @@ uuid = "76f85450-5226-5b5a-8eaa-529ad045b433" [[deps.LibSSH2_jll]] deps = ["Artifacts", "Libdl", "MbedTLS_jll"] uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" -version = "1.10.2+0" +version = "1.11.0+1" [[deps.Libdl]] uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" @@ -310,9 +344,9 @@ version = "0.6.3" [[deps.LogExpFunctions]] deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"] -git-tree-sha1 = "c3ce8e7420b3a6e071e0fe4745f5d4300e37b13f" +git-tree-sha1 = "7d6dd4e9212aebaeed356de34ccf262a3cd415aa" uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" -version = "0.3.24" +version = "0.3.26" [deps.LogExpFunctions.extensions] LogExpFunctionsChainRulesCoreExt = "ChainRulesCore" @@ -329,9 +363,9 @@ uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" [[deps.MacroTools]] deps = ["Markdown", "Random"] -git-tree-sha1 = "42324d08725e200c23d4dfb549e0d5d89dede2d2" +git-tree-sha1 = "9ee1618cbf5240e6d4e0371d6f24065083f60c48" uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" -version = "0.5.10" +version = "0.5.11" [[deps.Markdown]] deps = ["Base64"] @@ -395,20 +429,20 @@ version = "0.5.5+0" [[deps.Optim]] deps = ["Compat", "FillArrays", "ForwardDiff", "LineSearches", "LinearAlgebra", "NLSolversBase", "NaNMath", "Parameters", "PositiveFactorizations", "Printf", "SparseArrays", "StatsBase"] -git-tree-sha1 = "e3a6546c1577bfd701771b477b794a52949e7594" +git-tree-sha1 = "01f85d9269b13fedc61e63cc72ee2213565f7a72" uuid = "429524aa-4258-5aef-a3af-852621145aeb" -version = "1.7.6" +version = "1.7.8" [[deps.OrderedCollections]] -git-tree-sha1 = "2e73fe17cac3c62ad1aebe70d44c963c3cfdc3e3" +git-tree-sha1 = "dfdf5519f235516220579f949664f1bf44e741c5" uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" -version = "1.6.2" +version = "1.6.3" [[deps.PDMats]] deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "f6f85a2edb9c356b829934ad3caed2ad0ebbfc99" +git-tree-sha1 = "4e5be6bb265d33669f98eb55d2a57addd1eeb72c" uuid = "90014a1f-27ba-587c-ab20-58faa44d9150" -version = "0.11.29" +version = "0.11.30" [[deps.Parameters]] deps = ["OrderedCollections", "UnPack"] @@ -429,21 +463,21 @@ version = "0.2.4" [[deps.PrecompileTools]] deps = ["Preferences"] -git-tree-sha1 = "9673d39decc5feece56ef3940e5dafba15ba0f81" +git-tree-sha1 = "03b4c25b43cb84cee5c90aa9b5ea0a78fd848d2f" uuid = "aea7be01-6a6a-4083-8856-8a6e6704d82a" -version = "1.1.2" +version = "1.2.0" [[deps.Preferences]] deps = ["TOML"] -git-tree-sha1 = "7eb1686b4f04b82f96ed7a4ea5890a4f0c7a09f1" +git-tree-sha1 = "00805cd429dcb4870060ff49ef443486c262e38e" uuid = "21216c6a-2e73-6563-6e65-726566657250" -version = "1.4.0" +version = "1.4.1" [[deps.PrettyTables]] -deps = ["Crayons", "LaTeXStrings", "Markdown", "Printf", "Reexport", "StringManipulation", "Tables"] -git-tree-sha1 = "ee094908d720185ddbdc58dbe0c1cbe35453ec7a" +deps = ["Crayons", "LaTeXStrings", "Markdown", "PrecompileTools", "Printf", "Reexport", "StringManipulation", "Tables"] +git-tree-sha1 = "88b895d13d53b5577fd53379d913b9ab9ac82660" uuid = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" -version = "2.2.7" +version = "2.3.1" [[deps.Printf]] deps = ["Unicode"] @@ -504,6 +538,24 @@ git-tree-sha1 = "6ed52fdd3382cf21947b15e8870ac0ddbff736da" uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f" version = "0.4.0+0" +[[deps.Roots]] +deps = ["ChainRulesCore", "CommonSolve", "Printf", "Setfield"] +git-tree-sha1 = "0f1d92463a020321983d04c110f476c274bafe2e" +uuid = "f2b01f46-fcfa-551c-844a-d8ac1e96c665" +version = "2.0.22" + + [deps.Roots.extensions] + RootsForwardDiffExt = "ForwardDiff" + RootsIntervalRootFindingExt = "IntervalRootFinding" + RootsSymPyExt = "SymPy" + RootsSymPyPythonCallExt = "SymPyPythonCall" + + [deps.Roots.weakdeps] + ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" + IntervalRootFinding = "d2bf35a9-74e0-55ec-b149-d360ff49b807" + SymPy = "24249f21-da20-56a4-8eb1-6a02cf4ae2e6" + SymPyPythonCall = "bc8888f7-b21e-4b7c-a06a-5d9c9496438c" + [[deps.SHA]] uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" version = "0.7.0" @@ -526,9 +578,9 @@ uuid = "6462fe0b-24de-5631-8697-dd941f90decc" [[deps.SortingAlgorithms]] deps = ["DataStructures"] -git-tree-sha1 = "c60ec5c62180f27efea3ba2908480f8055e17cee" +git-tree-sha1 = "5165dfb9fd131cf0c6957a3a7605dede376e7b63" uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c" -version = "1.1.1" +version = "1.2.0" [[deps.SparseArrays]] deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"] @@ -536,19 +588,19 @@ uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [[deps.SpecialFunctions]] deps = ["IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] -git-tree-sha1 = "7beb031cf8145577fbccacd94b8a8f4ce78428d3" +git-tree-sha1 = "e2cfc4012a19088254b3950b85c3c1d8882d864d" uuid = "276daf66-3868-5448-9aa4-cd146d93841b" -version = "2.3.0" +version = "2.3.1" weakdeps = ["ChainRulesCore"] [deps.SpecialFunctions.extensions] SpecialFunctionsChainRulesCoreExt = "ChainRulesCore" [[deps.StaticArrays]] -deps = ["LinearAlgebra", "Random", "StaticArraysCore"] -git-tree-sha1 = "0adf069a2a490c47273727e029371b31d44b72b2" +deps = ["LinearAlgebra", "PrecompileTools", "Random", "StaticArraysCore"] +git-tree-sha1 = "5ef59aea6f18c25168842bded46b16662141ab87" uuid = "90137ffa-7385-5640-81b9-e52037218182" -version = "1.6.5" +version = "1.7.0" weakdeps = ["Statistics"] [deps.StaticArrays.extensions] @@ -566,15 +618,15 @@ version = "1.9.0" [[deps.StatsAPI]] deps = ["LinearAlgebra"] -git-tree-sha1 = "45a7769a04a3cf80da1c1c7c60caf932e6f4c9f7" +git-tree-sha1 = "1ff449ad350c9c4cbc756624d6f8a8c3ef56d3ed" uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0" -version = "1.6.0" +version = "1.7.0" [[deps.StatsBase]] deps = ["DataAPI", "DataStructures", "LinearAlgebra", "LogExpFunctions", "Missings", "Printf", "Random", "SortingAlgorithms", "SparseArrays", "Statistics", "StatsAPI"] -git-tree-sha1 = "75ebe04c5bed70b91614d684259b661c9e6274a4" +git-tree-sha1 = "1d77abd07f617c4868c33d4f5b9e1dbb2643c9cf" uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" -version = "0.34.0" +version = "0.34.2" [[deps.StatsFuns]] deps = ["HypergeometricFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] @@ -591,9 +643,10 @@ version = "1.3.0" InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" [[deps.StringManipulation]] -git-tree-sha1 = "46da2434b41f41ac3594ee9816ce5541c6096123" +deps = ["PrecompileTools"] +git-tree-sha1 = "a04cabe79c5f01f4d723cc6704070ada0b9d46d5" uuid = "892a3eda-7b42-436c-8928-eab12a02cf0e" -version = "0.3.0" +version = "0.3.4" [[deps.SuiteSparse]] deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"] @@ -616,10 +669,10 @@ uuid = "3783bdb8-4a98-5b6b-af9a-565f29a5fe9c" version = "1.0.1" [[deps.Tables]] -deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "OrderedCollections", "TableTraits", "Test"] -git-tree-sha1 = "1544b926975372da01227b382066ab70e574a3ec" +deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "OrderedCollections", "TableTraits"] +git-tree-sha1 = "cb76cf677714c095e535e3501ac7954732aeea2d" uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c" -version = "1.10.1" +version = "1.11.1" [[deps.Tar]] deps = ["ArgTools", "SHA"] @@ -631,10 +684,10 @@ deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [[deps.TimeSeries]] -deps = ["Dates", "DelimitedFiles", "DocStringExtensions", "RecipesBase", "Reexport", "Statistics", "Tables"] -git-tree-sha1 = "8b9288d84da88ea44693ca8cf9c236da1778f274" +deps = ["Dates", "DelimitedFiles", "DocStringExtensions", "PrettyTables", "RecipesBase", "Reexport", "Statistics", "Tables"] +git-tree-sha1 = "665545ab6c4d651d6ccbb00b08427a7b6eeb3da4" uuid = "9e3dc215-6440-5c97-bce1-76c03772f85e" -version = "0.23.2" +version = "0.24.0" [[deps.UUIDs]] deps = ["Random", "SHA"] @@ -650,9 +703,9 @@ uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" [[deps.WoodburyMatrices]] deps = ["LinearAlgebra", "SparseArrays"] -git-tree-sha1 = "de67fa59e33ad156a590055375a30b23c40299d3" +git-tree-sha1 = "5f24e158cf4cee437052371455fe361f526da062" uuid = "efce3f68-66dc-5838-9240-27a6d6f5f9b6" -version = "0.5.5" +version = "0.5.6" [[deps.Zlib_jll]] deps = ["Libdl"] @@ -667,7 +720,7 @@ version = "5.8.0+0" [[deps.nghttp2_jll]] deps = ["Artifacts", "Libdl"] uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" -version = "1.48.0+0" +version = "1.52.0+1" [[deps.p7zip_jll]] deps = ["Artifacts", "Libdl"] diff --git a/Project.toml b/Project.toml index 60868f7..2cc272c 100644 --- a/Project.toml +++ b/Project.toml @@ -5,21 +5,22 @@ version = "0.1.0" [deps] Base64 = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" +ChangePointDetection = "68f453bf-a161-50e8-aafe-8ad75c8d6944" Dates = "ade2ca70-3891-5945-98fb-dc099432e06a" DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab" -Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" Downloads = "f43a241f-c20a-4ad4-852c-f6b1247861c6" +HypothesisTests = "09f84164-cd44-5f33-b23f-e6b0d136a0d5" Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" Loess = "4345ca2d-374a-55d4-8d30-97f9976e7612" MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" Missings = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28" Optim = "429524aa-4258-5aef-a3af-852621145aeb" +OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" PrettyTables = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" TimeSeries = "9e3dc215-6440-5c97-bce1-76c03772f85e" [compat] -DelimitedFiles = "1" Downloads = "1" MacroTools = "0.5" Optim = "1.7" diff --git a/components/nitrogen.comp b/components/nitrogen.comp new file mode 100644 index 0000000..36596ee --- /dev/null +++ b/components/nitrogen.comp @@ -0,0 +1,5 @@ +Nitrogen +1 +N +2 +28.0134 \ No newline at end of file diff --git a/demo2.jl b/demo2.jl new file mode 100644 index 0000000..6fc47d1 --- /dev/null +++ b/demo2.jl @@ -0,0 +1,46 @@ +using FlowsheetTools + +syscomps = ComponentList() +count = readcomponentlist!(syscomps, "components", ["Ethylene", "Ethane", "Hydrogen"]) + + + + +# The next one should error +# sysstreams["Wrong"] = readstreamhistory(joinpath("streamhistories", "WrongStream.csv"), "Wrong", syscomps; ismoleflow=true) + +dummy = emptystream(sysstreams, "Dummy") +dummy = sysstreams["Feed"] +dummies = sysstreams[["Feed", "Product"]] + +count = length(sysstreams) + +dummy = dummies[1] + dummies[2] +dummy = 2 * dummy +dummy = dummy * 2 + +shortstreams = StreamList() +shortstreams["Feed"] = readstreamhistory(joinpath("streamhistories", "FeedStreamShort.csv"), "Feed", syscomps; ismoleflow=true) + +# Iterate through the CompList - returns pairs of names and Component objects +for (compname, comp) in syscomps + println(compname) +end + +# Iterate through the StreamList - returns pairs of names and Stream objects +for (streamname, stream) in sysstreams + println(streamname) +end + +sysunits = UnitOpList() + +@unitop begin + inlets --> ["Feed"] + outlets --> ["Product"] + calc --> mixer! +end "Mixer" sysstreams sysunits + +# Iterate through the UnitOpList - returns pairs of names and UnitOp objects +for (unitopname, unitop) in sysunits + println(unitopname) +end diff --git a/demos.jl b/demos.jl index edc4dd5..5bb7dae 100644 --- a/demos.jl +++ b/demos.jl @@ -1,34 +1,41 @@ -using FlowsheetTools - +using FlowsheetTools, Statistics #-Components--------------------------- +# We need a ComponentList to hold all the +# components, so we know where to find them later syscomps = ComponentList() -@comp begin - H --> 2 -end "Hydrogen" syscomps -@comp begin - C --> 2 - H --> 6 -end "Ethane" syscomps +# You can read them from a folder with saved components +count = readcomponentlist!(syscomps, "components", ["Ethylene", "Ethane", "Hydrogen"]) + +# Or define them directly @comp begin - C --> 2 - H --> 4 -end "Ethylene" syscomps + N --> 2 +end "Nitrogen" syscomps + +# And then save them to file + +writecomponent(joinpath("components/", "Nitrogen.comp"), syscomps["Nitrogen"]) + -writecomponent(joinpath("components/", "Hydrogen.comp"), syscomps["Hydrogen"]) -writecomponent(joinpath("components/", "Ethane.comp"), syscomps["Ethane"]) -writecomponent(joinpath("components/", "Ethylene.comp"), syscomps["Ethylene"]) #-Streams------------------------------ +# Like for the components, we need a container for streams +# so we have something to iterate through later sysstreams = StreamList() +# You can create the streams directly with instantaneous flows +# This can be in either mass or molar unitops +# The units are not specified - if you assume the mass flows are +# in kg/h, then the molar equivalent is kmol/hr, but this could +# as easily be lb/week and lbmole/week. + @stream mass begin "Ethylene" --> 2.8053 "Ethane" --> 27.06192 @@ -36,55 +43,84 @@ sysstreams = StreamList() end "Test" syscomps sysstreams @stream mole begin - "Ethylene" --> 0.1 "Ethane" --> 0.9 "Hydrogen" --> 1.1 + "Ethylene" --> 0.1 end "Product" syscomps sysstreams -sysstreams["Test"].moleflows ≈ sysstreams["Product"].moleflows +# One stream was specified as mass flows, the other as mass flows, +# but these streams are the same and the missing flows are calculated +# automatically + +sysstreams["Test"].moleflows .≈ sysstreams["Product"].moleflows + +# As are the atomic flows, in the same units as the molar flows + +all(getindex.(values(sysstreams["Test"].atomflows), "C") .== getindex.(values(sysstreams["Product"].atomflows), "C")) +all(getindex.(values(sysstreams["Test"].atomflows), "H") .== getindex.(values(sysstreams["Product"].atomflows), "H")) + +# When we want to deal with streams with multiple historic data points, +# the best option is to use readstreamhistory to read them from file + +sysstreams = StreamList() # Create a new container and dump the previous streams +sysstreams["Feed"] = readstreamhistory(joinpath("streamhistories", "FeedStream.csv"), "Feed", syscomps; ismoleflow=true) +sysstreams["Product"] = readstreamhistory(joinpath("streamhistories", "ProdStream.csv"), "Product", syscomps; ismoleflow=true) + +# In the files, we had data for Ethylene, Ethane and Hydrogen, but our list of +# components also includes nitrogen. Zero flows are automatically added for any +# components not in the file, so all streams contain all the components. +# We can also add components after reading the files, but then we need to +# call refreshcomplist(streamlist) to add the zero flows to all the new components +# to existing streams. + +@comp begin + Ar --> 1 +end "Argon" syscomps + +refreshcomplist(sysstreams) + +sysstreams["Feed"] + + +# ================================ + # Manipulate some streams +# Multiplication with a scalar +sysstreams["Prod2"] = 2.0*sysstreams["Product"] + +# Note the use of .≈ and .* - internally the data are stored in TimeArrays (TimeSeries.jl) +# and only the braodcasted operators are used on TimeArrays. Comparing to TimeArrays also +# returns a TimeArray and we extract the values using values() to get a BitVector +all(values(sysstreams["Prod2"].totalmassflow) .≈ values(2.0 .* sysstreams["Product"].totalmassflow)) + # Copy and copy with multiply -copystream!(sysstreams, "Product", "mystream") -copystream!(sysstreams, "Product", "mystream2"; factor=2.0) # double the flow! +copystream!(sysstreams, "Product", "MyStream") +copystream!(sysstreams, "Product", "MyStream2"; factor=2.0) # double the flow! -sysstreams["mystream2"].totalmassflow ≈ 2.0 * sysstreams["mystream"].totalmassflow +# Check that they are identical +all(values(sysstreams["MyStream2"].totalmassflow .≈ 2.0 .* sysstreams["MyStream"].totalmassflow)) # Different name, so not identical -sysstreams["Product"] == sysstreams["mystream"] +sysstreams["Product"] == sysstreams["MyStream"] # But all the flows are the same -sysstreams["Product"].atomflows == sysstreams["mystream"].atomflows +all(getindex.(values(sysstreams["Product"].atomflows), "C") .== getindex.(values(sysstreams["MyStream"].atomflows), "C")) +all(getindex.(values(sysstreams["Product"].atomflows), "H") .== getindex.(values(sysstreams["MyStream"].atomflows), "H")) +all(getindex.(values(sysstreams["Product"].atomflows), "N") .== getindex.(values(sysstreams["MyStream"].atomflows), "N")) # Do some more things with streams -renamestream!(sysstreams, "mystream", "dummy") -deletestream!(sysstreams, "dummy") - -# Multiplication with a scalar -sysstreams["Prod2"] = 2.0*sysstreams["Product"] -sysstreams["Prod2"].totalmassflow ≈ 2.0 * sysstreams["Product"].totalmassflow +renamestream!(sysstreams, "MyStream", "Dummy") +deletestream!(sysstreams, "Dummy") #-UnitOps, Boundaries and KPIs--------- - # Test unit ops and mass balance boundaries # Define some streams sysstreams = StreamList() -@stream mole begin - "Ethylene" --> 1.0 - "Hydrogen" --> 2.0 -end "Feed" syscomps sysstreams - -@stream mole begin - "Ethylene" --> 0.1 - "Ethane" --> 0.9 - "Hydrogen" --> 1.1 -end "Product" syscomps sysstreams - - @stream mole begin "Hydrogen" --> 1.1 end "H2" syscomps sysstreams @@ -94,23 +130,18 @@ end "H2" syscomps sysstreams "Ethane" --> 0.9 end "C2" syscomps sysstreams +sysstreams["Mixed"] = emptystream(sysstreams, "Mixed") + @stream mole begin - "Hydrogen" --> 0.0 -end "Mixed" syscomps sysstreams + "Ethylene" --> 0.0 + "Ethane" --> 1.0 + "Hydrogen" --> 1.0 +end "Product" syscomps sysstreams + # Define some unit ops sysunitops = UnitOpList() -@unitop begin - inlets --> ["Feed"] - outlets --> ["Product"] -end "Reactor" sysstreams sysunitops - -@unitop begin - inlets --> ["Product"] - outlets --> ["C2", "H2"] -end "Membrane" sysstreams sysunitops - @unitop begin inlets --> ["H2", "C2"] outlets --> ["Mixed"] @@ -118,9 +149,14 @@ end "Membrane" sysstreams sysunitops end "Mixer" sysstreams sysunitops sysunitops["Mixer"]() +@unitop begin + inlets --> ["Mixed"] + outlets --> ["Product"] +end "Reactor" sysstreams sysunitops + # Define a mass balance boundary @boundary begin - unitops --> ["Reactor", "Membrane"] + unitops --> ["Mixer", "Reactor"] end b sysunitops # Check the closures @@ -132,124 +168,95 @@ b.total_out.totalmassflow # And some defined KPIs conversion(b, "Ethane") conversion(b, "Ethylene") -selectivity(b, "Ethylene", "Ethane") +conversion(b, "Hydrogen") +molar_selectivity(b, "Ethylene", "Ethane") -# Now test the reconciliations -# Change some flows to simulate measurement errors -copystream!(sysstreams, "Feed", "Feed2"; factor=0.95) -copystream!(sysstreams, "Product", "Prod2"; factor=1.01) +# Let's do the same with some history attached to the streams + +sysstreams = StreamList() # Create a new container and dump the previous streams +sysstreams["C2"] = readstreamhistory(joinpath("streamhistories", "C2.csv"), "C2", syscomps; ismoleflow=true) +sysstreams["H2"] = readstreamhistory(joinpath("streamhistories", "Hydrogen.csv"), "H2", syscomps; ismoleflow=true) +sysstreams["Product"] = readstreamhistory(joinpath("streamhistories", "Product.csv"), "Product", syscomps; ismoleflow=true) +sysstreams["Mixed"] = emptystream(sysstreams, "Mixed") + +# Define some unit ops +sysunitops = UnitOpList() -# Define the unit ops and boundary @unitop begin - inlets --> ["Feed2"] - outlets --> ["Prod2"] -end "Reactor2" sysstreams sysunitops + inlets --> ["H2", "C2"] + outlets --> ["Mixed"] + calc --> mixer! +end "Mixer" sysstreams sysunitops +sysunitops["Mixer"]() @unitop begin - inlets --> ["Prod2"] - outlets --> ["C2", "H2"] -end "Membrane2" sysstreams sysunitops + inlets --> ["Mixed"] + outlets --> ["Product"] +end "Reactor" sysstreams sysunitops +# Define a mass balance boundary @boundary begin - unitops --> ["Reactor2", "Membrane2"] -end b2 sysunitops - -# Get the correction factors on the inlets and outlets -corrections = calccorrections(b2) -b2 = closemb(b2) -b2 - -# Test the pretty printing -print(syscomps["Ethane"]) -print(sysstreams["Feed"]) -print(sysunitops["Membrane"]) -print(b) - -# Read in a list of components -syscomps2 = ComponentList() -count = readcomponentlist!(syscomps2, "components", ["Ethylene", "Ethane", "Hydrogen"]) -count == 3 - + unitops --> ["Mixer", "Reactor"] +end b sysunitops -#-Streams with historical data--------- +# Check the closures +b.atomclosures +b.closure +b.total_in.totalmassflow +b.total_out.totalmassflow +# And some defined KPIs +c1 = conversion(b, "Ethane") +c2 = conversion(b, "Ethylene") +sc2 = molar_selectivity(b, "Ethylene", "Ethane") -histstreams = StreamHistoryList() -# Read in some stream histories -histstreams["Feed"] = readstreamhistory(joinpath("streamhistories", "FeedStream.csv"), "Feed", syscomps; ismoleflow=true) -histstreams["Product"] = readstreamhistory(joinpath("streamhistories", "ProdStream.csv"), "Product", syscomps; ismoleflow=true) -histstreams -# And mix them -histstreams["Comb"] = histstreams["Feed"] + histstreams["Product"] +mean(values(c1)) +mean(values(c2)) +mean(values(sc2)) -feeddata = showdata(histstreams["Feed"]); -println(feeddata) -# And scale the result -histstreams["Comb2"] = 2.0 * histstreams["Comb"] +# Lets introduce some errors and check our closure corrections -#-UnitOps etc with historical data----- +copystream!(sysstreams, "C2", "eC2", factor = 1.05) +copystream!(sysstreams, "H2", "eH2", factor = 0.95) +copystream!(sysstreams, "Product", "eProduct") +sysstreams["eMixed"] = emptystream(sysstreams, "eMixed") -# Test UnitOpHistory -histops = UnitOpHistoryList() -@unitophist begin - inlets --> ["Feed", "Product"] - outlets --> ["Comb2"] +# Define the unit ops and boundary +@unitop begin + inlets --> ["eH2", "eC2"] + outlets --> ["eMixed"] calc --> mixer! -end "Mixer" histstreams histops -histops["Mixer"]() +end "eMixer" sysstreams sysunitops +sysunitops["eMixer"]() -histstreams["Comb"].massflows == histstreams["Comb2"].massflows -histstreams["Comb"].comps == histstreams["Comb2"].comps -histstreams["Comb"].timestamps == histstreams["Comb2"].timestamps +@unitop begin + inlets --> ["eMixed"] + outlets --> ["eProduct"] +end "eReactor" sysstreams sysunitops +# Define a mass balance boundary +@boundary begin + unitops --> ["eMixer", "eReactor"] +end b sysunitops -# Test BalanceBoundaryHistory -@unitophist begin - inlets --> ["Feed"] - outlets --> ["Product"] -end "RX101" histstreams histops -@boundaryhist begin - unitops --> ["RX101"] -end bh histops +# Get the correction factors on the inlets and outlets +corrections = calccorrections(b, "eProduct") +b2 = closemb(b, anchor = "eProduct") +sysunitops["eMixer"]() -corrections = calccorrections(bh) -bh = closemb(bh) -conversion(bh, "Ethylene") -selectivity(bh, "Ethylene", "Ethane") +mean(values(b.closure)) +mean(values(b2.closure)) -print(showdata(bh)) +print(showdata(b2)) # Test Flowsheet -sysunitops = UnitOpList() - -@unitop begin - inlets --> ["Feed"] - outlets --> ["Product"] -end "Reactor" sysstreams sysunitops - -@unitop begin - inlets --> ["Product"] - outlets --> ["C2", "H2"] -end "Membrane" sysstreams sysunitops - -@unitop begin - inlets --> ["H2", "C2"] - outlets --> ["Mixed"] - calc --> mixer! -end "Mixer" sysstreams sysunitops - - fs = Flowsheet(sysunitops, ["Reactor"], [1]) -addunitop!(fs, ["Membrane", "Mixer"]) +addunitop!(fs, ["Mixer"]) fs() -sysstreams["Dummy"] = sysstreams["H2"] + sysstreams["C2"] -sysstreams["Dummy"].massflows == sysstreams["Mixed"].massflows -sysstreams["Dummy"].atomflows == sysstreams["Mixed"].atomflows - generateBFD(fs, "./myflowsheet.svg") diff --git a/src/FlowsheetTools.jl b/src/FlowsheetTools.jl index c2aa724..c2af5d9 100644 --- a/src/FlowsheetTools.jl +++ b/src/FlowsheetTools.jl @@ -1,21 +1,45 @@ +# TODO Should macros use ~ instead of ->?? +# TODO replace all @error and @assert with exceptions (where relevant) + module FlowsheetTools -export Component, ComponentList, @comp, writecomponent, readcomponent, readcomponentlist!, - Stream, StreamHistory, StreamList, StreamHistoryList, @stream, copystream!, deletestream!, renamestream!, - copystreamhistory!, deletestreamhistory!, renamestreamhistory!, readstreamhistory, showdata, - UnitOp, UnitOpHistory, UnitOpList, UnitOpHistoryList, - BalanceBoundary, BalanceBoundaryHistory, @boundary, @boundaryhist, +export Component, ComponentList, @comp, writecomponent, readcomponent, readcomponentlist!, names, + Stream, StreamList, @stream, copystream!, deletestream!, renamestream!, renamestream, emptystream, readstreamhistory, refreshcomplist, + UnitOp, UnitOpList, @unitop, mixer!, + BalanceBoundary, @boundary, showdata, calccorrections, closemb, - conversion, selectivity, - @unitop, @unitophist, mixer!, + conversion, molar_selectivity, Flowsheet, addunitop!, setorder!, generateBFD -using MacroTools, Optim, Statistics, Dates, DelimitedFiles, PrettyTables, Base64, Downloads +using Base64, # Used to encode strings to pass Mermaid diagram definition to server + ChangePointDetection, # Used to detect changepoints (step changes) + Dates, # Used with all DateTime types + DelimitedFiles, # Used for reading stream histories + Downloads, # Used to pass Mermaid diagram definition to server + HypothesisTests, # Used for Augmented Dickey-Fuller test to see if data has a slop between changepoints + Interpolations, # Used for extrapolation to endpoints in cleaning up data + Loess, # Used for smoothing and filling data sets + MacroTools, # Used for @comp, @stream etc + Missings, # Used for handling missing data in data cleanup + Optim, # Used to minimize objective function for mass balance reconciliations (replace with JuMP?) + OrderedCollections, # Used for OrderedDict + PrettyTables, # Used for pretty printing + Statistics, # Basic statistical functions + TimeSeries # Used for TimeArray in Streams + + import Base.setindex! import Base.getindex import Base.length import Base.copy +import Base.iterate +import Base.show +import Base.:+ +import Base.:* +# TODO Add these +# import Base.== +# import Base.≈ diff --git a/src/boundaries.jl b/src/boundaries.jl index 240d18a..0d85b21 100644 --- a/src/boundaries.jl +++ b/src/boundaries.jl @@ -9,28 +9,6 @@ struct BalanceBoundary unitlist::UnitOpList units::Vector{String} - # Streams crossing the boundary limits - # Calculate from inlets and outlets of included unit ops - inlets::Vector{String} - outlets::Vector{String} - - # Combined inlet and outlet streams - total_in::Stream - total_out::Stream - - # Total mass balance closure: (out - in)/in - closure::Float64 - - # Elemental closures, Dict{String, Float64} - atomclosures::Dict{String, Float64} -end - - -struct BalanceBoundaryHistory - # Units included in the boundary - unitlist::UnitOpHistoryList - units::Vector{String} - # Number of historic data points numdata::Integer @@ -40,18 +18,20 @@ struct BalanceBoundaryHistory outlets::Vector{String} # Combined inlet and outlet streams - total_in::StreamHistory - total_out::StreamHistory - - # Total mass balance closure: (out - in)/in - closure::Vector{Float64} + total_in::Stream + total_out::Stream - # Elemental closures, Dict{String, Float64} - atomclosures::Vector{Dict{String, Float64}} + # Total mass balance closure: out/in + closure::TimeArray + # Elemental closures: out/in + atomclosures::TimeArray - function BalanceBoundaryHistory(unitlist, units, numdata, inlets, outlets, total_in, total_out, closure, atomclosures) - total_in.numdata != total_out.numdata && error("All in/outlets must must have similar history lengths.") + # Internal constructor to ensurte that all inlets and outlets have the same number of historic data points + # and identical timestamps. It is only required here in case the user doesn't use the outer constructor + function BalanceBoundary(unitlist, units, numdata, inlets, outlets, total_in, total_out, closure, atomclosures) + total_in.numdata != total_out.numdata && error("all in/outlets must must have similar history lengths.") + !all(timestamp(total_in.massflows) .== timestamp(total_out.massflows)) && error("all in/outlets must must have identical timestamps.") new(unitlist, units, numdata, inlets, outlets, total_in, total_out, closure, atomclosures) end end @@ -67,14 +47,17 @@ end """ BalanceBoundary(unitlist, units) -Constructor for a BalanceBoundary. Inputs are a UnitOpList and and array of names of UnitOps in the list +Constructor for a BalanceBoundary. Inputs are a UnitOpList and and array of names of `UnitOp`s in the list that are inside the boundary. Inlet and outlet streams crossing the boundary are automatically calculated. -Fill error if there are either no inlets or outlets. +Will error if there are either no inlets or outlets. + +Since not all atoms referenced in the streams will be present, closures for atoms not present will be indicated by +setting the values to -1.0 """ function BalanceBoundary(unitlist::UnitOpList, units::Vector{String}) # Get the streams that cross the boundary - inlets, outlets = boundarystreams(unitlist, units) + inlets, outlets, _ = boundarystreams(unitlist, units) # Convert to list of names inletnames = String[] @@ -91,65 +74,28 @@ function BalanceBoundary(unitlist::UnitOpList, units::Vector{String}) total_in = sum(inlets) total_out = sum(outlets) - # Calculate the closures - closure = total_out.totalmassflow/total_in.totalmassflow - - atomclosures = Dict{String, Float64}() - for atom in keys(total_in.atomflows) - in = total_in.atomflows[atom] - out = total_out.atomflows[atom] - atomclosures[atom] = out/in - end - - BalanceBoundary(unitlist, units, inletnames, outletnames, total_in, total_out, closure, atomclosures) -end - - -""" - BalanceBoundaryHistory(unitlist, units) - -Constructor for a BalanceBoundaryHistory. Inputs are a UnitOpHistoryList and and array of names of `UnitOpHistory`s in the list -that are inside the boundary. Inlet and outlet streams crossing the boundary are automatically calculated. - -Fill error if there are either no inlets or outlets. -""" -function BalanceBoundaryHistory(unitlist::UnitOpHistoryList, units::Vector{String}) - # Get the streams that cross the boundary - inlets, outlets = boundarystreams(unitlist, units) - - # Convert to list of names - inletnames = String[] - outletnames = String[] - - for inlet in inlets - push!(inletnames, inlet.name) - end - for outlet in outlets - push!(outletnames, outlet.name) - end - - # Now add the inlets and outlets together to get one total inlet and one total outlet - total_in = sum(inlets) - total_out = sum(outlets) + # Ensure that all inlets and outlets have the same number of historic data points at the same timestamps + total_in.numdata != total_out.numdata && error("all in/outlets must must have similar history lengths.") + !all(timestamp(total_in.massflows) .== timestamp(total_out.massflows)) && error("all in/outlets must must have identical timestamps.") # Calculate the closures - closure = total_out.totalmassflow./total_in.totalmassflow + closure = total_out.totalmassflow ./ total_in.totalmassflow + rename!(closure, Symbol("Total Mass Closure")) numdata = length(closure) - atomclosurehistory = Vector{Dict{String, Float64}}(undef, numdata) - + atomclosures = Vector{Dict{String, Float64}}(undef, numdata) for datum in 1:numdata - atomclosures = Dict{String, Float64}() - for atom in keys(total_in.atomflows[datum]) - in = total_in.atomflows[datum][atom] - out = total_out.atomflows[datum][atom] - atomclosures[atom] = out/in + _atomclosures = Dict{String, Float64}() + for atom in keys(values(total_in.atomflows)[datum]) + in = values(total_in.atomflows)[datum][atom] + out = values(total_out.atomflows)[datum][atom] + _atomclosures[atom] = (in == 0.0) ? -1.0 : out/in end - atomclosurehistory[datum] = atomclosures + atomclosures[datum] = _atomclosures end - - BalanceBoundaryHistory(unitlist, units, numdata, inletnames, outletnames, total_in, total_out, closure, atomclosurehistory) + atomclosures_ta = TimeArray(timestamp(closure), atomclosures, [Symbol("Elemental Closures")]) + BalanceBoundary(unitlist, units, numdata, inletnames, outletnames, total_in, total_out, closure, atomclosures_ta) end @@ -164,31 +110,16 @@ end function Base.show(io::IO, b::BalanceBoundary) println(io, "Balance Boundary:\n") println(io, "Enclosed units: ", [u for u in b.units]) - println(io, "Closure: ", prettyround(b.closure)) println(io) - print(io, "Combined Feed ") - println(io, b.total_in) - print(io, "Combined Product ") - println(io, b.total_out) - println(io, "\nElemental closures:") - println(io, "="^19) - - atoms = collect(keys(b.atomclosures)) - closures = collect(values(b.atomclosures)) - for i in eachindex(atoms) - println(io, " ", rpad(atoms[i],2), lpad(prettyround(closures[i]), 14)) - end -end - - -# Pretty printing for BalanceBoundary objects -function Base.show(io::IO, b::BalanceBoundaryHistory) - println(io, "Balance Boundary:\n") - println(io, "Enclosed units: ", [u for u in b.units]) + println(io, "Closure:") + pretty_table(io, b.closure, display_size=(14, -1)) println(io) - println(io, "Data length:\t", b.numdata) - println(io, "Data starts:\t", b.total_in.timestamps[begin]) - println(io, "Data ends:\t", b.total_in.timestamps[begin]) + println(io, "Combined Feed Mass Flows:") + pretty_table(io, b.total_in.massflows, display_size=(14, -1)) + println(io, "Combined Product Mass Flows:") + pretty_table(io, b.total_out.massflows, display_size=(14, -1)) + println(io, "\nElemental closures (-1.0 if not present):") + pretty_table(io, b.atomclosures, display_size=(14, -1)) end @@ -201,46 +132,29 @@ end """ @boundary begin unitops --> ["Reactor", "Membrane"] - end b sysunitops + end b1 sysunitops -Create a boundary, b, that includes UnitOps "Reactor" amd "Membrane" from the UnitOpsList sysunitops. -""" -macro boundary(ex::Expr, name::Symbol, unitoplist::Symbol) - local unitops = String[] - - for line in ex.args - match_comp = @capture(line, unitops --> [ous__]) - if match_comp - for ou in ous - push!(unitops, ou) - end - end - end - - return :($(esc(name)) = BalanceBoundary($(esc(unitoplist)), $unitops)) -end +Create a boundary, b1, that includes UnitOps "Reactor" amd "Membrane" from the UnitOpsList sysunitops. - -""" @boundaryhist begin unitops --> ["RX101"] - end bh histops + end b2 sysunitstops -Create a boundary, bh, that includes UnitOps "RX101" from the UnitOpsList histops. +Create a boundary, b2, that includes UnitOps "RX101" from the UnitOpsList sysunitops. """ -macro boundaryhist(ex::Expr, name::Symbol, unitoplist::Symbol) +macro boundary(ex::Expr, name::Symbol, unitoplist::Symbol) local unitops = String[] for line in ex.args - match_comp = @capture(line, unitops --> [ous__]) + match_comp = @capture(line, unitops --> [uops__]) if match_comp - for ou in ous - push!(unitops, ou) + for uop in uops + push!(unitops, uop) end end end - return :($(esc(name)) = BalanceBoundaryHistory($(esc(unitoplist)), $unitops)) + return :($(esc(name)) = BalanceBoundary($(esc(unitoplist)), $unitops)) end @@ -252,32 +166,25 @@ end """ - showdata(bh::BalanceBoundaryHistory) + showdata(b::BalanceBoundary) -Returns a data table for a BalanceBoundaryHistory, as a String. +Returns a data table for a BalanceBoundary, as a String. Example: - print(showdata(bh)) + print(showdata(b)) """ -function showdata(bh::BalanceBoundaryHistory) - titles_in = copy(bh.total_in.comps) - titles_in .= "In:" .* titles_in - titles_out = copy(vcat(bh.total_out.comps)) - titles_out .= "Out:" .* titles_out - titles = vcat("Timestamp", titles_in, titles_out) +function showdata(b::BalanceBoundary) + invals = pretty_table(String, b.total_in.massflows) + outvals = pretty_table(String, b.total_out.massflows) - data_in = hcat(bh.total_in.timestamps, bh.total_in.massflows') - data_out = bh.total_out.massflows' - data = hcat(data_in, data_out) - - str = pretty_table(String, data, header=titles) + str = "Mass Balance Boundary:\n" * "-"^22 * "\n\nTotal Mass Flows In:\n" * invals * "\n\n" * "Total Mass Flows In:\n" * outvals return str end """ - inlets, outlets = boundarystreams(unitlist, units) + inlets, outlets, internals = boundarystreams(unitlist, units) Find the streams that enter and leave the boundary that includes the specified unitops. """ @@ -313,7 +220,7 @@ function boundarystreams(unitlist::UnitOpList, units::Vector{String}) for (i, stream) in enumerate(inlets) if stream in outlets keep_in[i] = false - push!(internals, stream) # only need this in this loop, as any internal stream is some block's inlet + push!(internals, stream) end end for (i, stream) in enumerate(outlets) @@ -330,54 +237,4 @@ function boundarystreams(unitlist::UnitOpList, units::Vector{String}) @assert length(outlets) > 0 "zero outlet streams from the boundary" return inlets, outlets, internals -end - - -function boundarystreams(unitlist::UnitOpHistoryList, units::Vector{String}) - #Figure out which streams cross the boundary: - # 1. Take all the input stream from the units - # 2. Subtract the ones that are product stream, as they will be internal streams - # Do the same for outlet streams, just the other way around - # 3. Calculate the elemental closures - - inlets = StreamHistory[] - outlets = StreamHistory[] - - # 1. Take all the input stream from the units - for unitname in units - unit = unitlist[unitname] - - for feedname in unit.inlets - feed = unit.streamlist[feedname] - push!(inlets, feed) - end - for prodname in unit.outlets - prod = unit.streamlist[prodname] - push!(outlets, prod) - end - end - - # 2a. Check which ones should be kept. - # Don't delete now, or the check on outlets will miss the ones already deleted from inlets! - keep_in = trues(length(inlets)) - keep_out = trues(length(outlets)) - for (i, stream) in enumerate(inlets) - if stream in outlets - keep_in[i] = false - end - end - for (i, stream) in enumerate(outlets) - if stream in inlets - keep_out[i] = false - end - end - - # 2b. Now delete the ones we don't want - inlets = inlets[keep_in] - outlets = outlets[keep_out] - - @assert length(inlets) > 0 "zero inlet streams to the boundary" - @assert length(outlets) > 0 "zero outlet streams from the boundary" - - return inlets, outlets end \ No newline at end of file diff --git a/src/closure.jl b/src/closure.jl index e412d58..42b6f36 100644 --- a/src/closure.jl +++ b/src/closure.jl @@ -5,62 +5,75 @@ #---------------------------------------------------------------------------- + """ - function calccorrections(boundary::BalanceBoundary; totalweight=1.0, elementweight=1.0) - function calccorrections(boundary::BalanceBoundaryHistory; totalweight=1.0, elementweight=1.0) + function calccorrections(boundary::BalanceBoundary, anchor::String; totalweight=1.0, elementweight=1.0) Basic mass balance reconciliation. Error in total mass closure and average element closures are weighted by *totalweight* and *elementweight* respectively and the squared weighted error is minimized. Results are returned as a dict of streams and corrections to their flows. """ -function calccorrections(boundary::BalanceBoundary; totalweight=1.0, elementweight=1.0) +function calccorrections(boundary::BalanceBoundary, anchor::String; totalweight=1.0, elementweight=1.0) # Pull the streamlist of the first unit op in the list. Since this is from a UnitOpList, # all of the unit ops must have the same stream list - streamlist = first(boundary.unitlist.list).second.streamlist - + streamlist = first(boundary.unitlist.list).second.streamlist + corrections = Dict{String, Float64}() - + inlets = boundary.inlets outlets = boundary.outlets ins = length(inlets) outs = length(outlets) - - factors = ones(ins + outs) # Correction factors for flows - # Check for viable closure, i.e. there is at least some mass entering and some mass leaving. - allzero = true - for stream in inlets - streamlist[stream].totalmassflow > 0.0 && (allzero = false) + if anchor in inlets + anchorinlet = true + elseif anchor in outlets + anchorinlet = false + else + error("anchor $anchor not in inlets or outlets") end - allzero && throw(DomainError(x, "at least one inlet must have non-zero flow")) - - allzero = true - for stream in inlets - streamlist[stream].totalmassflow > 0.0 && (allzero = false) + + if anchorinlet + anchorindex = findfirst(isequal(anchor), inlets) + else + anchorindex = findfirst(isequal(anchor), outlets) + ins end - allzero && throw(DomainError(x, "at least one inlet must have non-zero flow")) + + + factors = ones(ins + outs - 1) # No factor for the anchor stream + + numdata = streamlist[inlets[1]].numdata + + function f(_factors) + factors = vcat(_factors[1:anchorindex-1], 1.0, _factors[anchorindex:end]) - function f(factors) + # Since inlets and outlets are arrays of Stream, summing them produces Stream objects total_in = sum(factors[1:ins] .* streamlist[inlets]) total_out = sum(factors[ins+1:end] .* streamlist[outlets]) - masserror = abs2(total_out.totalmassflow/total_in.totalmassflow - 1.0) - - atomclosures = Dict{String, Float64}() - for atom in keys(total_in.atomflows) - in = total_in.atomflows[atom] - out = total_out.atomflows[atom] - atomclosures[atom] = out/in + + + masserrors = (total_out.totalmassflow ./ total_in.totalmassflow) .- 1.0 + masserror = sum(abs2, values(masserrors)) + + atomerror = 0.0 + for datum = 1:numdata + for atom in keys(values(total_in.atomflows)[datum]) + inflow = values(total_in.atomflows)[datum][atom] + if inflow > 0.0 + outflow = values(total_out.atomflows)[datum][atom] + atomerror += abs2(outflow/inflow - 1.0) + end + end end - numatoms = length(atomclosures) - atomerror = abs2(sum(values(atomclosures)) - numatoms) - return totalweight*masserror + elementweight*atomerror + totalerr = totalweight*masserror + elementweight*atomerror + return totalerr end - - res = optimize(f, factors) - optfactors = res.minimizer + res = optimize(f, factors, BFGS()) + _optfactors = res.minimizer + optfactors = vcat(_optfactors[1:anchorindex-1], 1.0, _optfactors[anchorindex:end]) i = 1 for stream in inlets @@ -75,11 +88,8 @@ function calccorrections(boundary::BalanceBoundary; totalweight=1.0, elementweig return corrections end - """ - function closemb(boundary::BalanceBoundary, [corrections::Dict{String, Float64}]) - function closemb(boundary::BalanceBoundaryHistory, [corrections = Dict{String, Float64}]) Apply the mass balance reconciliation. The corrections can first calculated using `calccorrections` and can then be applied to multiple boundaries using this function. If no corrections are passed, @@ -87,83 +97,11 @@ and can then be applied to multiple boundaries using this function. If no correc Since balance boundaries are immutable, a new boundary instance is returned. """ -function closemb(boundary::BalanceBoundary, corrections=nothing; totalweight=1.0, elementweight=1.0) - isnothing(corrections) && (corrections = calccorrections(boundary; totalweight, elementweight)) - - # Pull the streamlist of the first unit op in the list. Since this is from a UnitOpList, - # all of the unit ops must have the same stream list - streamlist = first(boundary.unitlist.list).second.streamlist - - # Apply the stream corrections - for stream in keys(corrections) - streamlist[stream] *= corrections[stream] - end - - #Recreate the boundary - newboundary = BalanceBoundary(boundary.unitlist, boundary.units) - - return newboundary -end - - -function calccorrections(boundary::BalanceBoundaryHistory; totalweight=1.0, elementweight=1.0) - # Pull the streamlist of the first unit op in the list. Since this is from a UnitOpList, - # all of the unit ops must have the same stream list - streamlist = first(boundary.unitlist.list).second.streamlist - - corrections = Dict{String, Float64}() - - inlets = boundary.inlets - outlets = boundary.outlets - - ins = length(inlets) - outs = length(outlets) - - factors = ones(ins + outs) - - numdata = streamlist[inlets[1]].numdata - - function f(factors) - # Since inlets and outlets are arrays of StreamHistory, summing them produces StreamHistory objects - total_in = sum(factors[1:ins] .* streamlist[inlets]) - total_out = sum(factors[ins+1:end] .* streamlist[outlets]) - - masserror = 0.0 - atomerror = 0.0 - for datum = 1:numdata - masserror += abs2(total_out.totalmassflow[datum]/total_in.totalmassflow[datum] - 1.0) - - atomclosures = Dict{String, Float64}() - for atom in keys(total_in.atomflows[datum]) - in = total_in.atomflows[datum][atom] - out = total_out.atomflows[datum][atom] - atomclosures[atom] = out/in - end - numatoms = length(atomclosures) - atomerror += abs2(sum(values(atomclosures)) - numatoms) - end - return totalweight*masserror + elementweight*atomerror - end - - res = optimize(f, factors) - optfactors = res.minimizer - - i = 1 - for stream in inlets - corrections[stream] = optfactors[i] - i += 1 +function closemb(boundary::BalanceBoundary; corrections=nothing, anchor=nothing, totalweight=1.0, elementweight=1.0) + if isnothing(corrections) + isnothing(anchor) && error("an anchor stream must be specified") + corrections = calccorrections(boundary, anchor; totalweight, elementweight) end - for stream in outlets - corrections[stream] = optfactors[i] - i += 1 - end - - return corrections -end - - -function closemb(boundary::BalanceBoundaryHistory, corrections=nothing; totalweight=1.0, elementweight=1.0) - isnothing(corrections) && (corrections = calccorrections(boundary; totalweight, elementweight)) # Pull the streamlist of the first unit op in the list. Since this is from a UnitOpList, # all of the unit ops must have the same stream list @@ -175,7 +113,7 @@ function closemb(boundary::BalanceBoundaryHistory, corrections=nothing; totalwei end #Recreate the boundary - newboundary = BalanceBoundaryHistory(boundary.unitlist, boundary.units) + newboundary = BalanceBoundary(boundary.unitlist, boundary.units) return newboundary end \ No newline at end of file diff --git a/src/components.jl b/src/components.jl index a2eb978..a5db69e 100644 --- a/src/components.jl +++ b/src/components.jl @@ -4,16 +4,15 @@ # #---------------------------------------------------------------------------- - struct Component - name - atoms - counts + name::String + atoms::Vector{String} + counts::Vector{Int} Mr::Float64 end struct ComponentList - list + list::OrderedDict{String, Component} end @@ -28,7 +27,7 @@ end Component(name, atoms, counts) Constructor for molecular species that defines the atomic make-up. -It is recommended to rather use the @comp macro to create components. +It is recommended to rather use the @comp macro to create components interactively """ function Component(name, atoms, counts) Mr = 0.0 @@ -38,13 +37,14 @@ function Component(name, atoms, counts) return Component(name, atoms, counts, Mr) end + """ ComponentList() Constructor for an empty component list. Components are added when created via the @comp macro """ function ComponentList() - l = Dict{String, Component}() + l = OrderedDict{String, Component}() return ComponentList(l) end @@ -70,8 +70,21 @@ function Base.length(A::ComponentList) return length(A.list) end +function Base.iterate(A::ComponentList) + return iterate(A.list) +end + +function Base.iterate(A::ComponentList, state) + return iterate(A.list, state) +end + +#---------------------------------------------------------------------------- +# +#----Pretty Printing--------------------------------------------------------- +# +#---------------------------------------------------------------------------- + -# Pretty printing for component objects function Base.show(io::IO, c::Component) println(io, "Component: $(c.name)\n") println(io, " Atom\t\tCount") @@ -210,4 +223,14 @@ function readcomponentlist!(complist::ComponentList, folder::String, filenames:: end return count +end + + +""" + names(complist::ComponentList) + +Return the list of names of components in the ComponenList. +""" +function names(complist::ComponentList) + return collect(keys(complist.list)) end \ No newline at end of file diff --git a/src/datacleaning.jl b/src/datacleaning.jl index 8abd7bc..9ce3d15 100644 --- a/src/datacleaning.jl +++ b/src/datacleaning.jl @@ -1,7 +1,7 @@ -using Dates, Loess, Interpolations, Missings, TimeSeries +# using Dates, Loess, Interpolations, Missings, TimeSeries -# These are used only during testing -using Plots, Distributions +# # These are used only during testing +# using Plots, Distributions """ calcHoL(timestamps) @@ -17,22 +17,22 @@ end # Generate dummy data for testing -function gendata(timestamps, period, fracfilled) - HoL = calcHoL(timestamps) - data = zeros(Union{Float64, Missing}, length(times)) - - for i in eachindex(HoL) - if rand() > fracfilled - data[i] = missing - else - data[i] = sin(π*period*(HoL[i]/(Hour(endtime - starttime)/Hour(1)))) - norm = Normal(0, 0.1*abs(data[i])) - data[i] += rand(norm) - end - end - - return data -end +# function gendata(timestamps, period, fracfilled) +# HoL = calcHoL(timestamps) +# data = zeros(Union{Float64, Missing}, length(times)) + +# for i in eachindex(HoL) +# if rand() > fracfilled +# data[i] = missing +# else +# data[i] = sin(π*period*(HoL[i]/(Hour(endtime - starttime)/Hour(1)))) +# norm = Normal(0, 0.1*abs(data[i])) +# data[i] += rand(norm) +# end +# end + +# return data +# end """ smooth_and_fill(raw; loessspan=0.3, suggest_start=false, startval=0.0, suggest_end=false, endval=0.0) @@ -79,15 +79,15 @@ function smooth_and_fill(raw; loessspan=0.3, suggest_start=false, startvals=Floa end -# Generate dummy data with missing values -starttime = DateTime(2023, 1, 1, 0, 0) -endtime = DateTime(2023, 1, 12, 24, 0) -times = starttime:Hour(6):endtime +# # Generate dummy data with missing values +# starttime = DateTime(2023, 1, 1, 0, 0) +# endtime = DateTime(2023, 1, 12, 24, 0) +# times = starttime:Hour(6):endtime -raw = TimeArray(times, [gendata(times, 2, 0.95) gendata(times, 2, 0.95)], [:raw1, :raw2]) -plot(raw, leg=:bottomleft) +# raw = TimeArray(times, [gendata(times, 2, 0.95) gendata(times, 2, 0.95)], [:raw1, :raw2]) +# plot(raw, leg=:bottomleft) -sf1 = smooth_and_fill(raw, suggest_start=true, suggest_end=true, startvals=[0.0, 0.0], endvals=[0.0, 0.0]) -sf2 = smooth_and_fill(raw) -plot!(sf1) -plot!(sf2) +# sf1 = smooth_and_fill(raw, suggest_start=true, suggest_end=true, startvals=[0.0, 0.0], endvals=[0.0, 0.0]) +# sf2 = smooth_and_fill(raw) +# plot!(sf1) +# plot!(sf2) diff --git a/src/flowsheets.jl b/src/flowsheets.jl index a01e796..34440e0 100644 --- a/src/flowsheets.jl +++ b/src/flowsheets.jl @@ -8,12 +8,13 @@ struct Flowsheet unitlist::UnitOpList units::Vector{String} - order::Vector{Int} + execution_order::Vector{Int} end + function (fs::Flowsheet)(neworder = nothing) if isnothing(neworder) - for i in fs.order + for i in fs.execution_order fs.unitlist[fs.units[i]]() end else @@ -27,6 +28,28 @@ function (fs::Flowsheet)(neworder = nothing) end +#---------------------------------------------------------------------------- +# +#----Base overloads---------------------------------------------------------- +# +#---------------------------------------------------------------------------- + + +function Base.show(io::IO, fs::Flowsheet) + d = OrderedDict{Int, String}() + for i in 1:length(fs.execution_order) + d[fs.execution_order[i]] = fs.units[i] + end + sort!(d) + + println(io, "Flowsheet:") + println(io) + println(io, "Execution Order:") + for (k, v) in d + println(k, " ", v) + end +end + #---------------------------------------------------------------------------- # #----Utilities--------------------------------------------------------------- @@ -37,7 +60,7 @@ end function addunitop!(fs::Flowsheet, u::String) if haskey(fs.unitlist.list, u) push!(fs.units, u) - push!(fs.order, length(fs.order) + 1) + push!(fs.execution_order, length(fs.execution_order) + 1) else error("unitop $u not defined in list $(fs.unitlist)") end @@ -50,7 +73,7 @@ function addunitop!(fs::Flowsheet, us::Vector{String}) for u in us if haskey(fs.unitlist.list, u) push!(fs.units, u) - push!(fs.order, length(fs.order) + 1) + push!(fs.execution_order, length(fs.execution_order) + 1) else error("unitop $u not defined in list $(fs.unitlist)") end @@ -61,26 +84,26 @@ end function setorder!(fs::Flowsheet, neworder) - numold = length(fs.order) + numold = length(fs.execution_order) numnew = length(neworder) if numold == numnew for i in 1:numold - fs.order[i] = neworder[i] + fs.execution_order[i] = neworder[i] end elseif numold < numnew for i in 1:numold - fs.order[i] = neworder[i] + fs.execution_order[i] = neworder[i] end for j = 1:(numnew - numold) - push!(fs.order, neworder[numold+j]) + push!(fs.execution_order, neworder[numold+j]) end else for i in 1:numold - fs.order[i] = neworder[i] + fs.execution_order[i] = neworder[i] end for j = 1:(numold - numnew) - pop!(fs.order) + pop!(fs.execution_order) end end diff --git a/src/kpis.jl b/src/kpis.jl index e023170..4b8c13b 100644 --- a/src/kpis.jl +++ b/src/kpis.jl @@ -7,142 +7,93 @@ """ conversion(b::BalanceBoundary, component::String) - conversion(b::BalanceBoundaryHistory, component::String) Calculate the fractional conversion of the component over the balance boundary """ function conversion(b::BalanceBoundary, component::String) # Find the component in the feed - i = findfirst(x->x==component, b.total_in.comps) - if isnothing(i) - feed = 0.0 + if component in keys(b.total_in.complist.list) + feeds = values(b.total_in.massflows[Symbol(component)]) else - feed = b.total_in.massflows[i] + feeds = zeros(length(b.total_in.massflows)) end # Find the component in the product - i = findfirst(x->x==component, b.total_out.comps) - if isnothing(i) - prod = 0 + if component in keys(b.total_in.complist.list) + prods = values(b.total_out.massflows[Symbol(component)]) else - prod = b.total_out.massflows[i] + prods = zeros(length(b.total_out.massflows)) end - return feed > 0.0 ? (feed - prod)/feed : 0.0 -end - - -function conversion(b::BalanceBoundaryHistory, component::String) - # Find the component in the feed - i = findfirst(x->x==component, b.total_in.comps) - - # Conversion zero if nothing in the feed - if isnothing(i) - results = zeros(b.total_in.numdata) - return results + conversions = zeros(length(b.total_in.massflows)) + for (i, (feed, prod)) in enumerate(zip(feeds, prods)) + if feed > 0.0 + conversions[i] = (feed - prod)/feed + else + conversions[i] = 0.0 + end end - # Find the component in the product - j = findfirst(x->x==component, b.total_out.comps) - - # If nothing left in the product, then return ones - if isnothing(j) - results = ones(b.total_in.numdata) - return results - end - - # Create place to put the results. - results = zeros(b.total_in.numdata) - - for datum = 1:b.total_in.numdata - feed = b.total_in.massflows[i, datum] - prod = b.total_out.massflows[j, datum] - results[datum] = feed > 0.0 ? (feed - prod)/feed : 0.0 + if length(conversions) == 1 + return conversions[1] + else + timestamps = timestamp(b.total_in.massflows) + conversions_ta = TimeArray(timestamps, conversions, [Symbol("Conversion: " * component)]) + return conversions_ta end - - return results end """ - function selectivity(b, reactant, product) - function selectivity(b::BalanceBoundaryHistory, reactant::String, product::String) + function molar_selectivity(b, reactant, product) Calculate the molar selectivity of the converted reactant to the product, over the balance boundary (b). """ -function selectivity(b::BalanceBoundary, reactant::String, product::String) +function molar_selectivity(b::BalanceBoundary, reactant::String, product::String) # Find the reactant inflow - i = findfirst(x->x==reactant, b.total_in.comps) - if isnothing(i) - return 0.0 + if reactant in keys(b.total_in.complist.list) + r_ins = values(b.total_in.moleflows[Symbol(reactant)]) else - r_in = b.total_in.moleflows[i] + r_ins = zeros(length(b.total_in.moleflows)) end # Find the reactant outflow - i = findfirst(x->x==reactant, b.total_out.comps) - if isnothing(i) - r_out = 0.0 + if reactant in keys(b.total_out.complist.list) + r_outs = values(b.total_out.moleflows[Symbol(reactant)]) else - r_out = b.total_out.moleflows[i] + r_ins = zeros(length(b.total_out.moleflows)) end # Find the product inflow - i = findfirst(x->x==product, b.total_in.comps) - if isnothing(i) - p_in = 0 + if product in keys(b.total_in.complist.list) + p_ins = values(b.total_in.moleflows[Symbol(product)]) else - p_in = b.total_in.moleflows[i] + p_ins = zeros(length(b.total_in.moleflows)) end # Find the product outflow - i = findfirst(x->x==product, b.total_out.comps) - if isnothing(i) - p_out = 0.0 + if product in keys(b.total_out.complist.list) + p_outs = values(b.total_out.moleflows[Symbol(product)]) else - p_out = b.total_out.moleflows[i] + p_ins = zeros(length(b.total_out.moleflows)) end - return r_in > 0.0 ? (p_out - p_in)/(r_in - r_out) : 0.0 -end - - -function selectivity(b::BalanceBoundaryHistory, reactant::String, product::String) - results = zeros(b.total_in.numdata) - - # Find the reactant in the feed - i = findfirst(x->x==reactant, b.total_in.comps) - # Selectivity zero if nothing in the feed - isnothing(i) && (return results) - r_in = b.total_in.moleflows[i, :] - - # Find the reactant in the product - i = findfirst(x->x==reactant, b.total_out.comps) - if isnothing(i) - r_out = zeros(b.total_out.numdata) - else - r_out = b.total_out.moleflows[i, :] + selectivities = zeros(length(b.total_in.massflows)) + for (i, (r_in, r_out, p_in, p_out)) in enumerate(zip(r_ins, r_outs, p_ins, p_outs)) + if r_in > 0.0 + selectivities[i] = (p_out - p_in)/(r_in - r_out) + else + selectivities[i] = 0.0 + end end - # Find the product in the feed - i = findfirst(x->x==product, b.total_in.comps) - if isnothing(i) - p_in = zeros(b.total_in.numdata) + if length(selectivities) == 1 + return selectivities[1] else - p_in = b.total_in.moleflows[i, :] - end - - # Find the product in the product - i = findfirst(x->x==product, b.total_out.comps) - if isnothing(i) - p_out = zeros(b.total_out.numdata) - else - p_out = b.total_out.moleflows[i, :] - end - - for datum = 1:b.total_in.numdata - results[datum] = r_in[datum] > 0.0 ? (p_out[datum] - p_in[datum])/(r_in[datum] - r_out[datum]) : 0.0 + timestamps = timestamp(b.total_in.massflows) + selectivities_ta = TimeArray(timestamps, selectivities, [Symbol("Molar selectivity: " * reactant * " -> " * product)]) + return selectivities_ta end - return results -end \ No newline at end of file + return selectivities +end diff --git a/src/stepchange.jl b/src/stepchange.jl index da2f8f9..adede57 100644 --- a/src/stepchange.jl +++ b/src/stepchange.jl @@ -1,3 +1,5 @@ +# TODO Kill this code and replace with ChangePointDetection.changepoints(ts; threshold = 0.5, window = 150)... + """ function findsteps(data) diff --git a/src/streams.jl b/src/streams.jl index 0b0dcea..24025d4 100644 --- a/src/streams.jl +++ b/src/streams.jl @@ -1,3 +1,5 @@ +# TODO Add uncertainty / variance to each stream for various flows + #---------------------------------------------------------------------------- # #----Definitions------------------------------------------------------------- @@ -8,47 +10,21 @@ struct Stream name::String - complist::ComponentList - comps::Vector{String} - - # Specifiy mass flows, calculate mole flows - massflows::Vector{Float64} - moleflows::Vector{Float64} - totalmassflow::Float64 - - # Molar flow of atoms, calculated from above - atomflows::Dict{String, Float64} -end - - -struct StreamHistory - name::String - # Number of historic data points numdata::Integer complist::ComponentList - comps::Vector{String} - - timestamps::Vector{DateTime} - # Specifiy mass flows, calculate mole flows - massflows::Matrix{Float64} - moleflows::Matrix{Float64} - totalmassflow::Vector{Float64} - - # Molar flow of atoms, calculated from above - atomflows::Vector{Dict{String, Float64}} + # Flow data, mass, molar and atoms + massflows::TimeArray + moleflows::TimeArray + totalmassflow::TimeArray + atomflows::TimeArray end struct StreamList - list -end - - -struct StreamHistoryList - list + list::OrderedDict{String, Stream} end @@ -60,128 +36,108 @@ end """ - function Stream(name, complist, comps, flows, ismoleflow=false) + function Stream(name, complist, comps, timestamps, flowdata, ismoleflow=false) -Constructor for a stream that defines the stream name and component flows. -The mass flows are specified and a molar composition and atomic molar flows calculated. +Constructor for a stream history object that defines the stream name and component flows +for various past measurements. + +The mass flows are specified and molar flows and atomic molar flows calculated + OR +The mole flows are specified and mass flows and atomic molar flows calculated + +The flow data is passed as a matrix where each column represents a component. Internally, +flows for all the components in the specified component list is stored, but only a sublist +need be specified in the constructor. Other components are assigned zero flows. This is to +expidite the addition of streams using the internal TimeArray objects. -It is recommeded to rather use the @stream macro to create streams. +It is recommended to rather create StreamHistory objects via readstreamhistory(). """ -function Stream(name, complist, comps, flows, ismoleflow=false) +function Stream(name, complist, comps, timestamps, flowdata, ismoleflow=false) numcomps = length(comps) - - length(flows) != numcomps && error("mismatch between number of components and available data.") - - massflows = zeros(numcomps) - moleflows = zeros(numcomps) + numdata = size(flowdata, 1) - atomflows = Dict{String, Float64}() - totalmassflow = 0.0 - - for (i, compname) in enumerate(comps) + # Sanity checks on the data + size(flowdata, 2) != numcomps && error("mismatch between number of components and available data.") + length(timestamps) != numdata && error("length mismatch between data and timestamps.") + + # Build the TimeArrays + allcomps = names(complist) # Also those not present in the stream + massflows = zeros(numdata, length(allcomps)) + moleflows = zeros(numdata, length(allcomps)) + + # Collect the atoms to create a blank dictionary - all values are zero + # This Dict only carries atoms actually present in the stream + # When adding streams, this list is recreated from the new list of components + emptyatomflows = Dict{String, Float64}() + for compname in allcomps comp = complist[compname] - - if !ismoleflow - massflows[i] = flows[i] - moleflows[i] = massflows[i]/comp.Mr - else - moleflows[i] = flows[i] - massflows[i] = moleflows[i]*comp.Mr - end - - totalmassflow += massflows[i] - - for (j, atom) in enumerate(comp.atoms) - if atom in keys(atomflows) - atomflows[atom] += moleflows[i]*comp.counts[j] - else - atomflows[atom] = moleflows[i]*comp.counts[j] + for atom in comp.atoms + if atom ∉ keys(emptyatomflows) + emptyatomflows[atom] = 0.0 end end - end + - return Stream(name, complist, comps, massflows, moleflows, totalmassflow, atomflows) -end - - -""" - function StreamHistory(name, complist, comps, timestamps, flowshistory, ismoleflow=false) - -Constructor for a stream history object that defines the stream name and component flows -for various past measurements. -The mass flows are specified and a molar composition and atomic molar flows calculated. -The mass flow history is passed as a matrix where each column is a datum - -It is recommended to rather create StreamHistory objects via readstreamhistory(filename, streamname, complist). -""" -function StreamHistory(name, complist, comps, timestamps, flowshistory, ismoleflow=false) - numcomps = length(comps) - numdata = size(flowshistory, 2) - - size(flowshistory, 1) != numcomps && error("mismatch between number of components and available data.") - length(timestamps) != numdata && error("length mismatch between data and timestamps.") - - massflowshistory = zeros(numcomps, numdata) - moleflowshistory = zeros(numcomps, numdata) - - atomflowshistory = Vector{Dict{String, Float64}}(undef, numdata) - totalmassflowhistory = zeros(numdata) + atomflows = Vector{Dict{String, Float64}}(undef, numdata) + totalmassflows = zeros(numdata) for datum in 1:numdata - atomflows = Dict{String, Float64}() - - for (i, compname) in enumerate(comps) - comp = complist[compname] - - if !ismoleflow - massflowshistory[i, datum] = flowshistory[i, datum] - moleflowshistory[i, datum] = flowshistory[i, datum]/comp.Mr - else - moleflowshistory[i, datum] = flowshistory[i, datum] - massflowshistory[i, datum] = flowshistory[i, datum]*comp.Mr - end - - totalmassflowhistory[datum] += massflowshistory[i, datum] - - for (j, atom) in enumerate(comp.atoms) - if atom in keys(atomflows) - atomflows[atom] += moleflowshistory[i, datum]*comp.counts[j] + # Get a copy of the empty atomic flows list + _atomflows = deepcopy(emptyatomflows) + + # Here we loop through all the components in the complist and + # use zeros when they are not present + for (i, compname) in enumerate(allcomps) + if compname in comps + comp = complist[compname] + idx = findfirst(x->x==compname, comps) + + if !ismoleflow + massflows[datum, i] = flowdata[datum, idx] + moleflows[datum, i] = flowdata[datum, idx]/comp.Mr else - atomflows[atom] = moleflowshistory[i, datum]*comp.counts[j] + moleflows[datum, i] = flowdata[datum, idx] + massflows[datum, i] = flowdata[datum, idx]*comp.Mr end + + totalmassflows[datum] += massflows[datum, i] + + for (j, atom) in enumerate(comp.atoms) + _atomflows[atom] += moleflows[datum, i]*comp.counts[j] + end + else + massflows[datum, i] = 0.0 + moleflows[datum, i] = 0.0 end end - atomflowshistory[datum] = atomflows + atomflows[datum] = _atomflows end - return StreamHistory(name, numdata, complist, comps, timestamps, massflowshistory, moleflowshistory, - totalmassflowhistory, atomflowshistory) + massflows_ta = TimeArray(timestamps, massflows, allcomps) + moleflows_ta = TimeArray(timestamps, moleflows, allcomps) + totalmassflows_ta = TimeArray(timestamps, totalmassflows, [:totalmassflows]) + atomflows_ta = TimeArray(timestamps, atomflows, [:atomflows]) + + + return Stream(name, numdata, complist, massflows_ta, moleflows_ta, totalmassflows_ta, atomflows_ta) end + """ Streamlist() -Constructor for an empty stream list. Streams are added when created via the @stream macro +Constructor for an empty stream list. Streams are added when created via the @stream macro or +when read from file with readstreamhistory() """ function StreamList() - l = Dict{String, Stream}() + l = OrderedDict{String, Stream}() return StreamList(l) end -""" - StreamHistorylist() - -Constructor for an empty stream list with historical data. Streams are added when read from file with readstreamhistory(). -""" -function StreamHistoryList() - l = Dict{String, StreamHistory}() - return StreamHistoryList(l) -end - #---------------------------------------------------------------------------- # #----Base overloads---------------------------------------------------------- @@ -189,110 +145,6 @@ end #---------------------------------------------------------------------------- -""" - Base.:+(a::Stream, b::Stream) - -Extend the addition operator to add to streams to each other - a mixer. -It is assumed that the streams will have different components in arbitrary order. -""" -function Base.:+(a::Stream, b::Stream) - # Make sure the streams use the same system components - a.complist != b.complist && error("cannot add streams with different system component lists") - - comps = copy(a.comps) - massflows = copy(a.massflows) - - compmap = indexin(b.comps, comps) - for i in eachindex(b.comps) - if isnothing(compmap[i]) - # Add the component - push!(comps, b.comps[i]) - push!(massflows, b.massflows[i]) - else - # Add the flow to the existing component - massflows[compmap[i]] += b.massflows[i] - end - end - - return Stream(a.name * "+" * b.name, a.complist, comps, massflows) -end - - -""" - Base.+(a::StreamHistory, b::StreamHistory) - -Extend the addition operator to add to streams histories to each other - a mixer. -It is assumed that the streams histories will have different components in arbitrary order. -""" -function Base.:+(a::StreamHistory, b::StreamHistory) - # Make sure the streams use the same system components - a.complist != b.complist && error("cannot add streams with different system component lists") - - # Check that the data peridos are the same. - # The constructor already checks that the timestamp length matches the data length. - # We don't check the comps, as the streams could have different compositions, which we mix. - a.numdata != b.numdata && error("history length of a not identical to that of b") - - # Check that the timestamps are identical! - for i in eachindex(a.timestamps) - a.timestamps[i] != b.timestamps[i] && error("timestamp values do not match at entry $i") - end - - comps = copy(a.comps) - massflows = copy(a.massflows) - - compmap = indexin(b.comps, comps) - for i in eachindex(b.comps) - if isnothing(compmap[i]) - # Add the component - push!(comps, b.comps[i]) - vcat(massflows, b.massflows[i, :]') - else - # Add the flow to the existing component - massflows[compmap[i], :] .+= b.massflows[i, :] - end - end - - return StreamHistory(a.name * "+" * b.name, a.complist, comps, a.timestamps, massflows) -end - - -""" - Base.*(a::T, b::Stream) where T <: Real - -Extend the multiplication operator to scale a stream's flows by a scalar value. -Used in mass balance reconciliations to apply flow corrections. -""" -function Base.:*(a::T, b::Stream) where T <: Real - return Stream(b.name, b.complist, b.comps, a .* b.massflows) -end - -function Base.:*(b::Stream, a::T) where T <: Real - return Stream(b.name, b.complist, b.comps, a .* b.massflows) -end - -""" - Base.*(a::T, b::StreamHistory) where T <: Real - -Extend the multiplication operator to scale a stream history's flows by a scalar value. -Used in mass balance reconciliations to apply flow corrections. -""" -function Base.:*(a::T, b::StreamHistory) where T <: Real - return StreamHistory(b.name, b.complist, b.comps, b.timestamps, a .* b.massflows) -end - - -""" - Base.*(b::StreamHistory, a::T) where T <: Real - -Extend the multiplication operator to scale a stream history's flows by a scalar value. -Used in mass balance reconciliations to apply flow corrections. -""" -function Base.:*(b::StreamHistory, a::T) where T <: Real - return StreamHistory(b.name, b.complist, b.comps, b.timestamps, a .* b.massflows) -end - - function Base.setindex!(A::StreamList, X::Stream, idx::String) if length(A.list) == 0 A.list[idx] = X @@ -302,7 +154,10 @@ function Base.setindex!(A::StreamList, X::Stream, idx::String) # We get the value entry using the `second` field of the `Pair`, which returns a `Stream`, # of which we get the `complist` field. currentlist = first(A.list).second.complist + current_ts = timestamp(first(A.list).second.massflows) + X.complist != currentlist && error("all streams in StreamList must reference the same ComponentList") + !all(timestamp(X.massflows) .== current_ts) && error("all streams in StreamList must have the same timestamps") A.list[idx] = X end @@ -330,109 +185,128 @@ function Base.length(A::StreamList) end -function Base.setindex!(A::StreamHistoryList, X::StreamHistory, idx::String) - if length(A.list) == 0 - A.list[idx] = X - else - # Verify that all the streams reference the same ComponentList. - # Get the first item in the `list` field. As `list` is a `Dict`, this returns a `Pair`. - # We get the value entry using the `second` field of the `Pair`, which returns a `StreamHistory`, - # of which we get the `complist` field. - currentlist = first(A.list).second.complist - X.complist != currentlist && error("all streams in StreamHistoryList must reference the same ComponentList") +""" + Base.:+(a::Stream, b::Stream) - A.list[idx] = X - end -end +Extend the addition operator to add to streams to each other - a mixer. +All streams must refer to the same component list . +""" +function Base.:+(a::Stream, b::Stream) + # Make sure the streams use the same system components and timestamps + # TimeArrays will add only matching timestamps, but this result in varying data lengths, which must be avoided + a.complist != b.complist && error("cannot add streams with different system component lists") + a.numdata != b.numdata && error("cannot add streams with different data lengths") + !all(timestamp(a.massflows) .== timestamp(b.massflows)) && error("cannot add streams with different timestamps") + comps = string.(colnames(a.massflows)) + timestamps = timestamp(a.massflows) + flowdata = values(a.massflows .+ b.massflows) -function Base.getindex(A::StreamHistoryList, idx::String) - return A.list[idx] + + return Stream(a.name * "+" * b.name, a.complist, comps, timestamps, flowdata) end -function Base.getindex(A::StreamHistoryList, idxs::Vector{String}) - res = StreamHistory[] - for idx in idxs - push!(res, A.list[idx]) - end - return res +""" + Base.*(a::T, b::Stream) where T <: Real + +Extend the multiplication operator to scale a stream's flows by a scalar value. +Used in mass balance reconciliations to apply flow corrections. +""" +function Base.:*(a::T, b::Stream) where T <: Real + comps = string.(colnames(b.massflows)) + timestamps = timestamp(b.massflows) + flowdata = values(a .* b.massflows) + + return Stream(b.name, b.complist, comps, timestamps, flowdata) end -function Base.length(A::StreamHistoryList) - return length(A.list) +function Base.:*(b::Stream, a::T) where T <: Real + comps = string.(colnames(b.massflows)) + timestamps = timestamp(b.massflows) + flowdata = values(a .* b.massflows) + + return Stream(b.name, b.complist, comps, timestamps, flowdata) end function Base.copy(A::Stream) - return Stream(A.name, A.complist, A.comps, A.massflows, A.moleflows, A.totalmassflow, A.atomflows) + return deepcopy(A) end -# Pretty printing for stream objects -function Base.show(io::IO, s::Stream) - println(io, "Stream: ", s.name, "\nTotal mass flow: ", prettyround(s.totalmassflow), "\n") - println(io, "Component\tMass Flow\tMolar Flow") - println(io, "-"^42) - - for (i, compname) in enumerate(s.comps) - comp = s.complist[compname] - println(io, " ", rpad(comp.name, 9), "\t", rpad(prettyround(s.massflows[i]), 9), "\t", prettyround(s.moleflows[i])) - end - +# # Pretty printing for stream objects +function Base.show(io::IO, stream::Stream) + println(io, "Stream: $(stream.name)") println(io) - println(io, "Atom\t\tFlow") - println(io, "-"^20) - atoms = collect(keys(s.atomflows)) - flows = collect(values(s.atomflows)) - for i in eachindex(atoms) - println(io, " ", rpad(atoms[i],2), lpad(prettyround(flows[i]), 17)) - end -end - - -# Pretty printing for stream history objects -function Base.show(io::IO, s::StreamHistory) - println(io, "StreamHistory: $(s.name)") - for compname in s.comps - comp = s.complist[compname] - println(io, " ", rpad(comp.name, 9)) + if stream.numdata == 1 + header = string.(colnames(stream.massflows)) + pushfirst!(header, " ") + massflows = ["Mass flows " (values(stream.massflows))] + moleflows = ["Molar flows" (values(stream.moleflows))] + data = vcat(massflows, moleflows) + pretty_table(io, data, header = header) + println(io) + println(io, "Total mass flow: $(prettyround(values(stream.totalmassflow)[begin]))") + println(io) + pretty_table(io, values(stream.atomflows)[1], header = ["Atom", "Molar flow"]) + else + println(io, "Mass flows:") + pretty_table(io, stream.massflows, display_size=(14, -1)) + println(io) + println(io, "Molar flows:") + pretty_table(io, stream.moleflows, display_size=(14, -1)) + println(io) + println(io) + println(io, "Atom Flows:") + pretty_table(io, stream.atomflows, display_size=(14, -1)) + println(io) + println(io, "Data length:\t$(stream.numdata)") + println(io, "Data starts:\t$(timestamp(stream.massflows)[begin])") + println(io, "Data ends:\t$(timestamp(stream.massflows)[end])") end - - println(io) - println(io, "Data length:$(length(s.totalmassflow))") - println(io, "Data starts:\t$(s.timestamps[begin])") - println(io, "Data ends:\t$(s.timestamps[end])") end function Base.show(io::IO, sl::StreamList) println(io, "Stream list:") - for strm in sl.list - println(io, " ", strm.first) - end -end - - -function Base.show(io::IO, sl::StreamHistoryList) - println(io, "Stream list:") + println(io) + println(io, "Streams:") if length(sl.list) > 0 - for strm in sl.list - println(io, " ", strm.first) + for (name, _) in sl + println(io, " ", name) + end + + stream = first(sl.list).second + + println(io) + println(io, "Components:") + for (name, _) in stream.complist + println(io, " ", name) end + println(io) - strm = first(sl.list) - println(io, "Data length:\t$(length(strm.second.totalmassflow))") - println(io, "Data starts:\t$(strm.second.timestamps[begin])") - println(io, "Data ends:\t$(strm.second.timestamps[end])") + println(io, "Data length:\t$(stream.numdata)") + println(io, "Data starts:\t$(timestamp(stream.massflows)[begin])") + println(io, "Data ends:\t$(timestamp(stream.massflows)[end])") else - println(io, "Empty list") + println(io, "\tEmpty list") end end +function Base.iterate(A::StreamList) + return iterate(A.list) +end + + +function Base.iterate(A::StreamList, state) + return iterate(A.list, state) +end + + # ---------------------------------------------------------------------------- # #----Macros------------------------------------------------------------------- @@ -448,29 +322,31 @@ The names of components must match those in the specified componentlist (`syscomps::ComponentList` in the example). - @stream "mass" begin + @stream mass begin "Ethylene" --> 2.8053 "Ethane" --> 27.06192 "Hydrogen" --> 2.21738 end "Test" syscomps sysstreams - @stream "mole" begin + @stream mole begin "Ethylene" --> 0.1 "Ethane" --> 0.9 "Hydrogen" --> 1.1 end "Product" syscomps sysstreams +The created Stream object will have data at a single timestamp of DateTime(0). """ macro stream(flowtype::Symbol, ex::Expr, name::String, complist::Symbol, streamlist::Symbol) local comps = String[] local flows = Float64[] + local timestamps = [DateTime(0)] if flowtype == :mass local ismoleflow = false elseif flowtype == :mole local ismoleflow = true else - error("flow basis specification must be mass or mole") + error("flow basis specification must be \"mass\" or \"mole\"") end for line in ex.args @@ -487,7 +363,7 @@ macro stream(flowtype::Symbol, ex::Expr, name::String, complist::Symbol, streaml end end end - return :($(esc(streamlist))[$name] = Stream($name, $(esc(complist)), $comps, $flows, $ismoleflow)) + return :($(esc(streamlist))[$name] = Stream($name, $(esc(complist)), $comps, $timestamps, $(flows'), $ismoleflow)) end @@ -499,13 +375,17 @@ end """ - function copystream(list::StreamList, from::String, to::String) + function copystream!(list::StreamList, from::String, to::String) Copy a stream in the stream list """ -function copystream!(list::StreamList, from::String, to::String; factor=1.0) - str = factor*list[from] - list[to] = Stream(to, str.complist, str.comps, str.massflows, str.moleflows, str.totalmassflow, str.atomflows) +function copystream!(list::StreamList, from::String, to::String; factor=1.0) + fromstream = list[from] + comps = string.(colnames(fromstream.massflows)) + timestamps = timestamp(fromstream.massflows) + flowdata = factor .* values(fromstream.massflows) + + list[to] = Stream(to, fromstream.complist, comps, timestamps, flowdata) return nothing end @@ -523,84 +403,50 @@ function deletestream!(list::StreamList, from::String) end + """ function renamestream!(list::StreamList, from::String, to::String) Rename a stream in the stream list """ function renamestream!(list::StreamList, from::String, to::String) - str = list[from] - list[to] = Stream(to, str.complist, str.comps, str.massflows, str.moleflows, str.totalmassflow, str.atomflows) - delete!(list.list, from) + fromstream = list[from] + comps = string.(colnames(fromstream.massflows)) + timestamps = timestamp(fromstream.massflows) + flowdata = values(fromstream.massflows) + list[to] = Stream(to, fromstream.complist, comps, timestamps, flowdata) + deletestream!(list, from) return nothing end - """ - function renamestream(str::Stream, newname::String) + renamestream(strm::Stream, to::String) -Returns a copy of the stream with the new name. Internal use only, not exported. +Return a copy of the stream with a new name """ -function renamestream(str::Stream, newname::String) - return Stream(newname, str.complist, str.comps, str.massflows, str.moleflows, str.totalmassflow, str.atomflows) -end - +function renamestream(fromstream::Stream, to::String) + comps = string.(colnames(fromstream.massflows)) + timestamps = timestamp(fromstream.massflows) + flowdata = values(fromstream.massflows) -function renamestream(str::StreamHistory, newname::String) - return StreamHistory(newname, str.numdata, str.complist, str.comps, str.timestamps, str.massflows, str.moleflows, str.totalmassflow, str.atomflows) + return Stream(to, fromstream.complist, comps, timestamps, flowdata) end """ - function copystreamhistory(list::StreamHistoryList, from::String, to::String) + emptystream(list::StreamList, name::String) -Copy a StreamHistory in the StreamHistoryList +Returns a stream with the same components and timestamp as the specified StreamList and all flows set to zero. """ -function copystreamhistory!(list::StreamHistory, from::String, to::String; factor=1.0) - str = factor*list[from] - list[to] = StreamHistory(to, str.numdata, str.complist, str.comps, str.timestamps, str.massflowshistory, - str.moleflowshistory, str.totalmassflowhistory, str.atomflowshistory) +function emptystream(list::StreamList, name::String) + # Take the first stream in the StreamList as a reference + refstrm = first(list).second - return nothing + return Stream(name, refstrm.complist, string.(colnames(refstrm.massflows)), timestamp(refstrm.massflows), zeros(size(refstrm.massflows))) end -""" - function deletestreamhistory!(list::StreamHistoryList, from::String) - -Delete a StreamHistory from the StreamHistoryList -""" -function deletestreamhistory!(list::StreamHistoryList, from::String) - delete!(list.list, from) - - return nothing -end - - -""" - function renamestream!(list::StreamHistoryList, from::String, to::String) - -Rename a StreamHistory in the StreamHistoryList -""" -function renamestreamhistory!(list::StreamHistoryList, from::String, to::String) - str = list[from] - list[to] = StreamHistory(to, str.numdata, str.complist, str.comps, str.timestamps, str.massflowshistory, - str.moleflowshistory, str.totalmassflowhistory, str.atomflowshistory) - delete!(list.list, from) - - return nothing -end - -""" - function renamestreamhistory(str::Stream, newname::String) - -Returns a copy of the stream with the new name. Internal use only, not exported. -""" -function renamestreamhistory(str::StreamHistory, newname::String) - return StreamHistory(newname, str.numdata, str.complist, str.comps, str.timestamps, str.massflows, str.moleflows, str.totalmassflow, str.atomflows) -end - """ readstreamhistory(filename, streamname, complist) @@ -621,26 +467,22 @@ function readstreamhistory(filename, streamname, complist; ismoleflow=false) flows = data[:, 2:end] timestamps = DateTime.(data[:, 1], "yyyy/mm/dd HH:MM") - return StreamHistory(streamname, complist, comps, timestamps, transpose(flows), ismoleflow) + return Stream(streamname, complist, comps, timestamps, flows, ismoleflow) end -""" - showdata(hs::StreamHistory) -Display a table of the component flow history from hs::StreamHistory. """ -function showdata(hs::StreamHistory) - titles = vcat("Timestamp", hs.comps) - data = hcat(hs.timestamps, hs.massflows') + refreshcomplist(streamlist::StreamList) - str = pretty_table(String, data, header=titles) +Refresh all the streams in the StreamList to reflect components added to the StreamList after the streams were created. +""" +function refreshcomplist(streamlist::StreamList) + for (name, stream) in streamlist + complist = stream.complist + comps = string.(colnames(stream.massflows)) + timestamps = timestamp(stream.massflows) + flowdata = values(stream.massflows) - return str + streamlist[name] = Stream(name, complist, comps, timestamps, flowdata) + end end - - - - - - - diff --git a/src/unitops.jl b/src/unitops.jl index 2f930ac..ab58584 100644 --- a/src/unitops.jl +++ b/src/unitops.jl @@ -1,3 +1,6 @@ +# TODO Add a splitter unit +# TODO Add a stoichiometric reactor with multiple reactions and conversions + #---------------------------------------------------------------------------- # #----Definitions------------------------------------------------------------- @@ -9,55 +12,31 @@ struct UnitOp name::String streamlist::StreamList - - # Connected streams - inlets::Vector{String} - outlets::Vector{String} - - f!::Function - params -end - -function (u::UnitOp)(newparams = nothing) - if isnothing(newparams) - u.f!(u.streamlist, u.outlets, u.inlets, u.params) - else - u.f!(u.streamlist, u.outlets, u.inlets, newparams) - end -end - - -struct UnitOpHistory - name::String - - streamlist::StreamHistoryList # Connected streams with history inlets::Vector{String} outlets::Vector{String} - f!::Function - params + f!::Function # An optional function that will write to outlets from inlets and params + params # Default paramaters can be added upon creation of the object # Inner constructor to validate input data - function UnitOpHistory(name, streamlist, inlets, outlets, f!, params) - numdata = length(streamlist[inlets[1]].totalmassflow) - timestamps = streamlist[inlets[1]].timestamps + function UnitOp(name, streamlist, inlets, outlets, f!, params) + numdata = streamlist[inlets[1]].numdata + timestamps = timestamp(streamlist[inlets[1]].massflows) for inletname in inlets inlet = streamlist[inletname] + inlet_ts = timestamp(inlet.massflows) - inlet.numdata != numdata && error("all streams must must have similar history lengths.") - for j in eachindex(timestamps) - timestamps[j] != inlet.timestamps[j] && error("all stream histories must have matching timestamps.") - end + inlet.numdata != numdata && error("all streams must must have identical history lengths.") + !all(inlet_ts .== timestamps) && error("all streams must must have identical timestamps.") end for outletname in outlets outlet = streamlist[outletname] + outlet_ts = timestamp(outlet.massflows) - outlet.numdata != numdata && error("all streams must must have similar history lengths.") - for j in eachindex(timestamps) - timestamps[j] != outlet.timestamps[j] && error("all stream histories must have matching timestamps.") - end + outlet.numdata != numdata && error("all streams must must have identical history lengths.") + !all(outlet_ts .== timestamps) && error("all streams must must have identical timestamps.") end new(name, streamlist, inlets, outlets, f!, params) @@ -65,7 +44,7 @@ struct UnitOpHistory end -function (u::UnitOpHistory)(newparams = nothing) +function (u::UnitOp)(newparams = nothing) if isnothing(newparams) u.f!(u.streamlist, u.outlets, u.inlets, u.params) else @@ -75,12 +54,7 @@ end struct UnitOpList - list::Dict{String, UnitOp} -end - - -struct UnitOpHistoryList - list::Dict{String, UnitOpHistory} + list::OrderedDict{String, UnitOp} end @@ -107,45 +81,17 @@ function UnitOp(name::String, streamlist::StreamList, inlets::Vector{String}, ou end -""" - UnitOpHistory(name::String, streamlist::StreamHistoryList, inlets::Vector{String}, outlets::Vector{String}) - UnitOpHistory(name::String, streamlist::StreamHistoryList, inlets::Vector{String}, outlets::Vector{String}, f!::Function) - UnitOpHistory(name::String, streamlist::StreamHistoryList, inlets::Vector{String}, outlets::Vector{String}, f!::Function, params) - -Constructor for a UnitOp. It is recommended to rather use the @unitop macro. -""" -function UnitOpHistory(name::String, streamlist::StreamHistoryList, inlets::Vector{String}, outlets::Vector{String}) - return UnitOpHistory(name, streamlist, inlets, outlets, passive, nothing) -end - - -function UnitOpHistory(name::String, streamlist::StreamHistoryList, inlets::Vector{String}, outlets::Vector{String}, f!::Function) - return UnitOpHistory(name, streamlist, inlets, outlets, f!, nothing) -end - - """ UnitOpList() Constructor for an empty UnitOpList. Unit operations are added when created with the @unitop macro. """ function UnitOpList() - l = Dict{String, UnitOp}() + l = OrderedDict{String, UnitOp}() return UnitOpList(l) end -""" - UnitOpHistoryList() - -Constructor for an empty UnitOpHistoryList. Unit operations are added when created with the @unitophist macro. -""" -function UnitOpHistoryList() - l = Dict{String, UnitOpHistory}() - return UnitOpHistoryList(l) -end - - #---------------------------------------------------------------------------- # #----Active UnitOps---------------------------------------------------------- @@ -157,14 +103,13 @@ end Default calculation for UnitOps. Simply a placeholder - does no calculation. """ -function passive(streamlist::StreamList, outlets::Vector{String}, inlets::Vector{String}, params) +@inline function passive(streamlist::StreamList, outlets::Vector{String}, inlets::Vector{String}, params) return nothing end """ mixer!(streamlist::StreamList, outlets::Vector{String}, inlets::Vector{String}, params) - mixer!(streamlist::StreamHistoryList, outlets::Vector{String}, inlets::Vector{String}, params) Calculation for mixer UnitOps. Combines all feed streams into a single outlet stream. Will error if more than one outlet stream is defined. Values in `params` are ignored. @@ -172,33 +117,11 @@ Will error if more than one outlet stream is defined. Values in `params` are ign function mixer!(streamlist::StreamList, outlets::Vector{String}, inlets::Vector{String}, params) @assert length(outlets) == 1 "mixers can have only one outlet stream" - old = streamlist[outlets[1]] - - tempstream = Stream(old.name, old.complist, String[], Float64[]) - for inlet in inlets[1:end] - tempstream += streamlist[inlet] - end - - tempstream = renamestream(tempstream, old.name) - streamlist[old.name] = tempstream - - return nothing -end - - -function mixer!(streamlist::StreamHistoryList, outlets::Vector{String}, inlets::Vector{String}, params) - @assert length(outlets) == 1 "mixers can have only one outlet stream" - - old = streamlist[inlets[1]] - tempstream = StreamHistory(old.name, old.numdata, old.complist, old.comps, old.timestamps, - copy(old.massflows), copy(old.moleflows), copy(old.totalmassflow), copy(old.atomflows)) - - for inlet in inlets[2:end] - tempstream += streamlist[inlet] - end + name = streamlist[outlets[1]].name + tempstream = sum(streamlist[inlets]) - tempstream = renamestream(tempstream, outlets[1]) - streamlist[outlets[1]] = tempstream + tempstream = renamestream(tempstream, name) + streamlist[name] = tempstream return nothing end @@ -232,24 +155,17 @@ function Base.getindex(A::UnitOpList, idx::String) end -function Base.setindex!(A::UnitOpHistoryList, X::UnitOpHistory, idx::String) - if length(A.list) == 0 - A.list[idx] = X - else - # Verify that all the streams reference the same ComponentList. - # Get the first item in the `list` field. As `list` is a `Dict`, this returns a `Pair`. - # We get the value entry using the `second` field of the `Pair`, which returns a `UnitOpHistory`, - # of which we get the `streamlist` field. - currentlist = first(A.list).second.streamlist - X.streamlist != currentlist && error("all streams in UnitOpHistoryList must reference the same StreamHistoryList") - - A.list[idx] = X +function Base.getindex(A::UnitOpList, idxs::Vector{String}) + res = UnitOp[] + for idx in idxs + push!(res, A.list[idx]) end + return res end -function Base.getindex(A::UnitOpHistoryList, idx::String) - return A.list[idx] +function Base.length(A::UnitOpList) + return length(A.list) end @@ -259,50 +175,56 @@ function Base.show(io::IO, u::UnitOp) println(io, "Feed streams: ", [s for s in u.inlets]) println(io) println(io, "Product streams: ", [s for s in u.outlets]) -end - - -# Pretty printing for UnitOp objects -function Base.show(io::IO, u::UnitOpHistory) - println(io, "Unit Operation: $(u.name)") - println(io, "Feed streams: ", [s for s in u.inlets]) - println(io) - println(io, "Product streams: ", [s for s in u.outlets]) - println(io) - strm = first(u.streamlist.list) - println(io, "Data length:\t$(strm.second.numdata)") - println(io, "Data starts:\t$(strm.second.timestamps[begin])") - println(io, "Data ends:\t$(strm.second.timestamps[end])") -end - - -function Base.show(io::IO, uol::UnitOpList) - println(io, "UnitOp list:") - for unitop in uol.list - println(io, " ", unitop.first) + + strm = first(u.streamlist).second + + # If streams hold only one data point, don't print this + if strm.numdata > 1 + ts = timestamp(strm.massflows) + println(io) + println(io, "Data length:\t$(strm.numdata)") + println(io, "Data starts:\t$(ts[begin])") + println(io, "Data ends:\t$(ts[end])") end end -function Base.show(io::IO, uol::UnitOpHistoryList) - println(io, "UnitOpHistory list:") +function Base.show(io::IO, uol::UnitOpList) + println(io, "UnitOperation list:") if length(uol.list) > 0 for unitop in uol.list println(io, " ", unitop.first) end - println(io) - sl = first(uol.list).second.streamlist - inlets = first(uol.list).second.inlets + + uo = first(uol.list).second + sl = uo.streamlist + inlets = uo.inlets strm = sl[inlets[1]] - println(io, "Data length:\t$(strm.numdata)") - println(io, "Data starts:\t$(strm.timestamps[begin])") - println(io, "Data ends:\t$(strm.timestamps[end])") + + # If streams hold only one data point, don't print this + if strm.numdata > 1 + println(io) + println(io, "Data length:\t$(strm.numdata)") + ts = timestamp(strm.massflows) + println(io, "Data starts:\t$(ts[begin])") + println(io, "Data ends:\t$(ts[end])") + end else println(io, "Empty list") end end +function Base.iterate(A::UnitOpList) + return iterate(A.list) +end + + +function Base.iterate(A::UnitOpList, state) + return iterate(A.list, state) +end + + #---------------------------------------------------------------------------- # #----Macros------------------------------------------------------------------ @@ -325,7 +247,9 @@ The specified streams refer to entries in sysstreams::StreamList. Two optional parameters may be specified, i.e. calc and params: calc is a function, e.g. + function mixer!(streamlist::StreamList, outlets::Vector{String}, inlets::Vector{String}, params) + params is an iterable passed to the function. Together, these specify a calculation to perform to calculate the outlet(s) from the inlet(s), when @@ -374,69 +298,3 @@ macro unitop(ex::Expr, name::String, streamlist::Symbol, unitoplist::Symbol) return :($(esc(unitoplist))[$name] = UnitOp($name, $(esc(streamlist)), $inlets, $outlets, $(esc(func)), $(esc(params)))) end end - - -""" -Defines a `UnitOpHist with` the specified name and connected streams. - -@unitophist begin - inlets --> ["H2", "C2"] - outlets --> ["Mixed"] - calc --> mixer! -end "Mixer" histstreams histunitops - -This will create a UnitOpHist named "Mixer", saved into histunitops["Mixer"]. -The specified streams refer to entries in histstreams::StreamHistoryList. - -Two optional parameters may be specified, i.e. calc and params: - -calc is a function, e.g. - function mixer!(streamlist::StreamHistoryList, outlets::Vector{String}, inlets::Vector{String}, params) -params is an iterable passed to the function. - -Together, these specify a calculation to perform to calculate the outlet(s) from the inlet(s), when -calling `histunitops["Mixer"]()` - -Alternative parameters may also be specified in the call, e.g. - histunitops["Mixer"]([1, 2, 3]) - -""" -macro unitophist(ex::Expr, name::String, streamlist::Symbol, unitoplist::Symbol) - local inlets = String[] - local outlets = String[] - local func = nothing - local params = nothing - - for line in ex.args - match_comp = @capture(line, inlets --> [ins__]) - if match_comp - for strm in ins - push!(inlets, strm) - end - end - - match_comp = @capture(line, outlets --> [outs__]) - if match_comp - for strm in outs - push!(outlets, strm) - end - end - - match_comp = @capture(line, calc --> unitopfunc_) - if match_comp - func = unitopfunc - end - - match_comp = @capture(line, params --> unitopparams_) - if match_comp - params = unitopparams - end - - end - - if isnothing(func) - return :($(esc(unitoplist))[$name] = UnitOpHistory($name, $(esc(streamlist)), $inlets, $outlets)) - else - return :($(esc(unitoplist))[$name] = UnitOpHistory($name, $(esc(streamlist)), $inlets, $outlets, $(esc(func)), $(esc(params)))) - end -end \ No newline at end of file diff --git a/streamhistories/C2.csv b/streamhistories/C2.csv new file mode 100644 index 0000000..1ebb40f --- /dev/null +++ b/streamhistories/C2.csv @@ -0,0 +1,28 @@ +TimeStamp,Ethylene,Ethane, +2020/01/01 00:00,0.107239036,0.89247781,0.999716846 +2020/01/01 06:00,0.107378006,0.896736538,1.004114544 +2020/01/02 00:00,0.092180185,0.985485177,1.077665362 +2020/01/02 06:00,0.106090983,0.970531008,1.076621991 +2020/01/03 00:00,0.103429382,0.923285307,1.026714689 +2020/01/03 06:00,0.102459451,0.949043158,1.051502609 +2020/01/04 00:00,0.109831854,0.885418833,0.995250687 +2020/01/04 06:00,0.10869931,0.904365099,1.013064409 +2020/01/05 00:00,0.097697742,0.830520181,0.928217923 +2020/01/05 06:00,0.102912242,0.978773788,1.08168603 +2020/01/06 00:00,0.103911053,0.870373916,0.974284969 +2020/01/06 06:00,0.105671277,0.868468211,0.974139488 +2020/01/07 00:00,0.09985486,0.82563294,0.9254878 +2020/01/07 06:00,0.100074529,0.9229478,1.023022329 +2020/01/08 00:00,0.094945703,0.825959831,0.920905534 +2020/01/08 06:00,0.096199477,0.987651733,1.08385121 +2020/01/09 00:00,0.092142514,0.889116851,0.981259365 +2020/01/09 06:00,0.095967638,0.963560661,1.059528299 +2020/01/10 00:00,0.108760194,0.967972888,1.076733082 +2020/01/10 06:00,0.096962919,0.930838356,1.027801275 +2020/01/11 00:00,0.095538445,0.824699673,0.920238118 +2020/01/11 06:00,0.101902136,0.941316514,1.04321865 +2020/01/12 00:00,0.107373055,0.934839923,1.042212978 +2020/01/12 06:00,0.102622132,0.816601807,0.919223939 +2020/01/13 00:00,0.097522962,0.813489671,0.911012633 +2020/01/13 06:00,0.109016434,0.839104313,0.948120747 +2020/01/14 00:00,0.106539973,0.8548783,0.961418273 diff --git a/streamhistories/FeedStreamShort.csv b/streamhistories/FeedStreamShort.csv new file mode 100644 index 0000000..3940698 --- /dev/null +++ b/streamhistories/FeedStreamShort.csv @@ -0,0 +1,2 @@ +TimeStamp,Ethylene,Ethane,Hydrogen +2020/01/01 00:00,2.106769468,1.003401256,1.982216977 diff --git a/streamhistories/Hydrogen.csv b/streamhistories/Hydrogen.csv new file mode 100644 index 0000000..b3372b7 --- /dev/null +++ b/streamhistories/Hydrogen.csv @@ -0,0 +1,28 @@ +TimeStamp,Hydrogen +2020/01/01 00:00,1.117239036 +2020/01/01 06:00,1.097378006 +2020/01/02 00:00,1.102180185 +2020/01/02 06:00,1.096090983 +2020/01/03 00:00,1.113429382 +2020/01/03 06:00,1.092459451 +2020/01/04 00:00,1.119831854 +2020/01/04 06:00,1.09869931 +2020/01/05 00:00,1.107697742 +2020/01/05 06:00,1.092912242 +2020/01/06 00:00,1.113911053 +2020/01/06 06:00,1.095671277 +2020/01/07 00:00,1.10985486 +2020/01/07 06:00,1.090074529 +2020/01/08 00:00,1.104945703 +2020/01/08 06:00,1.086199477 +2020/01/09 00:00,1.102142514 +2020/01/09 06:00,1.085967638 +2020/01/10 00:00,1.118760194 +2020/01/10 06:00,1.086962919 +2020/01/11 00:00,1.105538445 +2020/01/11 06:00,1.091902136 +2020/01/12 00:00,1.117373055 +2020/01/12 06:00,1.092622132 +2020/01/13 00:00,1.107522962 +2020/01/13 06:00,1.099016434 +2020/01/14 00:00,1.116539973 diff --git a/streamhistories/Product.csv b/streamhistories/Product.csv new file mode 100644 index 0000000..c5207b3 --- /dev/null +++ b/streamhistories/Product.csv @@ -0,0 +1,28 @@ +TimeStamp,Ethylene,Ethane,Hydrogen +2020/01/01 00:00,0,0.999716846,1.01 +2020/01/01 06:00,0,1.004114544,0.99 +2020/01/02 00:00,0,1.077665362,1.01 +2020/01/02 06:00,0,1.076621991,0.99 +2020/01/03 00:00,0,1.026714689,1.01 +2020/01/03 06:00,0,1.051502609,0.99 +2020/01/04 00:00,0,0.995250687,1.01 +2020/01/04 06:00,0,1.013064409,0.99 +2020/01/05 00:00,0.001,0.927217923,1.011 +2020/01/05 06:00,0,1.08168603,0.99 +2020/01/06 00:00,0,0.974284969,1.01 +2020/01/06 06:00,0,0.974139488,0.99 +2020/01/07 00:00,0.002,0.9234878,1.012 +2020/01/07 06:00,0.00001,1.023012329,0.99001 +2020/01/08 00:00,0,0.920905534,1.01 +2020/01/08 06:00,0,1.08385121,0.99 +2020/01/09 00:00,0,0.981259365,1.01 +2020/01/09 06:00,0.00003,1.059498299,0.99003 +2020/01/10 00:00,0,1.076733082,1.01 +2020/01/10 06:00,0,1.027801275,0.99 +2020/01/11 00:00,0,0.920238118,1.01 +2020/01/11 06:00,0,1.04321865,0.99 +2020/01/12 00:00,0,1.042212978,1.01 +2020/01/12 06:00,0,0.919223939,0.99 +2020/01/13 00:00,0,0.911012633,1.01 +2020/01/13 06:00,0,0.948120747,0.99 +2020/01/14 00:00,0,0.961418273,1.01 diff --git a/streamhistories/Wrongstream.csv b/streamhistories/Wrongstream.csv new file mode 100644 index 0000000..d81b3a3 --- /dev/null +++ b/streamhistories/Wrongstream.csv @@ -0,0 +1,28 @@ +TimeStamp,Ethylene,Ethane,Hydrogen +2020/01/01 00:10,2.106769468,1.003401256,1.982216977 +2020/01/01 06:00,2.131540904,1.013929486,2.060710012 +2020/01/02 00:00,2.185135078,0.900772132,2.089783814 +2020/01/02 06:00,2.121437851,0.996612314,1.949841967 +2020/01/03 00:00,1.870060192,0.908472049,2.133315741 +2020/01/03 06:00,2.162483792,1.017292333,1.981289862 +2020/01/04 00:00,1.874144316,1.016528818,1.869917964 +2020/01/04 06:00,1.930813671,0.987594037,2.057931188 +2020/01/05 00:00,2.12724289,0.923263342,1.942614287 +2020/01/05 06:00,2.166490061,1.09140564,2.103809748 +2020/01/06 00:00,2.171592849,0.938031976,2.082197957 +2020/01/06 06:00,1.804783711,1.009065835,2.174490206 +2020/01/07 00:00,1.957542321,0.968165676,2.0137447 +2020/01/07 06:00,2.090437725,0.959660933,1.981954587 +2020/01/08 00:00,1.826464608,1.032311353,1.970649718 +2020/01/08 06:00,2.0636046,0.929639599,1.981938883 +2020/01/09 00:00,1.939313926,1.021765505,1.83542421 +2020/01/09 06:00,1.812793332,1.014771253,2.025132329 +2020/01/10 00:00,2.154109894,0.94432342,1.807169447 +2020/01/10 06:00,2.046941449,1.072374428,2.04317426 +2020/01/11 00:00,2.123098553,1.039538123,2.153854047 +2020/01/11 06:00,2.164622176,0.924480721,1.844844328 +2020/01/12 00:00,2.029604398,1.096570384,2.008174739 +2020/01/12 06:00,1.812750071,0.93871181,2.080720078 +2020/01/13 00:00,2.040306101,1.054645559,1.832314174 +2020/01/13 06:00,1.825082427,1.086792417,1.870012827 +2020/01/14 00:00,2.071117783,0.923618539,1.946087723