# Out of memory error in \`viz\` of 3D \`RectilinearGrid\` from Meshes.jl

**URL:** <https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643>\
**Category:** Visualization\
**Tags:** hdf5, memory, glmakie, meshes\
**Created:** [May 2, 2025, 9:46pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643 "2025-05-02T21:46:04Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 2, 2025, 9:46pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/1 "2025-05-02T21:46:04Z")

</div>

I want to plot a grid using Meshes.jl and GLMakie.jl. I am reading the HDF5 data file output of [Athena++](https://github.com/PrincetonUniversity/athena/wiki/HDF5-Format) astrophysical code using HDF5.jl library. In Athena++ mesh grid is divided into smaller parts called _Block_s to use block level parallelization. In code below i create x ,y ad z direction coordinates points for 872 blocks using rand in matrices. I got _out of memory_ error. How can i plot the grid containing x, y and z coordinate points for every block?

```julia
using HDF5, Meshes, GLMakie
x = rand(Float32, 18,872)
y = rand(Float32, 4,872)
z = rand(Float32, 16,872)

grid = RectilinearGrid(vec(x), vec(y), vec(z))
viz(grid, alpha = 0.3)

```

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 2, 2025, 9:56pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/2 "2025-05-02T21:56:16Z")

</div>

The `RectilinearGrid` expects `x`, `y` and `z` vectors. The resulting grid is the Cartesian product of these vectors. In your example, you are defining matrices instead. Are you sure you don’t want `StructuredGrid`?

Please check the documentation:

> **[Meshes · Meshes.jl](https://juliageometry.github.io/MeshesDocs/stable/domains/meshes/#Meshes.RectilinearGrid)**
>
> Documentation for Meshes.jl.

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 2, 2025, 10:01pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/3 "2025-05-02T22:01:34Z")

</div>

> [@juliohm](#):
>
> In your example, you are defining matrices instead.

I am getting such matrices when i read Athena++ output HDF5 file. In those matrices 872 represents the number of blocks that have been used to create a grid during simulation.

 ![image](https://global.discourse-cdn.com/julialang/original/3X/3/8/38a13504f45e73be2a2a0335166bd2364e29b855.png)  
 ![image](https://global.discourse-cdn.com/julialang/original/3X/a/d/ad0244979dcaeb3d558f16321b94b1083c67937a.jpeg)

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 2, 2025, 10:13pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/4 "2025-05-02T22:13:20Z")

</div>

I have no idea what Athena++ is or does. You have to manipulate the data to feed the expected inputs to Meshes.jl

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 3, 2025, 8:56am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/5 "2025-05-03T08:56:27Z")

</div>

> [@juliohm](#):
>
> The `RectilinearGrid` expects `x`, `y` and `z` vectors.

I am using vec() to convert matrix to vector.

> [@juliohm](#):
>
> Are you sure you don’t want `StructuredGrid`?

I want to create a grid of those points for all blocks. I don’t mind whether it is ‘StructuredGrid’ or not.

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 3, 2025, 9:16am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/6 "2025-05-03T09:16:30Z")

</div>

The grid is successfully created. The out-of-memory error is thrown in the `viz` call.

Visualization of `RectilinearGrid` is optimized in the 2D case. If you can plot the grid with Makie.jl directly, please share the code used, and we can add the optimization.

Currently, the visualization of 3D grid with arbitrary spacing is using the fallback visualization recipe for grids, which is memory-hungry.

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 3, 2025, 10:18am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/7 "2025-05-03T10:18:07Z")

</div>

This code creates plot on using only some data points in GLMakie.jl.

```julia
using HDF5, Meshes, GLMakie

# Load Athena++ grid from file
filename = "athena_output.00000.h5"
h5 = h5open("sane00.prim.01800.athdf", "r")

@timev x1v = read(h5["x1v"]) # cell centers in x
x2v = read(h5["x2v"]) # cell centers in y
x3v = read(h5["x3v"]) # cell centers in z

close(h5)
@show size(x2v)
# Convert to RectilinearGrid (if the grid is structured and aligned)
#grid = RectilinearGrid(vec(x1v), vec(x2v), vec(x3v))
#return viz(grid, alpha = 0.3)

# Use stride to subsample (reduce memory usage)
stride = 100
x_sample = x1v[1:stride:end]
y_sample = x2v[1:stride:end]
z_sample = x3v[1:stride:end]

# Create Point3f0 objects
points = [GLMakie.Point(x, y, z) for x in x_sample, y in y_sample, z in z_sample]
flat_points = vec(points) # convert 3D array to 1D vector

# Plotting
fig = Figure(size = (800, 600))
ax = Axis3(fig[1, 1], title = "Grid Points", xlabel = "x", ylabel = "y", zlabel = "z")
scatter!(ax, flat_points, color = :blue, markersize = 9.5)
fig

```

If i use `stride` = 1 then also i get _out of memory_ error. My Athena++ data is at this [link](https://users.camk.edu.pl/miki/sane00.prim.01800.athdf). Please explain me solution in detail. 💬 Over this grid later i will plot density or Magnetic fields etc.

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 3, 2025, 10:30am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/8 "2025-05-03T10:30:44Z")

</div>

The problem is that you are instantiating the Cartesian product of x y z. You need to find the function in Makie.jl that plots the corresponding grid directiy without manual construction of points.

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 3, 2025, 11:26am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/9 "2025-05-03T11:26:59Z")

</div>

I am able to plot grid of single `block` using Meshes.jl. Can i merge all these block’s grids to get final grid?

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 3, 2025, 11:40am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/10 "2025-05-03T11:40:31Z")

