[cig-commits] r19920 - seismo/3D/SPECFEM3D_GLOBE/trunk

[dkomati1 at geodynamics.org]

Author: dkomati1
Date: 2012-04-04 15:22:43 -0700 (Wed, 04 Apr 2012)
New Revision: 19920

Modified:
   seismo/3D/SPECFEM3D_GLOBE/trunk/todo_list_please_dont_remove.txt
Log:
updated the todo list


Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/todo_list_please_dont_remove.txt
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/todo_list_please_dont_remove.txt	2012-04-04 16:09:11 UTC (rev 19919)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/todo_list_please_dont_remove.txt	2012-04-04 22:22:43 UTC (rev 19920)
@@ -4,117 +4,33 @@
 
 Things that could be done in a future version:
 
+- fix the 3D attenuation memory size problem by making it optional and off by default.
+To do that, the following command will be very useful:
+nm --print-size --size-sort --reverse-sort --radix=d ./bin/xspecfem3D | cut -f 2- -d ' ' | less 
+
 - check if sigma_xz is inverted with sigma_zx and so on in compute_forces; usually this does not matter because the stress tensor is symmetric, but for Cosser
 at or for some C-PML formulations this could matter
 
-- use a more accurate time scheme: either symplectic schemes (but there could be a problem in the viscoelastic case, see below) or else Runge Kutta RK4, which is non staggered and would thus solve the approximate implementation of the displacement gradient that I used to implement attenuation back in 1999. It is not clear if symplectic schemes could be suitable for viscoelastic media (Tarje implemented them for elastic media only), but RK4 is definitely suitable (Roland Martin used it successfully in 2D viscoelastic).
+- implement RK4 and LDDRK4-6 high-order time schemes, including for viscoelastic media: will be done by Zhinan Xie
 
-Subject: Re: SEM attenuation
-Date: Fri, 23 Jul 2010 12:16:40 +0200
-From: Tarje Nissen-Meyer <tarjen at ethz.ch>
-Organization: ETH Zurich
-To: Dimitri Komatitsch <dimitri.komatitsch at univ-pau.fr>
-CC: <yingz at vt.edu>, Jeroen Tromp <jtromp at princeton.edu>, Min Chen <mchen at gps.caltech.edu>, Vala Hjorleifsdottir <vala at ldeo.columbia.edu>, "Brian Savage" <savage at uri.edu>, Shiann-Jong Lee <sjlee at earth.sinica.edu.tw>, "Roland Martin" <roland.martin at univ-pau.fr>, Bernhard Schuberth   <mail at bernhard-schuberth.de>, Carl Tape <carltape at fas.harvard.edu>, "Anne Sieminski" <anne.sieminski at obs.ujf-grenoble.fr>, Paul Friberg   <p.friberg at isti.com>, Kasper van Wijk <kasper at cgiss.boisestate.edu>, "Dylan Mikesell" <dmikesell at cgiss.boisestate.edu>, Federica Magnoni  <federica.magnoni at ingv.it>, Ebru Bozdag <bozdag at princeton.edu>, Hejun Zhu     <hejunzhu at princeton.edu>, Pieyre Le Loher <pieyre.le_loher at inria.fr>,   Christina Morency <cmorency at Princeton.EDU>, Emanuele Casarotti    <emanuele.casarotti at gmail.com>, Piero Basini <basini at bo.ingv.it>, "Emmanuel Chaljub" <Emmanuel.Chaljub at obs.ujf-grenoble.fr>, Qinya Liu    <liuqy at physics.utoronto.ca>, Yang Luo <yangl at Princeton.EDU>, Daniel Peter     <dpeter at Princeton.EDU>
+- fix the RK4 attenuation implementation problem from my 1999 approximation: will be done by Zhinan Xie
 
-Hi all,
+- implement C-PML: will be done by Zhinan Xie
 
-Just to add that symplectic schemes are based on conservative systems
-(they approximate the Hamiltonian) and we don't know yet how they behave
-for dissipative media. But still of course something to be added and
-looked at. Should also give credit to Jean-Paul Ampuero who really first
-suggested trying these schemes for elastodynamics with the SEM.
+- use PT-Scotch to be able to switch from 6 N^2 MPI threads to any other number, and feed that to SPECFEM3D instead of SPECFEM3D_GLOBE: will be done by Joseph Charles
 
-Best regards,
-Tarje
+=============================
 
-On 23/07/10 02:21, Dimitri Komatitsch wrote:
->
-> Hi Ying,
->
-> What I know for sure if that back in 1999 when I developed the time
-> scheme for attenuation I used a trick that made implementation much
-> easier but that also makes the 4th order Runge Kutta (RK4) time scheme
-> become second order only, i.e. RK2 instead of RK4. This means that for
-> very long simulations (for instance multi-orbit surface waves)
-> attenuation can become inaccurate because of a lack of accuracy in the
-> time scheme (in which case reducing Delta_t purposely solves the
-> problem, but of course makes the simulation more expensive).
->
-> A way of solving this problem could be to switch to better time
-> schemes such as the symplectic time scheme that Tarje introduced in a
-> GJI paper a few years ago. Tarje and I should probably implement that
-> in the code one day... (therefore I cc him; we talked about this at
-> the EGU meeting in May)
->
-> Thank you,
-> Dimitri.
+Less important for now:
+----------------------
 
