Molecular Vibrations with Tako Script

tako.vibration(structure, settings) runs the molecular finite-difference vibration workflow and resolves to the worker result object. It is the Script interface for a finite molecule or non-periodic cluster. Periodic lattice dynamics belongs to tako.phonon; the scientific distinction is developed in Molecular Vibrations and IR/Raman Spectroscopy.

This chapter assumes that the input structure already represents the stationary point to be classified. tako.vibration does not optimize coordinates and does not automatically decide whether a signed-negative frequency is translational noise or an internal instability.

Signature and calculator boundary

const result = await tako.vibration(structure, {
  calculator,
  properties,
  vibrationDisplacementAngstrom,
  vibrationNfree,
  spectroscopy,
});

The required calculator is constructed with tako.calculator.gfn2(...) or tako.calculator.mlip(...). Validate it before submitting work:

const calculator = tako.calculator.gfn2({
  dispersion: 'd4',
  charge: 0,
  unpaired: 0,
});

const check = calculator.validate('vibration');
if (!check.ok) {
  throw new Error(check.reasons.join('\n'));
}

Both calculator families can return frequencies and normal modes. IR intensities and Raman activities require the GFN2 response path. An MLIP run still completes if infrared or raman is present in spectroscopy, but the worker logs the unsupported response and omits the corresponding result object. A script that requires spectra must enforce the calculator family rather than merely checking that the request array contains a name.

Settings

The current public call does not expose a subset of atom indices. Instead, the backend displaces every atom that is not completely fixed by the structure’s constraint mask. A partial Cartesian constraint does not remove the atom: all three coordinates remain in the active atom set while constraint-adjusted force components enter the derivative. Thus complete-atom constraints can change the Hessian dimension, and partial constraints can change its matrix elements; neither is an optional post-processing filter. Record them with the input structure.

For $M$ active atoms, a frequency pass contains one undisplaced evaluation and $6M$ displaced evaluations for vibrationNfree: 2, or one plus $12M$ for 4. IR and Raman use their own response finite-difference passes. Avoid requesting response channels that will not be used.

A guarded minimum calculation

const structure = await tako.input('input/optimized-water.extxyz');
const calculator = tako.calculator.gfn2({
  dispersion: 'd4',
  charge: 0,
  unpaired: 0,
});

const validation = calculator.validate('vibration');
if (!validation.ok) throw new Error(validation.reasons.join('\n'));

const result = await tako.vibration(structure, {
  calculator,
  vibrationDisplacementAngstrom: 0.01,
  vibrationNfree: 2,
  spectroscopy: ['frequencies', 'infrared'],
});

if (result.kind !== 'vibration') {
  throw new Error(`Unexpected result kind: ${result.kind}`);
}
if (!result.vibrations?.frequencies_cm_inv?.length) {
  throw new Error('The vibration result contains no frequencies.');
}
if (!result.infrared) {
  throw new Error('IR response was required but not returned.');
}

await tako.output('output/vibration-result.json', result);

The returned object is authoritative for the call. tako.output writes the object to the path chosen by the script; it does not reproduce the interface calculation folder’s derived trajectory files automatically. If another program needs a stable reduced schema, write one explicitly rather than relying on every backend detail:

const summary = {
  level_of_theory: result.level_of_theory,
  engine: result.engine,
  frequencies_cm_inv: result.vibrations.frequencies_cm_inv,
  zero_point_energy_ev: result.vibrations.zero_point_energy_ev,
  ir_intensities_au: result.infrared?.intensities_au ?? null,
  raman_activities: result.raman?.activities ?? null,
};

await tako.output('output/vibration-summary.json', summary);

Result hierarchy

Every successful call returns:

result
├── kind = "vibration"
├── level_of_theory
├── engine
└── vibrations
    ├── frequencies_cm_inv
    ├── frequency_components
    ├── energies_hartree
    ├── zero_point_energy_hartree
    ├── zero_point_energy_ev
    └── modes

When GFN2 IR response succeeds, result.infrared contains a second vibration payload under infrared.vibrations, plus intensities_au, static_dipole_au, and max_force_hartree_per_bohr. The first name is historical—the values are in (D/Å)² amu⁻¹—and static_dipole_au is a scalar norm, not a three-component vector. When GFN2 Raman response succeeds, result.raman contains raman.vibrations, activities, absolute_intensities, and the per-mode invariants alpha2, gamma2, and delta2. In the current backend the two intensity arrays are identical $45\alpha^2+7\gamma^2$ activities in Å$^4$/amu; no temperature or laser-frequency factor is included.

The complete field-by-field types, units, and availability conditions are generated from the runtime-owned documentation model in tako.vibration settings and result. The central distinctions are:

• frequencies_cm_inv is a signed real list for convenient classification; a negative entry represents the imaginary component with a minus sign;
• frequency_components preserves { real, imaginary } separately;
• energies_hartree also preserves complex components;
• modes[mode][atom][xyz] contains Cartesian displacements $\mathbf M^{-1/2}\mathbf e_k$ derived from normalized mass-weighted eigenvectors, in input atom order; fixed atoms are retained with zero displacement;
• zero_point_energy_* is the harmonic backend sum, not a temperature-dependent Gibbs correction;
• IR and Raman payloads are optional even when their names were requested, because calculator capability determines availability.

Define the magnitude threshold and scientific classification appropriate to the molecule, then preserve the full list for audit:

const imaginaryThresholdCmInv = 30;
const significantImaginary = result.vibrations.frequencies_cm_inv
  .map((frequency, modeIndex) => ({ modeIndex, frequency }))
  .filter(({ frequency }) => frequency < -imaginaryThresholdCmInv);

await tako.output('output/stationary-point-check.json', {
  threshold_cm_inv: imaginaryThresholdCmInv,
  significant_imaginary_modes: significantImaginary,
  classification:
    significantImaginary.length === 0 ? 'candidate minimum'
    : significantImaginary.length === 1 ? 'candidate first-order saddle'
    : 'higher-order instability or unconverged geometry',
});

That classification remains provisional until the corresponding eigenvectors are inspected. The numeric threshold is a reporting convention selected by the author, not a universal Tako acceptance criterion.

Displacement and stencil sensitivity

Quantitative use should establish that the relevant internal frequencies are stable with respect to the numerical derivative. Run a small sensitivity series without changing the structure or calculator:

const cases = [
  { displacement: 0.008, nfree: 2 },
  { displacement: 0.010, nfree: 2 },
  { displacement: 0.010, nfree: 4 },
];

const rows = [];
for (const entry of cases) {
  const value = await tako.vibration(structure, {
    calculator,
    vibrationDisplacementAngstrom: entry.displacement,
    vibrationNfree: entry.nfree,
    spectroscopy: ['frequencies'],
  });
  rows.push({
    displacement_angstrom: entry.displacement,
    nfree: entry.nfree,
    frequencies_cm_inv: value.vibrations.frequencies_cm_inv,
  });
}

await tako.output('output/displacement-sensitivity.json', rows);

Compare modes by displacement pattern as well as sorted frequency. Close or degenerate eigenvalues can exchange order between calculations, so array index alone is not always a stable mode identity.

Reproducibility record

A vibration report should preserve the input geometry, complete calculator construction, constraints, displacement, stencil, requested response channels, returned method identifiers, and raw result. It should also identify the optimization that supplied the stationary geometry. Frequencies from different geometries or calculators are not a numerical-convergence series.

The harmonic model has explicit boundaries: no anharmonic force constants, no thermal population model, no conformational ensemble, and no automatic standard-state thermochemistry are added by this call. Compute and report those layers separately if the scientific question requires them.