FUSE Introductory Tutorial

Download this tutorial from the FuseExamples repository


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.

FUSE.test_cases
ARC             ini, act = FUSE.case_parameters(:ARC)
CAT             ini, act = FUSE.case_parameters(:CAT)
D3D             ini, act = FUSE.case_parameters(:D3D; scenario=:default)
D3D_Hmode       ini, act = FUSE.case_parameters(:D3D; scenario_sources=true, scenario=:H_mode)
D3D_Lmode       ini, act = FUSE.case_parameters(:D3D; scenario_sources=false, scenario=:L_mode)
DTT             ini, act = FUSE.case_parameters(:DTT)
EXCITE          ini, act = FUSE.case_parameters(:EXCITE)
FPP             ini, act = FUSE.case_parameters(:FPP)
ITER_ods        ini, act = FUSE.case_parameters(:ITER; init_from=:ods)
ITER_scalars    ini, act = FUSE.case_parameters(:ITER; init_from=:scalars)
JET_HDB5        ini, act = FUSE.case_parameters(:HDB5; case=500, tokamak=:JET)
KDEMO           ini, act = FUSE.case_parameters(:KDEMO)
KDEMO_compact   ini, act = FUSE.case_parameters(:KDEMO_compact)
MANTA           ini, act = FUSE.case_parameters(:MANTA)
SPARC           ini, act = FUSE.case_parameters(:SPARC; init_from=:ods)

Get initial parameters (ini) and actions (act) for a given use-case

ini, act = FUSE.case_parameters(:KDEMO);

Modifying ini parameters.

ini.equilibrium.B0 = 7.8
ini.equilibrium.R0 = 6.5;

Modifying act parameters.

act.ActorCoreTransport.model = :FluxMatcher;

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: HCD
actors:  SimpleEC
actors:  SimpleIC
actors:  SimpleLH
actors:  SimpleNB
actors:  SimplePellet
actors: Current
actors:  SteadyStateCurrent
actors: CXbuild
actors: PassiveStructures

We can @checkin and @checkout variables with an associated tag.

This is handy to save and restore our progress (we'll use this later).

@checkin :init dd ini act

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[2].boundary
boundary
├─ elongation ➡ 1.99926
├─ elongation_lowerFunction
├─ elongation_upperFunction
├─ geometric_axis
│  ├─ r ➡ 6.49971 [m]
│  └─ z ➡ -0.00101717 [m]
├─ minor_radius ➡ 2.00728 [m]
├─ outline
│  ├─ r239-element Vector{Float64} [m]
│  │      min:4.51   avg:6.27   max:8.49 
│  └─ z239-element Vector{Float64} [m]min:-4.04   avg:-0.0701   max:3.87 
├─ ovality ➡ 0.00308297
├─ squareness ➡ -0.00753311
├─ squareness_lower_innerFunction
├─ squareness_lower_outerFunction
├─ squareness_upper_innerFunction
├─ squareness_upper_outerFunction
├─ strike_point
│  ├─ 1
│  │  ├─ r ➡ 5.50287 [m]
│  │  └─ z ➡ -5.01275 [m]
│  └─ 2
│     ├─ r ➡ 4.22213 [m]
│     └─ z ➡ -4.51762 [m]
├─ tilt ➡ -0.00645963
├─ triangularity ➡ 0.585259
├─ triangularity_lowerFunction
├─ triangularity_upperFunction
├─ twist ➡ 0.00262375
└─ x_point
   ├─ 1
   │  ├─ r ➡ 5.11438 [m]
   │  └─ z ➡ -4.04412 [m]
   └─ 2
      ├─ r ➡ 4.87066 [m]
      └─ z ➡ 4.32838 [m]

this can be done up to a certain depth with print_tree

print_tree(dd.equilibrium.time_slice[2].boundary; maxdepth=1)
boundary
├─ elongation ➡ 1.99926
├─ elongation_lower ➡ Function
├─ elongation_upper ➡ Function
├─ geometric_axis
│  ⋮
│
├─ minor_radius ➡ 2.00728 [m]
├─ outline
│  ⋮
│
├─ ovality ➡ 0.00308297
├─ squareness ➡ -0.00753311
├─ squareness_lower_inner ➡ Function
├─ squareness_lower_outer ➡ Function
├─ squareness_upper_inner ➡ Function
├─ squareness_upper_outer ➡ Function
├─ strike_point
│  ⋮
│
├─ tilt ➡ -0.00645963
├─ triangularity ➡ 0.585259
├─ triangularity_lower ➡ Function
├─ triangularity_upper ➡ Function
├─ twist ➡ 0.00262375
└─ 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)
Example block output

...or the core profiles

plot(dd.core_profiles)
Example block output

Whant to know what arguments can be passed? use help_plot() function

help_plot(dd.equilibrium; core_profiles_overlay=true, psi_levels_in=21, psi_levels_out=5, show_secondary_separatrix=true, coordinate=:psi_norm)
Example block output

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)
Example block output

