FUSE Introductory Tutorial
Download this tutorial from the FuseExamples repository
Basic concepts
To make sense of this tutorial, you'll need to know the following organization concepts of FUSE:
- 📂 Data storage: All data is stored in the
dd
structure, which follows the ITER IMAS ontology. - 🧠 Actors: The core components of FUSE simulations are physics and engineering actors.
- 🕹️ Control: Actor functionality is governed by
act
parameters. - 🚀 Initialization: The data structure can be initialized from 0D
ini
parameters. - 🔧 Use cases: FUSE includes templates for various machines (e.g., FPP, ITER, ARC).
- 🔄 Workflows: Self-contained studies and optimizations are conducted via workflows, typically involving multiple FUSE simulations.
- 🌍 Interoperability: FUSE interfaces with existing modeling tools like OMFIT/OMAS and the IMAS ecosystem.
A diagram illustrating these concepts is provided below:
Let's get started!
NOTE: Julia is a Just In Time (JIT) programming language. The first time something is executed it will take longer because of the compilation process. Subsequent calls the the same code will be blazingly fast.
Import the necessary packages
using Plots # for plotting
using FUSE # this will also import IMAS in the current namespace
Starting from a use-case
FUSE comes with some predefined use-cases, some of which are used for regression testing.
Note that some use cases are for non-nuclear experiments and certain Actors like Blankets or BalanceOfPlant will not perform any actions.
Here's the list of supported use-cases. These can be customized and you will also be able to build your own.
methods(FUSE.case_parameters)
# 22 methods for generic function case_parameters from [35mFUSE[39m:- case_parameters(::Type{Val{:HDB5}}; tokamak, case, database_case) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/HDB5.jl:9
- case_parameters(::Type{Val{:D3D_machine}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/D3D.jl:318
- case_parameters(::Type{Val{:FIRST}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/FIRST.jl:6
- case_parameters(::Type{Val{:SPARC}}; init_from, flux_matcher) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/SPARC.jl:6
- case_parameters(::Type{Val{:ARC}}; flux_matcher) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/ARC.jl:6
- case_parameters(::Type{Val{:CAT}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/CAT.jl:6
- case_parameters(::Type{Val{:UNIT}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/UNIT.jl:6
- case_parameters(::Type{Val{:DTT}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/DTT.jl:6
- case_parameters(::Type{Val{:EXCITE}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/EXCITE.jl:6
- case_parameters(::Type{Val{:KDEMO_compact}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/KDEMO.jl:96
- case_parameters(::Type{Val{:KDEMO}}) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/KDEMO.jl:4
- case_parameters(::Type{Val{:D3D}}, scenario::Symbol) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/D3D.jl:262
- case_parameters(::Type{Val{:D3D}}, shot::Int64; new_impurity_match_power_rad, EFIT_tree, PROFILES_tree, omega_user, omega_omfit_root, omega_omas_root, use_local_cache) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/D3D.jl:9
- case_parameters(::Type{Val{:D3D}}, dd::IMASdd.dd) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/D3D.jl:241
- case_parameters(::Type{Val{:D3D}}, ods_file::AbstractString) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/D3D.jl:224
- case_parameters(::Type{Val{:FPP}}; flux_matcher) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/FPP.jl:6
- case_parameters(::Type{Val{:MANTA}}; flux_matcher) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/MANTA.jl:10
- case_parameters(::Type{Val{:ITER}}; init_from, boundary_from, ne_setting, time_dependent) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/ITER.jl:11
- case_parameters(::Type{Val{:STEP}}; init_from, pf_from) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/STEP.jl:6
- case_parameters(::Type{Val{:baby_MANTA}}; flux_matcher) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/baby_MANTA.jl:4
- case_parameters(case::Symbol, args...; kw...) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/_cases.jl:36
- case_parameters(data_row::DataFrames.DataFrameRow) in FUSE at /home/runner/work/FUSE.jl/FUSE.jl/src/cases/HDB5.jl:20
Get initial parameters (ini
) and actions (act
) for a given use-case, let's use KDEMO for example
ini, act = FUSE.case_parameters(:KDEMO);
The ini
data structure contains 0D parameters that will be used to bootstrap the dd
with plausible data.
The ini
parameters can be modified.
ini.equilibrium.B0 = 7.8
ini.equilibrium.R0 = 6.5;
The act
data structure contains parameters that define how the actors (ie the models) will behave.
The act
parameters can also be modified.
act.ActorCoreTransport.model = :FluxMatcher;
ini
and act
can now be used to initialize the data dictionary (dd
) using the 0D parameters.
NOTE: init()
does not return a self-consistent solution, just a plausible starting point to initialize our simulations!
dd = IMAS.dd() # an empty dd
FUSE.init(dd, ini, act);
actors: Equilibrium
actors: TEQUILA
actors: CXbuild
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: PassiveStructures
We can @checkin
and @checkout
variables with an associated tag.
This is handy to save and restore (checkpoint) our progress without having to always start from scratch (we'll use this later).
@checkin :init dd ini act
Running Actors
Let's now run a series of actors and play around with plotting to get a sense of what each individual actor does.
Here's how we can restore things back to after the initialization stage (in case we did anything else in between)
@checkout :init dd ini act
Actors in FUSE can be executed by passing two arguments to them: dd
and act
.
Let's start by positioning the PF coils, so that we stand a chance to reproduce the desired plasma shape. This will be important to ensure the stability of the ActorStationaryPlasma
that we are going to run next.
FUSE.ActorPFdesign(dd, act; do_plot=true); # instead of setting `act.ActorPFdesign.do_plot=true` we can just pass `do_plot=true` as argument without chaning `act`
actors: PFdesign
GKS: cannot open display - headless operation mode active
The ActorStationaryPlasma
iterates between plasma transport, pedestal, equilibrium and sources to return a self-consistent plasma solution
peq = plot(dd.equilibrium; label="before")
pcp = plot(dd.core_profiles; color=:gray, label="before")
FUSE.ActorStationaryPlasma(dd, act);
actors: StationaryPlasma
actors: --------------- 1/5
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 1/5 @ 671.70%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 2/5 @ 652.18%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 3/5 @ 344.44%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 4/5 @ 229.85%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 5/5 @ 90.92%
we can compare equilibrium before and after the self-consistency loop
plot!(peq, dd.equilibrium; label="after")
we can compare core_profiles before and after the self-consistency loop
plot!(pcp, dd.core_profiles; label="after")
here are the sources
plot(dd.core_sources)
and the flux-matched transport
plot(dd.core_transport)
HFS sizing actor changes the thickness of the OH and TF layers on the high field side to satisfy current and stresses constraints
plot(dd.build)
FUSE.ActorHFSsizing(dd, act);
plot!(dd.build; cx=false)
The stresses on the center stack are stored in the solid_mechanics
IDS
plot(dd.solid_mechanics.center_stack.stress)
LFS sizing actors change location of the outer TF leg to meet ripple requirements
plot(dd.build)
FUSE.ActorLFSsizing(dd, act);
plot!(dd.build; cx=false)
A custom show()
method is defined to print the summary of dd.build.layer
dd.build.layer
23×10 DataFrame
Row │ group details type ΔR R_start R_end material area volume shape
│ String String String Float64 Float64 Float64 String Float64 Float64 String
─────┼───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
1 │ in 1.30697 0.0 1.30697 steel 20.151 99.383 rectangle
2 │ in oh 0.303946 1.30697 1.61092 nb3sn 6.99039 80.9118 rectangle
3 │ hfs tf 1.6871 1.61092 3.29802 nb3sn_kdemo 43.1068 893.261 convex hull
4 │ hfs gap tf vacuum vessel 0.0 3.29802 3.29802 vacuum 7.50531 494.647 double ellipse
5 │ hfs vacuum outer vessel 0.0983016 3.29802 3.39632 steel 3.07657 125.573 negative offset
6 │ hfs gap water 0.147452 3.39632 3.54377 water 4.50099 183.828 negative offset
7 │ hfs vacuum inner vessel 0.0983016 3.54377 3.64208 steel 2.92474 119.531 negative offset
8 │ hfs gap high temp shield vacuum vess… 0.00983016 3.64208 3.65191 vacuum 0.289134 11.8201 negative offset
9 │ hfs high temp shield 0.196603 3.65191 3.84851 steel 5.65514 231.326 negative offset
10 │ hfs blanket 0.373546 3.84851 4.22205 lithium_lead 21.3676 976.315 negative offset
11 │ hfs first wall 0.0196603 4.22205 4.24172 tungsten 1.02216 33.8526 offset
12 │ lhfs plasma 4.51651 4.24172 8.75822 plasma 32.1483 1253.13
13 │ lfs first wall 0.0196603 8.75822 8.77788 tungsten 1.02203 33.8512 offset
14 │ lfs blanket 1.17962 8.77788 9.9575 lithium_lead 21.3677 976.316 negative offset
15 │ lfs high temp shield 0.196603 9.9575 10.1541 steel 5.65514 231.326 negative offset
16 │ lfs gap high temp shield vacuum vess… 0.044173 10.1541 10.1983 vacuum 0.289134 11.8201 negative offset
17 │ lfs vacuum inner vessel 0.0983016 10.1983 10.2966 steel 2.92474 119.531 negative offset
18 │ lfs gap water 0.147452 10.2966 10.444 water 4.50099 183.828 negative offset
19 │ lfs vacuum outer vessel 0.0983016 10.444 10.5423 steel 3.07657 125.573 negative offset
20 │ lfs gap tf vacuum vessel 0.7756 10.5423 11.3179 vacuum 7.50531 494.647 double ellipse
21 │ lfs tf 1.6871 11.3179 13.005 nb3sn_kdemo 43.1068 893.261 convex hull
22 │ out 1.96603 0.0 14.9711 vacuum 164.962 8444.88
23 │ out cryostat 0.0983016 0.0 15.0694 steel 4.88694 311.924 silo
ActorHFSsizing and ActorLFSsizing only change the layer's thicknesses, so we then need to trigger a build of the 2D cross-sections after them:
FUSE.ActorCXbuild(dd, act);
plot(dd.build)
Generate passive structures information (for now the vacuum vessel)
FUSE.ActorPassiveStructures(dd, act)
plot(dd.pf_passive)
We can now give the PF coils their final position given the new build
actor = FUSE.ActorPFdesign(dd, act);
plot(actor) # some actors define their own plot
With information about both pfactive and pfpassive we can now evaluate vertical stability
FUSE.ActorVerticalStability(dd, act)
IMAS.freeze(dd.mhd_linear)
mhd_linear [DETACHED]
├─ time ➡ [0] [s]
└─ time_slice
└─ 1
├─ time ➡ 0 [s]
└─ toroidal_mode
├─ 1
│ ├─ n_tor ➡ 0
│ ├─ perturbation_type
│ │ ├─ description ➡ "Vertical stability margin > 0.15 for stability"
│ │ └─ name ➡ "m_s"
│ └─ stability_metric ➡ 0.0890334
└─ 2
├─ n_tor ➡ 0
├─ perturbation_type
│ ├─ description ➡ "Normalized vertical growth rate < 10 for stability"
│ └─ name ➡ "γτ"
└─ stability_metric ➡ 14.1594
The ActorNeutronics
calculates the heat flux on the first wall
FUSE.ActorNeutronics(dd, act);
p = plot(; layout=2, size=(900, 350))
plot!(p, dd.neutronics.time_slice[].wall_loading, subplot=1)
plot!(p, FUSE.define_neutrons(dd, 100000)[1], dd.equilibrium.time_slice[]; subplot=1, colorbar_entry=false)
plot!(p, dd.neutronics.time_slice[].wall_loading; cx=false, subplot=2, ylabel="")
The ActorBlanket
will change the thickess of the first wall, breeder, shield, and Li6 enrichment to achieve target TBR
FUSE.ActorBlanket(dd, act);
print_tree(IMAS.freeze(dd.blanket); maxdepth=5)
actors: Blanket
blanket [DETACHED]
├─ module
│ └─ 1
│ ├─ layer
│ │ ├─ 1
│ │ │ ├─ material ➡ "tungsten"
│ │ │ ├─ midplane_thickness ➡ 0.0199661 [m]
│ │ │ └─ name ➡ "lfs first wall"
│ │ ├─ 2
│ │ │ ├─ material ➡ "lithium-lead: Li6/7=90.000%"
│ │ │ ├─ midplane_thickness ➡ 1.28957 [m]
│ │ │ └─ name ➡ "lfs blanket"
│ │ └─ 3
│ │ ├─ material ➡ "steel"
│ │ ├─ midplane_thickness ➡ 0.0863486 [m]
│ │ └─ name ➡ "lfs high temp shield"
│ ├─ name ➡ "blanket"
│ └─ time_slice
│ └─ 1
│ ├─ peak_escape_flux ➡ 166364 [W/m^2]
│ ├─ peak_wall_flux ➡ 926975 [W/m^2]
│ ├─ power_incident_neutrons ➡ 8.96545e+06 [W]
│ ├─ power_incident_radiated ➡ 0 [W]
│ ├─ power_thermal_extracted ➡ 1.07585e+07 [W]
│ ├─ power_thermal_neutrons ➡ 1.07585e+07 [W]
│ ├─ power_thermal_radiated ➡ 0 [W]
│ ├─ time ➡ 0 [s]
│ └─ tritium_breeding_ratio ➡ 1.67697
├─ time ➡ [0] [s]
└─ tritium_breeding_ratio ➡ [0.0671234]
The ActorDivertors
actor calculates the divertors heat flux
FUSE.ActorDivertors(dd, act);
print_tree(IMAS.freeze(dd.divertors); maxdepth=4)
actors: Divertors
divertors [DETACHED]
├─ divertor
│ └─ 1
│ ├─ power_black_body
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_conducted
│ │ ├─ data ➡ [1.20932e+08] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_convected
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_currents
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_incident
│ │ ├─ data ➡ [3.41223e+07] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_neutrals
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_radiated
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_recombination_neutrals
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_recombination_plasma
│ │ ├─ data ➡ [0] [W]
│ │ └─ time ➡ [0] [s]
│ ├─ power_thermal_extracted
│ │ ├─ data ➡ [3.41223e+07] [W]
│ │ └─ time ➡ [0] [s]
│ └─ target
│ ├─ 1
│ │ ⋮
│ │
│ └─ 2
│ ⋮
│
└─ time ➡ [0] [s]
The ActorBalanceOfPlant
calculates the optimal cooling flow rates for the heat sources (breeder, divertor, and wall) and get an efficiency for the electricity conversion cycle
FUSE.ActorBalanceOfPlant(dd, act);
IMAS.freeze(dd.balance_of_plant)
balance_of_plant [DETACHED]
├─ Q_plant ➡ [0.126023]
├─ power_electric_net ➡ [-1.26727e+08] [W]
├─ power_electric_plant_operation
│ ├─ system
│ │ ├─ 1
│ │ │ ├─ index ➡ 1
│ │ │ ├─ name ➡ "HCD"
│ │ │ ├─ power ➡ [1e+08] [W]
│ │ │ └─ subsystem
│ │ │ ├─ 1
│ │ │ │ ├─ index ➡ 1
│ │ │ │ ├─ name ➡ "nbi"
│ │ │ │ └─ power ➡ [0] [W]
│ │ │ ├─ 2
│ │ │ │ ├─ index ➡ 2
│ │ │ │ ├─ name ➡ "ec_launchers"
│ │ │ │ └─ power ➡ [5e+07] [W]
│ │ │ ├─ 3
│ │ │ │ ├─ index ➡ 3
│ │ │ │ ├─ name ➡ "ic_antennas"
│ │ │ │ └─ power ➡ [5e+07] [W]
│ │ │ └─ 4
│ │ │ ├─ index ➡ 4
│ │ │ ├─ name ➡ "lh_antennas"
│ │ │ └─ power ➡ [0] [W]
│ │ ├─ 2
│ │ │ ├─ index ➡ 3
│ │ │ ├─ name ➡ "cryostat"
│ │ │ └─ power ➡ [3e+07] [W]
│ │ ├─ 3
│ │ │ ├─ index ➡ 4
│ │ │ ├─ name ➡ "tritium_handling"
│ │ │ └─ power ➡ [1.5e+07] [W]
│ │ └─ 4
│ │ ├─ index ➡ 6
│ │ ├─ name ➡ "pf_active"
│ │ └─ power ➡ [0] [W]
│ └─ total_power ➡ [1.45e+08] [W]
├─ power_plant
│ ├─ heat_load
│ │ ├─ breeder ➡ [1.07585e+07] [W]
│ │ ├─ divertor ➡ [3.41223e+07] [W]
│ │ └─ wall ➡ [3.55911e+07] [W]
│ ├─ power_cycle_type ➡ "rankine"
│ ├─ power_electric_generated ➡ [1.82733e+07] [W]
│ └─ total_heat_supplied ➡ [8.0472e+07] [W]
├─ thermal_efficiency_plant ➡ [0.227076]
└─ time ➡ [0] [s]
ActorCosting
will break down the capital and operational costs
FUSE.ActorCosting(dd, act)
plot(dd.costing)
Let's checkpoint our results
@checkin :manual dd ini act
Whole facility design
Here we restore the :init
checkpoint that we had previously stored. Resetting any changes to dd
, ini
, and act
that we did in the meantime.
@checkout :init dd ini act
Actors can call other actors, creating workflows. For example, the ActorWholeFacility
can be used to to get a self-consistent stationary whole facility design.
FUSE.ActorWholeFacility(dd, act);
actors: WholeFacility
actors: PFdesign
actors: StationaryPlasma
actors: --------------- 1/5
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 1/5 @ 671.70%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 2/5 @ 652.18%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 3/5 @ 344.44%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 4/5 @ 229.85%
actors: HCD
actors: SimpleEC
actors: SimpleIC
actors: NeutralFueling
actors: Current
actors: QED
actors: Pedestal
actors: EPED
actors: CoreTransport
actors: FluxMatcher
actors: Current
actors: QED
actors: Equilibrium
actors: TEQUILA
actors: --------------- 5/5 @ 90.92%
actors: HFSsizing
actors: FluxSwing
actors: Stresses
actors: LFSsizing
actors: CXbuild
actors: PFdesign
actors: Equilibrium
actors: TEQUILA
actors: CXbuild
actors: Neutronics
actors: Blanket
actors: CXbuild
actors: PassiveStructures
actors: Divertors
actors: PlasmaLimits
actors: VerticalStability
actors: TroyonBetaNN
actors: BalanceOfPlant
actors: ThermalPlant
actors: PowerNeeds
actors: Costing
actors: CostingARIES
Let's check what we got at a glance with the FUSE.digest(dd)
function:
FUSE.digest(dd)
GEOMETRY EQUILIBRIUM TEMPERATURES
──────────────────────────────────── ──────────────────────────────────── ────────────────────────────────────
R0 → 6.5 [m] B0 → 7.8 [T] Te0 → 19.3 [keV]
a → 2.01 [m] ip → 12.7 [MA] Ti0 → 18.5 [keV]
1/ϵ → 3.24 q95 → 6.11 <Te> → 8.66 [keV]
κ → 2 <Bpol> → 0.801 [T] <Ti> → 7.86 [keV]
δ → 0.588 βpol_MHD → 0.771 Te0/<Te> → 2.23
ζ → -0.0136 βtor_MHD → 0.00852 Ti0/<Ti> → 2.36
Volume → 928 [m³] βn_MHD → 1.04
Surface → 754 [m²]
DENSITIES PRESSURES TRANSPORT
──────────────────────────────────── ──────────────────────────────────── ────────────────────────────────────
ne0 → 9.19e+19 [m⁻³] P0 → 0.528 [MPa] τe → 2.37 [s]
ne_ped → 6.71e+19 [m⁻³] <P> → 0.206 [MPa] τe_exp → 1.83 [s]
ne_line → 8.32e+19 [m⁻³] P0/<P> → 2.56 H98y2 → 0.866
<ne> → 7.61e+19 [m⁻³] βn → 1.05 H98y2_exp → 0.8
ne0/<ne> → 1.21 βn_th → 1.05 Hds03 → 0.644
fGW → 0.829 Hds03_exp → 0.573
zeff_ped → 2 τα_thermalization → 0.91 [s]
<zeff> → 2 τα_slowing_down → 1.03 [s]
impurities → DT Ne20 He4
SOURCES EXHAUST CURRENTS
──────────────────────────────────── ──────────────────────────────────── ────────────────────────────────────
Pec → 50 [MW] Psol → 121 [MW] ip_bs_aux_ohm → 13 [MA]
rho0_ec → 0.56 [MW] PLH → 134 [MW] ip_ni → 6.47 [MA]
Pnbi → NaN [MW] Bpol_omp → 1.11 [T] ip_bs → 3.06 [MA]
Enbi1 → NaN [MeV] λq → 0.974 [mm] ip_aux → 3.41 [MA]
Pic → 50 [MW] qpol → 2.33e+03 [MW/m²] ip_ohm → 6.52 [MA]
Plh → NaN [MW] qpar → 1.27e+04 [MW/m²] ejima → 0.4
Paux_tot → 100 [MW] P/R0 → 18.6 [MW/m] flattop → 0.7 [Hours]
Pα → 56 [MW] PB/R0 → 145 [MW T/m]
Pohm → 0.529 [MW] PBp/R0 → 14.9 [MW T/m]
Pheat → 157 [MW] PBϵ/R0q95 → 7.33 [MW T/m]
Prad_tot → -35.6 [MW] neutrons_peak → 0.371 [MW/m²]
BOP BUILD COSTING
──────────────────────────────────── ──────────────────────────────────── ────────────────────────────────────
Pfusion → 280 [MW] PF_material → nb3sn capital_cost → 6.86 [$B]
Qfusion → 2.8 TF_material → nb3sn_kdemo levelized_CoE → Inf [$/kWh]
thermal_cycle_type → rankine OH_material → nb3sn TF_of_total → 17 [%]
thermal_efficiency_plant → 22.7 [%] TF_max_b → 15.4 [T] BOP_of_total → 1.81 [%]
thermal_efficiency_cycle → NaN [%] OH_max_b → 16 [T] blanket_of_total → 20.6 [%]
power_electric_generated → 19.1 [MW] TF_j_margin → 6.57 cryostat_of_total → 2.99 [%]
Pelectric_net → -126 [MW] OH_j_margin → 1.4
Qplant → 0.132 TF_stress_margin → 3.03
TBR → 0.0731 OH_stress_margin → 1.2
@ time = 0.0 [s]
GKS: could not find font middle.ttf
Like before we can checkpoint results for later use
@checkin :awf dd ini act
Getting into the weeds
Saving and loading data
tutorial_temp_dir = tempdir()
filename = joinpath(tutorial_temp_dir, "$(ini.general.casename).json")
"/tmp/K-DEMO.json"
When saving data to be shared outside of FUSE, one can set freeze=true
so that all expressions in the dd are evaluated and saved to file.
IMAS.imas2json(dd, filename; freeze=false, strict=false);
Load from JSON
dd1 = IMAS.json2imas(filename);
Exploring the data dictionary
- FUSE stores data following the IMAS data schema.
- The root of the data structure is
dd
, which stands for "Data Dictionary". - More details are available in the documentation.
Display part of the equilibrium data in dd
dd.equilibrium.time_slice[1].boundary
boundary
├─ elongation ➡ 1.999
├─ elongation_lower ➡ Function
├─ elongation_upper ➡ Function
├─ geometric_axis
│ ├─ r ➡ 6.49997 [m]
│ └─ z ➡ 0.0133353 [m]
├─ minor_radius ➡ 2.00734 [m]
├─ outline
│ ├─ r ➡ 243-element Vector{Float64} [m]
│ │ min:4.52 avg:6.25 max:8.49
│ └─ z ➡ 243-element Vector{Float64} [m]
│ min:-4.11 avg:-0.0562 max:3.9
├─ ovality ➡ 0.00572603
├─ squareness ➡ -0.0135544
├─ squareness_lower_inner ➡ Function
├─ squareness_lower_outer ➡ Function
├─ squareness_upper_inner ➡ Function
├─ squareness_upper_outer ➡ Function
├─ strike_point
│ ├─ 1
│ │ ├─ last_closed_flux_surface_gap ➡ 1.18089e-08 [m]
│ │ ├─ r ➡ 4.21832 [m]
│ │ └─ z ➡ -4.50887 [m]
│ └─ 2
│ ├─ last_closed_flux_surface_gap ➡ 1.18089e-08 [m]
│ ├─ r ➡ 5.50741 [m]
│ └─ z ➡ -5.06831 [m]
├─ tilt ➡ -0.00162958
├─ triangularity ➡ 0.587596
├─ triangularity_lower ➡ Function
├─ triangularity_upper ➡ Function
├─ twist ➡ 0.00165424
└─ x_point
├─ 1
│ ├─ r ➡ 5.10069 [m]
│ └─ z ➡ -4.11965 [m]
└─ 2
├─ r ➡ 4.86596 [m]
└─ z ➡ 4.34009 [m]
this can be done up to a certain depth with print_tree
print_tree(dd.equilibrium.time_slice[1].boundary; maxdepth=1)
boundary
├─ elongation ➡ 1.999
├─ elongation_lower ➡ Function
├─ elongation_upper ➡ Function
├─ geometric_axis
│ ⋮
│
├─ minor_radius ➡ 2.00734 [m]
├─ outline
│ ⋮
│
├─ ovality ➡ 0.00572603
├─ squareness ➡ -0.0135544
├─ squareness_lower_inner ➡ Function
├─ squareness_lower_outer ➡ Function
├─ squareness_upper_inner ➡ Function
├─ squareness_upper_outer ➡ Function
├─ strike_point
│ ⋮
│
├─ tilt ➡ -0.00162958
├─ triangularity ➡ 0.587596
├─ triangularity_lower ➡ Function
├─ triangularity_upper ➡ Function
├─ twist ➡ 0.00165424
└─ x_point
⋮
Plotting data from dd
FUSE uses Plots.jl
recipes for visualizing data from dd
.
This allows different plots to be shown when calling plot()
on different items in the data structure.
Learn more about Plots.jl here
For example plotting the equilibrium...
plot(dd.equilibrium)
...or the core profiles
plot(dd.core_profiles)
Whant to know what arguments can be passed? use help_plot()
function
help_plot(dd.equilibrium; core_profiles_overlay=true, levels_in=21, levels_out=5, show_secondary_separatrix=true, coordinate=:rho_tor_norm)
These plots can be composed by calling plot!()
instead of plot()
plot(dd.equilibrium; color=:gray, cx=true)
plot!(dd.build.layer)
plot!(dd.pf_active)
plot!(dd.pf_passive)
plot!(dd.pulse_schedule.position_control; color=:red)
Plotting an array...
plot(dd.core_profiles.profiles_1d[1].pressure_thermal)
...is different from plotting a field from the IDS (which plots the quantity against its coordinate and with units)
plot(dd.core_profiles.profiles_1d[1], :pressure_thermal)
Customizing plot attributes:
plot(dd.core_profiles.profiles_1d[1], :pressure_thermal; label="", linewidth=2, color=:red, labelfontsize=25)
Use findall(ids, r"...")
to search for certain fields. In Julia, string starting with r
are regular expressions.
findall(dd, r"\.psi")
[1] dd.core_profiles.profiles_1d[1].grid.psi_norm [101-element Vector{Float64}] (min:0, avg:0.468, max:1)
[2] dd.core_profiles.profiles_1d[1].grid.psi [Wb] [101-element Vector{Float64}] (min:-69.8, avg:-37.5, max:-0.883)
[3] dd.equilibrium.time_slice[1].global_quantities.psi_boundary [Wb] [Float64] (all:-0.798)
[4] dd.equilibrium.time_slice[1].global_quantities.psi_axis [Wb] [Float64] (all:-73.2)
[5] dd.equilibrium.time_slice[1].profiles_1d.psi_norm [129-element Vector{Float64}] (min:0, avg:0.5, max:1)
[6] dd.equilibrium.time_slice[1].profiles_1d.psi [Wb] [129-element Vector{Float64}] (min:-73.2, avg:-37, max:-0.798)
[7] dd.equilibrium.time_slice[1].profiles_2d[1].psi [Wb] [31×13 Matrix{Float64}] (min:-0.0395, avg:0.635, max:6.65)
[8] dd.equilibrium.time_slice[1].profiles_2d[2].psi [Wb] [66×129 Matrix{Float64}] (min:-73.2, avg:17.3, max:164)
findall(ids, r"...")
can be combined with plot()
to plot multiple fields
plot(findall(dd, r"\.psi"))
Working with time series
The IMAS data structure supports time-dependent data, and IMAS.jl provides ways to handle time data efficiently.
Each dd
has a global_time
attribute, which is used throughout FUSE and IMAS to indicate the time at which things should be operate.
dd.global_time
0.0
For the sake of demonstrating handling of time, let's add a new time_slice to the equilibrium.
NOTE: time dependent arrays of structures can be resized with resize!(ids, time0::Float64)
in addition to the usual resize!(ids, n::Int)
.
resize the time dependent array of structure
resize!(dd.equilibrium.time_slice, 1.0);
let's just populate it with the data from the previous time slice
dd.equilibrium.time_slice[2] = deepcopy(dd.equilibrium.time_slice[1]);
dd.equilibrium.time_slice[2].time = 1.0
1.0
Here we see that equilibrium has mulitiple time_slices
dd.equilibrium.time
2-element Vector{Float64}:
0.0
1.0
We can access time-dependent arrays of structures via integer index...
eqt = dd.equilibrium.time_slice[2]
eqt.time
1.0
...or at a given time, by passing the time as a floating point number (in seconds)
eqt = dd.equilibrium.time_slice[1.0]
eqt.time
1.0
NOTE: If we ask a time that is not exactly in the arrays of structures, we'll get the closest (causal!) time-slice
eqt = dd.equilibrium.time_slice[0.9]
eqt.time
eqt = dd.equilibrium.time_slice[1.1]
eqt.time
1.0
... or at the current dd.global_time
by leaving the square brackets empty []
NOTE: This is what you want to use in most situations that involve arrays of structures!
dd.global_time = 0.0
eqt = dd.equilibrium.time_slice[]
eqt.time
dd.global_time = 1.0
eqt = dd.equilibrium.time_slice[]
eqt.time
1.0
What we described above was for time-dependent arrays of structures.
The other place where time comes in, is when dealing with time-dependent arrays of data.
In this case, we can use the @ddtime
macro to manipulate these time-dependent arrays at dd.global_time
.
NOTE: Also in this case, @ddtime
will operate on the closest (causal!) time point
dd.equilibrium.vacuum_toroidal_field.b0
dd.global_time = 1.0
@ddtime(dd.equilibrium.vacuum_toroidal_field.b0 = 10.0)
dd.equilibrium.vacuum_toroidal_field.b0
dd.global_time = 0.0
@ddtime(dd.equilibrium.vacuum_toroidal_field.b0)
dd.global_time = 1.0
@ddtime(dd.equilibrium.vacuum_toroidal_field.b0)
10.0
Expressions in dd
Some fields in the data dictionary are expressions (ie. Functions). For example dd.core_profiles.profiles_1d[].pressure
is dynamically calculated as the product of thermal densities and temperature with addition of fast ions contributions
dd.global_time = 0.0
print_tree(dd.core_profiles.profiles_1d[]; maxdepth=1)
1
├─ conductivity_parallel ➡ Function [ohm^-1.m^-1]
├─ electrons
│ ⋮
│
├─ grid
│ ⋮
│
├─ ion
│ ⋮
│
├─ j_bootstrap ➡ 101-element Vector{Float64} [A/m^2]
│ min:433 avg:1.14e+05 max:5.48e+05
├─ j_non_inductive ➡ 101-element Vector{Float64} [A/m^2]
│ min:6.62e+03 avg:4.69e+05 max:1.21e+06
├─ j_ohmic ➡ 101-element Vector{Float64} [A/m^2]
│ min:982 avg:3.9e+05 max:6.66e+05
├─ j_tor ➡ 101-element Vector{Float64} [A/m^2]
│ min:7.83e+03 avg:8.52e+05 max:1.89e+06
├─ j_total ➡ 101-element Vector{Float64} [A/m^2]
│ min:8.69e+03 avg:8.59e+05 max:1.87e+06
├─ neutral
│ ⋮
│
├─ pressure ➡ Function [Pa]
├─ pressure_ion_total ➡ Function [Pa]
├─ pressure_parallel ➡ Function [Pa]
├─ pressure_perpendicular ➡ Function [Pa]
├─ pressure_thermal ➡ Function [Pa]
├─ rotation_frequency_tor_sonic ➡ 101-element Vector{Float64} [s^-1]
│ all:0
├─ t_i_average ➡ Function [eV]
├─ time ➡ 0 [s]
└─ zeff ➡ 101-element Vector{Float64}
all:2
accessing a dynamic expression, automatically evaluates it
dd.core_profiles.profiles_1d[].conductivity_parallel
101-element Vector{Float64}:
1.6902802295056438e9
1.6454277515304515e9
1.585301274775194e9
1.5227060693549705e9
1.4591628139007523e9
1.395091753049798e9
1.3405358363628967e9
1.2950311369621973e9
1.2516961922989542e9
1.214870334767901e9
⋮
9.309876449522336e7
7.780699010949866e7
6.16781928345333e7
4.6137807028042935e7
3.2350135798069518e7
2.0999535064624272e7
1.2157639079676557e7
5.413322584430922e6
643368.4050559333
In addition to evaluating expressions by accessing them, expressions in the tree can be evaluated using IMAS.freeze(ids)
NOTE: IMAS.freeze(ids, field::Symbol)
works on a single field and IMAS.refreeze!(ids, field)
forces re-evaluation of an expression. Also, IMAS.empty!(ids, field::Symbol)
can be used to revert a frozen field back into an expression.
print_tree(IMAS.freeze(dd.core_profiles.profiles_1d[1]); maxdepth=1)
profiles_1d [DETACHED]
├─ conductivity_parallel ➡ 101-element Vector{Float64} [ohm^-1.m^-1]
│ min:6.43e+05 avg:4.85e+08 max:1.69e+09
├─ electrons
│ ⋮
│
├─ grid
│ ⋮
│
├─ ion
│ ⋮
│
├─ j_bootstrap ➡ 101-element Vector{Float64} [A/m^2]
│ min:433 avg:1.14e+05 max:5.48e+05
├─ j_non_inductive ➡ 101-element Vector{Float64} [A/m^2]
│ min:6.62e+03 avg:4.69e+05 max:1.21e+06
├─ j_ohmic ➡ 101-element Vector{Float64} [A/m^2]
│ min:982 avg:3.9e+05 max:6.66e+05
├─ j_tor ➡ 101-element Vector{Float64} [A/m^2]
│ min:7.83e+03 avg:8.52e+05 max:1.89e+06
├─ j_total ➡ 101-element Vector{Float64} [A/m^2]
│ min:8.69e+03 avg:8.59e+05 max:1.87e+06
├─ neutral
│ ⋮
│
├─ pressure ➡ 101-element Vector{Float64} [Pa]
│ min:427 avg:2.76e+05 max:5.28e+05
├─ pressure_ion_total ➡ 101-element Vector{Float64} [Pa]
│ min:201 avg:1.24e+05 max:2.43e+05
├─ pressure_parallel ➡ 101-element Vector{Float64} [Pa]
│ min:142 avg:9.2e+04 max:1.76e+05
├─ pressure_perpendicular ➡ 101-element Vector{Float64} [Pa]
│ min:142 avg:9.2e+04 max:1.76e+05
├─ pressure_thermal ➡ 101-element Vector{Float64} [Pa]
│ min:427 avg:2.76e+05 max:5.28e+05
├─ rotation_frequency_tor_sonic ➡ 101-element Vector{Float64} [s^-1]
│ all:0
├─ t_i_average ➡ 101-element Vector{Float64} [eV]
│ min:83.7 avg:1.01e+04 max:1.85e+04
├─ time ➡ 0 [s]
└─ zeff ➡ 101-element Vector{Float64}
all:2
Comparing two IDSs
We can introduce a change in the dd1
and spot it with the diff
function
dd1.equilibrium.time_slice[1].time = -100.0
IMAS.diff(dd.equilibrium, dd1.equilibrium)
Dict{String, String} with 3 entries:
"time" => "length: 2 -- 1"
"vacuum_toroidal_field.b0" => "length: 2 -- 1"
"time_slice" => "length: 2 -- 1"
Summary
Snapshot of dd
in 0D quantities (evaluated at dd.global_time
).
Extract + plots saved to PDF (printed to screen if filename
is omitted). NOTE: For PDF creation to work, one may need to install of DejaVu Sans Mono
font.
filename = joinpath(tutorial_temp_dir, "$(ini.general.casename).pdf")
display(filename)
FUSE.digest(dd)#, filename)