[CIG-SEISMO] Mineos problem
Martin van Driel
vandriel at tomo.ig.erdw.ethz.ch
Fri Mar 4 00:45:43 PST 2016
Hi all,
that is more than four years ago, so I am scratching my head: the
problems I had where related to the solver not converging for shorter
periods (around 50s I think). I tried a few easy workarounds (higher
model sampling in the crust, ignoring these modes in the mode sum...),
but none of them really worked out for me, so I moved to using another
code. I was only interested in the synthetic seismograms though, not the
eigenfunctions.
See below for two threads from this discussion. Keep in mind that I had
just started in seismology back then ;)
Cheers,
Martin
On 03.03.2016 17:46, Matthew Knepley wrote:
> I think there are some posts from Martin van Driel in the archive of
> this mailing list that describe some of the problems in more detail.
>
>
> Are these archives searchable? The link from the Mineos webpage does not
> incorporate a search.
>
> Hey Evan,
>
> I got the script working, see attachment (this time compressed so its
> hopefully not filtered again).
>
> But what should I say, the results are pretty frustrating. Modes that do
> work with the same model when I call minos_bran on a whole set of modes
> suddenly have extremely large rayleigh coefficients. On the other hand,
> not a single mode in the frequency range I am interested in gets stuck
> up when calling minos bran on the single modes, but it does when
> minos_bran is called normally.
>
> This means, that minos_bran gives different results for the same model
> and mode depending on what it did before, you called it path dependent.
> That's a very strong indication, that sth is not physical.
>
> When I do the normal call, I could get seismograms pretty similar to the
> other method I am benchmarking (attached), at least for the Z component,
> not that good for E (which corresponds to the S Modes).
>
> When I do the same with my single mode script, the fit is a lot worse.
>
> Actually I planned to use Mineos as the good old working method for
> benchmarking my new code. For the moment it is the other way around.
> Very frustrating.
>
> Cheers,
> Martin
>
> On 11/07/2011 08:17 AM, Evan Welch wrote:
>> Hello Martin,
>>
>> Sorry for the delay. I was on a mountaineering trip actually and
>> wasn't able to re-respond to your email.
>>
>> I can definitely give you more information tomorrow when I get back to
>> campus but from what I can do remotely from my friend's computer I can
>> offer you these two scripts which are very useful for running mineos.
>> One, "run_minos_..." will run the code to compute the modes one by one
>> which is sort of slow-going but the only effective way i know to
>> prevent it from hanging. Although for most models like PREM I have
>> seen the code compute the first chunk of modes (for period >= 5
>> seconds) for n=0 up to about 40 before hanging. So you may want to
>> compute this first segment, and then use the script to make the rest
>> of the modes 1by1.
>>
>> The other script, sortONESBY.csh is a quick and dirty way to
>> concatenate all the modes you have computed with the previous modes
>> into larger binary files by N value and all in one folder. These files
>> you can further concatenate into a single binary file if you like
>> which will identically match the output of minos_bran.f if it were to
>> have worked properly. The one thing the script does not do, which I
>> will fix soon and can update you, is make a corresponding text output
>> file for these fully concatenated modes. When you have the
>> minos_bran.f binary output and the minos_bran text output. You can
>> freely continue with the MINEOS software by running EIGCON.f which
>> takes these two as inputs.
>>
>> Note that for both scripts you will have to change the variables that
>> specifically apply to my directory structure and naming convention of
>> modes.
>>
>> Sorry if this is a little scrambled -- I've been awake for over a day
>> and can't tell if I'm making sense!
>>
>> Best Wishes,
>> Evan
>>
>>
>> On Nov 4, 2011, at 10:48 AM, Martin van Driel wrote:
>>
>>> Hi Evan,
>>>
>>> thanks for the quick reply, sounds promising. It would be really
>>> kind if you could send me the scripts you are mentioning.
>>>
>>> As I am only interested in MINEOS to produce reference solutions,
>>> any dirty fix is fine for me, as long as I can rely on the results.
>>>
>>> Thanks again!
>>> Martin
>>>
>>> On 11/04/2011 03:42 PM, Evan Welch wrote:
>>> > Hi Martin,
>>> >
>>> > So I am still working out some kinks with MINEOS but I /have/
>>> > received
>>> > some advice since then. Just after I send this I will forward you my
>>> > correspondence with Guy Masters who actually wrote the problematical
>>> > piece of software, minos_bran.f.
>>> >
>>> > For computing a full mode catalogue I have found the only solution
>>> > that
>>> > works is computing the modes "one-by-one" in a shell script by
>>> > running
>>> > minos_bran.f in a while loop with the "timeout" function and then
>>> > concatenating all the resulting binary files. I have a simple
>>> > script if
>>> > you would like it as an example. Using the timeout function, one can
>>> > allow minos_bran.f only a few seconds to spend on a particular mode
>>> > and
>>> > then give up and move on if unsuccessful.
>>> >
>>> > The regular culprit for both of us seems to be the high frequency
>>> > Stoneley modes -- the integration of which apparently gets trapped on
>>> > the ICB. If you wish to debug, check out the subroutine match and
>>> > subroutine remedy in minos_bran.f, these are the codes attempts to
>>> > fix
>>> > these pathological eigenfunctions. If you insert a few print
>>> > statements
>>> > on the variables a1 and a2 or fcat in the match subroutine you'll
>>> > notice
>>> > that when the code hangs is when these variables become Infinite or
>>> > NaN.
>>> >
>>> > In the correspondence I will forward you, Guy provides a "fix" for
>>> > this.
>>> > Also, as for trying different models, the more model layers in the
>>> > upper
>>> > mantle the more successful is the cubic interpolation.
>>> >
>>> > Best of luck,
>>> > Evan
>>> > On Nov 4, 2011, at 10:27 AM, Martin van Driel wrote:
>>> >
>>> >> Hi Evan,
>>> >>
>>> >> I am having the very same problem you posted to the CIG-SEISMO
>>> >> mailing
>>> >> list in October. Did you finally find a solution to this?
>>> >>
>>> >> I am testing different models that I generate based on the software
>>> >> that I want to benchmark (AXISEM) and it seems that some models have
>>> >> less problems then others, but I don't see any systematics. As
>>> >> soon as
>>> >> I want to go to higher frequencies then 30mHz it seems to be more or
>>> >> less random whether it works or not.
>>> >>
>>> >> The other thing I am worrying about in this context is the "raylquo"
>>> >> output in the output of minos brans: according to the manual it
>>> >> should
>>> >> be on the order of "eps", but there are always modes where it goes
>>> >> up
>>> >> to the order of 1, no matter what eps I choose as input.
>>> >>
>>> >> Any hint would be apreciated,
>>> >> Thanks,
>>> >> Martin van Driel
>>> >>
>>> >>
>>> >> > Hello there,
>>> >> >
>>> >> > I am a student at Princeton working with Frederik Simons. We would
>>> >> > like to compute a full catalog of normal modes as alluded to in
>>> >> the
>>> >> > "Benchmark" section of the Mineos manual. That is nmin = 0, nmax =
>>> >> > 400, lmin = 0, lmax = 3000 from 0 Hz to .2 Hz. This is what the
>>> >> manual
>>> >> > says is done. However, after spending a week or so trying to get
>>> >> the
>>> >> > package to compute the modes, we are still unsuccessful. The
>>> >> problem
>>> >> > is that on various modes the code gets hung up.
>>> >> >
>>> >> > Is there a known strategy for computing a full mode catalog?
>>> >> Breaking
>>> >> > it up by n? by l? by frequency range?
>>> >> >
>>> >> > Also, how may we get access to the set of normal modes used in the
>>> >> > benchmarks section. This would be useful to do for
>>> >> > comparison. Unfortunately there are no error messages when the
>>> >> code
>>> >> > hangs indefinitely so we cannot diagnose the problem.
>>> >> >
>>> >> > It's possibly that you are not the person I am supposed to ask for
>>> >> > help but I hope, if not, you will lead me in the right direction.
>>> >> >
>>> >> > Thank You,
>>> >> > Evan Welch
>>> >> > (313) 303 9880
Martin, below is my correspondence with Guy Masters. Let me know if you
make progress or would like to talk more about how to compute modes.
Using MINEOS can be quite frustrating!
--Evan
Begin forwarded message:
> From: masters guy <tmasters at ucsd.edu>
> Date: November 2, 2011 8:11:07 PM EDT
> To: Evan Welch <efwelch at Princeton.EDU>
> Subject: Re: MINEOS QUESTION: minos_bran hanging up while making full
mode catalogue
>
> B1 only does first 4 branches which should be fine
>
> B2 turns gravity off -- i am not sure how many mode branches were
included but i don;t see any well-developed body waves
>
> I'm not sure how they did B3 but the seismograms look pretty
long-period to me -- i'm guessing they filtered out the bad stuff
>
> Guy
>
>
>
> On Nov 2, 2011, at 4:58 PM, Evan Welch wrote:
>
>> Do you know how the "Benchmarks" were computed then? All I am trying
to do is what was outlined in Appendix B. If someone has that mode
catalogue with some ~250,000 spheroidal modes lying around I would
happily take it as is.
>>
>> Thanks for this extra information below.
>>
>> Evan
>> On Nov 2, 2011, at 7:53 PM, masters guy wrote:
>>
>>> I refer you to page 10 section 1.3 of the handbook -- you are
trying to use mineos well out of its comfort range.
>>>
>>> match is used by remedy -- so yes, you are seeing problems with
stoneley modes -- these are incredibly locked on the core-mantle
boundary and are just extremely difficult to compute. My other code
which I use for high frequency modes cam manage up to about l=350 --
you may have to reduce this a bit to work with your version of mineos --
i am still running some tests but it seems that l=300 may be ok -- you
should check the rayleigh quotient for these modes as you go up in l --
at some point, the number will just go crazy -- another check is that
the group velocity is about 8.2km/s.
>>>
>>> I'll also run some tests to see when the code fails on the inner
core modes so the second part of the fix I gave you may also have to
adjusted.
>>>
>>> These bad modes must be removed from the catalog or they can screw
up your synthetics
>>>
>>> Guy
>>>
>>> On Nov 2, 2011, at 4:25 PM, Evan Welch wrote:
>>>
>>>> Thank you. I will try your fix. Just to update you after some more
debugging, the first time the code first generates "Infinity" values for
the variables a1 and a2 is within the subroutine "match". The statement
looks something like this:
>>>>
>>>> a2=(af(3,1)*afr(1)-af(1,1)*afr(3))/(af(1,2)*af(3,1)-af(1,1)
>>>> + *af(3,2))
>>>> a1=(af(3,2)*afr(1)-af(1,2)*afr(3))/(af(1,1)*af(3,2)-af(3,1)
>>>> + *af(1,2))
>>>>
>>>> Some of the inputs that go into a2 and a1 here are very large or
very small which cause a2 and a1 to become infinite. Later on this
causes NaNs to enter the mix. This specific example happened for n = 43,
and L = [560,570] and the code fails on L = 569.
>>>>
>>>> Is this more specific information helpful to identifying the
problem? Does it corroborate your suspicion below about the code failing
on Stoneley modes?
>>>>
>>>> We really appreciate your time helping us.
>>>>
>>>> -- Evan
>>>>
>>>>
>>>> On Nov 2, 2011, at 7:13 PM, masters guy wrote:
>>>>
>>>>> ok -- high frequency stoneley modes are really difficult to
compute and have zero surface displacement -- so if you are happy with
just ignoring then you should make the following fix. In subroutine
detqn, find the call toe subroutine remedy, it should look something like:
>>>>>
>>>>> 90 if(dabs(det).gt.1.d-3.and.ls.le.noc) call remedy(ls)
>>>>> call eifout(ls)
>>>>>
>>>>> then replace with
>>>>>
>>>>> 90 if(dabs(det).gt.1.d-3.and.ls.le.noc) then
>>>>> ccc**** catch high l stoneley modes which have no surface
displacement
>>>>> if(ls.gt.nic.and.l.gt.350) then
>>>>> cg=0.
>>>>> wray=0.
>>>>> qinv=0.
>>>>> ls=n+1
>>>>> return
>>>>> end if
>>>>> cc*** following is a conservative calculation of the frequency
>>>>> cc*** below which we know an IC mode has no surface displacement
>>>>> cc*** (note lower freq overtones are more concentrated on ICB)
>>>>> if(ls.lt.nic) then
>>>>> c fcut=.00242*l*l+0.798*l-8.64
>>>>> fcut=1.13*l-56.5
>>>>> fdim=500.*wdim/pi
>>>>> if(fdim.lt.fcut) then
>>>>> cg=0.
>>>>> wray=0.
>>>>> qinv=0.
>>>>> ls=n+1
>>>>> return
>>>>> end if
>>>>> end if
>>>>> call remedy(ls)
>>>>> end if
>>>>> call eifout(ls)
>>>>>
>>>>> hope that helps
>>>>>
>>>>> If you are trying to compute a complete catalog of modes for
prem_noocean, you will need to increase the number of knots near the
surface to get accurate values for the Q of the mode and the group
velocity (keep an eye on the rayleigh quotient at the end of each output
line). These are computed using a cubic interpolation of eigenfunctions
between model knots which can be poor for highly exponential modes. The
mode frequency is always accurate but you also need the Q to make
synthetics. Also, you will not be able to accurately interpolate the
mode excitation factors if the source depth is between model knots.
>>>>>
>>>>> I have a version of the code which uses RK to do the Q and U
integrals -- but this version doesn't include self-gravity (ie -- just
the cowling approx)
>>>>>
>>>>> Guy
>>>>>
>>>>> On Nov 2, 2011, at 1:37 PM, Evan Welch wrote:
>>>>>
>>>>>> Right so I am running minos_bran.f from the CIG package with the
following parameters:
>>>>>>
>>>>>> model: prem_noocean.txt [ the one i sent you ]
>>>>>> output file
>>>>>> output eigenfunction file
>>>>>>
>>>>>> EPS = 1e-12
>>>>>> wgrav = 200.0
>>>>>> jcom = 3 (for spheroidal modes)
>>>>>>
>>>>>> lmin = 0
>>>>>> lmax = 3000
>>>>>> fmin = 0.01
>>>>>> fmax = 200.0
>>>>>> nmin = 0
>>>>>> nmax = 400
>>>>>>
>>>>>> and the code hangs during the n=43 modes because the value of
"xsq" as the second argument of the subroutine bfs becomes NaN. However
if I find the mode for which it hangs and then try to run minos_bran on
that mode in isolation the code doesn't hang. In this sense the problem
is "path-dependent".
>>>>>>
>>>>>> Evan
>>>>>> On Nov 2, 2011, at 4:23 PM, masters guy wrote:
>>>>>>
>>>>>>> could you tell me exactly what control variables you are using
-- and where it first hangs for you
>>>>>>>
>>>>>>> Guy
>>>>>>>
>>>>>>> On Nov 2, 2011, at 11:58 AM, Evan Welch wrote:
>>>>>>>
>>>>>>>> Thank you for getting back to me Professor.
>>>>>>>>
>>>>>>>> I have attached the model file in this email. Essentially we are
>>>>>>>> trying to reproduce the "Benchmark" section of the attached
mineos
>>>>>>>> manual from the CIG site. This manual highlights a comparison
between
>>>>>>>> mineos and SPECFEM_3D.
>>>>>>>>
>>>>>>>> We would like to compute all modes with periods greater than 5
>>>>>>>> seconds, N = [0,400]. We are setting the integration precision to
>>>>>>>> 1e-12 as recommended in the manual and we are trying to not
ignore any
>>>>>>>> gravitation terms for the calculation.
>>>>>>>>
>>>>>>>> Since we are actually doing a lot of subsequent calculation of
>>>>>>>> seismograms in MATLAB and not relying on the other programs in
the
>>>>>>>> mineos package a different i/o format of a mode catalogue
won't be a
>>>>>>>> problem. As long as you can give me a hint at how to read it I
>>>>>>>> shouldn't have any trouble writing some quick code to parse the
>>>>>>>> binary. Do you have such a full catalogue computed?
>>>>>>>>
>>>>>>>> Thank you,
>>>>>>>> Evan
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> On Nov 2, 2011, at 2:17 PM, masters guy wrote:
>>>>>>>>
>>>>>>>> > hi evan
>>>>>>>> >
>>>>>>>> > sorry about not replying sooner
>>>>>>>> >
>>>>>>>> > please send me the model file as I do not use the CIG
version of
>>>>>>>> > mineos. This also means that the io format is different so
my mode
>>>>>>>> > catalogs might be hard for you to read.
>>>>>>>> >
>>>>>>>> > minos_bran makes an estimate of the frequency of the next
mode along
>>>>>>>> > the branch using the group velocity of the previous mode --
but each
>>>>>>>> > mode is essentially computed independently
>>>>>>>> >
>>>>>>>> > what value are you using for the integration precision?
>>>>>>>> >
>>>>>>>> > what frequency range are you computing
>>>>>>>> >
>>>>>>>> > do you want the catalog to compute long period synthetics or
short
>>>>>>>> > period synthetics? I use a different version of the code to
compute
>>>>>>>> > body wave synthetics that uses the cowling approximation. This
>>>>>>>> > handles difficult modes (inner core modes and stoneley
modes) better
>>>>>>>> >
>>>>>>>> > Guy
>>>>>>>> >
>>>>>>>> >
>>>>>>>> >
>>>>>>>> > On Nov 2, 2011, at 10:14 AM, Evan Welch wrote:
>>>>>>>> >
>>>>>>>> >> Professor Masters,
>>>>>>>> >>
>>>>>>>> >> Since it has been a while since my last email to you I
thought I'd
>>>>>>>> >> try again. Frederik and I are still working on trying to get
>>>>>>>> >> minos_bran.f to compute a full mode catalogue for
prem_noocean.txt
>>>>>>>> >> listed in the subpath DEMO/models/ of the mineos-1.0.2
release. Do
>>>>>>>> >> you have a robust strategy for computing these modes and
combat the
>>>>>>>> >> programs tendency to hang?
>>>>>>>> >>
>>>>>>>> >> Do you have a complete mode catalogue we could perhaps get
off your
>>>>>>>> >> ftp site?
>>>>>>>> >>
>>>>>>>> >> And one more question: Will there be any negative effect on
>>>>>>>> >> accuracy of the mode catalogue if the modes are computed
one at a
>>>>>>>> >> time? For example, does the minos_bran software use the
neighboring
>>>>>>>> >> modes to compute the next one?
>>>>>>>> >>
>>>>>>>> >> Thank you for your help. If you would prefer to skype or
talk on
>>>>>>>> >> the phone instead of reply via email that would be fine
with us.
>>>>>>>> >>
>>>>>>>> >> Best,
>>>>>>>> >> Evan Welch
>>>>>>>> >> (313) 303 - 9880
>>>>>>>> >> On Oct 17, 2011, at 2:33 PM, masters guy wrote:
>>>>>>>> >>
>>>>>>>> >>> Hi Evan -- please send me the model you are trying to
compute the
>>>>>>>> >>> mode catalog for
>>>>>>>> >>>
>>>>>>>> >>> Guy
>>>>>>>> >>>
>>>>>>>> >>>
>>>>>>>> >>> On Oct 17, 2011, at 9:04 AM, Evan Welch wrote:
>>>>>>>> >>>
>>>>>>>> >>>> Hello Professor Masters,
>>>>>>>> >>>>
>>>>>>>> >>>> I am a student at Princeton working with Frederik Simons. We
>>>>>>>> >>>> would like to compute a full catalogue of normal modes as
alluded
>>>>>>>> >>>> to in the "Benchmark" section of the Mineos manual. That
is nmin
>>>>>>>> >>>> = 0, nmax = 400, lmin = 0, lmax = 3000 from 0 Hz to .2
Hz. This
>>>>>>>> >>>> is what the manual says is done. However, after spending
a few
>>>>>>>> >>>> weeks or so trying to get the package to compute the
modes, we
>>>>>>>> >>>> are still unsuccessful. The problem is that on various
modes the
>>>>>>>> >>>> code gets hung up - seemingly due to the propagation of
NaNs -
>>>>>>>> >>>> the origin of which I have been able to trace as far back
as the
>>>>>>>> >>>> "gauslv" subroutine.
>>>>>>>> >>>>
>>>>>>>> >>>> Is there a known strategy for computing a full mode catalog?
>>>>>>>> >>>> Breaking it up by n? by l? by frequency range?
>>>>>>>> >>>>
>>>>>>>> >>>> Also, how may we get access to the set of normal modes
used in
>>>>>>>> >>>> the benchmarks section. This would be useful to do for
>>>>>>>> >>>> comparison. Unfortunately there are no error messages
when the
>>>>>>>> >>>> code hangs indefinitely so we cannot diagnose the problem.
>>>>>>>> >>>>
>>>>>>>> >>>> I have already extended my question to those on the CIG
seismo
>>>>>>>> >>>> list and they, in addition to Frederik, have pointed me
to you.
>>>>>>>> >>>> Please let me know if I can tell you anything more to
help you
>>>>>>>> >>>> answer my question.
>>>>>>>> >>>>
>>>>>>>> >>>> Thank You,
>>>>>>>> >>>> Evan Welch
>>>>>>>> >>>> (313) 303 9880
>>>>>>>> >>>
>>>>>>>> >>
>>>>>>>> >
>>>>>>>>
>>>>>>>> <prem_noocean.txt><mineos.pdf>
>>>>>>>
>>>>>>
>>>>>
>>>>
>>>
>>
>
>>> >
>>>
>>
>
More information about the CIG-SEISMO
mailing list