Open
Description
When creating a TreeMesh, if we access the cell_centers attribute before we finalize the mesh we end up getting an invalid mesh with an empty array for the cell_centers after we do finalize it.
The following example works as expected and a valid mesh is created:
import numpy as np
import discretize
# Minimum cell widths
dx, dy = 5, 5
# Mesh lengths along each direction
x_length, y_length = 300.0, 300.0
# Compute number of base mesh cells required in x and y
nbcx = 2 ** int(np.round(np.log2(x_length / dx)))
nbcy = 2 ** int(np.round(np.log2(y_length / dy)))
# Define the base mesh
hx = [(dx, nbcx)]
hy = [(dy, nbcy)]
mesh = discretize.TreeMesh([hx, hy], x0="CC")
# Refine near points
xx = np.array([0.0, 10.0, 0.0, -10.0])
yy = np.array([-20.0, -10.0, 0.0, -10])
pts = np.c_[discretize.utils.mkvc(xx), discretize.utils.mkvc(yy)]
mesh = discretize.utils.refine_tree_xyz(
mesh, pts, octree_levels=[2, 2], method="radial", finalize=False
)
mesh.finalize()
print(mesh)
print("\nSize of cell_centers:", mesh.cell_centers.size)QuadTreeMesh: 4.64% filled
Level : Number of cells Mesh Extent Cell Widths
----------------------- min , max min , max
2 : 4 --------------------------- --------------------
3 : 34 x: -160.0 , 160.0 5.0 , 80.0
4 : 36 y: -160.0 , 160.0 5.0 , 80.0
5 : 68
6 : 48
-----------------------
Total : 190
Size of cell_centers: 380
But if we access mesh.cell_centers before running mesh.finalize(), the mesh.cell_centers array is empty, making impossible to use this mesh in any further step (inversion, plotting, etc).
import numpy as np
import discretize
# Minimum cell widths
dx, dy = 5, 5
# Mesh lengths along each direction
x_length, y_length = 300.0, 300.0
# Compute number of base mesh cells required in x and y
nbcx = 2 ** int(np.round(np.log2(x_length / dx)))
nbcy = 2 ** int(np.round(np.log2(y_length / dy)))
# Define the base mesh
hx = [(dx, nbcx)]
hy = [(dy, nbcy)]
mesh = discretize.TreeMesh([hx, hy], x0="CC")
# Refine near points
xx = np.array([0.0, 10.0, 0.0, -10.0])
yy = np.array([-20.0, -10.0, 0.0, -10])
pts = np.c_[discretize.utils.mkvc(xx), discretize.utils.mkvc(yy)]
mesh = discretize.utils.refine_tree_xyz(
mesh, pts, octree_levels=[2, 2], method="radial", finalize=False
)
print(mesh.cell_centers)
mesh.finalize()
print(mesh)
print("\nSize of cell_centers:", mesh.cell_centers.size)[]
QuadTreeMesh: 4.64% filled
Level : Number of cells Mesh Extent Cell Widths
----------------------- min , max min , max
2 : 4 --------------------------- --------------------
3 : 34 x: -160.0 , 160.0 5.0 , 80.0
4 : 36 y: -160.0 , 160.0 5.0 , 80.0
5 : 68
6 : 48
-----------------------
Total : 190
Size of cell_centers: 0
cc @lheagy
Metadata
Metadata
Assignees
Labels
No labels