freesurface() error

[eikmeier-cn]

Hi there.

By running FWI with free surface, model.fs=True and model0.fs=True, we have a problem with freesurface().

model = demo_model('circle-isotropic', vp_circle=3.0, vp_background=2.5,
                    origin=origin, shape=shape, spacing=spacing, nbl=40)

model0 = demo_model('circle-isotropic', vp_circle=2.5, vp_background=2.5,
                     origin=origin, shape=shape, spacing=spacing, nbl=40,
                     grid = model.grid)

# Free surface setup
model.fs = True
model0.fs = True

# Compute synthetic data with forward operator 
from examples.seismic.acoustic import AcousticWaveSolver

solver = AcousticWaveSolver(model, geometry, space_order=4)
true_d, a_, b_ = solver.forward(vp=model.vp)

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
~/devito-4.2.2/devito/devito/tools/memoization.py in __call__(self, *args, **kw)
     88         try:
---> 89             res = cache[key]
     90         except KeyError:

KeyError: (<function AcousticWaveSolver.op_fwd at 0x7fc916ac7cb0>, (None,), frozenset())

During handling of the above exception, another exception occurred:

KeyError                                  Traceback (most recent call last)
<ipython-input-7-42bed3535510> in <module>
      3 
      4 solver = AcousticWaveSolver(model, geometry, space_order=4)
----> 5 true_d, a_, b_ = solver.forward(vp=model.vp)

~/devito-4.2.2/devito/examples/seismic/acoustic/wavesolver.py in forward(self, src, rec, u, vp, save, **kwargs)
    110 
    111         # Execute operator and return wavefield and receiver data
--> 112         summary = self.op_fwd(save).apply(src=src, rec=rec, u=u, vp=vp,
    113                                           dt=kwargs.pop('dt', self.dt), **kwargs)
    114         return rec, u, summary

~/devito-4.2.2/devito/devito/tools/memoization.py in __call__(self, *args, **kw)
     89             res = cache[key]
     90         except KeyError:
---> 91             res = cache[key] = self.func(*args, **kw)
     92         return res
     93 

~/devito-4.2.2/devito/examples/seismic/acoustic/wavesolver.py in op_fwd(self, save)
     51         return ForwardOperator(self.model, save=save, geometry=self.geometry,
     52                                kernel=self.kernel, space_order=self.space_order,
---> 53                                **self._kwargs)
     54 
     55     @memoized_meth

~/devito-4.2.2/devito/examples/seismic/acoustic/operators.py in ForwardOperator(model, geometry, space_order, save, kernel, **kwargs)
    129 
    130     s = model.grid.stepping_dim.spacing
--> 131     eqn = iso_stencil(u, model, kernel)
    132 
    133     # Construct expression to inject source values

~/devito-4.2.2/devito/examples/seismic/acoustic/operators.py in iso_stencil(field, model, kernel, **kwargs)
     92     # Add free surface
     93     if model.fs:
---> 94         eqns.append(freesurface(model, Eq(unext, eq_time)))
     95     return eqns
     96 

~/devito-4.2.2/devito/examples/seismic/acoustic/operators.py in freesurface(model, eq)
     20     lhs, rhs = eq.evaluate.args
     21     # Get vertical dimension and corresponding subdimension
---> 22     zfs = model.grid.subdomains['fsdomain'].dimensions[-1]
     23     z = zfs.parent
     24 

KeyError: 'fsdomain'

[mloubout]

You cannot change model.fs after creating the model because the subdomain associated with it won't exist. If you want to use a freesurface you have to create the model explicitly with it Model(...., fs=True)

[eikmeier-cn]

Thank you for the explanation @mloubout.