</div>

Plotting individual blocks one by one will take an enormous amount of time with multiple calls to the GPU. I repeat: find the Makie.jl function that can plot the grid directly with the x, y and z vectors. If this function exists, we can call it in our `viz` recipes.

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 3, 2025, 11:49am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/11 "2025-05-03T11:49:07Z")

</div>

This code gives correct scatter plot but i want 3D volume plot. My grid is non-uniform so i have trouble getting volume plot. Please guide me how to get volume plot?

```julia
using HDF5, GLMakie, Colors, LinearAlgebra
GLMakie.activate!(ssao=true)
GLMakie.closeall() # close any open screen

function plothdf(h5file::String)

    h5 = h5open(h5file, "r")
    num_blocks = read(HDF5.attributes(h5), "NumMeshBlocks")
    ρ = read(h5["prim"])[:,:,:,:,1] # shape: (Nx, Ny, Nz, Nvars, Nblocks)
    x_all = read(h5["x1v"]) # shape: (Nx, Nblocks)
    y_all = read(h5["x2v"])
    z_all = read(h5["x3v"])
    close(h5)

    all_xs, all_ys, all_zs, all_density = Float32[], Float32[], Float32[], Float32[]

    for b in 1:num_blocks
        x = x_all[:, b]
        y = y_all[:, b]
        z = z_all[:, b]
        density = ρ[:,:,:,b] # 3D array

        for i in 1:length(x), j in 1:length(y), k in 1:length(z)
            r, θ, ϕ = x[i], y[j], z[k]
            X = r * sin(θ) * cos(ϕ)
            Y = r * sin(θ) * sin(ϕ)
            Z = r * cos(θ)

            push!(all_xs, X)
            push!(all_ys, Y)
            push!(all_zs, Z)
            push!(all_density, density[i, j, k])
        end
    end

    # Normalize density
    norm_density = log10.(abs.(normalize!(all_density)))
    cmap = cgrad(:plasma)
    colors = get.(Ref(cmap), norm_density)
    fig = Figure(size = (1000, 800))
    ssao = Makie.SSAO(radius = 5.0, blur = 3)
    ax = LScene(fig[1, 1], scenekw = (ssao=ssao,))
    ax.scene.ssao.bias[] = 1.95
    plt = scatter!(ax, all_xs, all_ys, all_zs, color = norm_density, markersize = 5)
    Colorbar(fig[1, 2], plt, label = "log₁₀(norm(ρ))")
    display(fig)
end
plothdf("sane00.prim.01800.athdf")

```

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 24, 2025, 9:58pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/12 "2025-05-24T21:58:55Z")

</div>

