From c22a906817e7a7b95bb27bbb511ddfc779e71974 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Thu, 20 Jul 2023 15:59:45 -0700
Subject: [PATCH 01/13] Add NEO package to FUSE and integrate with
 ActorNeoclassical

---
 Project.toml                               |  1 +
 src/actors/transport/neoclassical_actor.jl | 36 ++++++++++++++++++++--
 2 files changed, 35 insertions(+), 2 deletions(-)

diff --git a/Project.toml b/Project.toml
index 0c531114a..8bd8059b9 100644
--- a/Project.toml
+++ b/Project.toml
@@ -39,6 +39,7 @@ Memoize = "c03570c3-d221-55d1-a50c-7939bbd78826"
 MeshTools = "804828d9-b4e9-4eca-b847-748d6c7bafa9"
 Metaheuristics = "bcdb8e00-2c21-11e9-3065-2b553b22f898"
 MillerExtendedHarmonic = "c82744c2-dc08-461a-8c37-87ab04d0f9b8"
+NEO = "081d0416-0f87-42c0-8034-5895805e76d8"
 NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56"
 NNeutronics = "a9424c20-d414-11ec-167b-9106c24d956c"
 NumericalIntegration = "e7bfaba1-d571-5449-8927-abc22e82249b"
diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 2bf9f1a32..181095b94 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -1,5 +1,6 @@
 import TGLFNN
 import TAUENN
+import NEO
 
 #= ===================== =#
 #  ActorNeoclassical      #
@@ -7,7 +8,7 @@ import TAUENN
 Base.@kwdef mutable struct FUSEparameters__ActorNeoclassical{T} <: ParametersActor where {T<:Real}
     _parent::WeakRef = WeakRef(nothing)
     _name::Symbol = :not_set
-    neoclassical_model::Switch{Symbol} = Switch{Symbol}([:changhinton], "-", "Neoclassical model to run"; default=:changhinton)
+    neoclassical_model::Switch{Symbol} = Switch{Symbol}([:changhinton, :neo], "-", "Neoclassical model to run"; default=:changhinton)
     rho_transport::Entry{AbstractVector{<:T}} = Entry{AbstractVector{<:T}}("-", "rho_tor_norm values to compute neoclassical fluxes on"; default=0.2:0.1:0.8)
 end
 
@@ -45,12 +46,43 @@ function _step(actor::ActorNeoclassical)
     model = resize!(dd.core_transport.model, :neoclassical)
     m1d = resize!(model.profiles_1d)
     m1d.grid_flux.rho_tor_norm = par.rho_transport
+    cp1d = dd.core_profiles.profiles_1d[]
 
     if par.neoclassical_model == :changhinton
         model.identifier.name = "Chang-Hinton"
         eqt = dd.equilibrium.time_slice[]
-        cp1d = dd.core_profiles.profiles_1d[]
         actor.flux_solutions = [TAUENN.neoclassical_changhinton(eqt, cp1d, rho, 1) for rho in par.rho_transport]
+    elseif par.neoclassical_model == :neo 
+        model.identifier.name = "NEO"
+        rho_cp = cp1d.grid.rho_tor_norm
+        gridpoint_cp = [argmin(abs.(rho_cp .- rho)) for rho in par.rho_transport]
+        
+        for i in gridpoint_cp
+            neo_solution = NEO.run_neo(NEO.InputNEO(dd, i))
+            total_ion_energy_flux = 0.0 
+           
+            for j in 1:11
+                if ismissing(getfield(neo_solution, Symbol("ENERGY_FLUX_$j")))
+                    energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(j-1)")) #subtract out electron energy flux, which is the last energy flux in the list
+                    total_ion_energy_flux -= energy_flux_electrons
+                    break
+                end
+            end
+
+            for field in fieldnames(NEO.flux_solution)
+                if ismissing(getfield(neo_solution, field))
+                    setfield!(neo_solution, field, 0.0)
+                end
+
+                if occursin("ENERGY", string(field))
+                    total_ion_energy_flux += getfield(neo_solution, field) 
+                end
+            end
+
+            solution = TGLFNN.flux_solution(0.0, 0.0, 0.0, total_ion_energy_flux)
+            push!(actor.flux_solutions, solution)
+        end
+  
     else
         error("$(par.neoclassical_model) is not implemented")
     end

