diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..94f1119 --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +.DS_Store +.vscode diff --git a/Manifest.toml b/Manifest.toml index b1f5d1f..0f15f87 100644 --- a/Manifest.toml +++ b/Manifest.toml @@ -12,6 +12,9 @@ git-tree-sha1 = "84918055d15b3114ede17ac6a7182f68870c16f7" uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" version = "3.3.1" +[[ArgTools]] +uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" + [[ArnoldiMethod]] deps = ["LinearAlgebra", "Random", "StaticArrays"] git-tree-sha1 = "f87e559f87a45bece9c9ed97458d3afe98b1ebb9" @@ -19,10 +22,7 @@ uuid = "ec485272-7323-5ecc-a04f-4719b315124d" version = "0.1.0" [[Artifacts]] -deps = ["Pkg"] -git-tree-sha1 = "c30985d8821e0cd73870b17b0ed0ce6dc44cb744" uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" -version = "1.3.0" [[AxisAlgorithms]] deps = ["LinearAlgebra", "Random", "SparseArrays", "WoodburyMatrices"] @@ -41,9 +41,9 @@ version = "1.0.1" [[ChainRulesCore]] deps = ["Compat", "LinearAlgebra", "SparseArrays"] -git-tree-sha1 = "5d64be50ea9b43a89b476be773e125cef03c7cd5" +git-tree-sha1 = "be770c08881f7bb928dfd86d1ba83798f76cf62a" uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" -version = "0.10.1" +version = "0.10.9" [[Classes]] deps = ["DataStructures", "InteractiveUtils", "MacroTools"] @@ -76,15 +76,13 @@ version = "1.0.2" [[Compat]] deps = ["Base64", "Dates", "DelimitedFiles", "Distributed", "InteractiveUtils", "LibGit2", "Libdl", "LinearAlgebra", "Markdown", "Mmap", "Pkg", "Printf", "REPL", "Random", "SHA", "Serialization", "SharedArrays", "Sockets", "SparseArrays", "Statistics", "Test", "UUIDs", "Unicode"] -git-tree-sha1 = "e4e2b39db08f967cc1360951f01e8a75ec441cab" +git-tree-sha1 = "dc7dedc2c2aa9faf59a55c622760a25cbefbe941" uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" -version = "3.30.0" +version = "3.31.0" [[CompilerSupportLibraries_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "8e695f735fca77e9708e795eda62afdb869cbb70" +deps = ["Artifacts", "Libdl"] uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" -version = "0.3.4+0" [[Compose]] deps = ["Base64", "Colors", "DataStructures", "Dates", "IterTools", "JSON", "LinearAlgebra", "Measures", "Printf", "Random", "Requires", "Statistics", "UUIDs"] @@ -94,9 +92,9 @@ version = "0.9.2" [[ConstructionBase]] deps = ["LinearAlgebra"] -git-tree-sha1 = "1dc43957fb9a1574fa1b7a449e101bd1fd3a9fb7" +git-tree-sha1 = "f74e9d5388b8620b4cee35d4c5a618dd4dc547f4" uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" -version = "1.2.1" +version = "1.3.0" [[Crayons]] git-tree-sha1 = "3f71217b538d7aaee0b69ab47d9b7724ca8afa0d" @@ -104,9 +102,9 @@ uuid = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f" version = "4.0.4" [[DataAPI]] -git-tree-sha1 = "dfb3b7e89e395be1e25c2ad6d7690dc29cc53b1d" +git-tree-sha1 = "ee400abb2298bd13bfc3df1c412ed228061a2385" uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a" -version = "1.6.0" +version = "1.7.0" [[DataFrames]] deps = ["Compat", "DataAPI", "Future", "InvertedIndices", "IteratorInterfaceExtensions", "LinearAlgebra", "Markdown", "Missings", "PooledArrays", "PrettyTables", "Printf", "REPL", "Reexport", "SortingAlgorithms", "Statistics", "TableTraits", "Tables", "Unicode"] @@ -150,10 +148,10 @@ uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" version = "0.24.18" [[DocStringExtensions]] -deps = ["LibGit2", "Markdown", "Pkg", "Test"] -git-tree-sha1 = "9d4f64f79012636741cf01133158a54b24924c32" +deps = ["LibGit2"] +git-tree-sha1 = "a32185f5428d3986f47c2ab78b1f216d5e6cc96f" uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" -version = "0.8.4" +version = "0.8.5" [[DoubleFloats]] deps = ["GenericLinearAlgebra", "GenericSchur", "LinearAlgebra", "Polynomials", "Printf", "Quadmath", "Random", "Requires", "SpecialFunctions"] @@ -161,6 +159,10 @@ git-tree-sha1 = "cfc5657a37c2881a728d76bd14ad808ca096d601" uuid = "497a8b3b-efae-58df-a0af-a86822472b78" version = "1.1.22" +[[Downloads]] +deps = ["ArgTools", "LibCURL", "NetworkOptions"] +uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" + [[Electron]] deps = ["Base64", "FilePaths", "JSON", "Pkg", "Sockets", "URIParser", "UUIDs"] git-tree-sha1 = "a53025d3eabe23659065b3c5bba7b4ffb1327aa0" @@ -179,22 +181,22 @@ uuid = "8f5d6c58-4d21-5cfd-889c-e3ad7ee6a615" version = "1.1.0" [[FFTW]] -deps = ["AbstractFFTs", "FFTW_jll", "IntelOpenMP_jll", "Libdl", "LinearAlgebra", "MKL_jll", "Reexport"] -git-tree-sha1 = "1b48dbde42f307e48685fa9213d8b9f8c0d87594" +deps = ["AbstractFFTs", "FFTW_jll", "LinearAlgebra", "MKL_jll", "Preferences", "Reexport"] +git-tree-sha1 = "f985af3b9f4e278b1d24434cbb546d6092fca661" uuid = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341" -version = "1.3.2" +version = "1.4.3" [[FFTW_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "5a0d4b6a22a34d17d53543bd124f4b08ed78e8b0" +git-tree-sha1 = "3676abafff7e4ff07bbd2c42b3d8201f31653dcc" uuid = "f5851436-0d7a-5f13-b9de-f02708fd171a" -version = "3.3.9+7" +version = "3.3.9+8" [[FileIO]] deps = ["Pkg", "Requires", "UUIDs"] -git-tree-sha1 = "ca36405da56db2f4730b29cd2ad2bf5869baea3c" +git-tree-sha1 = "256d8e6188f3f1ebfa1a5d17e072a0efafa8c5bf" uuid = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549" -version = "1.9.1" +version = "1.10.1" [[FilePaths]] deps = ["FilePathsBase", "MacroTools", "Reexport", "Requires"] @@ -210,9 +212,9 @@ version = "0.9.10" [[FillArrays]] deps = ["LinearAlgebra", "Random", "SparseArrays"] -git-tree-sha1 = "31939159aeb8ffad1d4d8ee44d07f8558273120a" +git-tree-sha1 = "a603e79b71bb3c1efdb58f0ee32286efe2d1a255" uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" -version = "0.11.7" +version = "0.11.8" [[FixedPointNumbers]] deps = ["Statistics"] @@ -255,10 +257,10 @@ uuid = "a2cc645c-3eea-5389-862e-a155d0052231" version = "0.4.4" [[HTTP]] -deps = ["Base64", "Dates", "IniFile", "MbedTLS", "Sockets"] -git-tree-sha1 = "c7ec02c4c6a039a98a15f955462cd7aea5df4508" +deps = ["Base64", "Dates", "IniFile", "MbedTLS", "NetworkOptions", "Sockets", "URIs"] +git-tree-sha1 = "86ed84701fbfd1142c9786f8e53c595ff5a4def9" uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3" -version = "0.8.19" +version = "0.9.10" [[Inflate]] git-tree-sha1 = "f5fc07d4e706b84f72d54eedcc1c13d92fb0871c" @@ -282,10 +284,10 @@ deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" [[Interpolations]] -deps = ["AxisAlgorithms", "LinearAlgebra", "OffsetArrays", "Random", "Ratios", "SharedArrays", "SparseArrays", "StaticArrays", "WoodburyMatrices"] -git-tree-sha1 = "1e0e51692a3a77f1eeb51bf741bdd0439ed210e7" +deps = ["AxisAlgorithms", "ChainRulesCore", "LinearAlgebra", "OffsetArrays", "Random", "Ratios", "Requires", "SharedArrays", "SparseArrays", "StaticArrays", "WoodburyMatrices"] +git-tree-sha1 = "1470c80592cf1f0a35566ee5e93c5f8221ebc33a" uuid = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" -version = "0.13.2" +version = "0.13.3" [[Intervals]] deps = ["Dates", "Printf", "RecipesBase", "Serialization", "TimeZones"] @@ -317,9 +319,9 @@ version = "1.0.0" [[JLD2]] deps = ["DataStructures", "FileIO", "MacroTools", "Mmap", "Pkg", "Printf", "Reexport", "TranscodingStreams", "UUIDs"] -git-tree-sha1 = "8ee0507881f9fbd2d6063c8befe1df4d95966d78" +git-tree-sha1 = "4813826871754cf52607e76ad37acb36ccf52719" uuid = "033835bb-8acc-5ee8-8aae-3f567f8a3819" -version = "0.4.8" +version = "0.4.11" [[JLLWrappers]] deps = ["Preferences"] @@ -346,23 +348,33 @@ uuid = "5ab0869b-81aa-558d-bb23-cbf5423bbe9b" version = "0.6.3" [[LazyArtifacts]] -deps = ["Pkg"] -git-tree-sha1 = "4bb5499a1fc437342ea9ab7e319ede5a457c0968" +deps = ["Artifacts", "Pkg"] uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3" -version = "1.3.0" + +[[LibCURL]] +deps = ["LibCURL_jll", "MozillaCACerts_jll"] +uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" + +[[LibCURL_jll]] +deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"] +uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" [[LibGit2]] -deps = ["Printf"] +deps = ["Base64", "NetworkOptions", "Printf", "SHA"] uuid = "76f85450-5226-5b5a-8eaa-529ad045b433" +[[LibSSH2_jll]] +deps = ["Artifacts", "Libdl", "MbedTLS_jll"] +uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" + [[Libdl]] uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" [[Libiconv_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "8e924324b2e9275a51407a4e06deb3455b1e359f" +git-tree-sha1 = "42b62845d70a619f063a7da093d995ec8e15e778" uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531" -version = "1.16.0+7" +version = "1.16.1+1" [[LightGraphs]] deps = ["ArnoldiMethod", "DataStructures", "Distributed", "Inflate", "LinearAlgebra", "Random", "SharedArrays", "SimpleTraits", "SparseArrays", "Statistics"] @@ -406,10 +418,8 @@ uuid = "739be429-bea8-5141-9913-cc70e7f3736d" version = "1.0.3" [[MbedTLS_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "0eef589dd1c26a3ac9d753fe1a8bcad63f956fa6" +deps = ["Artifacts", "Libdl"] uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" -version = "2.16.8+1" [[Measures]] git-tree-sha1 = "e498ddeee6f9fdb4551ce855a46f54dbd900245f" @@ -424,9 +434,9 @@ version = "0.6.7" [[Mimi]] deps = ["CSVFiles", "Classes", "DataFrames", "DataStructures", "Dates", "DelimitedFiles", "Distributions", "Electron", "FileIO", "FilePaths", "GlobalSensitivityAnalysis", "GraphPlot", "IterTools", "IteratorInterfaceExtensions", "JSON", "LightGraphs", "LinearAlgebra", "Logging", "MacroTools", "MetaGraphs", "NamedArrays", "Pkg", "ProgressMeter", "Random", "Serialization", "Statistics", "StatsBase", "StringBuilders", "TableTraits", "VegaLite"] -git-tree-sha1 = "b41925a57b0d0f59654546ed6fc5463de5596cdc" +git-tree-sha1 = "265766488a3c7c4e45a2dbeaea0c2de9e02208f8" uuid = "e4e893b0-ee5e-52ea-8111-44b3bdec128c" -version = "1.2.3" +version = "1.3.0" [[Missings]] deps = ["DataAPI"] @@ -443,11 +453,14 @@ git-tree-sha1 = "916b850daad0d46b8c71f65f719c49957e9513ed" uuid = "78c3b35d-d492-501b-9361-3d52fe80e533" version = "0.7.1" +[[MozillaCACerts_jll]] +uuid = "14a3606d-f60d-562e-9121-12d972cd8159" + [[MutableArithmetics]] deps = ["LinearAlgebra", "SparseArrays", "Test"] -git-tree-sha1 = "ad9b2bce6021631e0e20706d361972343a03e642" +git-tree-sha1 = "3927848ccebcc165952dc0d9ac9aa274a87bfe01" uuid = "d8a4904e-b15c-11e9-3269-09a3773c0cb0" -version = "0.2.19" +version = "0.2.20" [[NamedArrays]] deps = ["Combinatorics", "DataStructures", "DelimitedFiles", "InvertedIndices", "LinearAlgebra", "Random", "Requires", "SparseArrays", "Statistics"] @@ -455,6 +468,9 @@ git-tree-sha1 = "9ba8ddb0c06a08b1bad81b7120d13288e5d766fa" uuid = "86f7a689-2022-50b4-a561-43c23ac3c673" version = "0.9.5" +[[NetworkOptions]] +uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" + [[NodeJS]] deps = ["Pkg"] git-tree-sha1 = "905224bbdd4b555c69bb964514cfa387616f0d3a" @@ -474,15 +490,15 @@ version = "0.3.3" [[OffsetArrays]] deps = ["Adapt"] -git-tree-sha1 = "1381a7142eefd4cd12f052a4d2d790fe21bd1d55" +git-tree-sha1 = "93e9f0b571c1ddaabebe5a6a52e16be53a97fe25" uuid = "6fe1bfb0-de20-5000-8ca7-80f57d26f881" -version = "1.9.2" +version = "1.10.1" [[OpenSpecFun_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "9db77584158d0ab52307f8c04f8e7c08ca76b5b3" +git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1" uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e" -version = "0.5.3+4" +version = "0.5.5+0" [[OrderedCollections]] git-tree-sha1 = "85f8e6578bf1f9ee0d11e7bb1b1456435479d47c" @@ -491,9 +507,9 @@ version = "1.4.1" [[PDMats]] deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "f82a0e71f222199de8e9eb9a09977bd0767d52a0" +git-tree-sha1 = "4dd403333bcf0909341cfe57ec115152f937d7d8" uuid = "90014a1f-27ba-587c-ab20-58faa44d9150" -version = "0.11.0" +version = "0.11.1" [[Parsers]] deps = ["Dates"] @@ -502,7 +518,7 @@ uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0" version = "1.1.0" [[Pkg]] -deps = ["Dates", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "UUIDs"] +deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" [[Polynomials]] @@ -525,9 +541,9 @@ version = "1.2.2" [[PrettyTables]] deps = ["Crayons", "Formatting", "Markdown", "Reexport", "Tables"] -git-tree-sha1 = "b60494adf99652d220cdef46f8a32232182cc22d" +git-tree-sha1 = "0d1245a357cc61c8cd61934c07447aa569ff22e6" uuid = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" -version = "1.0.1" +version = "1.1.0" [[Printf]] deps = ["Unicode"] @@ -552,7 +568,7 @@ uuid = "be4d8f0f-7fa4-5f49-b795-2f01399ab2dd" version = "0.5.5" [[REPL]] -deps = ["InteractiveUtils", "Markdown", "Sockets"] +deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" [[Random]] @@ -570,9 +586,9 @@ uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" version = "1.1.1" [[Reexport]] -git-tree-sha1 = "57d8440b0c7d98fc4f889e478e80f268d534c9d5" +git-tree-sha1 = "5f6c21241f0f655da3952fd60aa18477cf96c220" uuid = "189a3867-3050-52da-a836-e630ba90ab69" -version = "1.0.0" +version = "1.1.0" [[Requires]] deps = ["UUIDs"] @@ -582,15 +598,15 @@ version = "1.1.3" [[Rmath]] deps = ["Random", "Rmath_jll"] -git-tree-sha1 = "86c5647b565873641538d8f812c04e4c9dbeb370" +git-tree-sha1 = "bf3188feca147ce108c76ad82c2792c57abe7b1f" uuid = "79098fc4-a85e-5d69-aa6a-4863f24498fa" -version = "0.6.1" +version = "0.7.0" [[Rmath_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "1b7bf41258f6c5c9c31df8c1ba34c1fc88674957" +git-tree-sha1 = "68db32dff12bb6127bac73c209881191bf0efbb7" uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f" -version = "0.2.2+2" +version = "0.3.0+0" [[SHA]] uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" @@ -641,9 +657,9 @@ version = "1.5.1" [[StaticArrays]] deps = ["LinearAlgebra", "Random", "Statistics"] -git-tree-sha1 = "42378d3bab8b4f57aa1ca443821b752850592668" +git-tree-sha1 = "745914ebcd610da69f3cb6bf76cb7bb83dcb8c9a" uuid = "90137ffa-7385-5640-81b9-e52037218182" -version = "1.2.2" +version = "1.2.4" [[Statistics]] deps = ["LinearAlgebra", "SparseArrays"] @@ -677,9 +693,7 @@ uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" [[TOML]] deps = ["Dates"] -git-tree-sha1 = "44aaac2d2aec4a850302f9aa69127c74f0c3787e" uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" -version = "1.0.3" [[TableShowUtils]] deps = ["DataValues", "Dates", "JSON", "Markdown", "Test"] @@ -701,12 +715,16 @@ version = "1.0.1" [[Tables]] deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "TableTraits", "Test"] -git-tree-sha1 = "aa30f8bb63f9ff3f8303a06c604c8500a69aa791" +git-tree-sha1 = "8ed4a3ea724dac32670b062be3ef1c1de6773ae8" uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c" -version = "1.4.3" +version = "1.4.4" + +[[Tar]] +deps = ["ArgTools", "SHA"] +uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" [[Test]] -deps = ["Distributed", "InteractiveUtils", "Logging", "Random"] +deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [[TextParse]] @@ -733,6 +751,11 @@ git-tree-sha1 = "53a9f49546b8d2dd2e688d216421d050c9a31d0d" uuid = "30578b45-9adc-5946-b283-645ec420af67" version = "0.4.1" +[[URIs]] +git-tree-sha1 = "97bbe755a53fe859669cd907f2d96aee8d2c1355" +uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4" +version = "1.3.0" + [[UUIDs]] deps = ["Random", "SHA"] uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" @@ -766,9 +789,9 @@ version = "0.5.3" [[XML2_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"] -git-tree-sha1 = "be0db24f70aae7e2b89f2f3092e93b8606d659a6" +git-tree-sha1 = "1acf5bdf07aa0907e0a37d3718bb88d4b687b74a" uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a" -version = "2.9.10+3" +version = "2.9.12+0" [[ZipFile]] deps = ["Libdl", "Printf", "Zlib_jll"] @@ -777,7 +800,13 @@ uuid = "a5390f91-8eb1-5f08-bee0-b1d1ffed6cea" version = "0.9.3" [[Zlib_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "320228915c8debb12cb434c59057290f0834dbf6" +deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" -version = "1.2.11+18" + +[[nghttp2_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" + +[[p7zip_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0" diff --git a/Project.toml b/Project.toml index 54f4559..757c5b0 100644 --- a/Project.toml +++ b/Project.toml @@ -1,4 +1,20 @@ +name = "MimiFAIRv2" +uui = "db6a9ac4-572f-40b4-b774-4040310d200a" +version = "0.0.1-DEV" + [deps] CSVFiles = "5d742f6a-9f54-50ce-8119-2520741973ca" DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0" Mimi = "e4e893b0-ee5e-52ea-8111-44b3bdec128c" + +[compat] +CSVFiles = "1.0" +DataFrames = "1.0" +Mimi = "1.3" +julia = "1.4" + +[extras] +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" + +[targets] +test = ["Test"] diff --git a/README.md b/README.md index a32efa1..f13d0cb 100644 --- a/README.md +++ b/README.md @@ -1,15 +1,17 @@ -# MimiFAIRv2.0.jl +# MimiFAIRv2.jl This is a work-in-progress respository for a Julia-Mimi implementation of the FAIRv2.0 simple climate model. The model description paper can be found at [FaIRv2.0.0: A Generalized Impulse Response Model for Climate Uncertainty and Future Scenario Exploration](https://gmd.copernicus.org/articles/14/3007/2021/gmd-14-3007-2021.html). -## Running Mimi-FAIRv2.0 +## Running Mimi-FAIRv2 To run the model, execute the following code: ```julia + # Load the model code. -include("src/MimiFAIRv2.0.jl") +include("src/MimiFAIRv2.jl") # load the MimiFAIRv2 module +using Main.MimiFAIRv2 # bring the module into your namespace -# Create an instance of Mimi-FAIRv2.0. +# Create an instance of Mimi-FAIRv2. m = get_model() # Run the model. @@ -33,4 +35,4 @@ The `get_model()` function currently has the following keyword arguments: \ \ -![Python vs. Julia temperature comparison](https://github.com/FrankErrickson/MimiFAIRv2.0.jl/blob/main/data/python_replication_data/Python_Mimi_FAIR2_temperature_comparison.png) +![Python vs. Julia temperature comparison](https://github.com/FrankErrickson/MimiFAIRv2.jl/blob/main/data/python_replication_data/Python_Mimi_FAIR2_temperature_comparison.png) \ No newline at end of file diff --git a/src/MimiFAIRv2.0.jl b/src/MimiFAIRv2.0.jl deleted file mode 100644 index 300bdf3..0000000 --- a/src/MimiFAIRv2.0.jl +++ /dev/null @@ -1,291 +0,0 @@ -# Load required packages. -using CSVFiles, DataFrames, Mimi - -# Load helper functions and Mimi-FAIR v2.0 model commponent files. -include("helper_functions.jl") -include("components/co2_cycle.jl") -include("components/ch4_cycle.jl") -include("components/n2o_cycle.jl") -include("components/montreal_gas_cycles.jl") -include("components/flourinated_gas_cycles.jl") -include("components/aerosol_plus_gas_cycles.jl") -include("components/radiative_forcing.jl") -include("components/temperature.jl") - - -function get_model(;emissions_forcing_scenario::String="ssp585", start_year::Int=1750, end_year::Int=2500, TCR::Float64=1.79, RWF::Float64=0.552, F2x::Float64=3.759) - - # --------------------------------------------- - # --------------------------------------------- - # Set Up Data and Parameter Values - # --------------------------------------------- - # --------------------------------------------- - - # Load emissions and forcing data and crop to appropriate model years (scenarios span 1750-2500 by default). - scenario_indices = indexin(start_year:end_year, 1750:2500) - forcing_data = DataFrame(load(joinpath(@__DIR__, "..", "data", "rcmip_"*emissions_forcing_scenario*"_effective_radiative_forcing_1750_to_2500.csv"), skiplines_begin=6))[scenario_indices,:] - emissions_data = DataFrame(load(joinpath(@__DIR__, "..", "data", "rcmip_"*emissions_forcing_scenario*"_emissions_1750_to_2500.csv"), skiplines_begin=6))[scenario_indices,:] - - # Load FAIR default gas cycle (gas) and indirect radiative forcing (irf_p) parameters. - gas_p = DataFrame(load(joinpath(@__DIR__, "..", "data", "default_gas_cycle_parameters.csv"), skiplines_begin=6)) - irf_p = DataFrame(load(joinpath(@__DIR__, "..", "data", "default_indirect_radiative_forcing_parameters.csv"), skiplines_begin=7)) - - # Load initial conditions for 1750 under the RCMIP emissions & forcing scenario. - init_gas_vals = DataFrame(load(joinpath(@__DIR__, "..", "data", "fair_initial_gas_cycle_conditions_1750.csv"), skiplines_begin=7)) - init_thermal_vals = DataFrame(load(joinpath(@__DIR__, "..", "data", "fair_initial_thermal_conditions_1750.csv"), skiplines_begin=7)) - - # Calculate exogenous effective radiative forcing values (sum of land use, solar, and volcanic forcings). - exogenous_forcing = sum(Array(forcing_data[:, [:land_use,:volcanic,:solar]]), dims=2)[:,1] - - # Isolate specific gas parameters for convenience. - co2_p = filter(:gas_name => ==("carbon_dioxide"), gas_p) - ch4_p = filter(:gas_name => ==("methane"), gas_p) - n2o_p = filter(:gas_name => ==("nitrous_oxide"), gas_p) - montreal_gas_p = filter(:gas_group => ==("montreal"), gas_p) - flourinated_gas_p = filter(:gas_group => ==("flourinated"), gas_p) - aerosol_plus_gas_p = filter(:gas_group => ==("aerosol_plus"), gas_p) - - # Sort arrays of parameters for multiple gases so they are listed alphabetically. - sort!(flourinated_gas_p, :gas_name) - sort!(montreal_gas_p, :gas_name) - sort!(aerosol_plus_gas_p, :gas_name) - - #Isolate Montreal indirect forcing parameters and sort alphabetically. - montreal_irf_p = filter(:gas_group => ==("montreal"), irf_p) - sort!(montreal_irf_p, :gas_name) - - # Isolate initial model conditions (initial conditions currently default to year 1750). - co2_init = filter(:gas_name => ==("carbon_dioxide"), init_gas_vals) - ch4_init = filter(:gas_name => ==("methane"), init_gas_vals) - n2o_init = filter(:gas_name => ==("nitrous_oxide"), init_gas_vals) - montreal_init = filter(:gas_group => ==("montreal"), init_gas_vals) - flourinated_init = filter(:gas_group => ==("flourinated"), init_gas_vals) - aerosol_plus_init = filter(:gas_group => ==("aerosol_plus"), init_gas_vals) - - # Sort arrays of initial conditions for multiple gases so they are listed alphabetically. - sort!(montreal_init, :gas_name) - sort!(flourinated_init, :gas_name) - sort!(aerosol_plus_init, :gas_name) - - # Extract emissions arrays for multi-gas groupings. - montreal_emissions = emissions_data[:, Symbol.(montreal_init.gas_name)] - flourinated_emissions = emissions_data[:, Symbol.(flourinated_init.gas_name)] - aerosol_plus_emissions = emissions_data[:, Symbol.(aerosol_plus_init.gas_name)] - - # Create arrays for 'a' parameter groups. - co2_a = vec(Array(co2_p[:,[:a1,:a2,:a3,:a4]])) - ch4_a = vec(Array(ch4_p[:,[:a1,:a2,:a3,:a4]])) - n2o_a = vec(Array(n2o_p[:,[:a1,:a2,:a3,:a4]])) - montreal_a = Array(montreal_gas_p[:,[:a1,:a2,:a3,:a4]]) - flourinated_a = Array(flourinated_gas_p[:,[:a1,:a2,:a3,:a4]]) - aerosol_plus_a = Array(aerosol_plus_gas_p[:,[:a1,:a2,:a3,:a4]]) - - # Create arrays for 'τ' parameter groups. - co2_τ = vec(Array(co2_p[:,[:tau1,:tau2,:tau3,:tau4]])) - ch4_τ = vec(Array(ch4_p[:,[:tau1,:tau2,:tau3,:tau4]])) - n2o_τ = vec(Array(n2o_p[:,[:tau1,:tau2,:tau3,:tau4]])) - montreal_τ = Array(montreal_gas_p[:,[:tau1,:tau2,:tau3,:tau4]]) - flourinated_τ = Array(flourinated_gas_p[:,[:tau1,:tau2,:tau3,:tau4]]) - aerosol_plus_τ = Array(aerosol_plus_gas_p[:,[:tau1,:tau2,:tau3,:tau4]]) - - # Calculate constants to approximate numerical solution for state-dependent timescale adjustment factor from Millar et al. (2017). - g0_co2, g1_co2 = calculate_g0_g1(co2_a, co2_τ) - g0_ch4, g1_ch4 = calculate_g0_g1(ch4_a, ch4_τ) - g0_n2o, g1_n2o = calculate_g0_g1(n2o_a, n2o_τ) - g0_montreal, g1_montreal = calculate_g0_g1(montreal_a, montreal_τ) - g0_flourinated, g1_flourinated = calculate_g0_g1(flourinated_a, flourinated_τ) - g0_aerosol_plus, g1_aerosol_plus = calculate_g0_g1(aerosol_plus_a, aerosol_plus_τ) - - # Calculate default thermal parameter values (defaults: transient climate response=1.79, realized warming fraction=0.552, forcing from a doubling of CO₂=3.759) - thermal_p = get_thermal_parameter_defaults(TCR=1.79, RWF=0.552, F2x=3.759) - - # Calculate thermal decay factors, defined as exp(-1/d). - thermal_decay_factors = exp.(-1.0 ./ vec(Array(thermal_p[:,[:d1,:d2,:d3]]))) - - - # --------------------------------------------- - # --------------------------------------------- - # Initialize Mimi model. - # --------------------------------------------- - # --------------------------------------------- - - # Create a Mimi model. - m = Model() - - # Set time and gas-grouping indices. - set_dimension!(m, :time, start_year:end_year) - set_dimension!(m, :montreal_gases, montreal_gas_p[:,:gas_name]) - set_dimension!(m, :flourinated_gases, flourinated_gas_p[:,:gas_name]) - set_dimension!(m, :aerosol_plus_gases, aerosol_plus_gas_p[:,:gas_name]) - - # --------------------------------------------- - # Add components to model - # --------------------------------------------- - - add_comp!(m, co2_cycle) - add_comp!(m, ch4_cycle) - add_comp!(m, n2o_cycle) - add_comp!(m, montreal_cycles) - add_comp!(m, flourinated_cycles) - add_comp!(m, aerosol_plus_cycles) - add_comp!(m, radiative_forcing) - add_comp!(m, temperature) - - # --------------------------------------------- - # Set component-specific parameters - # --------------------------------------------- - - # ---- Carbon Cycle ---- # - set_param!(m, :co2_cycle, :co2_0, co2_init.concentration[1]) - set_param!(m, :co2_cycle, :r0_co2, co2_p.r0[1]) - set_param!(m, :co2_cycle, :rU_co2, co2_p.rC[1]) - set_param!(m, :co2_cycle, :rT_co2, co2_p.rT[1]) - set_param!(m, :co2_cycle, :rA_co2, co2_p.rA[1]) - set_param!(m, :co2_cycle, :GU_co2_0, co2_init.cumulative_uptake[1]) - set_param!(m, :co2_cycle, :g0_co2, g0_co2) - set_param!(m, :co2_cycle, :g1_co2, g1_co2) - set_param!(m, :co2_cycle, :R0_co2, vec(Array(co2_init[:,[:pool1,:pool2,:pool3,:pool4]]))) - set_param!(m, :co2_cycle, :emiss2conc_co2, co2_p.emis2conc[1]) - set_param!(m, :co2_cycle, :a_co2, co2_a) - set_param!(m, :co2_cycle, :τ_co2, co2_τ) - set_param!(m, :co2_cycle, :E_co2, emissions_data.carbon_dioxide) - - # ---- Methane Cycle ---- # - set_param!(m, :ch4_cycle, :ch4_0, ch4_init.concentration[1]) - set_param!(m, :ch4_cycle, :r0_ch4, ch4_p.r0[1]) - set_param!(m, :ch4_cycle, :rU_ch4, ch4_p.rC[1]) - set_param!(m, :ch4_cycle, :rT_ch4, ch4_p.rT[1]) - set_param!(m, :ch4_cycle, :rA_ch4, ch4_p.rA[1]) - set_param!(m, :ch4_cycle, :GU_ch4_0, ch4_init.cumulative_uptake[1]) - set_param!(m, :ch4_cycle, :g0_ch4, g0_ch4) - set_param!(m, :ch4_cycle, :g1_ch4, g1_ch4) - set_param!(m, :ch4_cycle, :R0_ch4, vec(Array(ch4_init[:,[:pool1,:pool2,:pool3,:pool4]]))) - set_param!(m, :ch4_cycle, :emiss2conc_ch4, ch4_p.emis2conc[1]) - set_param!(m, :ch4_cycle, :a_ch4, ch4_a) - set_param!(m, :ch4_cycle, :τ_ch4, ch4_τ) - set_param!(m, :ch4_cycle, :E_ch4, emissions_data.methane) - - # ---- Nitrous Oxide Cycle ---- # - set_param!(m, :n2o_cycle, :n2o_0, n2o_init.concentration[1]) - set_param!(m, :n2o_cycle, :r0_n2o, n2o_p.r0[1]) - set_param!(m, :n2o_cycle, :rU_n2o, n2o_p.rC[1]) - set_param!(m, :n2o_cycle, :rT_n2o, n2o_p.rT[1]) - set_param!(m, :n2o_cycle, :rA_n2o, n2o_p.rA[1]) - set_param!(m, :n2o_cycle, :GU_n2o_0, n2o_init.cumulative_uptake[1]) - set_param!(m, :n2o_cycle, :g0_n2o, g0_n2o) - set_param!(m, :n2o_cycle, :g1_n2o, g1_n2o) - set_param!(m, :n2o_cycle, :R0_n2o, vec(Array(n2o_init[:,[:pool1,:pool2,:pool3,:pool4]]))) - set_param!(m, :n2o_cycle, :emiss2conc_n2o, n2o_p.emis2conc[1]) - set_param!(m, :n2o_cycle, :a_n2o, n2o_a) - set_param!(m, :n2o_cycle, :τ_n2o, n2o_τ) - set_param!(m, :n2o_cycle, :E_n2o, emissions_data.nitrous_oxide) - - # ---- Flourinated Gas Cycles ---- # - set_param!(m, :flourinated_cycles, :flourinated_0, flourinated_init[:, :concentration]) - set_param!(m, :flourinated_cycles, :r0_flourinated, flourinated_gas_p[:,:r0]) - set_param!(m, :flourinated_cycles, :rU_flourinated, flourinated_gas_p[:,:rC]) - set_param!(m, :flourinated_cycles, :rT_flourinated, flourinated_gas_p[:,:rT]) - set_param!(m, :flourinated_cycles, :rA_flourinated, flourinated_gas_p[:,:rA]) - set_param!(m, :flourinated_cycles, :GU_flourinated_0, flourinated_init[:, :cumulative_uptake]) - set_param!(m, :flourinated_cycles, :g0_flourinated, g0_flourinated) - set_param!(m, :flourinated_cycles, :g1_flourinated, g1_flourinated) - set_param!(m, :flourinated_cycles, :R0_flourinated, Array(flourinated_init[:,[:pool1,:pool2,:pool3,:pool4]])) - set_param!(m, :flourinated_cycles, :emiss2conc_flourinated, flourinated_gas_p[:,:emis2conc]) - set_param!(m, :flourinated_cycles, :a_flourinated, flourinated_a) - set_param!(m, :flourinated_cycles, :τ_flourinated, flourinated_τ) - set_param!(m, :flourinated_cycles, :E_flourinated, Array(flourinated_emissions)) - - # ---- Montreal Protocol Gas Cycles ---- # - set_param!(m, :montreal_cycles, :montreal_0, montreal_init[:, :concentration]) - set_param!(m, :montreal_cycles, :r0_montreal, montreal_gas_p[:,:r0]) - set_param!(m, :montreal_cycles, :rU_montreal, montreal_gas_p[:,:rC]) - set_param!(m, :montreal_cycles, :rT_montreal, montreal_gas_p[:,:rT]) - set_param!(m, :montreal_cycles, :rA_montreal, montreal_gas_p[:,:rA]) - set_param!(m, :montreal_cycles, :GU_montreal_0, montreal_init[:, :cumulative_uptake]) - set_param!(m, :montreal_cycles, :g0_montreal, g0_montreal) - set_param!(m, :montreal_cycles, :g1_montreal, g1_montreal) - set_param!(m, :montreal_cycles, :R0_montreal, Array(montreal_init[:,[:pool1,:pool2,:pool3,:pool4]])) - set_param!(m, :montreal_cycles, :emiss2conc_montreal, montreal_gas_p[:,:emis2conc]) - set_param!(m, :montreal_cycles, :a_montreal, montreal_a) - set_param!(m, :montreal_cycles, :τ_montreal, montreal_τ) - set_param!(m, :montreal_cycles, :E_montreal, Array(montreal_emissions)) - - # ---- Tropospheric Ozone Precursors, Aerosols, & Reactive Gas Cycles (Aerosol+) ---- # - set_param!(m, :aerosol_plus_cycles, :aerosol_plus_0, aerosol_plus_init[:, :concentration]) - set_param!(m, :aerosol_plus_cycles, :r0_aerosol_plus, aerosol_plus_gas_p[:,:r0]) - set_param!(m, :aerosol_plus_cycles, :rU_aerosol_plus, aerosol_plus_gas_p[:,:rC]) - set_param!(m, :aerosol_plus_cycles, :rT_aerosol_plus, aerosol_plus_gas_p[:,:rT]) - set_param!(m, :aerosol_plus_cycles, :rA_aerosol_plus, aerosol_plus_gas_p[:,:rA]) - set_param!(m, :aerosol_plus_cycles, :GU_aerosol_plus_0, aerosol_plus_init[:, :cumulative_uptake]) - set_param!(m, :aerosol_plus_cycles, :g0_aerosol_plus, g0_aerosol_plus) - set_param!(m, :aerosol_plus_cycles, :g1_aerosol_plus, g1_aerosol_plus) - set_param!(m, :aerosol_plus_cycles, :R0_aerosol_plus, Array(aerosol_plus_init[:,[:pool1,:pool2,:pool3,:pool4]])) - set_param!(m, :aerosol_plus_cycles, :emiss2conc_aerosol_plus, aerosol_plus_gas_p[:,:emis2conc]) - set_param!(m, :aerosol_plus_cycles, :a_aerosol_plus, aerosol_plus_a) - set_param!(m, :aerosol_plus_cycles, :τ_aerosol_plus, aerosol_plus_τ) - set_param!(m, :aerosol_plus_cycles, :E_aerosol_plus, Array(aerosol_plus_emissions)) - - # ---- Radiative Forcing ---- # - set_param!(m, :radiative_forcing, :co2_f, vec(Array(co2_p[:, [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :ch4_f, vec(Array(ch4_p[:, [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :ch4_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="methane|o3"), [:f1, :f2, :f3]])) ) - set_param!(m, :radiative_forcing, :ch4_h2o_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="methane|strat_h2o"), [:f1, :f2, :f3]])) ) - set_param!(m, :radiative_forcing, :n2o_f, vec(Array(n2o_p[:, [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :n2o_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nitrous_oxide|o3"), [:f1, :f2, :f3]])) ) - set_param!(m, :radiative_forcing, :montreal_f, Array(montreal_gas_p[:,[:f1, :f2, :f3]])) - set_param!(m, :radiative_forcing, :montreal_ind_f, Array(montreal_irf_p[:,[:f1, :f2, :f3]])) - set_param!(m, :radiative_forcing, :flourinated_f, Array(flourinated_gas_p[:,[:f1, :f2, :f3]])) - set_param!(m, :radiative_forcing, :bc_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="bc"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :bc_snow_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="bc|bc_on_snow"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :bc_aci_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="bc|aci"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :co_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="co"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :co_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="co|o3"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nh3_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nh3"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nmvoc_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nmvoc"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nmvoc_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nmvoc|o3"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nox_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nox"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nox_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nox|o3"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nox_avi_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nox_avi"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :nox_avi_contrails_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nox_avi|contrails"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :oc_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="oc"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :oc_aci_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="oc|aci"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :so2_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="so2"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :so2_aci_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="so2|aci"), [:f1, :f2, :f3]]))) - set_param!(m, :radiative_forcing, :exogenous_RF, exogenous_forcing) - - # ---- Global Temperature Anomaly ---- # - set_param!(m, :temperature, :Tj_0, vec(Array(init_thermal_vals[:,[:thermal_box_1,:thermal_box_2,:thermal_box_3]]))) - set_param!(m, :temperature, :T_0, init_thermal_vals.global_temp_anomaly[1]) - set_param!(m, :temperature, :q, vec(Array(thermal_p[:,[:q1,:q2,:q3]]))) - set_param!(m, :temperature, :decay_factor, thermal_decay_factors) - - # --- Set Parameters Common to Multiple Components ---- # - set_param!(m, :co2_pi, co2_p.PI_conc[1]) - set_param!(m, :ch4_pi, ch4_p.PI_conc[1]) - set_param!(m, :n2o_pi, n2o_p.PI_conc[1]) - set_param!(m, :montreal_pi, montreal_gas_p[:,:PI_conc]) - set_param!(m, :flourinated_pi, flourinated_gas_p[:,:PI_conc]) - set_param!(m, :aerosol_plus_pi, aerosol_plus_gas_p[:,:PI_conc]) - - # --------------------------------------------- - # Create Connections Between Mimi Components - # --------------------------------------------- - - # Syntax is :component_needing_a_parameter_input => :name_of_that_parameter, :component_calculating_required_values => :name_of_variable_output - connect_param!(m, :montreal_cycles => :Tj, :temperature => :Tj) - connect_param!(m, :flourinated_cycles => :Tj, :temperature => :Tj) - connect_param!(m, :aerosol_plus_cycles => :Tj, :temperature => :Tj) - connect_param!(m, :co2_cycle => :Tj, :temperature => :Tj) - connect_param!(m, :ch4_cycle => :Tj, :temperature => :Tj) - connect_param!(m, :n2o_cycle => :Tj, :temperature => :Tj) - connect_param!(m, :radiative_forcing => :co2_conc, :co2_cycle => :co2) - connect_param!(m, :radiative_forcing => :ch4_conc, :ch4_cycle => :ch4) - connect_param!(m, :radiative_forcing => :n2o_conc, :n2o_cycle => :n2o) - connect_param!(m, :radiative_forcing => :montreal_conc, :montreal_cycles => :montreal_conc) - connect_param!(m, :radiative_forcing => :flourinated_conc, :flourinated_cycles => :flourinated_conc) - connect_param!(m, :radiative_forcing => :aerosol_plus_conc, :aerosol_plus_cycles => :aerosol_plus_conc) - connect_param!(m, :temperature => :F, :radiative_forcing => :total_RF) - - # Return model. - return m -end diff --git a/src/MimiFAIRv2.jl b/src/MimiFAIRv2.jl new file mode 100644 index 0000000..de77963 --- /dev/null +++ b/src/MimiFAIRv2.jl @@ -0,0 +1,336 @@ +module MimiFAIRv2 + +# Load required packages. +using CSVFiles, DataFrames, Mimi + +# Load helper functions and MimiFAIRv2 model commponent files. +include("helper_functions.jl") +include("components/co2_cycle.jl") +include("components/ch4_cycle.jl") +include("components/n2o_cycle.jl") +include("components/montreal_gas_cycles.jl") +include("components/flourinated_gas_cycles.jl") +include("components/aerosol_plus_gas_cycles.jl") +include("components/radiative_forcing.jl") +include("components/temperature.jl") + +export get_model + +""" + get_model(;emissions_forcing_scenario::String="ssp585", start_year::Int=1750, + end_year::Int=2500, TCR::Float64=1.79, RWF::Float64=0.552, + F2x::Float64=3.759) + +Return a constructed model with default settings for emissions forcing scenario +(default to ssp585), start year (default to 1750), end year (default to 2500), +transient climate response (default to 1.79), realized warming fraction (default +to 5.22), and forcing from a doubling of CO₂ (default to 3.759). +""" +function get_model(;emissions_forcing_scenario::String="ssp585", start_year::Int=1750, end_year::Int=2500, TCR::Float64=1.79, RWF::Float64=0.552, F2x::Float64=3.759) + + # TODO - should this warn or error? + if start_year !== 1750 + @warn "Model should not be set to start with a year differing from 1750." + end + + # --------------------------------------------- + # --------------------------------------------- + # Set Up Data and Parameter Values + # --------------------------------------------- + # --------------------------------------------- + + # Load emissions and forcing data and crop to appropriate model years (scenarios span 1750-2500 by default). + scenario_indices = indexin(start_year:end_year, 1750:2500) + forcing_data = DataFrame(load(joinpath(@__DIR__, "..", "data", "rcmip_"*emissions_forcing_scenario*"_effective_radiative_forcing_1750_to_2500.csv"), skiplines_begin=6))[scenario_indices,:] + emissions_data = DataFrame(load(joinpath(@__DIR__, "..", "data", "rcmip_"*emissions_forcing_scenario*"_emissions_1750_to_2500.csv"), skiplines_begin=6))[scenario_indices,:] + + # Load FAIR default gas cycle (gas) and indirect radiative forcing (irf_p) parameters. + gas_p = DataFrame(load(joinpath(@__DIR__, "..", "data", "default_gas_cycle_parameters.csv"), skiplines_begin=6)) + irf_p = DataFrame(load(joinpath(@__DIR__, "..", "data", "default_indirect_radiative_forcing_parameters.csv"), skiplines_begin=7)) + + # Load initial conditions for 1750 under the RCMIP emissions & forcing scenario. + init_gas_vals = DataFrame(load(joinpath(@__DIR__, "..", "data", "fair_initial_gas_cycle_conditions_1750.csv"), skiplines_begin=7)) + init_thermal_vals = DataFrame(load(joinpath(@__DIR__, "..", "data", "fair_initial_thermal_conditions_1750.csv"), skiplines_begin=7)) + + # Calculate exogenous effective radiative forcing values (sum of land use, solar, and volcanic forcings). + exogenous_forcing = sum(Array(forcing_data[:, [:land_use,:volcanic,:solar]]), dims=2)[:,1] + + # Isolate specific gas parameters for convenience. + co2_p = filter(:gas_name => ==("carbon_dioxide"), gas_p) + ch4_p = filter(:gas_name => ==("methane"), gas_p) + n2o_p = filter(:gas_name => ==("nitrous_oxide"), gas_p) + montreal_gas_p = filter(:gas_group => ==("montreal"), gas_p) + flourinated_gas_p = filter(:gas_group => ==("flourinated"), gas_p) + aerosol_plus_gas_p = filter(:gas_group => ==("aerosol_plus"), gas_p) + + # Sort arrays of parameters for multiple gases so they are listed alphabetically. + sort!(flourinated_gas_p, :gas_name) + sort!(montreal_gas_p, :gas_name) + sort!(aerosol_plus_gas_p, :gas_name) + + #Isolate Montreal indirect forcing parameters and sort alphabetically. + montreal_irf_p = filter(:gas_group => ==("montreal"), irf_p) + sort!(montreal_irf_p, :gas_name) + + # Isolate initial model conditions (initial conditions currently default to year 1750). + co2_init = filter(:gas_name => ==("carbon_dioxide"), init_gas_vals) + ch4_init = filter(:gas_name => ==("methane"), init_gas_vals) + n2o_init = filter(:gas_name => ==("nitrous_oxide"), init_gas_vals) + montreal_init = filter(:gas_group => ==("montreal"), init_gas_vals) + flourinated_init = filter(:gas_group => ==("flourinated"), init_gas_vals) + aerosol_plus_init = filter(:gas_group => ==("aerosol_plus"), init_gas_vals) + + # Sort arrays of initial conditions for multiple gases so they are listed alphabetically. + sort!(montreal_init, :gas_name) + sort!(flourinated_init, :gas_name) + sort!(aerosol_plus_init, :gas_name) + + # Extract emissions arrays for multi-gas groupings. + montreal_emissions = emissions_data[:, Symbol.(montreal_init.gas_name)] + flourinated_emissions = emissions_data[:, Symbol.(flourinated_init.gas_name)] + aerosol_plus_emissions = emissions_data[:, Symbol.(aerosol_plus_init.gas_name)] + + # Create helper arrays of indices to use for pulling parameter groups. + a_idxs = [:a1,:a2,:a3,:a4] + τ_idxs = [:tau1,:tau2,:tau3,:tau4] + + # Create arrays for 'a' parameter groups. + co2_a = vec(Array(co2_p[:, a_idxs])) + ch4_a = vec(Array(ch4_p[:, a_idxs])) + n2o_a = vec(Array(n2o_p[:, a_idxs])) + montreal_a = Array(montreal_gas_p[:, a_idxs]) + flourinated_a = Array(flourinated_gas_p[:, a_idxs]) + aerosol_plus_a = Array(aerosol_plus_gas_p[:, a_idxs]) + + # Create arrays for 'τ' parameter groups. + co2_τ = vec(Array(co2_p[:, τ_idxs])) + ch4_τ = vec(Array(ch4_p[:, τ_idxs])) + n2o_τ = vec(Array(n2o_p[:, τ_idxs])) + montreal_τ = Array(montreal_gas_p[:, τ_idxs]) + flourinated_τ = Array(flourinated_gas_p[:, τ_idxs]) + aerosol_plus_τ = Array(aerosol_plus_gas_p[:, τ_idxs]) + + # Calculate constants to approximate numerical solution for state-dependent timescale adjustment factor from Millar et al. (2017). + g0_co2, g1_co2 = calculate_g0_g1(co2_a, co2_τ) + g0_ch4, g1_ch4 = calculate_g0_g1(ch4_a, ch4_τ) + g0_n2o, g1_n2o = calculate_g0_g1(n2o_a, n2o_τ) + g0_montreal, g1_montreal = calculate_g0_g1(montreal_a, montreal_τ) + g0_flourinated, g1_flourinated = calculate_g0_g1(flourinated_a, flourinated_τ) + g0_aerosol_plus, g1_aerosol_plus = calculate_g0_g1(aerosol_plus_a, aerosol_plus_τ) + + # Calculate default thermal parameter values + # Defaults pulled from get_model function defaults are are as follows: + # transient climate response = 1.79 + # realized warming fraction = 0.552 + # forcing from a doubling of CO₂ = 3.759 + thermal_p = get_thermal_parameter_defaults(TCR = TCR, RWF= RWF, F2x = F2x) + + # Calculate thermal decay factors, defined as exp(-1/d). + thermal_decay_factors = exp.(-1.0 ./ vec(Array(thermal_p[:,[:d1,:d2,:d3]]))) + + + # --------------------------------------------- + # --------------------------------------------- + # Initialize Mimi model. + # --------------------------------------------- + # --------------------------------------------- + + # Create a Mimi model. + m = Model() + + # Set time and gas-grouping indices. + set_dimension!(m, :time, start_year:end_year) + set_dimension!(m, :montreal_gases, montreal_gas_p[:,:gas_name]) + set_dimension!(m, :flourinated_gases, flourinated_gas_p[:,:gas_name]) + set_dimension!(m, :aerosol_plus_gases, aerosol_plus_gas_p[:,:gas_name]) + + # --------------------------------------------- + # Add components to model + # --------------------------------------------- + + add_comp!(m, co2_cycle) + add_comp!(m, ch4_cycle) + add_comp!(m, n2o_cycle) + add_comp!(m, montreal_cycles) + add_comp!(m, flourinated_cycles) + add_comp!(m, aerosol_plus_cycles) + add_comp!(m, radiative_forcing) + add_comp!(m, temperature) + + # --------------------------------------------- + # Set component-specific parameters + # --------------------------------------------- + + # ---- Carbon Cycle ---- # + update_param!(m, :co2_cycle, :co2_0, co2_init.concentration[1]) + update_param!(m, :co2_cycle, :r0_co2, co2_p.r0[1]) + update_param!(m, :co2_cycle, :rU_co2, co2_p.rC[1]) + update_param!(m, :co2_cycle, :rT_co2, co2_p.rT[1]) + update_param!(m, :co2_cycle, :rA_co2, co2_p.rA[1]) + update_param!(m, :co2_cycle, :GU_co2_0, co2_init.cumulative_uptake[1]) + update_param!(m, :co2_cycle, :g0_co2, g0_co2) + update_param!(m, :co2_cycle, :g1_co2, g1_co2) + update_param!(m, :co2_cycle, :R0_co2, vec(Array(co2_init[:,[:pool1,:pool2,:pool3,:pool4]]))) + update_param!(m, :co2_cycle, :emiss2conc_co2, co2_p.emis2conc[1]) + update_param!(m, :co2_cycle, :a_co2, co2_a) + update_param!(m, :co2_cycle, :τ_co2, co2_τ) + update_param!(m, :co2_cycle, :E_co2, emissions_data.carbon_dioxide) + + # ---- Methane Cycle ---- # + update_param!(m, :ch4_cycle, :ch4_0, ch4_init.concentration[1]) + update_param!(m, :ch4_cycle, :r0_ch4, ch4_p.r0[1]) + update_param!(m, :ch4_cycle, :rU_ch4, ch4_p.rC[1]) + update_param!(m, :ch4_cycle, :rT_ch4, ch4_p.rT[1]) + update_param!(m, :ch4_cycle, :rA_ch4, ch4_p.rA[1]) + update_param!(m, :ch4_cycle, :GU_ch4_0, ch4_init.cumulative_uptake[1]) + update_param!(m, :ch4_cycle, :g0_ch4, g0_ch4) + update_param!(m, :ch4_cycle, :g1_ch4, g1_ch4) + update_param!(m, :ch4_cycle, :R0_ch4, vec(Array(ch4_init[:,[:pool1,:pool2,:pool3,:pool4]]))) + update_param!(m, :ch4_cycle, :emiss2conc_ch4, ch4_p.emis2conc[1]) + update_param!(m, :ch4_cycle, :a_ch4, ch4_a) + update_param!(m, :ch4_cycle, :τ_ch4, ch4_τ) + update_param!(m, :ch4_cycle, :E_ch4, emissions_data.methane) + + # ---- Nitrous Oxide Cycle ---- # + update_param!(m, :n2o_cycle, :n2o_0, n2o_init.concentration[1]) + update_param!(m, :n2o_cycle, :r0_n2o, n2o_p.r0[1]) + update_param!(m, :n2o_cycle, :rU_n2o, n2o_p.rC[1]) + update_param!(m, :n2o_cycle, :rT_n2o, n2o_p.rT[1]) + update_param!(m, :n2o_cycle, :rA_n2o, n2o_p.rA[1]) + update_param!(m, :n2o_cycle, :GU_n2o_0, n2o_init.cumulative_uptake[1]) + update_param!(m, :n2o_cycle, :g0_n2o, g0_n2o) + update_param!(m, :n2o_cycle, :g1_n2o, g1_n2o) + update_param!(m, :n2o_cycle, :R0_n2o, vec(Array(n2o_init[:,[:pool1,:pool2,:pool3,:pool4]]))) + update_param!(m, :n2o_cycle, :emiss2conc_n2o, n2o_p.emis2conc[1]) + update_param!(m, :n2o_cycle, :a_n2o, n2o_a) + update_param!(m, :n2o_cycle, :τ_n2o, n2o_τ) + update_param!(m, :n2o_cycle, :E_n2o, emissions_data.nitrous_oxide) + + # ---- Flourinated Gas Cycles ---- # + update_param!(m, :flourinated_cycles, :flourinated_0, flourinated_init[:, :concentration]) + update_param!(m, :flourinated_cycles, :r0_flourinated, flourinated_gas_p[:,:r0]) + update_param!(m, :flourinated_cycles, :rU_flourinated, flourinated_gas_p[:,:rC]) + update_param!(m, :flourinated_cycles, :rT_flourinated, flourinated_gas_p[:,:rT]) + update_param!(m, :flourinated_cycles, :rA_flourinated, flourinated_gas_p[:,:rA]) + update_param!(m, :flourinated_cycles, :GU_flourinated_0, flourinated_init[:, :cumulative_uptake]) + update_param!(m, :flourinated_cycles, :g0_flourinated, g0_flourinated) + update_param!(m, :flourinated_cycles, :g1_flourinated, g1_flourinated) + update_param!(m, :flourinated_cycles, :R0_flourinated, Array(flourinated_init[:,[:pool1,:pool2,:pool3,:pool4]])) + update_param!(m, :flourinated_cycles, :emiss2conc_flourinated, flourinated_gas_p[:,:emis2conc]) + update_param!(m, :flourinated_cycles, :a_flourinated, flourinated_a) + update_param!(m, :flourinated_cycles, :τ_flourinated, flourinated_τ) + update_param!(m, :flourinated_cycles, :E_flourinated, Array(flourinated_emissions)) + + # ---- Montreal Protocol Gas Cycles ---- # + update_param!(m, :montreal_cycles, :montreal_0, montreal_init[:, :concentration]) + update_param!(m, :montreal_cycles, :r0_montreal, montreal_gas_p[:,:r0]) + update_param!(m, :montreal_cycles, :rU_montreal, montreal_gas_p[:,:rC]) + update_param!(m, :montreal_cycles, :rT_montreal, montreal_gas_p[:,:rT]) + update_param!(m, :montreal_cycles, :rA_montreal, montreal_gas_p[:,:rA]) + update_param!(m, :montreal_cycles, :GU_montreal_0, montreal_init[:, :cumulative_uptake]) + update_param!(m, :montreal_cycles, :g0_montreal, g0_montreal) + update_param!(m, :montreal_cycles, :g1_montreal, g1_montreal) + update_param!(m, :montreal_cycles, :R0_montreal, Array(montreal_init[:,[:pool1,:pool2,:pool3,:pool4]])) + update_param!(m, :montreal_cycles, :emiss2conc_montreal, montreal_gas_p[:,:emis2conc]) + update_param!(m, :montreal_cycles, :a_montreal, montreal_a) + update_param!(m, :montreal_cycles, :τ_montreal, montreal_τ) + update_param!(m, :montreal_cycles, :E_montreal, Array(montreal_emissions)) + + # ---- Tropospheric Ozone Precursors, Aerosols, & Reactive Gas Cycles (Aerosol+) ---- # + update_param!(m, :aerosol_plus_cycles, :aerosol_plus_0, aerosol_plus_init[:, :concentration]) + update_param!(m, :aerosol_plus_cycles, :r0_aerosol_plus, aerosol_plus_gas_p[:,:r0]) + update_param!(m, :aerosol_plus_cycles, :rU_aerosol_plus, aerosol_plus_gas_p[:,:rC]) + update_param!(m, :aerosol_plus_cycles, :rT_aerosol_plus, aerosol_plus_gas_p[:,:rT]) + update_param!(m, :aerosol_plus_cycles, :rA_aerosol_plus, aerosol_plus_gas_p[:,:rA]) + update_param!(m, :aerosol_plus_cycles, :GU_aerosol_plus_0, aerosol_plus_init[:, :cumulative_uptake]) + update_param!(m, :aerosol_plus_cycles, :g0_aerosol_plus, g0_aerosol_plus) + update_param!(m, :aerosol_plus_cycles, :g1_aerosol_plus, g1_aerosol_plus) + update_param!(m, :aerosol_plus_cycles, :R0_aerosol_plus, Array(aerosol_plus_init[:,[:pool1,:pool2,:pool3,:pool4]])) + update_param!(m, :aerosol_plus_cycles, :emiss2conc_aerosol_plus, aerosol_plus_gas_p[:,:emis2conc]) + update_param!(m, :aerosol_plus_cycles, :a_aerosol_plus, aerosol_plus_a) + update_param!(m, :aerosol_plus_cycles, :τ_aerosol_plus, aerosol_plus_τ) + update_param!(m, :aerosol_plus_cycles, :E_aerosol_plus, Array(aerosol_plus_emissions)) + + # ---- Radiative Forcing ---- # + update_param!(m, :radiative_forcing, :co2_f, vec(Array(co2_p[:, [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :ch4_f, vec(Array(ch4_p[:, [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :ch4_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="methane|o3"), [:f1, :f2, :f3]])) ) + update_param!(m, :radiative_forcing, :ch4_h2o_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="methane|strat_h2o"), [:f1, :f2, :f3]])) ) + update_param!(m, :radiative_forcing, :n2o_f, vec(Array(n2o_p[:, [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :n2o_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nitrous_oxide|o3"), [:f1, :f2, :f3]])) ) + update_param!(m, :radiative_forcing, :montreal_f, Array(montreal_gas_p[:,[:f1, :f2, :f3]])) + update_param!(m, :radiative_forcing, :montreal_ind_f, Array(montreal_irf_p[:,[:f1, :f2, :f3]])) + update_param!(m, :radiative_forcing, :flourinated_f, Array(flourinated_gas_p[:,[:f1, :f2, :f3]])) + update_param!(m, :radiative_forcing, :bc_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="bc"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :bc_snow_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="bc|bc_on_snow"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :bc_aci_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="bc|aci"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :co_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="co"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :co_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="co|o3"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nh3_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nh3"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nmvoc_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nmvoc"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nmvoc_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nmvoc|o3"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nox_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nox"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nox_o3_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nox|o3"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nox_avi_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="nox_avi"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :nox_avi_contrails_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="nox_avi|contrails"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :oc_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="oc"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :oc_aci_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="oc|aci"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :so2_f, vec(Array(aerosol_plus_gas_p[(aerosol_plus_gas_p.gas_name.=="so2"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :so2_aci_f, vec(Array(irf_p[(irf_p.indirect_forcing_effect.=="so2|aci"), [:f1, :f2, :f3]]))) + update_param!(m, :radiative_forcing, :exogenous_RF, exogenous_forcing) + + # ---- Global Temperature Anomaly ---- # + update_param!(m, :temperature, :Tj_0, vec(Array(init_thermal_vals[:,[:thermal_box_1,:thermal_box_2,:thermal_box_3]]))) + update_param!(m, :temperature, :T_0, init_thermal_vals.global_temp_anomaly[1]) + update_param!(m, :temperature, :q, vec(Array(thermal_p[:,[:q1,:q2,:q3]]))) + update_param!(m, :temperature, :decay_factor, thermal_decay_factors) + + # --- Set Parameters Common to Multiple Components ---- # + add_shared_param!(m, :shared_co2_pi, co2_p.PI_conc[1]) + connect_param!(m, :co2_cycle, :co2_pi, :shared_co2_pi) + connect_param!(m, :radiative_forcing, :co2_pi, :shared_co2_pi) + + add_shared_param!(m, :shared_ch4_pi, ch4_p.PI_conc[1]) + connect_param!(m, :ch4_cycle, :ch4_pi, :shared_ch4_pi) + connect_param!(m, :radiative_forcing, :ch4_pi, :shared_ch4_pi) + + add_shared_param!(m, :shared_n2o_pi, n2o_p.PI_conc[1]) + connect_param!(m, :n2o_cycle, :n2o_pi, :shared_n2o_pi) + connect_param!(m, :radiative_forcing, :n2o_pi, :shared_n2o_pi) + + add_shared_param!(m, :shared_montreal_pi, montreal_gas_p[:,:PI_conc], dims = [:montreal_gases]) + connect_param!(m, :montreal_cycles, :montreal_pi, :shared_montreal_pi) + connect_param!(m, :radiative_forcing, :montreal_pi, :shared_montreal_pi) + + add_shared_param!(m, :shared_flourinated_pi, flourinated_gas_p[:,:PI_conc], dims = [:flourinated_gases]) + connect_param!(m, :flourinated_cycles, :flourinated_pi, :shared_flourinated_pi) + connect_param!(m, :radiative_forcing, :flourinated_pi, :shared_flourinated_pi) + + add_shared_param!(m, :shared_aerosol_plus_pi, aerosol_plus_gas_p[:,:PI_conc], dims = [:aerosol_plus_gases]) + connect_param!(m, :aerosol_plus_cycles, :aerosol_plus_pi, :shared_aerosol_plus_pi) + connect_param!(m, :radiative_forcing, :aerosol_plus_pi, :shared_aerosol_plus_pi) + + # --------------------------------------------- + # Create Connections Between Mimi Components + # --------------------------------------------- + + # Syntax is :component_needing_a_parameter_input => :name_of_that_parameter, :component_calculating_required_values => :name_of_variable_output + connect_param!(m, :montreal_cycles => :Tj, :temperature => :Tj) + connect_param!(m, :flourinated_cycles => :Tj, :temperature => :Tj) + connect_param!(m, :aerosol_plus_cycles => :Tj, :temperature => :Tj) + connect_param!(m, :co2_cycle => :Tj, :temperature => :Tj) + connect_param!(m, :ch4_cycle => :Tj, :temperature => :Tj) + connect_param!(m, :n2o_cycle => :Tj, :temperature => :Tj) + connect_param!(m, :radiative_forcing => :co2_conc, :co2_cycle => :co2) + connect_param!(m, :radiative_forcing => :ch4_conc, :ch4_cycle => :ch4) + connect_param!(m, :radiative_forcing => :n2o_conc, :n2o_cycle => :n2o) + connect_param!(m, :radiative_forcing => :montreal_conc, :montreal_cycles => :montreal_conc) + connect_param!(m, :radiative_forcing => :flourinated_conc, :flourinated_cycles => :flourinated_conc) + connect_param!(m, :radiative_forcing => :aerosol_plus_conc, :aerosol_plus_cycles => :aerosol_plus_conc) + connect_param!(m, :temperature => :F, :radiative_forcing => :total_RF) + + # Return model. + return m +end + +end # module diff --git a/src/components/aerosol_plus_gas_cycles.jl b/src/components/aerosol_plus_gas_cycles.jl index 1505d26..e5ff6d7 100644 --- a/src/components/aerosol_plus_gas_cycles.jl +++ b/src/components/aerosol_plus_gas_cycles.jl @@ -4,31 +4,31 @@ @defcomp aerosol_plus_cycles begin - aerosol_plus_gases = Index() # Index for tropospheric ozone precursors, aerosols and reactive gas - - aerosol_plus_0 = Parameter(index=[aerosol_plus_gases]) # Aerosol+ gas concentration in initial model period (ppb). - aerosol_plus_pi = Parameter(index=[aerosol_plus_gases]) # Pre-industrial Aerosol+ concentrations (ppb). - g0_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 Aerosol+ gas cycle parameterization at α=1. - g1_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 Aerosol+ gas cycle parameterization at α=1. - emiss2conc_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Conversion factor between emissions and concentrations. - GU_aerosol_plus_0 = Parameter(index=[aerosol_plus_gases]) # Initial model period value for cumulative uptake of agent (unit of E⁻¹). - r0_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Strength of pre-industrial uptake from atmosphere. - rA_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). - rT_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). - rU_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). - τ_aerosol_plus = Parameter(index=[aerosol_plus_gases, 4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). - a_aerosol_plus = Parameter(index=[aerosol_plus_gases, 4]) # Fraction of emissions entering iᵗʰ Aerosol+ gas reservior. - R0_aerosol_plus = Parameter(index=[aerosol_plus_gases, 4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). + aerosol_plus_gases = Index() # Index for tropospheric ozone precursors, aerosols and reactive gas + + aerosol_plus_0 = Parameter(index=[aerosol_plus_gases]) # Aerosol+ gas concentration in initial model period (ppb). + aerosol_plus_pi = Parameter(index=[aerosol_plus_gases]) # Pre-industrial Aerosol+ concentrations (ppb). + g0_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 Aerosol+ gas cycle parameterization at α=1. + g1_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 Aerosol+ gas cycle parameterization at α=1. + emiss2conc_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Conversion factor between emissions and concentrations. + GU_aerosol_plus_0 = Parameter(index=[aerosol_plus_gases]) # Initial model period value for cumulative uptake of agent (unit of E⁻¹). + r0_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Strength of pre-industrial uptake from atmosphere. + rA_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). + rT_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). + rU_aerosol_plus = Parameter(index=[aerosol_plus_gases]) # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). + τ_aerosol_plus = Parameter(index=[aerosol_plus_gases, 4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). + a_aerosol_plus = Parameter(index=[aerosol_plus_gases, 4]) # Fraction of emissions entering iᵗʰ Aerosol+ gas reservior. + R0_aerosol_plus = Parameter(index=[aerosol_plus_gases, 4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). E_aerosol_plus = Parameter(index=[time, aerosol_plus_gases]) # Annual Aerosol+ emissions (Tg yr⁻¹). - Tj = Parameter(index=[time, 3]) # Temperature change for three thermal pools (K). + Tj = Parameter(index=[time, 3]) # Temperature change for three thermal pools (K). α_aerosol_plus = Variable(index=[time, aerosol_plus_gases]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. aerosol_plus_conc = Variable(index=[time, aerosol_plus_gases]) # Total atmospheric Aerosol+ concentrations (ppb). GA_aerosol_plus = Variable(index=[time, aerosol_plus_gases]) # Atmospheric burden of agent above pre-industrial levels (unit of E). GU_aerosol_plus = Variable(index=[time, aerosol_plus_gases]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). iIRFT100_aerosol_plus = Variable(index=[time, aerosol_plus_gases]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). - decay_rate_aerosol_plus = Variable(index=[time, aerosol_plus_gases, 4]) - R_aerosol_plus = Variable(index=[time, aerosol_plus_gases, 4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). + decay_rate_aerosol_plus = Variable(index=[time, aerosol_plus_gases, 4]) # Decay rates + R_aerosol_plus = Variable(index=[time, aerosol_plus_gases, 4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). function run_timestep(p, v, d, t) diff --git a/src/components/ch4_cycle.jl b/src/components/ch4_cycle.jl index cae4ae6..83fcd12 100644 --- a/src/components/ch4_cycle.jl +++ b/src/components/ch4_cycle.jl @@ -4,29 +4,29 @@ @defcomp ch4_cycle begin - ch4_0 = Parameter() # Methane concentration in initial model period (ppb). - ch4_pi = Parameter() # Pre-industrial methane concentrations (ppb). - g0_ch4 = Parameter() # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 methane cycle parameterization at α=1. - g1_ch4 = Parameter() # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 methane cycle parameterization at α=1. - emiss2conc_ch4 = Parameter() # Conversion factor between emissions and concentrations. - GU_ch4_0 = Parameter() # Initial model period value for cumulative uptake of agent (unit of E⁻¹). - r0_ch4 = Parameter() # Strength of pre-industrial uptake from atmosphere. - rA_ch4 = Parameter() # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). - rT_ch4 = Parameter() # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). - rU_ch4 = Parameter() # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). - τ_ch4 = Parameter(index=[4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). - a_ch4 = Parameter(index=[4]) # Fraction of emissions entering iᵗʰ methane reservior. - R0_ch4 = Parameter(index=[4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). - E_ch4 = Parameter(index=[time]) # Annual methane emissions (TgCH₄ yr⁻¹). - Tj = Parameter(index=[time,3]) # Temperature change for three thermal pools (K). - - α_ch4 = Variable(index=[time]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. - ch4 = Variable(index=[time]) # Total atmospheric methane concentrations (ppb). - GA_ch4 = Variable(index=[time]) # Atmospheric burden of agent above pre-industrial levels (unit of E). - GU_ch4 = Variable(index=[time]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). - iIRFT100_ch4 = Variable(index=[time]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). - decay_rate_ch4 = Variable(index=[time,4]) - R_ch4 = Variable(index=[time,4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). + ch4_0 = Parameter() # Methane concentration in initial model period (ppb). + ch4_pi = Parameter() # Pre-industrial methane concentrations (ppb). + g0_ch4 = Parameter() # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 methane cycle parameterization at α=1. + g1_ch4 = Parameter() # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 methane cycle parameterization at α=1. + emiss2conc_ch4 = Parameter() # Conversion factor between emissions and concentrations. + GU_ch4_0 = Parameter() # Initial model period value for cumulative uptake of agent (unit of E⁻¹). + r0_ch4 = Parameter() # Strength of pre-industrial uptake from atmosphere. + rA_ch4 = Parameter() # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). + rT_ch4 = Parameter() # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). + rU_ch4 = Parameter() # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). + τ_ch4 = Parameter(index=[4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). + a_ch4 = Parameter(index=[4]) # Fraction of emissions entering iᵗʰ methane reservior. + R0_ch4 = Parameter(index=[4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). + E_ch4 = Parameter(index=[time]) # Annual methane emissions (TgCH₄ yr⁻¹). + Tj = Parameter(index=[time,3]) # Temperature change for three thermal pools (K). + + α_ch4 = Variable(index=[time]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. + ch4 = Variable(index=[time]) # Total atmospheric methane concentrations (ppb). + GA_ch4 = Variable(index=[time]) # Atmospheric burden of agent above pre-industrial levels (unit of E). + GU_ch4 = Variable(index=[time]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). + iIRFT100_ch4 = Variable(index=[time]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). + decay_rate_ch4 = Variable(index=[time,4]) # Decay rates + R_ch4 = Variable(index=[time,4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). function run_timestep(p, v, d, t) diff --git a/src/components/co2_cycle.jl b/src/components/co2_cycle.jl index 477ecda..7416d47 100644 --- a/src/components/co2_cycle.jl +++ b/src/components/co2_cycle.jl @@ -4,29 +4,29 @@ @defcomp co2_cycle begin - co2_0 = Parameter() # Carbon dioxide concentration in initial model period (ppm). - co2_pi = Parameter() # Pre-industrial carbon dioxide concentrations (ppm). - g0_co2 = Parameter() # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 carbon cycle parameterization at α=1. - g1_co2 = Parameter() # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 carbon cycle parameterization at α=1. - emiss2conc_co2 = Parameter() # Conversion factor between emissions (GtC) and concentrations (ppm). - GU_co2_0 = Parameter() # Initial model period value for cumulative uptake of agent (unit of E⁻¹). - r0_co2 = Parameter() # Strength of pre-industrial uptake from atmosphere. - rA_co2 = Parameter() # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). - rT_co2 = Parameter() # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). - rU_co2 = Parameter() # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). - τ_co2 = Parameter(index=[4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). - a_co2 = Parameter(index=[4]) # Fraction of emissions entering iᵗʰ carbon reservior. - R0_co2 = Parameter(index=[4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). - E_co2 = Parameter(index=[time]) # Annual carbon dioxide emissions (GtC yr⁻¹). - Tj = Parameter(index=[time,3]) # Temperature change for three thermal pools (K). - - α_co2 = Variable(index=[time]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. - co2 = Variable(index=[time]) # Total atmospheric carbon dioxide concentrations (ppm). - GA_co2 = Variable(index=[time]) # Atmospheric burden of agent above pre-industrial levels (unit of E). - GU_co2 = Variable(index=[time]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). - iIRFT100_co2 = Variable(index=[time]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). - decay_rate_co2 = Variable(index=[time,4]) - R_co2 = Variable(index=[time,4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). + co2_0 = Parameter() # Carbon dioxide concentration in initial model period (ppm). + co2_pi = Parameter() # Pre-industrial carbon dioxide concentrations (ppm). + g0_co2 = Parameter() # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 carbon cycle parameterization at α=1. + g1_co2 = Parameter() # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 carbon cycle parameterization at α=1. + emiss2conc_co2 = Parameter() # Conversion factor between emissions (GtC) and concentrations (ppm). + GU_co2_0 = Parameter() # Initial model period value for cumulative uptake of agent (unit of E⁻¹). + r0_co2 = Parameter() # Strength of pre-industrial uptake from atmosphere. + rA_co2 = Parameter() # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). + rT_co2 = Parameter() # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). + rU_co2 = Parameter() # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). + τ_co2 = Parameter(index=[4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). + a_co2 = Parameter(index=[4]) # Fraction of emissions entering iᵗʰ carbon reservior. + R0_co2 = Parameter(index=[4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). + E_co2 = Parameter(index=[time]) # Annual carbon dioxide emissions (GtC yr⁻¹). + Tj = Parameter(index=[time,3]) # Temperature change for three thermal pools (K). + + α_co2 = Variable(index=[time]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. + co2 = Variable(index=[time]) # Total atmospheric carbon dioxide concentrations (ppm). + GA_co2 = Variable(index=[time]) # Atmospheric burden of agent above pre-industrial levels (unit of E). + GU_co2 = Variable(index=[time]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). + iIRFT100_co2 = Variable(index=[time]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). + decay_rate_co2 = Variable(index=[time,4]) # Decay rates + R_co2 = Variable(index=[time,4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). function run_timestep(p, v, d, t) diff --git a/src/components/flourinated_gas_cycles.jl b/src/components/flourinated_gas_cycles.jl index db93aa7..1769941 100644 --- a/src/components/flourinated_gas_cycles.jl +++ b/src/components/flourinated_gas_cycles.jl @@ -4,31 +4,31 @@ @defcomp flourinated_cycles begin - flourinated_gases = Index() # Index for flourinated gases controlled under the Kyoto Protocol. - - flourinated_0 = Parameter(index=[flourinated_gases]) # flourinated gas concentration in initial model period (ppb). - flourinated_pi = Parameter(index=[flourinated_gases]) # Pre-industrial flourinated gas concentrations (ppb). - g0_flourinated = Parameter(index=[flourinated_gases]) # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 flourinated gas cycle parameterization at α=1. - g1_flourinated = Parameter(index=[flourinated_gases]) # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 flourinated gas cycle parameterization at α=1. - emiss2conc_flourinated = Parameter(index=[flourinated_gases]) # Conversion factor between emissions and concentrations. - GU_flourinated_0 = Parameter(index=[flourinated_gases]) # Initial model period value for cumulative uptake of agent (unit of E⁻¹). - r0_flourinated = Parameter(index=[flourinated_gases]) # Strength of pre-industrial uptake from atmosphere. - rA_flourinated = Parameter(index=[flourinated_gases]) # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). - rT_flourinated = Parameter(index=[flourinated_gases]) # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). - rU_flourinated = Parameter(index=[flourinated_gases]) # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). - τ_flourinated = Parameter(index=[flourinated_gases, 4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). - a_flourinated = Parameter(index=[flourinated_gases, 4]) # Fraction of emissions entering iᵗʰ flourinated gas reservior. - R0_flourinated = Parameter(index=[flourinated_gases, 4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). - E_flourinated = Parameter(index=[time, flourinated_gases]) # Annual flourinated gas emissions (Tg yr⁻¹). - Tj = Parameter(index=[time, 3]) # Temperature change for three thermal pools (K). - - α_flourinated = Variable(index=[time, flourinated_gases]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. - flourinated_conc = Variable(index=[time, flourinated_gases]) # Total atmospheric flourinated gas concentrations (ppb). - GA_flourinated = Variable(index=[time, flourinated_gases]) # Atmospheric burden of agent above pre-industrial levels (unit of E). - GU_flourinated = Variable(index=[time, flourinated_gases]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). - iIRFT100_flourinated = Variable(index=[time, flourinated_gases]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). - decay_rate_flourinated = Variable(index=[time, flourinated_gases, 4]) - R_flourinated = Variable(index=[time, flourinated_gases, 4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). + flourinated_gases = Index() # Index for flourinated gases controlled under the Kyoto Protocol. + + flourinated_0 = Parameter(index=[flourinated_gases]) # flourinated gas concentration in initial model period (ppb). + flourinated_pi = Parameter(index=[flourinated_gases]) # Pre-industrial flourinated gas concentrations (ppb). + g0_flourinated = Parameter(index=[flourinated_gases]) # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 flourinated gas cycle parameterization at α=1. + g1_flourinated = Parameter(index=[flourinated_gases]) # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 flourinated gas cycle parameterization at α=1. + emiss2conc_flourinated = Parameter(index=[flourinated_gases]) # Conversion factor between emissions and concentrations. + GU_flourinated_0 = Parameter(index=[flourinated_gases]) # Initial model period value for cumulative uptake of agent (unit of E⁻¹). + r0_flourinated = Parameter(index=[flourinated_gases]) # Strength of pre-industrial uptake from atmosphere. + rA_flourinated = Parameter(index=[flourinated_gases]) # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). + rT_flourinated = Parameter(index=[flourinated_gases]) # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). + rU_flourinated = Parameter(index=[flourinated_gases]) # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). + τ_flourinated = Parameter(index=[flourinated_gases, 4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). + a_flourinated = Parameter(index=[flourinated_gases, 4]) # Fraction of emissions entering iᵗʰ flourinated gas reservior. + R0_flourinated = Parameter(index=[flourinated_gases, 4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). + E_flourinated = Parameter(index=[time, flourinated_gases]) # Annual flourinated gas emissions (Tg yr⁻¹). + Tj = Parameter(index=[time, 3]) # Temperature change for three thermal pools (K). + + α_flourinated = Variable(index=[time, flourinated_gases]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. + flourinated_conc = Variable(index=[time, flourinated_gases]) # Total atmospheric flourinated gas concentrations (ppb). + GA_flourinated = Variable(index=[time, flourinated_gases]) # Atmospheric burden of agent above pre-industrial levels (unit of E). + GU_flourinated = Variable(index=[time, flourinated_gases]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). + iIRFT100_flourinated = Variable(index=[time, flourinated_gases]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). + decay_rate_flourinated = Variable(index=[time, flourinated_gases, 4]) # Decay rates + R_flourinated = Variable(index=[time, flourinated_gases, 4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). function run_timestep(p, v, d, t) diff --git a/src/components/montreal_gas_cycles.jl b/src/components/montreal_gas_cycles.jl index c555190..d0ab6af 100644 --- a/src/components/montreal_gas_cycles.jl +++ b/src/components/montreal_gas_cycles.jl @@ -4,31 +4,31 @@ @defcomp montreal_cycles begin - montreal_gases = Index() # Index for Ozone Depleting Substances controlled under the Montreal Protocol. + montreal_gases = Index() # Index for Ozone Depleting Substances controlled under the Montreal Protocol. montreal_0 = Parameter(index=[montreal_gases]) # Montreal gas concentration in initial model period (ppb). - montreal_pi = Parameter(index=[montreal_gases]) # Pre-industrial Montreal gas concentrations (ppb). + montreal_pi = Parameter(index=[montreal_gases]) # Pre-industrial Montreal gas concentrations (ppb). g0_montreal = Parameter(index=[montreal_gases]) # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 Montreal gas cycle parameterization at α=1. g1_montreal = Parameter(index=[montreal_gases]) # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 Montreal gas cycle parameterization at α=1. - emiss2conc_montreal = Parameter(index=[montreal_gases]) # Conversion factor between emissions and concentrations. + emiss2conc_montreal = Parameter(index=[montreal_gases]) # Conversion factor between emissions and concentrations. GU_montreal_0 = Parameter(index=[montreal_gases]) # Initial model period value for cumulative uptake of agent (unit of E⁻¹). - r0_montreal = Parameter(index=[montreal_gases]) # Strength of pre-industrial uptake from atmosphere. + r0_montreal = Parameter(index=[montreal_gases]) # Strength of pre-industrial uptake from atmosphere. rA_montreal = Parameter(index=[montreal_gases]) # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). rT_montreal = Parameter(index=[montreal_gases]) # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). rU_montreal = Parameter(index=[montreal_gases]) # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). - τ_montreal = Parameter(index=[montreal_gases, 4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). - a_montreal = Parameter(index=[montreal_gases, 4]) # Fraction of emissions entering iᵗʰ Montreal gas reservior. - R0_montreal = Parameter(index=[montreal_gases, 4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). - E_montreal = Parameter(index=[time, montreal_gases]) # Annual Montreal gas emissions (Tg yr⁻¹). - Tj = Parameter(index=[time, 3]) # Temperature change for three thermal pools (K). + τ_montreal = Parameter(index=[montreal_gases, 4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). + a_montreal = Parameter(index=[montreal_gases, 4]) # Fraction of emissions entering iᵗʰ Montreal gas reservior. + R0_montreal = Parameter(index=[montreal_gases, 4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). + E_montreal = Parameter(index=[time, montreal_gases]) # Annual Montreal gas emissions (Tg yr⁻¹). + Tj = Parameter(index=[time, 3]) # Temperature change for three thermal pools (K). α_montreal = Variable(index=[time, montreal_gases]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. montreal_conc = Variable(index=[time, montreal_gases]) # Total atmospheric Montreal gas concentrations (ppb). GA_montreal = Variable(index=[time, montreal_gases]) # Atmospheric burden of agent above pre-industrial levels (unit of E). GU_montreal = Variable(index=[time, montreal_gases]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). iIRFT100_montreal = Variable(index=[time, montreal_gases]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). - decay_rate_montreal = Variable(index=[time, montreal_gases, 4]) - R_montreal = Variable(index=[time, montreal_gases, 4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). + decay_rate_montreal = Variable(index=[time, montreal_gases, 4]) # Decay rates + R_montreal = Variable(index=[time, montreal_gases, 4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). function run_timestep(p, v, d, t) diff --git a/src/components/n2o_cycle.jl b/src/components/n2o_cycle.jl index 28a29c7..60b2547 100644 --- a/src/components/n2o_cycle.jl +++ b/src/components/n2o_cycle.jl @@ -4,29 +4,29 @@ @defcomp n2o_cycle begin - n2o_0 = Parameter() # Nitrous oxide concentration in initial model period (ppb). - n2o_pi = Parameter() # Pre-industrial nitrous oxide concentrations (ppb). - g0_n2o = Parameter() # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 nitrous oxide cycle parameterization at α=1. - g1_n2o = Parameter() # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 nitrous oxide cycle parameterization at α=1. - emiss2conc_n2o = Parameter() # Conversion factor between emissions and concentrations. - GU_n2o_0 = Parameter() # Initial model period value for cumulative uptake of agent (unit of E⁻¹). - r0_n2o = Parameter() # Strength of pre-industrial uptake from atmosphere. - rA_n2o = Parameter() # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). - rT_n2o = Parameter() # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). - rU_n2o = Parameter() # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). - τ_n2o = Parameter(index=[4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). - a_n2o = Parameter(index=[4]) # Fraction of emissions entering iᵗʰ nitrous oxide reservior. - R0_n2o = Parameter(index=[4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). - E_n2o = Parameter(index=[time]) # Annual nitrous oxide emissions (TgN yr⁻¹). - Tj = Parameter(index=[time,3]) # Temperature change for three thermal pools (K). - - α_n2o = Variable(index=[time]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. - n2o = Variable(index=[time]) # Total atmospheric nitrous oxide concentrations (ppb). - GA_n2o = Variable(index=[time]) # Atmospheric burden of agent above pre-industrial levels (unit of E). - GU_n2o = Variable(index=[time]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). - iIRFT100_n2o = Variable(index=[time]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). - decay_rate_n2o = Variable(index=[time,4]) - R_n2o = Variable(index=[time,4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). + n2o_0 = Parameter() # Nitrous oxide concentration in initial model period (ppb). + n2o_pi = Parameter() # Pre-industrial nitrous oxide concentrations (ppb). + g0_n2o = Parameter() # Constant to set approximation of value for α equal to Millar et al. (2017) numerical solution for iIRF100 nitrous oxide cycle parameterization at α=1. + g1_n2o = Parameter() # Constant to set approximation of gradient for α equal to Millar et al. (2017) numerical solution for iIRF100 nitrous oxide cycle parameterization at α=1. + emiss2conc_n2o = Parameter() # Conversion factor between emissions and concentrations. + GU_n2o_0 = Parameter() # Initial model period value for cumulative uptake of agent (unit of E⁻¹). + r0_n2o = Parameter() # Strength of pre-industrial uptake from atmosphere. + rA_n2o = Parameter() # Sensitivity of uptake from atmosphere to current atmospheric burden of agent (unit of E⁻¹). + rT_n2o = Parameter() # Sensitivity of uptake from atmosphere to model temperature change since initialization (K⁻¹). + rU_n2o = Parameter() # Sensitivity of uptake from atmosphere to cumulative uptake of agent since model initialization (unit of E⁻¹). + τ_n2o = Parameter(index=[4]) # Atmospheric lifetime of gas in iᵗʰ reservior (years). + a_n2o = Parameter(index=[4]) # Fraction of emissions entering iᵗʰ nitrous oxide reservior. + R0_n2o = Parameter(index=[4]) # Initial model period value for quantity of agent in iᵗʰ atmospheric reservior (unit of E). + E_n2o = Parameter(index=[time]) # Annual nitrous oxide emissions (TgN yr⁻¹). + Tj = Parameter(index=[time,3]) # Temperature change for three thermal pools (K). + + α_n2o = Variable(index=[time]) # State-dependent multiplicative adjustment coefficient of reservior lifetimes. + n2o = Variable(index=[time]) # Total atmospheric nitrous oxide concentrations (ppb). + GA_n2o = Variable(index=[time]) # Atmospheric burden of agent above pre-industrial levels (unit of E). + GU_n2o = Variable(index=[time]) # Cumulative uptake of agent since model initialization (unit of E⁻¹). + iIRFT100_n2o = Variable(index=[time]) # 100-year integrated impulse response function (the average airborne fraction over a 100-year period). + decay_rate_n2o = Variable(index=[time,4]) # Decay rates + R_n2o = Variable(index=[time,4]) # Quantity of agent in iᵗʰ atmospheric reservior (unit of E). function run_timestep(p, v, d, t) diff --git a/src/components/radiative_forcing.jl b/src/components/radiative_forcing.jl index 676f704..5ee7755 100644 --- a/src/components/radiative_forcing.jl +++ b/src/components/radiative_forcing.jl @@ -4,64 +4,80 @@ @defcomp radiative_forcing begin - montreal_gases = Index() # Index for Ozone Depleting Substances controlled under the Montreal Protocol. - aerosol_plus_gases = Index() # Index for tropospheric ozone precursors, aerosols and reactive gas - flourinated_gases = Index() # Index for flourinated gases controlled under the Kyoto Protocol. - - # CO2, CH4, and N2O parametes - co2_conc = Parameter(index=[time]) - ch4_conc = Parameter(index=[time]) - n2o_conc = Parameter(index=[time]) - co2_pi = Parameter() - ch4_pi = Parameter() - n2o_pi = Parameter() - co2_f = Parameter(index=[3]) # f coefficient for - ch4_f = Parameter(index=[3]) # f coefficient for - ch4_o3_f = Parameter(index=[3]) # f coefficient for - ch4_h2o_f = Parameter(index=[3]) # f coefficient for - n2o_f = Parameter(index=[3]) # f coefficient for - n2o_o3_f = Parameter(index=[3]) # f coefficient for + montreal_gases = Index() # Index for Ozone Depleting Substances controlled under the Montreal Protocol. + aerosol_plus_gases = Index() # Index for tropospheric ozone precursors, aerosols and reactive gas + flourinated_gases = Index() # Index for flourinated gases controlled under the Kyoto Protocol. + + ## + ## Parameters + ## + + # CO2, CH4, and N2O parameters + co2_conc = Parameter(index=[time]) + ch4_conc = Parameter(index=[time]) + n2o_conc = Parameter(index=[time]) + co2_pi = Parameter() + ch4_pi = Parameter() + n2o_pi = Parameter() + co2_f = Parameter(index=[3]) # f coefficient for co2 + ch4_f = Parameter(index=[3]) # f coefficient for ch4 + ch4_o3_f = Parameter(index=[3]) # f coefficient for ch4_o3 + ch4_h2o_f = Parameter(index=[3]) # f coefficient for ch4_h2o + n2o_f = Parameter(index=[3]) # f coefficient for n2o_f + n2o_o3_f = Parameter(index=[3]) # f coefficient for n2o_o3 # Aerosol+ parameters - aerosol_plus_conc = Parameter(index=[time, aerosol_plus_gases]) - aerosol_plus_pi = Parameter(index=[aerosol_plus_gases]) - bc_f = Parameter(index=[3]) # f coefficient for - bc_snow_f = Parameter(index=[3]) # f coefficient for - bc_aci_f = Parameter(index=[3]) # f coefficient for - co_f = Parameter(index=[3]) # f coefficient for - co_o3_f = Parameter(index=[3]) # f coefficient for - nh3_f = Parameter(index=[3]) # f coefficient for - nmvoc_f = Parameter(index=[3]) # f coefficient for - nmvoc_o3_f = Parameter(index=[3]) # f coefficient for - nox_f = Parameter(index=[3]) # f coefficient for - nox_o3_f = Parameter(index=[3]) # f coefficient for - nox_avi_f = Parameter(index=[3]) # f coefficient for - nox_avi_contrails_f = Parameter(index=[3]) # f coefficient for - oc_f = Parameter(index=[3]) # f coefficient for - oc_aci_f = Parameter(index=[3]) # f coefficient for - so2_f = Parameter(index=[3]) # f coefficient for - so2_aci_f = Parameter(index=[3]) # f coefficient for + aerosol_plus_conc = Parameter(index=[time, aerosol_plus_gases]) + aerosol_plus_pi = Parameter(index=[aerosol_plus_gases]) + bc_f = Parameter(index=[3]) # f coefficient for bc + bc_snow_f = Parameter(index=[3]) # f coefficient for bc_snow + bc_aci_f = Parameter(index=[3]) # f coefficient for bc_aci + co_f = Parameter(index=[3]) # f coefficient for co + co_o3_f = Parameter(index=[3]) # f coefficient for co_o3 + nh3_f = Parameter(index=[3]) # f coefficient for nh3 + nmvoc_f = Parameter(index=[3]) # f coefficient for nmvoc + nmvoc_o3_f = Parameter(index=[3]) # f coefficient for nmvoc_o3 + nox_f = Parameter(index=[3]) # f coefficient for nox + nox_o3_f = Parameter(index=[3]) # f coefficient for nox_o3 + nox_avi_f = Parameter(index=[3]) # f coefficient for nox_avi + nox_avi_contrails_f = Parameter(index=[3]) # f coefficient for nox_avi_contrails + oc_f = Parameter(index=[3]) # f coefficient for oc + oc_aci_f = Parameter(index=[3]) # f coefficient for oc_aci + so2_f = Parameter(index=[3]) # f coefficient for so2 + so2_aci_f = Parameter(index=[3]) # f coefficient for so2_aci # Montreal gas parameters - montreal_conc = Parameter(index=[time, montreal_gases]) - montreal_pi = Parameter(index=[montreal_gases]) - montreal_f = Parameter(index=[montreal_gases, 3]) # f coefficient - montreal_ind_f = Parameter(index=[montreal_gases, 3]) # f coefficient for indirect forcing effect on O3 + montreal_conc = Parameter(index=[time, montreal_gases]) + montreal_pi = Parameter(index=[montreal_gases]) + montreal_f = Parameter(index=[montreal_gases, 3]) # f coefficient + montreal_ind_f = Parameter(index=[montreal_gases, 3]) # f coefficient for indirect forcing effect on O3 # Flourinated gas parameters - flourinated_conc = Parameter(index=[time, flourinated_gases]) - flourinated_pi = Parameter(index=[flourinated_gases]) - flourinated_f = Parameter(index=[flourinated_gases, 3]) # f coefficient + flourinated_conc = Parameter(index=[time, flourinated_gases]) + flourinated_pi = Parameter(index=[flourinated_gases]) + flourinated_f = Parameter(index=[flourinated_gases, 3]) # f coefficient - exogenous_RF = Parameter(index=[time]) + # Other parameters + exogenous_RF = Parameter(index=[time]) + ## + ## Variables + ## + + # CO2, CH4, and N2O variables + co2_RF = Variable(index=[time]) + ch4_RF = Variable(index=[time]) + ch4_o3_RF = Variable(index=[time]) + ch4_h2o_RF = Variable(index=[time]) + n2o_RF = Variable(index=[time]) + n2o_o3_RF = Variable(index=[time]) + + # Aerosol+ variables bc_RF = Variable(index=[time]) bc_snow_RF = Variable(index=[time]) bc_aci_RF = Variable(index=[time]) - co_RF = Variable(index=[time]) co_o3_RF = Variable(index=[time]) - nh3_RF = Variable(index=[time]) nmvoc_RF = Variable(index=[time]) nmvoc_o3_RF = Variable(index=[time]) @@ -74,20 +90,14 @@ so2_RF = Variable(index=[time]) so2_aci_RF = Variable(index=[time]) - co2_RF = Variable(index=[time]) - - ch4_RF = Variable(index=[time]) - ch4_o3_RF = Variable(index=[time]) - ch4_h2o_RF = Variable(index=[time]) - - n2o_RF = Variable(index=[time]) - n2o_o3_RF = Variable(index=[time]) - + # Montreal gas variables montreal_RF = Variable(index=[time, montreal_gases]) montreal_ind_RF = Variable(index=[time, montreal_gases]) - + + # Flourinated gas variables flourinated_RF = Variable(index=[time, flourinated_gases]) + # Other variables total_RF = Variable(index=[time]) function run_timestep(p, v, d, t) diff --git a/src/components/temperature.jl b/src/components/temperature.jl index 57b3605..208deb3 100644 --- a/src/components/temperature.jl +++ b/src/components/temperature.jl @@ -6,10 +6,10 @@ #d = Parameter(index=[3]) # Thermal response timescales: [1] Thermal equilibration of deep ocean & [2] Thermal admustment of upper ocean (years). decay_factor = Parameter(index=[3]) # Thermal response decay factor, calculated as exp(-1/d) where d represents thermal response timescale of jth thermal box. - q = Parameter(index=[3]) # Raditive forcing coefficient: [1] Thermal equilibration of deep ocean & [2] Thermal admustment of upper ocean (K W⁻¹m²). - F = Parameter(index=[time]) # Total radiative forcing (Wm⁻²). - Tj_0 = Parameter(index=[3]) - T_0 = Parameter() + q = Parameter(index=[3]) # Raditive forcing coefficient: [1] Thermal equilibration of deep ocean & [2] Thermal admustment of upper ocean (K W⁻¹m²). + F = Parameter(index=[time]) # Total radiative forcing (Wm⁻²). + Tj_0 = Parameter(index=[3]) + T_0 = Parameter() T = Variable(index=[time]) # Global mean surface temperature anomaly (K). Tj = Variable(index=[time,3]) # Temperature change for three thermal pools (K). diff --git a/src/helper_functions.jl b/src/helper_functions.jl index 7da2811..41cb731 100644 --- a/src/helper_functions.jl +++ b/src/helper_functions.jl @@ -1,6 +1,6 @@ # #---------------------------------------------------------------------------------------------------------------------- # #---------------------------------------------------------------------------------------------------------------------- -# This file contains functions and other snippets of code that are used in various calculations for Mimi-FAIRv2.0. +# This file contains functions and other snippets of code that are used in various calculations for MimiFAIRv2. # #---------------------------------------------------------------------------------------------------------------------- # #---------------------------------------------------------------------------------------------------------------------- @@ -21,6 +21,20 @@ # F2x: Radiative forcing from a doubling of carbon dioxide concentrations (Wm⁻²). #---------------------------------------------------------------------------------------------------------------------- +""" + get_thermal_parameter_defaults(;TCR::Float64=1.79, RWF::Float64=0.552, F2x::Float64=3.759) + +Return a dataframe of the FAIRv2.0.0-alpha default climate parameters using the +optional keyword arguments for transient climate response (K) (TCR; defaults to 1.79), +realized warming fraction in ratio of TCR/ECS (RWF; defaults to 0.552), and radiative +forcing from a double of carbon dioxide concentrations (Wm⁻²) (F2x; defaults to +3.759). + +The response timescales d1-3 (and the shortest-timescale coefficient, q1) are set +to the central estimate of a CMIP6 inferred distribution constrained with observational +warming. The constraint does not significantly affect the central estimates of +the prior (ie. raw CMIP6 inference) distribution. +""" function get_thermal_parameter_defaults(;TCR::Float64=1.79, RWF::Float64=0.552, F2x::Float64=3.759) # Set default values for d1, d2, d3, and q1 parameters. @@ -63,6 +77,15 @@ end #---------------------------------------------------------------------------------------------------------------------- # Function version for a single gas (where a and τ are vectors). +""" + calculate_g0_g1(a::Array{Float64,1}, τ::Array{Float64,1}) + +For a SINGLE gas, calculate and return constants for estimating the state-dependent +adjustment coefficient of a reservior's lifetime (α) such that it approximates +the Millar et al. (2017) numerical solution for a iIRF100 carbon cycle parameterization +at α=1. The two required arguments are a, the fraction of emissions entering iᵗʰ +atmospheric pool, and τ, the atmospheric lifetime of gas in iᵗʰ pool. +""" function calculate_g0_g1(a::Array{Float64,1}, τ::Array{Float64,1}) g1 = sum(a .* τ .* (1.0 .- (1.0 .+ 100.0 ./ τ) .* exp.(-100.0 ./ τ))) g0 = exp(-1.0 * sum(a .* τ .* (1.0 .- exp.(-100.0 ./ τ))) / g1) @@ -70,6 +93,15 @@ function calculate_g0_g1(a::Array{Float64,1}, τ::Array{Float64,1}) end # Function version for multiple gases (where a and τ are 2-d arrays). +""" + calculate_g0_g1(a::Array{Float64,1}, τ::Array{Float64,1}) + +For MULITPLE gases, calculate and return constants for estimating the state-dependent +adjustment coefficient of a reservior's lifetime (α) such that it approximates +the Millar et al. (2017) numerical solution for a iIRF100 carbon cycle parameterization +at α=1. The two required arguments are a, the fraction of emissions entering iᵗʰ +atmospheric pool, and τ, the atmospheric lifetime of gas in iᵗʰ pool. +""" function calculate_g0_g1(a::Array{Float64,2}, τ::Array{Float64,2}) g1 = vec(sum(a .* τ .* (1.0 .- (1.0 .+ 100.0 ./ τ) .* exp.(-100.0 ./ τ)), dims=2)) g0 = vec(exp.(-1.0 .* sum(a .* τ .* (1.0 .- exp.(-100.0 ./ τ)), dims=2) ./ g1)) @@ -91,6 +123,18 @@ end # C_pi: Pre-industrial concentration of forcing agent in atmosphere. #---------------------------------------------------------------------------------------------------------------------- +""" + calculate_RF(f, C, C_pi) + +Calculate and return the effective radiative forcing based on the atmospheric +concentration of a greenhouse gas or other forcing agent following Equation (6) +in Leach et al. (2021). + +The required arguments are as follows: + f: A vector of three concentration-forcing coefficients (1 = logarithmic, 2 = linear, 3 = square root) + C: Concentration of forcing agent in atmosphere + C_pi: Pre-industrial concentration of forcing agent in atmosphere +""" function calculate_RF(f, C, C_pi) if C <= 0.0 diff --git a/test/runtests.jl b/test/runtests.jl new file mode 100644 index 0000000..0b870db --- /dev/null +++ b/test/runtests.jl @@ -0,0 +1,106 @@ +using DataFrames +using CSVFiles +using Test + +include(joinpath(@__DIR__, "..", "src", "MimiFAIRv2.jl")) +using Main.MimiFAIRv2 + +@testset "UI" begin + + # run basic steps presented in README + m = get_model() + run(m) + fair_temps = m[:temperature, :T] + + # Run model with altered keyword args and make sure results reflect the + # appropriate changes, or at least a change at all + # - start and end years + # - emissions scenario + # - TCR, RWF, and F2x + + # end year + m = get_model(end_year = 2100) + run(m) + fair_temps = m[:temperature, :T] + @test length(fair_temps) == length(1750:2100) + + # start year - should throw a warning if use a start year other than 1750 + @test_logs (:warn, "Model should not be set to start with a year differing from 1750.") get_model(start_year = 1760) + m = get_model(start_year = 1760) + run(m) + fair_temps = m[:temperature, :T] + @test length(fair_temps) == length(1760:2500) + + # emissions scenario + m1 = get_model(emissions_forcing_scenario = "ssp126") + m2 = get_model(emissions_forcing_scenario = "ssp585") + run(m1) + run(m2) + m1_fair_temps = m1[:temperature, :T] + m2_fair_temps = m2[:temperature, :T] + # Mimi.plot(m1, :temperature, :T) + # Mimi.plot(m2, :temperature, :T) + @test m1_fair_temps !== m2_fair_temps + @test maximum(m1_fair_temps) < maximum(m2_fair_temps) + + # TCR - transient climate response (default to 1.79) + m1 = get_model() + m2 = get_model(TCR = 2.0) + run(m1) + run(m2) + m1_fair_temps = m1[:temperature, :T] + m2_fair_temps = m2[:temperature, :T] + # Mimi.plot(m1, :temperature, :T) + # Mimi.plot(m2, :temperature, :T) + @test m1_fair_temps !== m2_fair_temps + @test maximum(m1_fair_temps) < maximum(m2_fair_temps) + + # RWF - realized warming fraction (default to 0.552) + m1 = get_model() + m2 = get_model(RWF = 0.5) + run(m1) + run(m2) + m1_fair_temps = m1[:temperature, :T] + m2_fair_temps = m2[:temperature, :T] + # Mimi.plot(m1, :temperature, :T) + # Mimi.plot(m2, :temperature, :T) + @test m1_fair_temps !== m2_fair_temps + @test maximum(m1_fair_temps) < maximum(m2_fair_temps) + + # F2x - forcing from a doubling of CO₂ (default to 3.759) + m1 = get_model() + m2 = get_model(F2x = 3.5) + run(m1) + run(m2) + m1_fair_temps = m1[:temperature, :T] + m2_fair_temps = m2[:temperature, :T] + # Mimi.plot(m1, :temperature, :T) + # Mimi.plot(m2, :temperature, :T) + @test m1_fair_temps !== m2_fair_temps + @test maximum(m1_fair_temps) < maximum(m2_fair_temps) + +end + +@testset "Python Comparison" begin + + # run model for each possible SSP and compare variable values between julia + # and python versons + + for SSP in ["ssp119", "ssp126", "ssp245", "ssp370", "ssp585"] + + m = get_model(;emissions_forcing_scenario = SSP, end_year = 2100) + run(m) + + # compare temperatures (TODO - this is failing) + julia_temps = m[:temperature, :T] + + data = load(joinpath(@__DIR__, "..", "data", "python_replication_data", string(SSP, "_temperature.csv"))) |> DataFrame + rename!(data, [:year, :T]) + python_temps = data[!, :T] + + @test maximum(abs, julia_temps .- python_temps) ≈ 0. atol = 1e-3 + + # TODO - what other variables do we want to compare? + + end +end