.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "_build/auto_examples/convergence/Integration_on_segment.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr__build_auto_examples_convergence_Integration_on_segment.py: Convergence analysis of integration methods on segment ====================================================== .. GENERATED FROM PYTHON SOURCE LINES 6-7 The following plots show the relative error made by different quadrature methods available as :class:`choreo.segm.quad.QuadTable` on the integration of a smooth function (a scaled exponential). .. GENERATED FROM PYTHON SOURCE LINES 7-74 .. code-block:: Python all_integrands = [ "exp", ] methods = [ 'Gauss' , 'Radau_I' , 'Radau_II' , 'Lobatto_III' , 'Cheb_I' , 'Cheb_II' , 'ClenshawCurtis' , 'NewtonCotesOpen' , 'NewtonCotesClosed' , ] all_nsteps = [5,6] refinement_lvl = np.array(range(1,100)) def setup(fun_name, quad_method, quad_nsteps, nint, DoEFT = False): return {'fun_name': fun_name, 'quad_method': quad_method, 'quad_nsteps': quad_nsteps, 'nint': nint, 'DoEFT':DoEFT} all_args = { 'integrand': all_integrands , 'quad_method': methods , 'quad_nsteps': all_nsteps , 'nint': refinement_lvl , } all_funs = [choreo.segm.test.Quad_cpte_error_on_test] bench_filename = os.path.join(bench_folder, 'quad_cvg_bench.npz') all_errors = pyquickbench.run_benchmark( all_args , all_funs , setup = setup , mode = "scalar_output" , filename = bench_filename , ForceBenchmark = ForceBenchmark , StopOnExcept = True , pooltype = "process" , ) plot_intent = { 'integrand': "subplot_grid_x" , 'quad_method': "curve_color" , 'quad_nsteps': "subplot_grid_y" , 'nint': "points" , pyquickbench.fun_ax_name : "same" , } fig, ax = pyquickbench.plot_benchmark( all_errors , all_args , all_funs , plot_intent = plot_intent , title = 'Relative error of quadrature' , ylabel = "Relative error" , ) plot_xlim = ax[0,0].get_xlim() fig.tight_layout() fig.show() .. image-sg:: /_build/auto_examples/convergence/images/sphx_glr_Integration_on_segment_001.png :alt: Relative error of quadrature, integrand : exp, quad_nsteps : 5, integrand : exp, quad_nsteps : 6 :srcset: /_build/auto_examples/convergence/images/sphx_glr_Integration_on_segment_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 118-120 The following plots give the measured convergence rate as a function of the number of quadrature subintervals. The dotted lines are theoretical convergence rates. .. GENERATED FROM PYTHON SOURCE LINES 120-150 .. code-block:: Python plot_ylim = [0,15] fig, ax = pyquickbench.plot_benchmark( all_errors , all_args , all_funs , transform = "pol_cvgence_order" , plot_xlim = plot_xlim , plot_ylim = plot_ylim , logx_plot = True , clip_vals = True , stop_after_first_clip = True , plot_intent = plot_intent , title = 'Approximate convergence rate' , ) for iy, nstep in enumerate(all_nsteps): for method in methods: quad = choreo.segm.multiprec_tables.ComputeQuadrature(nstep, method=method) th_order = quad.th_cvg_rate xlim = ax[iy,0].get_xlim() ax[iy,0].plot(xlim, [th_order, th_order], linestyle='dotted') fig.tight_layout() fig.show() .. image-sg:: /_build/auto_examples/convergence/images/sphx_glr_Integration_on_segment_002.png :alt: Approximate convergence rate, integrand : exp, quad_nsteps : 5, integrand : exp, quad_nsteps : 6 :srcset: /_build/auto_examples/convergence/images/sphx_glr_Integration_on_segment_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 151-158 We can see 3 distinct phases on these plots: * A first pre-convergence phase, where the convergence rate is growing towards its theoretical value. the end of the pre-convergence phase occurs for a number of sub-intervals roughtly independent of the convergence order of the quadrature method. * A steady convergence phase where the convergence remains close to the theoretical value * A final phase, where the relative error stagnates arround 1e-15. The value of the integral is computed with maximal accuracy given floating point precision. The approximation of the convergence rate is dominated by seemingly random floating point errors. Now, if the number of integration steps is increased beyond what is necessary, these floating points error will start to steadily grow. This effect can be mitigated using a compensated summation algorithm. This is enabled by setting the argument ``DoEFT = True`` as shown below. .. GENERATED FROM PYTHON SOURCE LINES 158-200 .. code-block:: Python bench_filename = os.path.join(bench_folder,'quad_cvg_bench_EFT.npz') refinement_lvl = np.array([2**n for n in range(30)]) all_args = { 'integrand': all_integrands , 'quad_method': ['Gauss'] , 'quad_nsteps': [1, 2, 3, 4] , 'nint': refinement_lvl , 'DoEFT': [False, True] , } plot_intent = { 'integrand': "subplot_grid_x" , 'quad_method': "subplot_grid_y" , 'quad_nsteps': "curve_color" , 'nint': "points" , pyquickbench.fun_ax_name : "same" , 'DoEFT' : "curve_linestyle" , } all_errors = pyquickbench.run_benchmark( all_args , all_funs , setup = setup , mode = "scalar_output" , filename = bench_filename , ForceBenchmark = ForceBenchmark , StopOnExcept = True , pooltype = "process" , ) fig, ax = pyquickbench.plot_benchmark( all_errors , all_args , all_funs , plot_intent = plot_intent , title = 'Relative error of quadrature' , ylabel = "Relative error" , ) .. image-sg:: /_build/auto_examples/convergence/images/sphx_glr_Integration_on_segment_003.png :alt: Relative error of quadrature :srcset: /_build/auto_examples/convergence/images/sphx_glr_Integration_on_segment_003.png :class: sphx-glr-single-img .. _sphx_glr_download__build_auto_examples_convergence_Integration_on_segment.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: Integration_on_segment.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: Integration_on_segment.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: Integration_on_segment.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_