Montecarlo Methods to Eficently Compute Volume of Solids

Abstract

In this work we present and compare chosen Monte Carlo methods, used to calculate solid's volume. Such solids are defined with a certain surface (, ) wich results in classical methods being hard or impossible to apply. In designating those volume, due to the desire to shorten calculation time, we will apply parallel calculations.

Keywords

• eol>Double integral
• Monte Carlo methods
• parallel computations

1. Introduction

above and below respectively. If the domain D s normal (e.g) with respect to the x-axis, i.e. can be described by dependencies: In classic task of calculating solid’s volume we use double integral (or triple integral), often utilizing variable’s conversion. Most often they are polar, cylindrical or spherical coordinates. The use of this approach most often occurs, when surfaces limiting given solid are given in open form ( 2 ) = (, )en however, whether for curves in a plane or for surfaces in space, such limitation are given diferent hen the volume ( 1 ) can be described using ( 2 ), in the (non-classical) way. For curve this could be for example: form: parametric or polar form. In case of a surface from 3 surface, the situation may look similar, but then there is an additional variable. Wanting to calculate area of a given flat surface or volume of a surface in space, appropriate formulas should be designated. Authors discussed their use in [1] for flat areas and in [ 2] for spatial areas.

Especially in the 3 space, where surface limiting given volume is dependent from (, ) suitable integral de- ( 3 ) scribing this volume may be hard or even impossible to determine analytically. In this work we propose the use of certain approximate methods on calculating numerically given volume, independent from the form of the funcion (, ) 2. Domains in 3 space In which the volume V is described by the integral:

where continuous in D domain functions 1 and 2 are limiting the space above and below D domain, from ( 1 )

To describe such domain also can be used polar coordinates, cylindrical or spherical. Remembering the Jacobian determinant of given variables swap and describing right Ddomain (for flat surfaces) or Space (in case of describing solid, e.g. using spherical coordinates) we get appropriate, known to us equations for calculating volume of the solid with given coordinates. We want to focus on the case, in which area limiting given solid is defined by a function (, ) this case any point of this surface is deifned by three variables. The first of them is the distance from this point to the beginning of the coordinate system. The second is an angle between positive Ox is and a view, on a Oxy face of a segment joining the beginning of a coordinate system and given point. The third one is an angle between positive Oz is and a segment, joining the beginning of a coordinate system and given point (view Fig. 2). If for certain angles and for certain angles e define non-negative function (, ) then all points received this way will define a certain surface in 3 space.

It turns out that then the volume of a solid defined this way can e described by integral [2]:

In this connection that the form of the function can lead to integrals that do not express themselves with elementary functions, thereby the values of them cannot be determined precisely, therefore we will use some numeric methods to calculate approximate volume.

3. Monte Carlo methods

Monte Carlo methods appeared (most frequently this is assumed, although the theory describing these methods appeared the most frequently) in 1949, after the publication of the work [3], and allowed their creators to solve a computationally complex problem, that emerged during the work on the construction of the atomic bomb. The Monte Carlo methods family provides an approximate solution to a wide class of computationally complex problems, including numerical integration. We will focus on this group of problems. The computational complexity and error estimation obtained by this method show [4] that as the dimension increases, Monte Carlo methods become more and more useful in relation to classical methods (e.g. to quadrature based methods). We will deal with the double integral, but the ideas of the Monte Carlo methods, for the sake of clarity, will be shown on the example of the integral.

3.1. Monte Carlo methods no. 1

In classical terms, the domain between the graphs of the functions f(x) and g(x)we assume that () ⩾ () for ⩾ ⩾ can be approximated as follows: enter the ifrst graphs of the functions f and g into the rectangle R: 1–9 ( 5 ) ( 6 ) ( 7 ) ( 8 )

