diff --git a/examples/example_3dq8BMS.ipynb b/examples/example_3dq8BMS.ipynb new file mode 100644 index 0000000..5caf359 --- /dev/null +++ b/examples/example_3dq8BMS.ipynb @@ -0,0 +1,263 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Example usage of NRSur3dq8BMSRemnant fit." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import surfinBH" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "fit_name = 'NRSur3dq8BMSRemnant'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load the fit, this only needs to be done once at the start of a script" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Loaded NRSur3dq8BMSRemnant fit.\n" + ] + } + ], + "source": [ + "fit = surfinBH.LoadFits(fit_name)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Read the documentation" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Help on Fit3dq8BMS in module surfinBH._fit_evaluators.fit_3dq8BMS object:\n", + "\n", + "class Fit3dq8BMS(surfinBH.surfinBH.SurFinBH)\n", + " | Fit3dq8BMS(name)\n", + " |\n", + " | A class for the NRSur3dq8BMSRemnant model presented in Da Re et al.,\n", + " | arxiv:????.?????. This model is referred to as surfinBH3dq8BMS in the paper.\n", + " |\n", + " | This model predicts the supertransation modes up to ell = 8 and the\n", + " | 3 components of the boost velocity of the BMS transformation from the\n", + " | inspiral (PN) BMS frame to the remnant black hole BMS frame. The boost\n", + " | velocity coincides with the remnant black hole kick velocity as observed\n", + " | from the inspiral (PN) BMS frame.\n", + " |\n", + " | The model was trained on nonprecessing binary black hole systems. The fits\n", + " | are done using Gaussian Process Regression (GPR) and also provide an error\n", + " | estimate along with the fit value.\n", + " |\n", + " | This model has been trained in the parameter space:\n", + " | q <= 8, |chiAz| <= 0.8, |chiBz| <= 0.8\n", + " |\n", + " | =========================================================================\n", + " | Usage:\n", + " |\n", + " | import surfinBH\n", + " |\n", + " | # Load the fit\n", + " | fit = surfinBH.LoadFits('NRSur3dq8BMSRemnant')\n", + " |\n", + " | #To evaluate the supertranslation parameter \"alpha\" and the boost velocity\n", + " | \"boost\" together with their respective 1-sigma error estimates call:\n", + " | alpha, boost, alpha_err, boost_err = fit.all(q, chiA, chiB, **kwargs)\n", + " |\n", + " | # NOTE: The ell=0 mode, corresponding to a time translation,\n", + " | is set to be identically zero.\n", + " |\n", + " | # alpha is expressed as a complex array of spherical harmonics modes in the\n", + " | order\n", + " | (0,0),(1,-1),(1,0),(1,1),(2,-2),(2,-1),(2,0),...\n", + " | (same convention as the spherical_functions and scri packages)\n", + " |\n", + " | The arguments for each of these call methods are as follows:\n", + " | Arguments:\n", + " | q: Mass ratio (q>=1)\n", + " |\n", + " | chiA: Dimensionless spin of the larger BH (array of size 3).\n", + " |\n", + " | chiB: Dimensionless spin of the smaller BH (array of size 3).\n", + " | This model allows only nonprecessing spins, so only the\n", + " | z-components of these arrays should be non-zero.\n", + " |\n", + " | Optional arguments:\n", + " | allow_extrap:\n", + " | If False, raises a warning when q > 8.1 or |chiA|,|chiB| > 0.81,\n", + " | and raises an error when q > 10.1 or |chiA|,|chiB| > 1.\n", + " | If True, allows extrapolation to any q and |chiA|,|chiB| <= 1.\n", + " | Use at your own risk.\n", + " | Default: False.\n", + " |\n", + " | The inspiral frame can be defined up to an arbitrary U(1) rotation about\n", + " | the z-axis. We fix this gauge freedom using the same frame alignment choice\n", + " | used in NrSur3dq8Remnant. The inspiral frame at -100M from the peak of the\n", + " | waveform is defined as:\n", + " | The z-axis is along the orbital angular momentum direction of the binary.\n", + " | The x-axis is along the line of separation from the smaller BH to\n", + " | the larger BH at this time.\n", + " | The y-axis completes the triad.\n", + " | The rotation leading to this frame choice is included in the BMS transformation\n", + " | as outlined in arxiv:????.?????.\n", + " |\n", + " | Method resolution order:\n", + " | Fit3dq8BMS\n", + " | surfinBH.surfinBH.SurFinBH\n", + " | builtins.object\n", + " |\n", + " | Methods defined here:\n", + " |\n", + " | __init__(self, name)\n", + " | name: Name of the fit excluding the surfinBH prefix. Ex: 7dq2.\n", + " | soft_param_lims: param limits beyond which to raise a warning.\n", + " | hard_param_lims: param limits beyond which to raise an error.\n", + " | aligned_spin_only: raise an error if given precessing spins.\n", + " | See _fit_evaluators.fit_7dq2.py for an example.\n", + " |\n", + " | all(self, *args, **kwargs)\n", + " | Evaluates fit and 1-sigma error estimate for supertranslation parameter\n", + " | alpha, and boost velocity of the BMS transformation from the inspiral BMS\n", + " | frame to the remnant BMS frame.\n", + " | Returns:\n", + " | alpha, boost, alpha_err, boost_err\n", + " |\n", + " | chif(self, *args, **kwargs)\n", + " | chif is not implemented in this model. Will return (None, None),\n", + " |\n", + " | mf(self, *args, **kwargs)\n", + " | mf is not implemented in this model. Will return (None, None),\n", + " |\n", + " | vf(self, *args, **kwargs)\n", + " | vf is not implemented in this model. Will return (None, None),\n", + " |\n", + " | ----------------------------------------------------------------------\n", + " | Data descriptors inherited from surfinBH.surfinBH.SurFinBH:\n", + " |\n", + " | __dict__\n", + " | dictionary for instance variables\n", + " |\n", + " | __weakref__\n", + " | list of weak references to the object\n", + "\n" + ] + } + ], + "source": [ + "help(fit)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Evaluate the fits" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "q = 4.3 # Mass ratio q>=1\n", + "chiA = [0,0,0.6] # Spin of larger BH (z-direction only)\n", + "chiB = [0,0,-0.7] # Spin of smaller BH (z-direction only)\n", + "\n", + "## Evaluate the fits and GPR error estimate.\n", + "\n", + "#this model fits the supertranslation parameter alpha and the boost \n", + "#velocity v associated to the BMS transformation from the PN BMS frame to\n", + "#the remnant BMS frame\n", + "\n", + "# Outputs:\n", + "# - supertranslation parameter alpha and 1-sigma error estimate alpha_err\n", + "# NOTE: alpha is a complex array. Each component is a spherical harmonic mode of alpha\n", + "# the order of the modes in the array is (0,0),(1,-1),(1,0),(1,1),(2,-2),(2,-1),(2,0),(2,1),(2,2),...\n", + "# (same order used in spherical_functions and scri package)\n", + "# - boost velocity boost and 1-sigma error estimate boost_err\n", + "\n", + "alpha, boost, alpha_err, boost_err = fit.all(q, chiA, chiB)\n", + "\n", + "# to get a specific supertranslation mode from alpha you can also use the spherical_function package function LM_index\n", + "def LM_index(ell,m):\n", + " \"\"\" Returns the index in the alpha array corresponding to the ell,m mode\n", + " \"\"\"\n", + " return ell * (ell + 1) + m\n", + "\n", + "\n", + "# access a specific ell,m mode of alpha\n", + "ell,m = 2,0\n", + "alpha[LM_index(ell,m)]\n", + "# cartesian components of boost velocity\n", + "boost[0] # x-component\n", + "boost[1] # y-component\n", + "boost[2] # z-component" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.12.3" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/surfinBH/_fit_evaluators/__init__.py b/surfinBH/_fit_evaluators/__init__.py index f34a214..5864a4e 100644 --- a/surfinBH/_fit_evaluators/__init__.py +++ b/surfinBH/_fit_evaluators/__init__.py @@ -3,3 +3,4 @@ from .fit_7dq4 import Fit7dq4 from .fit_7dq4Emri import Fit7dq4Emri from .fit_3dq8_RD import Fit3dq8_RD +from .fit_3dq8BMS import Fit3dq8BMS diff --git a/surfinBH/_fit_evaluators/fit_3dq8BMS.py b/surfinBH/_fit_evaluators/fit_3dq8BMS.py new file mode 100644 index 0000000..2666052 --- /dev/null +++ b/surfinBH/_fit_evaluators/fit_3dq8BMS.py @@ -0,0 +1,204 @@ +import numpy as np +from surfinBH import surfinBH +import warnings + + +#============================================================================= +class Fit3dq8BMS(surfinBH.SurFinBH): + """ A class for the NRSur3dq8BMSRemnant model presented in Da Re et al., + arxiv:????.?????. This model is referred to as surfinBH3dq8BMS in the paper. + + This model predicts the supertransation modes up to ell = 8 and the + 3 components of the boost velocity of the BMS transformation from the + inspiral (PN) BMS frame to the remnant black hole BMS frame. The boost + velocity coincides with the remnant black hole kick velocity as observed + from the inspiral (PN) BMS frame. + + The model was trained on nonprecessing binary black hole systems. The fits + are done using Gaussian Process Regression (GPR) and also provide an error + estimate along with the fit value. + + This model has been trained in the parameter space: + q <= 8, |chiAz| <= 0.8, |chiBz| <= 0.8 + + ========================================================================= + Usage: + + import surfinBH + + # Load the fit + fit = surfinBH.LoadFits('NRSur3dq8BMSRemnant') + + #To evaluate the supertranslation parameter "alpha" and the boost velocity + "boost" together with their respective 1-sigma error estimates call: + alpha, boost, alpha_err, boost_err = fit.all(q, chiA, chiB, **kwargs) + + # NOTE: The ell=0 mode, corresponding to a time translation, + is set to be identically zero. + + # alpha is expressed as a complex array of spherical harmonics modes in the + order + (0,0),(1,-1),(1,0),(1,1),(2,-2),(2,-1),(2,0),... + (same convention as the spherical_functions and scri packages) + + The arguments for each of these call methods are as follows: + Arguments: + q: Mass ratio (q>=1) + + chiA: Dimensionless spin of the larger BH (array of size 3). + + chiB: Dimensionless spin of the smaller BH (array of size 3). + This model allows only nonprecessing spins, so only the + z-components of these arrays should be non-zero. + + Optional arguments: + allow_extrap: + If False, raises a warning when q > 8.1 or |chiA|,|chiB| > 0.81, + and raises an error when q > 10.1 or |chiA|,|chiB| > 1. + If True, allows extrapolation to any q and |chiA|,|chiB| <= 1. + Use at your own risk. + Default: False. + + The inspiral frame can be defined up to an arbitrary U(1) rotation about + the z-axis. We fix this gauge freedom using the same frame alignment choice + used in NrSur3dq8Remnant. The inspiral frame at -100M from the peak of the + waveform is defined as: + The z-axis is along the orbital angular momentum direction of the binary. + The x-axis is along the line of separation from the smaller BH to + the larger BH at this time. + The y-axis completes the triad. + The rotation leading to this frame choice is included in the BMS transformation + as outlined in arxiv:????.?????. + """ + + #------------------------------------------------------------------------- + def __init__(self, name): + + # Param limits beyond which to raise a warning + soft_param_lims = { + 'q': 8.1, + 'chiAmag': 0.81, + 'chiBmag': 0.81, + } + + # Param limits beyond which to raise an error + hard_param_lims = { + 'q': 10.1, + 'chiAmag': 1, + 'chiBmag': 1, + } + super(Fit3dq8BMS, self).__init__(name, soft_param_lims, hard_param_lims, + aligned_spin_only=True) + + #------------------------------------------------------------------------- + def _load_fits(self, h5file): + """ Loads fits from h5file and returns a dictionary of fits. """ + fits = {} + ell_max = 8 + keys_sup_real = [f'S_{ell}{m}_real' for ell in range(1,ell_max + 1) + for m in range(0,ell + 1)] + keys_sup_imag = [f'S_{ell}{m}_imag' for ell in range(1,ell_max + 1) + for m in range(1,ell + 1)] + keys_vel = ['boost_x','boost_y','boost_z'] + for key in keys_vel + keys_sup_real + keys_sup_imag: + fits[key] = self._load_scalar_fit(fit_key=key, h5file=h5file) + return fits + + #------------------------------------------------------------------------- + def _get_fit_params(self, x, fit_key): + """ Transforms the input parameter to fit parameters for the 3dq8BMS model. + That is, maps from [q, chiAz, chiBz] to [np.log(q), chiHat, chi_a] + chiHat is defined in Eq.(3) of 1508.07253. + chi_a = (chiAz - chiBz)/2. + """ + q, chiAz, chiBz = x + eta = q/(1.+q)**2 + chi_wtAvg = (q*chiAz+chiBz)/(1.+q) + chiHat = (chi_wtAvg - 38.*eta/113.*(chiAz + chiBz))/(1. - 76.*eta/113.) + chi_a = (chiAz - chiBz)/2. + fit_params = [np.log(q), chiHat, chi_a] + return fit_params + + #------------------------------------------------------------------------- + def mf(self, *args, **kwargs): + """mf is not implemented in this model. Will return (None, None),""" + return self._eval_wrapper('mf', *args, **kwargs) + + def chif(self, *args, **kwargs): + """chif is not implemented in this model. Will return (None, None),""" + return self._eval_wrapper('chif', *args, **kwargs) + + def vf(self, *args, **kwargs): + """vf is not implemented in this model. Will return (None, None),""" + return self._eval_wrapper('vf', *args, **kwargs) + + def all(self, *args, **kwargs): + """ Evaluates fit and 1-sigma error estimate for supertranslation parameter + alpha, and boost velocity of the BMS transformation from the inspiral BMS + frame to the remnant BMS frame. + Returns: + alpha, boost, alpha_err, boost_err + """ + return self._eval_wrapper('all', *args, **kwargs) + + #------------------------------------------------------------------------- + def _eval_wrapper(self, fit_key, q, chiA, chiB, **kwargs): + """ Evaluates the NRSur3dq8BMSRemnant model. + """ + + # mf,chif,vf not implemented for this model + if fit_key == 'mf': + mf, mf_err = None, None + return mf, mf_err + + if fit_key == 'chif': + chif, chif_err = None, None + return chif, chif_err + + if fit_key == 'vf': + vf, vf_err = None, None + return vf, vf_err + + chiA = np.array(chiA) + chiB = np.array(chiB) + + # Warn/Exit if extrapolating + allow_extrap = kwargs.pop('allow_extrap', False) + self._check_param_limits(q, chiA, chiB, allow_extrap) + + self._check_unused_kwargs(kwargs) + + x = [q, chiA[2], chiB[2]] + + ell_max = 8 + LM_total_size = ell_max * (ell_max + 2) + 1 + alpha = np.zeros((LM_total_size,), dtype=complex) + alpha_err = np.zeros((LM_total_size,), dtype=complex) + LMs = [(ell,m) for ell in range(1,ell_max + 1) for m in range(-ell,ell + 1)] + for i, (ell, m) in enumerate(LMs): + k = i + 1 #first slot is reserved for ell = 0 mode + sgn = int(np.sign(m)) + # the fits are trained with the variable S_{ell |m|} = 1/2(alpha_{ell |m|} + (-1)^m alpha_{ell -|m|}) + # the variable S_{ell |m|} coincides with alpha_{ell |m|} because alpha is real + S_real, S_real_err = self._evaluate_fits(x, f'S_{ell}{abs(m)}_real') + if m == 0: + alpha[k] = S_real + alpha_err[k] = S_real_err + # avoid the else to be faster + continue + S_imag, S_imag_err = self._evaluate_fits(x, f'S_{ell}{abs(m)}_imag') + # general formula alpha_{ell m} = (Re{alpha_{ell |m|}} + i sgn(m) Im{alpha_{ell |m|}}) sgn(m)^m for all ell,m + alpha[k] = (S_real + sgn * 1j * S_imag) * sgn**int(m) + alpha_err[k] = (S_real_err + sgn * 1j * S_imag_err) * sgn**int(m) + + boost_x, boost_x_err = self._evaluate_fits(x, 'boost_x') + boost_y, boost_y_err = self._evaluate_fits(x, 'boost_y') + boost_z, boost_z_err = self._evaluate_fits(x, 'boost_z') + boost = np.array([boost_x, boost_y, boost_z]) + boost_err = np.array([boost_x_err, boost_y_err, boost_z_err]) + + if fit_key == 'all': + return alpha, boost, alpha_err, boost_err + + + diff --git a/surfinBH/_loadFits.py b/surfinBH/_loadFits.py index 68fd0ac..7ad3efe 100644 --- a/surfinBH/_loadFits.py +++ b/surfinBH/_loadFits.py @@ -113,3 +113,13 @@ def DownloadData(name='all', data_dir=DataPath()): data_url = 'https://zenodo.org/records/13307895/files/fit_3dq8_RD.h5', refs = 'arxiv:2408.05300', ) + +#update arxiv number and url (zenodo) when available +fits_collection['NRSur3dq8BMSRemnant'] = FitAttributes( \ + fit_class = _fit_evaluators.Fit3dq8BMS, + desc = 'Fits for supertranslation parameter modes and boost velocity of the BMS' + ' transformation from the PN BMS frame to the remnant BMS frame ' + ' for nonprecessing BBH systems.', + data_url = 'https://raw.githubusercontent.com/gda98re/files/surfinBH/fit_3dq8BMS.h5', + refs = 'arxiv:????.?????', + ) diff --git a/surfinBH/surfinBH.py b/surfinBH/surfinBH.py index 8cb29e2..3c580b5 100644 --- a/surfinBH/surfinBH.py +++ b/surfinBH/surfinBH.py @@ -310,6 +310,18 @@ def all(self, *args, **kwargs): mf, chif, vf, mf_err_est, chif_err_est, vf_err_est chif, vf, chif_err_est and vf_err_est are arrays of size 3. + + In the case of NRsur3dq8BMSRemnant evaluates fit and 1-sigma error estimate for the + supertranslation parameter's modes and boost velocity of the BMS transformation from the PN + BMS frame to the remnant BMS frame. + Returns: + alpha, boost_vel, alpha_err, boost_vel_err + + alpha, alpha_err complex arrays of spherical harmonics modes from ell = 0 to ell = 8 with spherical_function package mode ordering. + [(0,0),(1,-1),(1,0),(1,1),(2,-2),...,(8,8)] + + boost_vel, boost_vel_err are arrays of size 3. """ return self._eval_wrapper('all', *args, **kwargs) + diff --git a/test/generate_regression_data.py b/test/generate_regression_data.py index 0f882f1..cfa7681 100644 --- a/test/generate_regression_data.py +++ b/test/generate_regression_data.py @@ -5,6 +5,35 @@ import surfinBH +def save_data_3dq8BMS(h5grp, fit, num_tests, kwargs={}): + + for i in range(num_tests): + # save each test evaluation as group + test_h5grp = h5grp.create_group('test_%d'%i) + + # Generate params randomly within allowed values + q, chiA, chiB = fit._generate_random_params_for_tests() + + # save params + test_h5grp.create_dataset('q', data=q) + test_h5grp.create_dataset('chiA', data=chiA) + test_h5grp.create_dataset('chiB', data=chiB) + + # save evaluations as a group + y_h5grp = test_h5grp.create_group('y') + + alpha, boost, alpha_err, boost_err = fit.all(q, chiA, chiB, **kwargs) + + # supertranslation + y = alpha, alpha_err + y_h5grp.create_dataset('alpha', data=y) + + # boost velocity + y = boost, boost_err + y_h5grp.create_dataset('boost', data=y) + + + def save_data(h5grp, fit, num_tests, kwargs={}): for i in range(num_tests): @@ -62,7 +91,10 @@ def save_data(h5grp, fit, num_tests, kwargs={}): # save regression without optional kwargs h5grp = h5file.create_group('No_kwargs') - save_data(h5grp, fit, num_tests) + if name_tag == '3dq8BMS': + save_data_3dq8BMS(h5grp, fit, num_tests) + else: + save_data(h5grp, fit, num_tests) # save regression with optional kwargs. # NOTE: this would be empty unless _extra_regression_kwargs() is diff --git a/test/regression_data/fit_3dq8BMS.h5 b/test/regression_data/fit_3dq8BMS.h5 new file mode 100644 index 0000000..900fb87 Binary files /dev/null and b/test/regression_data/fit_3dq8BMS.h5 differ diff --git a/test/test_interface.py b/test/test_interface.py index 2035efe..d5b4389 100644 --- a/test/test_interface.py +++ b/test/test_interface.py @@ -64,6 +64,11 @@ def test_interface(): # List of all available fits fit_names = surfinBH.fits_collection.keys() for name in fit_names: + + # since there are no methods for the fit parameters of 3dq8BMS yet, these tests + # cannot be applied to this model + if name == 'NRSur3dq8BMSRemnant': continue + # Load fit fit = surfinBH.LoadFits(name) diff --git a/test/test_regression.py b/test/test_regression.py index c5ec973..c1f104f 100644 --- a/test/test_regression.py +++ b/test/test_regression.py @@ -70,22 +70,38 @@ def test_fit_regression(): q = test_h5grp['q'][()] chiA = test_h5grp['chiA'][()] chiB = test_h5grp['chiB'][()] - - # remnant mass - y_reg = test_h5grp['y/mf'][()] - y_fit = fit.mf(q, chiA, chiB, **kwargs) - np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) - - # remnant spin - y_reg = test_h5grp['y/chif'][()] - y_fit = fit.chif(q, chiA, chiB, **kwargs) - np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) - # remnant kick - # Needed for NRSur7dq4EmriRemnant - if 'vf' in test_h5grp['y'].keys(): - y_reg = test_h5grp['y/vf'][()] - y_fit = fit.vf(q, chiA, chiB, **kwargs) + if name == 'NRSur3dq8BMSRemnant': + #compute fits + alpha, boost, alpha_err, boost_err = fit.all(q, chiA, chiB, **kwargs) + + # supertranslation + y_reg = test_h5grp['y/alpha'][()] + y_fit = alpha, alpha_err np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) + # boost velocity + y_reg = test_h5grp['y/boost'][()] + y_fit = boost, boost_err + np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) + + else: + # remnant mass + y_reg = test_h5grp['y/mf'][()] + y_fit = fit.mf(q, chiA, chiB, **kwargs) + np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) + + # remnant spin + y_reg = test_h5grp['y/chif'][()] + y_fit = fit.chif(q, chiA, chiB, **kwargs) + np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) + + # remnant kick + # Needed for NRSur7dq4EmriRemnant + if 'vf' in test_h5grp['y'].keys(): + y_reg = test_h5grp['y/vf'][()] + y_fit = fit.vf(q, chiA, chiB, **kwargs) + np.testing.assert_allclose(y_fit, y_reg, rtol=rtol, atol=atol) + + regression_h5file.close()