Over-prediction of fast He ion pressure

[daveweisberg]

I am concerned that FUSE is over-predicting the fast He ion pressure, resulting in overly conservative FPP solutions (i.e. fusion power is too low relative to plasma beta). Nan Shi has re-analyzed our 8MA STEP solution and sees that / = 1.07. Compare this to the many FPP FUSE optimization scans I've run over the past week, all of which have a betaN / betaN_th ratio greater than this:

Notice that most FUSE solutions cluster around the betaN / betaN_th = 1.33 line, well above the STEP result of 1.07.

This is a critical issue, because fusion power scales like the the square of P_thermal. Since we are optimizing our FPP solutions to minimize the total betaN, increasing the betaN / betaN_th ratio will dramatically decrease the fusion produced at a given betaN. For two solutions that have equal betaN but different betaN / betaN_th ratios of 1.33 and 1.07, the 1.33 case will have (1.07/1.33)^2 = 65% of the fusion power of the 1.07 case!

I am asking the FUSE team to do two things:

1. Look at Nan's STEP solution (located on saturn here: /fusion/ga/projects/ird/ptp/shinan/omfit/toCol/toDave/8MA_fastion_1160.json) and run it in FUSE to see what the predicted fast ion pressure is.
2. Determine whether the FUSE vs. STEP fast ion pressure calculations are different enough to warrant a change to the FUSE model, and understand what needs to be done to bring FUSE into agreement with STEP.

[daveweisberg]

Update: it appears this is partially a collisionality effect. Plotting vs. ne0 / Te0^1.5 shows a distinct inverse relationship, with additional dependency on major radius:

[orso82]

@daveweisberg can you please try to do one more run and set approximate=true here: `https://github.com/ProjectTorreyPines/IMAS.jl/blob/d4fed452ed2372df7b949d785e3385a62d151167/src/physics/fast.jl#L144

[daveweisberg]

@orso82 Ok I finished an optimization run with approximate=true, and I don't see any real difference from the corresponding approximate=false run.

Here are the the Pareto-optimal solutions:

And here are all the (filtered) optimization solutions:

[daveweisberg]

@orso82 Any progress on comparing FUSE fast ion calculations on Nan's ODS file to STEP results?

[orso82]

@daveweisberg I can confirm that FUSE is calculating a fast ion pressure that is higher than the one calculated in STEP. This is how I have checked it:

dd = IMAS.json2imas(dirname(pathof(FUSE)) * "/../sample/highbatap_fpp_8MA_adhoc_EC.json")

p = plot(dd.core_profiles.profiles_1d[], :pressure, label="Pressure ODS")
plot!(dd.core_profiles.profiles_1d[], :pressure_thermal, label="Pressure_thermal ODS")

IMAS.fusion_source!(dd.core_sources, dd.core_profiles; only_DT=true)
IMAS.fast_particles!(dd.core_sources, dd.core_profiles.profiles_1d[]; verbose=true)

plot!(dd.core_profiles.profiles_1d[], :pressure, label="Pressure FUSE", ls=:dash, lw=3)
plot!(dd.core_profiles.profiles_1d[], :pressure_thermal, label="Pressure_thermal FUSE", ls=:dash, lw=3)

display(p)

There seems to be a missing a factor of 1/3 here: https://github.com/ProjectTorreyPines/IMAS.jl/blob/d59319c31805cd87664b5b3fb920f38ac7b0249e/src/physics/fast.jl#L285

with pressa = taus .* 2.0 ./ 3.0 .* (sion.particles .* particle_energy .* constants.e) ./ 3.0 I get near perfect agreement:

Can @jmcclena and @lstagner please take a look and see why we're missing such factor?

[daveweisberg]

@orso82 Great detective work! Let's get this fix merged ASAP :)

[bclyons12]

I would urge figuring out where this factor of three comes from before pushing any changes. Are we certain STEP isn't wrong?

[bclyons12]

@orso82 unrelated, but this line needs an @views https://github.com/ProjectTorreyPines/IMAS.jl/blob/d59319c31805cd87664b5b3fb920f38ac7b0249e/src/physics/fast.jl#L282

[jmcclena]

From freya, it looks like the equation should be:

p_f = 1/2 * taus * g_e* q_f * 2/3 e

where g_e is the fraction of heating to the electrons. I haven't figured it comes from yet. Looking at what Jeff started to put in TGYRO makes much more sense to me:

 ! JC 2022:
     ! Adding alpha-particle density/temperature calculation [Estrada 2006]
     ! Eqs (9),(10), where a=v_cross/v_alpha
     a = 1/sqrt(x_a)
     i2 = (1/3.0)*log((1+a**3)/a**3)
     i4 = 0.5-a**2*((1/6.0)*log((1-a+a**2)/(1+a)**2)+ &
          1/sqrt(3.0)*(atan((2-a)/(a*sqrt(3.0)))+pi/6))

Then p_f can be calculated with p_f = q_f * taus* i4 * 2/3 e. So it's like freya is assuming that i4/2 =g_e

[daveweisberg]

@jmcclena Does this reasoning agree with Orso's observation that FUSE over-predicts the fast ion pressure by a factor of 3?

[orso82]

Thank you for tracking this down @jmcclena !

One minor correction, rewriting the equations so they have the factors in the same order

p_f = q_f * taus * 2/3 * e * 1/2 * g_e       # [FREYA]
p_f = q_f * taus * 2/3 * e * i4              # [Estrada 2006]

it looks like freya is assuming that g_e/2 = i4

Does this factor have no radial dependence, and comes out to be ~1/3? is it 1/3 only for the case that we are looking at or that's just ~1/3 constant?

What do you think is the best way forward?

• Should we decide to just trust something that we can track the source of? namely Estrada 2006 ?
• Should we look into what TRANSP calculates for some benchmark case?

[orso82]

Reporting here an email by @bgrierson21

---

TRANSP run 20100W53. Look at time 162 +/- 1 second before the density ramp. Most of the species are equilibrated at this time. The multi-plot “pdens” shows all species, and lists the variables in the terminal.

[bclyons12]

A chat with Mike Van Zeeland led me to Waltz 2014, which again refers to Estrada 2006, so I think we're good to use that. The factor I_4 that is currently missing in FUSE varies with the electron temperature, but is typically in the range of 0.3-0.4. Here's a table, as well as the formula defining the fast-ion pressure from a slowing-down distribution (Eq. 14)

where those quantities are defined by

[bclyons12]

Appears to have been fixed by https://github.com/ProjectTorreyPines/IMAS.jl/pull/83. I'll let @daveweisberg close this when he's satisfied.

[orso82]

I think we can close this. @daveweisberg please re-open if that's not the case