[CIG-SEISMO] [SPECFEM3D_GLOBE]: Solver error when CUT_SUPERBRICK_XI/ETA are true

[Elliott Sales de Andrade]

Hi Dimitri,

In fact, I am also running with SUPPRESS_CRUSTAL_MESH = .false. and 
ADD_4TH_DOUBLING = .false..

But CUT_SUPERBRICK_XI/ETA can still be set to .true. by this code in
read_compute_parameters when NEX/16 is not a multiple of NPROC:

if(mod(NEX_XI,16) /= 0) stop 'NEX_XI must be a multiple of 16'
if(mod(NEX_ETA,16) /= 0) stop 'NEX_ETA must be a multiple of 16'
if(mod(NEX_XI/8,NPROC_XI) /= 0) stop 'NEX_XI must be a multiple of
8*NPROC_XI'
if(mod(NEX_ETA/8,NPROC_ETA) /= 0) stop 'NEX_ETA must be a multiple of
8*NPROC_ETA'
if(mod(NEX_XI/16,NPROC_XI) /=0) CUT_SUPERBRICK_XI = .true.
if(mod(NEX_ETA/16,NPROC_ETA) /=0) CUT_SUPERBRICK_ETA = .true.

which means that NGLOB1D_RADIAL_CORNER is changed in meshfem3D.f90's
meshfem3D_setup_counters.

However, it appears that NGLOB1D_RADIAL_CM, NGLOB1D_RADIAL_OC, and
NGLOB1D_RADIAL_IC in the values_from_mesher.h file are correct. Perhaps
the solver should use those numbers instead?

On Wed, 2012-12-19 at 18:50 +0100, Dimitri Komatitsch wrote:
> Dear Elliott,
> 
> This was written by David Michéa in my group back in 2008 and used to be 
> OK for sure in 2008, 2009 and 2010. I am not sure if side effects from 
> other changes could make it fail in the current version; in principle 
> there is no reason for npoin1D_corner to be different from NGLOB1D_RADIAL.
> This will be a bit difficult to debug because David left around 2009 (he 
> got a permanent job elsewhere).
> 
> Note that in most cases we use SUPPRESS_CRUSTAL_MESH = .false. and 
> ADD_4TH_DOUBLING = .false.; that case should work fine; other 
> combinations are not so well tested. (they were fully tested when we 
> developed them in 2008, but these days I think nobody uses them any more).
> 
> Best,
> Dimitri.
> 
> On 12/19/2012 06:55 AM, Elliott Sales de Andrade wrote:
> > Hi,
> >
> > I have been attempting to run SPECFEM3D_GLOBE on low-resolution settings
> > so that I do not have to use up much of the cluster when testing.
> >
> > I have attached the Par_file I am using, but I believe the relevant
> > parameters are:
> > NCHUNKS                         = 2
> > NEX_XI                          = 112
> > NEX_ETA                         = 112
> > NPROC_XI                        = 2
> > NPROC_ETA                       = 2
> >
> > The mesher runs without error, but the solver prints out the error:
> >   error rank           0 :
> >   incorrect nb of points in corner buffer
> >   Error detected, aborting MPI... proc            0
> >
> > This error arises in `read_arrays_buffers_2` when `npoin1D_corner` is
> > greater than `NGLOB1D_RADIAL`. Looking through the code, it appears that
> > `NGLOB1D_RADIAL` is computed "theoretically" and `npoin1D_corner` is
> > read from the files from the mesher.
> >
> > The mesher calculates `npoin1D_corner` in `meshfem3D_setup_counters`
> > which starts with the values from `NGLOB1D_RADIAL`, and then adds some
> > additional "DIFF" values when either `CUT_SUPERBRICK_XI` or
> > `CUT_SUPERBRICK_ETA` are true. These logicals are set to true when
> > NEX/16 is not a multiple of NPROC (for the default values of
> > SUPPRESS_CRUSTAL_MESH and ADD_4TH_DOUBLING), which occurs for NPROC=2
> > when NEX=80,112,144,176,208 (all of which are "correct" values from
> > Table 3.2).
> >
> > Given that the "DIFF" values are all positive, I'm not sure how the
> > condition in the solver could ever be true when CUT_SUPERBRICK is true.
> > Am I just reading this code incorrectly, or is the check in the solver
> > not right?
> >

-- 
Elliott Sales de Andrade <esalesde at physics.utoronto.ca>
University of Toronto (Physics)