From c5fa69eb7a24858f7ccfa9157dbc01383fc4efc3 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Thu, 10 Aug 2023 14:46:19 -0700
Subject: [PATCH 02/13] Include electron energy and particle fluxes in NEO
 solution

---
 src/actors/transport/neoclassical_actor.jl | 22 ++++++++++++++++++++--
 1 file changed, 20 insertions(+), 2 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 181095b94..63fd286dc 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -60,15 +60,24 @@ function _step(actor::ActorNeoclassical)
         for i in gridpoint_cp
             neo_solution = NEO.run_neo(NEO.InputNEO(dd, i))
             total_ion_energy_flux = 0.0 
+            particle_flux_electrons = 0.0
+            energy_flux_electrons = 0.0
            
             for j in 1:11
                 if ismissing(getfield(neo_solution, Symbol("ENERGY_FLUX_$j")))
-                    energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(j-1)")) #subtract out electron energy flux, which is the last energy flux in the list
+                    energy_flux_electrons += getfield(neo_solution, Symbol("ENERGY_FLUX_$(j-1)")) #subtract out electron energy flux, which is the last energy flux in the list
                     total_ion_energy_flux -= energy_flux_electrons
                     break
                 end
             end
 
+            for j in 1:11
+                if ismissing(getfield(neo_solution, Symbol("PARTICLE_FLUX_$j")))
+                    particle_flux_electrons += getfield(neo_solution, Symbol("PARTICLE_FLUX_$(j-1)")) #subtract out electron particle flux, which is the last particle flux in the list
+                    break
+                end
+            end
+
             for field in fieldnames(NEO.flux_solution)
                 if ismissing(getfield(neo_solution, field))
                     setfield!(neo_solution, field, 0.0)
@@ -79,7 +88,7 @@ function _step(actor::ActorNeoclassical)
                 end
             end
 
-            solution = TGLFNN.flux_solution(0.0, 0.0, 0.0, total_ion_energy_flux)
+            solution = TGLFNN.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
             push!(actor.flux_solutions, solution)
         end
   
