<p>You can partition your domain for subsequent parallel processing as follows. There are a few possible options. You can try a static geometric decomposition, including recursive bisection, quad or oct-trees and a space-filling curve.</p>

<p>Recursive bisection techniques are used to partition a domain (e.g., a finite element grid) into subdomains of approximately equal computational cost while attempting to minimize communication costs, that is, the number of channels crossing task boundaries
A divide and conquer approach is taken. The domain is first cut in one dimension to yield two subdomains. Cuts are then made recursively in the new subdomains until we have as many subdomains as we require tasks. Notice that this recursive strategy allows the partitioning algorithm itself to be executed in parallel. Advanced variants include recursive graph bisection and recursive spectral bisection.</p>

<p>Here is an example. In the following picture, assume that magenta points requires twice as much work as cyan ones.</p>

<p><img src="https://i.sstatic.net/GVhjR.png" alt="enter image description here"></p>

<p>Then, you continue applying recursive bisection by cycling through the axes until # pieces = # processors, obtaining something like the following:</p>

<p><img src="https://i.sstatic.net/cZddS.png" alt="enter image description here"></p>

<p>The idea behind quad-trees is the following one. Partition all of space into quadrants (octants in 3-space). If any piece contains too much work, partition it into quadrants, and so on. See the next picture for an example.</p>

<p><img src="https://i.sstatic.net/DhRPy.png" alt="enter image description here"></p>

<p>The best general-purpose geometric load-balancing comes from space-ﬁlling curves.
The order in which points are visited in the space-ﬁlling curve determines how the geometric objects are grouped together to be assigned to the processors. The following image shows an Hilbert space-filling curve.</p>

<p><img src="https://i.sstatic.net/QfrdS.png" alt="enter image description here"></p>

<p>You problem may also be amenable to a solution using what is known as a scattered decomposition. Used when there is a structured domain space (e.g., an image) and the processing requirements are clustered, such as modeling a crash or processing an image with only a few items of interest. Suppose there are P processors. Cover the problem domain with non-overlapping copies of a grid of size P and assign each processor a cell in each of the grids. If you use more pieces:</p>

<ul>
<li>Less load imbalance, i.e., computation time </li>
<li>More overhead and/or communication time</li>
</ul>

<p>Finding a good tradeoff may require experimenting running the application and measuring how parallel time is affected by different scattering strategies.</p>

<p>Here is an example. The first picture show an example image that must be processed. Say it contains sea water and an island, and you may want to process the vegetation on the island. The image is partitioned using a very simple domain decomposition technique. The second picture shows an example of scattered decomposition.</p>

<p><img src="https://i.sstatic.net/LOXDN.png" alt="enter image description here"></p>

<p><img src="https://i.sstatic.net/Ca0rZ.png" alt="enter image description here"></p>

<p>Finally, you can use an advanced technique known as Adaptive Mesh refinement (AMR): grids are broken into blocks of ﬁxed extents, when needed blocks are reﬁned into children with same extents. Block size are chosen to optimize balance of cache efﬁciency (optimization possible since loops operate on ﬁxed size array), communication, and load balance.
Often balancing blocks per processor sufﬁces; if communication is excessive then use a space-ﬁlling curve. Rebalancing requires only simple collective communication operations to decide where blocks go. The technique can be used for arbitrary dimensions. If block size is 1, adaptive blocks are quad or oct-trees.
See the next picture for an example.</p>

<p><img src="https://i.sstatic.net/bYNZA.png" alt="enter image description here"></p>

<p>Hope this helps.</p>