where c ⩾ min g(x), d ≤ max f(x) The next step is to draw n points belonging to the rectangle R. If this choice is random (from uniform distribution) and independent, then these points should cover the rectangle R evenly. Therefore, the ratio of the number of m points lying between the graphs of the functions f and g to the number of all points n is the same as the ratio of the fields: S– between the graphs of the functions f and g |R| – area of the rectangle R We have:

In Fig. 2 we present these dependencies graphically, with green color marking the points that went to the area of interest, and red color those that did not go to this area. 3.1.1. Example 1 We will show how it is possible to use this method to ifnd the approximate value of the number . For this we consider the integral:

Knowing its exact value, we can estimate the value of the number rising the formula ( 4 ). Let’s take a rectangle: is divided evenly, that is Δ = Δ = − . Thanks to this, formula ( 5 ) will take the form: conducted 50 draws (we want to check the stability of the method) and determined the following values: min – smallest area value, max – largest area value, S- average area value (as the arithmetic mean of the 50 results obtained), – standard deviation of the obtained results, Δ – calculation error in the formula Δ = | − ¯| the area S is the exact value of the area, and thus the number |e| and ¯ s the obtained value.

The conducted experiments show that this method is stable and gives a good approximation. The Monte Carlo error theory says that the error obtained is proportional to the number √ , where < 0 is a constant and s the number of points drawn. Figure 3 shows that this relationship is preserved (in this figure, 3 has the value 0.5)

3.2. Monte Carlo methods no. 2

In this method, we have a slightly diferent approach to determining the value of the integral previously approximated by the formula ( 4 ). This time we refer to the definition of the Riemann definite integral:

where , = 1, 2, ... are any points up to the subinterval [, +1], = 0, 1, . . . of the interval [, ] These sub-intervals are determined by a certain natural divi- 3.3. Other Monte Carlo methods sion, Δ is the length of such a sub-interval. If we run ndraws again (the assumptions for the drawing are the same as before), we can assume that the interval [, ] In addition to the many advantages of Monte Carlo methods, such as simplicity of implementation, regardless of the dimension of the task, or the non-growing of errors