@@ -96,6 +105,7 @@ end
 Writes ActorNeoclassical results to dd.core_transport
 """
 function _finalize(actor::ActorNeoclassical)
+    par = actor.par
     dd = actor.dd
     cp1d = dd.core_profiles.profiles_1d[]
     eqt = dd.equilibrium.time_slice[]
@@ -103,10 +113,18 @@ function _finalize(actor::ActorNeoclassical)
     model = findfirst(:neoclassical, actor.dd.core_transport.model)
     m1d = model.profiles_1d[]
     m1d.total_ion_energy.flux = zeros(length(actor.par.rho_transport))
+    m1d.electrons.particles.flux = zeros(length(actor.par.rho_transport))
+    m1d.electrons.energy.flux = zeros(length(actor.par.rho_transport))
+
     for (neoclassical_idx, rho) in enumerate(actor.par.rho_transport)
         rho_transp_idx = findfirst(i -> i == rho, m1d.grid_flux.rho_tor_norm)
         rho_cp_idx = argmin(abs.(cp1d.grid.rho_tor_norm .- rho))
         m1d.total_ion_energy.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].ENERGY_FLUX_i * IMAS.gyrobohm_energy_flux(cp1d, eqt)[rho_cp_idx] # W / m^2
+
+        if par.neoclassical_model == :neo
+            m1d.electrons.particles.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].PARTICLE_FLUX_e * IMAS.gyrobohm_particle_flux(cp1d, eqt)[rho_cp_idx]
+            m1d.electrons.energy.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].ENERGY_FLUX_e * IMAS.gyrobohm_energy_flux(cp1d, eqt)[rho_cp_idx]
+        end
     end
 
     return actor

From 514dad60bfccc5f8f11910143ec75f2ec3a166d1 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Mon, 21 Aug 2023 11:10:59 -0700
Subject: [PATCH 03/13] Include suggested changes to parsing neo_solution

Co-authored-by: Orso Meneghini <orso82@gmail.com>
---
 src/actors/transport/neoclassical_actor.jl | 8 +++++---
 1 file changed, 5 insertions(+), 3 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 63fd286dc..b59bd2f6f 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -63,17 +63,19 @@ function _step(actor::ActorNeoclassical)
             particle_flux_electrons = 0.0
             energy_flux_electrons = 0.0
            
+            #subtract from ions energy the electron energy flux, which is the last energy flux in the list (this will be added again later)
             for j in 1:11
                 if ismissing(getfield(neo_solution, Symbol("ENERGY_FLUX_$j")))
-                    energy_flux_electrons += getfield(neo_solution, Symbol("ENERGY_FLUX_$(j-1)")) #subtract out electron energy flux, which is the last energy flux in the list
-                    total_ion_energy_flux -= energy_flux_electrons
+                    energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(j-1)"))
+                    total_ion_energy_flux = -energy_flux_electrons
                     break
                 end
             end
 
+            #electron particle flux is the last particle flux in the list
             for j in 1:11
                 if ismissing(getfield(neo_solution, Symbol("PARTICLE_FLUX_$j")))
-                    particle_flux_electrons += getfield(neo_solution, Symbol("PARTICLE_FLUX_$(j-1)")) #subtract out electron particle flux, which is the last particle flux in the list
+                    particle_flux_electrons = getfield(neo_solution, Symbol("PARTICLE_FLUX_$(j-1)"))
                     break
                 end
             end

From fc3d1ddd12d54fdd50bb3c2a158c860bd704e489 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Mon, 21 Aug 2023 11:39:16 -0700
Subject: [PATCH 04/13] actor.par --> par

---
 src/actors/transport/neoclassical_actor.jl | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 63fd286dc..47dea4c62 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -112,11 +112,11 @@ function _finalize(actor::ActorNeoclassical)
 
     model = findfirst(:neoclassical, actor.dd.core_transport.model)
     m1d = model.profiles_1d[]
-    m1d.total_ion_energy.flux = zeros(length(actor.par.rho_transport))
-    m1d.electrons.particles.flux = zeros(length(actor.par.rho_transport))
-    m1d.electrons.energy.flux = zeros(length(actor.par.rho_transport))
+    m1d.total_ion_energy.flux = zeros(length(par.rho_transport))
+    m1d.electrons.particles.flux = zeros(length(par.rho_transport))
+    m1d.electrons.energy.flux = zeros(length(par.rho_transport))
 
-    for (neoclassical_idx, rho) in enumerate(actor.par.rho_transport)
+    for (neoclassical_idx, rho) in enumerate(par.rho_transport)
         rho_transp_idx = findfirst(i -> i == rho, m1d.grid_flux.rho_tor_norm)
         rho_cp_idx = argmin(abs.(cp1d.grid.rho_tor_norm .- rho))
         m1d.total_ion_energy.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].ENERGY_FLUX_i * IMAS.gyrobohm_energy_flux(cp1d, eqt)[rho_cp_idx] # W / m^2

From 016a2cb8928d9b35d513fb96fba83b757bb948b8 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Mon, 21 Aug 2023 11:53:55 -0700
Subject: [PATCH 05/13] Simpler way of parsing neo_solution outputs

---
 src/actors/transport/neoclassical_actor.jl | 24 +++++-----------------
 1 file changed, 5 insertions(+), 19 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index d4bcfef05..0cd30ce0c 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -56,29 +56,15 @@ function _step(actor::ActorNeoclassical)
         model.identifier.name = "NEO"
         rho_cp = cp1d.grid.rho_tor_norm
         gridpoint_cp = [argmin(abs.(rho_cp .- rho)) for rho in par.rho_transport]
+
+        n_species = length(cp1d.ion)+1
         
         for i in gridpoint_cp
             neo_solution = NEO.run_neo(NEO.InputNEO(dd, i))
-            total_ion_energy_flux = 0.0 
-            particle_flux_electrons = 0.0
-            energy_flux_electrons = 0.0
            
-            #subtract from ions energy the electron energy flux, which is the last energy flux in the list (this will be added again later)
-            for j in 1:11
-                if ismissing(getfield(neo_solution, Symbol("ENERGY_FLUX_$j")))
-                    energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(j-1)"))
-                    total_ion_energy_flux = -energy_flux_electrons
-                    break
-                end
-            end
-
-            #electron particle flux is the last particle flux in the list
-            for j in 1:11
-                if ismissing(getfield(neo_solution, Symbol("PARTICLE_FLUX_$j")))
-                    particle_flux_electrons = getfield(neo_solution, Symbol("PARTICLE_FLUX_$(j-1)"))
-                    break
-                end
-            end
+            energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(n_species)")) # electrons are always the last species in the list
+            total_ion_energy_flux = -energy_flux_electrons # exclude electron energy flux from total ion energy flux
+            particle_flux_electrons = getfield(neo_solution, Symbol("PARTICLE_FLUX_$(n_species)"))
 
             for field in fieldnames(NEO.flux_solution)
                 if ismissing(getfield(neo_solution, field))

From 8f444bf52514a1b86a0d64ea23a1a51aa5e57944 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Mon, 21 Aug 2023 11:58:48 -0700
Subject: [PATCH 06/13] Add NEO to Makefile

---
 Makefile | 5 ++++-
 1 file changed, 4 insertions(+), 1 deletion(-)

diff --git a/Makefile b/Makefile
index e6fcb844e..d04b98062 100644
--- a/Makefile
+++ b/Makefile
@@ -29,7 +29,7 @@ CURRENTDIR := $(shell (pwd -P))
 TODAY := $(shell date +'%Y-%m-%d')
 export JULIA_NUM_THREADS ?= $(shell julia -e "println(length(Sys.cpu_info()))")
 
-FUSE_PACKAGES_MAKEFILE := ADAS BoundaryPlasmaModels CHEASE CoordinateConventions EPEDNN FiniteElementHermite Fortran90Namelists FusionMaterials IMAS IMASDD MXHEquilibrium MeshTools MillerExtendedHarmonic NNeutronics QED SimulationParameters TAUENN TEQUILA TGLFNN VacuumFields 
+FUSE_PACKAGES_MAKEFILE := ADAS BoundaryPlasmaModels CHEASE CoordinateConventions EPEDNN FiniteElementHermite Fortran90Namelists FusionMaterials IMAS IMASDD MXHEquilibrium MeshTools MillerExtendedHarmonic NEO NNeutronics QED SimulationParameters TAUENN TEQUILA TGLFNN VacuumFields 
 FUSE_PACKAGES_MAKEFILE := $(sort $(FUSE_PACKAGES_MAKEFILE))
 FUSE_PACKAGES := $(shell echo '$(FUSE_PACKAGES_MAKEFILE)' | awk '{printf("[\"%s\"", $$1); for (i=2; i<=NF; i++) printf(", \"%s\"", $$i); print "]"}')
 DEV_PACKAGES := $(shell find ../*/.git/config -exec grep ProjectTorreyPines \{\} \; | cut -d'/' -f 2 | cut -d'.' -f 1 | tr '\n' ' ')
@@ -270,6 +270,9 @@ SimulationParameters:
 BoundaryPlasmaModels:
 	$(call clone_pull_repo,$@)
 
+NEO:
+	$(call clone_pull_repo,$@)
+
 # Install IJulia
 IJulia:
 	julia -e '\

From 8d9de003e8a86b4669ccdcd8713d795e9648225c Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Wed, 6 Sep 2023 11:45:11 -0700
Subject: [PATCH 07/13] Call changhinton from NEO instead of Tauenn

---
 src/actors/transport/neoclassical_actor.jl | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 1b565f6df..6a33cfa80 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -51,7 +51,7 @@ function _step(actor::ActorNeoclassical)
     if par.neoclassical_model == :changhinton
         model.identifier.name = "Chang-Hinton"
         eqt = dd.equilibrium.time_slice[]
-        actor.flux_solutions = [TAUENN.neoclassical_changhinton(eqt, cp1d, rho, 1) for rho in par.rho_transport]
+        actor.flux_solutions = [NEO.changhinton(eqt, cp1d, rho, 1) for rho in par.rho_transport]
     elseif par.neoclassical_model == :neo 
         model.identifier.name = "NEO"
         rho_cp = cp1d.grid.rho_tor_norm
@@ -76,7 +76,7 @@ function _step(actor::ActorNeoclassical)
                 end
             end
 
-            solution = TGLFNN.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
+            solution = IMAS.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
             push!(actor.flux_solutions, solution)
         end
   

From 71a4e5c19492d1f8d581f673499690238fceb78e Mon Sep 17 00:00:00 2001
From: TimSlendebroek <32385057+TimSlendebroek@users.noreply.github.com>
Date: Wed, 6 Sep 2023 12:11:04 -0700
Subject: [PATCH 08/13] universal flux_solution inside IMAS

---
 src/actors/transport/neoclassical_actor.jl | 8 +++-----
 src/actors/transport/tglf_actor.jl         | 4 ++--
 2 files changed, 5 insertions(+), 7 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 1b565f6df..4d57ea530 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -1,7 +1,5 @@
-import TGLFNN
 import TAUENN
 import NEO
-
 #= ===================== =#
 #  ActorNeoclassical      #
 #= ===================== =#
@@ -15,7 +13,7 @@ end
 mutable struct ActorNeoclassical{D,P} <: PlasmaAbstractActor
     dd::IMAS.dd{D}
     par::FUSEparameters__ActorNeoclassical{P}
-    flux_solutions::Vector{<:TGLFNN.flux_solution}
+    flux_solutions::Vector{<:IMAS.flux_solution}
 end
 
 """