Plotting an array...

plot(dd.core_profiles.profiles_1d[1].pressure_thermal)
Example block output

...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)
Example block output

Customizing plot attributes:

plot(dd.core_profiles.profiles_1d[1], :pressure_thermal; label="", linewidth=2, color=:red, labelfontsize=25)
Example block output

Use findall(ids, r"...") to search for certain fields. In Julia string starting with r are regular expressions.

findall(dd, r"pressure")
[1] dd.core_profiles.profiles_1d[1].ion[1].pressure_fast_perpendicular [Pa] [101-element Vector{Float64}] (all:0)
[2] dd.core_profiles.profiles_1d[1].ion[1].pressure_fast_parallel [Pa] [101-element Vector{Float64}] (all:0)
[3] dd.core_profiles.profiles_1d[1].ion[2].pressure_fast_perpendicular [Pa] [101-element Vector{Float64}] (all:0)
[4] dd.core_profiles.profiles_1d[1].ion[2].pressure_fast_parallel [Pa] [101-element Vector{Float64}] (all:0)
[5] dd.core_profiles.profiles_1d[1].ion[3].pressure_fast_perpendicular [Pa] [101-element Vector{Float64}] (min:7.84e-16, avg:1.8e+04, max:5.37e+04)
[6] dd.core_profiles.profiles_1d[1].ion[3].pressure_fast_parallel [Pa] [101-element Vector{Float64}] (min:7.84e-16, avg:1.8e+04, max:5.37e+04)
[7] dd.equilibrium.time_slice[2].profiles_1d.pressure [Pa] [129-element Vector{Float64}] (min:0, avg:3.24e+05, max:7.33e+05)
[8] dd.equilibrium.time_slice[2].profiles_1d.dpressure_dpsi [Pa.Wb^-1] [129-element Vector{Float64}] (min:-1.38e+04, avg:-1.15e+04, max:-7.92)

findall(ids, r"...") can be combined with plot() to plot multiple fields

plot(findall(dd, r"\.pressure"))
Example block output

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, and actors operate at such time

dd.global_time
0.0

Here we see that equilibrium has mulitiple time_slices

dd.equilibrium.time
2-element Vector{Float64}:
 -Inf
   0.0

Accessing time-dependent arrays of structures, via integer index

eqt = dd.equilibrium.time_slice[2]
eqt.time
0.0

At a given time, by passing the time as a floating point number (in seconds)

eqt = dd.equilibrium.time_slice[0.0]
eqt.time
0.0

At the global time, leaving the square brackets empty

eqt = dd.equilibrium.time_slice[]
eqt.time
0.0

Using the @ddtime macro to access and modify time-dependent arrays at dd.global_time:

dd.equilibrium.vacuum_toroidal_field.b0
2-element Vector{Float64}:
 7.80034989842101
 7.80034989842101

Accessing data at dd.global_time

my_b0 = @ddtime(dd.equilibrium.vacuum_toroidal_field.b0)
7.80034989842101

Writin data at dd.global_time

@ddtime(dd.equilibrium.vacuum_toroidal_field.b0 = my_b0 + 1)

dd.equilibrium.vacuum_toroidal_field.b0
2-element Vector{Float64}:
 7.80034989842101
 8.800349898421011

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

print_tree(dd.core_profiles.profiles_1d[1]; maxdepth=1)
1
├─ conductivity_parallel ➡ Function [ohm^-1.m^-1]
├─ electrons
│  ⋮
│
├─ grid
│  ⋮
│
├─ ion
│  ⋮
│
├─ j_bootstrap ➡ 101-element Vector{Float64} [A/m^2]
│                min:3.33e+03   avg:1.56e+05   max:2e+05
├─ j_non_inductive ➡ 101-element Vector{Float64} [A/m^2]
│                    min:3.33e+03   avg:8.38e+05   max:4.98e+06
├─ j_ohmic ➡ 101-element Vector{Float64} [A/m^2]
│            min:535   avg:4.54e+05   max:1.53e+06
├─ j_tor ➡ 101-element Vector{Float64} [A/m^2]
│          min:3.25e+03   avg:1.31e+06   max:6.65e+06
├─ j_total ➡ 101-element Vector{Float64} [A/m^2]
│            min:3.87e+03   avg:1.29e+06   max:6.48e+06
├─ 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 (in the pressure example, we get an array with data)

dd.core_profiles.profiles_1d[1].electrons.pressure
101-element Vector{Float64}:
 382381.5183953427
 381049.9804102575
 379610.10311554925
 378000.1428571452
 376220.778286027
 374276.7775090226
 372173.5159291897
 369916.29381371004
 367510.18551253923
 364960.023529924
      ⋮
  19923.746255848662
  18510.116050993194
  16824.74656249592
  14486.796954063524
  11100.203362314784
   6820.574494367021
   3122.7927779947017
   1038.701917412065
    218.50372479733514

