judiTopmute says: ERROR: MethodError: Cannot convert an object of type Nothing to an object of type Integer

[tjk9501]

Using a modified LSRTM 2D code

I will post my second trial first

I am defining right-hand preconditioners (model topmute)

Mr = judiTopmute(model0; taperwidth=10)

The output error is

ERROR: MethodError: Cannot `convert` an object of type Nothing to an object of type Integer

Closest candidates are:
  convert(::Type{Integer}, ::PyCall.PyObject)
   @ PyCall ~/.julia/packages/PyCall/1gn3u/src/conversions.jl:47
  convert(::Type{T}, ::T) where T<:Number
   @ Base number.jl:6
  convert(::Type{T}, ::T) where T
   @ Base Base.jl:84
  ...

Stacktrace:
 [1] setindex!(A::Vector{Integer}, x::Nothing, i1::Int64)
   @ Base ./array.jl:1021
 [2] setindex!
   @ ./multidimensional.jl:698 [inlined]
 [3] find_water_bottom(m::Matrix{Float32}; eps::Float64)
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/utils.jl:31
 [4] find_water_bottom(m::Matrix{Float32})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/utils.jl:25
 [5] judiTopmute(model::JUDI.IsoModel{Float32, 2}; taperwidth::Int64)
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:62
 [6] top-level scope
   @ REPL[72]:1

So, I am using a different Topmute definition method

wb_mask = ones(Float32,size(m))
wb_mask[1:Int(Xsrc[1])+20, :] .= 0;
wb_mask[Int(Xrec[1])-20:end, :] .= 0;
Mr = judiTopmute(model0.n, wb_mask, 10)

Then I begin inversion process

Mr = judiTopmute(model0.n, wb_mask, 10)
Ml = judiDataMute(q.geometry, d_obs.geometry)

x_iter = convert(Array{Float32,2}, model0.m)
fval = zeros(Float32, maxiter)

for j=1:maxiter
println("Iteration: ", j)

i = randperm(d_obs.nsrc)[1:batchsize]

r = Ml[i]*J[i]*Mr*x_iter - Ml[i]*d_obs[i]
g = adjoint(Mr)*adjoint(J[i])*adjoint(Ml[i])*r

fval[j] = .5f0*norm(r)^2
t = norm(r)^2/norm(g)^2
global x_iter -= t*g

The error begins at

r = Ml[i]*J[i]*Mr*x_iter - Ml[i]*d_obs[i]

Output message\

ERROR: joLinearOperatorException("shape mismatch")
Stacktrace:
 [1] *(A::joLinearOperator{Float32, Float32}, mv::Matrix{Float32})
   @ JOLI ~/.julia/packages/JOLI/S3LJy/src/joAbstractLinearOperator/base_functions.jl:155
 [2] afoldl
   @ ./operators.jl:544 [inlined]
 [3] *(a::DataMute{…}, b::judiJacobian{…}, c::TopMute{…}, xs::Matrix{…})
   @ Base ./operators.jl:587
 [4] top-level scope
   @ REPL[69]:1
Some type information was truncated. Use `show(err)` to see complete types.

So, I check the dimension of preconditioners and Jacobian

julia> Ml[i]
(DataMute{Float32, :reflection}, (8755776, 8755776), "DataMute{Float32, :reflection}")

julia> J[i]
JUDI born{Float32} propagator (z * x) -> (src * rec * time)

julia> Mr
(TopMute{Float32, 3, 2}, (8194613, 8194613), "TopMute{Float32, 3, 2}")

Please point out possible errors in my operations and instruct some solutions.

[tjk9501]

The above trial is originated from a first trial in testing LSRTM code

ERROR: MethodError: no method matching matvec(::TopMute{Float32, 3, 2}, ::Matrix{Float32})

Closest candidates are:
  matvec(::JUDI.TimeGain{T, K}, ::Array{T}) where {T, K}
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/DataPreconditioners.jl:374
  matvec(::TopMute{T, N}, ::Vector{T}) where {T, N}
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:82
  matvec(::TopMute{T, N}, ::judiWeights{T}) where {T, N}
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:81
  ...