CC BY-SA 3.0
2011-10-25T08:11:23.880Z
thanks a lot for your really impressive and exhaustive answer
Open the way
CC BY-SA 3.0
2011-10-25T08:12:58.337Z
you gave me a lot of new ideas I did not know the exact names for which, I will have a look at them in detail. Last thing, would you recommend any GPU implementation (if any) for AMR and quad tree decomposition?
Open the way
CC BY-SA 3.0
2011-10-25T11:32:44.670Z
Unfortunately, I am not aware of public domain, open source GPU implementations, unless recently released. I know for sure that someone implemented AMR and Quad-trees on GPUs, since I remember reading a few papers about.
Massimo Cafaro
CC BY-SA 3.0
2011-10-25T13:36:06.270Z
ok thanks a lot; you deserve 1000 points for the answer! last thing, which tools did you use for the figures?
Open the way
CC BY-SA 3.0
2011-10-25T13:46:03.647Z
The figures are not mine; I use them in my course on Parallel Algorithms in order to explain possible options for domain decompositions.
Massimo Cafaro
CC BY-SA 3.0
2011-10-24T19:21:56.653Z
2011-10-24T19:35:09.763Z
2011-10-24T19:35:09.763Z
Massimo Cafaro
<p>I posted already <a href="https://stackoverflow.com/questions/7805260/best-parallel-method-for-calculating-the-integral-of-a-2d-function">this question</a> to SO but got no answer so I try it now here:</p>

<p>In some crunching number program, I have a function which can be just 1 or 0 in three dimensions. I do not know in advance the function, but I need to know the total "surface" of the function which is equal to zero. In a similar problem I could draw a rectangle over the 2D representation of the map of United Kingdom. The function is equal to 0 at sea, and 1 at the earth. I need to know the total water surface. I wonder what is the best parallel algorithm or method for doing this.</p>

<p>I thought first about the following approach; a) divide 2D map area into a rectangular grid. For each point that belongs to the center of each cell, check whether it is earth of water. This can be done in parallel. At the end of the procedure I will have a matrix with ones and zeroes. I will get the area with some precision. Now I want to increase this precision, so b) choose the cells that are in the border regions between zeroes and ones (what is the best criterion for doing this?) and in those cells, divide them again into successive cells and repeat the process until one gets the desired accuracy. I guess that in this process, the critical parameters are the grid size for each new stage, and how to store and check the cells that belong to the border area. Finally the most optimal method, from the computational point of view, is the one that performs the minimal number of checks in order to get the value of the total surface with the desired accuracy.</p>

CC BY-SA 3.0
2011-10-24T07:26:54.363Z
yes, some heuristic would be great
Open the way
CC BY-SA 3.0
2011-10-24T07:27:53.940Z
what do you exactly mean about the boundary? the function to integrate or the country to calculate the area to is different on each different run of the program, if that is what you mean
Open the way
<dc.parallel-comp><na.numerical-analysis><integer-lattice>
Most optimal parallel method for calculating the integral of a 2D function
Open the way
# Most optimal parallel method for calculating the integral of a 2D function

(Asked by: username_0 on 2011-10-23)

I posted already [this question](https://stackoverflow.com/questions/7805260/best-parallel-method-for-calculating-the-integral-of-a-2d-function) to SO but got no answer so I try it now here:

In some crunching number program, I have a function which can be just 1 or 0 in three dimensions. I do not know in advance the function, but I need to know the total "surface" of the function which is equal to zero. In a similar problem I could draw a rectangle over the 2D representation of the map of United Kingdom. The function is equal to 0 at sea, and 1 at the earth. I need to know the total water surface. I wonder what is the best parallel algorithm or method for doing this.

I thought first about the following approach; a) divide 2D map area into a rectangular grid. For each point that belongs to the center of each cell, check whether it is earth of water. This can be done in parallel. At the end of the procedure I will have a matrix with ones and zeroes. I will get the area with some precision. Now I want to increase this precision, so b) choose the cells that are in the border regions between zeroes and ones (what is the best criterion for doing this?) and in those cells, divide them again into successive cells and repeat the process until one gets the desired accuracy. I guess that in this process, the critical parameters are the grid size for each new stage, and how to store and check the cells that belong to the border area. Finally the most optimal method, from the computational point of view, is the one that performs the minimal number of checks in order to get the value of the total surface with the desired accuracy.

**Tags:** dc.parallel-comp, na.numerical-analysis, integer-lattice
> **Comments:**
> **username_0**: yes, some heuristic would be great  
> **username_0**: what do you exactly mean about the boundary? the function to integrate or the country to calculate the area to is different on each different run of the program, if that is what you mean  

