Geobody Upscaling
Upscaling to Structured and Unstructured Grids
Andreas Lusth and Bjorn Lonn
Report in Scienti c Computing Advanced Course June 2010
PROJECT REPORT
Abstract
Upscaling to a grid enables the use of different sources of information in a given model. The combination of using a grid constructed from geological data and a geobody filtered from seismic data yields the grid cells of the model which represent the geobody voxels. Earlier upscaling methods took either very long time or skipped cells which should have been activated. One of the problems dealt with was finding which grid cells intersect which geobody voxels in space and determining the intersections accurately as fast as possible. The final upscaling method use a similar approach to divide and conquer, by splitting the grid into thousands of cell pillars bound by cuboid boxes. To get an accurate upscaling the method uses volume checking and several point tests. Plane intersection calculations between cells and voxels are possible but costly, in this project they have not been pursued.
Keywords: Upscaling, Geobody and Grid cell activation
Contents
1 Introduction 3
2 Background 3
3 Development 4
3.1 Petrel . . . 4
3.2 The API . . . 4
3.3 Choice of Language . . . 4
3.4 Development Platform . . . 4
3.5 Repository . . . 4
3.6 Profiling tool . . . 5
3.7 Visualizing of data . . . 5
4 Theory 5 4.1 Geobody data structure . . . 5
4.2 Grid data structure . . . 5
4.3 The connection between geobody and grid . . . 6
4.4 Different algorithm approaches . . . 6
4.4.1 Voxel inside a cell . . . 6
4.4.2 Cell inside a voxel . . . 7
4.4.3 Algorithm Complexity . . . 7
5 Upscaling 8 5.1 Methods . . . 10
5.2 Optimizations . . . 11
6 Downscaling 13 7 Results 14 7.1 Different Datamodels . . . 14
7.2 Upscaling Methods . . . 15
7.3 Upscaling Results . . . 15
7.3.1 Model 1 . . . 16
7.3.2 Model 2 . . . 17
7.3.3 Model 3 . . . 19
7.3.4 Model 4 . . . 21
8 Conclusions 22
9 Acknowledgements 23
1 Introduction
In geological research, seismic data of the Earth crust is often collected using seismic streamers and an acoustic source. Due to holes, irregularities and the fine resolution of the data received it is not suited for numerical calculations such as fluid dynamics. Geobody upscaling is a way of representing the seismic data in a predefined coarser grid making it accessible for numerical calculations.
Since the geobody often consist of a huge amount of fine elements (voxels) the upscaling process may be time consuming. Therefore, the task in this project is to upscale a general geobody as fast as possible while also making the upscaled grid cells cover the geobody as well as possible. Apart from the limitations on the upscaling methods complexity, the issue of numerically checking if a certain element of the geobody is inside a specific cell needs to be constructed to detect as many intersections as possible. Our final method, constructed for this project, produces results which are significantly better compared to the results produced by earlier methods. They were either too time consuming or skipping a lot of cells. The method handles both the case with cells much larger than geobody elements and the case with cells smaller than the geobody elements. In order to determine what case should be handled our method used a simple average volume checking to determine the volume ratio between cells and geobody elements. Depending on the case a geobody element was said to be inside a cell either by checking if any of the element’s describing points was inside all of the cell’s planes or if any of the cell’s describing points was inside a geobody element. To determine which cells should be handled by a certain geobody element we used a boundingbox concept, splitting the grid into thousands of smaller boxes, depending on data size.
2 Background
Two important information sources in geoscience are the geologist’s model repre- sented by a certain, usually unconstructed, grid and the seismic data represented in the form of a geobody. The information from the two sources is connected by performing a geobody upscaling. The upscaling produces a mix of the infor- mation from the geobody and the grid by finding which cells in the grid that should represent the geobody.
The geobody is a voxel representation of filtered seismic data. The seismic data is collected by towing seismic streamers equipped with hydrophones behind an acoustic source. The acoustic source transmits pulses of sound waves which are reflected by the sea bottom. They are then recorded by the hydrophones on the streamers. Once the signals from the hydrophones have been processed a seismic cube can be selected. Using different settings a geobody representing certain wanted properties can be extracted. These can for example be a specific material property or a certain shape in the geobody.[1]
The grid which is used in the upscaling is constructed to fit a geologist’s model of a specific area. The model can be totally independent of the seismic data. The difference between a model and a geobody is that the geobody is a representation of seismic signals while the model tries to keep information about rock type, porosity permeability etc. The grid is constructed by first selecting the upper and lower horizons. Horizons are distinct changes of material in the
crust. This is similar to the situation where the atmosphere meets the ground and therefore they are refered to as horizons. A 3D grid can be constructed by defining a certain amount of layers between the horizons and interpolating between each layer. Once we have both the geobody and the grid we can perform an upscaling, activating cells which are intersected by voxels.[2]
3 Development
3.1 Petrel
Petrel is a commercial multipurpose software developed by Schlumberger which is for instance used for seismic and reservoir modeling, seismic analysis, flow analysis and well planning. The program is versatile, it can for example interpret seismic data and perform calculations and visualize the results. The core of Petrel is written in C++, however, plugins for the software are often done in C#.
3.2 The API
We did not have access to Petrel since the Uppsala University does not have any licenses, and therefore an API was used to be able to work with data from Petrel, the interface is developed by our supervisor Bjarte Dysvik at Schlumberger.
The API is written in C# and contains several classes for easy handling of the geobody and the grid. There are two separate classes for the geobody and the grid containing all the voxels in the geobody and the cells in the grid. Those classes also include several methods for easy access and manipulation of the data. The API also contains classes for exporting and importing data for Petrel such as the grid and the geobody. Furthermore, there are sub-classes which build up the geobody- and grid classes.
3.3 Choice of Language
The code for this project is written in C#, since Petrel uses C# for plugins.
Furthermore, the tools provided were originally written in C#, and instead of porting the API to another language, it was easier to continue working with it in C#.
3.4 Development Platform
Visual Studio 2008 Professional Edition was chosen as the development environ- ment. The reason for this is that Visual Studio is a comprehensive environment and very good programming platform. Visual Studio has both a good compiler and an excellent debugger for testing. Moreover, the API is done in Visual Studio, which made it easy to continue to work and use it.
3.5 Repository
To efficiently work collaborative and always have the code up to date we used Subversion as a version tracking software. The system was provided to us by
Uppsala University. The repository was also used to exchange data with our supervisor Bjarte Dysvik.
3.6 Profiling tool
A profiling tool can be used to examine the program’s performance. A profiler can measure the number of calls to a method and the amount of times a single line is executed. The measurements give information about the program’s bottle necks. This analysis makes optimizations of the code easier and the profiling tool used was used dotTRACE Performance 4.0 from JetBrains.
3.7 Visualizing of data
The images of the geobodies and the grids in the report are produced in Petrel, however, since we did not have access to the software all images are rendered by our supervisor.
4 Theory
There is not much theory behind Geobody upscaling even though the concept has been in use for a long time. However, some knowledge of search algorithms for data structures together with information of computer architecture are very useful in making the upscaling faster while also producing a more accurate result. To make use of these theories, knowledge of the data structures used in the project is needed.
4.1 Geobody data structure
The geobody is a very sparse data structure consisting of voxels. All voxels are assigned an integer value and each value represents a material, for example sand or shale. The geobody is sparse since most voxels have the material value zero which represents a transparent object, that is a material in the geobody which the user have no interest in. It is often the case that only a few percent of the voxels represent interesting materials. Voxels are cuboids which means they all have the same measurements, but they are not necessarily cubic. This makes the geobody structured and the connection between index- and position space is straight forward. Geobodies are always have the same resolution as the seismic data which is huge, therefore, it is very time consuming to access ran- dom elements of the geobody. However, every geobody have an attribute called nonempty subcubes. They are defined by their minimum and maximum in- dices. The subcubes are constructed for better memory locality and are usually 120x120x120 voxels large. They are also nonempty which means they contain voxels with material values of interest. Since the geobody often is sparse the subcubes only cover a fraction of the geobody.
4.2 Grid data structure
The grid consists of cells which are described by eight corners in space and six planes representing the cell faces, constructed from the corners. The grid is constructed by requesting a certain number of layers between the upper and
lower horizons of geological data. The cells are then defined by interpolating between the layers. In an unstructured grid every cell extent can be unique in shape and size. However, sometimes the grid is structured in the x and y directions presenting an easy connection between position and index for those directions. All cells in the grid have a value like the voxels. This value is set to represent the most common material of the voxels that the cell intersects. A typical grid is of the order 128x86x30 up to 106x286x50 number of cells, but it can consist of smaller or larger numbers.
4.3 The connection between geobody and grid
There are numerous ways of determining the material value of a cell when upscaling. Checking through all voxels in some way is essential since we want the majority of the voxel materials in a cell to represent the cells value. To check if a voxel is intersecting a cell can be done in many ways. Note that the cells are made up from six planes, this gives an easy way to determine if any point is inside a cell. The distance between a point and a plane can be found by performing a scalar product between the normal to the plane and the vector between a point in the plane and the observed point. This works since all faces of the cells have their normals pointing out of the cell. Therefore, when a point is inside a cell it will create a vector pointing the opposite direction of the normal which results in a negative scalar product. In contrast to that, when we need to determine if a cell is inside a voxel, the check is straight forward since the voxels are cubiods.
4.4 Different algorithm approaches
The purpose of the upscaling is to map the information from the geobody to the grid as accurate as possible, and the main problem is to balance the importance of an accurate result against the runtime of the algorithm. Furthermore, there are basically two different ways of performing the upscaling. Both techniques involve looping over all elements in either the grid or the geobody and finding the correspondent element’s position in the other structure. Moreover, both are similar in structure but they produce widely different results, and both of them have their pros and cons.
4.4.1 Voxel inside a cell
The technique of looking whether a voxel is inside a cell starts out with looping over all voxels in the geobody. There can exist more than one hundred thousand voxels or even millions. However, many of the voxels are not relevant in the current geobody, they are referred to as transparent. This makes the upscaling process faster since only opaque voxels are interesting. Even though it is only necessary to loop over the voxels of the subcubes it is still a significant amount, and to reduce the amount of work one technique is to skip voxels according to some scheme. However, this can have bad effects on the results. The main problem is to find the cells with similar coordinates as a specific voxel, this because the grid is often unstructured. The worst case scenario would be to loop over the whole grid to find the cells of interest. Furthermore, we have to check if the given voxel is inside any of the cells. This can be done in several
ways, either by checking if some points in and on the voxel are inside the cell or by using intersection tests between the sides of the voxels and the planes of the cells. The first technique is much faster and less complicated, however, it does not yield as much information unless a lot of points of the voxel are tested.
Moreover, since this technique is fast, easy to test and gives good enough results it is what majority of performed testing is based on.
4.4.2 Cell inside a voxel
Instead of looping over all voxels it is possible to loop over all cells in the grid. However, it is not known if the corresponding voxels in the geobody contains any interesting material, and the only way is to check this for each cell. Furthermore, until confirmed there is no way of knowing if a certain cell is even inside the current nonempty subcube which is required, otherwise all locality is lost. However, since the geobody is structured it makes the coordinate transformation from the grid to the geobody easy. Furthermore, the test for intersection with a cell can be calculated in the same ways as earlier.
4.4.3 Algorithm Complexity
Investigating algorithm complexity is useful since it gives a good insight to how a certain change effects the performance of the algorithm. The brute force upscaling method goes through all voxels and for each cell, and checks every grid cell. The brute force method has a complexity of m· n where m is the number of voxels and n is the number of cells. Considering the complexity of the data sets this is huge. Another brute force method neglects memory locality, loops through all grid cells and for each cell it iterate through all the voxels. Since there is no easy way of telling if a certain cell is inside the space coordinates of a specific nonempty subcube the method often miss the subcubes. Accessing random elements of the geobody instead of using the subcubes is not an option since the memory overhead will be huge. Even though the latter method has the same complexity, the time consumption is 100 times more, due to poor usage of nonempty subcubes.
Sometimes the grid is structured in the x- and y-directions presenting an opportunity to replace the searches in these directions with easy calculations thus reducing the complexity to m· k where k < n, is the number of cells in the z direction. However, when the grid is unstructured the calculations of indices become impossible.
There is a more advanced way of reducing the complexity by splitting the grid into a number of boxes with given index size. The box concept is used by first checking if a box intersects the geobody subcube. Then continue with checking which voxels in the subcube that intersects the box and finally check which cells of the box are intersected by a specific voxel. The complexity of this method is the number of boxes Y times the number of subcubes V times the number of voxels in a subcube P times the number of cells in a box a· b· k. By comparing this method with its predecessors it is possible to analyze of how its variables effect the overall complexity. If a and b both are selected to have the value 1, that is all boxes represent a pillar in the grid with k cells, this will make the method very similar to the case with a structured grid. The only difference is that now we allow each voxel to be checked for more than one pillar. Increasing
the value of a and b will result in less pillars thus producing less pillar tests.
However, we have to check more cells each time when a voxel intersects a pillar.
Furthermore, since the one pillar case is hitting several pillars for each voxel, it is similar to using a nine pillar case. This means that the nine pillar case will not produce that many more intersection tests. Increasing a and b further will ultimately hit a break-even point depending on the particular geobody. The tipping point occurs when the time gain from fewer box tests equals the time loss from unnecessary cell intersection tests. If a and b are set to the extent of the grid it will result in one box with the same size as the whole grid. This forces all cells to be tested for all voxels thus the complexity is the same as the brute force method.
5 Upscaling
This section covers our insight and strategies to produce accurate upscalings while keeping the computational load as small as possible.
The most important thing regarding the connection between cells and voxels is that a voxel can intersect more than one cell. This is done by simply adding the index of an intersected cell to a list. Determining if a voxel intersects a cell can require a lot of computations. The simplest way is to use the normals of the planes of a cell to determine if a certain point is contained by the planes of the cell. However, in many cases a single point representation of a voxel is not sufficient. Luckily, to determine if a point is inside a cell is fast, if this is performed only when intersections are probable. Therefore, it is a possibility to use more points when checking for intersections.
In some cases it is not a point representation of the voxel which produces a better result but a point representation of the cell which is the key. In Figure 1 the voxel has its center in a corner of a cell. In the worst case only slightly more than 12% of the voxels volume is inside the cell, and only checking the center point will result in one upscaled cell.
Figure 1: Two cells intersected by a voxel which is only represented by its center point. The cells are blue and red and the voxel is black. The voxel center (orange dot) is inside the blue cell thus only upscaling to the blue cell even though the voxel also intersect the red cell as well as possible cells above.
In order to activate the neighboring cells which the voxel also intersects, the
strategy was to introduce more points, for example the corners of the voxel.
Note that as soon as an intersection between cell and voxel has been found the upscaling can break and continue with another cell. Figure 2 illustrates a special case of an unstructured grid, when one cell fails to be upscaled although all corner points are tested. Using more points to check for intersections will make the upscaling more accurate but will soon become time consuming as all checks are performed for cells which are not upscaled but still have the voxel within its bounding box.
Figure 2: Unstructured cells make it harder to find all intersections between cells and voxels. Even though the black voxel is represented by its center and corner points it fails to upscale the left cell since all the points lie outside the cell.
Two techniques has been used to minimize the amount of intersection tests between cells and voxels which do not lie inside the cells. When unnecessary tests are reduced the detection can be improved by using more points to de- scribe the voxels. One technique is to use bounding boxes for the cells, and the other is to create bounding boxes for groups of cells as well. The gain when using bounding boxes is that they are cuboids which makes it easy to detect intersections with voxels which are also cuboids. However, as in Figure 3 the bounding boxes of a group of neighboring cells produce only a rough filtering in the sense that when the bounding box is approximating a large group it will cover a much larger volume than the group itself. In these cases smaller groups may be considered at the cost of increasing the amount of voxel intersection checks with bounding boxes.
Figure 3: Pillars of cells from an unstructured grid have bounding boxes which cover more area than the cells in the pillar. Thus the bounding boxes overlap and it happens that a voxel is inside a box but not the cells of the box. Pillars of cells in red and blue. Voxel in black.
In all considered cases so far a voxel have been assumed to be smaller than a cell in any direction. When this assumption does not hold, i.e, the voxel is larger than the cell in some direction, the intersection tests become reversed in that direction. In such cases it is no longer useful to search for voxels intersecting cells but instead search for cells intersecting voxels as can be seen in Figure 4
Figure 4: The two red cells are thin and long compared to the black voxel.
Describing the cells with points instead of the voxel and perform inside checks for those is a much better approach
5.1 Methods
The structure and flow of three similar upscaling methods implementing differ- ent levels of complexity in the upscaling can be seen in Figure 5. The implemen- tation of a basic method is shown in black with the two other implementations branching off of it. The basic implementation which is called original method is
obtained by only following the black lines. The implementation which assumes a structured grid is obtained by prioritizing green over black. The third imple- mentation uses bounding boxes to bound pillar groups of cells and is obtained by prioritizing red over green and black.
By branching off from the original method we show the difference in flow of the methods and provide insight into the machinery of a upscaling method.
Following the black lines you see that this method loops over all voxels of a nonempty subcube and for each voxel with a positive material value, (which means it is not transparent), it loops over all cells of the grid. If prioritizing the green lines you see that instead of looping over all cells for each voxel this method now calculates the I and J index for the cells using the center point of the voxel. Since the I and J indices have been calculated we only need to loop over the K index of the cells thus this is a one pillar method.
A method, similar to our final version for unstructured grids is obtained by following the priority order red, green , black. When upscaling to an unstruc- tured grid the method first create bounding boxes for groups of cells then it continues by looping over the voxels just like the other methods. For each voxel the method perform an intersection test with the bounding boxes. Furthermore, for all cells of the bounding boxes which are intersected by a voxel, the method continues with further intersection tests now between cell and voxel. Once a cell has been found intersected by a voxel the material value represented by the voxel is added to the cell and the method can start over with another cell or voxel.
5.2 Optimizations
Apart from regular code optimizations such as minimizing method calls, moving unnecessary calculations out of loop structure and memory locality optimiza- tions, the profiling of the code pointed out that some parts of the code were doing a lot of unnecessary work.
From studying cell construction we found that since the cells are defined by their eight corners, each time an intersection test must be performed the planes representing the faces of the a cell must first be constructed. The plane construction itself is in fact relatively fast, however, if done every time when intersection tests are encountered it will be done so often that it is one of the most time consuming parts of the method. The solution which reduces plane creation to an absolute minimum came out of modifications of the cell construc- tion process. When investigating intersections the method first investigates if the voxel intersect the bounding box of a cell, therefore, in this level of inves- tigations the planes are not needed. By separating the plane creation from the cell construction and inserting a variable which tracks if the planes have been created for the cell or not, the method can create the planes precisely when needed. Furthermore, the cells are now stored in the grid data structure instead of deleted, which enables the method to reuse cells which already have been created, thus saving run time while the memory consumption is increased.
When searching through a pillar of cells, that is keeping the I and J indices fixed as the K index is changed, the method often encounters a cell where no intersection is found while the previous cell was upscaled. It is then safe to assume that for this specific voxel, no other cells in the same pillar can be upscaled. Thus the method can instead continue with another pillar. A small
Figure 5: Flowchart of a brute force upscaling method in black. Branching of the brute force method is, in green the method which assumes a structured grid and in red a method using bounding boxes to reduce the amount of tests for unstructured grids
change concerning when the final value of a cell is set, enabled the method to only set the value once, instead of every time the cell is found to be intersected by a voxel. This was done by only recording the index of newly created cells which are intersected by a voxel. However, this change gave no significant speed up but was worth the effort since it makes the indices of upscaled cells easily accessible.
Several optimizations were also done to the initial API. One of the most im- portant changes was in the geobody class, in order to to simplify the calculation of the center point of a voxel. The center point is calculated from the index and dimensions of the voxel. When this is done, the point has to be rotated into the correct position according to the rotation of the geobody. However, the rotation is not needed unless the geobody itself is rotated. By using this insight one transformation, which is a 4x4 matrix multiplication can be saved for each voxel lookup. An other improvement of the API was the distance calculation between a point and a plane which is used when checking if a voxel is inside the cell.
6 Downscaling
In some rare occasions, the model consists of a grid where the cells are smaller than the voxels. This yields problems when checking if a voxel is inside a cell because the voxel will span more cells. Therefore, the voxels have to be checked against more cells. Furthermore, to effectively hit all cells the voxel has to be represented by more points. One way of dealing with the problem is to check if the cells are inside the voxels, instead of the other way around. However, the method has to effectively be able to upscale and downscale a model without the user knowing about it. Both the structured and the unstructured methods do this by calculating the volume of the voxels and the cells, and based upon the result, then decide which method to use. If the cells volume are larger then an upscale is done, however, if the voxels are larger then a downscaling is performed.
It is trivial to calculate the volume of the voxels, however, calculating the volume of the cells is harder since they may have different sizes and shapes. To be able to calculate the volume the method assumes the cells are parallelepipeds, see Figure 6 and calculates the volume by |a·(b×c)|. To get a representative volume value for the cells, an averaged value of 100 cells is used.
Figure 6: A cell in the form of a parallelepiped, the opposite sides are parallel.
7 Results
7.1 Different Datamodels
The algorithms has been tested on four different data models with a geobody and a grid. Two models contains the same geobody but with different material properties and different grids. The geobody of the first model can be seen in Figure 7a and it contains about 4 million voxels, but only 7343 of them are not transparent. Furthermore, the geobody also contains only two material properties out of the 256 which geobodies can hold. The grid for model 1 is structured in the x and y directions. Since there are such a low amount of voxels representing a material the runtimes were much faster than many of the others.
The second model’s geobody has 1524755 voxels with material values, the grid is also bigger. The third model is a more realistic model in the sense that the grid is unstructured, the geobody contains about half the amount of voxels compared to the previous geobodies, but approximately the same amount of voxels containing materials. Moreover, the grid is completely unstructured, as seen in the pictures below. The last model is similar to the second one, however, the cells of the grid are smaller than the voxels. Complete information about all the models can be found in Table 1.
Data models Number of voxels Number of voxels with materials
Number of nonempty subcubes
Number of cells
Model 1 40085010 7343 2 401100
Model 2 40085010 1524755 32 101760
Model 3 40085010 7343 22 22440
Model 4 25622625 1607294 2 1496250
Table 1: The table contains numbers of the total voxels, the voxels with material properties, how many nonempty subcubes these are spread out on as well as the number of cells in the grid.
(a) Model 1 (b) Model 2
(c) Model 3 (d) Model 4
Figure 7: The geobodies for the different models, the colors represent different materials
7.2 Upscaling Methods
These were described in 5.1 and the upscalings have been done with the original method, the structured method and the unstructured method. All these have been tested on all data models and they produce different results. Important to note is that the original method only can handle one material property and it needs a grid that is structured in two directions, hence that method does less computations in this regard compared to the other methods.
7.3 Upscaling Results
In the following section we will present the runtimes and the upscaled grids for the different models. Each model will be presented in its own section. The tests are done on a HP Compaq dc7900 with a Intel CoreR TM2 Quad Processor Q9400 at a clock frequency of 2.66 GHz and a 6 MB L2 cache. Furthermore, the computer is equipped with 4096 MB of RAM at 800 MHz, and the operating system is Windows 7. The execution times presented are an average value of five runs.
7.3.1 Model 1
The geobody of model 1 can be seen in Figure 7a, and the different runtimes with the number of activated cells in each grid for the methods are presented in table 2.
The result of the original method can be seen in Figure 8. In the figure a lot of voxels along the outer edges haven’t been upscaled. This means that this result can be improved on several places. The structured method can be seen in Figure 9, and from the picture it is clear that the result has improved, the complete geobody has been covered. Furthermore, the execution time has also decreased much, compared to the original method, with some optimizing discussed previously in 5.2 together with some overall improvements to the API.
The result of the unstructured method in Figure 10 is the same result as for the structured method. However, the execution time was slower since the the method does not utilize the structured properties in the grid.
Method Cells upscaled Execution time (ms)
Original 1150 627
Structured 1534 919
Unstructured 1534 1484
Table 2: Execution time and the number of activated cells for the different methods on Model 1.
Figure 8: The upscaled grid produced by the original method on model 1. The geobody is purple and the grid has a green color. Note that there is only one material in the figure.
Figure 9: The upscaled grid produced by the structured method on model 1.
The geobody is purple and the grid is turquoise and green. Note that the complete geobody has been covered by the grid.
Figure 10: The upscaled grid produced by the unstructured method on model 1. The geobody is purple and the grid is turquoise and green. Note that the complete geobody has been covered by the grid.
7.3.2 Model 2
The geobody of model 2 can be seen in Figure 7b, and the same methods as before were tested. The results of the different methods are found in Table 3 and Figures 11-13. A quick glance at all pictures shows that there is not much difference between the methods in the number of upscaled cells, the unstructured method upscales 7% more voxels than the original method. However, Table 3 reveals that the structured methods and the unstructured methods are a lot faster than the original. The structured method can be up to ten times faster on this grid.
Method Cells upscaled Execution time (ms)
Original 38442 68814
Structured with 1 pillar 40883 12904
Structured with 9 pillars 41300 118873
Unstructured 41300 59143
Table 3: Execution time and the number of activated cells for the different methods on Model 2.
Figure 11: The upscaled grid produced by the original method on model 2. Note that there is only one material in the figure.
(a) 1 Pillar (b) 9 Pillars
Figure 12: The upscaled grid produced by the structured method on model 2.
Figure 13: The upscaled grid produced by the unstructured method on model 2.
7.3.3 Model 3
The geobody of model 3 can be seen in Figure 7c and the results of the different methods are found in Table 4 and Figures 14-16. For this model there is a clear difference between the methods and a quick look at Figures 14-15 shows that the structured methods perform very poorly on a unstructured grid. It is only the unstructured method, Figure 16 which produces a nice result, there are no visible voxels which have not been upscaled. Furthermore, Table 4 reveals that the original method is the slowest, and that looping over more pillars costs a lot of time. Moreover, the unstructured method is almost as fast as structured looping over nine pillars.
Method Cells upscaled Execution time (ms)
Original 1359 23467
Structured with 1 pillar 493 5363
Structured with 9 pillars 1414 22782
Unstructured 6534 25537
Table 4: Execution time and the number of activated cells for the different methods on Model 3.
Figure 14: The upscaled grid produced by the original method on model 3. Note that there is only one material in the figure.
(a) 1 Pillar (b) 9 Pillars
Figure 15: The upscaled grid produced by the structured method on model 3.
Figure 16: The upscaled grid produced by the unstructured method on model 3.
7.3.4 Model 4
The geobody of model 4 can be seen in Figure 7d and the results of the different methods are found in Table 5 and Figures 17-19. The original method preforms poorly when the cells are smaller than the voxels, and Figure 17 shows that only a handful of cells has been upscaled. Furthermore, there is a large difference in the results for the structured method depending on how many pillars are chosen. A test run with 9 pillars upscales more than 6 times the cells upscaled by a 1 pillar attempt. However, the unstructured method seen in Figure 19 upscales the most cells, and it produces results which are 30 times better than the original method. A look at Table 5 shows that the structured method with 9 pillars and the unstructured method are considerably slower than the other two. However, the other two produce extremely poor results.
Method Cells upscaled Execution time (ms)
Original 1228 977
Structured with 1 pillar 6119 445
Structured with 9 pillars 39314 1630
Unstructured 40689 4922
Table 5: Execution time and the number of activated cells for the different methods on Model 4.
Figure 17: The upscaled grid produced by the original method on model 2. The geobody is purple and the grid has a green color. Note that there is only one material in the figure.
(a) 1 Pillar (b) 9 Pillars
Figure 18: The upscaled grid produced by the structured method on model 4.
Figure 19: The upscaled grid produced by the unstructured method on model 4.
8 Conclusions
Analyzing the results of the performance of the considered methods on the benchmark problems, we can conlude the following. It is clearly seen that the original method produces a worse upscaling compared to the structured method.
The only difference effecting the accuracy between these methods when the cells are smaller than the voxels is that the structured method uses a total of 35 points to describe a voxel compared to the original method which use only the center point. Furthermore, even though the structured method uses 35 times more intersection tests to detect intersections it is 6 times faster than the original method. This is the effect of using bounding boxes for the cells to rule out cells not in the vicinity of the voxel together with the optimizations described in section 5.2. In most cases the 9 pillar version of the structured method proved not to be worth the longer run times given the small increase in upscaled cells. However, it was crucial in the case when the voxels are lager than the cells. There, six times more cells were upscaled while the time only increased by a factor of four. The idea behind the 9-pillar version is that voxels may be shared among the pillars of cells surrounding the pillar which indices was calculated. However, for a rather homogeneous geobody such as model 1, the effect is very low. It is also instructive to see how the method scales from model 1 to model 2, the amount of voxels between the models are the same, but the second model contains 350 times more voxels with intresting material properties. The runtime for the original method increase with a factor of 100 while the execution time for the structured method over one pillar only increases 50 times. The unstructured method runs 40 times slower than on model 2. The structured- and unstructured methods scale better because they are designed to minimize the number of intersection tests.
Upscaling to unstructured grids with methods which assume that the grid is structured fail as expected. However it is interesting to note the very low amount of upscaled cells produced by the structured method on model 3 compared to the original method. This was at first not expected since they are based on the same idea and the structured method over one pillar is always at least suppose to reproduce the result of the original method. On model 3 the unstructured method really shines, it is the only method which produces a good result and it even maintains a low run time. Even though the grid is unstructured, the relations for the run time between the structured and unstructured methods are the same.
Upscaling a model where the cells are smaller than the voxels is a very similar procedure, the issue is that, since the cells are much smaller then the voxels, it is hard to hit them accurately with the voxel point tests used in all other methods. To solve this, it is natural to reverse the way of thinking, that is describe the cells with points instead and check if any of those points are inside the voxel. Checking if a point is inside a voxel is straight forward since the voxels are cuboids the method can check if the point is bounded by the voxels min and max coordinates. Determining if this or the previous methods should be used can be done by comparing the mean volume of a random selection of cells, with the voxel volume.
Even though the structured method and unstructured method produce a better upscaling while preserving a low runtime there are still some issues to be resolved. The main problem is that a cell is always upscaled by any voxel which is found to intersect the cell with no consideration of its volume inside the cell.
This may in some cases lead to more upscaled cells than necessary, and it may also give some cells a undesired material type. The best way to handle this is to implement some volume checking where the user can determine how much of the voxel has to lie inside the cell. However, since volume checking requires more advanced calculations, it slows down the performance of the upscaling procedure.
9 Acknowledgements
The authors would like to thank Bjarte Dysvik from Schlumberger for his active participation and Maya Neytcheva the course coordinator.
References
[1] Munday Brink and Jan Gateman, Method for simultaneous collec- tion of seismic data from shallow and deep targets. Patent 5148406, http://www.freepatentsonline.com/5148406.html, 1992.
[2] Alexis Carrilant and Brice Vall`es, From 3D Seismic Facies to Reservoir Simulation: An Example From the Grane Field. Mathematical Methods and Modelling in Hydrocarbon Exploration and Production, Volume 7, 301-338, Springer Berlin Heidelberg, 2005.