@@ -32,7 +30,7 @@ end
 
 function ActorNeoclassical(dd::IMAS.dd, par::FUSEparameters__ActorNeoclassical; kw...)
     par = par(kw...)
-    return ActorNeoclassical(dd, par, TGLFNN.flux_solution[])
+    return ActorNeoclassical(dd, par, IMAS.flux_solution[])
 end
 
 """
@@ -76,7 +74,7 @@ function _step(actor::ActorNeoclassical)
                 end
             end
 
-            solution = TGLFNN.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
+            solution = IMAS.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
             push!(actor.flux_solutions, solution)
         end
   
diff --git a/src/actors/transport/tglf_actor.jl b/src/actors/transport/tglf_actor.jl
index daefc237b..53e582ff0 100644
--- a/src/actors/transport/tglf_actor.jl
+++ b/src/actors/transport/tglf_actor.jl
@@ -17,7 +17,7 @@ mutable struct ActorTGLF{D,P} <: PlasmaAbstractActor
     dd::IMAS.dd{D}
     par::FUSEparameters__ActorTGLF{P}
     input_tglfs::Vector{<:TGLFNN.InputTGLF}
-    flux_solutions::Vector{<:TGLFNN.flux_solution}
+    flux_solutions::Vector{<:IMAS.flux_solution}
 end
 
 """
