diff --git a/.github/workflows/make_docs.yml b/.github/workflows/make_docs.yml new file mode 100644 index 0000000..e14ab77 --- /dev/null +++ b/.github/workflows/make_docs.yml @@ -0,0 +1,35 @@ +name: Make Docs +on: + pull_request: + branches: ["master"] + push: + branches: + - master + - dev + - docs + paths: + - '.github/workflows/make_docs.yml' + - 'src/' + - 'docs/**' + tags: '*' + workflow_dispatch: + +jobs: + make_docs: + permissions: + contents: write + statuses: write + name: Documentation + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v4 + - uses: julia-actions/setup-julia@latest + - uses: julia-actions/cache@v1 + - name: Install dependencies + run: | + julia --project=docs/ -e 'using Pkg; Pkg.add(; url="https://github.com/ProjectTorreyPines/IMASDD.jl.git"); Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()' + - name: Build and deploy + env: + GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} + DOCUMENTER_KEY: ${{ secrets.DOCUMENTER_KEY }} + run: julia --project=docs/ docs/make.jl \ No newline at end of file diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml new file mode 100644 index 0000000..df9cbb2 --- /dev/null +++ b/.github/workflows/test.yml @@ -0,0 +1,43 @@ +name: Test + +on: + push: + branches: ["master", "dev", "autotest"] + pull_request: + branches: ["master", "dev"] +jobs: + test: + runs-on: ${{ matrix.os }} + strategy: + matrix: + julia-version: [1.9.3] + julia-arch: [x86] + os: [ubuntu-latest] + steps: + - uses: julia-actions/setup-julia@latest + with: + version: ${{ matrix.julia-version }} + - uses: actions/checkout@v4 + - uses: webfactory/ssh-agent@v0.9.0 + with: + ssh-private-key: | + ${{ secrets.SOLPSTESTSAMPLES_SSH_KEY}} + ${{ secrets.DVC_SSH_KEY }} + - name: Configure ssh + run: | + echo "${{ secrets.DVC_KNOWN_HOSTS }}" >> ~/.ssh/known_hosts + echo "${{ secrets.DVC_SSH_CONFIG }}" >> ~/.ssh/config + - uses: iterative/setup-dvc@v1 + - name: DVC Pull + run: | + dvc pull + - uses: julia-actions/cache@v1 + - name: Install dependencies + run: | + julia --project=. -e 'using Pkg; Pkg.rm(["IMASDD"]); Pkg.add(url="git@github.com:ProjectTorreyPines/IMASDD.jl.git", rev="master")' + - uses: julia-actions/julia-runtest@v1 + # Not set up yet + # - uses: julia-actions/julia-processcoverage@v1 + # - uses: codecov/codecov-action@v4 + # with: + # files: lcov.info \ No newline at end of file diff --git a/.gitignore b/.gitignore index d31be6e..b69d3d3 100644 --- a/.gitignore +++ b/.gitignore @@ -4,3 +4,4 @@ # environment. Manifest.toml examples/Project.toml +docs/build/ \ No newline at end of file diff --git a/LICENSE b/LICENSE index 28143ce..171d76c 100644 --- a/LICENSE +++ b/LICENSE @@ -1,21 +1,201 @@ -MIT License - -Copyright (c) 2023 Project Torrey Pines - -Permission is hereby granted, free of charge, to any person obtaining a copy -of this software and associated documentation files (the "Software"), to deal -in the Software without restriction, including without limitation the rights -to use, copy, modify, merge, publish, distribute, sublicense, and/or sell -copies of the Software, and to permit persons to whom the Software is -furnished to do so, subject to the following conditions: - -The above copyright notice and this permission notice shall be included in all -copies or substantial portions of the Software. - -THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE -SOFTWARE. + Apache License + Version 2.0, January 2004 + http://www.apache.org/licenses/ + + TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION + + 1. Definitions. + + "License" shall mean the terms and conditions for use, reproduction, + and distribution as defined by Sections 1 through 9 of this document. + + "Licensor" shall mean the copyright owner or entity authorized by + the copyright owner that is granting the License. + + "Legal Entity" shall mean the union of the acting entity and all + other entities that control, are controlled by, or are under common + control with that entity. For the purposes of this definition, + "control" means (i) the power, direct or indirect, to cause the + direction or management of such entity, whether by contract or + otherwise, or (ii) ownership of fifty percent (50%) or more of the + outstanding shares, or (iii) beneficial ownership of such entity. + + "You" (or "Your") shall mean an individual or Legal Entity + exercising permissions granted by this License. + + "Source" form shall mean the preferred form for making modifications, + including but not limited to software source code, documentation + source, and configuration files. + + "Object" form shall mean any form resulting from mechanical + transformation or translation of a Source form, including but + not limited to compiled object code, generated documentation, + and conversions to other media types. + + "Work" shall mean the work of authorship, whether in Source or + Object form, made available under the License, as indicated by a + copyright notice that is included in or attached to the work + (an example is provided in the Appendix below). + + "Derivative Works" shall mean any work, whether in Source or Object + form, that is based on (or derived from) the Work and for which the + editorial revisions, annotations, elaborations, or other modifications + represent, as a whole, an original work of authorship. For the purposes + of this License, Derivative Works shall not include works that remain + separable from, or merely link (or bind by name) to the interfaces of, + the Work and Derivative Works thereof. + + "Contribution" shall mean any work of authorship, including + the original version of the Work and any modifications or additions + to that Work or Derivative Works thereof, that is intentionally + submitted to Licensor for inclusion in the Work by the copyright owner + or by an individual or Legal Entity authorized to submit on behalf of + the copyright owner. For the purposes of this definition, "submitted" + means any form of electronic, verbal, or written communication sent + to the Licensor or its representatives, including but not limited to + communication on electronic mailing lists, source code control systems, + and issue tracking systems that are managed by, or on behalf of, the + Licensor for the purpose of discussing and improving the Work, but + excluding communication that is conspicuously marked or otherwise + designated in writing by the copyright owner as "Not a Contribution." + + "Contributor" shall mean Licensor and any individual or Legal Entity + on behalf of whom a Contribution has been received by Licensor and + subsequently incorporated within the Work. + + 2. Grant of Copyright License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + copyright license to reproduce, prepare Derivative Works of, + publicly display, publicly perform, sublicense, and distribute the + Work and such Derivative Works in Source or Object form. + + 3. Grant of Patent License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + (except as stated in this section) patent license to make, have made, + use, offer to sell, sell, import, and otherwise transfer the Work, + where such license applies only to those patent claims licensable + by such Contributor that are necessarily infringed by their + Contribution(s) alone or by combination of their Contribution(s) + with the Work to which such Contribution(s) was submitted. If You + institute patent litigation against any entity (including a + cross-claim or counterclaim in a lawsuit) alleging that the Work + or a Contribution incorporated within the Work constitutes direct + or contributory patent infringement, then any patent licenses + granted to You under this License for that Work shall terminate + as of the date such litigation is filed. + + 4. Redistribution. You may reproduce and distribute copies of the + Work or Derivative Works thereof in any medium, with or without + modifications, and in Source or Object form, provided that You + meet the following conditions: + + (a) You must give any other recipients of the Work or + Derivative Works a copy of this License; and + + (b) You must cause any modified files to carry prominent notices + stating that You changed the files; and + + (c) You must retain, in the Source form of any Derivative Works + that You distribute, all copyright, patent, trademark, and + attribution notices from the Source form of the Work, + excluding those notices that do not pertain to any part of + the Derivative Works; and + + (d) If the Work includes a "NOTICE" text file as part of its + distribution, then any Derivative Works that You distribute must + include a readable copy of the attribution notices contained + within such NOTICE file, excluding those notices that do not + pertain to any part of the Derivative Works, in at least one + of the following places: within a NOTICE text file distributed + as part of the Derivative Works; within the Source form or + documentation, if provided along with the Derivative Works; or, + within a display generated by the Derivative Works, if and + wherever such third-party notices normally appear. The contents + of the NOTICE file are for informational purposes only and + do not modify the License. You may add Your own attribution + notices within Derivative Works that You distribute, alongside + or as an addendum to the NOTICE text from the Work, provided + that such additional attribution notices cannot be construed + as modifying the License. + + You may add Your own copyright statement to Your modifications and + may provide additional or different license terms and conditions + for use, reproduction, or distribution of Your modifications, or + for any such Derivative Works as a whole, provided Your use, + reproduction, and distribution of the Work otherwise complies with + the conditions stated in this License. + + 5. Submission of Contributions. Unless You explicitly state otherwise, + any Contribution intentionally submitted for inclusion in the Work + by You to the Licensor shall be under the terms and conditions of + this License, without any additional terms or conditions. + Notwithstanding the above, nothing herein shall supersede or modify + the terms of any separate license agreement you may have executed + with Licensor regarding such Contributions. + + 6. Trademarks. This License does not grant permission to use the trade + names, trademarks, service marks, or product names of the Licensor, + except as required for reasonable and customary use in describing the + origin of the Work and reproducing the content of the NOTICE file. + + 7. Disclaimer of Warranty. Unless required by applicable law or + agreed to in writing, Licensor provides the Work (and each + Contributor provides its Contributions) on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or + implied, including, without limitation, any warranties or conditions + of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A + PARTICULAR PURPOSE. You are solely responsible for determining the + appropriateness of using or redistributing the Work and assume any + risks associated with Your exercise of permissions under this License. + + 8. Limitation of Liability. In no event and under no legal theory, + whether in tort (including negligence), contract, or otherwise, + unless required by applicable law (such as deliberate and grossly + negligent acts) or agreed to in writing, shall any Contributor be + liable to You for damages, including any direct, indirect, special, + incidental, or consequential damages of any character arising as a + result of this License or out of the use or inability to use the + Work (including but not limited to damages for loss of goodwill, + work stoppage, computer failure or malfunction, or any and all + other commercial damages or losses), even if such Contributor + has been advised of the possibility of such damages. + + 9. Accepting Warranty or Additional Liability. While redistributing + the Work or Derivative Works thereof, You may choose to offer, + and charge a fee for, acceptance of support, warranty, indemnity, + or other liability obligations and/or rights consistent with this + License. However, in accepting such obligations, You may act only + on Your own behalf and on Your sole responsibility, not on behalf + of any other Contributor, and only if You agree to indemnify, + defend, and hold each Contributor harmless for any liability + incurred by, or claims asserted against, such Contributor by reason + of your accepting any such warranty or additional liability. + + END OF TERMS AND CONDITIONS + + APPENDIX: How to apply the Apache License to your work. + + To apply the Apache License to your work, attach the following + boilerplate notice, with the fields enclosed by brackets "[]" + replaced with your own identifying information. (Don't include + the brackets!) The text should be enclosed in the appropriate + comment syntax for the file format. We also recommend that a + file or class name and description of purpose be included on the + same "printed page" as the copyright notice for easier + identification within third-party archives. + + Copyright 2024 General Atomics + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License. diff --git a/NOTICE.md b/NOTICE.md new file mode 100644 index 0000000..01d03d0 --- /dev/null +++ b/NOTICE.md @@ -0,0 +1,15 @@ +GGDUtils.jl +========= + +The purpose of this NOTICE file is to provide legal notices and acknowledgments that must be displayed to users in any derivative works or distributions. This file does not alter the terms of the Apache 2.0 license that governs the use and distribution of the GGDUtils.jl package. + +GGDUtils.jl was originally developed under the ProjectTorreyPines by the Magnetic Fusion Energy group at General Atomics. + +If this software contributes to an academic publication, please cite it as follows: +A. Gupta, D. Eldon, H. Anand, A. Dautt-Silva, S. De Pascuale, J. Lore, O. Meneghini, and J.S. Park, "Plasma boundary control development using a time-dependent scrape-off layer model in closed-loop simulations", Nucl. Mater. Energy, manuscript in preparation for PSI conference (2024). + +The names "General Atomics", and any associated logos or images, are trademarks of General Atomics. Use of these trademarks without prior written consent from General Atomics is strictly prohibited. Users cannot imply endorsement by General Atomics or contributors to the project simply because the project is part of their work. + +Copyright (c) 2024 General Atomics + +Version: v1.0 \ No newline at end of file diff --git a/Project.toml b/Project.toml index 073e5ca..5711386 100644 --- a/Project.toml +++ b/Project.toml @@ -1,14 +1,14 @@ name = "GGDUtils" uuid = "b7b5e640-9b39-4803-84eb-376048795def" authors = ["Anchal Gupta "] -version = "0.2.1" +version = "1.0.0" [deps] ArgParse = "c7e460c6-2fb9-53a9-8c5b-16f535851c63" ColorSchemes = "35d6a980-a343-548e-a6ea-1d62b119f2f4" +IMASDD = "06b86afa-9f21-11ec-2ef8-e51b8960cfc5" Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce" -OMAS = "91cfaa06-6526-4804-8666-b540b3feef2f" RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" diff --git a/README.md b/README.md index af5e098..0f3b390 100644 --- a/README.md +++ b/README.md @@ -1,5 +1,28 @@ # GGDUtils ![Format Check](https://github.com/ProjectTorreyPines/GGDUtils.jl/actions/workflows/format_check.yml/badge.svg) +![Docs](https://github.com/ProjectTorreyPines/GGDUtils.jl/actions/workflows/make_docs.yml/badge.svg) +![Tests](https://github.com/ProjectTorreyPines/GGDUtils.jl/actions/workflows/test.yml/badge.svg) -Package holding utilities for Generalized Grid Description (GGD) objects in IMAS datastructure. Primary goals are interpolation and core profile extrapolation. +Package holding utilities for Generalized Grid Description (GGD) objects in IMAS datastructure. Primary goals are interpolation and core profile extrapolation. For installation and usage instructions, see the [online documentation](https://projecttorreypines.github.io/GGDUtils.jl/dev) + +## Building julia environment (installation) + +After cloning this repo, check the make menu: +``` +GGDUtils.jl % make help +Help Menu + +make env_with_cloned_repo (or make r): Creates a Julia environment with the cloned repositories +make env_with_git_url (or make u): Creates a Julia environment with the git urls without creating local clones +make clean: Deletes Project.toml and Manifest.toml for a fresh start +``` + +### make r +This option creates local copies of required private repositories at the same level as current repository and uses them in develop mode to create a Manifest.toml + +### make u +This option uses url of required private repositories to create a static Manifest.toml attached to current master branches of these repositories. + +### make clean +Deletes Manifest.toml so that environment can be recreated, to update or change the last used method. \ No newline at end of file diff --git a/docs/Project.toml b/docs/Project.toml new file mode 100644 index 0000000..b519e7f --- /dev/null +++ b/docs/Project.toml @@ -0,0 +1,3 @@ +[deps] +Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" +RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" diff --git a/docs/make.jl b/docs/make.jl new file mode 100644 index 0000000..1a176b4 --- /dev/null +++ b/docs/make.jl @@ -0,0 +1,23 @@ +# using Pkg +# Pkg.add("Documenter") +# Pkg.add(; url="https://github.com/ProjectTorreyPines/IMASDD.jl.git"); +# Pkg.develop(PackageSpec(; path="../")); +# Pkg.instantiate() +using Documenter +using GGDUtils +using RecipesBase: RecipesBase + +makedocs(; + modules=[GGDUtils], + format=Documenter.HTML(), + sitename="GGDUtils", + checkdocs=:none, +) + +deploydocs(; + repo="github.com/ProjectTorreyPines/GGDUtils.jl.git", + target="build", + branch="gh-pages", + devbranch="dev", + versions=["stable" => "v^", "v#.#"], +) diff --git a/docs/src/index.md b/docs/src/index.md new file mode 100644 index 0000000..d9db9e2 --- /dev/null +++ b/docs/src/index.md @@ -0,0 +1,86 @@ + +# GGDUtils.jl + +```@contents +Pages = ["index.md"] +Depth = 5 +``` + +## Installation + +### Using make: +After cloning this repo, check the make menu: +``` +GGDUtils.jl % make help +Help Menu + +make env_with_cloned_repo (or make r): Creates a Julia environment with the cloned repositories +make env_with_git_url (or make u): Creates a Julia environment with the git urls without creating local clones +make clean: Deletes Project.toml and Manifest.toml for a fresh start +``` + +#### make r +This option creates local copies of required private repositories at the same level as current repository and uses them in develop mode to create a Manifest.toml + +#### make u +This option uses url of required private repositories to create a static Manifest.toml attached to current master branches of these repositories. + +#### make clean +Deletes Manifest.toml so that environment can be recreated, to update or change the last used method. + +### Using Julia REPL and installing using Github url + +Or, in julia REPL: +```julia +julia> using Pkg; +julia> Pkg.add(; url="https://github.com/ProjectTorreyPines/IMASDD.jl.git"); +julia> Pkg.add(; url="https://github.com/ProjectTorreyPines/GGDUtils.jl.git"); +julia> Pkg.instantiate() +``` + +## Interpolations + +Several interpolation functions are available to create interpolaiton functions for data present in a GGD represented over a particular grid subset: + +```@docs +interp +get_TPS_mats +get_kdtree +``` + +## Subset Tools + +```@docs +add_subset_element! +get_subset_space +get_grid_subset +get_subset_boundary_inds +get_subset_boundary +subset_do +get_subset_centers +project_prop_on_subset! +deepcopy_subset +Base.:∈ +get_prop_with_grid_subset_index +``` + +## Types + +```@docs +get_types_with +``` + +## Plot recipes + +Several plot recipes have been defined for easy visualization. +```@docs +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.edge_profiles__grid_ggd___space) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.edge_profiles__grid_ggd___space, ::GGDUtils.IMASDD.edge_profiles__grid_ggd___grid_subset) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.edge_profiles__grid_ggd, ::GGDUtils.IMASDD.IDSvectorElement) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::AbstractVector{<:GGDUtils.IMASDD.edge_profiles__grid_ggd}, ::GGDUtils.IMASDD.IDSvectorElement) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.interferometer) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.interferometer__channel) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.interferometer__channel___line_of_sight) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.interferometer__channel___n_e_line) +RecipesBase.apply_recipe(::Dict{Symbol,Any}, ::GGDUtils.IMASDD.interferometer__channel___n_e_line_average) +``` \ No newline at end of file diff --git a/examples/plotting.ipynb b/examples/plotting.ipynb index eabcd48..5c2f533 100644 --- a/examples/plotting.ipynb +++ b/examples/plotting.ipynb @@ -1,157 +1,163 @@ { - "cells": [ - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Plotting examples using GGDUtils\n", - " \n", - " For running this notebook, you need to install package IJulia in your home environment (that is messy, but that is the only way I know right now). So in your terminal:\n", - " ```\n", - " % julia\n", - " julia > ]\n", - " (@v1.9 pkg) pkg> add IJulia\n", - " ```\n", - "\n", - " After this, Julia kernel would appear in your jupyter notebooks as an option. This also works for julia notebooks directly opened on VSCode. Select the Julia kernel to run this notebook." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "using Pkg\n", - "Pkg.activate(\"./\")\n", - "Pkg.add(url=\"git@github.com:ProjectTorreyPines/OMAS.jl.git\")\n", - "Pkg.develop(path=\"../\")\n", - "Pkg.add(\"Plots\")\n", - "Pkg.add(\"LaTeXStrings\")" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "using OMAS\n", - "using GGDUtils\n", - "using Plots\n", - "using LaTeXStrings" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "dd = OMAS.json2imas(\"../samples/time_dep_edge_profiles_last_step_only.json\")\n", - "grid_ggd = dd.edge_profiles.grid_ggd[1]\n", - "space = grid_ggd.space[1]" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## Plotting grid and subsets" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Choose backend\n", - "gr() # Fast and can save pdf\n", - "# plotlyjs() # Use for interactive plot, can only save png\n", - "\n", - "plot(space) # Simply plot the grid described in space, all common arguments to plot can be given here\n", - "\n", - "# You can overlay any subset by giving a second argument\n", - "# Labels \n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 6), markercolor=:chocolate1)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 7), linecolor=:red, linewidth=2)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 8), linecolor=:darkred, linewidth=2)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 9), linecolor=:limegreen, linewidth=2)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 10), linecolor=:darkgreen, linewidth=2)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 11), linecolor=:cyan, linewidth=2)\n", - "# plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 12), linecolor=:teal, linewidth=1)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 13), linecolor=:royalblue1, linewidth=2)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 14), linecolor=:navyblue, linewidth=2)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 15), linecolor=:fuchsia, linewidth=2, linestyle=:dash)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 16), linecolor=:purple4, linewidth=2, linestyle=:dash)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 101), markershape=:rect, markercolor=:royalblue1)\n", - "# plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 102), markershape=:rect, markercolor=:maroon)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 103), markershape=:diamond, markercolor=:fuchsia)\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 104), markershape=:diamond, markercolor=:purple4)\n", - "\n", - "# Legend is supressed unless asked for specifically\n", - "plot!(legend=true)\n", - "# Default labels are subset.identifier.name but can be changed by providing a label argument" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## Plotting 2D quantities as heatmaps" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Choose backend\n", - "gr() # Fast and can save pdf\n", - "# plotlyjs() # Use for interactive plot, can only save png\n", - "\n", - "n_e = GGDUtils.get_prop_with_grid_subset_index(dd.edge_profiles.ggd[1].electrons.density, 5)\n", - "plot(dd.edge_profiles.grid_ggd, n_e, colorbar_title=\"Electrons density / \" * L\"m^{-3}\")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### You can overlap any grid on top of a quantity" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Choose backend\n", - "gr() # Fast and can save pdf\n", - "# plotlyjs() # Use for interactive plot, can only save png\n", - "\n", - "plot(dd.edge_profiles.grid_ggd, n_e) # Note default label in colorbar\n", - "plot!(space, GGDUtils.get_grid_subset_with_index(grid_ggd, 16), linecolor=:red, linewidth=2, linestyle=:solid, label=\"Separatix\", legend=true)" - ] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Julia 1.9.2", - "language": "julia", - "name": "julia-1.9" - }, - "language_info": { - "file_extension": ".jl", - "mimetype": "application/julia", - "name": "julia", - "version": "1.9.2" - }, - "orig_nbformat": 4 - }, - "nbformat": 4, - "nbformat_minor": 2 + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Plotting examples using GGDUtils\n", + " \n", + " For running this notebook, you need to install package IJulia in your home environment (that is messy, but that is the only way I know right now). So in your terminal:\n", + " ```\n", + " % julia\n", + " julia > ]\n", + " (@v1.9 pkg) pkg> add IJulia\n", + " ```\n", + "\n", + " After this, Julia kernel would appear in your jupyter notebooks as an option. This also works for julia notebooks directly opened on VSCode. Select the Julia kernel to run this notebook." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "using Pkg\n", + "Pkg.activate(\"./\")\n", + "Pkg.add(url=\"git@github.com:ProjectTorreyPines/IMASDD.jl.git\")\n", + "Pkg.develop(path=\"../\")\n", + "Pkg.add(\"Plots\")\n", + "Pkg.add(\"LaTeXStrings\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "using IMASDD\n", + "using GGDUtils\n", + "using Plots\n", + "using LaTeXStrings" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dd = IMASDD.json2imas(\"../samples/time_dep_edge_profiles_last_step_only.json\");\n", + "grid_ggd = dd.edge_profiles.grid_ggd[1];\n", + "space = grid_ggd.space[1];" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plotting grid and subsets" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Choose backend\n", + "gr() # Fast and can save pdf\n", + "# plotlyjs() # Use for interactive plot, can only save png\n", + "\n", + "plot(space) # Simply plot the grid described in space, all common arguments to plot can be given here\n", + "\n", + "# You can overlay any subset by giving a second argument\n", + "# Labels \n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, \"x_points\"), markercolor=:chocolate1)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, \"core_cut\"), linecolor=:red, linewidth=2)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 8), linecolor=:darkred, linewidth=2)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 9), linecolor=:limegreen, linewidth=2)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 10), linecolor=:darkgreen, linewidth=2)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 11), linecolor=:cyan, linewidth=2)\n", + "# plot!(space, GGDUtils.get_grid_subset(grid_ggd, 12), linecolor=:teal, linewidth=1)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 13), linecolor=:royalblue1, linewidth=2)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 14), linecolor=:navyblue, linewidth=2)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 15), linecolor=:fuchsia, linewidth=2, linestyle=:dash)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 16), linecolor=:purple4, linewidth=2, linestyle=:dash)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 101), markershape=:rect, markercolor=:royalblue1)\n", + "# plot!(space, GGDUtils.get_grid_subset(grid_ggd, 102), markershape=:rect, markercolor=:maroon)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 103), markershape=:diamond, markercolor=:fuchsia)\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 104), markershape=:diamond, markercolor=:purple4)\n", + "\n", + "# Legend is supressed unless asked for specifically\n", + "plot!(legend=true)\n", + "# Default labels are subset.identifier.name but can be changed by providing a label argument\n", + "\n", + "# It is recommended to add atleast 10 pt margin to the left to accomodate the y-axis label\n", + "plot!(left_margin=10Plots.pt)\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plotting 2D quantities as heatmaps" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Choose backend\n", + "gr() # Fast and can save pdf\n", + "# plotlyjs() # Use for interactive plot, can only save png\n", + "\n", + "n_e = GGDUtils.get_prop_with_grid_subset_index(dd.edge_profiles.ggd[1].electrons.density, 5)\n", + "plot(dd.edge_profiles.grid_ggd, n_e, colorbar_title=\"Electrons density / \" * L\"m^{-3}\",\n", + " left_margin=10Plots.pt)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### You can overlap any grid on top of a quantity" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Choose backend\n", + "gr() # Fast and can save pdf\n", + "# plotlyjs() # Use for interactive plot, can only save png\n", + "\n", + "plot(dd.edge_profiles.grid_ggd, n_e) # Note default label in colorbar\n", + "plot!(space, GGDUtils.get_grid_subset(grid_ggd, 16), linecolor=:red,\n", + " linewidth=2, linestyle=:solid, label=\"Separatix\", legend=true,\n", + " left_margin=10Plots.pt)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.9.2", + "language": "julia", + "name": "julia-1.9" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.9.2" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 } diff --git a/examples/plotting_interferometer.ipynb b/examples/plotting_interferometer.ipynb index 917d0a5..2df0e0d 100644 --- a/examples/plotting_interferometer.ipynb +++ b/examples/plotting_interferometer.ipynb @@ -24,7 +24,7 @@ "source": [ "using Pkg\n", "Pkg.activate(\"./\")\n", - "Pkg.add(url=\"git@github.com:ProjectTorreyPines/OMAS.jl.git\")\n", + "Pkg.add(url=\"git@github.com:ProjectTorreyPines/IMASDD.jl.git\")\n", "Pkg.develop(path=\"../\")\n", "Pkg.add(\"Plots\")" ] @@ -35,7 +35,7 @@ "metadata": {}, "outputs": [], "source": [ - "using OMAS\n", + "using IMASDD\n", "using GGDUtils\n", "using Plots" ] @@ -53,9 +53,9 @@ "metadata": {}, "outputs": [], "source": [ - "ids = OMAS.json2imas(\"$(@__DIR__)/../samples/time_dep_edge_profiles_with_interferometer.json\")\n", - "grid_ggd = ids.edge_profiles.grid_ggd[1]\n", - "space = grid_ggd.space[1]" + "ids = IMASDD.json2imas(\"$(@__DIR__)/../samples/time_dep_edge_profiles_with_interferometer.json\");\n", + "grid_ggd = ids.edge_profiles.grid_ggd[1];\n", + "space = grid_ggd.space[1];" ] }, { @@ -77,7 +77,7 @@ "\n", "plot(space)\n", "plot!(ids.interferometer) # Default plot_type is :los \n", - "plot!(legend=true)" + "plot!(legend=true, left_margin=10Plots.pt)" ] }, { @@ -99,7 +99,7 @@ "\n", "plot(space)\n", "plot!(ids.interferometer, mirror_length=0.7, linewidth=4, mirror_thickness=0.2)\n", - "plot!(legend=true)" + "plot!(legend=true, left_margin=10Plots.pt)" ] }, { @@ -121,7 +121,7 @@ "\n", "plot(space)\n", "plot!(ids.interferometer, mirror=false)\n", - "plot!(legend=true)" + "plot!(legend=true, left_margin=10Plots.pt)" ] }, { @@ -143,7 +143,7 @@ "\n", "plot(space)\n", "plot!(ids.interferometer.channel[1], label=\"Channel 1\")\n", - "plot!(legend=true)" + "plot!(legend=true, left_margin=10Plots.pt)" ] }, { @@ -166,7 +166,7 @@ "gr() # Fast and can save pdf\n", "# plotlyjs() # Use for interactive plot, can only save png\n", "\n", - "plot(ids.interferometer, plot_type=:n_e)" + "plot(ids.interferometer, plot_type=:n_e, left_margin=10Plots.pt)" ] }, { @@ -179,7 +179,7 @@ "gr() # Fast and can save pdf\n", "# plotlyjs() # Use for interactive plot, can only save png\n", "\n", - "plot(ids.interferometer, plot_type=:n_e_average)" + "plot(ids.interferometer, plot_type=:n_e_average, left_margin=10Plots.pt)" ] }, { @@ -199,7 +199,28 @@ "gr() # Fast and can save pdf\n", "# plotlyjs() # Use for interactive plot, can only save png\n", "\n", - "plot(ids.interferometer.channel[1], plot_type=:n_e_average)" + "plot(ids.interferometer.channel[1], plot_type=:n_e_average, left_margin=10Plots.pt)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plotting the interferometer geometry on top of 2D property data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "gr()\n", + "\n", + "n_e = GGDUtils.get_prop_with_grid_subset_index(ids.edge_profiles.ggd[1].electrons.density, 5)\n", + "plot(ids.edge_profiles.grid_ggd, n_e, colorbar_title=\"Electrons density / m^(-3)\", left_margin=10Plots.pt)\n", + "plot!(space)\n", + "plot!(ids.interferometer, legend=true, size=[635, 900]) # Adding a size comment to make plot aspect ratio better" ] } ], diff --git a/makefile b/makefile new file mode 100644 index 0000000..1124392 --- /dev/null +++ b/makefile @@ -0,0 +1,27 @@ +SHELL := /bin/zsh +help: + @echo "Help Menu" + @echo + @echo "make env_with_cloned_repo (or make r): Creates a Julia environment with the cloned repositories" + @echo "make env_with_git_url (or make u): Creates a Julia environment with the git urls without creating local clones" + @echo "make clean: Deletes Project.toml and Manifest.toml for a fresh start" + @echo + +env_with_cloned_repo r: + @echo "Pulling sample files using dvc" + -dvc pull + @echo "Creating Julia environment by creating local clones of dependent repositories" + @echo "Cloning the repositories and generating Manifest.toml" + -git clone "git@github.com:ProjectTorreyPines/IMASDD.jl.git" ../IMASDD; \ + julia --project=. -e 'using Pkg; Pkg.rm(["IMASDD"]); Pkg.develop(path="../IMASDD"); Pkg.instantiate()' + +env_with_git_url u: + @echo "Pulling sample files using dvc" + -dvc pull + @echo "Creating Julia environment with the git urls without creating local clones" + @echo "Generating Project.toml and Manifest.toml" + julia --project=. -e 'using Pkg; Pkg.rm(["IMASDD"]); Pkg.add(url="git@github.com:ProjectTorreyPines/IMASDD.jl.git", rev="master"); Pkg.instantiate()' + +clean: + @echo "Deleting Manifest.toml" + -rm Manifest.toml diff --git a/samples/time_dep_edge_profiles_last_step_only.json.dvc b/samples/time_dep_edge_profiles_last_step_only.json.dvc index 67f20fd..299351c 100644 --- a/samples/time_dep_edge_profiles_last_step_only.json.dvc +++ b/samples/time_dep_edge_profiles_last_step_only.json.dvc @@ -1,12 +1,12 @@ -md5: 3ea764e4095ca2b8a5359cd76259e26a +md5: cdd59460686d5d4ee901e1f9339f480d frozen: true deps: - path: ITER_Lore_2296_00000/IMAS/time_dep_edge_profiles_last_step_only.json repo: url: git@github.com:ProjectTorreyPines/SOLPSTestSamples.git - rev_lock: 88f14955a52cd4d4fd7e86e2274436a34f4203ec + rev_lock: cf392d9f0292d4f4da4c027f5cc9b2a640fa7fb8 outs: -- md5: db47dee590648770a5ce49c0e0706acf - size: 4087084 +- md5: 3c54b821ce167c41f270999e77a6602c + size: 11793595 hash: md5 path: time_dep_edge_profiles_last_step_only.json diff --git a/samples/time_dep_edge_profiles_with_interferometer.json.dvc b/samples/time_dep_edge_profiles_with_interferometer.json.dvc index c45c30b..72d3c7a 100644 --- a/samples/time_dep_edge_profiles_with_interferometer.json.dvc +++ b/samples/time_dep_edge_profiles_with_interferometer.json.dvc @@ -1,12 +1,12 @@ -md5: fb9792d06bb5d13f3a8803654b78624c +md5: 6aa2ce8a23a582ee8cc1aa3d030070f2 frozen: true deps: - path: ITER_Lore_2296_00000/IMAS/time_dep_edge_profiles_with_interferometer.json repo: url: git@github.com:ProjectTorreyPines/SOLPSTestSamples.git - rev_lock: 80688afd131b17455c3320a4485bbf2b7a55c80c + rev_lock: 790f37427d8e6809449684e4f1d518eab626e780 outs: -- md5: 610272ceabe4a62dd776e0f3c08d2036 - size: 40662828 +- md5: bba50e6afd5cfba436a0b0927a9d5bc9 + size: 41544711 hash: md5 path: time_dep_edge_profiles_with_interferometer.json diff --git a/src/GGDUtils.jl b/src/GGDUtils.jl index 6c7fb25..3837305 100644 --- a/src/GGDUtils.jl +++ b/src/GGDUtils.jl @@ -1,13 +1,9 @@ module GGDUtils -using OMAS: OMAS +using IMASDD: IMASDD const inv_16pi = 1.0 / (16π) -export project_prop_on_subset! -export get_subset_centers -export get_prop_with_grid_subset_index - include("types.jl") include("subset_tools.jl") diff --git a/src/interpolations.jl b/src/interpolations.jl index 1d4f9c2..d9f465f 100644 --- a/src/interpolations.jl +++ b/src/interpolations.jl @@ -7,7 +7,12 @@ export interp export get_kdtree export get_TPS_mats -function get_kdtree(space::OMAS.edge_profiles__grid_ggd___space) +""" + get_kdtree(space::IMASDD.edge_profiles__grid_ggd___space) + +Get a KDTree for all the cells in the space for search for nearest neighbours. +""" +function get_kdtree(space::IMASDD.edge_profiles__grid_ggd___space) grid_nodes = space.objects_per_dimension[1].object grid_faces = space.objects_per_dimension[3].object grid_faces = [cell for cell ∈ grid_faces if length(cell.nodes) == 4] @@ -18,9 +23,17 @@ function get_kdtree(space::OMAS.edge_profiles__grid_ggd___space) return KDTree(grid_centers; leafsize=10) end +""" + get_kdtree( + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + ) + +Get a KDTree for a subset of the space for search for nearest neighbours. +""" function get_kdtree( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) subset_centers = get_subset_centers(space, subset) return KDTree([SVector{2}(sc) for sc ∈ subset_centers]; leafsize=10) @@ -28,12 +41,11 @@ end """ interp( - prop_values::Vector{T}, - kdtree::KDTree; - use_nearest_n::Int=4, - weighing::Function=(d) -> 1 / d, - -) where {T <: Real} + prop_values::Vector{T}, + kdtree::KDTree; + use_nearest_n::Int=4, + weighing::Function=(d) -> 1 / d, + ) where {T <: Real} Lowest level interpolation function. It takes a vector of property values and a KDTree defined over a 2D space with the same number of nodes as the property values. It returns @@ -129,10 +141,9 @@ end """ interp( - y::Vector{T}, - TPS_mats::Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, - -) where {T <: Real, U <: Real} + y::Vector{T}, + TPS_mats::Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, + ) where {T <: Real, U <: Real} Lowest level function for Thin Plate Spline method """ @@ -167,17 +178,16 @@ function interp(y::Vector{T}, x::Vector{Tuple{U, U}}) where {T <: Real, U <: Rea return interp(y, get_TPS_mats(x)) end -function get_TPS_mats(space::OMAS.edge_profiles__grid_ggd___space) +function get_TPS_mats(space::IMASDD.edge_profiles__grid_ggd___space) nodes = [Tuple(node.geometry) for node ∈ space.objects_per_dimension[1].object] return get_TPS_mats(nodes) end """ interp( - prop_values::Vector{T}, - space::OMAS.edge_profiles__grid_ggd___space - -) where {T <: Real} + prop_values::Vector{T}, + space::IMASDD.edge_profiles__grid_ggd___space + ) where {T <: Real} If the whole space is provided instead of a kdtree, calculate the kdtree for whole space. Again, here it is assumed that the property values are porvided for each node @@ -185,77 +195,80 @@ of the space. """ function interp( prop_values::Vector{T}, - space::OMAS.edge_profiles__grid_ggd___space, + space::IMASDD.edge_profiles__grid_ggd___space, ) where {T <: Real} return interp(prop_values, get_TPS_mats(space)) end function get_TPS_mats( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) return get_TPS_mats(get_subset_centers(space, subset)) end """ interp( - prop_values::Vector{Real}, - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset - -) + prop_values::Vector{Real}, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset + ) If a subset of the space is provided, calculate the kdtree for the subset. In this case it is assumed that the property values are provided for each element of the subset. """ function interp( prop_values::Vector{T}, - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) where {T <: Real} return interp(prop_values, get_TPS_mats(space, subset)) end """ interp( - prop::edge_profiles__prop_on_subset, - grid_ggd::OMAS.edge_profiles__grid_ggd, - value_field::Symbol=:values - -) + prop::edge_profiles__prop_on_subset, + grid_ggd::IMASDD.edge_profiles__grid_ggd, + value_field::Symbol=:values + ) Example: + +```julia grid_ggd = dd.edge_profiles.grid_ggd[1] get_electron_density = interp(dd.edge_profiles.ggd[1].electrons.density[1], grid_ggd) get_e_field_par = interp(dd.edge_profiles.ggd[1].e_field[1], grid_ggd, :parallel) +``` """ function interp( prop::edge_profiles__prop_on_subset, - grid_ggd::OMAS.edge_profiles__grid_ggd, + grid_ggd::IMASDD.edge_profiles__grid_ggd, value_field::Symbol=:values, ) - subset = get_grid_subset_with_index(grid_ggd, prop.grid_subset_index) + subset = get_grid_subset(grid_ggd, prop.grid_subset_index) space = grid_ggd.space[subset.element[1].object[1].space] return interp(getfield(prop, value_field), space, subset) end """ interp( - prop_arr::Vector{T}, - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, - value_field::Symbol=:values - -) where {T <: edge_profiles__prop_on_subset} + prop_arr::AbstractVector{T}, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + value_field::Symbol=:values + ) where {T <: edge_profiles__prop_on_subset} Example: + +```julia sol = get_grid_subset_with_index(dd.edge_profiles.grid_ggd[1], 23) get_electron_density = interp(dd.edge_profiles.ggd[1].electrons.density, space, sol) +``` """ function interp( - prop_arr::Vector{T}, - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + prop_arr::AbstractVector{T}, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, value_field::Symbol=:values, ) where {T <: edge_profiles__prop_on_subset} prop = get_prop_with_grid_subset_index(prop_arr, subset.identifier.index) @@ -264,59 +277,61 @@ end """ interp( - prop_arr::Vector{T}, - grid_ggd::OMAS.edge_profiles__grid_ggd, - grid_subset_index::Int, - value_field::Symbol=:values - -) where {T <: edge_profiles__prop_on_subset} + prop_arr::AbstractVector{T}, + grid_ggd::IMASDD.edge_profiles__grid_ggd, + grid_subset_index::Int, + value_field::Symbol=:values + ) where {T <: edge_profiles__prop_on_subset} Example: + +```julia get_n_e_sep = interp(dd.edge_profiles.ggd[1].electrons.density, grid_ggd, 16) +``` """ function interp( - prop_arr::Vector{T}, - grid_ggd::OMAS.edge_profiles__grid_ggd, + prop_arr::AbstractVector{T}, + grid_ggd::IMASDD.edge_profiles__grid_ggd, grid_subset_index::Int, value_field::Symbol=:values, ) where {T <: edge_profiles__prop_on_subset} prop = get_prop_with_grid_subset_index(prop_arr, grid_subset_index) - subset = get_grid_subset_with_index(grid_ggd, grid_subset_index) + subset = get_grid_subset(grid_ggd, grid_subset_index) space = grid_ggd.space[subset.element[1].object[1].space] return interp(getfield(prop, value_field), space, subset) end -function get_TPS_mats(grid_ggd::OMAS.edge_profiles__grid_ggd, grid_subset_index::Int) - subset = get_grid_subset_with_index(grid_ggd, grid_subset_index) +function get_TPS_mats(grid_ggd::IMASDD.edge_profiles__grid_ggd, grid_subset_index::Int) + subset = get_grid_subset(grid_ggd, grid_subset_index) space = grid_ggd.space[subset.element[1].object[1].space] return get_TPS_mats(space, subset) end -#! format off """ interp( - prop_arr::Vector{T}, - TPS_mats::Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, - grid_subset_index::Int, - value_field::Symbol=:values, - -) where {T <: edge_profiles__prop_on_subset} + prop_arr::AbstractVector{T}, + TPS_mats::Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, + grid_subset_index::Int, + value_field::Val{V}=Val(:values), + ) where {T <: edge_profiles__prop_on_subset, U <: Real, V} Same use case as above but allows one to reuse previously calculated TPS matrices. Example: -TPS_mat = get_TPS_mats(dd.edge_profiles.grid_ggd[1], 5) + +```julia +TPS_mat = get_TPS_mats(dd.edge_profiles.grid_ggd[1], 5); for it ∈ eachindex(dd.edge_profiles.ggd) get_n_e = interp(dd.edge_profiles.ggd[it].electrons.density, TPS_mat_sep, 5) - println("This time step has n_e at (0, 0) = ", get_n_e(, 0)) + println("This time step has n_e at (0, 0) = ", get_n_e(0, 0)) end +``` -This will run faster has heavy matrix calculations will happen only once. +This will run faster as heavy matrix calculations will happen only once. """ -#! format on function interp( - prop_arr::Vector{T}, + prop_arr::AbstractVector{T}, TPS_mats::Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, grid_subset_index::Int, value_field::Val{V}=Val(:values), @@ -330,16 +345,19 @@ const RHO_EXT_POS = [1.0001, 1.1, 5] const RHO_EXT_NEG = [-5, -0.0001] # I guess this would be at a PF coil or something? """ - interp(eqt::OMAS.equilibrium__time_slice) + interp(eqt::IMASDD.equilibrium__time_slice) For a given equilibrium time slice, return a function that can be used to interpolate from (r, z) space to rho (normalized toroidal flux coordinate)space. Example: + +```julia rz2rho = interp(dd.equilibrium.time_slice[1]) -rho = rz2rho.([(r, z) for r in 3:0.01:9, for z in -5:0.01:5]) +rho = rz2rho.([(r, z) for r ∈ 3.0:0.01:9.0, z ∈ -5.0:0.01:5.0]) +``` """ -function interp(eqt::OMAS.equilibrium__time_slice) +function interp(eqt::IMASDD.equilibrium__time_slice) p1 = eqt.profiles_1d p2 = eqt.profiles_2d[1] gq = eqt.global_quantities @@ -369,22 +387,24 @@ end """ interp( - prop::Vector{T}, - prof::OMAS.core_profiles__profiles_1d, - -) where {T <: Real} + prop::Vector{T}, + prof::IMASDD.core_profiles__profiles_1d, + ) where {T <: Real} Returns an inteprolation function for the core profile property values defined on normalized toroidal flux coordinate rho. Example: + +```julia core_profile_n_e = dd.core_profiles.profiles_1d[1].electrons.density get_n_e = interp(core_profile_n_e, dd.core_profiles.profiles_1d[1]) get_n_e(1) # Returns electron density at rho = 1 (separatix) +``` """ function interp( prop::Vector{T}, - prof::OMAS.core_profiles__profiles_1d, + prof::IMASDD.core_profiles__profiles_1d, ) where {T <: Real} rho_prof = copy(prof.grid.rho_tor_norm) length(prop) == length(rho_prof) || @@ -399,11 +419,10 @@ end """ interp( - prop::Vector{T}, - prof::OMAS.core_profiles__profiles_1d, - rz2rho::Function, - -) + prop::Vector{T}, + prof::IMASDD.core_profiles__profiles_1d, + rz2rho::Function, + ) Returns an inteprolation function in (R, Z) domain for the core profile property values defined on normalized toroidal flux coordinate rho and with a provided function to @@ -411,14 +430,16 @@ convert (R,Z) to rho. Example: +```julia rz2rho = interp(dd.equilibrium.time_slice[1]) core_profile_n_e = dd.core_profiles.profiles_1d[1].electrons.density get_n_e = interp(core_profile_n_e, dd.core_profiles.profiles_1d[1], rz2rho) get_n_e(5.0, 3.5) # Returns electron density at (R, Z) = (5.0, 3.5) +``` """ function interp( prop::Vector{T}, - prof::OMAS.core_profiles__profiles_1d, + prof::IMASDD.core_profiles__profiles_1d, rz2rho::Function, ) where {T <: Real} itp = interp(prop, prof) @@ -430,11 +451,10 @@ end """ interp( - prop::Vector{T}, - prof::OMAS.core_profiles__profiles_1d, - eqt::OMAS.equilibrium__time_slice, - -) where {T <: Real} + prop::Vector{T}, + prof::IMASDD.core_profiles__profiles_1d, + eqt::IMASDD.equilibrium__time_slice, + ) where {T <: Real} Returns an inteprolation function in (R, Z) domain for the core profile property values defined on normalized toroidal flux coordinate rho and with a provided equilibrium time @@ -442,15 +462,17 @@ slice to get (R, Z) to rho conversion. Example: +```julia eqt = dd.equilibrium.time_slice[1] core_profile_n_e = dd.core_profiles.profiles_1d[1].electrons.density get_n_e = interp(core_profile_n_e, dd.core_profiles.profiles_1d[1], eqt) get_n_e(5.0, 3.5) # Returns electron density at (R, Z) = (5.0, 3.5) +``` """ function interp( prop::Vector{T}, - prof::OMAS.core_profiles__profiles_1d, - eqt::OMAS.equilibrium__time_slice, + prof::IMASDD.core_profiles__profiles_1d, + eqt::IMASDD.equilibrium__time_slice, ) where {T <: Real} rz2rho = interp(eqt) return interp(prop, prof, rz2rho) diff --git a/src/recipes.jl b/src/recipes.jl index 3e4ef79..c85140b 100644 --- a/src/recipes.jl +++ b/src/recipes.jl @@ -2,10 +2,17 @@ using RecipesBase using ColorSchemes: ColorSchemes import Statistics: norm, dot -@recipe function f(space::OMAS.edge_profiles__grid_ggd___space) +""" + plot(space::IMASDD.edge_profiles__grid_ggd___space) + +Plot the grid_ggd space object. Defaults to size of [600, 900] and linecolor of :black, +linewidth of 0.2, and no legend. +""" +@recipe function f(space::IMASDD.edge_profiles__grid_ggd___space) nodes = space.objects_per_dimension[1].object edges = space.objects_per_dimension[2].object legend --> false + subplot --> 1 linewidth --> 0.2 linecolor --> :black size --> [600, 900] @@ -32,14 +39,24 @@ import Statistics: norm, dot end end +""" + plot( + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + ) + +Plot the a subset of a space. Defaults to size of [600, 900] and linecolor of :black, +linewidth of 0.2, and no legend. +""" @recipe function f( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) nodes = space.objects_per_dimension[1].object edges = space.objects_per_dimension[2].object cells = space.objects_per_dimension[3].object legend --> false + subplot --> 1 linewidth --> 0.2 linecolor --> :black size --> [600, 900] @@ -49,8 +66,8 @@ end label_assigned = false if subset.element[1].object[1].dimension == 3 for ele ∈ subset.element - for bnd_ind ∈ cells[ele.object[1].index].boundary - union!(subset_edge_inds, bnd_ind) + for bnd ∈ cells[ele.object[1].index].boundary + union!(subset_edge_inds, bnd.index) end end elseif subset.element[1].object[1].dimension == 2 @@ -97,13 +114,28 @@ end end end -@recipe function f(grid_ggd::OMAS.edge_profiles__grid_ggd, prop::OMAS.IDSvectorElement) - subset = get_grid_subset_with_index(grid_ggd, prop.grid_subset_index) +""" + plot( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + prop::IMASDD.IDSvectorElement, + ) + +Plot 2D heatmap of edge_profiles_ggd property on a grid_ggd space object. Defaults to +size of [635, 900], xaxis of "R / m", yaxis of "Z / m", and no legend. If :seriescolor +is not provided, :inferno color scheme is used. If :colorbar_title is not provided, the +property name is used. This function creates a plot with layout [a{0.95w} b] where a +is the heatmap and b is the colorbar. +""" +@recipe function f( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + prop::IMASDD.IDSvectorElement, +) + subset = get_grid_subset(grid_ggd, prop.grid_subset_index) space = grid_ggd.space[subset.element[1].object[1].space] nodes = space.objects_per_dimension[1].object cells = space.objects_per_dimension[3].object legend --> false - size --> [600, 900] + size --> [635, 900] xaxis --> "R / m" yaxis --> "Z / m" layout := @layout [a{0.95w} b] @@ -157,9 +189,21 @@ end end end +""" + plot( + grid_ggd_arr::AbstractVector{<:IMASDD.edge_profiles__grid_ggd}, + prop::IMASDD.IDSvectorElement, + ) + +Plot 2D heatmap of edge_profiles_ggd property on a grid_ggd space object. Defaults to +size of [635, 900], xaxis of "R / m", yaxis of "Z / m", and no legend. If :seriescolor +is not provided, :inferno color scheme is used. If :colorbar_title is not provided, the +property name is used. This function creates a plot with layout [a{0.95w} b] where a +is the heatmap and b is the colorbar. +""" @recipe function f( - grid_ggd_arr::Vector{OMAS.edge_profiles__grid_ggd}, - prop::OMAS.IDSvectorElement, + grid_ggd_arr::AbstractVector{<:IMASDD.edge_profiles__grid_ggd}, + prop::IMASDD.IDSvectorElement, ) found = false for grid_ggd ∈ grid_ggd_arr @@ -177,8 +221,23 @@ end end end +""" + plot( + ifo::IMASDD.interferometer, + ) + +Plot all the channels of interferometer object. +Optional keywords: + + - :plot_type: :los(default), :n_e, or :n_e_average. :los plots the line of sight of + the channel in a 2D plot, :n_e plots the integrated n_e along the line of sight vs + time, and :n_e_average plots the average n_e vs time. + - :mirror: true(default) or false. + - :mirror_length: 0.5(default). + - :mirror_thickness: 0.1(default). +""" @recipe f( - ifo::OMAS.interferometer, + ifo::IMASDD.interferometer, ) = for ch ∈ ifo.channel @series begin @@ -186,8 +245,23 @@ end end end +""" + plot( + ifo_ch::IMASDD.interferometer__channel, + ) + +Plot individual channel of interferometer. +Optional keywords: + + - :plot_type: :los(default), :n_e, or :n_e_average. :los plots the line of sight of + the channel in a 2D plot, :n_e plots the integrated n_e along the line of sight vs + time, and :n_e_average plots the average n_e vs time. + - :mirror: true(default) or false. + - :mirror_length: 0.5(default). + - :mirror_thickness: 0.1(default). +""" @recipe function f( - ifo_ch::OMAS.interferometer__channel, + ifo_ch::IMASDD.interferometer__channel, ) if :plot_type ∈ keys(plotattributes) plot_type = plotattributes[:plot_type] @@ -216,16 +290,37 @@ end end end +""" + plot( + ifo_ch_los::IMASDD.interferometer__channel___line_of_sight, + ) + +Plot line of sight of a channel of interferometer. + +Default plot settings: + + - subplot: 1 + - size: [600, 900] + - xaxis: "R / m" + - yaxis: "Z / m" + +Optional keywords: + + - :mirror: true(default) or false. + - :mirror_length: 0.5(default). + - :mirror_thickness: 0.1(default). +""" @recipe function f( - ifo_ch_los::OMAS.interferometer__channel___line_of_sight, + ifo_ch_los::IMASDD.interferometer__channel___line_of_sight, ) + subplot --> 1 size --> [600, 900] xaxis --> "R / m" yaxis --> "Z / m" fp = ifo_ch_los.first_point sp = ifo_ch_los.second_point tp = ifo_ch_los.third_point - if tp == OMAS.interferometer__channel___line_of_sight__third_point() + if tp == IMASDD.interferometer__channel___line_of_sight__third_point() tp = fp end @@ -295,16 +390,34 @@ end end end +""" + plot( + ifo_ch_n_e_line::IMASDD.interferometer__channel___n_e_line, + ) + +Plot line integrated electron density of a channel of interferometer. + +Default plot settings: + + - subplot: 1 + - xaxis: "time / s" + - yaxis: "Integerated n_e / m^-2" + +Optional keywords: + + - :average: true or false(default). If true, plot the average n_e vs time. +""" @recipe function f( - ifo_ch_n_e_line::OMAS.interferometer__channel___n_e_line, + ifo_ch_n_e_line::IMASDD.interferometer__channel___n_e_line, ) if :average ∈ keys(plotattributes) @series begin ifo_ch.n_e_line_average end else + subplot --> 1 xaxis --> "time / s" - yaxis --> "Integrrated n_e / m^-2" + yaxis --> "Integerated n_e / m^-2" @series begin seriestype := :path linewidth --> 2 @@ -313,9 +426,23 @@ end end end +""" + plot( + ifo_ch_n_e_line_average::IMASDD.interferometer__channel___n_e_line_average, + ) + +Plot average electron density of a channel of interferometer. + +Default plot settings: + + - subplot: 1 + - xaxis: "time / s" + - yaxis: "Average n_e / m^-3" +""" @recipe function f( - ifo_ch_n_e_line_average::OMAS.interferometer__channel___n_e_line_average, + ifo_ch_n_e_line_average::IMASDD.interferometer__channel___n_e_line_average, ) + subplot --> 1 xaxis --> "time / s" yaxis --> "Average n_e / m^-3" @series begin diff --git a/src/subset_tools.jl b/src/subset_tools.jl index 1638574..c755b18 100644 --- a/src/subset_tools.jl +++ b/src/subset_tools.jl @@ -1,19 +1,29 @@ +export add_subset_element! +export get_subset_space +export get_grid_subset +export get_subset_boundary_inds +export get_subset_boundary +export subset_do +export get_subset_centers +export project_prop_on_subset! +export deepcopy_subset +export get_prop_with_grid_subset_index + """ add_subset_element!( - subset::OMAS.edge_profiles__grid_ggd___grid_subset, - sn::Int, - dim::Int, - index::Int, - in_subset=(x...) -> true; - kwargs..., - -) + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + sn::Int, + dim::Int, + index::Int, + in_subset=(x...) -> true; + kwargs..., + ) Adds a new element to gird_subset with properties space number (sn), dimension (dim), and index (index). The element is added only if the function in_subset returns true. """ function add_subset_element!( - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, sn::Int, dim::Int, index::Int, @@ -32,19 +42,18 @@ end """ add_subset_element!( - subset, - sn, - dim, - index::Vector{Int}, - in_subset=(x...) -> true; - kwargs..., - -) + subset, + sn, + dim, + index::Vector{Int}, + in_subset=(x...) -> true; + kwargs..., + ) Overloaded to work differently (faster) with list of indices to be added. """ function add_subset_element!( - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, sn::Int, dim::Int, index::Vector{Int}, @@ -65,30 +74,33 @@ function add_subset_element!( end """ - get_subset_space(space::OMAS.edge_profiles__grid_ggd___space, - elements::Vector{OMAS.edge_profiles__grid_ggd___grid_subset___element}) + get_subset_space( + space::IMASDD.edge_profiles__grid_ggd___space, + elements::AbstractVector{<:IMASDD.edge_profiles__grid_ggd___grid_subset___element} + ) Returns an array of space object indices corresponding to the correct objects_per_dimension (nodes, edges or cells) for the subset elements. """ -function get_subset_space(space::OMAS.edge_profiles__grid_ggd___space, - elements::Vector{OMAS.edge_profiles__grid_ggd___grid_subset___element}) +function get_subset_space( + space::IMASDD.edge_profiles__grid_ggd___space, + elements::AbstractVector{<:IMASDD.edge_profiles__grid_ggd___grid_subset___element}, +) nD = elements[1].object[1].dimension nD_objects = space.objects_per_dimension[nD].object return [nD_objects[ele.object[1].index] for ele ∈ elements] end """ - get_grid_subset_with_index( - grid_ggd::OMAS.edge_profiles__grid_ggd, - grid_subset_index::Int, - -) + get_grid_subset( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + grid_subset_index::Int, + ) Returns the grid_subset in a grid_ggd with the matching grid_subset_index """ -function get_grid_subset_with_index( - grid_ggd::OMAS.edge_profiles__grid_ggd, +function get_grid_subset( + grid_ggd::IMASDD.edge_profiles__grid_ggd, grid_subset_index::Int, ) for subset ∈ grid_ggd.grid_subset @@ -96,19 +108,36 @@ function get_grid_subset_with_index( return subset end end - # BCL 12/8: Creates type instability, but maybe okay since it's a "simple" Union - # subset::Union{Int64, OMAS.edge_profiles__grid_ggd___grid_subset} - # - # Better would be to immediately throw an error or return nothing - return 0 # Indicates failure + + return error("Subset ", grid_subset_index, " not found.") end """ - get_subset_boundary_inds( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + get_grid_subset( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + grid_subset_name::String, + ) +Returns the grid_subset in a grid_ggd with the matching grid_subset_name +""" +function get_grid_subset( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + grid_subset_name::String, ) + for subset ∈ grid_ggd.grid_subset + if subset.identifier.name == grid_subset_name + return subset + end + end + + return error("Subset ", grid_subset_name, " not found.") +end + +""" + get_subset_boundary_inds( + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + ) Returns an array of space object indices corresponding to the boundary of the subset. That means, it returns indices of nodes that are at the end of open edge subset or @@ -116,8 +145,8 @@ it returns the indices of edges that are the the boundary of a cell subset. Returns an empty array if the subset is 1D (nodes). """ function get_subset_boundary_inds( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) nD = subset.element[1].object[1].dimension if nD > 1 # Only 2D (edges) and 3D (cells) subsets have boundaries @@ -134,19 +163,18 @@ end """ get_subset_boundary( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, - -) + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + ) Returns an array of elements of grid_subset generated from the boundary of the subset provided. The dimension of these elments is reduced by 1. """ function get_subset_boundary( - space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) - ret_subset = OMAS.edge_profiles__grid_ggd___grid_subset() + ret_subset = IMASDD.edge_profiles__grid_ggd___grid_subset() boundary_inds = get_subset_boundary_inds(space, subset) bnd_dim = subset.element[1].object[1].dimension - 1 space_number = subset.element[1].object[1].space @@ -155,22 +183,30 @@ function get_subset_boundary( end """ - subset_do(set_operator, - itrs::Vararg{Vector{OMAS.edge_profiles__grid_ggd___grid_subset___element}}; - space::OMAS.edge_profiles__grid_ggd___space=OMAS.edge_profiles__grid_ggd___space(), - use_nodes=false) + subset_do( + set_operator, + itrs::Vararg{ + AbstractVector{<:IMASDD.edge_profiles__grid_ggd___grid_subset___element}, + }; + space::IMASDD.edge_profiles__grid_ggd___space=IMASDD.edge_profiles__grid_ggd___space(), + use_nodes=false + ) Function to perform any set operation (intersect, union, setdiff etc.) on subset.element to generate a list of elements to go to subset object. If use_nodes is true, the set operation will be applied on the set of nodes from subset.element, space argument is required for this. Note: that the arguments are subset.element (not the subset itself). Similarly, the -return object is a list of OMAS.edge_profiles__grid_ggd___grid_subset___element. +return object is a list of IMASDD.edge_profiles__grid_ggd___grid_subset___element. """ -function subset_do(set_operator, - itrs::Vararg{Vector{OMAS.edge_profiles__grid_ggd___grid_subset___element}}; - space::OMAS.edge_profiles__grid_ggd___space=OMAS.edge_profiles__grid_ggd___space(), - use_nodes=false) +function subset_do( + set_operator, + itrs::Vararg{ + AbstractVector{<:IMASDD.edge_profiles__grid_ggd___grid_subset___element}, + }; + space::IMASDD.edge_profiles__grid_ggd___space=IMASDD.edge_profiles__grid_ggd___space(), + use_nodes=false, +) if use_nodes ele_inds = set_operator( [ @@ -187,23 +223,29 @@ function subset_do(set_operator, ) dim = itrs[1][1].object[1].dimension end - ret_subset = OMAS.edge_profiles__grid_ggd___grid_subset() + ret_subset = IMASDD.edge_profiles__grid_ggd___grid_subset() space_number = itrs[1][1].object[1].space add_subset_element!(ret_subset, space_number, dim, ele_inds) return ret_subset.element end """ - get_subset_centers(space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset) + get_subset_centers( + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + ) Returns an array of tuples corresponding to (r,z) coordinates of the center of cells or the center of edges in the subset space. """ -function get_subset_centers(space::OMAS.edge_profiles__grid_ggd___space, - subset::OMAS.edge_profiles__grid_ggd___grid_subset, +function get_subset_centers( + space::IMASDD.edge_profiles__grid_ggd___space, + subset::IMASDD.edge_profiles__grid_ggd___grid_subset, ) subset_space = get_subset_space(space, subset.element) + if subset.element[1].object[1].dimension == 1 + return [Tuple(obj.geometry) for obj ∈ subset_space] + end grid_nodes = space.objects_per_dimension[1].object return [ Tuple(mean(SVector{2}(grid_nodes[node].geometry) for node ∈ obj.nodes)) for @@ -211,56 +253,67 @@ function get_subset_centers(space::OMAS.edge_profiles__grid_ggd___space, ] end -#! format: off """ - project_prop_on_subset!(prop_arr::Vector{T}, - from_subset::OMAS.edge_profiles__grid_ggd___grid_subset, - to_subset::OMAS.edge_profiles__grid_ggd___grid_subset, - space::OMAS.edge_profiles__grid_ggd___space, -) + project_prop_on_subset!( + prop_arr::AbstractVector{T}, + from_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + to_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, + value_field::Symbol=:values; + TPS_mats::Union{ + Nothing, + Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, + }=nothing, + ) where {T <: edge_profiles__prop_on_subset, U <: Real} This function can be used to add another instance on a property vector representing the value in a new subset that can be taken as a projection from an existing larger subset. -Arguments: -prop: A property like electrons.density that is a vector of objects with fields - coefficients, grid_index, grid_subset_index, and values. The different instances - in the vector correspond to different grid_subset for which the property is - provided. -from_subset: grid_subset object which is already represented in the property instance. - grid subset with index 5 is populated in electrons.density already if the - values for all cells are present. -to_subset: grid_subset which is either a smaller part of from_subset (core, sol, idr, - odr) but has same dimensions as from_subset - OR - is smaller in dimension that goes through the from_subset (core_boundary, - separatix etc.) -space: (optional) space object in grid_ggd is required only when from_subset is - higher dimensional than to_subset. + +Input Arguments: + + - prop: A property like electrons.density that is a vector of objects with fields + coefficients, grid_index, grid_subset_index, and values. The different instances + in the vector correspond to different grid_subset for which the property is + provided. + - from_subset: grid_subset object which is already represented in the property instance. + grid subset with index 5 is populated in electrons.density already if the + values for all cells are present. + - to_subset: grid_subset which is either a smaller part of from_subset (core, sol, idr, + odr) but has same dimensions as from_subset + OR + is smaller in dimension that goes through the from_subset (core_boundary, + separatix etc.) + - space: (optional) space object in grid_ggd is required only when from_subset is + higher dimensional than to_subset. + Returns: NOTE: This function ends in ! which means it updates prop argument in place. But for the additional utility, this function also returns a tuple (to_subset_centers, to_prop_values) when from_subset dimension is greater than - to_subset dimension +to_subset dimension OR (to_subset_ele_obj_inds, to_prop_values) when from_subset dimension is same as - to_subset dimension) +to_subset dimension) + Descriptions: to_subset_centers: center of cells or center of edges of the to_subset where property - values are defined and stored +values are defined and stored to_subset_ele_obj_inds: Indices of the elements of to_subset where property values are - defined and stored +defined and stored to_prop_values: The projected values of the properties added to prop object in a new - instance +instance """ -#! format: on -function project_prop_on_subset!(prop_arr::Vector{T}, - from_subset::OMAS.edge_profiles__grid_ggd___grid_subset, - to_subset::OMAS.edge_profiles__grid_ggd___grid_subset, - space::OMAS.edge_profiles__grid_ggd___space, +function project_prop_on_subset!( + prop_arr::AbstractVector{T}, + from_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + to_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + space::IMASDD.edge_profiles__grid_ggd___space, value_field::Symbol=:values; - interp_method::Symbol=:thin_plate_spline, - interp_kwargs=Dict(), -) where {T <: edge_profiles__prop_on_subset} + TPS_mats::Union{ + Nothing, + Tuple{Matrix{U}, Matrix{U}, Matrix{U}, Vector{Tuple{U, U}}}, + }=nothing, +) where {T <: edge_profiles__prop_on_subset, U <: Real} if from_subset.element[1].object[1].dimension == to_subset.element[1].object[1].dimension return project_prop_on_subset!(prop, from_subset, to_subset) @@ -279,26 +332,36 @@ function project_prop_on_subset!(prop_arr::Vector{T}, to_prop_values = getfield(to_prop, value_field) from_prop_values = getfield(from_prop, value_field) resize!(to_prop_values, length(to_subset.element)) - if interp_method == :thin_plate_spline + if isnothing(TPS_mats) prop_interp = interp(prop_arr, space, from_subset) - elseif interp_method == :KDTree - prop_interp = interp( - from_prop_values, - get_kdtree(space, from_subset; interp_kwargs...), - ) else - error("Supported interpolation methods are :thin_plate_spline and :KDTree") + prop_interp = interp(from_prop_values, TPS_mats) end to_prop_values = prop_interp.(to_subset_centers) + setproperty!(to_prop, value_field, to_prop_values) return to_subset_centers, to_prop_values else error("to_subset is higher dimensional than from_subset") end end -function project_prop_on_subset!(prop_arr::Vector{T}, - from_subset::OMAS.edge_profiles__grid_ggd___grid_subset, - to_subset::OMAS.edge_profiles__grid_ggd___grid_subset, +""" + project_prop_on_subset!( + prop_arr::AbstractVector{T}, + from_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + to_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + value_field::Symbol=:values, + ) where {T <: edge_profiles__prop_on_subset} + +If the dimensions of from_subset and to_subset are the same, this function can be used +to add another instance on a property vector representing the value in to_subset without +any interpolation or use of space object. The function returns a tuple of indices of +elements of to_subset and the values of the property in to_subset. +""" +function project_prop_on_subset!( + prop_arr::AbstractVector{T}, + from_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, + to_subset::IMASDD.edge_profiles__grid_ggd___grid_subset, value_field::Symbol=:values, ) where {T <: edge_profiles__prop_on_subset} from_prop = get_prop_with_grid_subset_index(prop_arr, from_subset.identifier.index) @@ -344,14 +407,89 @@ function project_prop_on_subset!(prop_arr::Vector{T}, end """ - Base.:∈( - point::Tuple{Real, Real}, - subset_of_space::Tuple{ - OMAS.edge_profiles__grid_ggd___grid_subset, - OMAS.edge_profiles__grid_ggd___space, - }, + deepcopy_subset(subset::IMASDD.edge_profiles__grid_ggd___grid_subset) -) +Faster deepcopy function for grid_subset object. This function is used to create a deep +copy of a grid_subset object bypassing several checks performed by IMASDD. +""" +function deepcopy_subset(subset::IMASDD.edge_profiles__grid_ggd___grid_subset) + new_subset = IMASDD.edge_profiles__grid_ggd___grid_subset() + + base = getfield(subset, :base) + new_base = getfield(new_subset, :base) + resize!(new_base, length(base)) + for (ii, b) ∈ enumerate(base) + if !IMASDD.ismissing(b, :jacobian) + new_base[ii].jacobian = deepcopy(getfield(b, :jacobian)) + end + if !IMASDD.ismissing(b, :tensor_contravariant) + new_base[ii].tensor_contravariant = + deepcopy(getfield(b, :tensor_contravariant)) + end + if !IMASDD.ismissing(b, :tensor_covariant) + new_base[ii].tensor_covariant = deepcopy(getfield(b, :tensor_covariant)) + end + end + + if !IMASDD.ismissing(subset, :dimension) + new_subset.dimension = deepcopy(getfield(subset, :dimension)) + end + + element = getfield(subset, :element) + new_element = getfield(new_subset, :element) + resize!(new_element, length(element)) + for (ii, ele) ∈ enumerate(element) + object = getfield(ele, :object) + new_object = getfield(new_element[ii], :object) + resize!(new_object, length(object)) + for (jj, obj) ∈ enumerate(object) + if !IMASDD.ismissing(obj, :dimension) + new_object[jj].dimension = deepcopy(getfield(obj, :dimension)) + end + if !IMASDD.ismissing(obj, :index) + new_object[jj].index = deepcopy(getfield(obj, :index)) + end + if !IMASDD.ismissing(obj, :space) + new_object[jj].space = deepcopy(getfield(obj, :space)) + end + end + end + + identifier = getfield(subset, :identifier) + new_identifier = getfield(new_subset, :identifier) + if !IMASDD.ismissing(identifier, :index) + new_identifier.index = deepcopy(getfield(identifier, :index)) + end + if !IMASDD.ismissing(identifier, :name) + new_identifier.name = deepcopy(getfield(identifier, :name)) + end + if !IMASDD.ismissing(identifier, :description) + new_identifier.description = deepcopy(getfield(identifier, :description)) + end + + metric = getfield(subset, :metric) + new_metric = getfield(new_subset, :metric) + if !IMASDD.ismissing(metric, :jacobian) + new_metric.jacobian = deepcopy(getfield(metric, :jacobian)) + end + if !IMASDD.ismissing(metric, :tensor_contravariant) + new_metric.tensor_contravariant = + deepcopy(getfield(metric, :tensor_contravariant)) + end + if !IMASDD.ismissing(metric, :tensor_covariant) + new_metric.tensor_covariant = deepcopy(getfield(metric, :tensor_covariant)) + end + return new_subset +end + +""" + Base.:∈( + point::Tuple{Real, Real}, + subset_of_space::Tuple{ + IMASDD.edge_profiles__grid_ggd___grid_subset, + IMASDD.edge_profiles__grid_ggd___space, + }, + ) Overloading ∈ operator to check if a point is inside a subset of space. @@ -360,27 +498,38 @@ it is checked if the point is within the enclosed area. It is assumed that a 2-dimensional subset used in such a context will form a closed area. If the subset is 3-dimensional, its boundary is calculated on the fly. If used multiple times, it is recommended to calculate the boundary once and store it in a variable. + +Example: + +```julia +if (5.5, 0.0) ∈ (subset_sol, space) + println("Point (5.5, 0.0) is inside the SOL subset.") +else + println("Point (5.5, 0.0) is outside the SOL subset.") +end +``` """ function Base.:∈( point::Tuple{Real, Real}, subset_of_space::Tuple{ - OMAS.edge_profiles__grid_ggd___grid_subset, - OMAS.edge_profiles__grid_ggd___space, + IMASDD.edge_profiles__grid_ggd___grid_subset, + IMASDD.edge_profiles__grid_ggd___space, }, ) r, z = point subset, space = subset_of_space - dim = subset.element[1].object[1].dimension - nodes = space.objects_per_dimension[1].object - edges = space.objects_per_dimension[2].object + dim = getfield(getfield(getfield(subset, :element)[1], :object)[1], :dimension) + opd = getfield(space, :objects_per_dimension) + nodes = getfield(opd[1], :object) + edges = getfield(opd[2], :object) if dim == 3 - subset_bnd = OMAS.edge_profiles__grid_ggd___grid_subset() + subset_bnd = IMASDD.edge_profiles__grid_ggd___grid_subset() subset_bnd.element = get_subset_boundary(space, subset) elseif dim == 2 subset_bnd = subset elseif dim == 1 - for ele ∈ subset.element - node = nodes[ele.object[1].index] + for ele ∈ getfield(subset, :element) + node = nodes[getfield(getfield(ele, :object)[1], :index)] if node.geometry[1] == r && node.geometry[2] == z return true end @@ -391,12 +540,18 @@ function Base.:∈( end # Count number of times an upward going ray from (r,z) intersects the boundary count = 0 - for ele ∈ subset_bnd.element - edge = edges[ele.object[1].index] - r_max = maximum(nodes[node].geometry[1] for node ∈ edge.nodes) - r_min = minimum(nodes[node].geometry[1] for node ∈ edge.nodes) + for ele ∈ getfield(subset_bnd, :element) + edge = edges[getfield(getfield(ele, :object)[1], :index)] + edge_nodes_r = zeros(2) + edge_nodes_z = zeros(2) + for (ii, node) ∈ enumerate(getfield(edge, :nodes)) + edge_nodes_r[ii] = getfield(nodes[node], :geometry)[1] + edge_nodes_z[ii] = getfield(nodes[node], :geometry)[2] + end + r_max = maximum(edge_nodes_r) + r_min = minimum(edge_nodes_r) if r_min <= r < r_max - z_max = maximum(nodes[node].geometry[2] for node ∈ edge.nodes) + z_max = maximum(edge_nodes_z) if z < z_max count += 1 end @@ -406,8 +561,17 @@ function Base.:∈( return count % 2 == 1 end +""" + get_prop_with_grid_subset_index( + prop::AbstractVector{T}, + grid_subset_index::Int, + ) where {T <: edge_profiles__prop_on_subset} + +Find the edge_profiles property instance in an array of properties that corresponds to +the grid_subset_index provided. +""" function get_prop_with_grid_subset_index( - prop::Vector{T}, + prop::AbstractVector{T}, grid_subset_index::Int, ) where {T <: edge_profiles__prop_on_subset} for p ∈ prop diff --git a/src/types.jl b/src/types.jl index b9850df..7f0e54d 100644 --- a/src/types.jl +++ b/src/types.jl @@ -1,13 +1,17 @@ +export get_types_with + """ get_types_with(parent::Type, field::Symbol) -A type creation utility meant for searching types in OMAS database. This function +A type creation utility meant for searching types in IMAS database. This function returns a list of types that are fields at any level below the parent data type which have a particular field present in it. Example: -get_types_with(OMAS.edge_profiles, :grid_subset_index) +```julia +get_types_with(IMASDD.edge_profiles, :grid_subset_index) +``` returns all edge_profiles types that have a subfield named grid_subset_index. """ @@ -37,4 +41,4 @@ function get_types_with(parent::Type, field::Symbol) end edge_profiles__prop_on_subset = - Union{get_types_with(OMAS.edge_profiles, :grid_subset_index)...} + Union{get_types_with(IMASDD.edge_profiles, :grid_subset_index)...} diff --git a/test/runtests.jl b/test/runtests.jl index c0c2b74..ce8c153 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,5 +1,5 @@ -import GGDUtils: interp, get_kdtree, project_prop_on_subset!, get_grid_subset_with_index -import OMAS: json2imas +import GGDUtils: interp, get_kdtree, project_prop_on_subset!, get_grid_subset +using IMASDD: IMASDD import Statistics: mean using Test using ArgParse: ArgParse @@ -34,7 +34,7 @@ end args = parse_commandline() print("json2imas() time: ") -@time ids = json2imas( +@time ids = IMASDD.json2imas( "$(@__DIR__)/../samples/time_dep_edge_profiles_last_step_only.json", ) @@ -66,7 +66,7 @@ if args["interp"] @test abs.((grid_val .- searched_val) ./ grid_val) < allowed_rtol # test interp(prop_arr, space, subset) - subset = get_grid_subset_with_index(grid_ggd, 5) + subset = get_grid_subset(grid_ggd, -5) print("interp(prop_arr, space, subset) time: ") @time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, space, subset) @@ -75,7 +75,7 @@ if args["interp"] # test interp(prop_arr, grid_ggd, grid_subset_index) print("interp(prop_arr, grid_ggd, grid_subset_index) time: ") - @time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, grid_ggd, 5) + @time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, grid_ggd, -5) searched_val = get_n_e(cell_center...) @test abs.((grid_val .- searched_val) ./ grid_val) < allowed_rtol @@ -128,9 +128,9 @@ if args["projection"] space = ids.edge_profiles.grid_ggd[1].space[1] prop = ids.edge_profiles.ggd[1].electrons.density # All cells - from_subset = get_grid_subset_with_index(ids.edge_profiles.grid_ggd[1], 5) + from_subset = get_grid_subset(ids.edge_profiles.grid_ggd[1], -5) # separatix - to_subset = get_grid_subset_with_index(ids.edge_profiles.grid_ggd[1], 16) + to_subset = get_grid_subset(ids.edge_profiles.grid_ggd[1], 16) print("project_prop_on_subset!(prop, from_subset, to_subset, space) time: ") @time separatix_centers, values_at_separatix = project_prop_on_subset!(prop, from_subset, to_subset, space) @@ -140,7 +140,7 @@ if args["projection"] # end subset_core = - get_grid_subset_with_index(ids.edge_profiles.grid_ggd[1], 22) + get_grid_subset(ids.edge_profiles.grid_ggd[1], 22) print("project_prop_on_subset!(prop, from_subset, subset_core) time: ") @time core_element_inds, values_at_core = project_prop_on_subset!(prop, from_subset, subset_core) @@ -161,9 +161,9 @@ if args["in"] @testset "test ∈" begin grid_ggd = ids.edge_profiles.grid_ggd[1] space = grid_ggd.space[1] - subset_corebnd = get_grid_subset_with_index(grid_ggd, 15) - subset_sol = get_grid_subset_with_index(grid_ggd, 23) - subset_odr = get_grid_subset_with_index(grid_ggd, 24) + subset_corebnd = get_grid_subset(grid_ggd, 15) + subset_sol = get_grid_subset(grid_ggd, 23) + subset_odr = get_grid_subset(grid_ggd, 24) @test (6.0, 0.0) ∈ (subset_corebnd, space) @test (5.0, -2.5) ∉ (subset_corebnd, space) diff --git a/test/speedtest.jl b/test/speedtest.jl new file mode 100644 index 0000000..51293a1 --- /dev/null +++ b/test/speedtest.jl @@ -0,0 +1,83 @@ +import GGDUtils: interp, get_kdtree, project_prop_on_subset!, get_grid_subset, IMASDD + +println("-----------------------------------------------------------------------------") +print("json2imas() time with compilation: ") +@time ids = IMASDD.json2imas( + "$(@__DIR__)/../samples/time_dep_edge_profiles_last_step_only.json", +) +print("json2imas() time (true runtime): ") +@time ids = IMASDD.json2imas( + "$(@__DIR__)/../samples/time_dep_edge_profiles_last_step_only.json", +) + +println("-----------------------------------------------------------------------------") +n_e = ids.edge_profiles.ggd[1].electrons.density[1]; +grid_ggd = ids.edge_profiles.grid_ggd[1]; +space = grid_ggd.space[1] +print("interp(prop, grid_ggd) time with compilation: ") +@time get_n_e = interp(n_e, grid_ggd) +print("interp(prop, grid_ggd) (true runtime): ") +@time get_n_e = interp(n_e, grid_ggd) + +println("-----------------------------------------------------------------------------") +subset = get_grid_subset(grid_ggd, 5) + +print("interp(prop_arr, space, subset) time with compilation: ") +@time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, space, subset) +print("interp(prop_arr, space, subset) time (true runtime): ") +@time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, space, subset) + +println("-----------------------------------------------------------------------------") +print("interp(prop_arr, grid_ggd, grid_subset_index) time with compilation: ") +@time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, grid_ggd, 5) +print("interp(prop_arr, grid_ggd, grid_subset_index) time (true runtime): ") +@time get_n_e = interp(ids.edge_profiles.ggd[1].electrons.density, grid_ggd, 5) + +println("-----------------------------------------------------------------------------") +kdtree = get_kdtree(space) +print("interp(prop_arr, kdtree) time with compilation: ") +@time get_T_e = interp(ids.edge_profiles.ggd[1].electrons.temperature[1].values, kdtree) +print("interp(prop_arr, kdtree) time (true runtime): ") +@time get_T_e = interp(ids.edge_profiles.ggd[1].electrons.temperature[1].values, kdtree) + +println("-----------------------------------------------------------------------------") +prop = ids.edge_profiles.ggd[1].electrons.density +# All cells +from_subset = get_grid_subset(ids.edge_profiles.grid_ggd[1], 5) +# separatix +to_subset = get_grid_subset(ids.edge_profiles.grid_ggd[1], 16) +print( + "project_prop_on_subset!(prop, from_subset, to_subset, space) time with compilation: ", +) +@time separatix_centers, values_at_separatix = + project_prop_on_subset!(prop, from_subset, to_subset, space) +print( + "project_prop_on_subset!(prop, from_subset, to_subset, space) time (true runtime): ", +) +@time separatix_centers, values_at_separatix = + project_prop_on_subset!(prop, from_subset, to_subset, space) + +println("-----------------------------------------------------------------------------") +subset_core = + get_grid_subset(ids.edge_profiles.grid_ggd[1], 22) +print("project_prop_on_subset!(prop, from_subset, subset_core) time with compilation: ") +@time core_element_inds, values_at_core = + project_prop_on_subset!(prop, from_subset, subset_core) +print("project_prop_on_subset!(prop, from_subset, subset_core) time (true runtime): ") +@time core_element_inds, values_at_core = + project_prop_on_subset!(prop, from_subset, subset_core) + +println("-----------------------------------------------------------------------------") +subset_corebnd = get_grid_subset(grid_ggd, 15) +subset_sol = get_grid_subset(grid_ggd, 23) +subset_odr = get_grid_subset(grid_ggd, 24) + +print("test ∈ (edges) time with compilation: ") +@time (6.0, 0.0) ∈ (subset_corebnd, space) +print("test ∈ (edges) time (true runtime): ") +@time (6.0, 0.0) ∈ (subset_corebnd, space) + +print("test ∈ (cells) time with compilation: ") +@time (6.0, 4.0) ∈ (subset_sol, space) +print("test ∈ (cells) time (true runtime): ") +@time (6.0, 4.0) ∈ (subset_sol, space)