Repository navigation
Derive plasma beta explicitly from the current plasma state - #4640
grmtrkngtn wants to merge 31 commits into
Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #4640 +/- ##
==========================================
+ Coverage 49.92% 49.98% +0.06%
==========================================
Files 151 151
Lines 30236 30451 +215
==========================================
+ Hits 15094 15220 +126
- Misses 15142 15231 +89 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
|
Hmm.. this only happens when beta is used to calculate things that relate to how beta is calculated. As far as I am aware, this was only happening in the Saamuli pedestal temperature scaling Where is the fixed point problem now for the open-source PROCESS LAR point? I just think it's really important to identify what these "circular loops" are. More importantly, are you sure this is a complete closure of the plasma calculation chain, or are there other fixed point issues? |
|
Going to invite comments from as many @ukaea/process-model-review on this one, thanks @CoronelBuendia for already commenting. @grmtrkngtn I will do a code review once modeller questions/concerns have been answered/addressed. |
Yes, that is the primary motivation.
I'm not sure what you mean, it doesn't need another PR - it's independent of anything else. |
As in, why does this PR not remove the beta consistency equation from the code and input files? |
| 2: IterationVariable("b_plasma_toroidal_on_axis", "physics", 0.010, 30.00), | ||
| 3: IterationVariable("rmajor", "physics", 0.1, 50.00), | ||
| 4: IterationVariable("temp_plasma_electron_vol_avg_kev", "physics", 5.0, 150.0), | ||
| 5: IterationVariable("beta_total_vol_avg", "physics", 0.001, 1.0), |
There was a problem hiding this comment.
Why are we removing it as an iteration variable in this case? How would a user constraint the beta value then. We can still have this as an iteration but now without the constraint
There was a problem hiding this comment.
Now that beta is calculated directly from the plasma state, I don't think it makes sense for it to remain an iteration variable.
If a user wants to restrict beta, that should be done by applying a constraint to the derived beta rather than allowing the optimiser to vary beta independently of the plasma state. We already have beta upper and lower-limit constraints for this purpose.
There was a problem hiding this comment.
I would agree. Now that beta_total_vol_avg is calculated in the code (here) having it as an iteration variable would be incorrect at best and dangerous at worse.
The solver would set beta_total_vol_avg = x and then the code would update it to beta_total_vol_avg. = y and all of the gradients would be calculated wrt y but the solver would think they related to x, which would cause the solver to move in incorrect ways.
| constraint_registration, | ||
| ) | ||
|
|
||
|
|
There was a problem hiding this comment.
Removing this breaks every previous PROCESS IN.DAT. If we want to progress with this change, could we have a grace period where both methods of determining beta are valid.
There was a problem hiding this comment.
I'm indifferent to this. Happy to allow both approaches although I'd argue it's better to rip the band-aid off in on go.
It is straight forward to modify any IN.DAT by removing icc = 1, ixc = 5 and reducing n_equality_constraints by 1.
grmtrkngtn
changed the title
|
Prompted by Matti's comment - I've made a significant change to the proposed solution in this PR although it achieves the same effect. Updated description. |
|
|
||
| """ | ||
| # Thermal beta is derived directly from the current plasma state. | ||
| self.data.physics.beta_thermal_vol_avg = self.calculate_plasma_beta( |
There was a problem hiding this comment.
Just wondering if we could do this all tidy and just put this in the run method of the beta class?
There was a problem hiding this comment.
Yes, I've tidied this up in the latest commits
There was a problem hiding this comment.
I cant see the thermal calc inside the beat run method its still showing as a random calc for me
I've added this into the latest version (as in removed beta consistency from code and inputs). Some may prefer a grace period although I think it's better to just remove it if it's not needed any more. |
|
So at present, @grmtrkngtn, none of the regression tests actually have a fixed point issue in beta, meaning it is not possible to verify that removing beta as an iteration variable is successfully addressed by fixed point iteration (a.k.a. idempotence loop in PROCESS). I don't actually doubt that this is the case. Just that in future (as we have had in the past e.g. for EU-DEMO with the Saamuli T_e_ped scaling) we may need beta as a variable in the calculation chain of beta; i.e. beta = f(beta). Note that this will (or already does happen with other variables). For example the Giacomin pedestal density is f(P_sol), so if this is calculated before P_sol, we would have P_sol = g(P_sol). Again, this should be treated by the existing Picard iteration (idempotence loop), but if these variables ever grow in number (and I think it's fair to say it is not known how many exist at present, see e.g. #4638), it's possible Picard iteration will not cut it. |
…fter fast particle calculation
This reverts commit 8dc8e88.
|
After some discussion we have agreed to leave in icc = 1 tool for checking potential fpp issues which may arise in future. IXC = 5 will be removed. |
| ... | ||
| ``` | ||
| The beta and plasma power balance equality constraints are recommended to ensure model consistency, and the `temp_plasma_electron_vol_avg_kev` and `nd_plasma_electrons_vol_avg` iteration variables used to solve them. | ||
| The plasma power balance equality constraint is recommended to ensure model consistency, and the `temp_plasma_electron_vol_avg_kev` and `nd_plasma_electrons_vol_avg` may be used as iteration variables used to solve them. |
There was a problem hiding this comment.
Now that there is only one consistency equation enabled only one active iteration variable could be used
There was a problem hiding this comment.
Yes - I've changed this to tell the user to select one of these to solve. This should be a temporary change while we add the fuel-ion equilibrium into code to solve with temperature and density.
There was a problem hiding this comment.
This file needs an iteration variable removing since it now only has one consistency constraint
There was a problem hiding this comment.
I've removed density. This should be a temporary change pending #4639 which will give us two equations for which we can solve with te and ne.
| ) | ||
|
|
||
| # Total beta is derived from the thermal and fast-particle components. | ||
| self.data.physics.beta_total_vol_avg = ( |
There was a problem hiding this comment.
Now that beta_total_vol_avg is calculated we should also remove it as an input
There was a problem hiding this comment.
Yes, done - also removed from regression input files accordingly.
There was a problem hiding this comment.
Again will need an iteration variable removed because it now only has 2 consistency constraints
There was a problem hiding this comment.
I've removed density. This should be a temporary change pending #4639 which will give us two equations for which we can solve with te and ne.
Co-authored-by: Timothy <75321887+timothy-nunn@users.noreply.github.com>
… equality contraint
|
|
||
| """ | ||
| # Thermal beta is derived directly from the current plasma state. | ||
| self.data.physics.beta_thermal_vol_avg = self.calculate_plasma_beta( |
There was a problem hiding this comment.
I cant see the thermal calc inside the beat run method its still showing as a random calc for me
Co-authored-by: Christopher Ashe <91618944+chris-ashe@users.noreply.github.com>
ab551cb to
4508b31
Compare
6 participants
Description
This PR changes plasma beta from an independently optimised quantity to a quantity derived explicitly from the current plasma state. This has been motivated by UQ work, where we want to define solutions to PROCESS in a determined system of equations.
Previously,
beta_total_vol_avgcould be included as iteration variableixc=5and constrained using the beta consistency equationicc=1. Investigation of the non-idempotent beta behaviour showed that the current model no longer contains a meaningfulbeta = f(beta)feedback loop.Instead, the remaining non-idempotence came from calculation ordering:
beta_total_vol_avgcould be evaluated beforebeta_fast_alphahad been updated for the current model evaluation.This PR therefore:
calculates the fast-particle beta contributions before finalising plasma beta
derives thermal beta directly from the current plasma pressure and magnetic field
derives total beta as the sum of thermal, fast-alpha and beam contributions
removes beta_total_vol_avg as iteration variable ixc=5
removes the beta consistency constraint icc=1
Investigation
Thanks to Matti for questioning whether this was still genuinely a fixed-point problem - it appears that it is not.
To test this, I took the converged DEMO LAR point on main, held all other optimisation variables fixed, and varied only the supplied
beta_total_vol_avgChanging the supplied beta therefore does not change the upstream plasma quantities from which beta is calculated. This shows that the current model no longer contains a meaningful beta -> plasma state -> beta fixed-point loop.
New solution
There is stil some kind of non-idempotent behaviour that Jon picked up on while developing UQ.
Inspection of the model ordering showed that
PlasmaBeta.run()could be called beforebeta_fast_alphawas recalculated for the current model evaluation.Total beta could therefore use the fast-alpha contribution left from the previous evaluation, with a subsequent evaluation then using the updated value.
The issue is therefore calculation ordering/stale state rather than a physical beta fixed point.
New calculation
Thermal beta is now calculated directly from the current thermal plasma pressure and the total magnetic field
beta_thermal = 2 μ0 p_thermal / B²Total beta is then:
beta_total = beta_thermal + beta_fast_alpha + beta_beamThis makes
beta_total_vol_avga derived state quantity rather than an independent optimisation variable.Independent consistency check
Before removing
icc=1, I tested both regression cases with:ixc=5 OFFicc=1 ONThis leaves the old consistency equation in place as an independent check while preventing the optimiser from varying beta to satisfy it.
For both cases the beta consistency residual was exactly zero:
DEMO LAR: beta consistency residual = 0.0ST: beta consistency residual = 0.0Adding or removing
icc=1did not otherwise change either converged solution or the number of VMCON iterations.This confirms that once beta is derived from the current plasma state, the old consistency equation is redundant.
Regression results
DEMO LAR
The revised formulation reproduces the previous LAR solution very closely.
mainThe physical and engineering design is virtually unchanged. The main difference is the optimisation trajectory: removing the redundant beta variable/equality increases the number of VMCON iterations for this case.
ST regression
The ST case also recovers effectively the same plasma optimum:
mainThe plasma temperature, density, current and
q95are effectively equal..The ST optimiser does settle on a different feasible engineering layout, including a change in net-electric output from approximately
106 MWto170 MW.This does not result from a change in the plasma optimum: the optimisation objective is
Q.This is not caused by a change in the plasma optimum: the optimisation objective is Q, which is effectively identical between the two runs. Removing the redundant beta optimisation degree of freedom changes the optimisation path and leads to a different feasible engineering point with the same plasma optimum.
OUT.DATs for inspection
beta_fix_low_aspect_ratio_DEMO.OUT.DAT.txt
icc1_beta_fix_low_aspect_ratio_DEMO.OUT.DAT.txt
beta_fix_st_regression.OUT.DAT.txt
icc1_beta_fix_st_regression.OUT.DAT.txt