In addition to evaluating expressions by accessing them, expressions in the tree can be evaluated using IMAS.freeze()

print_tree(IMAS.freeze(dd.core_profiles.profiles_1d[1]); maxdepth=1)
profiles_1d
├─ conductivity_parallel ➡ 101-element Vector{Float64} [ohm^-1.m^-1]
│                          min:2.88e+06   avg:3.62e+09   max:1.25e+10
├─ electrons
│  ⋮
│
├─ grid
│  ⋮
│
├─ ion
│  ⋮
│
├─ j_bootstrap ➡ 101-element Vector{Float64} [A/m^2]
│                min:3.33e+03   avg:1.56e+05   max:2e+05
├─ j_non_inductive ➡ 101-element Vector{Float64} [A/m^2]
│                    min:3.33e+03   avg:8.38e+05   max:4.98e+06
├─ j_ohmic ➡ 101-element Vector{Float64} [A/m^2]
│            min:535   avg:4.54e+05   max:1.53e+06
├─ j_tor ➡ 101-element Vector{Float64} [A/m^2]
│          min:3.25e+03   avg:1.31e+06   max:6.65e+06
├─ j_total ➡ 101-element Vector{Float64} [A/m^2]
│            min:3.87e+03   avg:1.29e+06   max:6.48e+06
├─ pressure ➡ 101-element Vector{Float64} [Pa]
│             min:413   avg:4.09e+05   max:8.85e+05
├─ pressure_ion_total ➡ 101-element Vector{Float64} [Pa]
│                       min:195   avg:1.67e+05   max:3.41e+05
├─ pressure_parallel ➡ 101-element Vector{Float64} [Pa]
│                      min:138   avg:1.36e+05   max:2.95e+05
├─ pressure_perpendicular ➡ 101-element Vector{Float64} [Pa]
│                           min:138   avg:1.36e+05   max:2.95e+05
├─ pressure_thermal ➡ 101-element Vector{Float64} [Pa]
│                     min:413   avg:3.55e+05   max:7.23e+05
├─ rotation_frequency_tor_sonic ➡ 101-element Vector{Float64} [s^-1]
│                                 all:0
├─ t_i_average ➡ 101-element Vector{Float64} [eV]
│                min:80   avg:1.32e+04   max:2.5e+04
├─ time ➡ 0 [s]
└─ zeff ➡ 101-element Vector{Float64}
          all:2

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 in FUSE can be executed by passing two arguments to them: dd and act. Internally, actors can call other actors, creating workflows. For example, the ActorWholeFacility can be used to to get a self-consistent stationary whole facility design. The actors: print statements with their nested output tell us what actors are calling other actors.

FUSE.ActorWholeFacility(dd, act);
actors: WholeFacility
actors:  PFdesign
actors:  StationaryPlasma
actors:   --------------- 1/5
actors:   HCD
actors:    SimpleEC
actors:    SimpleIC
actors:    SimpleLH
actors:    SimpleNB
actors:    SimplePellet
actors:   Current
actors:    SteadyStateCurrent
actors:   Pedestal
actors:    EPED
actors:   CoreTransport
actors:    FluxMatcher
actors:   Current
actors:    SteadyStateCurrent
actors:   Equilibrium
actors:    TEQUILA
actors:   --------------- 1/5 @ 456.12%
actors:   HCD
actors:    SimpleEC
actors:    SimpleIC
actors:    SimpleLH
actors:    SimpleNB
actors:    SimplePellet
actors:   Current
actors:    SteadyStateCurrent
actors:   Pedestal
actors:    EPED
actors:   CoreTransport
actors:    FluxMatcher
actors:   Current
actors:    SteadyStateCurrent
actors:   Equilibrium
actors:    TEQUILA
actors:   --------------- 2/5 @ 247.35%
actors:   HCD
actors:    SimpleEC
actors:    SimpleIC
actors:    SimpleLH
actors:    SimpleNB
actors:    SimplePellet
actors:   Current
actors:    SteadyStateCurrent
actors:   Pedestal
actors:    EPED
actors:   CoreTransport
actors:    FluxMatcher
actors:   Current
actors:    SteadyStateCurrent
actors:   Equilibrium
actors:    TEQUILA
actors:   --------------- 3/5 @ 102.32%
actors:   HCD
actors:    SimpleEC
actors:    SimpleIC
actors:    SimpleLH
actors:    SimpleNB
actors:    SimplePellet
actors:   Current
actors:    SteadyStateCurrent
actors:   Pedestal
actors:    EPED
actors:   CoreTransport
actors:    FluxMatcher
actors:   Current
actors:    SteadyStateCurrent
actors:   Equilibrium
actors:    TEQUILA
actors:   --------------- 4/5 @ 80.31%
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
┌ Warning: act.ActorPlasmaLimits.models = [:vertical_stability, :beta_troyon_nn, :q95_gt_2, :gw_density, :κ_controllability]
│
│ Some stability rules exceed their limit threshold:
│  193% of Vertical stability margin > 0.15 for stability
│  163% of Normalized vertical growth rate < 10 for stability
│  109% of βn < BetaTroyonNN n=1
│
│ Some stability rules satisfy their limit threshold:
│  32% of βn < BetaTroyonNN n=2
│  33% of βn < BetaTroyonNN n=3
│  31% of q(rho=0.95) > 2.0
│  92% of greenwald_fraction < 1.0
│  82% of elongation < IMAS.elongation_limit
└ @ FUSE ~/work/FUSE.jl/FUSE.jl/src/actors/stability/limits_actor.jl:92
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 → 22.2 [keV]
a → 2.01 [m]                            ip → 12.7 [MA]                          Ti0 → 21.5 [keV]
1/ϵ → 3.24                              q95 → 6.46                              <Te> → 9.2 [keV]
κ → 2                                   <Bpol> → 0.799 [T]                      <Ti> → 8.34 [keV]
δ → 0.585                               βpol_MHD → 0.945                        Te0/<Te> → 2.42
ζ → -0.00753                            βtor_MHD → 0.0102                       Ti0/<Ti> → 2.58
Volume → 927 [m³]                       βn_MHD → 1.25
Surface → 760 [m²]

