Regarding Parallelization in JUDI

[asidosaputra]

I’m experiencing an issue where the per-shot per-worker forward modeling process becomes significantly slower when executed in parallel.

My setup:

Ubuntu 22.04 running inside Windows Subsystem for Linux (WSL)

Julia version 1.10 LTS, JUDI v4.1.3

PC workstation with 32 cores, 64 threads, 256 GB RAM

I attempted to parallelize data generation from the Overthrust example provided [https://github.com/slimgroup/JUDI.jl/blob/master/examples/software_paper/generate_data_overthrust_3D.jl], with some minor modifications (code attached below).

using Distributed
addprocs(16)

@everywhere begin
using JUDI, HDF5
end

# Load overthrust model
if ~isfile("$(JUDI.JUDI_DATA)/overthrust_3D_true_model.h5")
    ftp_data("ftp://slim.gatech.edu/data/SoftwareRelease/WaveformInversion.jl/3DFWI/overthrust_3D_true_model.h5")
end
n, d, o, m = read(h5open("$(JUDI.JUDI_DATA)/overthrust_3D_true_model.h5","r"), "n", "d", "o", "m")

# Set up model structure
model = Model((n[1], n[2], n[3]), (d[1], d[2], d[3]), (o[1], o[2], o[3]), m)

# Set up source grid

xsrc = 400f0: 600f0: 19600f0	# 400 m space from boundary
ysrc = 400f0: 600f0: 19600f0
zsrc = 50f0

nsrc = length(xsrc)^2
(xsrc, ysrc, zsrc) = setup_3D_grid(xsrc, ysrc, zsrc)

xsrc = convertToCell(xsrc)
ysrc = convertToCell(ysrc)
zsrc = convertToCell(zsrc)

timeS = 3000f0
dtS = 4f0
src_geometry = Geometry(xsrc, ysrc, zsrc; dt = dtS, t = timeS)

## Set up receivers
xmin = 100f0; xmax = 19900f0
ymin = 100f0; ymax = 19900f0
offset = 6000f0
nxrec = 240
nyrec = 240
xrec = Array{Any,1}(undef, nsrc)
yrec = Array{Any,1}(undef, nsrc)
zrec = Array{Any}(undef, nsrc)

for j = 1: nsrc
	# 6 km max. offset around current source
    xlocal1 = minimum(clamp.(xsrc[j] .- offset, xmin, xmax))  # Lower bound
    xlocal2 = maximum(clamp.(xsrc[j] .+ offset, xmin, xmax))  # Upper bound

    ylocal1 = minimum(clamp.(ysrc[j] .- offset, ymin, ymax))  # Lower bound
    ylocal2 = maximum(clamp.(ysrc[j] .+ offset, ymin, ymax))  # Upper bound

	xlocal = xlocal1: 50f0: xlocal2
	ylocal = ylocal1: 50f0: ylocal2
    
	zlocal = 500f0
	(xlocal,ylocal,zlocal) = setup_3D_grid(xlocal, ylocal, zlocal)
	xrec[j] = xlocal
	yrec[j] = ylocal
	zrec[j] = zlocal
end



timeR = 3000f0	# receiver recording time [ms]
dtR = 4f0	# receiver sampling interval
rec_geometry = Geometry(xrec, yrec, zrec; dt = dtR, t = timeR)

# Set up source
f0 = 0.008
wavelet = ricker_wavelet(timeS, dtS, f0)

# #######################################################################

# Set up modeling options
opt = Options(limit_m = true,
			  save_data_to_disk = true,
			  file_path = "/mnt/e/Asido/Data/Overthrust/Gathers",   # replace w/ path to directory in which data will be stored (~total of 1.7 TB)
 			  file_name = "overthrust_3D_shot_"   # adds x-y coordinates to name
              )

# Set up operators
Pr = judiProjection(rec_geometry)
Ps = judiProjection(src_geometry)
A_inv = judiModeling(model; options=opt)
q = judiVector(src_geometry, wavelet)

# Generate data and save as individual SEG-Y files to disk
d_obs = Pr*A_inv*adjoint(Ps)*q

The main problem is that the per-shot runtime becomes extremely long when running in parallel — in my case, about 1400 seconds per worker (with 16 workers).
In contrast, running the same task serially only takes about 60–110 seconds per shot.
This is almost 10× slower in parallel mode.

Is this expected behavior, or am I making a mistake in how I execute the code?

I have also attached some data for reference.
Any explanation or guidance would be greatly appreciated. Thank you.