---

## 1 answers

### Answer by username_1 on 2011-10-24. Score: 4 | Accepted answer

You can partition your domain for subsequent parallel processing as follows. There are a few possible options. You can try a static geometric decomposition, including recursive bisection, quad or oct-trees and a space-filling curve.

Recursive bisection techniques are used to partition a domain (e.g., a finite element grid) into subdomains of approximately equal computational cost while attempting to minimize communication costs, that is, the number of channels crossing task boundaries
A divide and conquer approach is taken. The domain is first cut in one dimension to yield two subdomains. Cuts are then made recursively in the new subdomains until we have as many subdomains as we require tasks. Notice that this recursive strategy allows the partitioning algorithm itself to be executed in parallel. Advanced variants include recursive graph bisection and recursive spectral bisection.

Here is an example. In the following picture, assume that magenta points requires twice as much work as cyan ones.

![enter image description here](https://i.sstatic.net/GVhjR.png)

Then, you continue applying recursive bisection by cycling through the axes until # pieces = # processors, obtaining something like the following:

![enter image description here](https://i.sstatic.net/cZddS.png)

The idea behind quad-trees is the following one. Partition all of space into quadrants (octants in 3-space). If any piece contains too much work, partition it into quadrants, and so on. See the next picture for an example.

![enter image description here](https://i.sstatic.net/DhRPy.png)

The best general-purpose geometric load-balancing comes from space-ﬁlling curves.
The order in which points are visited in the space-ﬁlling curve determines how the geometric objects are grouped together to be assigned to the processors. The following image shows an Hilbert space-filling curve.

![enter image description here](https://i.sstatic.net/QfrdS.png)

You problem may also be amenable to a solution using what is known as a scattered decomposition. Used when there is a structured domain space (e.g., an image) and the processing requirements are clustered, such as modeling a crash or processing an image with only a few items of interest. Suppose there are P processors. Cover the problem domain with non-overlapping copies of a grid of size P and assign each processor a cell in each of the grids. If you use more pieces:

* Less load imbalance, i.e., computation time
* More overhead and/or communication time

Finding a good tradeoff may require experimenting running the application and measuring how parallel time is affected by different scattering strategies.

Here is an example. The first picture show an example image that must be processed. Say it contains sea water and an island, and you may want to process the vegetation on the island. The image is partitioned using a very simple domain decomposition technique. The second picture shows an example of scattered decomposition.

![enter image description here](https://i.sstatic.net/LOXDN.png)

![enter image description here](https://i.sstatic.net/Ca0rZ.png)

Finally, you can use an advanced technique known as Adaptive Mesh refinement (AMR): grids are broken into blocks of ﬁxed extents, when needed blocks are reﬁned into children with same extents. Block size are chosen to optimize balance of cache efﬁciency (optimization possible since loops operate on ﬁxed size array), communication, and load balance.
Often balancing blocks per processor sufﬁces; if communication is excessive then use a space-ﬁlling curve. Rebalancing requires only simple collective communication operations to decide where blocks go. The technique can be used for arbitrary dimensions. If block size is 1, adaptive blocks are quad or oct-trees.
See the next picture for an example.

![enter image description here](https://i.sstatic.net/bYNZA.png)

Hope this helps.

> **Comments:**
> **username_0**: thanks a lot for your really impressive and exhaustive answer  
> **username_0**: you gave me a lot of new ideas I did not know the exact names for which, I will have a look at them in detail. Last thing, would you recommend any GPU implementation (if any) for AMR and quad tree decomposition?  
> **username_1**: Unfortunately, I am not aware of public domain, open source GPU implementations, unless recently released. I know for sure that someone implemented AMR and Quad-trees on GPUs, since I remember reading a few papers about.  
> **username_0**: ok thanks a lot; you deserve 1000 points for the answer! last thing, which tools did you use for the figures?  
> **username_1**: The figures are not mine; I use them in my course on Parallel Algorithms in order to explain possible options for domain decompositions.  

---

