def jacobian(self, coords, system, strip_dummies=None, converter_options=None, analytic_deriv_order=None, **kwargs):
    if (converter_options is None):
        converter_options = {}
    merged_convert_options = dict(self.converter_options, **converter_options)
    try:
        strip_dummies = (merged_convert_options['strip_dummies'] if (strip_dummies is None) else strip_dummies)
    except KeyError:
        strip_dummies = False
    try:
        analytic_deriv_order = (merged_convert_options['analytic_deriv_order'] if (analytic_deriv_order is None) else analytic_deriv_order)
    except KeyError:
        analytic_deriv_order = 0
    if strip_dummies:
        dummies = self.molecule.dummy_positions
        if (len(dummies) == 0):
            dummies = None
    else:
        dummies = None
    if (dummies is not None):
        main_excludes = np.setdiff1d(np.arange(self.molecule.num_atoms), dummies)
    else:
        main_excludes = None
    jacs = super().jacobian(coords, system, analytic_deriv_order=analytic_deriv_order, converter_options=converter_options, **kwargs)
    raw_jacs = []
    for (n, j) in enumerate(jacs):
        baseline = ((2 * analytic_deriv_order) + len(coords.shape))
        ext_dim = (j.ndim - baseline)
        shp = (sum((((j.shape[i] // 3), 3) for i in range(ext_dim)), ()) + j.shape[(- baseline):])
        j = j.reshape(shp)
        if (dummies is not None):
            for i in range(ext_dim):
                j = np.take(j, main_excludes, axis=(2 * i))
            for i in range(analytic_deriv_order):
                j = np.take(j, main_excludes, axis=((- 2) * (i + 2)))
        if (len(coords.shape) > 2):
            j = np.moveaxis(j, (- 3), 0)
        raw_jacs.append(j)
    jacs = raw_jacs
    return jacs
def jacobian(self, coords, system, strip_dummies=None, converter_options=None, analytic_deriv_order=None, **kwargs):
    if (converter_options is None):
        converter_options = {}
    merged_convert_options = dict(self.converter_options, **converter_options)
    try:
        strip_dummies = (merged_convert_options['strip_dummies'] if (strip_dummies is None) else strip_dummies)
    except KeyError:
        strip_dummies = False
    try:
        analytic_deriv_order = (merged_convert_options['analytic_deriv_order'] if (analytic_deriv_order is None) else analytic_deriv_order)
    except KeyError:
        analytic_deriv_order = 0
    if strip_dummies:
        dummies = self.molecule.dummy_positions
        if (len(dummies) == 0):
            dummies = None
    else:
        dummies = None
    if (dummies is not None):
        main_excludes = np.setdiff1d(np.arange(self.molecule.num_atoms), dummies)
    else:
        main_excludes = None
    jacs = super().jacobian(coords, system, analytic_deriv_order=analytic_deriv_order, converter_options=converter_options, **kwargs)
    raw_jacs = []
    for (n, j) in enumerate(jacs):
        baseline = ((2 * analytic_deriv_order) + len(coords.shape))
        ext_dim = (j.ndim - baseline)
        shp = (sum((((j.shape[i] // 3), 3) for i in range(ext_dim)), ()) + j.shape[(- baseline):])
        j = j.reshape(shp)
        if (dummies is not None):
            for i in range(ext_dim):
                j = np.take(j, main_excludes, axis=(2 * i))
            for i in range(analytic_deriv_order):
                j = np.take(j, main_excludes, axis=((- 2) * (i + 2)))
        if (len(coords.shape) > 2):
            j = np.moveaxis(j, (- 3), 0)
        raw_jacs.append(j)
    jacs = raw_jacs
    return jacs