Stacktrace:
 [1] matvec(D::TopMute{Float32, 3, 2}, x::PhysicalParameter{Float32, 2})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:80
 [2] *(J::TopMute{Float32, 3, 2}, ms::PhysicalParameter{Float32, 2})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/base.jl:29
 [3] (::JOLI.var"#175#179"{joLinearOperator{Float32, Float32}, TopMute{Float32, 3, 2}})(v3::Vector{Float32})
   @ JOLI ~/.julia/packages/JOLI/S3LJy/src/joAbstractLinearOperator/base_functions.jl:141
 [4] *(A::joLinearOperator{Float32, Float32}, v::Vector{Float32})
   @ JOLI ~/.julia/packages/JOLI/S3LJy/src/joAbstractLinearOperator/base_functions.jl:191
 [5] *(J::joLinearOperator{Float32, Float32}, x::judiVector{Float32, Matrix{Float32}})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Types/abstract.jl:86
 [6] mul!(out::Vector{Float32}, F::joLinearOperator{Float32, Float32}, q::judiVector{Float32, Matrix{Float32}})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/LinearOperators/operators.jl:183
 [7] lsqr_method!(log::ConvergenceHistory{…}, x::Vector{…}, A::joLinearOperator{…}, b::judiVector{…}; damp::Int64, atol::Float32, btol::Float32, conlim::Float32, maxiter::Int64, verbose::Bool)
   @ IterativeSolvers ~/.julia/packages/IterativeSolvers/bjiMt/src/lsqr.jl:130
 [8] lsqr!(x::Vector{…}, A::joLinearOperator{…}, b::judiVector{…}; maxiter::Int64, log::Bool, kwargs::@Kwargs{})
   @ IterativeSolvers ~/.julia/packages/IterativeSolvers/bjiMt/src/lsqr.jl:74
 [9] top-level scope
   @ REPL[65]:1

caused by: MethodError: no method matching matvec(::TopMute{Float32, 3, 2}, ::Matrix{Float32})

Closest candidates are:
  matvec(::JUDI.TimeGain{T, K}, ::Array{T}) where {T, K}
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/DataPreconditioners.jl:374
  matvec(::TopMute{T, N}, ::Vector{T}) where {T, N}
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:82
  matvec(::TopMute{T, N}, ::judiWeights{T}) where {T, N}
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:81
  ...

Stacktrace:
 [1] matvec(D::TopMute{Float32, 3, 2}, x::PhysicalParameter{Float32, 2})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:80
 [2] *(J::TopMute{Float32, 3, 2}, ms::PhysicalParameter{Float32, 2})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/base.jl:29
 [3] (::JOLI.var"#175#179"{joLinearOperator{…}, TopMute{…}})(v3::judiVector{Float32, Matrix{…}})
   @ JOLI ~/.julia/packages/JOLI/S3LJy/src/joAbstractLinearOperator/base_functions.jl:141
 [4] *(J::joLinearOperator{Float32, Float32}, x::judiVector{Float32, Matrix{Float32}})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Types/abstract.jl:86
 [5] mul!(out::Vector{Float32}, F::joLinearOperator{Float32, Float32}, q::judiVector{Float32, Matrix{Float32}})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/LinearOperators/operators.jl:183
 [6] lsqr_method!(log::ConvergenceHistory{…}, x::Vector{…}, A::joLinearOperator{…}, b::judiVector{…}; damp::Int64, atol::Float32, btol::Float32, conlim::Float32, maxiter::Int64, verbose::Bool)
   @ IterativeSolvers ~/.julia/packages/IterativeSolvers/bjiMt/src/lsqr.jl:130
 [7] lsqr!(x::Vector{…}, A::joLinearOperator{…}, b::judiVector{…}; maxiter::Int64, log::Bool, kwargs::@Kwargs{})
   @ IterativeSolvers ~/.julia/packages/IterativeSolvers/bjiMt/src/lsqr.jl:74
 [8] top-level scope
   @ REPL[65]:1
Some type information was truncated. Use `show(err)` to see complete types.

I am using a homogeneous model as starting model

model0 = Model(n, d, o, m, nb, rho)

Since source and receivers are located at top and bottom of the model, I manually created water mask between transducer locations and edge of the model

wb_mask = ones(Float32,size(m))
wb_mask[1:Int(Xsrc[1])+20, :] .= 0;
wb_mask[Int(Xrec[1])-20:end, :] .= 0;

Now I am defining right-hand preconditioners (model topmute)

Mr = judiTopmute(model0; taperwidth=10)

The output error is

ERROR: MethodError: Cannot `convert` an object of type Nothing to an object of type Integer

Closest candidates are:
  convert(::Type{Integer}, ::PyCall.PyObject)
   @ PyCall ~/.julia/packages/PyCall/1gn3u/src/conversions.jl:47
  convert(::Type{T}, ::T) where T<:Number
   @ Base number.jl:6
  convert(::Type{T}, ::T) where T
   @ Base Base.jl:84
  ...

Stacktrace:
 [1] setindex!(A::Vector{Integer}, x::Nothing, i1::Int64)
   @ Base ./array.jl:1021
 [2] setindex!
   @ ./multidimensional.jl:698 [inlined]
 [3] find_water_bottom(m::Matrix{Float32}; eps::Float64)
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/utils.jl:31
 [4] find_water_bottom(m::Matrix{Float32})
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/utils.jl:25
 [5] judiTopmute(model::JUDI.IsoModel{Float32, 2}; taperwidth::Int64)
   @ JUDI ~/.julia/packages/JUDI/tagnw/src/TimeModeling/Preconditioners/ModelPreconditioners.jl:62
 [6] top-level scope
   @ REPL[72]:1

So, I am using a different Topmute definition method

Mr = judiTopmute(model0.n, wb_mask, 10)

Then define

Ml = judiDataMute(q.geometry, d_obs.geometry)
indsrc = randperm(d_obs.nsrc)[1:d_obs.nsrc]
lsqr_sol = zeros(Float32, prod(size(model0)))

The optimization step is

dinv = d_obs[indsrc]
Jinv = J[indsrc]

Jp = Ml[indsrc]*Jinv*Mr
dinvp = Ml[indsrc]*dinv

lsqr!(lsqr_sol, Jp, dinvp; maxiter=maxiter)

The error at top of the post occurred at

lsqr!(lsqr_sol, Jp, dinvp; maxiter=maxiter)

[mloubout]

What does your model look like?

[mloubout]

I mean the values

[tjk9501]

a homogeneous model with m[i,j]=0.0818f0

[mloubout]

Ok that's hwat I thought, if your model is constant then it cannot detect the "bottom" of the first/water layer so it will error because it will not be able to setup correctly.

Your only option in this case is to define the top mute with the indices directly judiTopmute(n::NTuple{N, Integer}, wb::Array{T, Nw}, taperwidth::Integer)
https://github.com/slimgroup/JUDI.jl/blob/395793c9fd98f48922926f82cc2ab1d2a3bc3ff0/src/TimeModeling/Preconditioners/ModelPreconditioners.jl#L58C1-L58C74

[tjk9501]

Yes, I have tried that before, with

wb_mask = ones(Float32,size(m))
wb_mask[1:Int(Xsrc[1])+20, :] .= 0;
wb_mask[Int(Xrec[1])-20:end, :] .= 0;
Mr = judiTopmute(model0.n, wb_mask, 10)

but it will prompt error with

ERROR: joLinearOperatorException("shape mismatch")
Stacktrace:
 [1] *(A::joLinearOperator{Float32, Float32}, mv::Matrix{Float32})
   @ JOLI ~/.julia/packages/JOLI/S3LJy/src/joAbstractLinearOperator/base_functions.jl:155
 [2] afoldl
   @ ./operators.jl:544 [inlined]
 [3] *(a::DataMute{…}, b::judiJacobian{…}, c::TopMute{…}, xs::Matrix{…})
   @ Base ./operators.jl:587
 [4] top-level scope
   @ REPL[69]:1
Some type information was truncated. Use `show(err)` to see complete types.

in

r = Ml[i]*J[i]*Mr*x_iter - Ml[i]*d_obs[i]

[mloubout]

No you don't give a mask you give the indices of the bottom, i.e it's an array of one less dimension than the model.

If you want to use directly a mask you can just make a plain JOLI operator, i.e joDiag(wb_mask[:])

[tjk9501]

What will happen to the mute part? Will it be excluded in simulations? I also see there are data mute operators, e.g.,

Ml = judiDatMute(q.geometry, dD.geometry; t0=.120)	# data topmute starting at time 120ms

So, are the mute grids in judiTopmute never received observed data after 120ms?

[mloubout]

Then I need to provide a 800 by 1000 ones array to wb_mask, is that right?

No, judiTopMute is not a mask and does not accept a mask as an input. It is a "first x = 0" operation that accept indices, so for an x, z model it accept a list of zi for each x so that you set x[xi, 1:wb[xi]] = 0. If you want to use a mask you need to use joDiag or just apply your mask directly wb_mask .* m