DENSITIES                               PRESSURES                               TRANSPORT
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
ne0 → 1.06e+20 [m⁻³]                    P0 → 0.839 [MPa]                        τe → 1.88 [s]
ne_ped → 6.93e+19 [m⁻³]                 <P> → 0.248 [MPa]                       τe_exp → 2.5 [s]
ne_line → 9.19e+19 [m⁻³]                P0/<P> → 3.38                           H98y2 → 0.849
<ne> → 8.28e+19 [m⁻³]                   βn → 1.23                               H98y2_exp → 0.927
ne0/<ne> → 1.28                         βn_th → 1.17                            Hds03 → 0.597
fGW → 0.918                                                                     Hds03_exp → 0.678
zeff_ped → 2                                                                    τα_thermalization → 0.805 [s]
<zeff> → 2                                                                      τα_slowing_down → 1.1 [s]
impurities → DT Ne20 He4

SOURCES                                 EXHAUST                                 CURRENTS
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
Pec → 50 [MW]                           Psol → 131 [MW]                         ip_bs_aux_ohm → 13 [MA]
rho0_ec → 0.5 [MW]                      PLH → 71.8 [MW]                         ip_ni → 6.64 [MA]
Pnbi → NaN [MW]                         Bpol_omp → 1.2 [T]                      ip_bs → 3.2 [MA]
Enbi1 → NaN [MeV]                       λq → 0.904 [mm]                         ip_aux → 3.44 [MA]
Pic → 50 [MW]                           qpol → 2.71e+03 [MW/m²]                 ip_ohm → 6.37 [MA]
Plh → NaN [MW]                          qpar → 1.37e+04 [MW/m²]                 ejima → 0.4
Paux_tot → 100 [MW]                     P/R0 → 20.1 [MW/m]                      flattop → 0.7 [Hours]
Pα → 73.7 [MW]                          PB/R0 → 157 [MW T/m]
Pohm → 0.0942 [MW]                      PBp/R0 → 16.1 [MW T/m]
Pheat → 174 [MW]                        PBϵ/R0q95 → 7.51 [MW T/m]
Prad_tot → -42.9 [MW]                   neutrons_peak → 0.502 [MW/m²]

BOP                                     BUILD                                   COSTING
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
Pfusion → 368 [MW]                      PF_material → nb3sn                     capital_cost → 6.43 [$B]
Qfusion → 3.68                          TF_material → nb3sn_kdemo               levelized_CoE → Inf [$/kWh]
thermal_cycle_type → rankine            OH_material → nb3sn                     TF_of_total → 15.1 [%]
thermal_efficiency_plant → 22.7 [%]     TF_max_b → 15.4 [T]                     BOP_of_total → 2.08 [%]
thermal_efficiency_cycle → NaN [%]      OH_max_b → 19.9 [T]                     blanket_of_total → 23.6 [%]
power_electric_generated → 21.6 [MW]    TF_j_margin → 5.9                       cryostat_of_total → 3.07 [%]
Pelectric_net → -123 [MW]               OH_j_margin → 1.71
Qplant → 0.149                          TF_stress_margin → 3.15
TBR → 0.0662                            OH_stress_margin → 1.37

@ 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

Running a custom workflow

Let's now run a series of actors similar to what ActorWholeFacility does and play around with plotting to get a sense of what each individual actor does.

Let's start again from after the initialization stage

