From 054e270fd6e61f2f00783bcc05e441b4e3f86092 Mon Sep 17 00:00:00 2001 From: vsnever Date: Fri, 2 Aug 2024 10:29:34 +0200 Subject: [PATCH 1/4] Replace deprecated scipy.integrate.quadrature with GaussianQuadrature in test_lineshapes.py --- cherab/core/tests/test_lineshapes.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/cherab/core/tests/test_lineshapes.py b/cherab/core/tests/test_lineshapes.py index b6827256..c245e0c3 100644 --- a/cherab/core/tests/test_lineshapes.py +++ b/cherab/core/tests/test_lineshapes.py @@ -20,7 +20,6 @@ import numpy as np from scipy.special import erf, hyp2f1 -from scipy.integrate import quadrature from raysect.core import Point3D, Vector3D from raysect.core.math.function.float import Arg1D, Constant1D @@ -372,13 +371,14 @@ def stark_lineshape_sigma_minus(x): weight_poly_coeff = [5.14820e-04, 1.38821e+00, -9.60424e-02, -3.83995e-02, -7.40042e-03, -5.47626e-04] lorentz_weight = np.exp(np.poly1d(weight_poly_coeff[::-1])(np.log(fwhm_lorentz / fwhm_full))) + integrator_pi = GaussianQuadrature(stark_lineshape_pi, relative_tolerance=integrator.relative_tolerance) + integrator_sigma_plus = GaussianQuadrature(stark_lineshape_sigma_plus, relative_tolerance=integrator.relative_tolerance) + integrator_sigma_minus = GaussianQuadrature(stark_lineshape_sigma_minus, relative_tolerance=integrator.relative_tolerance) + for i in range(bins): - lorentz_bin = 0.5 * sin_sqr * quadrature(stark_lineshape_pi, wavelengths[i], wavelengths[i + 1], - rtol=integrator.relative_tolerance)[0] / delta - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quadrature(stark_lineshape_sigma_plus, wavelengths[i], wavelengths[i + 1], - rtol=integrator.relative_tolerance)[0] / delta - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quadrature(stark_lineshape_sigma_minus, wavelengths[i], wavelengths[i + 1], - rtol=integrator.relative_tolerance)[0] / delta + lorentz_bin = 0.5 * sin_sqr * integrator_pi(wavelengths[i], wavelengths[i + 1]) / delta + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * integrator_sigma_plus(wavelengths[i], wavelengths[i + 1]) / delta + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * integrator_sigma_minus(wavelengths[i], wavelengths[i + 1]) / delta ref_value = lorentz_bin * lorentz_weight + gaussian[i] * (1. - lorentz_weight) self.assertAlmostEqual(ref_value, spectrum.samples[i], delta=1e-9, msg='StarkBroadenedLine.add_line() method gives a wrong value at {} nm.'.format(wavelengths[i])) From 08a14b7214c2156412ccec2f62eb034d6bd66511 Mon Sep 17 00:00:00 2001 From: vsnever Date: Fri, 2 Aug 2024 13:02:25 +0200 Subject: [PATCH 2/4] Replace GaussianQuadrature with scipy.integrate.quad in test_lineshapes.py --- cherab/core/tests/test_lineshapes.py | 17 +++++++---------- 1 file changed, 7 insertions(+), 10 deletions(-) diff --git a/cherab/core/tests/test_lineshapes.py b/cherab/core/tests/test_lineshapes.py index c245e0c3..4e51e38a 100644 --- a/cherab/core/tests/test_lineshapes.py +++ b/cherab/core/tests/test_lineshapes.py @@ -20,6 +20,7 @@ import numpy as np from scipy.special import erf, hyp2f1 +from scipy.integrate import quad from raysect.core import Point3D, Vector3D from raysect.core.math.function.float import Arg1D, Constant1D @@ -295,7 +296,7 @@ def test_stark_broadened_line(self): line = Line(deuterium, 0, (6, 2)) # D-delta line target_species = self.plasma.composition.get(line.element, line.charge) wavelength = 656.104 - integrator = GaussianQuadrature(relative_tolerance=1.e-5) + integrator = GaussianQuadrature(relative_tolerance=1.e-8) stark_line = StarkBroadenedLine(line, wavelength, target_species, self.plasma, self.atomic_data, integrator=integrator) # spectrum parameters @@ -371,16 +372,12 @@ def stark_lineshape_sigma_minus(x): weight_poly_coeff = [5.14820e-04, 1.38821e+00, -9.60424e-02, -3.83995e-02, -7.40042e-03, -5.47626e-04] lorentz_weight = np.exp(np.poly1d(weight_poly_coeff[::-1])(np.log(fwhm_lorentz / fwhm_full))) - integrator_pi = GaussianQuadrature(stark_lineshape_pi, relative_tolerance=integrator.relative_tolerance) - integrator_sigma_plus = GaussianQuadrature(stark_lineshape_sigma_plus, relative_tolerance=integrator.relative_tolerance) - integrator_sigma_minus = GaussianQuadrature(stark_lineshape_sigma_minus, relative_tolerance=integrator.relative_tolerance) - for i in range(bins): - lorentz_bin = 0.5 * sin_sqr * integrator_pi(wavelengths[i], wavelengths[i + 1]) / delta - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * integrator_sigma_plus(wavelengths[i], wavelengths[i + 1]) / delta - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * integrator_sigma_minus(wavelengths[i], wavelengths[i + 1]) / delta - ref_value = lorentz_bin * lorentz_weight + gaussian[i] * (1. - lorentz_weight) - self.assertAlmostEqual(ref_value, spectrum.samples[i], delta=1e-9, + lorentz_bin = 0.5 * sin_sqr * quad(stark_lineshape_pi, wavelengths[i], wavelengths[i + 1], epsrel=1.e-8)[0] + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_plus, wavelengths[i], wavelengths[i + 1], epsrel=1.e-8)[0] + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_minus, wavelengths[i], wavelengths[i + 1], epsrel=1.e-8)[0] + ref_value = lorentz_bin / delta * lorentz_weight + gaussian[i] * (1. - lorentz_weight) + self.assertAlmostEqual(ref_value, spectrum.samples[i], delta=1e-8, msg='StarkBroadenedLine.add_line() method gives a wrong value at {} nm.'.format(wavelengths[i])) def test_beam_emission_multiplet(self): From 95197e51fd8eb7a7e9adb9c3fd8ecf0d70e01efe Mon Sep 17 00:00:00 2001 From: vsnever Date: Fri, 2 Aug 2024 14:52:35 +0200 Subject: [PATCH 3/4] Define relative_tolerance in test_stark_broadened_line() and compare test/true ratio to 1 instead of comparing absolute values. --- cherab/core/tests/test_lineshapes.py | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/cherab/core/tests/test_lineshapes.py b/cherab/core/tests/test_lineshapes.py index 4e51e38a..4f2ef142 100644 --- a/cherab/core/tests/test_lineshapes.py +++ b/cherab/core/tests/test_lineshapes.py @@ -296,7 +296,8 @@ def test_stark_broadened_line(self): line = Line(deuterium, 0, (6, 2)) # D-delta line target_species = self.plasma.composition.get(line.element, line.charge) wavelength = 656.104 - integrator = GaussianQuadrature(relative_tolerance=1.e-8) + relative_tolerance = 1.e-8 + integrator = GaussianQuadrature(relative_tolerance=relative_tolerance) stark_line = StarkBroadenedLine(line, wavelength, target_species, self.plasma, self.atomic_data, integrator=integrator) # spectrum parameters @@ -373,12 +374,17 @@ def stark_lineshape_sigma_minus(x): lorentz_weight = np.exp(np.poly1d(weight_poly_coeff[::-1])(np.log(fwhm_lorentz / fwhm_full))) for i in range(bins): - lorentz_bin = 0.5 * sin_sqr * quad(stark_lineshape_pi, wavelengths[i], wavelengths[i + 1], epsrel=1.e-8)[0] - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_plus, wavelengths[i], wavelengths[i + 1], epsrel=1.e-8)[0] - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_minus, wavelengths[i], wavelengths[i + 1], epsrel=1.e-8)[0] + lorentz_bin = 0.5 * sin_sqr * quad(stark_lineshape_pi, wavelengths[i], wavelengths[i + 1], epsrel=relative_tolerance)[0] + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_plus, wavelengths[i], wavelengths[i + 1], epsrel=relative_tolerance)[0] + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_minus, wavelengths[i], wavelengths[i + 1], epsrel=relative_tolerance)[0] ref_value = lorentz_bin / delta * lorentz_weight + gaussian[i] * (1. - lorentz_weight) - self.assertAlmostEqual(ref_value, spectrum.samples[i], delta=1e-8, - msg='StarkBroadenedLine.add_line() method gives a wrong value at {} nm.'.format(wavelengths[i])) + if ref_value: + print(ref_value) + self.assertAlmostEqual(spectrum.samples[i] / ref_value, 1., delta=relative_tolerance, + msg='StarkBroadenedLine.add_line() method gives a wrong value at {} nm.'.format(wavelengths[i])) + else: + self.assertAlmostEqual(ref_value, spectrum.samples[i], delta=relative_tolerance, + msg='StarkBroadenedLine.add_line() method gives a wrong value at {} nm.'.format(wavelengths[i])) def test_beam_emission_multiplet(self): # Test MSE line shape From 1e08f1e2ab7254b60b7fb85a9aa636ed4902cf0a Mon Sep 17 00:00:00 2001 From: vsnever Date: Fri, 2 Aug 2024 15:20:08 +0200 Subject: [PATCH 4/4] Remove print() in test_lineshapes.py. --- cherab/core/tests/test_lineshapes.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/cherab/core/tests/test_lineshapes.py b/cherab/core/tests/test_lineshapes.py index 4f2ef142..e2122d41 100644 --- a/cherab/core/tests/test_lineshapes.py +++ b/cherab/core/tests/test_lineshapes.py @@ -374,12 +374,14 @@ def stark_lineshape_sigma_minus(x): lorentz_weight = np.exp(np.poly1d(weight_poly_coeff[::-1])(np.log(fwhm_lorentz / fwhm_full))) for i in range(bins): - lorentz_bin = 0.5 * sin_sqr * quad(stark_lineshape_pi, wavelengths[i], wavelengths[i + 1], epsrel=relative_tolerance)[0] - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_plus, wavelengths[i], wavelengths[i + 1], epsrel=relative_tolerance)[0] - lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_minus, wavelengths[i], wavelengths[i + 1], epsrel=relative_tolerance)[0] + lorentz_bin = 0.5 * sin_sqr * quad(stark_lineshape_pi, wavelengths[i], wavelengths[i + 1], + epsrel=relative_tolerance)[0] + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_plus, wavelengths[i], wavelengths[i + 1], + epsrel=relative_tolerance)[0] + lorentz_bin += (0.25 * sin_sqr + 0.5 * cos_sqr) * quad(stark_lineshape_sigma_minus, wavelengths[i], wavelengths[i + 1], + epsrel=relative_tolerance)[0] ref_value = lorentz_bin / delta * lorentz_weight + gaussian[i] * (1. - lorentz_weight) if ref_value: - print(ref_value) self.assertAlmostEqual(spectrum.samples[i] / ref_value, 1., delta=relative_tolerance, msg='StarkBroadenedLine.add_line() method gives a wrong value at {} nm.'.format(wavelengths[i])) else: