-
Notifications
You must be signed in to change notification settings - Fork 0
Inversion
Important
The latest update has updated the command-line usage so that df only needs to be defined in the function file if using Newton-Raphson inversion or hybrid inversion. Otherwise, the inversion will run as expected without needing df.
Note
This code was integrated from another repository located here.
This code is for inversion steps on astrophysical functions (but can also be used generally).
Currently, the user must provide a single-variable inversion function f and its derivative function df.
The module path and filename are then passed as a command line argument and the inversion is done.
This module can also be used normally, by important functions in another python file and passing in values.
The code is original but based off of multiple pseudo-algorithms and mathematical formulae.
The hybrid inversion is a more robust way of doing inversion. The code is original but the idea comes from Matt Coleman, Research Scientist at Princeton University, who utilized the concept in the Athena++ astrophysical MHD codebase: https://github.com/PrincetonUniversity/athena.
The steps are as follows:
- Check relative error.
- If error > 0.01, use secant method for one step.
- Else use Newton-Raphson for one step.
- If step goes outside bracketing values, revert to Brent-Dekker (assuming bracketing values are given and exist).
- Update bracketing values. Loop untils error < tolerance or n iterations > n max iterations.
The inversion formulae used are (typically) from Wikipedia.
Newton-Raphson: https://en.wikipedia.org/wiki/Newton%27s_method.
Brent-Dekker: https://en.wikipedia.org/wiki/Brent%27s_method.
Secant: https://en.wikipedia.org/wiki/Secant_method.
Bisection: https://en.wikipedia.org/wiki/Bisection_method.
On general root-finding algorithms: https://en.wikipedia.org/wiki/Root-finding_algorithm.
Copyright © 2025 RandomKiddo
- Newton-Raphson:
newton_raphson(f, df, y0, tol, max_iter, verbose) - Brent-Dekker:
brent_dekker(f, a, b, tol, max_iter, verbose) - Secant:
secant(f, y0, y1, tol, max_iter, stopping_condition, verbose) - Bisection:
bisection(f, a, b, tol, max_iter, verbose) - Hybrid (as defined above):
hybrid_inversion(f, df, y0, a, b, tol, max_iter, verbose)
Common parameters:
f: The inversion callable function, taking a float and returning a float.
tol: The precision or tolerance for the relative error to use for the inversion, defaults to 1e-5.
max_iter: The max number of iterations to use in the inversion, defaults to 100.
verbose: If convergence warnings should be returned to the console, defaults to True.
Note
Convergence warnings means warnings that the function could not converge due to issues with precision, iterations, etc. Warnings will still be outputted when caused by swapped bracketed values, faulty initial values, and more.
Other parameters (depending on function):
df: The inversion callable derivative function, taking a float and returning a float.
y0: The float initial guess value to use. Pick a sensible one for your task.
y1: The float second initial guess value to use, where required. Pick a sensible one for your task.
a: The left bracketing value, where required. Pick a sensible one for your task.
b: The right bracketing value, where required. Pick a sensible one for your task.
stopping_condition: Int value representing which stopping condition to use (1, 2, or 3), defaults to 1.
Important
On the stopping condition (which is utilized for the secant method), we have three different possible conditions:
- Calculate error using
$|y_0-y_1| < \varepsilon$ . - Calculate error using
$|y_0/y_1 - 1| < \varepsilon$ . - Calculate error using
$|f(y_1)| < \varepsilon$ .
For a given tolerance
If any other integer is given, the default behavior is then 1.
Importing the subpackage can be done (if installed through PyPI or equiv.):
import astrocore.inversionor
from astrocore.inversion import *Using the command line leverages an ArgumentParser. The only required argument is the .py file with the required functions.
To then use in the command line, you can use any Python installation (Python, Conda, etc.) to run the inversion.py file. The below code line may need to be changed depending on the Python instance or the alias for Python (e.g., python3 instead of python).
python inversion.py "/path/to/inversion/function/file"The "inversion function file" referenced in the above path is a Python file with the inversion functions to be used by the inversion code. Below is a simple example:
def f(x: float) -> float:
return x**2 + 5
def df(x: float) -> float:
return 2*xUsing the command line python call as above would use default values for all values. Instead, we can use flags to specify values. To see these flags in the command line itself, you can use:
python inversion.py -hor
python inversion.py --helpHere are a list of the inversion flags you can use (their meanings are the same as the parameters for the inversion functions itself, see the arguments in this section for more info.):
-y0 Y0 Initial guess for inversion solver. Defaults to 1.0.
-tol TOL Convergence tolerance. Defaults to 1e-5.
-max_iter MAX_ITER Maximum number of iterations. Defaults to 100.
-verbose If convergence warnings should be returned to the console.
-y1 Y1 Initial second guess for inversion solver. Defaults to 0.0.
-a A The left bracketing value. Defaults to 0.0.
-b B The right bracketing value. Defaults to 1.0.
-stopcon STOPCON The int stopping condition number for secant inversion. Defaults to 1.
Note
For those who have never used command line flags, the capital letters following the flag means that it is asking for a value after a flag. So for a flag like -verbose, it requires no values after it, so you can use it as follows (or simply don't include it to exclude its behavior):
python inversion.py -verboseWhereas, flags like -y0 expect a value (in this case a float) Y0 to override the default value of 1.0 (again, in this case). So one uses it as follows (or excludes it if it's not needed or the default value is ok):
python inversion.py -y0 5.0After inputting all the necessary flags (if you know which method you want beforehand, then you can just include the ones you need instead of thinking about all the values), an example command line command will look like:
python inversion.py test_inversion_funcs.py -a -4 -b 1.3333After running this, the module will prompt which inversion method you like. It will use default arguments if you pick a method that requires values not provided in the shell command. You must enter the integer next to the inversion method you wish to use, and the prompt will continuously repeat if faulty input values are given. Then the inversion will occur, and the root, error, and iterations used will be returned.
Tip
Command line usage is generally best for using static functions (meaning you aren't trying many different f and df), where the root is the only important step. The values cannot be used in Python after running in the command line. If you wish to do analysis on the results in Python itself (like using the root to then estimate some other value), you should use this code in module form (below).
Using this code like a module means utilizing the source code in Python itself, after importing the module.
You can define functions and import inversion methods to then use in Python directly.
from astrocore.inversion import brent_dekker
def f(x: float) -> float:
return (x+3)*(x-1)**2
brent_dekker(f, -4, 4/3)Warning
Creating f and df for use in inversion is assumed to only take one positional argument, and return one value (a float of the value of the function). When using the command line, the method signature is validated for both f and df so no problems will occur. If, however, you are using this code as a module, the functions are not checked for validity, which could cause issues when doing inversion. Create functions that only take one positional argument, which should be a float like the following:
def df(x: float) -> float:
return 2*xThe actual parameter name does not matter. The functions must be called f and df for command-line usage. For module usage, the function name does not matter.
This page was last edited on 09.16.2025