> [@juliohm](#):
>
> find the [Makie.jl](https://juliaregistries.github.io/General/packages/redirect_to_repo/Makie) function that can plot the grid directly with the x, y and z vectors.

In `Makie.jl`, the `volume` function does **not** natively support non-uniform grids via direct `x`, `y`, `z` vectors — it assumes a uniform grid internally. My x, y and z vectors contain non-uniform data. How can i get uniform grid using interpolation? My x, y and z vectors are points of Mesh with many level of refinement.

 ![image](https://global.discourse-cdn.com/julialang/original/3X/e/3/e348c2262852c0da39504c956fdf4d5a04c49457.jpeg)

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 25, 2025, 9:54am UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/13 "2025-05-25T09:54:58Z")

</div>

Yes, that is what I remember from the last time I looked into it. If Makie.jl doesn’t offer a function for this curently, we can’t use it to optimize the vizualization of 3D irregular grids.

---

<div class="post-metadata">

**Author:** ![jules](https://avatars.discourse-cdn.com/v4/letter/j/41988e/32.png) [@jules](https://discourse.julialang.org/u/jules)\
**Post date:** [May 25, 2025, 2:06pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/14 "2025-05-25T14:06:49Z")

</div>

The way to go in these scenarios is usually to interpolate onto a regular grid and use that

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 25, 2025, 3:15pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/15 "2025-05-25T15:15:51Z")

</div>

That is a good point. You can use the `InterpolateNeighbors` transform to interpolate your `RectilinearGrid` values into a `RegularGrid` and then call `viz`. The transform is documented in GeoStats.jl and the GDSJL book.

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 30, 2025, 2:36pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/16 "2025-05-30T14:36:49Z")

</div>

My code in [above post](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/11) is able to get 3D scatter plot as shown earlier. Please can you see my grid points and tell how to interpolate?

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 30, 2025, 2:51pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/17 "2025-05-30T14:51:51Z")

</div>

Please fill in the gaps:

```julia
using GeoStats

# your original irregular grid with x, y, z vectors
grid1 = RectilinearGrid(...)

# new regular grid for visualization
grid2 = CartesianGrid(minimum(grid1), maximum(grid1), dims=(100,100,100))

# table of values in grid1
temperature = rand(nelements(grid1))
pressure = rand(nelements(grid1))
geotable = georef((; temperature, pressure), grid1)

# interpolate over grid2 with nearest neighbors
newgeotable = geotable |> InterpolateNeighbors(grid2)

# extract interpolated values from geotable
newgeotable.temperature
newgeotable.pressure

# alternatively, visualize directly
import GLMakie as Mke

newgeotable |> viewer # option 1

viz(newgeotable.geometry, color = newgeotable.temperature) # option 2

```

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 30, 2025, 4:25pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/18 "2025-05-30T16:25:02Z")

</div>

I got this error.

> **Summary**
>
> julia\> include(“Density.jl”)  
> size(ρ) = (18, 4, 16, 872)  
> size(z\_all) = (16, 872)  
> size(all\_zs) = (1004544,)  
> grid1 = RectilinearGrid(all\_xs, all\_ys, all\_zs) = 1004543×1004543×1004543 RectilinearGrid  
> ERROR: LoadError: ArgumentError: spacing must be positive  
> Stacktrace:  
> [1] RegularGrid  
> @ ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/regulargrid.jl:63 [inlined]  
> [2] RegularGrid  
> @ ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/regulargrid.jl:88 [inlined]  
> [3] RegularGrid(dims::Tuple{…}, origin::Meshes.Point{…}, spacing::Tuple{…}, offset::Tuple{…})  
> @ Meshes ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/regulargrid.jl:100  
> [4] RegularGrid  
> @ ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/regulargrid.jl:97 [inlined]  
> [5] #RegularGrid#195  
> @ ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/regulargrid.jl:125 [inlined]  
> [6] RegularGrid  
> @ ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/regulargrid.jl:121 [inlined]  
> [7] RegularGrid  
> @ ~/.julia/packages/Meshes/Pbr8O/src/domains/meshes/cartesiangrid.jl:17 [inlined]  
> [8] plothdf(h5file::String)  
> @ Main ~/Videos/Density.jl:41  
> [9] macro expansion  
> @ ./timing.jl:581 [inlined]  
> [10] top-level scope  
> @ ~/Videos/Density.jl:387  
> [11] include(fname::String)  
> @ Main ./sysimg.jl:38  
> [12] top-level scope  
> @ REPL[4]:1  
> in expression starting at /home/raman/Videos/Density.jl:53  
> Some type information was truncated. Use `show(err)` to see complete types.

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [May 30, 2025, 4:45pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/19 "2025-05-30T16:45:10Z")

</div>

Can you please confirm you have unique coordinates in x y z? Are they sorted?

---

<div class="post-metadata">