where , 1, 2, . . . n are the drawn points of the interval [a, b. On this basis, the field described by formula ( 4 ) will take the form: 3.2.1. Example Now let’s do the same experiment as in Example 3.2.1.

Table 2 shows the results of these experiments, and Fig. 4 we present the equivalent of Fig. 3, where the constant c as the value 0.1.

Intuitively, it is dificult to estimate which of the formulas ( 4 ) or ( 6 ) is represented by the greater error. It ( 9 ) turns out that the formula ( 6 ) is more useful, and the results summarized in Tables 1 and 2 confirm this theory ( 10 ) ( 11 ) as the dimension increases, there are also disadvantages of this method. The most frequently mentioned is the 1 value of the obtained errors, which is proportional to √ , where nis the number of points drawn. To increase the error by one row, e.g. from 1/10 to 1/100 , the number of points drawn should be increased from 100 to 104 . The methods discussed by us in this part of the work want to counteract this problem. There are many diferent approaches to such a task, we will cover the three most commonly used.

3.3.1. Separating of the main part

This method consists in finding a certain function which is close to the integral function, but for which we can find the exact value of the integral. To refer to our example of finding the approximate value of e where we integrate the function ( = + 1), we can use the well-known Maclaurin series expansion of the function Also note that:

3.3.2. Importance sampling

This method [5] is based on the assumption that more points should be drawn in more significant parts of the interval [, ]. To do this, one needs to find a function () which in the interval[, ] is the probability density, i.e. () > 0

Of course, there are infinitely many such functions, and its choice has an impact on the resulting error. Most often, choosing the best such function is impossible, because it is related to the value of the integral (which is unknown). In practice, the approach is to select the function ()in such a way that, similar to the previous example, it is close to the function () i.e. that (()) ≈ where c is a constant. After finding such a function, we can use the formula:

which leads to the formula for the approximate value of this integral: 1–9 (17) (18) (19) (20) and in the discussed example, this formula takes the form: which leads to the formula for the approximate value of this integral:

Since the diference between the functions f(x and g(x) is determined, so:

where = 1, 2, . . . are random interval points [, ] from uniform distribution.

where = 1, 2, . . . are random interval points [, ] In the discussed example, using the previously found function () = 2 + + 0.52 we can take: () = 38 ()

So we will get the following formula which allows us to find the approximate value of the number ( 12 ) ( 13 ) ( 14 ) (15) (16) where = 1, 2, . . . are random distribution () (You can read about the drawing of points from a given distribution in e.g. [6, 7, 8]).

3.3.3. Stratified sampling

This approach [9] is similar to the previous one. This time we divide the integration interval [a,b into k subintervals and randomize points in each of them. The numbers e chosen so that their sum is and that their size is proportional to the rate of change of the integral value of the function () Based on this approach, the sought integral can be approximated by the formula: () is the i-th number drawn from the j-th where subinterval, and is the length of the jth sub-interval. 4. Monte Carlo methods for

If the interval [0, 1] we will divide it into 5 equal parts double integrals and choose the next sub-ranges accordingly 1/25 , 3/25 , 5/25 , 7/25 and 9/25 the number of all drawn points, The advantages of the Monte Carlo method, highlighted then the formula ( 9 ) adjusted to approximate the value by us previously, are perfect for double integrals, both of the number e will take the form: the simplicity of implementation (for both classical methods and their improvements) and the estimated error are transferred to double integrals (in general, to multiple integrals). As we mentioned in the introduction, we want to use Monte Carlo methods for integrals ( 3 ). This integral is formed when the surface bounding a given spatial domain is given by the function(, ) [2] authors discussed the formula ( 3 ) and applied it to several surfaces, incl. peanut, bluebell, shell and exotic fruit. In this paper, (22) we will use three approaches to the problem of determining these integrals using the Monte Carlo method I, the 3.3.4. Result for other Monte Carlo methods Monte Carlo method II and parallel calculations used in the Monte Carlo method II, the aim of which is to shorten As for the previous methods, for the same values of n we the time needed to determine a given integral. Moreover, carry out the process of approximating the value of the as some integrals ( 3 ) are reduced to a single integral in number e using the formulas ( 7 ), ( 8 ) and ( 9 ). In tables the above-mentioned surfaces, we will consider domains 3, 4 and 5 we will collect the appropriate results and in both on the plane and in space. Parallel computations drawing 5 we illustrate the absolute error Δ similarly as will be performed using a multithreaded process, and the we did in the drawings 3 and 4. aaabbbccc process itself is very simple. As the points

The conducted research confirms that, compared to drawn are independent and do not afect each other, as the best known Monte Carlo method I, all other methods well as the calculations of the function values in these are more efective. For the same number n of drawn points, it is enough to divide the points, they give a better (statistically) approximation of number of draws n into k equal parts and each thread the value of e will perform its calculations in parallel. After making the calculations by all threads, the results obtained are summed up to the final result. Note that the formula ( 6 ) can be implemented in parallel – the number n of all points is divided into k equal parts, where k is the number of threads. Each thread will do its own sum:

where ℎ() = () − () are random points from the interval [, ] And then we will create the final sum:

This process is shown in Fig. 6. (27) The calculations were carried out on a computer with a 12-thread Intel Core i7-3930k CPU 3.20 GHz processor, where, like before, m is the number of points that hit 16 Gb RAM and a 64 bit Windows 7 system. Programs the inside of the bell, and n is the number of all drawn that implement all Monte Carlo methods were written points of the cuboid However, using the spaces perin Mathematica 12.2 (including those for parallel com- formed in the work [2] the volume of the bell is described puting) [10]. We will now present the application of the by the single integral: proposed Monte Carlo methods on the example of two surfaces introduced and discussed in [2] – bell and exotic lfower (28) (23) (24) (25) (26)

This time, the equivalent of formula ( 4 ) should be a formula adapted to the size, so the volume V of the sought solid would be approximate by the formula:

So we have a situation where (without taking into account the unit) the volume of the bell equals (numerically) the area under the curve:

for x in [0, ], which means that it is enough to draw n points of the rectangle R.

(the number 12 can be determined on the basis in Fig. 8 which shows the graph of the function () along with the drawn = 103 points) and use the formula ( 4 ).

After conducting the experiment consisting in drawing = 9! points of the rectangle e obtained the following approximate result of the bell volume value: = 39.3611 what the error corresponds to Δ = 0.0214 This result, which will be important when comparing the (29) (30)

4.1. Bell Area

This area is represented by the relationship: (, ) = + (2 ) + (3 ) and its graph, for a= 2 is shown in Fig. 7.

As we showed in [2] the volume of the solid limited by this area is: which for a= 2 gives the value 39.382572945...

4.1.1. Monte Carlo method no. I

According to the idea of this method, the surface of the bell will be entered in a cuboid methods, was achieved in = 3.19802 seconds. Proceeding similarly to chapter 3 we then conducted a series of 50 such experiments, and their results (designations analogous to the previous ones with an additional value of the calculation time ¯ are presented in Table 6.

4.1.2. Monte Carlo method no.II

As before, we will go straight to the calculation of the area under the curve ( 10 ) and use the formula ( 6 ). Therefore, we will select = 9! points from the interval [0, ]and a searched area (and thus the volume of the bell) we will determine approximately. This time we got an area (volume) equal to 39.3762, which corresponds to error Δ = 0.0064. We obtained this result in = 1.903 seconds. After conducting 50 such experiments, we obtained the results summarized in Tab. 7.

4.1.3. Parallel Monte Carlo method no.II

Monte Carlo methods are perfect for parallelization. We can do it because the independent draws of n points from the interval [, ] is equivalent to k draws of n/k points 4.2.1. Monte Carlo method no. I from this interval (we sum the k set obtained in this way, This time we will refer to the formula ( 4 ) which says that obtaining one n-elementary set). Let us now carry out the surface of an exotic fruit is limited by an area whose experiments analogous to the above experiments, except what for a=1 gives value of 7.120943348. that we will use parallel computations with the use of = 2, 3, ..., 10 threads. The results of these experiments are collected in Table 8, we will pay special attention to the average time ¯ – its change is to prove the usefulness of parallel computations. As indicated by the data in table 8, the introduction of calculations, on the one hand does not worsen the obtained volume results, and on the other hand reduces the time needed to obtain this result.

The comparison of the dependence of this time t to the number of threads k is shown in Fig. 9.

4.2. Exotic fruit surface

This area will be represented by the relationship: 0 ≤ ≤ 2, 0 ≤ ≤ and ≥ 1, and its graph for a= 1 is shown in Fig. 10.

As we showed in [2], the volume of a solid limited by this area is: (31) (32) volume corresponds to the volume of the area limited from above by the area (, ) above the rectangle where which means that it is enough to draw n points of the cuboid

(the number 8/3 can be determined on the basis in Fig. 11, which shows the graph of the function (, ) with the drawn = 104 points) and use the formula ( 6 ), which this time will have to change by one form dimension: where m is the number of points that hit the inside of the fruit (and thus are in the cuboid C below the surface (, ))

After conducting a series of 50 experiments described above, assuming n= 9!, we obtained the results presented in table Tab. 9.

4.2.2. Monte Carlo method no. II

Proceeding as before, i.e. drawing n= 9! points and using the formula ( 11 ), after adjusting it to the dimension of the task: (33) (34) (35) (36) (37) where , = 1, 2, ... are the points of the rectangle on which the integration takes place, we will get the results (we have carried out 50 such experiments again) which are summarized in Tab. 10.

4.2.3. Parallel Monte Carlo method

In this method we will act in a similar way to what we did in subsection 4.1.3. Again we will draw parallel points, using k threads, each of which will draw / points. As before, we will use = 2, 3, ..., 10 threads, and the results of these experiments are presented in Tab. 11. Again, it turned out that multi-threading does not deteriorate the quality of the solution, while reducing the computation time. A graphical interpretation of the dependence of this time on the number of threads used is shown in Fig. 12.

5. Conclusions

The authors in [2] introduced a certain formula for determining the volume of solids. They also calculated the volume of the selected solids.Since the formula is non-standard (in this work it is formula ( 4 )) the authors decided to verify it. for this purpose, we used selected Monte Carlo methods that are well suited for this purpose. We discussed 5 such methods, and to verify the formula ( 4 ) we used 2 of them and an additional sixth method, which is a parallel version of the second method. The conducted research has shown that Monte Carlo methods are efective in this type of tasks, that their modifications are more useful, that formula ( 4 ) is correct and that the disadvantage of Monte Carlo methods can be efectively dealt with. The introduced parallel method does not have a negative impact on the quality, and at the same time shortens the calculation time.

References

[1] K. Czapla , M. Pleszczyński , Some formulas for determining areas and volumes, part I: areas in space 2, MINUT 3 ( 2021 ) 208 - 228 . In Polish

[2] K. Czapla , M. Pleszczyński , Some formulas for determining areas and volumes, part II: areas in space 3, MINUT 3 ( 2021 ) 229 - 244 . In Polish.

[3] N. Metropolis , S. Ulam, The monte carlo method , Journal of the American Statistical Association 44 ( 1949 ) 335 - 341 .

[4] N. P. Buslenko , D. I. Golenko , I. M. Sobol , W. G. Sragowicz , J. A. Szrejder , The Monte Carlo Method, PWN , Warsaw, 1967 . In Polish.

[5] Capizzi G. , Sciuto G.L. , Napoli C. , Shikler R. , Wozniak M. , Optimizing the organic solar cell manufacturing process by means of AFM measurements and neural networks ( 2018 ) Energies, 11 ( 5 ), doi: 10.3390/en11051221.

[6] Kloek , T. , and Van Dijk, H. K. ( 1978 ). Bayesian estimates of equation system parameters: an application of integration by Monte Carlo . Econometrica: Journal of the Econometric Society , 1 - 19 .

[7] Devroye , Luc. "Sample-based non-uniform random variate generation . " Proceedings of the 18th conference on Winter simulation . 1986 .

[8] Wróbel M. , Starczewski J.T. , Napoli C. , Handwriting recognition with extraction of letter fragments ( 2017 ) Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics) , 10246 LNAI, pp. 183 - 192 , doi: 10.1007/978-3- 319 -59060-8_ 18 .

[9] The stochastic collocation monte carlo sampler: Highly eficient sampling from 'expensive' distributions ,Quantitative Finance, 19 ( 2 ), 339 - 356 .

[10] Połap D. , Wózniak M. , Napoli C. , Tramontana E. , Is Swarm Intelligence Able to Create Mazes? (2015) International Journal of Electronics and Telecommunications , 61 ( 4 ), pp. 305 - 310 , doi: 10.1515/eletel2015- 0039 .

[11] Połap D. , Wózniak M. , Borowik G. , Napoli C. , A First Attempt to Cloud-Based User Verification in Distributed System ( 2015 ) Proceedings - 2015 AsiaPacific Conference on Computer-Aided System Engineering, APCASE 2015 , art . no. 7287024 , pp. 226 - 231 , doi: 10.1109/APCASE. 2015 . 47 .

[12] Steinbrecher , G. , and Shaw , W. T. ( 2008 ). Quantile mechanics . European journal of applied mathematics , 19 ( 2 ), 87 - 112 .

[13] Botev , Z. , and Ridder , A. ( 2017 ). Variance reduction . Wiley statsRef: Statistics reference online, 1-6.

[14] Morrison , Michael, Clif Hastings , and Kelvin Mischo . "Hands-on start to wolfram mathematica . " ( 2015 ).