@checkout :init dd ini 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

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:   SimpleLH
actors:   SimpleNB
actors:   SimplePellet
actors:  Current
actors:   SteadyStateCurrent
actors:  Pedestal
actors:   EPED
actors:  CoreTransport
actors:   FluxMatcher
actors:  Current
actors:   SteadyStateCurrent
actors:  Equilibrium
actors:   TEQUILA
actors:  --------------- 1/5 @ 456.12%
actors:  HCD
actors:   SimpleEC
actors:   SimpleIC
actors:   SimpleLH
actors:   SimpleNB
actors:   SimplePellet
actors:  Current
actors:   SteadyStateCurrent
actors:  Pedestal
actors:   EPED
actors:  CoreTransport
actors:   FluxMatcher
actors:  Current
actors:   SteadyStateCurrent
actors:  Equilibrium
actors:   TEQUILA
actors:  --------------- 2/5 @ 247.35%
actors:  HCD
actors:   SimpleEC
actors:   SimpleIC
actors:   SimpleLH
actors:   SimpleNB
actors:   SimplePellet
actors:  Current
actors:   SteadyStateCurrent
actors:  Pedestal
actors:   EPED
actors:  CoreTransport
actors:   FluxMatcher
actors:  Current
actors:   SteadyStateCurrent
actors:  Equilibrium
actors:   TEQUILA
actors:  --------------- 3/5 @ 102.32%
actors:  HCD
actors:   SimpleEC
actors:   SimpleIC
actors:   SimpleLH
actors:   SimpleNB
actors:   SimplePellet
actors:  Current
actors:   SteadyStateCurrent
actors:  Pedestal
actors:   EPED
actors:  CoreTransport
actors:   FluxMatcher
actors:  Current
actors:   SteadyStateCurrent
actors:  Equilibrium
actors:   TEQUILA
actors:  --------------- 4/5 @ 80.31%

we can compare equilibrium before and after the self-consistency loop

plot!(peq, dd.equilibrium; label="after")
Example block output

we can compare core_profiles before and after the self-consistency loop

plot!(pcp, dd.core_profiles; label="after")
Example block output

here are the sources

plot(dd.core_sources)
Example block output

and the flux-matched transport

plot(dd.core_transport)
Example block output

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)
Example block output

The stresses on the center stack are stored in the solid_mechanics IDS

plot(dd.solid_mechanics.center_stack.stress)
Example block output

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)
Example block output

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                                                   0.185414    0.0        0.185414  steel          19.9166      98.2221
   2 │ in                                         oh        1.66873     0.185414   1.85414   nb3sn           6.90907     79.9666
   3 │ hfs                                        tf        1.44372     1.85414    3.29786   nb3sn_kdemo    42.7593     886.02    convex hull
   4 │ hfs     gap tf vacuum vessel                         0.0         3.29786    3.29786   vacuum          7.40954    486.553   double ellipse
   5 │ hfs     vacuum  outer                      vessel    0.098297    3.29786    3.39616   steel           3.04733    123.995   negative offset
   6 │ hfs     gap water                                    0.147445    3.39616    3.54361   water           4.45714    181.478   negative offset
   7 │ hfs     vacuum  inner                      vessel    0.098297    3.54361    3.6419    steel           2.89552    117.976   negative offset
   8 │ hfs     gap high temp shield vacuum vess…            0.0098297   3.6419     3.65173   vacuum          0.286212    11.6651  negative offset
   9 │ hfs     high temp                          shield    0.196594    3.65173    3.84833   steel           5.5967     228.246   negative offset
  10 │ hfs                                        blanket   0.373528    3.84833    4.22185   lithium_lead   20.5289     935.506   negative offset
  11 │ hfs     first                              wall      0.0196594   4.22185    4.24151   tungsten        0.978743    32.4571  offset
  12 │ lhfs                                       plasma    4.51639     4.24151    8.7579    plasma         32.0006    1247.4
  13 │ lfs     first                              wall      0.0196594   8.7579     8.77756   tungsten        0.978736    32.4574  offset
  14 │ lfs                                        blanket   1.17956     8.77756    9.95713   lithium_lead   20.5289     935.505   negative offset
  15 │ lfs     high temp                          shield    0.196594    9.95713   10.1537    steel           5.5967     228.246   negative offset
  16 │ lfs     gap high temp shield vacuum vess…            0.0441988  10.1537    10.1979    vacuum          0.286212    11.6651  negative offset
  17 │ lfs     vacuum  inner                      vessel    0.098297   10.1979    10.2962    steel           2.89552    117.976   negative offset
  18 │ lfs     gap water                                    0.147445   10.2962    10.4437    water           4.45714    181.478   negative offset
  19 │ lfs     vacuum  outer                      vessel    0.098297   10.4437    10.542     steel           3.04733    123.995   negative offset
  20 │ lfs     gap tf vacuum vessel                         0.775563   10.542     11.3175    vacuum          7.40954    486.553   double ellipse
  21 │ lfs                                        tf        1.44372    11.3175    12.7612    nb3sn_kdemo    42.7593     886.02    convex hull
  22 │ out                                                  1.96594     0.0       14.7272    vacuum        163.122     8312.78
  23 │ out                                        cryostat  0.098297    0.0       14.8255    steel           4.8503     308.042   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)
