[CIG-SEISMO] How to choose time step and other parameters to avoid floating-point error
Chunyang Xiao
u6488884 at anu.edu.au
Mon May 7 06:55:30 PDT 2018
Hi Dimitri,
Thank you for your reply. I checked the line of 'Max stability for P wave velocity...' as you told me to, and found it was over 0.5.
After I lowered my time step, the problem was solved! 😊 Your suggestions are really helpful! Thank you for your assistance.
Best regards,
Charles
________________________________
From: Dimitri Komatitsch <komatitsch at lma.cnrs-mrs.fr>
Sent: Friday, 4 May 2018 5:14:25 AM
To: cig-seismo at geodynamics.org; Chunyang Xiao
Subject: Re: [CIG-SEISMO] How to choose time step and other parameters to avoid floating-point error
Hi Charles,
Thanks! That is good news. Now it is probably simply that your time step
is too high and you are above the CFL stability limit. Could you check
what the following line says:
Max stability for P wave velocity (must be below about 0.50 or so) = ...
and check that the number is below 0.5?
Best regards,
Dimitri.
On 05/03/2018 05:18 AM, Chunyang Xiao wrote:
> Hi Dimitri,
>
>
> Thank you so much for identifying the problem in my code. I didn't
> realize that the error message '0x4192BA in read_regions_ at
> read_regions.f90:99' means that I have problems in line 99 of
> read_regions.f90, and thank you for telling me.
>
>
> I have exactly the error you said: I set Vp = 1, and Vs = 1 in the
> example Rayleigh_wave_no_crack. After I set Vp = 2, the error was gone.
>
>
> Now as I changed my materials to be
>
>
> 1 1 2000.d0 3000.d0 2000.d0 0 0 9999.d0 9999.d0 0 0 0 0 0 0
> 2 1 1.d0 3000.d0 0.d0 0 0 9999.d0 9999.d0 0 0 0 0 0 0
>
>
> An error came again, which said
>
>
> #0 0x7F5BFED14E08
> #1 0x7F5BFED13F90
> #2 0x7F5BFE4454AF
> #3 0x4489B6 in compute_forces_viscoelastic_ at
> compute_forces_viscoelastic.F90:601
> #4 0x44E240 in compute_forces_viscoelastic_main_ at
> compute_forces_viscoelastic_calling_routine.F90:67
> #5 0x4936C7 in iterate_time_ at iterate_time.F90:165
> #6 0x40230A in specfem2d at specfem2D.F90:381
> #7 0x7F5BFE43082F
>
>
> Line 601 of compute_forces_viscoelastic.F90 says 'tempx1(i,j) =
> jacobianl * (sigma_xx*xixl + sigma_zx*xizl) ! this goes to accel_x'. I
> couldn't figure out what that means. I am thankful for any input or
> suggestion about it.
>
>
> I am also curious about how to find max stability and output_solver.txt
> because floating point errors sometimes also happened when I changed the
> time step.
>
>
> Thank you again for your assistance.
>
> Sincerely,
> Charles Xiao
> ------------------------------------------------------------------------
> *From:* Dimitri Komatitsch <komatitsch at lma.cnrs-mrs.fr>
> *Sent:* Thursday, 3 May 2018 7:24 AM
> *To:* cig-seismo at geodynamics.org; Chunyang Xiao
> *Subject:* Re: [CIG-SEISMO] How to choose time step and other parameters
> to avoid floating-point error
>
> Hi all,
>
> I have added the stop statement.
>
> Best regards,
> Dimitri.
>
> On 05/02/2018 11:20 PM, Dimitri Komatitsch wrote:
>>
>> Hi Charles,
>>
>> Thanks for your message!
>> Line 99 of src/meshfem2D/read_regions.f90 is this:
>>
>> poisson_ratio = 0.5d0*(vpregion*vpregion -
>> 2.d0*vsregion*vsregion) / (vpregion*vpregion - vsregion*vsregion)
>>
>> thus it is just an error in your input file (Par_file), you have at
>> least one material for which Vs = Vp and thus you get a division by zero.
>>
>> Materials always have a value of Poisson's ratio > -1 (in practice, for
>> natural materials, >= 0), thus for all natural materials you must have
>> Vs < Vp / sqrt(2).
>>
>> I will add a stop statement in the code, right before that line, to
>> detect that.
>>
>> Best regards,
>> Dimitri.
>>
>> On 05/01/2018 02:33 PM, Chunyang Xiao wrote:
>>> Hi there,
>>>
>>>
>>> I am an undergraduate student majoring in mechanical engineering.
>>> Recently when I changed the speed of P-wave (Vp from the section
>>> 'velocity and density models' in Par_file) in SPECFEM2D, an error
>>> occurred.
>>>
>>>
>>> Program received signal SIGFPE: Floating-point exception - erroneous
>>> arithmetic operation.
>>>
>>> Backtrace for this error:
>>>
>>> #0 0x7F209C4FDE08
>>>
>>> #1 0x7F209C4FCF90
>>>
>>> #2 0x7F209C14D4AF
>>>
>>> #3 0x4192BA in read_regions_ at read_regions.f90:99
>>>
>>> #4 0x418CEE in read_parameter_file_ at read_parameter_file.F90:88
>>>
>>> #5 0x42521C in MAIN__ at meshfem2D.F90:403
>>>
>>> Floating point exception (core dumped)
>>>
>>>
>>> This happened a few times also when changed other parameters, such as
>>> time step. A similar error happened to Mr. Moritz Wehler before, and
>>> Dimitri's answer on Wed May 27 07:40:46 PDT 2015 was
>>>
>>>
>>> It seems you simply forgot to also divide the time step by two to
>>> remain
>>>
>>> stable when you increased mesh resolution by two:
>>>
>>> *** Max stability for P wave velocity (must be below about 0.50
>>> or so)
>>>
>>> = 0.972936358376379
>>>
>>>
>>> I was curious where we could see the 'Max stability for P wave
>>> velocity' for a given case. As I turned to the SPECFEM2D user manual
>>> (version 7.0) for help, I saw chapter 4.6 (How to choose the time
>>> step) say:
>>>
>>>
>>> You can directly compare these values with the value given in
>>> sentence ‘_Max stability for P wave velocity_’ in file
>>> _output_solver.txt _to see whether you set the correct ∆t in
>>> Par_file or not. For elastic simulation, the CFL value given in
>>> _output_solver.txt _does not consider the V p /V s ratio......
>>>
>>>
>>> That makes me more confused because I can find neither the sentence
>>> ‘Max stability for P wave velocity’ nor the file output_solver.txt.
>>> The file named output_solver.txt only exists in the folders of a few
>>> examples I did not use or run.
>>>
>>> Also, in Ms. Parvathy Sindhu Nair's question on Mon Feb 15 14:25:50
>>> PST 2016, she and Dimitri mentioned CFL, which seems to matter but
>>> does not appear in my Par_file.
>>>
>>> this is either the time step that is a bit too high (unlikely, since
>>> the CFL you used (0.44) seems fine), or if the amplitude of the
>>> source you used is very high......
>>>
>>>
>>> Would anyone kindly let me know how I can choose appropriate
>>> parameters (time step, speed, etc.) to avoid such floating-point
>>> errors and where I can see my CFL, output_solver.txt, and max
>>> stability for P wave?
>>>
>>> Sincerely,
>>> Charles Xiao
>>>
>>>
>>> _______________________________________________
>>> CIG-SEISMO mailing list
>>> CIG-SEISMO at geodynamics.org
>>> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-seismo
>>>
>>
>
> --
> Dimitri Komatitsch, CNRS Research Director (DR CNRS)
> Laboratory of Mechanics and Acoustics, Marseille, France
> http://komatitsch.free.fr
>
>
> _______________________________________________
> CIG-SEISMO mailing list
> CIG-SEISMO at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-seismo
>
--
Dimitri Komatitsch, CNRS Research Director (DR CNRS)
Laboratory of Mechanics and Acoustics, Marseille, France
http://komatitsch.free.fr
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/cig-seismo/attachments/20180507/f1e00e8c/attachment-0001.html>
More information about the CIG-SEISMO
mailing list