diff --git a/.github/workflows/make_docs.yml b/.github/workflows/make_docs.yml new file mode 100644 index 0000000..a305e69 --- /dev/null +++ b/.github/workflows/make_docs.yml @@ -0,0 +1,40 @@ +name: Make Docs +on: + pull_request: + branches: ["master", "dev"] + push: + branches: + - master + - dev + - docs + paths: + - '.github/workflows/make_docs.yml' + - 'src/' + - 'docs/**' + tags: '*' + workflow_dispatch: + +jobs: + make_docs: + permissions: + actions: write + 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: Extract branch name + shell: bash + run: echo "branch=${GITHUB_BASE_REF:-${GITHUB_REF#refs/heads/}}" >> $GITHUB_OUTPUT + id: extract_branch + - name: Install dependencies + run: | + julia --project=docs/ -e 'using Pkg; Pkg.add(; url="https://github.com/ProjectTorreyPines/IMASDD.jl.git"); Pkg.add(; url="https://github.com/ProjectTorreyPines/GGDUtils.jl.git", rev="${{ steps.extract_branch.outputs.branch }}"); Pkg.add(; url="https://github.com/ProjectTorreyPines/SOLPS2IMAS.jl.git", rev="${{ steps.extract_branch.outputs.branch }}"); 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 diff --git a/.gitignore b/.gitignore index 47e8665..d5946fc 100644 --- a/.gitignore +++ b/.gitignore @@ -4,4 +4,5 @@ # environment. Manifest.toml sd_input_data.json -example/Project.toml +example/*.toml +docs/build diff --git a/docs/Project.toml b/docs/Project.toml new file mode 100644 index 0000000..9311ade --- /dev/null +++ b/docs/Project.toml @@ -0,0 +1,2 @@ +[deps] +Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" \ No newline at end of file diff --git a/docs/make.jl b/docs/make.jl new file mode 100644 index 0000000..4740ba7 --- /dev/null +++ b/docs/make.jl @@ -0,0 +1,17 @@ +using Documenter +using SD4SOLPS + +makedocs(; + modules=[SD4SOLPS], + format=Documenter.HTML(), + sitename="SD4SOLPS", + checkdocs=:none, +) + +deploydocs(; + repo="github.com/ProjectTorreyPines/SD4SOLPS.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..6340d95 --- /dev/null +++ b/docs/src/index.md @@ -0,0 +1,99 @@ + +# SD4SOLPS.jl + +```@contents +Pages = ["index.md"] +Depth = 3 +``` + +This repository serves as the top most workflow manager with helpful utilities to use other repositories in this project. + +## Documentation of other repositories in this project + +### [GGDUtils.jl](https://projecttorreypines.github.io/GGDUtils.jl/stable) + +### [SOLPS2IMAS.jl](https://projecttorreypines.github.io/SOLPS2IMAS.jl/stable) + +### [SynthDiag.jl](https://projecttorreypines.github.io/SynthDiag.jl/stable) + +## Installation + +### Using make: +After cloning this repo, check the make menu: +``` +SD4SOLPS.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.add(; url="https://github.com/ProjectTorreyPines/SOLPS2IMAS.jl.git"); +julia> Pkg.add(; url="https://github.com/JuliaFusion/EFIT.jl.git"); +julia> Pkg.add(; url="https://github.com/ProjectTorreyPines/SD4SOLPS.jl.git"); +julia> Pkg.instantiate() +``` + +## Top file handling functions + +```@docs +find_files_in_allowed_folders +geqdsk_to_imas! +preparation +``` + +## Repairing/filling out partial equilibrium files + +Tools for repairing/filling out partial equilibrium files. + +Some of the added fields may not be totally accurate, so it is recommended to +use this tool mainly for test cases, as a utility. For a real equilibrium, +problems should be fixed properly. + +```@docs +add_rho_to_equilibrium! +check_rho_1d +``` + +## Extrapolations + +Utilities for extrapolating profiles + +### Core profile extrapolations + +```@docs +extrapolate_core +fill_in_extrapolated_core_profile! +``` + +### Edge profiles extrapolations + +These functions have not been fully tested and/or supported yet. + +```@docs +mesh_psi_spacing +cached_mesh_extension! +``` + +## Unit conversion utilities + +```@docs +gas_unit_converter +``` \ No newline at end of file diff --git a/example/demo.ipynb b/example/demo.ipynb index 8be154c..45420d1 100644 --- a/example/demo.ipynb +++ b/example/demo.ipynb @@ -8,6 +8,15 @@ "source": [ "using Pkg;\n", "Pkg.activate(\"./\")\n", + "\n", + "# Can comment these after the first run\n", + "Pkg.add(; url=\"git@github.com:ProjectTorreyPines/IMASDD.jl.git\");\n", + "Pkg.add(; url=\"git@github.com:ProjectTorreyPines/GGDUtils.jl.git\");\n", + "Pkg.add(; url=\"git@github.com:ProjectTorreyPines/SOLPS2IMAS.jl.git\");\n", + "Pkg.add(; url=\"git@github.com:ProjectTorreyPines/SynthDiag.jl.git\");\n", + "Pkg.add(; url=\"git@github.com:JuliaFusion/EFIT.jl.git\");\n", + "Pkg.add(; url=\"git@github.com:ProjectTorreyPines/SD4SOLPS.jl.git\");\n", + "\n", "Pkg.instantiate()" ] }, diff --git a/example/makefile b/example/makefile deleted file mode 100644 index 72a84bc..0000000 --- a/example/makefile +++ /dev/null @@ -1,43 +0,0 @@ -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; \ - git clone "git@github.com:ProjectTorreyPines/GGDUtils.jl.git" ../../GGDUtils; \ - git clone "git@github.com:ProjectTorreyPines/SOLPS2IMAS.jl.git" ../../SOLPS2IMAS; \ - git clone "git@github.com:ProjectTorreyPines/SynthDiag.jl.git" ../../SynthDiag; \ - git clone "git@github.com:ProjectTorreyPines/SD4SOLPS.jl.git" ../../SD4SOLPS; - -julia --project=. -e 'using Pkg; Pkg.rm(["IMASDD"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["GGDUtils"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["SOLPS2IMAS"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["SynthDiag"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["SD4SOLPS"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["EFIT"]);' - -julia --project=. -e 'using Pkg; Pkg.develop(path="../../IMASDD"); Pkg.develop(path="../../GGDUtils"); Pkg.develop(path="../../SOLPS2IMAS"); Pkg.develop(path="../../SynthDiag"); Pkg.develop(path="../../SD4SOLPS"); Pkg.add(url="git@github.com:JuliaFusion/EFIT.jl.git", rev="master"); 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"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["GGDUtils"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["SOLPS2IMAS"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["SynthDiag"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["SD4SOLPS"]);' - -julia --project=. -e 'using Pkg; Pkg.rm(["EFIT"]);' - -julia --project=. -e 'using Pkg; Pkg.add(url="git@github.com:ProjectTorreyPines/IMASDD.jl.git", rev="master"); Pkg.add(url="git@github.com:ProjectTorreyPines/GGDUtils.jl.git", rev="master"); Pkg.add(url="git@github.com:ProjectTorreyPines/SOLPS2IMAS.jl.git", rev="master"); Pkg.add(url="git@github.com:ProjectTorreyPines/SynthDiag.jl.git", rev="master"); Pkg.add(url="git@github.com:ProjectTorreyPines/SD4SOLPS.jl.git", rev="master"); Pkg.add(url="git@github.com:JuliaFusion/EFIT.jl.git", rev="master"); Pkg.instantiate()' - -clean: - @echo "Deleting Manifest.toml" - - rm Manifest.toml diff --git a/src/SD4SOLPS.jl b/src/SD4SOLPS.jl index b9df33d..fea8b06 100644 --- a/src/SD4SOLPS.jl +++ b/src/SD4SOLPS.jl @@ -4,35 +4,37 @@ using IMASDD: IMASDD using SOLPS2IMAS: SOLPS2IMAS using EFIT: EFIT using Interpolations: Interpolations -#import GGDUtils -export find_files_in_allowed_folders -export geqdsk_to_imas! +export find_files_in_allowed_folders, geqdsk_to_imas!, preparation include("$(@__DIR__)/supersize_profile.jl") include("$(@__DIR__)/repair_eq.jl") -include("$(@__DIR__)/actuator_model.jl") - -greet() = print("Hello World!") +include("$(@__DIR__)/unit_utils.jl") """ - find_files_in_allowed_folders() + find_files_in_allowed_folders( + input_dirs::String...; + eqdsk_file::String, + recursive::Bool=true, + ) Searches a list of allowed folders for a set of filenames that will provide information about the SOLPS case. Returns a list of filenames with complete paths. Example: + +```julia SD4SOLPS.find_files_in_allowed_folders( -"/D3D_Ma_184833_03600", -eqdsk_file="g184833.03600", + "/D3D_Ma_184833_03600"; + eqdsk_file="g184833.03600", ) +``` """ function find_files_in_allowed_folders( - input_dirs...; - eqdsk_file, - recursive=true, - allow_reduced_versions=false, -) + input_dirs::String...; + eqdsk_file::String, + recursive::Bool=true, +)::Vector{String} files = ["b2fgmtry", "b2time.nc", "b2mn.dat", eqdsk_file] reduced_files = ["b2fgmtry_red", "b2time_red.nc", "b2mn.dat", eqdsk_file] @@ -46,7 +48,7 @@ function find_files_in_allowed_folders( else dirs = input_dirs end - for i ∈ 1:length(files) + for i ∈ eachindex(files) for dir ∈ dirs file = dir * "/" * files[i] reduced_file = dir * "/" * reduced_files[i] @@ -63,12 +65,22 @@ function find_files_in_allowed_folders( end """ - geqdsk_to_imas!() - -Transfers the equilibrium reconstruction in an EFIT-style gEQDSK file into + geqdsk_to_imas!( + eqdsk_file::String, + dd::IMASDD.dd; + set_time::Union{Nothing, Float64}=nothing, + time_index::Int=1, + ) + +Transfers the equilibrium reconstruction from an EFIT-style gEQDSK file into the IMAS DD structure. """ -function geqdsk_to_imas!(eqdsk_file, dd; set_time=nothing, time_index=1) +function geqdsk_to_imas!( + eqdsk_file::String, + dd::IMASDD.dd; + set_time::Union{Nothing, Float64}=nothing, + time_index::Int=1, +) # https://github.com/JuliaFusion/EFIT.jl/blob/master/src/io.jl g = EFIT.readg(eqdsk_file; set_time=set_time) # Copying ideas from OMFIT: omfit/omfit_classes/omfit_eqdsk.py / to_omas() @@ -134,7 +146,7 @@ function geqdsk_to_imas!(eqdsk_file, dd; set_time=nothing, time_index=1) if length(xrs) > 0 bx = eqt.boundary.x_point resize!(bx, length(xrs)) - for i ∈ length(xrs) + for i ∈ eachindex(xrs) bx[i].r = xrs[i] bx[i].z = xzs[i] end @@ -179,98 +191,28 @@ function geqdsk_to_imas!(eqdsk_file, dd; set_time=nothing, time_index=1) end """ - core_profile_2d(dd, prof_time_idx, eq_time_idx, quantity) - -Reads a 1D core profile and flux map and returns a quantity at requested R,Z points -dd: a data dictionary instance with equilibrium and core profile data loaded -quantity: a string specifying the quantity to fetch - -Returns a callable function that takes (r, z) as arguments and returns the quantity -""" -function core_profile_2d(dd, prof_time_idx, eq_time_idx, quantity) - if !check_rho_1d(dd; time_slice=eq_time_idx) - throw(ArgumentError("Equilibrium rho profile in input DD was missing.")) - end - prof = dd.core_profiles.profiles_1d[prof_time_idx] - rho_prof = prof.grid.rho_tor_norm - quantity_fields = split(quantity, ".") - p = prof - for field ∈ quantity_fields - p = getproperty(p, Symbol(field)) - end - eqt = dd.equilibrium.time_slice[eq_time_idx] - p1 = eqt.profiles_1d - p2 = eqt.profiles_2d[1] - gq = eqt.global_quantities - psi_a = gq.psi_axis - psi_b = gq.psi_boundary - rhon_eq = p1.rho_tor_norm - psi_eq = p1.psi - psin_eq = (psi_eq .- psi_a) ./ (psi_b - psi_a) - psirz = p2.psi - psinrz = (psirz .- psi_a) ./ (psi_b - psi_a) - r_eq = p2.grid.dim1 - z_eq = p2.grid.dim2 - extension = [1.0001, 1.1, 5] - # rho_N isn't defined on open flux surfaces, so it is extended by copying psi_N - psin_eq_ext = copy(psin_eq) - append!(psin_eq_ext, extension) - rhon_eq_ext = copy(rhon_eq) - append!(rhon_eq_ext, extension) - neg_extension = [-5, -0.0001] # I guess this would be at a PF coil or something? - prepend!(psin_eq_ext, neg_extension) - prepend!(rhon_eq_ext, neg_extension) - - rho_prof_ext = vcat(rho_prof, extension) - p_ext = vcat(p, zeros(size(extension))) - rho_prof_ext = prepend!(rho_prof_ext, neg_extension) - p_ext = prepend!(p_ext, zeros(size(neg_extension))) - - # Data are set up and ready to process - - # EDGE PROFILES (the input data): - # rho_prof_ext: rho_N values associated with the profile of some quantity - # p_ext : profile of some quantity vs. rho_prof - - # EQUILBRIUM (the connection from rho to R,Z via psi): - # psin_eq_ext : 1D array of psi_N values that corresponds to rhon_eq_ext - # rhon_eq_ext : 1D array of rho_N values that corresponds to psin_eq_ext - # --> connects rho to everything else via psi - # psinrz : 2D array of psi_N values that corresponds to r_eq and z_eq - # r_eq and z_eq : 1D coordinate arrays that correspond to psinrz - - # OUTPUT INSTRUCTIONS: - # r and z : coordinates of output points where values of p are desired - - # psi_at_requested_points = - # Interpolations.linear_interpolation((r_eq, z_eq), psinrz).(r, z) - # rhonpsi = Interpolations.linear_interpolation(psin_eq_ext, rhon_eq_ext) - # rho_at_requested_points = rhonpsi.(psi_at_requested_points) - # itp = Interpolations.linear_interpolation(rho_prof_ext, p_ext) - # p_at_requested_points = itp.(rho_at_requested_points) - # return p_at_requested_points - rz2psin = Interpolations.linear_interpolation((r_eq, z_eq), psinrz) - psin2rhon = Interpolations.linear_interpolation(psin_eq_ext, rhon_eq_ext) - rhon2prof = Interpolations.linear_interpolation(rho_prof_ext, p_ext) - return (r, z) -> rhon2prof(psin2rhon(rz2psin(r, z))) -end - -""" - preparation() + preparation( + eqdsk_file::String, + dirs::String...; + core_method::String="simple", + filename::String="sd_input_data", + output_format::String="json", + eqdsk_set_time::Union{Nothing, Float64}=nothing, + eq_time_index::Int64=1, + )::IMASDD.dd Gathers SOLPS and EFIT files and loads them into IMAS structure. Extrapolates profiles as needed to get a complete picture. """ function preparation( - eqdsk_file, - dirs...; + eqdsk_file::String, + dirs::Vector{String}; core_method::String="simple", - edge_method::String="simple", filename::String="sd_input_data", output_format::String="json", - eqdsk_set_time=nothing, - eq_time_index=1, -) + eqdsk_set_time::Union{Nothing, Float64}=nothing, + eq_time_index::Int64=1, +)::IMASDD.dd b2fgmtry, b2time, b2mn, eqdsk = find_files_in_allowed_folders(dirs...; eqdsk_file=eqdsk_file) println("Found source files:") diff --git a/src/repair_eq.jl b/src/repair_eq.jl index 48e2100..3d9bae7 100644 --- a/src/repair_eq.jl +++ b/src/repair_eq.jl @@ -13,11 +13,19 @@ export add_rho_to_equilibrium! export check_rho_1d """ - function check_rho_1d() + check_rho_1d( + dd::IMASDD.dd; + time_slice::Int64=1, + throw_on_fail::Bool=false, + )::Bool Checks to see if rho exists and is valid in the equilibrium 1d profiles """ -function check_rho_1d(dd::IMASDD.dd; time_slice::Int64=1, throw_on_fail::Bool=false) +function check_rho_1d( + dd::IMASDD.dd; + time_slice::Int64=1, + throw_on_fail::Bool=false, +)::Bool rho = dd.equilibrium.time_slice[time_slice].profiles_1d.rho_tor_norm if length(rho) < 1 rho_okay = false diff --git a/src/supersize_profile.jl b/src/supersize_profile.jl index cb9f070..9a60a00 100644 --- a/src/supersize_profile.jl +++ b/src/supersize_profile.jl @@ -2,7 +2,6 @@ Utilities for extrapolating profiles """ -# import CalculusWithJulia using IMASDD: IMASDD using Interpolations: Interpolations using GGDUtils: @@ -11,13 +10,14 @@ using GGDUtils: using PolygonOps: PolygonOps using JSON: JSON -export extrapolate_core -export fill_in_extrapolated_core_profile -export mesh_psi_spacing -export find_x_points! +export extrapolate_core, + fill_in_extrapolated_core_profile!, mesh_psi_spacing, cached_mesh_extension! """ - cumul_integrate(x::AbstractVector, y::AbstractVector) + cumul_integrate( + x::AbstractVector, + y::AbstractVector, + )::AbstractVector Computes cumulative integral of y(x) using trapezoidal rule. Source code from obsolete NumericalIntegration.jl package. @@ -25,7 +25,10 @@ https://github.com/dextorious/NumericalIntegration.jl/blob/540b6c4bee089dfef7b9a Modified the code to remove use of @inbounds, @fastmath macros that Julia documentation recommends to use with caution. """ -function cumul_integrate(x::AbstractVector, y::AbstractVector) +function cumul_integrate( + x::AbstractVector, + y::AbstractVector, +)::AbstractVector init = (x[2] - x[1]) * (y[1] + y[2]) n = length(x) retarr = Vector{typeof(init)}(undef, n) @@ -38,7 +41,11 @@ function cumul_integrate(x::AbstractVector, y::AbstractVector) end """ - extrapolate_core(edge_rho, edge_quantity, rho_output) + extrapolate_core( + edge_rho::Vector{Float64}, + edge_quantity::Vector{Float64}, + rho_output::Vector{Float64}, + )::Vector{Float64} Function for assuming a core profile when given edge profile data. @@ -54,7 +61,11 @@ Concept: 5. after making up a very simple gradient profile out of a few line segments, integrate it to get the profile of the quantity in question """ -function extrapolate_core(edge_rho, edge_quantity, rho_output) +function extrapolate_core( + edge_rho::Vector{Float64}, + edge_quantity::Vector{Float64}, + rho_output::Vector{Float64}, +)::Vector{Float64} grad = IMASDD.gradient(edge_rho, edge_quantity) gf = grad[1] rf = edge_rho[1] @@ -89,46 +100,48 @@ function extrapolate_core(edge_rho, edge_quantity, rho_output) ).(rho_output[rho_output.>=rf]) return output_profile end -#!format off + """ fill_in_extrapolated_core_profile!( - dd::IMASDD.dd, - quantity_name::String; - method::String="simple", - eq_time_idx::Int64=1, - eq_profiles_2d_idx::Int64=1, - grid_ggd_idx::Int64=1, - space_idx::Int64=1, -) - + dd::IMASDD.dd, + quantity_name::String; + method::String="simple", + eq_time_idx::Int64=1, + eq_profiles_2d_idx::Int64=1, + grid_ggd_idx::Int64=1, + space_idx::Int64=1, + ) -This function accepts a DD that should be populated with equilibrium and edge_profiles -as well as a request for a quantity to extrapolate into the core. It then maps -edge_profiles data to rho, calls the function that performs the extrapolation (which is +This function accepts a DD that should be populated with `equilibrium` and +`edge_profiles` as well as a request for a quantity to extrapolate into the core. It +then maps `edge_profiles` data to rho, calls the function that performs the extrapolation (which is not a simple linear extrapolation but has some trickery to attempt to make a somewhat convincing profile shape), and writes the result to core_profiles. This involves a bunch of interpolations and stuff. -dd: an IMAS data dictionary -quantity_name: the name of a quantity in edge_profiles.profiles_2d and - core_profiles.profiles_1d, such as "electrons.density" -method: Extrapolation method. -eq_time_idx: index of the equilibrium time slice to use. For a typical SOLPS run, - the SOLPS mesh will be based on the equilibrium reconstruction at a single - time, so the DD associated with the SOLPS run only needs one equilibrium - time slice to be loaded. However, one could combine the complete - equilibrium time series with the SOLPS run and then have to specify which - slice of the equilibrium corresponds to the SOLPS mesh. -grid_ggd_idx: index of the grid_ggd to use. For a typical SOLPS run, the SOLPS grid is - fixed, so this index defaults to 1. But in future, if a time varying grid - is used, then this index will need to be specified. -space_idx: index of the space to use. For a typical SOLPS run, there will be only one - space so this index will mostly remain at 1. -cell_subset_idx: index of the subset of cells to use for the extrapolation. The default - is 5, which is the subset of all cells. If edge_profiles data is - instead present for a different subset, for instance, -5, which are - b2.5 cells only, then this index should be set to -5. + +Input arguments: + + - `dd`: an IMAS data dictionary + - `quantity_name`: the name of a quantity in `edge_profiles.profiles_2d` and + `core_profiles.profiles_1d`, such as "electrons.density" + - `method`: Extrapolation method. + - `eq_time_id`x: index of the equilibrium time slice to use. For a typical SOLPS run, + the SOLPS mesh will be based on the equilibrium reconstruction at a single + time, so the DD associated with the SOLPS run only needs one equilibrium + time slice to be loaded. However, one could combine the complete + equilibrium time series with the SOLPS run and then have to specify which + slice of the equilibrium corresponds to the SOLPS mesh. + - `eq_profiles_2d_idx`: index of the `profiles_2D` in equilibrium `time_slice`. + - `grid_ggd_idx`: index of the `grid_ggd` to use. For a typical SOLPS run, the SOLPS + grid is fixed, so this index defaults to 1. But in future, if a time varying grid is + used, then this index will need to be specified. + - `space_id`x: index of the space to use. For a typical SOLPS run, there will be only + one space so this index will mostly remain at 1. + - `cell_subset_idx`: index of the subset of cells to use for the extrapolation. The + default is 5, which is the subset of all cells. If `edge_profiles` data is instead + present for a different subset, for instance, -5, which are b2.5 cells only, then + this index should be set to -5. """ -#!format on function fill_in_extrapolated_core_profile!( dd::IMASDD.dd, quantity_name::String; @@ -142,10 +155,8 @@ function fill_in_extrapolated_core_profile!( check_rho_1d(dd; time_slice=eq_time_idx, throw_on_fail=true) grid_ggd = dd.edge_profiles.grid_ggd[grid_ggd_idx] space = grid_ggd.space[space_idx] - cell_subset = - get_grid_subset(grid_ggd, cell_subset_idx) - midplane_subset = - get_grid_subset(grid_ggd, 11) + cell_subset = get_grid_subset(grid_ggd, cell_subset_idx) + midplane_subset = get_grid_subset(grid_ggd, 11) if length(midplane_subset.element) < 1 throw( @@ -258,38 +269,43 @@ function fill_in_extrapolated_core_profile!( end """ - function extrapolate_edge_exp( + extrapolate_edge_exp( quantity_edge::Vector{Float64}, dqdpsi::Vector{Float64}, psin_out::Vector{Float64}, - ) + )::Matrix{Float64} Exponential decay version of edge profile extrapolation. Should work well for many quantities, including Te and ne. If the exponential profile has no background offset, then its amplitude and scale length can be completely defined by matching the quantity value and first derivative at the edge of the mesh. -quantity_edge: values of some physics quantity in cells along the outer edge of the mesh -dqdpsi: Gradient of the quantity vs. psi, aligned perpendicular to the row of cells -being used. -psin_out: Normalized psi values along a vector orthogonal to the row of cells along the -edge. These psi_N values should be outside of the mesh (because the quantity is -already known in the mesh). +Input Arguments: + + - `quantity_edge`: values of some physics quantity in cells along the outer edge of + the mesh. + + - `dqdpsi`: Gradient of the quantity vs. psi, aligned perpendicular to the row of + cells being used. + - `psin_out`: Normalized psi values along a vector orthogonal to the row of cells + along the edge. These `psi_N` values should be outside of the mesh (because the + quantity is already known in the mesh). + The output will be a matrix. """ function extrapolate_edge_exp( quantity_edge::Vector{Float64}, dqdpsi::Vector{Float64}, psin_out::Vector{Float64}, -) +)::Matrix{Float64} x = psin_out - 1.0 lambda = -quantity_edge / dqdpsi q0 = quantity_edge ./ exp(-x' ./ lambda) - return y0 * exp(-x ./ lambda) + return q0 * exp(-x ./ lambda) end """ - prep_flux_map() + prep_flux_map(dd::IMASDD.dd; eq_time_idx::Int64=1, eq_profiles_2d_idx::Int64=1) Reads equilibrium data and extracts/derives some useful quantities. This is very basic, but it was being repeated and that's a no-no. @@ -313,38 +329,40 @@ function prep_flux_map(dd::IMASDD.dd; eq_time_idx::Int64=1, eq_profiles_2d_idx:: return r_eq, z_eq, psin_eq, rzpi end -#! format off """ mesh_psi_spacing( - dd::IMASDD.dd; - eq_time_idx::Int64=1, - eq_profiles_2d_idx::Int64=1, - grid_ggd_idx::Int64=1, - space_idx::Int64=1, - -) + dd::IMASDD.dd; + eq_time_idx::Int64=1, + eq_profiles_2d_idx::Int64=1, + grid_ggd_idx::Int64=1, + space_idx::Int64=1, + avoid_guard_cell::Bool=true, + spacing_rule="mean", + ) Inspects the mesh to see how far apart faces are in psi_N. Requires that GGD and equilibrium are populated. -dd: a data dictionary instance with required data loaded into it -eq_time_idx: index of the equilibrium time slice to use. For a typical SOLPS run, - the SOLPS mesh will be based on the equilibrium reconstruction at a single - time, so the DD associated with the SOLPS run only needs one equilibrium - time slice to be loaded. However, one could combine the complete - equilibrium time series with the SOLPS run and then have to specify which - slice of the equilibrium corresponds to the SOLPS mesh. -grid_ggd_idx: index of the grid_ggd to use. For a typical SOLPS run, the SOLPS grid is - fixed, so this index defaults to 1. But in future, if a time varying grid - is used, then this index will need to be specified. -space_idx: index of the space to use. For a typical SOLPS run, there will be only one - space so this index will mostly remain at 1. -avoid_guard_cell: assume that the last cell is a guard cell so take end-2 and end-1 - instead of end and end-1 -spacing_rule: "edge" or "mean" to make spacing of new cells (in psi_N) be the same - as the spacing at the edge of the mesh, or the same as the average spacing +Input Arguments: + + - `dd`: a data dictionary instance with required data loaded into it + - `eq_time_idx`: index of the equilibrium time slice to use. For a typical SOLPS run, + the SOLPS mesh will be based on the equilibrium reconstruction at a single + time, so the DD associated with the SOLPS run only needs one equilibrium + time slice to be loaded. However, one could combine the complete + equilibrium time series with the SOLPS run and then have to specify which + slice of the equilibrium corresponds to the SOLPS mesh. + - `eq_profiles_2d_id`x: index of the `profiles_2D` in equilibrium `time_slice`. + - `grid_ggd_idx`: index of the `grid_ggd` to use. For a typical SOLPS run, the SOLPS + grid is fixed, so this index defaults to 1. But in future, if a time varying grid + is used, then this index will need to be specified. + - `space_idx`: index of the space to use. For a typical SOLPS run, there will be only + one space so this index will mostly remain at 1. + - `avoid_guard_cell`: assume that the last cell is a guard cell so take `end-2` and + `end-1` instead of `end` and `end-1` + - `spacing_rule`: "edge" or "mean" to make spacing of new cells (in `psi_N`) be the + same as the spacing at the edge of the mesh, or the same as the average spacing """ -#! format on function mesh_psi_spacing( dd::IMASDD.dd; eq_time_idx::Int64=1, @@ -382,8 +400,7 @@ function mesh_psi_spacing( # weird. So use the outboard midplane. That's always a solid choice. grid_ggd = dd.edge_profiles.grid_ggd[grid_ggd_idx] space = grid_ggd.space[space_idx] - midplane_subset = - get_grid_subset(grid_ggd, 11) + midplane_subset = get_grid_subset(grid_ggd, 11) midplane_cell_centers = get_subset_centers(space, midplane_subset) r_mesh = [midplane_cell_centers[i][1] for i ∈ eachindex(midplane_cell_centers)] z_mesh = [midplane_cell_centers[i][2] for i ∈ eachindex(midplane_cell_centers)] @@ -409,12 +426,18 @@ function mesh_psi_spacing( end """ - pick_extension_psi_range() - -Defines the psi_N levels for an extended mesh. The range of psi_N levels starts -at the outer edge of the existing edge_profiles mesh at the midplane and goes + pick_extension_psi_range( + dd::IMASDD.dd; + eq_time_idx::Int64=1, + eq_profiles_2d_idx::Int64=1, + grid_ggd_idx::Int64=1, + space_idx::Int64=1, + )::Vector{Float64} + +Defines the `psi_N` levels for an extended mesh. The range of `psi_N` levels starts +at the outer edge of the existing `edge_profiles` mesh at the midplane and goes out to the most distant (in flux space) point on the limiting surface. -Returns a vector of psi_N levels. +Returns a vector of `psi_N` levels. """ function pick_extension_psi_range( dd::IMASDD.dd; @@ -422,7 +445,7 @@ function pick_extension_psi_range( eq_profiles_2d_idx::Int64=1, grid_ggd_idx::Int64=1, space_idx::Int64=1, -) +)::Vector{Float64} r_eq, z_eq, psin_eq, rzpi = prep_flux_map(dd; eq_time_idx, eq_profiles_2d_idx) # Use wall to find maximum extent of contouring project @@ -465,17 +488,20 @@ function pick_extension_psi_range( end """ - pick_mesh_ext_starting_points(dd; grid_ggd_idx, space_idx) + pick_mesh_ext_starting_points( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + space::IMASDD.edge_profiles__grid_ggd___space, + )::Tuple{Vector{Float64}, Vector{Float64}} Picks starting points for the radial lines of the mesh extension. The strategy is to start from the outer edge of the existing mesh and follow the steepest gradient (of psi_N) to extend these gridlines outward. -dd: a data dictionary instance with edge_profiles ggd and equilibrium loaded -grid_ggd_idx: index within ggd -space_idx: space number / index of the space to work with within edge_profiles Returns a tuple with vectors of R and Z starting points. """ -function pick_mesh_ext_starting_points(grid_ggd, space) +function pick_mesh_ext_starting_points( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + space::IMASDD.edge_profiles__grid_ggd___space, +)::Tuple{Vector{Float64}, Vector{Float64}} # Choose starting points for the orthogonal (to the contour) gridlines # Use the existing cells of the standard mesh all_cell_subset = get_grid_subset(grid_ggd, 5) @@ -508,28 +534,39 @@ function pick_mesh_ext_starting_points(grid_ggd, space) return r, z end -#!format off """ - mesh_ext_follow_grad() + mesh_ext_follow_grad( + r_eq::Vector{Float64}, + z_eq::Vector{Float64}, + psin_eq::Matrix, + rstart::Vector{Float64}, + zstart::Vector{Float64}, + nlvl::Int64, + dpsin::Float64, + rzpi=nothing, + )::Tuple{Matrix{Float64}, Matrix{Float64}} Follows the steepest gradient from a set of starting points, dropping nodes at approximately regular intervals in psi_N. Due to the numerical techniques used, the node spacing may be imperfect (especially prone to error in regions where curvature -of psi_N is large compared to its gradient). -r_eq: Equilibrium reconstruction's grid, R coordinates -z_eq: Equilibrium reconstruction's grid, Z coordinates -psin_eq: Normalized poloidal flux in the equilibrium reconstruction as a function of R and Z -rstart: R coordinates of starting points for the gradient following. -zstart: Z coordinates of starting points -nlvl: number of nodes to drop while following the gradient -dpsin: node spacing in delta psi_N -rzpi: linear interpolation of psin_eq() as a function of r_eq and z_eq - This was probably already computed and I think time would be saved by reusing it. - If you don't already have it, you can pass in nothing and let this function calculate it. +of `psi_N` is large compared to its gradient). + +Input Arguments: + + - `r_eq`: Equilibrium reconstruction's grid, R coordinates + - `z_eq`: Equilibrium reconstruction's grid, Z coordinates + - `psin_eq`: Normalized poloidal flux in the equilibrium reconstruction as a function + of R and Z + - `rstar`t: R coordinates of starting points for the gradient following. + - `zstart`: Z coordinates of starting points + - `nlvl`: number of nodes to drop while following the gradient + - `dpsin`: node spacing in delta psi_N + - `rzpi`: linear interpolation of psin_eq() as a function of r_eq and z_eq. This was + probably already computed and I think time would be saved by reusing it. If you + don't already have it, you can pass in nothing and let this function calculate it. Returns two matrices with R and Z coordinates of the mesh extension """ -#!format on function mesh_ext_follow_grad( r_eq::Vector{Float64}, z_eq::Vector{Float64}, @@ -539,7 +576,7 @@ function mesh_ext_follow_grad( nlvl::Int64, dpsin::Float64, rzpi=nothing, -) +)::Tuple{Matrix{Float64}, Matrix{Float64}} npol = length(rstart) mesh_r = zeros((npol, nlvl)) mesh_z = zeros((npol, nlvl)) @@ -587,13 +624,20 @@ function mesh_ext_follow_grad( end """ - modify_mesh_ext_near_x! + modify_mesh_ext_near_x!( + eqt::IMASDD.equilibrium__time_slice, + mesh_r::Matrix{Float64}, + mesh_z::Matrix{Float64}, + ) Modifies an extended mesh near a secondary X-point to compensate for the tendency of the mesh to go nuts near the X-point. -eqt: equilibrium.time_slice information -mesh_r: matrix of R values for the extended mesh -mesh_z: matrix of Z values for the extended mesh + +Input Arguments: + + - `eqt`: equilibrium.time_slice information + - `mesh_r`: matrix of R values for the extended mesh + - `mesh_z`: matrix of Z values for the extended mesh """ function modify_mesh_ext_near_x!( eqt::IMASDD.equilibrium__time_slice, @@ -677,24 +721,29 @@ function modify_mesh_ext_near_x!( end end -#!format off """ - record_regular_mesh!() + record_regular_mesh!( + grid_ggd::IMASDD.edge_profiles__grid_ggd, + space::IMASDD.edge_profiles__grid_ggd___space, + mesh_r::Matrix{Float64}, + mesh_z::Matrix{Float64}, + cut::Int64, + ) Records arrays of mesh data from regular 2D arrays into the DD -grid_ggd: grid_ggd within edge_profiles -space: space in edge_profiles -mesh_r: Matrix of R values along a mesh. Should be 2D. The two dimensions are in - the radial and poloidal directions. -mesh_z: Z values to go with mesh_r. -cut: Poloidal index of a cut between two groups of disconnected cells. Poloidal - connections (faces, cells) will not be added between this poloidal index - and the next index. + + - `grid_ggd`: `grid_ggd` within edge_profiles + - `space`: space in edge_profiles + - `mesh_r`: Matrix of R values along a mesh. Should be 2D. The two dimensions are in + the radial and poloidal directions. + - `mesh_z`: Z values to go with mesh_r. + - `cut`: Poloidal index of a cut between two groups of disconnected cells. Poloidal + connections (faces, cells) will not be added between this poloidal index + and the next index. """ -#!format on function record_regular_mesh!( - grid_ggd, - space, + grid_ggd::IMASDD.edge_profiles__grid_ggd, + space::IMASDD.edge_profiles__grid_ggd___space, mesh_r::Matrix{Float64}, mesh_z::Matrix{Float64}, cut::Int64, @@ -826,7 +875,7 @@ function record_regular_mesh!( end """ - convert_filename(filename::String) + convert_filename(filename::String)::String Converts a filename into a string that doesn't have illegal characters. The main application is removing the path separator from source files with full @@ -835,7 +884,7 @@ files used to form some data can be part of the cache name, allowing quick looku the cache filename is defined by the input files, and if it doesn't exist, it needs to be generated. """ -function convert_filename(filename::String) +function convert_filename(filename::String)::String filename_mod = replace(filename, "/" => "__") # Illegal on *nix bc it's the path separator filename_mod = replace(filename_mod, ":" => "--") # Illegal on mac and windows filename_mod = replace(filename_mod, "\\" => "__") # Illegal on windows @@ -848,24 +897,34 @@ function convert_filename(filename::String) return filename_mod end -#!format off """ - cached_mesh_extension!() + cached_mesh_extension!( + dd::IMASDD.dd, + eqdsk_file::String, + b2fgmtry::String; + eq_time_idx::Int64=1, + eq_profiles_2d_idx::Int64=1, + grid_ggd_idx::Int64=1, + space_idx::Int64=1, + clear_cache::Bool=false, + )::String Adds an extended mesh to a data dictionary, possibly from a cached result. -dd: The data dictionary. It will be modified in place. -eqdsk_file: the name of the EQDSK file that was used to get equilibrium data in - the dd. -b2fgmtry: the name of the SOLPS geometry file that was used to get GGD info in - edge_profiles in the dd. -eq_time_idx: Index of the time slice in equilibrium -eq_profiles_2d_idx: Index of the 2D profile set in equilibrium - (there is usually only one) -grid_ggd_idx: Index of the grid_ggd set in edge_profiles -space_idx: Index of the space -clear_cache: delete any existing cache file (for use in testing) + +Input Arguments: + + - `dd`: The data dictionary. It will be modified in place. + - `eqdsk_file`: the name of the EQDSK file that was used to get equilibrium data in + the dd. + - `b2fgmtry`: the name of the SOLPS geometry file that was used to get GGD info in + `edge_profiles` in the dd. + - `eq_time_idx`: Index of the time slice in equilibrium + - `eq_profiles_2d_idx`: Index of the 2D profile set in equilibrium + (there is usually only one) + - `grid_ggd_idx`: Index of the `grid_ggd` set in edge_profiles + - `space_idx`: Index of the space + - `clear_cache`: delete any existing cache file (for use in testing) """ -#!format on function cached_mesh_extension!( dd::IMASDD.dd, eqdsk_file::String, @@ -874,8 +933,8 @@ function cached_mesh_extension!( eq_profiles_2d_idx::Int64=1, grid_ggd_idx::Int64=1, space_idx::Int64=1, - clear_cache=false, -) + clear_cache::Bool=false, +)::String path = "$(@__DIR__)/../data/" cached_ext_name = path * string(hash(eqdsk_file * b2fgmtry)) * ".mesh_ext.json" if clear_cache @@ -931,7 +990,13 @@ function cached_mesh_extension!( end """ - function mesh_extension_sol() + mesh_extension_sol!( + dd::IMASDD.dd; + eq_time_idx::Int64=1, + eq_profiles_2d_idx::Int64=1, + grid_ggd_idx::Int64=1, + space_idx::Int64=1, + ) Extends the mesh out into the SOL """ @@ -983,7 +1048,10 @@ end """ fill_in_extrapolated_edge_profile!( - dd::IMASDD.dd, quantity_name::String; method::String="simple", + dd::IMASDD.dd, + quantity_name::String; + method::String="simple", + eq_time_idx::Int64=1, ) JUST A PLACEHOLDER FOR NOW. DOESN'T ACTUALLY WORK YET. diff --git a/src/actuator_model.jl b/src/unit_utils.jl similarity index 59% rename from src/actuator_model.jl rename to src/unit_utils.jl index ec0a0f0..739a55d 100644 --- a/src/actuator_model.jl +++ b/src/unit_utils.jl @@ -1,9 +1,7 @@ -# Actuator models to translate commands (probably in V) into gas flows - -using PhysicalConstants.CODATA2018 +using PhysicalConstants.CODATA2018: BoltzmannConstant, ElementaryCharge using Unitful: Unitful -using Interpolations: Interpolations -using YAML: YAML + +export gas_unit_converter torrl_per_pam3 = Float64( Unitful.upreferred((1 * Unitful.Torr * Unitful.L) / (Unitful.Pa * Unitful.m^3)), @@ -28,7 +26,13 @@ electrons_per_molecule = Dict( ) """ -gas_unit_converter() + gas_unit_converter( + value_in::Float64, + units_in::String, + units_out::String; + species::String="H", + temperature::Float64=293.15, + ) Converts gas flows between different units. Uses ideal gas law to convert between Pressure * volume type flows / quantities and count / current types of units. @@ -36,8 +40,11 @@ There is a version that accepts floats in and outputs floats, and another that deals in Unitful quantities. """ function gas_unit_converter( - value_in::Float64, units_in::String, units_out::String; species::String="H", - temperature=293.15, + value_in::Float64, + units_in::String, + units_out::String; + species::String="H", + temperature::Float64=293.15, ) if units_in == units_out return value_in @@ -76,26 +83,38 @@ function gas_unit_converter( end """ -gas_unit_converter() + gas_unit_converter( + value_in::Unitful.Quantity, + units_in::String, + units_out::String; + species::String="H", + temperature=293.15 * Unitful.K, + ) Converts gas flows between different units. Uses ideal gas law to convert between Pressure * volume type flows / quantities and count / current types of units. This is the Unitful version. + Output will be unitful, but the units are not simplified automatically. You can perform operations such as -(output |> Unitful.upreferred).val -Unitful.uconvert(Unitful.whatever, output).val + + - `(output |> Unitful.upreferred).val` + - `Unitful.uconvert(Unitful.whatever, output).val` + to handle simplification or conversion of units. -Although this function pretends torr L s^-1 and Pa m^3 s^-1 are different, use of -Unitful should cause them to behave the same way as long as you simplify or convert -units at the end. This means that you can use other pressure*volume type gas units -and call them torr L s^-1 and the script will deal with them up to having messy -units in the output. +Although this function pretends torr L s``^{-1}`` and Pa m``^3`` s``^{-1}`` are +different, use of Unitful should cause them to behave the same way as long as you +simplify or convert units at the end. This means that you can use other pressure*volume +type gas units and call them torr L s``^{-1}`` and the script will deal with them up to +having messy units in the output. """ function gas_unit_converter( - value_in::Unitful.Quantity, units_in::String, units_out::String; - species::String="H", temperature=293.15 * Unitful.K, + value_in::Unitful.Quantity, + units_in::String, + units_out::String; + species::String="H", + temperature=293.15 * Unitful.K, ) if units_in == units_out return value_in @@ -144,79 +163,3 @@ function gas_unit_converter( end return value_in .* conversion_factor end - -function select_default_config(model::String) - froot = model - if model == "instant" - froot = "simple" - end - filename = froot * "_gas_valve.yml" - return filename -end - -""" - model_gas_valve() - -The main function for handling a gas valve model. Has logic for selecting models -and configurations. -""" -function model_gas_valve( - model::String; configuration_file::String="auto", species::String="D2", -) - # Select configuration file - if configuration_file == "auto" - configuration_file = select_default_config(model) - end - default_config_path = "$(@__DIR__)/../config/" - if !occursin("/", configuration_file) - configuration_file = default_config_path * configuration_file - end - if !isfile(configuration_file) - throw(ArgumentError("Configuration file not found: " * configuration_file)) - end - config = YAML.load_file(configuration_file) - if (species in ["H", "H2", "D", "T", "T2"]) & !(species in keys(config)) - config = config["D2"] # They should all be the same - elseif (species in ["CH"]) & !(species in keys(config)) - config = config["CD"] - elseif (species in ["CH4", "CD4"]) & !(species in keys(config)) - config = config["CD4"] - else - config = config[species] - end - - if model == "simple" - function simple_gas_model(t, command) - flow0 = instant_gas_model(command, config) - t_ext = copy(t) - prepend!(t_ext, minimum(t) - config["delay"]) - flow0_ext = copy(flow0) - prepend!(flow0_ext, flow0[1]) - interp = Interpolations.linear_interpolation(t_ext, flow0_ext) - delayed_flow = interp.(t .- config["delay"]) - return lowpass_filter(t, delayed_flow, config["tau"]) - end - return simple_gas_model - elseif model == "instant" - instant_gas_model_(t, command) = instant_gas_model(command, config) - return instant_gas_model_ - else - throw(ArgumentError("Unrecognized model: " * model)) - end -end - -function instant_gas_model(command, config) - return config["p1"] .* (sqrt.(((command * config["p2"]) .^ 2.0 .+ 1) .- 1)) -end - -function lowpass_filter_(raw, previous_smooth, dt, tau) - return previous_smooth + dt / (dt + tau) * (raw - previous_smooth) -end - -function lowpass_filter(t, x, tau) - xs = zeros(length(t)) - for i ∈ 2:length(t) - xs[i] = lowpass_filter_(x[i], xs[i-1], t[i] - t[i-1], tau) - end - return xs -end diff --git a/test/Project.toml b/test/Project.toml deleted file mode 100644 index 495d8a9..0000000 --- a/test/Project.toml +++ /dev/null @@ -1,3 +0,0 @@ -[deps] -ArgParse = "c7e460c6-2fb9-53a9-8c5b-16f535851c63" -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/test/runtests.jl b/test/runtests.jl index 3fcb573..4cc73ce 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -13,11 +13,8 @@ function parse_commandline() s = ArgParse.ArgParseSettings(; description="Run tests. Default is all tests.") ArgParse.add_arg_table!(s, - ["--lightweight_utilities"], - Dict(:help => "Test only lightweight utilities", - :action => :store_true), - ["--actuator"], - Dict(:help => "Test only actuator model", + ["--unit_utils"], + Dict(:help => "Test only unit conversion utilities", :action => :store_true), ["--core_profile_extension"], Dict(:help => "Test only core profile extension", @@ -104,7 +101,6 @@ function define_default_sample_set() file_list = SD4SOLPS.find_files_in_allowed_folders( sample_path, sample_path2, sample_path3, sample_path4; eqdsk_file="thereisntoneyet", - allow_reduced_versions=false, ) b2fgmtry, b2time, b2mn, eqdsk = file_list eqdsk = @@ -113,8 +109,8 @@ function define_default_sample_set() return b2fgmtry, b2time, b2mn, eqdsk end -if args["lightweight_utilities"] - @testset "lightweight_utilities" begin +if args["unit_utils"] + @testset "Unit conversion utilities" begin # Gas unit converter flow_tls = 40.63 * Unitful.Torr * Unitful.L / Unitful.s flow_pam3 = SD4SOLPS.gas_unit_converter(flow_tls, "torr L s^-1", "Pa m^3 s^-1") @@ -147,21 +143,6 @@ if args["lightweight_utilities"] end end -if args["actuator"] - @testset "actuator" begin - t = collect(LinRange(0, 2.0, 1001)) - cmd = (t .> 1.0) * 1.55 + (t .> 1.5) * 0.93 + (t .> 1.8) * 0.25 - - instant_flow_function = SD4SOLPS.model_gas_valve("instant") - instant_flow = instant_flow_function(t, cmd) - @test length(instant_flow) == length(cmd) - - simple_flow_function = SD4SOLPS.model_gas_valve("simple") - simple_flow = simple_flow_function(t, cmd) - @test length(simple_flow) == length(cmd) - end -end - if args["core_profile_extension"] @testset "core_profile_extension" begin # Just the basic profile extrapolator ------------------ @@ -278,7 +259,7 @@ if args["heavy_utilities"] # Test for finding files in allowed folders sample_path = splitdir(pathof(SOLPS2IMAS))[1] * "/../samples" file_list = SD4SOLPS.find_files_in_allowed_folders( - sample_path; eqdsk_file="thereisntoneyet", allow_reduced_versions=true, + sample_path; eqdsk_file="thereisntoneyet", ) @test length(file_list) == 4 b2fgmtry, b2time, b2mn, eqdsk = file_list @@ -315,7 +296,11 @@ if args["heavy_utilities"] println("DD was repaired (rho added) for core 2d utility test") end density_on_grid = - SD4SOLPS.core_profile_2d(dd, prof_time_idx, eq_time_idx, quantity).(r, z) + GGDUtils.interp( + dd.core_profiles.profiles_1d[1].electrons.density, + dd.core_profiles.profiles_1d[1], + dd.equilibrium.time_slice[1], + ).(r, z) @test size(density_on_grid) == (length(rg), length(zg)) end end @@ -449,9 +434,8 @@ if args["preparation"] output_format = "json" dd = SD4SOLPS.preparation( eqdsk_file, - sample_paths...; + sample_paths; core_method=core_method, - edge_method=edge_method, filename=filename, output_format=output_format, )