adjtomo / seisflows

An automated workflow tool for full waveform inversion and adjoint tomography
http://seisflows.readthedocs.org
BSD 2-Clause "Simplified" License
178 stars 122 forks source link

Issue with Adjoint Waveform Calculation in the Latest Devel Branch. #219

Open trinitite271 opened 3 months ago

trinitite271 commented 3 months ago

Hello,

I have been using the latest devel branch of seisflows to set up an example. The real and initial models are smaller, closer to exploration-sized models. I generated these models using specfem, incorporating a tomo layer as done in the specfem marmousi example (OUTPUT_FILES_INIT and OUTPUT_FILES_TRUE).

The issue I'm encountering is that the observed data (obs data) and synthetic data (syn data) appear identical, resulting in zero adjoint data (adj data), which prevents the correct calculation of the adjoint waveform. However, the output for both the initial model and the observed model seems to be correct.

I have documented my process in a Jupyter notebook named test.ipynb available at my GitHub repository: https://github.com/trinitite271/seisflows2. This notebook includes steps for reading the Matlab models, writing models, setting up models, running the program, etc. Could you please review the results in this notebook and identify any possible errors in my steps?

Thank you for your assistance.

test1 test2

bch0w commented 3 months ago

Hey @trinitite271, thanks for opening up this issue. And thanks for providing the GitHub directory, that was really helpful for digging deeper into this issue.

Looks like for some reason SPECFEM2D is not using the gll option for defining the starting model, and ends up re-using your target model which would explain why your data and synthetics are the exact same and that your adjoint sources are zero.

Incorrect starting model (this should read MODEL: gll): https://github.com/trinitite271/seisflows2/blob/master/scratch/solver/001/fwd_solver.log#L233-L251 Target model showing same parameters: https://github.com/trinitite271/seisflows2/blob/master/specfem2d_workdir/OUTPUT_FILES_TRUE/solver_log.txt#L233-L251

Weirdly though the model parameter in the solver's Par_file is set correctly to gll so I'm not sure why SPECFEM2D is ignoring this and using the tomo option in the solver run instead.

Par_file: https://github.com/trinitite271/seisflows2/blob/master/scratch/solver/001/DATA/Par_file#L65

It's also set correctly in the master file which means it was never accidentally set to tomo at a previous step. Also pretty strange!: https://github.com/trinitite271/seisflows2/blob/master/specfem2d_workdir/DATA/Par_file#L65

I'll need to do a little digging, at first glance I'm wondering if this is an issue with SPECFEM2D overwriting the model parameter when it recognizes a tomography model? Not sure, hopefully I can get back to you with some useful information

@scott314159 this may also be related to the issue you're seeing? If you could take a look at the relevant log files and see what model type specfem is using during the solver run, that should help confirm/deny.

bch0w commented 3 months ago

Aha! And here is the magic line:

https://github.com/SPECFEM/specfem2d/blob/devel/src/specfem2D/setup_mesh.F90#L582-L585

if (tomo_material > 0) MODEL = 'tomo'

Corresponding to this line in the Par_file: https://github.com/trinitite271/seisflows2/blob/master/specfem2d_workdir/DATA/Par_file#L289

1 -1   0.d0    0.d0    1.d0    0 0    0    0 0 0 0 0 0 0     # solid from tomo file

Which is forcing SPECFEM to read your tomo file rather than use the required GLL model. I am pretty confident this is what is causing these issues. A quick fix would be to change that leading '1' to '-1' and I think things should work. For a longer term fix we'll have to fix that in SPECFEM2D, or atleast throw a warning that the Par_file parameter is being ignored.

Thanks for sticking with me there, hopefully this fixes your current issue. @scott314159 I think this is probably also the issue you are facing as well but let me know if you see a similar trend.

trinitite271 commented 3 months ago

I tried implementing your approach, and it seems to have resolved the issue. I reviewed the log file and understood your idea. The problem was that specfem failed to use the GLL mode when using the tomo files. The format in specfem for nbmodels is as follows: model_number material rho Vp 0 0 0 QKappa 9999 0 0 0 0 0 0. I suspect the actual bug occurs when material = -1, which forces the use of tomo, ignoring the model.

Based on my understanding of seisflows, during the inversion process, seisflows calls binary model files (gll) from OUTPUT_TRUE and OUTPUT_INIT. When material = -1, tomo is forcibly used; when material = 1, the nbmodel is ignored and seisflows use the gll model.

A quick fix would be to change that leading '1' to '-1' and I think things should work.

Therefore, the quick fix is to change -1 to 1 in seisflows before submitting, to avoid using tomo.

Besides, there is one small problem. The small one is, base on the solution above, I try to change the nbmodel before 'seisflows submit',but:

seisflows sempar velocity_model '1 1 1960.d0 1600.d0 800.0d0 0 0 9999 9999 0 0 0 0 0 0'

Traceback (most recent call last): File "/home/zhangchang/miniconda3/envs/seisflow_devel/bin/seisflows", line 8, in sys.exit(main()) ^^^^^^ File "/csim1/zhangchang/seisflow_devel/seisflows/seisflows.py", line 1456, in main sf() File "/csim1/zhangchang/seisflow_devel/seisflows/seisflows.py", line 456, in call getattr(self, self._args.command)(**vars(self._args)) File "/csim1/zhangchang/seisflow_devel/seisflows/seisflows.py", line 972, in sempar setpar_vel_model(file=par_file, model=value.split("+")) File "/csim1/zhangchang/seisflow_devel/seisflows/tools/specfem.py", line 620, in setpar_vel_model model_idx_start = model_lines[0]


IndexError: list index out of range