Example block output

Generate passive structures information (for now the vacuum vessel)

FUSE.ActorPassiveStructures(dd, act)
plot(dd.pf_passive)
Example block output

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
Example block output

With information about both pfactive and pfpassive we can now evaluate vertical stability

FUSE.ActorVerticalStability(dd, act)
IMAS.freeze(dd.mhd_linear)
mhd_linear
├─ time[0] [s]
└─ time_slice
   └─ 1
      ├─ time ➡ 0 [s]
      └─ toroidal_mode
         ├─ 1
         │  ├─ n_tor0
         │  ├─ perturbation_type
         │  │  ├─ description"Vertical stability margin > 0.15 for stability"
         │  │  └─ name"m_s"
         │  └─ stability_metric ➡ 0.081577
         └─ 2
            ├─ n_tor0
            ├─ perturbation_type
            │  ├─ description"Normalized vertical growth rate < 10 for stability"
            │  └─ name"γτ"
            └─ stability_metric ➡ 15.4671

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="")
Example block output

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
├─ module
│  └─ 1
│     ├─ layer
│     │  ├─ 1
│     │  │  ├─ material ➡ "tungsten"
│     │  │  ├─ midplane_thickness ➡ 0.0199664 [m]
│     │  │  └─ name ➡ "lfs first wall"
│     │  ├─ 2
│     │  │  ├─ material ➡ "lithium-lead: Li6/7=90.000%"
│     │  │  ├─ midplane_thickness ➡ 1.35585 [m]
│     │  │  └─ name ➡ "lfs blanket"
│     │  └─ 3
│     │     ├─ material ➡ "steel"
│     │     ├─ midplane_thickness ➡ 0.0200027 [m]
│     │     └─ name ➡ "lfs high temp shield"
│     ├─ name ➡ "blanket"
│     └─ time_slice
│        └─ 1
│           ├─ peak_escape_flux ➡ 215662 [W/m^2]
│           ├─ peak_wall_flux ➡ 1.22796e+06 [W/m^2]
│           ├─ power_incident_neutrons ➡ 1.25268e+07 [W]
│           ├─ power_incident_radiated ➡ 0 [W]
│           ├─ power_thermal_extracted ➡ 1.50321e+07 [W]
│           ├─ power_thermal_neutrons ➡ 1.50321e+07 [W]
│           ├─ power_thermal_radiated ➡ 0 [W]
│           ├─ time ➡ 0 [s]
│           └─ tritium_breeding_ratio ➡ 1.54136
├─ time ➡ [0] [s]
└─ tritium_breeding_ratio ➡ [0.0655252]

The ActorDivertors actor calculates the divertors heat flux

FUSE.ActorDivertors(dd, act);
print_tree(IMAS.freeze(dd.divertors); maxdepth=4)
actors: Divertors
divertors
├─ divertor
│  └─ 1
│     ├─ power_conducted
│     │  ├─ data ➡ [1.30866e+08] [W]
│     │  └─ time ➡ [0] [s]
│     ├─ power_convected
│     │  ├─ data ➡ [0] [W]
│     │  └─ time ➡ [0] [s]
│     ├─ power_incident
│     │  ├─ data ➡ [3.21524e+07] [W]
│     │  └─ time ➡ [0] [s]
│     ├─ power_thermal_extracted
│     │  ├─ data ➡ [3.21524e+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
├─ power_electric_plant_operation
│  ├─ system
│  │  ├─ 1
│  │  │  ├─ index1
│  │  │  ├─ name"HCD"
│  │  │  ├─ power[1e+08] [W]
│  │  │  └─ subsystem
│  │  │     ├─ 1
│  │  │     │  ├─ index1
│  │  │     │  ├─ name"nbi"
│  │  │     │  └─ power[0] [W]
│  │  │     ├─ 2
│  │  │     │  ├─ index2
│  │  │     │  ├─ name"ec_launchers"
│  │  │     │  └─ power[5e+07] [W]
│  │  │     ├─ 3
│  │  │     │  ├─ index3
│  │  │     │  ├─ name"ic_antennas"
│  │  │     │  └─ power[5e+07] [W]
│  │  │     └─ 4
│  │  │        ├─ index4
│  │  │        ├─ name"lh_antennas"
│  │  │        └─ power[0] [W]
│  │  ├─ 2
│  │  │  ├─ index3
│  │  │  ├─ name"cryostat"
│  │  │  └─ power[3e+07] [W]
│  │  ├─ 3
│  │  │  ├─ index4
│  │  │  ├─ name"tritium_handling"
│  │  │  └─ power[1.5e+07] [W]
│  │  └─ 4
│  │     ├─ index6
│  │     ├─ name"pf_active"
│  │     └─ power[0] [W]
│  └─ total_power[1.45e+08] [W]
├─ power_plant
│  ├─ heat_load
│  │  ├─ breeder[1.50321e+07] [W]
│  │  ├─ divertor[3.21524e+07] [W]
│  │  └─ wall[4.28941e+07] [W]
│  ├─ power_cycle_type"rankine"
│  ├─ power_electric_generated[2.04501e+07] [W]
│  └─ total_heat_supplied[9.00786e+07] [W]
├─ thermal_efficiency_plant[0.227025]
└─ time[0] [s]

ActorCosting will break down the capital and operational costs

FUSE.ActorCosting(dd, act)
plot(dd.costing)
Example block output

Let's checkpoint our results

@checkin :manual dd ini act

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);

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 1 entry:
  "time_slice[1].time" => "value:  -Inf --  -100.0"