-- use RK4 or symplectic time scheme (will be done by Tarje Nissen-Meyer and/or Jean-Paul Ampuero) (would be useful in the 2D version of the code as well):
-Hi Jeroen, Perfect. I think talking to Jean-Paul Ampuero would be useful
-as well because in Utrecht last year he had told us that
-he had implemented some nice 4th-order symplectic schemes
-in his version of SEM2D. Dimitri.
-Jeroen Tromp wrote:
-> Hi Dimitri:
-> This is one of the first things Tarje and I plan to work on after he
-> arrives.
-> Jeroen
-> Dimitri Komatitsch wrote:
->> Hi Jeroen,
->> I think the last important thing that is missing
->> in SPECFEM3D (and SPECFEM2D) is a fourth-order time scheme.
->> I think both Jean-Paul Ampuero and Tarje Nissen-Meyer
->> have worked on this, they will both be at Caltech
->> soon therefore maybe they could take care of adding it?
->> This would definitely increase the accuracy of very long
->> simulations (e.g. multi-orbit surface waves).
->> Dimitri.
-
-- From Dimitri: there is something in SPECFEM3D_GLOBE that we noticed a few years ago
-regarding attenuation but never fixed: on page 813 of our 1999 paper
-I used a trick suggested by Robertsson et al. (1994)
-to use a non-staggered Runge-Kutta (RK4) scheme for the attenuation
-equations while using a staggered finite-difference (Newmark)
-time scheme for all the other equations.
-(combining staggered and non-staggered formulations being difficult)
-The trick works fine in most cases but there is a hidden problem for very long
-simulations, and in particular for multi-orbit surface waves:
-because of that trick, the RK4 is not really a real Runge-Kutta scheme
-(because grad(displacement) is used as a source but not known half way between
-time steps, therefore I replaced it with an average between
-t and t + Delta_t, which implies that the approximation is not fourth-order
-accurate any more because of that smoothing/interpolation). Therefore
-dispersion due to attenuation is not very accurately computed in the case of
-very long runs; that matters mostly for surface waves.
-We should fix that one day. The only thing to do would be to design a better time
-integration scheme for the attenuation equation, without that trick; the rest
-(Newmark etc) is fine and does not need to change.
-Daniel Peter did not fix that in 2011; he only fixed the old implementation in SPECFEM3D which used a
-fixed period absorption band to an adaptive one, which uses the resolved
-periods of the mesh (as is done in the global version). The time stepping of the memory variables did not change.
-Thus this IMPORTANT KNOWN BUG remains: in compute_forces_viscoelastic.f90 the calculation of
-attenuation (viscoelasticity) is slightly incorrect because the gradients
-are computed twice but *at the same time step* instead of at different
-(staggered) time steps (t_{n-1} and t_n or something like that).
-That's easy to fix but I have no time for now. Let us fix that in the future.
-It will only make a very small difference in the final seismograms (except for very long runs,
-e.g., for multi-orbit surface waves) therefore the current code with the bug can be used without any major problem.
-
-Reply from Daniel Peter (Princeton University, USA) on February 23, 2010:
-for higher orbit arrivals, we do might want to consider the symplectic time
-schemes as well. Tarje made some tests and showed some "shocking" results
-for longer distance paths. However, I am not sure if his scheme currently
-considers purely elastic simulations, without incorporating attenuation.
-It would be nice having both improved formulations included for those cases.
-This seems to be a nice outlook for a future version update.
-
 - use a better checkpointing/restarting system with CRC-32 as in SPECFEM3D_GLOBE/tags/v4.1.0_beta_merged_mesher_solver_non_blocking_MPI/src
 
 - use a potential of (rho * u) instead of u in the fluid, in case of fluid-fluid discontinuities
 
 - merge SPECFEM3D_GLOBE into SPECFEM3D, i.e. create a single package and see the mesh of the Earth as a general mesh that we would decompose in parallel using PT-Scotch. Problem: would this scale for huge global Earth simulations, e.g. at a seismic period of 1 second on 200,000 processor cores?
 
-- make the code compatible with helioseismology / general Cowling formulation (could be done with Tarje Nissen-Meyer)
+- make the code compatible with helioseismology / general Cowling formulation
 
 - could be done by Vala:
 we could use heuristic rules to make source and receiver detection much
@@ -140,26 +56,6 @@
 ---------------------------------------------------------------------------------
 ---------------------------------------------------------------------------------
 
-Things to do later:
--------------------
-
-- afficher le modele PREM,les courbes de nombre de points par longueur d'onde et de stabilite et les worst elements (stability et nb pts per lambda) en OpenDX en serial sur le maitre dans le nouveau mesher, en copiant/collant ce qu'il faut depuis display_prem_sampling_doubling.f90 sur le master seulement (rank == 0)
-
----------------------------------------------------------------------------------
----------------------------------------------------------------------------------
----------------------------------------------------------------------------------
----------------------------------------------------------------------------------
-
-Things to add to the manual:
-----------------------------
-
-- ajouter les courbes de dispersion et stabilite au manuel pour tous les NEX classiques (160, 256, 380 etc) dans un Appendix special la-dessus; dire dans cet appendix que dans le outer core ce n'est bien sur pas Vs qui est represente mais Vp / 1.25 (on a pris 25 % de marge car on travaille en velocity potential dans le outer core, ce qui demande un echantillonnage un peu meilleur)
-
----------------------------------------------------------------------------------
----------------------------------------------------------------------------------
----------------------------------------------------------------------------------
----------------------------------------------------------------------------------
-
 Done:
 -----