@@ -35,7 +35,7 @@ end
 function ActorTGLF(dd::IMAS.dd, par::FUSEparameters__ActorTGLF; kw...)
     par = par(kw...)
     input_tglfs = Vector{TGLFNN.InputTGLF}(undef, length(par.rho_transport))
-    return ActorTGLF(dd, par, input_tglfs, TGLFNN.flux_solution[])
+    return ActorTGLF(dd, par, input_tglfs, IMAS.flux_solution[])
 end
 
 """

From 4cc2c34f07dc00aaaeb0f61704a82395255e86e0 Mon Sep 17 00:00:00 2001
From: TimSlendebroek <32385057+TimSlendebroek@users.noreply.github.com>
Date: Wed, 6 Sep 2023 13:11:16 -0700
Subject: [PATCH 09/13] change neoclassical_model to model

---
 src/actors/transport/neoclassical_actor.jl | 19 ++++++++++---------
 1 file changed, 10 insertions(+), 9 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 2b2d290af..3c30f3491 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -6,7 +6,7 @@ import NEO
 Base.@kwdef mutable struct FUSEparameters__ActorNeoclassical{T} <: ParametersActor where {T<:Real}
     _parent::WeakRef = WeakRef(nothing)
     _name::Symbol = :not_set
-    neoclassical_model::Switch{Symbol} = Switch{Symbol}([:changhinton, :neo], "-", "Neoclassical model to run"; default=:changhinton)
+    model::Switch{Symbol} = Switch{Symbol}([:changhinton, :neo], "-", "Neoclassical model to run"; default=:changhinton)
     rho_transport::Entry{AbstractVector{<:T}} = Entry{AbstractVector{<:T}}("-", "rho_tor_norm values to compute neoclassical fluxes on"; default=0.2:0.1:0.8)
 end
 
@@ -46,11 +46,11 @@ function _step(actor::ActorNeoclassical)
     m1d.grid_flux.rho_tor_norm = par.rho_transport
     cp1d = dd.core_profiles.profiles_1d[]
 
-    if par.neoclassical_model == :changhinton
+    if par.model == :changhinton
         model.identifier.name = "Chang-Hinton"
         eqt = dd.equilibrium.time_slice[]
         actor.flux_solutions = [NEO.changhinton(eqt, cp1d, rho, 1) for rho in par.rho_transport]
-    elseif par.neoclassical_model == :neo 
+    elseif par.model == :neo 
         model.identifier.name = "NEO"
         rho_cp = cp1d.grid.rho_tor_norm
         gridpoint_cp = [argmin(abs.(rho_cp .- rho)) for rho in par.rho_transport]
@@ -77,9 +77,6 @@ function _step(actor::ActorNeoclassical)
             solution = IMAS.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
             push!(actor.flux_solutions, solution)
         end
-  
-    else
-        error("$(par.neoclassical_model) is not implemented")
     end
 
     return actor
@@ -98,16 +95,20 @@ function _finalize(actor::ActorNeoclassical)
 
     model = findfirst(:neoclassical, actor.dd.core_transport.model)
     m1d = model.profiles_1d[]
+
+    if par.model == :neo
+        m1d.electrons.particles.flux = zeros(length(par.rho_transport))
+        m1d.electrons.energy.flux = zeros(length(par.rho_transport))
+    end
+
     m1d.total_ion_energy.flux = zeros(length(par.rho_transport))
-    m1d.electrons.particles.flux = zeros(length(par.rho_transport))
-    m1d.electrons.energy.flux = zeros(length(par.rho_transport))
 
     for (neoclassical_idx, rho) in enumerate(par.rho_transport)
         rho_transp_idx = findfirst(i -> i == rho, m1d.grid_flux.rho_tor_norm)
         rho_cp_idx = argmin(abs.(cp1d.grid.rho_tor_norm .- rho))
         m1d.total_ion_energy.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].ENERGY_FLUX_i * IMAS.gyrobohm_energy_flux(cp1d, eqt)[rho_cp_idx] # W / m^2
 
-        if par.neoclassical_model == :neo
+        if par.model == :neo
             m1d.electrons.particles.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].PARTICLE_FLUX_e * IMAS.gyrobohm_particle_flux(cp1d, eqt)[rho_cp_idx]
             m1d.electrons.energy.flux[rho_transp_idx] = actor.flux_solutions[neoclassical_idx].ENERGY_FLUX_e * IMAS.gyrobohm_energy_flux(cp1d, eqt)[rho_cp_idx]
         end

From 3d58ba07d8a957bcc6c21dd28258b30fc6b3f23b Mon Sep 17 00:00:00 2001
From: TimSlendebroek <32385057+TimSlendebroek@users.noreply.github.com>
Date: Wed, 6 Sep 2023 13:14:20 -0700
Subject: [PATCH 10/13] remove unces error handling

---
 src/actors/transport/flux_calculator_actor.jl | 4 ----
 1 file changed, 4 deletions(-)

diff --git a/src/actors/transport/flux_calculator_actor.jl b/src/actors/transport/flux_calculator_actor.jl
index 6b1911976..c0ed7ef1a 100644
--- a/src/actors/transport/flux_calculator_actor.jl
+++ b/src/actors/transport/flux_calculator_actor.jl
@@ -36,8 +36,6 @@ function ActorFluxCalculator(dd::IMAS.dd, par::FUSEparameters__ActorFluxCalculat
     elseif par.turbulence_model == :TGLF
         act.ActorTGLF.rho_transport = par.rho_transport
         turb_actor = ActorTGLF(dd, act.ActorTGLF)
-    else
-        error("turbulence_model `$(par.turbulence_model)` is not supported yet")
     end
 
     if par.neoclassical_model == :none
@@ -45,8 +43,6 @@ function ActorFluxCalculator(dd::IMAS.dd, par::FUSEparameters__ActorFluxCalculat
     elseif par.neoclassical_model == :neoclassical
         act.ActorNeoclassical.rho_transport = par.rho_transport
         neoc_actor = ActorNeoclassical(dd, act.ActorNeoclassical)
-    else
-        error("neoclassical_model `$(par.neoclassical_model)` is not supported yet")
     end
 
     return ActorFluxCalculator(dd, par, turb_actor, neoc_actor)

From 328a46a634ac7c642ed90cc11bd82106fa44cd46 Mon Sep 17 00:00:00 2001
From: TimSlendebroek <32385057+TimSlendebroek@users.noreply.github.com>
Date: Wed, 6 Sep 2023 14:15:44 -0700
Subject: [PATCH 11/13] store input.neo in the actor when full neo is ran

---
 src/actors/transport/neoclassical_actor.jl | 11 +++++++++--
 1 file changed, 9 insertions(+), 2 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index 3c30f3491..a02a43e0d 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -13,6 +13,7 @@ end
 mutable struct ActorNeoclassical{D,P} <: PlasmaAbstractActor
     dd::IMAS.dd{D}
     par::FUSEparameters__ActorNeoclassical{P}
+    input_neos::Union{Vector{<:NEO.InputNEO},Missing}
     flux_solutions::Vector{<:IMAS.flux_solution}
 end
 
@@ -30,7 +31,12 @@ end
 
 function ActorNeoclassical(dd::IMAS.dd, par::FUSEparameters__ActorNeoclassical; kw...)
     par = par(kw...)
-    return ActorNeoclassical(dd, par, IMAS.flux_solution[])
+    if par.model == :neo
+        input_neos = Vector{NEO.InputNEO}(undef, length(par.rho_transport))
+    else
+        input_neos = missing
+    end
+    return ActorNeoclassical(dd, par, input_neos, IMAS.flux_solution[])
 end
 
 """
