Controlling the order of Gaussian quadrature

1

Hello,

I have written a code for solving Laplace’s equation which works well on some grids but not so well when there is a short distance from an element to another node point. I would like to try increasing the order of the Gaussian quadrature scheme. How do I do this?

In the BEM++ project there were commands like:

bempp.api.global_parameters.quadrature.medium.double_order = 4
bempp.api.global_parameters.quadrature.far.double_order = 4

What are the equivalent commands in bempp-cl?

To a more general question: Is it possible to print out the results of the command:

bempp.api.utils.parameters._Quadrature

and get the default values used in bempp-cl?

Thanks very much,

Peter.

2

Hi Peter,

Welcome to this forum.

The quadrature order is part of the global parameters in bempp-cl. Default parameters are set when the library is imported.

You can print the quadrature order by:

print(bempp.api.GLOBAL_PARAMETERS.quadrature.regular)
print(bempp.api.GLOBAL_PARAMETERS.quadrature.singular)

This prints the quadrature order for regular and singular integrals. Bempp-cl does not make a distinction between near, medium, and far interactions anymore.
Notice that these are the order of quadrature, not the number of points. This information can be accessed by

print(bempp.api.integration.triangle_gauss.get_number_of_quad_points(4))
print(bempp.api.integration.triangle_gauss.rule(4))

The second command gives the location and weights of the quadrature point on a unit triangle.

To change the order of quadrature, the easiest way is to set the global parameters:

bempp.api.GLOBAL_PARAMETERS.quadrature.regular = 6

By default, all operators take the global parameters when assembling the matrix. So set this before defining the operators.

If op is an operator, you can check the quadrature order by

print(op.parameters.quadrature.regular)

If you need to define different operators with different quadrature order in the same script, you can pass an operator object to the parameters argument of an operator.

I hope this helps and good luck with using BEMPP.

Best,
Elwin