def METHOD_NAME(self, atoms, outdir, grid_dimensions, **kwargs):
    """Perform actual bispectrum calculation."""
    use_fp64 = kwargs.get("use_fp64", False)
    lammps_format = "lammps-data"
    ase_out_path = os.path.join(outdir, "lammps_input.tmp")
    ase.io.write(ase_out_path, atoms, format=lammps_format)
    nx = grid_dimensions[0]
    ny = grid_dimensions[1]
    nz = grid_dimensions[2]
    # Create LAMMPS instance.
    lammps_dict = {}
    lammps_dict["twojmax"] = self.parameters.bispectrum_twojmax
    lammps_dict["rcutfac"] = self.parameters.bispectrum_cutoff
    lammps_dict["atom_config_fname"] = ase_out_path
    lmp = self._setup_lammps(nx, ny, nz, outdir, lammps_dict,
                             log_file_name="lammps_bgrid_log.tmp")
    # An empty string means that the user wants to use the standard input.
    # What that is differs depending on serial/parallel execution.
    if self.parameters.lammps_compute_file == "":
        filepath = __file__.split("bispectrum")[0]
        if self.parameters._configuration["mpi"]:
            if self.parameters.use_z_splitting:
                self.parameters.lammps_compute_file = \
                    os.path.join(filepath, "in.bgridlocal.python")
            else:
                self.parameters.lammps_compute_file = \
                    os.path.join(filepath,
                                 "in.bgridlocal_defaultproc.python")
        else:
            self.parameters.lammps_compute_file = \
                os.path.join(filepath, "in.bgrid.python")
    # Do the LAMMPS calculation.
    lmp.file(self.parameters.lammps_compute_file)
    # Set things not accessible from LAMMPS
    # First 3 cols are x, y, z, coords
    ncols0 = 3
    # Analytical relation for fingerprint length
    ncoeff = (self.parameters.bispectrum_twojmax + 2) * \
             (self.parameters.bispectrum_twojmax + 3) * (self.parameters.bispectrum_twojmax + 4)
    ncoeff = ncoeff // 24   # integer division
    self.fingerprint_length = ncols0+ncoeff
    # Extract data from LAMMPS calculation.
    # This is different for the parallel and the serial case.
    # In the serial case we can expect to have a full bispectrum array at
    # the end of this function.
    # This is not necessarily true for the parallel case.
    if self.parameters._configuration["mpi"]:
        nrows_local = extract_compute_np(lmp, "bgridlocal",
                                         lammps_constants.LMP_STYLE_LOCAL,
                                         lammps_constants.LMP_SIZE_ROWS)
        ncols_local = extract_compute_np(lmp, "bgridlocal",
                                         lammps_constants.LMP_STYLE_LOCAL,
                                         lammps_constants.LMP_SIZE_COLS)
        if ncols_local != self.fingerprint_length + 3:
            raise Exception("Inconsistent number of features.")
        snap_descriptors_np = \
            extract_compute_np(lmp, "bgridlocal",
                               lammps_constants.LMP_STYLE_LOCAL, 2,
                               array_shape=(nrows_local, ncols_local),
                               use_fp64=use_fp64)
        lmp.close()
        # Copy the grid dimensions only at the end.
        self.grid_dimensions = [nx, ny, nz]
        return snap_descriptors_np, nrows_local
    else:
        # Extract data from LAMMPS calculation.
        snap_descriptors_np = \
            extract_compute_np(lmp, "bgrid", 0, 2,
                               (nz, ny, nx, self.fingerprint_length),
                               use_fp64=use_fp64)
        lmp.close()
        # switch from x-fastest to z-fastest order (swaps 0th and 2nd
        # dimension)
        snap_descriptors_np = snap_descriptors_np.transpose([2, 1, 0, 3])
        # Copy the grid dimensions only at the end.
        self.grid_dimensions = [nx, ny, nz]
        if self.parameters.descriptors_contain_xyz:
            return snap_descriptors_np, nx*ny*nz
        else:
            return snap_descriptors_np[:, :, :, 3:], nx*ny*nz