@@ -57,7 +63,8 @@ function _step(actor::ActorNeoclassical)
 
         n_species = length(cp1d.ion)+1
         
-        for i in gridpoint_cp
+        for (idx,i) in enumerate(gridpoint_cp)
+            actor.input_neos[idx] = NEO.InputNEO(dd, i)
             neo_solution = NEO.run_neo(NEO.InputNEO(dd, i))
            
             energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(n_species)")) # electrons are always the last species in the list

From 4e480630a5d9548c637d7405a9aa7176445fad18 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Wed, 6 Sep 2023 16:06:19 -0700
Subject: [PATCH 12/13] Move conversion from NEO.flux_solution format to
 IMAS.flux_solution format out of ActorNeoclassical

This conversion is now done within the NEO.run_neo function and the outputs are returned in the form of the IMAS.flux_solution object
---
 src/actors/transport/neoclassical_actor.jl | 24 +++-------------------
 1 file changed, 3 insertions(+), 21 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index a02a43e0d..d963203ce 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -61,29 +61,11 @@ function _step(actor::ActorNeoclassical)
         rho_cp = cp1d.grid.rho_tor_norm
         gridpoint_cp = [argmin(abs.(rho_cp .- rho)) for rho in par.rho_transport]
 
-        n_species = length(cp1d.ion)+1
-        
         for (idx,i) in enumerate(gridpoint_cp)
