slimgroup / JUDI.jl

Julia Devito inversion.
https://slimgroup.github.io/JUDI.jl
MIT License
97 stars 31 forks source link

Forward modeling error #76

Closed ziyiyin97 closed 2 years ago

ziyiyin97 commented 2 years ago

MFE is

using JUDI
using LinearAlgebra

n = (50, 50)
d = (10f0, 10f0)
o = (0f0, 0f0)
v = 1.5f0 * ones(Float32, n)
m = 1f0./v.^2f0

model = Model(n,d,o,m; nb=80)

nsrc = 3
nrec = 50

dtS = dtR = 4f0
timeS = timeR = 1000f0

xsrc = convertToCell([100f0, 200f0, 300f0])
ysrc = convertToCell(range(0f0,stop=0f0,length=nsrc))
zsrc = convertToCell(range(10f0,stop=10f0,length=nrec))

srcGeometry = Geometry(xsrc, ysrc, zsrc; dt=dtS, t=timeS)
wavelet = ricker_wavelet(timeS,dtS,0.015f0)

q = judiVector(srcGeometry, wavelet)

xrec = range(d[1],stop=(n[1]-1)*d[1],length=nrec)
yrec = 0f0
zrec = range(10f0, stop=10f0, length=nrec)

recGeometry = Geometry(xrec,yrec,zrec; dt=dtR, t=timeR, nsrc=nsrc)

ntComp = get_computational_nt(srcGeometry, recGeometry, model)
info = Info(prod(n), nsrc, ntComp)

opt = Options()

F = judiProjection(info, recGeometry)*judiModeling(info, model; options=opt)*adjoint(judiProjection(info, srcGeometry))

dobs = F*q

println(norm(dobs.data[2]-dobs.data[3])) # this is 0

It works fine with 2 sources but when it goes to 3 sources, the 2nd and 3rd shot records are exactly the same. I'll look into how the new parallel strategy affected this (possibly the error comes from #71 )

ziyiyin97 commented 2 years ago

It seems to work fine w/ 4 sources as well. So I would guess it works well w/ 2^n sources

mloubout commented 2 years ago

Then it's probably an issue with the reduction