REVIEW 4 major objections 5 minor 10 references
Delayed photoisomerisation of the trans-PSB3 retinal toy model using on-the-fly quantum dynamics
T0 review · 4 major / 5 minor · reviewed 2026-08-07 · deepseek-v4-flash
Pith's one-line read The paper claims that a full-dimensional quantum wavepacket method, DD-vMCG, makes the benchmark retinal model trans-PSB3 isomerise roughly 200-250 femtoseconds slower than trajectory methods, with the delay caused by a barrier before the…
desk verdict First DD-vMCG run on PSB3 shows a genuinely delayed isomerization, but the headline method-vs-trajectory comparison is confounded by different electronic structure levels; the abstract overclaims. read the letter →
The pith
A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.
The reading
What carries the argument
The load-bearing object is the DD-vMCG wavefunction ansatz, a variational multiconfigurational Gaussian wavepacket whose centers and momenta are propagated by the Dirac-Frenkel principle instead of by classical trajectories. On-the-fly electronic structure supplies energies, gradients, Hessians, and nonadiabatic couplings; propagation diabatisation converts the adiabatic surfaces into smooth diabatic ones to avoid geometric phase issues, and Shepard interpolation over a growing database of 48,857 electronic-structure points makes the propagation affordable. The local harmonic approximation around each Gaussian center supplies the potential integrals needed for a variational propagation without a precomputed global grid.
What would settle it
Repeat the benchmark surface-hopping and ab initio multiple spawning simulations using exactly the same interpolated two-state surface stored in the paper's database; if those trajectory methods also show a delay of roughly 200 fs, the delay is a property of the surface, not of quantum wavepacket dynamics.
Extended reading notes
Core claim
The central discovery is that DD-vMCG with SA(2)-CAS(6,6)/6-31G gives a qualitatively different deactivation timescale for PSB3 than independent-trajectory methods: the $S_1$ population starts transferring to $S_0$ only after about 200 fs and sits at a 60:40 mixture of $S_1$ and $S_0$ at 300 fs, even though the torsional angle has moved from about 180 degrees to 0 degrees and the BLA parameter has stabilised near 1.4 \AA. A database of 48,857 electronic-structure points shows the wavepacket lingering at the $S_1$-trans minimum, with only 0.4% of stored structures near the lowest-energy conical intersection, and the energy profile along $S_1$ shows a barrier before that intersection. The dynamics follow a sequential route: bond-length stretching in the first 25 fs, torsional activation only after about 150 fs, population transfer after 200 fs, and CH stretch modes dominating after 220 fs and driving the system to the cis minimum. Reduced-dimensionality runs that exclude CH stretches never reach the cis geometry, which the paper reads as evidence that those modes are required for successful isomerisation.
Load-bearing premise
The paper's central comparison assumes that the simpler two-state electronic-structure surface used here is close enough to the three-state-averaged and MRCI-corrected surfaces used in the benchmark trajectories that the observed delay is attributable to the quantum dynamics method, not to the missing second excited state or missing dynamic electron correlation.
Editorial extensions
If this is right
- PSB3 photoisomerisation timescales are method-dependent: trajectory-based values near 60 fs and wavepacket values near 250-300 fs are both possible outputs, so benchmark comparisons need to specify the nuclear-dynamics method.
- A barrier before the conical intersection becomes a concrete dynamical feature: the wavepacket lingers at the $S_1$-trans minimum, so transient trapping should be observable in time-resolved measurements if this picture holds.
- CH stretch modes gate the isomerisation: reduced-dimensional models that omit them will not reproduce the cis product, so effective-mode models of retinal chromophores should keep the high-frequency stretches.
- Two deactivation routes are available, one through a 120-degree conical intersection and one through an approx 80-degree near-degenerate region, and the wavepacket chooses between them depending on the shape of the surface.
- Full-dimensional on-the-fly quantum dynamics is now feasible for a 36-mode benchmark system, so similar wavepacket analyses can be applied to other photoactive molecules.
Reading between the lines
- Beyond the paper: because the earlier trajectory benchmarks used three-state averaging or added MRCI corrections, part of the reported delay could come from the simpler two-state CASSCF surface rather than from quantum wavepacket dynamics; a method-only comparison cannot separate these contributions.
- Beyond the paper: the variational coupling between Gaussian centres maintains quantum coherence across all 36 modes, so the lingering $S_1$ population may reflect coherent wavepacket motion that trajectory methods average away; the paper shows time-dependent mode widths but does not quantify coherence loss.
- Beyond the paper: an isotopic test follows immediately, deuterating the CH bonds should shift the timing or yield of isomerisation if those stretches are the gating modes, and the same dynamics could be rerun with deuterium masses.
Editorial analysis
A structured set of objections, weighed in public.
Referee Report
Summary. The manuscript reports on-the-fly direct dynamics variational multiconfigurational Gaussian (DD-vMCG) simulations of the trans-PSB3 retinal model in full dimensionality (36 normal modes), using a two-state-averaged CAS(6,6)SCF/6-31G electronic structure method. The central result is that the S1-to-S0 population transfer occurs only after about 200 fs and that the average torsional angle reaches cis-like values by 250-300 fs, which is hundreds of femtoseconds slower than previous AIMS and surface-hopping benchmarks. The authors attribute this delay to trapping in the S1 trans minimum and to a barrier en route to the conical intersection, and they use reduced-dimensionality runs to argue that CH stretch modes are essential for successful isomerisation.
Significance. If the central claim holds, this is a useful contribution to the method-comparison literature for nonadiabatic dynamics: it is, to my knowledge, the first full-dimensional on-the-fly vMCG study of PSB3, and the delayed timescale is an output of the propagation rather than a fitted parameter. The paper also ships a large electronic-structure database and Zenodo-hosted data, and it explicitly compares its population and geometry time series with well-known benchmarks, which is valuable for reproducibility. The main weakness is that the headline comparison is not method-only: the benchmark trajectory calculations used different electronic structure treatments, so the delayed isomerisation may be a property of the SA(2)-CASSCF surface rather than of quantum wavepacket dynamics. This confound must be resolved or carefully qualified before the paper's central attribution can be accepted.
major comments (4)
- [Abstract; Section 2.1] The central claim that DD-vMCG shows isomerisation 'hundreds of femtoseconds slower using the same electronic structure method' is not supported by the text. Section 2.1 states that the Martínez AIMS benchmark used a three-state average and that Szymczak et al. used SA(2)-CASSCF followed by MRCI, while the present work uses plain SA(2)-CAS(6,6)/6-31G without post-CASSCF correlation. The delayed timescale could therefore reflect the different potential energy surfaces (missing S2 state, missing dynamic correlation) rather than the nuclear dynamics method. Please either qualify the claim, run a trajectory method on the identical SA(2)-CAS(6,6)/6-31G surface, or provide quantitative evidence that the relevant surface features (S1 minimum, barrier height, CI accessibility) are comparable across these electronic-structure choices. The barrier in Figure 5 is computed on the present PES and cannot independently resolve this confound.
- [Section 3.1; Conclusions] The conclusion that PSB3 'has fully isomerised' by 250-300 fs is in tension with Figure 3(a), which shows that about 60% of the population remains in S1 at 300 fs. The average torsional angle and BLA in Figures 3(c) and 3(d) are weighted averages over Gaussian centers and can mask a bimodal wavepacket (for example, trans-like population on S1 and cis-like population on S0). The reported standard deviations describe the spread of Gaussian centers, not a state-resolved nuclear density. Please report state-resolved torsional and BLA distributions, or at least diabatic-population-weighted histograms, and rephrase the 'fully isomerised' claim to reflect that the average structure becomes cis-like while substantial S1 population remains.
- [Section 3.2; Figure 4] The database histogram is used to infer that the wavepacket is 'trapped' in the S1_trans minimum and that two deactivation routes exist near MECI_db and MECI_opt. However, the density of stored electronic-structure points is not a probability distribution of the nuclear wavepacket; it reflects on-the-fly sampling, database interpolation queries, and Hessian-update decisions. The observation that 3.6% of the 48,857 structures lie near S1_trans does not by itself establish trapping unless it is shown that these points carry significant Gaussian weight in the vMCG expansion. Please support the mechanistic conclusions with GBF-population-weighted geometrical distributions, and similarly qualify the two-route interpretation based on Figures S13 and S14.
- [Section 3.3; Figure 6] The claim that CH stretch modes are crucial for successful isomerisation is based on the observation that the 27-mode calculation (which excludes CH stretches) does not relax to the cis structure. This comparison is not fully controlled: the reduced-dimensionality runs use independent databases, freeze several modes that are active in the full run (for example q16, q23, and q27 appear in the active-mode analysis of Section 3.1), and removing modes also removes their zero-point energy and changes the local harmonic approximation. The 13- and 19-mode calculations exclude many other modes as well, so they do not isolate CH stretches. A more direct test would be to freeze only the CH-stretch modes while keeping all other modes active, ideally using the same database, before drawing the mechanistic conclusion.
minor comments (5)
- [Section 2.2] The text refers to 'equation 3' when introducing the vMCG ansatz, but the ansatz is equation (1); also, the local harmonic expansion in equation (3) omits the remainder term and does not define the notation for the Hessian.
- [Figure 5] The x-axis of Figure 5 is labeled in arbitrary units and the 'vector that connects the S0_trans and MECI_db geometries' is not defined; please describe how the one-dimensional cut was constructed and report the barrier height in a well-defined coordinate.
- [References] Reference 39 is incomplete: it gives only 'J. Chem. Phys., 2016, 144, 024111' without authors or title, which makes it difficult for readers to identify the coupled-coherent-states paper being cited.
- [Section 3.3] The text mentions '19 and 24mode' calculations, but the preceding sentence defines only 13, 19, and 27 modes; please correct this typo and verify the number of modes in each reduced-dimensionality run.
- [Data availability] The data-availability statement lists the 20GWP and 16GWP databases but does not specify the versions of Molpro and Quantics used or the convergence thresholds for the SA(2)-CASSCF energies, gradients, and nonadiabatic couplings; adding these details would improve reproducibility.
Circularity Check
No significant circularity: the delayed isomerisation is an output of an ab initio on-the-fly propagation, not a fitted or self-defined result; the main caveat is a method-comparison confound, which is a correctness issue, not circularity.
full rationale
The paper's central claim—that DD-vMCG dynamics isomerises trans-PSB3 hundreds of femtoseconds slower than trajectory-based methods—is presented as a direct output of the on-the-fly quantum dynamics simulation. No parameter is fitted to reproduce the delay; the timescale emerges from propagating a nuclear wavepacket on SA(2)-CAS(6,6)SCF/6-31G surfaces, with energies, gradients, couplings, and Hessians evaluated on the fly. The S1-to-S0 population transfer and the 250–300 fs cis formation are read off from the propagation (diabatic populations, torsional-angle expectation values, and BLA), not imposed by construction. The explanation that an energy barrier on the S1 surface delays access to the conical intersection is derived from the same simulation database (Figure 5), but this is a post-hoc mechanistic rationalisation of an observed dynamical result, not an input that predetermined it. The reduced-dimensionality runs that identify CH stretches as crucial are likewise new simulations whose outcomes (no isomerisation when those modes are frozen) are not built into the method. The only self-citations are to previous MCTDH/vMCG studies of other molecules, used as contextual examples that quantum wavepacket methods can differ from trajectory methods; these citations are not load-bearing for the PSB3 result itself. The skeptical concern that the comparison to AIMS and surface-hopping benchmarks is not strictly method-only (the benchmarks used SA(3)-CASSCF or MRCI corrections, while the present work uses plain SA(2)-CASSCF) is a legitimate threat to the attribution of the delay to the nuclear-dynamics method, but it is a question of accuracy and comparability, not circularity: the delay is not defined as, nor fitted to, the benchmark timescales. No self-definitional step, renamed fit, or uniqueness theorem invoked from the authors' prior work was found.
Assumptions & free parameters
free parameters (1)
- Width parameters of the frozen Gaussian basis functions =
not reported in main text
assumptions (6)
- domain assumption Born-Oppenheimer separation and restriction to two electronic states (S0 and S1) is sufficient for the dynamics.
- domain assumption The Local Harmonic Approximation with a second-order Taylor expansion of the potential around each Gaussian center is accurate for the sampled geometries.
- domain assumption The propagation diabatisation and Shepard interpolation between stored database points give a smooth and accurate diabatic representation near conical intersections.
- domain assumption Hessians from the initial structure, updated by gradient information, remain accurate throughout the 300 fs propagation.
- domain assumption The initial nuclear wavefunction is the harmonic ground state of S0 projected vertically onto S1.
- ad hoc to paper The literature trajectory results (AIMS, TSH) can be directly compared with the DD-vMCG simulation despite differences in electronic structure treatment (SA(3) or MRCI).
Cite this review
Pith. "Pith review of Delayed photoisomerisation of the trans-PSB3 retinal toy model using on-the-fly quantum dynamics." pith.science (2026). https://pith.science/paper/O4AQ4UE4
@misc{pith2026250603811,
author = {Pith},
title = {Pith review of: Delayed photoisomerisation of the trans-PSB3 retinal toy model using on-the-fly quantum dynamics},
year = {2026},
howpublished = {\url{https://pith.science/paper/O4AQ4UE4}},
note = {Machine review of arXiv:2506.03811}
}
read the original abstract
We explore the trans-cis photoisomerisation process in a representative retinal protonated Schiff base known as trans-PSB3, employing the quantum dynamics method direct dynamics variational multiconfigurational gaussian -- DD-vMCG -- in full dimensionality, i.e., 36 degrees of freedom on potential energy surfaces computed on-the-fly using the SA(2)-CAS(6,6)SCF electronic structure method with the 6-31G basis set. Although the toy molecule has been extensively studied using trajectory methods such as Tully Surface Hopping and Ab Initio Multiple Spawning, the on-the-fly quantum dynamics method DD-vMCG shows a trans-cis isomerisation hundreds of femtoseconds slower using the same electronic structure method, which can be explained in terms of the accesibility to the conical intersection connecting the ground and the excited state.
Reference graph
Works this paper leans on
-
[1]
is one of the simplest minimal models for 11-cis-retinylidene, 31 capturing key aspects of the isomerisation mechanism and enabling qualitative a Departamento de Química Física, Universidad de Salamanca, https://ror.org/02f40zc51, 37008, Spain; b Departamento de Química, Módulo 13, Universidad Autónoma de Madrid, https://ror.org/01cby8j38, 28049 Madrid, S...
-
[2]
Notably, around the top right yellow cross representing the minimum of the excited state for the trans-PSB3 geometry, we observe a high density of structures with increasing BLA values, consistent with the dynamical trends shown in Figure 3(d). The MECI_db identified in our simulation lies near the center of the most frequently visited regions of the datab...
work page 1939
-
[6]
14 O. Kobayashi and S. Nanbu,Chem. Phys., 2015,461, 47–57. 15 A. Valentini, D. Rivero, F. Zapata, C. García-Iriepa, M. Marazzi, R. Palmeiro, I. Fdez. Galván, D. Sampedro, M. Olivucci and L. M. Frutos,Angew. Chem., Int. Ed., 2017,56, 3842–3846. 16 J. K. Yu, R. Liang, F. Liu and T. Martínez,J. Am. Chem. Soc.,
work page 2015
-
[7]
Programa Propio C1
In conclusion, in this work we tested the DD-vMCG method on the trans-PSB3 system, which has been widely studied using mixed quantum–classical trajectory methods. The tested method is able to predict the photoisomerisation of the retinal toy model, although discrepancies in timescales and pathway exploration between quantum and mixed quantum–classical met...
-
[120]
23 D. Grabarek, E. Walczak and T. Andruniów,J. Chem. Theory Comput., 2016,12 5, 2346–56. 1–11 | 10 24 D. Tuna, D. Lefrancois,Łukasz Wolański, S. Gozem, I. Schapiro, T. Andruniów, A. Dreuw and M. Olivucci,J. Chem. Theory Comput., 2015,11 12, 5758–81. 25 L. Liu, J. Liu and T. J. Martinez,J. Phys. Chem. B, 2016,120, 1940–1949. 26 M. Filatov, S. K. Min and K....
work page 2016
-
[2012]
36 G. A. Worth and B. Lasorne,Quantum Chemistry and Dynamics of Excited States: Methods and Applications, 2020, 413–433. 37 G. A. Worth, M. A. Robb and I. Burghardt,Faraday Discuss., 2004,127, 307–323. 38 T. J. Martinez, M. Ben-Nun and R. D. Levine,J. Phys. Chem., 1996,100, 7884–7895. 39J. Chem. Phys., 2016,144, 024111. 40 A. Raab,Chem. Phys. Lett., 2000,...
work page 2020
-
[2015]
48 A. Viel, R. P. Krawczyk, U. Manthe and W. Domcke,J. Chem. Phys., 2004,120, 11000–11010. 49 S. Gómez, E. Spinlove and G. Worth,Phys. Chem. Chem. Phys., 2024,26, 1829–1844. 50 S. Gómez, L. M. Ibele and L. González,Phys. Chem. Chem. Phys., 2019,21, 4871–4878. 51 S. Gómez, N. K. Singer, L. González and G. A. Worth,Can. J. Chem., 2022,101, 745–757. 1–11 | 11
work page 2004
-
[2019]
17 L. E. Cigrang, B. F. E. Curchod, R. A. Ingle, A. Kelly, J. R. Mannouch, D. Accomasso, A. Alijah, M. Barbatti, W. Chebbi, N. Došlić, E. C. Eklund, S. Fernandez-Alberti, A. Freibert, L. González, G. Granucci, F. J. Hernández, J. Hernández- Rodríguez, A. Jain, J. Janoš, I. Kassal, A. Kirrander, Z. Lan, H. R. Larsson, D. Lauvergnat, B. L. Dé, Y. Lee, N. T....
work page Pith review arXiv 2025
Show all 10 references
-
[2020]
Vreven, F
11 T. Vreven, F. Bernardi, M. Garavelli, M. Olivucci, M. Robb and H. Schlegel,J. Am. Chem. Soc., 1997,119, 12687–12688. 12 Y. Liu and C. Zhu,Phys. Chem. Chem. Phys.,
1997
-
[2022]
Liang, F
7 R. Liang, F. Liu and T. J. Martínez,J. Phys. Chem. Lett., 2019,10, 2862–2868. 8 B. G. Levine and T. J. Martínez,Annu. Rev. Phys. Chem., 2007,58, 613–634. 9 M. Ruckenbauer, M. Barbatti, T. Muüller and H. Lischka,J. Phys. Chem. A, 2013,117, 2790–2799. 10 S. Gozem, P. J. M. Joh...
2019
Reviewed August 7, 2026 · model on record in the stance chip above.
Discussion (0). Continue with ORCID to comment.