-            actor.input_neos[idx] = NEO.InputNEO(dd, i)
-            neo_solution = NEO.run_neo(NEO.InputNEO(dd, i))
-           
-            energy_flux_electrons = getfield(neo_solution, Symbol("ENERGY_FLUX_$(n_species)")) # electrons are always the last species in the list
-            total_ion_energy_flux = -energy_flux_electrons # exclude electron energy flux from total ion energy flux
-            particle_flux_electrons = getfield(neo_solution, Symbol("PARTICLE_FLUX_$(n_species)"))
-
-            for field in fieldnames(NEO.flux_solution)
-                if ismissing(getfield(neo_solution, field))
-                    setfield!(neo_solution, field, 0.0)
-                end
-
-                if occursin("ENERGY", string(field))
-                    total_ion_energy_flux += getfield(neo_solution, field) 
-                end
-            end
-
-            solution = IMAS.flux_solution(particle_flux_electrons, 0.0, energy_flux_electrons, total_ion_energy_flux)
-            push!(actor.flux_solutions, solution)
+                actor.input_neos[idx] = NEO.InputNEO(dd, i)
         end
+
+        actor.flux_solutions = [NEO.run_neo(actor.input_neos[idx]) for idx in 1:length(par.rho_transport)]
     end
 
     return actor

