8000 736 warnings for function ncore by chris-ashe · Pull Request #3060 · ukaea/PROCESS · GitHub
[go: up one dir, main page]
More Web Proxy on the site http://driver.im/
Skip to content

736 warnings for function ncore #3060

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 10 commits into from
Feb 20, 2024
Merged

736 warnings for function ncore #3060

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
167b650
ncore warning
chris-ashe Feb 14, 2024
355e43c
Merge branch 'main' into 736-warnings-for-function-ncore
chris-ashe Feb 14, 2024
3b33e4d
syntax fix
chris-ashe Feb 14, 2024
2ba076e
return ncore
chris-ashe Feb 14, 2024
bf25b7c
change error type
chris-ashe Feb 14, 2024
52fd0de
tidy up
chris-ashe Feb 15, 2024
cc821c0
update docs
chris-ashe Feb 15, 2024
9f319f4
updated bounds of pedestal params
chris-ashe Feb 20, 2024
1567752
add ncore = 1E-6 back in
chris-ashe Feb 20, 2024
229fdda
re-update pedestal bounds
chris-ashe Feb 20, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions documentation/proc-pages/physics-models/plasma.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -267,6 +267,11 @@ simulations.
If `ipedestal` = 1 or 2 then the pedestal density `neped` is set as a fraction `fgwped` of the
Greenwald density (providing `fgwped` >= 0). The default value of `fgwped` is 0.8[^7].

!!! warning " Un-realistic profiles"

If `ipedestal >= 1` it is highly recommended to use constraint equation 81 (icc=81). This enforces solutions in which $n_0$ has to be greater than $n_{ped}$.
Negative $n_0$ values can also arise during iteration, so it is important to be weary on how low the lower bound for $n_e (\mathtt{dene})$ is set.

## Beta Limit

The plasma beta limit[^7] is given by
Expand Down
24 changes: 9 additions & 15 deletions process/profiles.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,10 +3,7 @@
from scipy import integrate
from abc import ABC, abstractmethod

from process.fortran import (
maths_library,
physics_variables,
)
from process.fortran import maths_library, physics_variables, error_handling

logger = logging.getLogger(__name__)
# Logging handler for console output
Expand Down Expand Up @@ -127,7 +124,9 @@ def calculate_profile_y(self, rho, rhopedn, n0, nped, nsep, alphan):
)

@staticmethod
def ncore(rhopedn, nped, nsep, nav, alphan):
def ncore(
rhopedn: float, nped: float, nsep: float, nav: float, alphan: float
) -> float:
"""This routine calculates the core denesity of a pedestalised profile.

:param rhopedn: normalised minor radius pedestal position
Expand All @@ -138,11 +137,10 @@ def ncore(rhopedn, nped, nsep, nav, alphan):
:type nsep: float
:param nav: electron density (/m3)
:type nav: float
:param alphan: _description_
:param alphan: density peaking parameter
:type alphan: float
:return: Core density
:rtype: numpy.array
:type: float
"""

ncore = (
Expand All @@ -155,15 +153,11 @@ def ncore(rhopedn, nped, nsep, nav, alphan):
)
)

if ncore < 0:
# Prevent ncore from going negative (and terminating the optimisation) by
# kludging to small positive value. Allows solver to continue and
# hopefully be constrained away from this point (e.g. constraint 81, ne0 > neped)
logger.exception("ncore is negative. Kludging to 1e-6.")
if ncore < 0.0:
# Allows solver to continue and
# warns the user to raise the lower bound on dene if the run did not converge
error_handling.report_error(282)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What was the reason this was not done before? What was the point of the kludge if the solver can deal with a negative ncore? I always thought this was to stop NaNs further down the code

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Running the files it dosent seem to make a difference as they all fail eventually anyway. It seems better to have it fail quickly and the user updates the lower bound to prevent it going negative again in future runs

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I probably need @jonmaddock to comment on this, as he is the one that wrote the kludge. For now, since this is a level 2 error, I will approve it but will make a follow-up issue to check that this fix is ok with the UQ tools (which is why the kludge was added)

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If needs be we could put the kludge back in without adding back the logger output as it just slows down the terminal output

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Probably a good idea, and then #3059 can address making this model more robust

ncore = 1.0e-6
# raise ValueError(
# f"Core density is negative: {ncore}. {nped = }, {nsep = }, {nav = }"
# )
return ncore

def set_physics_variables(self):
Expand Down
7 changes: 6 additions & 1 deletion process/utilities/errorlist.json
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@
"comment2": [
"Increment n_errortypes if an error is added to this list"
],
"n_errortypes": 281,
"n_errortypes": 282,
"errors": [
{
"no": 1,
Expand Down Expand Up @@ -1414,6 +1414,11 @@
"no": 281,
"level": 3,
"message": "CHECK: Cannot have i_tf_bucking >= 2 when tf_in_cs = 1"
},
{
"no": 282,
"level": 2,
"message": "CHECK: ncore is going negative when solving. Please raise the value of dene and or its lower limit."
}
]
}
8 changes: 4 additions & 4 deletions source/fortran/iteration_variables.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,7 @@ subroutine init_itv_6
use numerics, only: lablxc, boundl, boundu
implicit none
lablxc(6) = 'dene '
boundl(6) = 1.00D19
boundl(6) = 2.00D19
boundu(6) = 1.00D21
end subroutine init_itv_6

Expand Down Expand Up @@ -3244,8 +3244,8 @@ subroutine init_itv_145
use numerics, only: lablxc, boundl, boundu
implicit none
lablxc(145) = 'fgwped '
boundl(145) = 0.500D0
boundu(145) = 1.000D0
boundl(145) = 0.100D0
boundu(145) = 0.9D0
end subroutine init_itv_145

real(kind(1.d0)) function itv_145()
Expand Down Expand Up @@ -3377,7 +3377,7 @@ subroutine init_itv_152
implicit none
lablxc(152) = 'fgwsep '
boundl(152) = 0.001D0
boundu(152) = 1.000D0
boundu(152) = 0.5D0
end subroutine init_itv_152
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am a little concerned about changing the default bounds as this could lead to existing input files (that use the default bounds) starting outside of the feasible region. In this instance, the 3 changes are minor so I think I am happy with it.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The output values in the large tokamak regression are 0.85 and 0.5 for pedestal and separatrix respectively. @mkovari
Which is outside the bounds recommended

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Output values? Have you calculated neped/nGW and nesep/nGW from the data in the output file?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Electron density at separatrix / nGW (fgwsep_out) 5.000E-01
Electron density at pedestal / nGW (fgwped_out) 8.500E-01

From the output file of the large tokamak regression test

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would still like the warning to be explicit about what sort of lower and upper bounds people should set.

For example:
dene Recommended value of boundl(6) > 80% of expected volume-averaged density.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this true for all cases? Even for a stellerator?
The output warning are limited to one line so I will need to amend the current one as I cant add it under unless we change their length

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe point people to a section in the documentation with the warning? And have a longer explanation in the docs.


real(kind(1.d0)) function itv_152()
Expand Down
0