Summary

Snapshot of dd in 0D quantities (evaluated at dd.global_time)

FUSE.extract(dd)
GEOMETRY                                EQUILIBRIUM                             TEMPERATURES                            
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────    
R0 → 6.5 [m]                            B0 → 7.8 [T]                            Te0 → 22.2 [keV]                        
a → 2.01 [m]                            ip → 12.7 [MA]                          Ti0 → 21.5 [keV]                        
1/ϵ → 3.24                              q95 → 6.64                              <Te> → 9.2 [keV]                        
κ → 2                                   <Bpol> → 0.798 [T]                      <Ti> → 8.34 [keV]                       
δ → 0.585                               βpol_MHD → 0.95                         Te0/<Te> → 2.42                         
ζ → -0.00753                            βtor_MHD → 0.0102                       Ti0/<Ti> → 2.58                         
Volume → 929 [m³]                       βn_MHD → 1.25                                                                   
Surface → 761 [m²]                                                                                                      
                                                                                                                        
DENSITIES                               PRESSURES                               TRANSPORT                               
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────    
ne0 → 1.06e+20 [m⁻³]                    P0 → 0.839 [MPa]                        τe → 1.88 [s]                           
ne_ped → 6.93e+19 [m⁻³]                 <P> → 0.248 [MPa]                       τe_exp → 2.5 [s]                        
ne_line → 9.19e+19 [m⁻³]                P0/<P> → 3.38                           H98y2 → 0.85                            
<ne> → 8.28e+19 [m⁻³]                   βn → 1.23                               H98y2_exp → 0.928                       
ne0/<ne> → 1.28                         βn_th → 1.17                            Hds03 → 0.597                           
fGW → 0.918                                                                     Hds03_exp → 0.679                       
zeff_ped → 2                                                                    τα_thermalization → 0.805 [s]           
<zeff> → 2                                                                      τα_slowing_down → 1.1 [s]               
impurities → DT Ne20 He4                                                                                                
                                                                                                                        
SOURCES                                 EXHAUST                                 CURRENTS                                
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────    
Pec → 50 [MW]                           Psol → 131 [MW]                         ip_bs_aux_ohm → 13 [MA]                 
rho0_ec → 0.5 [MW]                      PLH → 71.9 [MW]                         ip_ni → 6.64 [MA]                       
PnbiNaN [MW]                         Bpol_omp → 1.18 [T]                     ip_bs → 3.2 [MA]                        
Enbi1NaN [MeV]                       λq → 0.922 [mm]                         ip_aux → 3.44 [MA]                      
Pic → 50 [MW]                           qpol → 2.66e+03 [MW/m²]                 ip_ohm → 6.36 [MA]                      
PlhNaN [MW]                          qpar → 1.37e+04 [MW/m²]                 ejima → 0.4                             
Paux_tot → 100 [MW]                     P/R0 → 20.1 [MW/m]                      flattop → 0.7 [Hours]                   
 → 73.7 [MW]                          PB/R0 → 157 [MW T/m]                                                            
Pohm → 0.0939 [MW]                      PBp/R0 → 16.1 [MW T/m]                                                          
Pheat → 174 [MW]                        PBϵ/R0q95 → 7.31 [MW T/m]                                                       
Prad_tot → -42.9 [MW]                   neutrons_peak → 0.501 [MW/m²]                                                   
                                                                                                                        
BOP                                     BUILD                                   COSTING                                 
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────    
Pfusion → 368 [MW]                      PF_material → nb3sn                     capital_cost → 6.33 [$B]                
Qfusion → 3.68                          TF_material → nb3sn_kdemo               levelized_CoE → Inf [$/kWh]             
thermal_cycle_type → rankine            OH_material → nb3sn                     TF_of_total → 15.5 [%]                  
thermal_efficiency_plant → 22.7 [%]     TF_max_b → 15.4 [T]                     BOP_of_total → 2.04 [%]                 
thermal_efficiency_cycleNaN [%]      OH_max_b → 19.9 [T]                     blanket_of_total → 20.3 [%]             
power_electric_generated → 20.5 [MW]    TF_j_margin → 5.9                       cryostat_of_total → 3.13 [%]            
Pelectric_net → -125 [MW]               OH_j_margin → 1.71                                                              
Qplant → 0.141                          TF_stress_margin → 3.15                                                         
TBR → 0.0655                            OH_stress_margin → 1.37                                                         
                                                                                                                        