From a5cef08dbe02f94904a21cb83c40390fbf889a82 Mon Sep 17 00:00:00 2001
From: adrianaghiozzi <67669644+adrianaghiozzi@users.noreply.github.com>
Date: Wed, 6 Sep 2023 16:16:54 -0700
Subject: [PATCH 13/13] Remove TAUENN + fix formatting

---
 src/actors/transport/neoclassical_actor.jl | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/src/actors/transport/neoclassical_actor.jl b/src/actors/transport/neoclassical_actor.jl
index d963203ce..1309c44ce 100644
--- a/src/actors/transport/neoclassical_actor.jl
+++ b/src/actors/transport/neoclassical_actor.jl
@@ -1,4 +1,3 @@
-import TAUENN
 import NEO
 #= ===================== =#
 #  ActorNeoclassical      #
@@ -62,7 +61,7 @@ function _step(actor::ActorNeoclassical)
         gridpoint_cp = [argmin(abs.(rho_cp .- rho)) for rho in par.rho_transport]
 
         for (idx,i) in enumerate(gridpoint_cp)
-                actor.input_neos[idx] = NEO.InputNEO(dd, i)
+            actor.input_neos[idx] = NEO.InputNEO(dd, i)
         end
 
         actor.flux_solutions = [NEO.run_neo(actor.input_neos[idx]) for idx in 1:length(par.rho_transport)]