fdWaveModel simulation methods¶
- void
fdWaveModel
::
forward_simulate
(int i_shot, bool store_fields, bool verbose)¶ Method to forward simulate wavefields for a specific shot.
Forward simulate wavefields of shot i_shot based on currently loaded models. Storage of wavefields and verbosity of run can be toggled.
- Parameters
-
i_shot
: Integer controlling which shot to simulate.store_fields
: Boolean to control storage of wavefields. If storage is not required (i.e. no adjoint modeling), forward simulation should be faster without storage.verbose
: Boolean controlling if modelling should be verbose.
- void
fdWaveModel
::
adjoint_simulate
(int i_shot, bool verbose)¶ Method to adjoint simulate wavefields for a specific shot.
Adjoint simulate wavefields of shot i_shot based on currently loaded models and calculate adjoint sources. Verbosity of run can be toggled.
- Parameters
-
i_shot
: Integer controlling which shot to simulate.verbose
: Boolean controlling if modelling should be verbose.
- void
fdWaveModel
::
write_receivers
()¶ Method to write out synthetic seismograms to plaintext.
This method writes out the synthetic seismograms to a plaintext file. Allows one to e.g. subsequently import these files as observed seismograms later using load_receivers(). Every shot generates a separate ux and uz receiver file (rtf_ux/rtf_uz), with every receiver being a single line in these files.
- void
fdWaveModel
::
write_sources
()¶ Method to write out source signals to plaintext.
This method writes out the source time function (without moment tensor) to plaintext file. Useful for e.g. visualizing the source staggering.
- void
fdWaveModel
::
load_receivers
(bool verbose)¶ Method to load receiver files.
This method loads receiver data from observed_data_folder folder into the object. The data has to exactly match the set-up, and be named according to component and shot (as generated by write_receivers() ).
- Parameters
-
verbose
: Controls the verbosity of the method during loading.
- void
fdWaveModel
::
load_model
(const std::string &de_path, const std::string &vp_path, const std::string &vs_path, bool verbose)¶ Method to load models from plaintext into the model.
This methods loads any appropriate model (expressed in density, P-wave velocity, and S-wave velocity) into the class and updates the Lamé fields accordingly.
- Parameters
-
de_path
: Relative path to plaintext de file.vp_path
: Relative path to plaintext vp file.vs_path
: Relative path to plaintext vs file.verbose
: Boolean controlling the verbosity of the method.
- void
fdWaveModel
::
map_kernels_to_velocity
()¶ Method to map kernels to velocity parameter set.
This method takes the kernels (lambda, mu, rho) as originally calculated on the Lamé’s parameter set and maps them to the velocity parameter set (vp, vs, rho).
- void
fdWaveModel
::
reset_kernels
()¶ Method to reset all Lamé sensitivity kernels to zero.
This method resets all sensitivity kernels to zero. Essential before performing new adjoint simulations, as otherwise the kernels of subsequent simulations would stack.
- void
fdWaveModel
::
update_from_velocity
()¶ Method to map velocities into Lamé’s parameter set.
This method updates Lamé’s parameters of the current model based on the velocity parameters of the current model. Typically has to be done every time after updating velocity.
- void
fdWaveModel
::
calculate_l2_misfit
()¶ Method to calculate L2 misfit.
This method calculates L2 misfit between observed seismograms and synthetic seismograms and stores it in the misfit field.
- void
fdWaveModel
::
calculate_l2_adjoint_sources
()¶ Method to calculate L2 adjoint sources This method calculates L2 misfit between observed seismograms and synthetic seismograms and stores it in the misfit field.
- void
fdWaveModel
::
run_model
(bool verbose) Method to perform all steps necessary for FWI.
This method performs all necessary steps in FWI; forward modelling, misfit calculation adjoint source calculation, adjoint modelling and kernel projection. Delegates to run_model(bool verbose, bool simulate_adjoint).
- Parameters
-
verbose
: Boolean controlling the verbosity of the method.
- void
fdWaveModel
::
run_model
(bool verbose, bool simulate_adjoint) Method to perform all steps necessary for FWI, with additional control over adjoint simulation.
This method performs all necessary steps in FWI; forward modelling, misfit calculation and optionally adjoint source calculation, adjoint modelling and kernel projection.
- Parameters
-
verbose
: Boolean controlling the verbosity of the method.simulate_adjoint
: Boolean controlling the execution of the adjoint simulation and kernel computation.