Fails to simulate an acoustic model

[lpk-py]

I replaced with zeros all the vs values in <model>_axisem.bm output by AxiSEM in /MESHER/Diags/, to make an acoustic Earth model. Then assign the file as external model in inparam_mesh.

However, the SOLVER exited due to some errors related to source location. Below is the OUTPUT_<RUNNAME> file on my PC:

    Initialized run for nproc =    2

 MAIN: Welcome to AxiSEM!

     Simulation started on 04/01/2017 at 20h 42min

 MAIN: Reading parameters..................................

        ///////////////////////////////////////////////////////////////
        //                                                           //
        //                  A   x   i   S   E   M                    //
        //                                                           //
        //                                                           //
        //         Parallel spectral-element solution to             //
        //                                                           //
        //           3-D seismic wave propagation for                //
        //                                                           //
        //          axially symmetric background models              //
        //                                                           //
        //               in a spherical 2-D domain                   //
        //                                                           //
        //                                                           //
        //  Authors : Tarje Nissen-Meyer (Oxford University)         //
        //              Martin van Driel (ETH Zurich)                //
        //                 Simon Stähler (LMU Munich)                //
        //                Kasra Hosseini (LMU Munich)                //
        //               Stefanie Hempel (University of Muenster)    //
        //            Alexandre Fournier (IPG Paris)                 //
        //                   Tony Dahlen (Princeton University)      //
        //                                                           //
        //   Contact:     info@axisem.info                           //
        //   Information: www.axisem.info                            //
        //                                                           //

        //                                                           //
        //     If you are publishing results obtained with this      //
        //          code, please cite this paper:                    //
        //                                                           //
        // (1) T. Nissen-Meyer, M. van Driel, S. C. Staehler,        //
        //     K. Hosseini, S. Hempel, L. Auer, A. Colombi           //
        //     and A. Fournier:                                      //
        //     "AxiSEM: broadband 3-D seismic wavefields in          //
        //              axisymmetric media"                          //
        //     Solid Earth, 5, 425-445, 2014                         //
        //     doi:10.5194/se-5-425-2014                             //
        //                                                           //
        //       Comprehensive description of the underlying         //
        //           numerical analysis can be found in:             //
        //                                                           //
        // (2) Tarje Nissen-Meyer, F. A. Dahlen, A Fournier (2007)   //
        //     "Spherical-earth Frechet sensitivity kernels"         //
        //     Geophysical Journal International 168(3),1051-1066.   //
        //     doi:10.1111/j.1365-246X.2006.03123.x                  //
        //                                                           //
        // (3) Tarje Nissen-Meyer, A Fournier, F. A. Dahlen (2007)   //
        //     "A two-dimensional spectral-element method for        //
        //        computing spherical-earth seismograms -            //
        //        I. Moment-tensor source"                           //
        //     Geophysical Journal International 168(3), 1067-1092.  //
        //     doi:10.1111/j.1365-246X.2006.03121.x                  //
        //                                                           //
        // (4) Tarje Nissen-Meyer, A Fournier, F. A. Dahlen (2007)   //
        //     "A two-dimensional spectral-element method for        //
        //        computing spherical-earth seismograms -            //
        //        II.  Waves in solid-fluid media"                   //
        //     Geophysical Journal International 174(3), 873-888.    //
        //     doi:10.1111/j.1365-246X.2008.03813.x                  //
        //                                                           //
        // (5) Martin van Driel and Tarje Nissen-Meyer (2014)        //
        //     "Seismic wave propagation in fully anisotropic        //
        //        axisymmetric media"                                //
        //      Geophysical Journal International 199 (2): 880-893.  //
        //      doi: 10.1093/gji/ggu269                              //
        //                                                           //
        // (6) Martin van Driel and Tarje Nissen-Meyer (2014)        //
        //     "Optimized visco-elastic wave propagation for         //
        //        weakly dissipative media"                          //
        //      Geophysical Journal International 199 (2): 1078-1093.//
        //      doi: 10.1093/gji/ggu314                              //
        //                                                           //
        //                                                           //
        //  April 2016: version 1.3                                  //
        //                                                           //
        ///////////////////////////////////////////////////////////////

        =============  I N P U T    P A R A M E T E R S ===============
            Data I/O path:                      ./Data              
            Info I/O path:                      ./Info              
            Simulation length [s]:               3600.000
            Enforced time step [s]:               0.000
            Enforced source period [s]:           0.000
            Simulation type:                    single
            Receiver file type:                 colatlon
            Sum seismograms?                     F
            Sum wavefields?                      F
            Time extrapolation scheme:          symplec4
            Seismogram sampling rate [s]:         0.000
            Dump kin./pot. energy?               F
            Dump global snaps?                   F
            Dump strain?                         F
            Wavefield dumping type:             displ_only  
            First GLL to save in strains:        1
            Last GLL to save in strains:         3
            First GLL to save in strains:        1
            Last GLL to save in strains:         3
            Samples per period for strains:       4.000
            Source dumping type:                igno
            Add heterogeneous region?            F
            Perform extensive mesh tests?        F
            Output format (seism., wavefields): netcdf
        ===============================================================
 MAIN: Reading mesh database...............................

   General numerical input/output parameters================
     grid pts/wavelngth =   1.5000000000000000     
     source period [s]  =   50.000000000000000     
     courant number     =  0.60000002384185791     
     time step [s]      =  0.14321715989718831     
 Checking value order in external model
  Value ' rad'  is in column  1
  Value ' vpv'  is in column  3
  Value ' vsv'  is in column  4
  Value ' rho'  is in column  2
  Value ' qka'  is in column  5
  Value ' qmu'  is in column  6
 Model in file ./external_model.bm has   160 layers and is isotropic and anelastic...
 Depths/radii are assumed to be defined in meters
 Layers in file ./external_model.bm start at surface
 Checking for discontinuities in the external velocity model
  1st order disc. at radius   6356000.0, layer:     2
  1st order disc. at radius   6346600.0, layer:     4
  1st order disc. at radius   6291000.0, layer:    10
  1st order disc. at radius   6151000.0, layer:    16
  1st order disc. at radius   5971000.0, layer:    22
  1st order disc. at radius   5771000.0, layer:    28
  1st order disc. at radius   5701000.0, layer:    34
  1st order disc. at radius   5600000.0, layer:    40
  1st order disc. at radius   3630000.0, layer:    81
  1st order disc. at radius   3480000.0, layer:    87
  1st order disc. at radius   1221500.0, layer:   134
 External model has  12 discontinuities
 
 Creating interpolation objects
    idom, upper_layer, lower_layer,      r(ul),      r(ll)
       1            1            2   6371000.0   6356000.0
       2            3            4   6356000.0   6346600.0
       3            5           10   6346600.0   6291000.0
       4           11           16   6291000.0   6151000.0
       5           17           22   6151000.0   5971000.0
       6           23           28   5971000.0   5771000.0
       7           29           34   5771000.0   5701000.0
       8           35           40   5701000.0   5600000.0
       9           41           81   5600000.0   3630000.0
      10           82           87   3630000.0   3480000.0
      11           88          134   3480000.0   1221500.0
      12          135          160   1221500.0         0.0
 

   Background model=========================================
     bkgrdmodel = external
     radius [m] =    6371000.0000000000     
     have_fluid =  T

   Min/max grid spacing=====================================
     hmin (global) [m]   :    1623.1277476708324     
     hmax (global) [m]   :    92925.020380043366     
     min_distance_dim [m]:    162.31277476708325     

   Axialogy=================================================

     Global total axial elements:          84
     Global solid axial elements:           0
     Global fluid axial elements:          84


 MAIN: Initializing grid...................................
 MAIN: Starting wave preparation...........................
 
 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 ++++++++    SEISMIC WAVE PROPAGATION: SOLID-FLUID CASE  ++++++++
 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
   *****************GIVEN SOURCE PARAMETERS*****************
      Magnitude [Nm]:             1.000E+20
      Excitation type:        monopole    mrr       
      Depth [km]:                 0.000E+00
      Colat. [deg]:               0.000E+00
      Long. [deg]:                0.000E+00
      Source time function:       gauss_0
      Dom. period mesh [s]:         50.0000
   *********************************************************
 
    construction of mapping for kwf output...
    ...solid part...
    ...fluid part...
 local point number:               38400
 after removing duplicates:        24937
 compression:                 0.649401069    
    .... constructing midpoint grid for kwf output
    .... constructing finite element grid for kwf output
    .... constructing spectral element grid for kwf output
    .... finished construction of mapping for kwf output

     Using period of the mesh:   50.000000000000000     

     desired simulation length  :  3600.00   seconds
     offered simulation length  :  3600.05   seconds
     number time loop iterations:   16758

     desired seismogram sampling:     0.00   seconds
     offered seismogram sampling:     0.21   seconds
     ...that is, every          :       1 timesteps
     number of samples          :   16759

     Number of snapshots        :      61
     ...approximately every     :   60.00   seconds
     ...that is, every          :     279 timesteps
 
   *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-
   SOURCE TIME FUNCTION: gauss_0
    coarsest dump every           1 th time step, dt:  0.21482573984578246     
 
   SHIFT FACTOR of source time function [s]:   75.1890106    
    # SEM, seis, coarse points per shift factor:   350.000000       350.000000       350.000000    
    # simul. half widths per shift factor:   1.50378025    
   *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-
 

  ::::::::: BACKGROUND MODEL & PRECOMPUTED MATRICES:::::::
 WARNING: boundary term not all that precise!
  Term should equal four for 2 boundaries (i.e. int (sin) )
  Actual numerical value:   0.0000000000000000     
  :::::::DONE BACKGROUND MODEL & PRECOMPUTED MATRICES:::::
  *****************************************
     Welcome to the source term calculator 
  *****************************************
  locating the source....