**Author:** ![raman\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raman_kumar/32/26782_2.png) [@raman\_kumar](https://discourse.julialang.org/u/raman_kumar)\
**Post date:** [May 30, 2025, 6:07pm UTC](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643/20 "2025-05-30T18:07:36Z")

</div>

On using `unique( )` i got x, y and z arrays with different length so i didn’t use it finally. But `sort( )` helped me proceed further.

> [@raman\_kumar](#):
>
> My Athena++ data is at this [link](https://users.camk.edu.pl/miki/sane00.prim.01800.athdf).

Now my this

> **updated Code**
>
> using HDF5, GLMakie, Colors, LinearAlgebra, GeoStats  
> GLMakie.activate!(ssao=true)  
> GLMakie.closeall() # close any open screen
> 
> function plothdf(h5file::String)
> 
> ```
> h5 = h5open(h5file, "r")
> num_blocks = read(HDF5.attributes(h5), "NumMeshBlocks")
> ρ = read(h5["prim"])[:,:,:,:,1] # shape: (Nx, Ny, Nz, Nvars, Nblocks)
> @show size(ρ)
> x_all = read(h5["x1v"]) # shape: (Nx, Nblocks)
> y_all = read(h5["x2v"])
> z_all = read(h5["x3v"])
> close(h5)
> @show size(z_all)
> 
> all_xs, all_ys, all_zs = Float32[], Float32[], Float32[]
> 
> for b in 1:num_blocks
> x = x_all[:, b]
> y = y_all[:, b]
> z = z_all[:, b]
> 
> for i in 1:length(x), j in 1:length(y), k in 1:length(z)
> r, θ, ϕ = x[i], y[j], z[k]
> X = r * sin(θ) * cos(ϕ)
> Y = r * sin(θ) * sin(ϕ)
> Z = r * cos(θ)
> 
> push!(all_xs, X)
> push!(all_ys, Y)
> push!(all_zs, Z)
> end
> end
> @show size(unique(all_ys))
> @show size(all_ys)
> @show size(unique(all_zs))
> all_xss = sort(all_xs)
> all_yss = sort(all_ys)
> all_zss = sort(all_zs)
> 
> # your original irregular grid with x, y, z vectors
> grid1 = RectilinearGrid(all_xss, all_yss, all_zss)
> 
> # new regular grid for visualization
> @show grid2 = CartesianGrid(minimum(grid1), maximum(grid1), dims=(1000,1000,1000))
> @show nelements(grid1)
> @show length(ρ)
> temperature = rand(nelements(grid1))
> pressure = rand(nelements(grid1))
> geotable = georef((; temperature, pressure), grid1)
> newgeotable = geotable |> InterpolateNeighbors(grid2)
> # Normalize density
> norm_density = log10.(abs.(normalize!(vec(ρ))))
> fig = Figure(size = (1000, 800))
> ssao = Makie.SSAO(radius = 5.0, blur = 3)
> ax = LScene(fig[1, 1], scenekw = (ssao=ssao,))
> ax.scene.ssao.bias[] = 1.95
> plt = scatter!(ax, all_xs, all_ys, all_zs, color = norm_density, markersize = 5, transparency =true)
> Colorbar(fig[1, 2], plt, label = "log₁₀(norm(ρ))")
> display(fig)
> 
> ```
> 
> end  
> @timev plothdf(“sane00.prim.01800.athdf”)

gives this error as `nelements(grid1) = 1013691010309291007`. I need volume plot of **ρ** values on these grid points.

> julia\> include(“Density.jl”)  
> size(ρ) = (18, 4, 16, 872)  
> size(z\_all) = (16, 872)  
> size(unique(all\_ys)) = (515483,)  
> size(all\_ys) = (1004544,)  
> size(unique(all\_zs)) = (17856,)  
> grid2 = CartesianGrid(minimum(grid1), maximum(grid1), dims = (1000, 1000, 1000)) = 1000×1000×1000 CartesianGrid  
> nelements(grid1) = 1013691010309291007  
> length(ρ) = 1004544  
> ERROR: LoadError: OutOfMemoryError()  
> Stacktrace:  
> [1] GenericMemory  
> @ ./boot.jl:516 [inlined]  
> [2] Array  
> @ ./boot.jl:578 [inlined]

[Next page](https://discourse.julialang.org/t/out-of-memory-error-in-viz-of-3d-rectilineargrid-from-meshes-jl/128643.md?page=2)