double and triple integrals kind of convergence rate estimate could be used to verify the implementationbut this topic is beyond the scope of this book test function for function with random numbers to make test functionwe need unit test that has identical behavior each time we run the test this seems difficult when random numbers are involvedbecause these numbers are different every time we run the algorithmand each run hence produces (slightlydifferent result standard way to test algorithms involving random numbers is to fix the seed of the random number generator then the sequence of numbers is the same every time we run the algorithm assuming that the montecarlo_double function workswe fix the seedobserve certain resultand take this result as the correct result provided the test function always uses this seedwe should get exactly this result every time the montecarlo_double function is called our test function can then be written as shown below def test_montecarlo_double_rectangle_area()"""check the area of rectangle ""def (xy)return ( if ( < < and < < else - embedded rectangle np random seed( must fix the seedi_expected computed with this seed i_computed montecarlo_doublelambda xy gx nassert abs(i_expected i_computed - (see the file mc_double py integral over circle the test above involves trivial function xyd we should also test non-constant function and more complicated domain let be circle at the origin with radius and let ry this choice makes it possible to compute an exact resultin polar coordinatesf xy/dxdy simplir fies to dr = we must be prepared for quite crude approximations that fluctuate around this exact result as in the test case abovewe experience better results with larger number of points when we have such evidence for working implementationwe can turn the test into proper test function here is an exampledef test_montecarlo_double_circle_r()"""check the integral of over circle with radius ""def (xy)xcyc center radius return ** (( -xc)** ( -yc)** exactintegral of * *dr over circle with radius becomes *pi* / * ** import sympy sympy symbols(' 'i_exact sympy integrate( *sympy pi* * ( )