Extract + plots saved to PDF (printed to screen it filename is omitted)

filename = joinpath(tutorial_temp_dir, "$(ini.general.casename).pdf")
FUSE.digest(dd)#, filename)
GEOMETRY                                EQUILIBRIUM                             TEMPERATURES
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
R0 → 6.5 [m]                            B0 → 7.8 [T]                            Te0 → 22.2 [keV]
a → 2.01 [m]                            ip → 12.7 [MA]                          Ti0 → 21.5 [keV]
1/ϵ → 3.24                              q95 → 6.64                              <Te> → 9.2 [keV]
κ → 2                                   <Bpol> → 0.798 [T]                      <Ti> → 8.34 [keV]
δ → 0.585                               βpol_MHD → 0.95                         Te0/<Te> → 2.42
ζ → -0.00753                            βtor_MHD → 0.0102                       Ti0/<Ti> → 2.58
Volume → 929 [m³]                       βn_MHD → 1.25
Surface → 761 [m²]

DENSITIES                               PRESSURES                               TRANSPORT
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
ne0 → 1.06e+20 [m⁻³]                    P0 → 0.839 [MPa]                        τe → 1.88 [s]
ne_ped → 6.93e+19 [m⁻³]                 <P> → 0.248 [MPa]                       τe_exp → 2.5 [s]
ne_line → 9.19e+19 [m⁻³]                P0/<P> → 3.38                           H98y2 → 0.85
<ne> → 8.28e+19 [m⁻³]                   βn → 1.23                               H98y2_exp → 0.928
ne0/<ne> → 1.28                         βn_th → 1.17                            Hds03 → 0.597
fGW → 0.918                                                                     Hds03_exp → 0.679
zeff_ped → 2                                                                    τα_thermalization → 0.805 [s]
<zeff> → 2                                                                      τα_slowing_down → 1.1 [s]
impurities → DT Ne20 He4

SOURCES                                 EXHAUST                                 CURRENTS
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
Pec → 50 [MW]                           Psol → 131 [MW]                         ip_bs_aux_ohm → 13 [MA]
rho0_ec → 0.5 [MW]                      PLH → 71.9 [MW]                         ip_ni → 6.64 [MA]
Pnbi → NaN [MW]                         Bpol_omp → 1.18 [T]                     ip_bs → 3.2 [MA]
Enbi1 → NaN [MeV]                       λq → 0.922 [mm]                         ip_aux → 3.44 [MA]
Pic → 50 [MW]                           qpol → 2.66e+03 [MW/m²]                 ip_ohm → 6.36 [MA]
Plh → NaN [MW]                          qpar → 1.37e+04 [MW/m²]                 ejima → 0.4
Paux_tot → 100 [MW]                     P/R0 → 20.1 [MW/m]                      flattop → 0.7 [Hours]
Pα → 73.7 [MW]                          PB/R0 → 157 [MW T/m]
Pohm → 0.0939 [MW]                      PBp/R0 → 16.1 [MW T/m]
Pheat → 174 [MW]                        PBϵ/R0q95 → 7.31 [MW T/m]
Prad_tot → -42.9 [MW]                   neutrons_peak → 0.501 [MW/m²]

BOP                                     BUILD                                   COSTING
────────────────────────────────────    ────────────────────────────────────    ────────────────────────────────────
Pfusion → 368 [MW]                      PF_material → nb3sn                     capital_cost → 6.33 [$B]
Qfusion → 3.68                          TF_material → nb3sn_kdemo               levelized_CoE → Inf [$/kWh]
thermal_cycle_type → rankine            OH_material → nb3sn                     TF_of_total → 15.5 [%]
thermal_efficiency_plant → 22.7 [%]     TF_max_b → 15.4 [T]                     BOP_of_total → 2.04 [%]
thermal_efficiency_cycle → NaN [%]      OH_max_b → 19.9 [T]                     blanket_of_total → 20.3 [%]
power_electric_generated → 20.5 [MW]    TF_j_margin → 5.9                       cryostat_of_total → 3.13 [%]
Pelectric_net → -125 [MW]               OH_j_margin → 1.71
Qplant → 0.141                          TF_stress_margin → 3.15
TBR → 0.0655                            OH_stress_margin → 1.37

@ time = 0.0 [s]
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​
​

More things

disable or enable printing of what actor is executed

FUSE.actor_logging(dd, false) # or true
true