ERROR: Source should be on the axis, hence theta = 0, but:
ERROR: Source should be on the axis, hence theta = 0, but:
At line 494 of file source.f90 (unit = 6, file = 'stdout')
Fortran runtime error: Expected INTEGER for item 2 in formatted transfer, got CHARACTER
(a,/,i7,2i2,/,a,3e10.2)
     ^
At line 494 of file source.f90 (unit = 6, file = 'stdout')
Fortran runtime error: Expected INTEGER for item 2 in formatted transfer, got CHARACTER
(a,/,i7,2i2,/,a,3e10.2)
     ^
-------------------------------------------------------
Primary job  terminated normally, but 1 process returned
a non-zero exit code.. Per user-direction, the job has been aborted.
-------------------------------------------------------
--------------------------------------------------------------------------
mpirun detected that one or more processes exited with non-zero status, thus causing
the job to be terminated. The first process to do so was:

  Process name: [[6601,1],0]
  Exit code:    2
--------------------------------------------------------------------------

On the cluster, the SOLVER also raise errors due to source (OUTPUT tells that multiple processes have the source).

[lpk-py]

Is this community still active? I noticed the last post was in Nov. 2015.

[martinvandriel]

I think fully fluid planets are currently not supported, because we did not yet implement the free surface (which is implicit for the elastic case). I guess we also assume the source is in the elastic part of the planet, hence the error you are getting. Both should not be terribly hard to add though.

[lpk-py]

@martinvandriel Thanks for the clarification. Look forward to the new features in next versions. If you are planning the support to the fluid Earth, I suggest to add a controlling parameter in inparam_mesh to make the builtin models to be acoustic automatically.

BTW, is there a limit for the parallel subprocesses, related to the dominant period of meshes? When I tried to divide the model into slices, there could be an error in the OUTPUT of MESHER

   define search sflobal index....
   partition sflobal index....
 ERORR: having more then 8 neighbours (+myself)
        check mesh decomposition in the solid)
        iproct =            0

or sometimes

   define search sflobal index....
 ERROR: too many neighbours in solid