| ECCC-3 Contents Page | THEOCHEM Home Page | Elsevier Chemistry Home Page |
Jeremy R. Greenwood, Hugh R. Capper, Robin D. Allan, and Graham A.R. Johnston
We present a comprehensive theoretical study of the structure and tautomerism of the four isomeric hydroxy-pyridazine N-oxides, as well as pyridazine-1,2-dioxide. Gas phase properties are modelled with high-level ab initio calculations employing large basis sets (6-311++G(3df,3pd)) and quadratic configuration interaction treatment of electron correlation (QCISD(T)). Since these acidic heterocycles are of interest as novel carboxylate bio-isosteres, their anionic conjugate bases are also examined. Aqueous solution-phase properties are investigated using the isodensity polarised continuum model (IPCM) , and the semi-empirical AM1-SM2 and PM3-SM3 models, and relative acidities compared. The calculated properties are generally in good agreement with existing experimental data, indicating that the 1-hydroxy tautomer predominates both in the gas phase and in solution in the case of the 6- substituted system, and that hydroxy- 1-oxide tautomers predominate in the 3- and 5- substituted systems. The behaviour of the 4-substituted isomer is less clear, with non-planar 1-hydroxy and planar 4-hydroxy tautomers being similar in stability.
Table 1. hydroxy conformers
Table 2. all tautomers at HF/6-31G(d)
Table 3. compound energies and free energies at 298.15K
Table 4. higher level calculations on 4-hydroxy-pyridazine 1-oxide
Table 5. free energies of solvation
Table 6. tautomer ratios in solution
Table 7. relative acidities
Table 8. literature pKa values
Figure 1. pyridazine 1,2-dioxide
Figure 2. 3-hydroxy-pyridazine 1-oxide
Figure 3. 4-hydroxy-pyridazine 1-oxide
Figure 4. 5-hydroxy-pyridazine 1-oxide
Figure 5. 6-hydroxy-pyridazine 1-oxide
Chart 1. structures referred to in text
As a class, the pyridazines are the least studied and understood of the six-membered diazine heterocycles, both experimentally, and theoretically [1]. With the resurgence of interest in pyridazine derivatives as novel bio-active molecules [2,3,4] and with the close structural similarity of some hydroxy-pyridazines to better known systems such as uracil, we expect these compounds to be the subject of further synthetic and medicinal study. In particular, being relatively strong acids and having comparatively large dipoles, these heterocycles show promise as novel carboxylate bio-isosteres.
Pyridazine-1,2-dioxide 1, and 3-, 4-, and 5-hydroxy-pyridazine 1-oxides 2, 3, 4 , are known synthetically [5,6,7,8]. Various ring-substituted derivatives of 6-hydroxy-pyridazine 1-oxide 5 are also known [9], although the parent compound appears to have eluded synthesis, despite the fact that its immediate precursors 6-chloro- and 6-methoxy-pyridazine 1-oxide are known [10].
Experimental studies on the tautomerism of these compounds are few, and apparently the subject has been neglected since the 1960s. Itai [11] suggests that 2c is favoured over 2b, 3a over 3d, 4d over 4b, and 5a over 5c, i.e. that hydroxy N-oxides are preferred for 2 and 4, and N-hydroxy ketones for 3 and 5, citing [6,8,9,12,13] as evidence. Igeta [6,14] has made a thorough study of 2, and based on NMR and UV spectroscopic data has concluded the preference for the lactim form 2c, and found a pKa value of 4.1 for this system. 4 is also well studied by UV spectroscopy by Okusa and Kamiya [8] who concluded that the phenolic 4d is preferred. By contrast, 3 and 5 are poorly studied. Aspects of the chemistry of derivatives of 5, such as facile demethylation of 6-methoxy-pyridazine 1-oxides or methyl migration to form 1-methoxy 6(1H)pyridazinones [15], along with comparisons with other better-studied N-hydroxy alpha-oxo heterocyles (e.g. 9, [37]), may be taken to support the hypothesis of the predominance of 1-hydroxy tautomers in the absence of direct studies. The tautomerism of 3 and its derivatives is even less well studied; certainly alkylation and acylation yield 1-OR derivatives [7] suggesting 3a, but this is not conclusive, as kinetics or thermodynamics may render alkylation favourable or unfavourable regardless of tautomer prevalence. The conclusions of Itai [11] would seem valid in the case of the solution-phase behaviour of 2 and 4, and open to further investigation in the case of 5 and particularly 3.
The hydroxy-pyridazine 1-oxides are apparently previously unstudied theoretically. Pyridazine 1,2-dioxide 1 has received limited attention as the subject of semi-empirical MINDO/2, CNDO/2 and other calculations [16,17]. The closely related 3- and 4-hydroxy-pyridazines 12, 13 are better studied. Fabian [18,19] reports a semi-empirical (AM1, PM3, MNDO, MINDO/3) and ab initio (HF/3-21G) study of these compounds, in conjunction with 3,6-pyridazinedione 15 (maleic hydrazide) and 3,5-pyridazinedione 14, concluding inter alia that PM3 handles these systems better than other low-cost methods. La Manna et al. [20,21] report calculations to HF/6-31G(d)//3-21G, finding 3(2H)-pyridazinone and 4(1H)-pyridazinone to be favoured tautomers. Lapinski et al. [22,23] studied the 3-hydroxy-pyridazine 12 system at HF/3-21g and SCF-MBPT(2)/DZ(d,p) theories, in conjunction with IR experiments, and predicted 3(2H)-pyridazinone to be heavily favoured, in agreement with experiment. Maleic hydrazide 15, being perhaps the best known pyridazine derivative and a commercial plant-growth inhibitor, and having been determined as preferring the somewhat unusual lactam-lactim structure in numerous experiments [24], has been more thoroughly studied, and attempts made to include solvation effects. In addition to the work of Fabian [18,19], Hofmann et al. [25] studied maleic hydrazide tautomerism up to HF/6-31G(d), and calculated the free-energy of solvation of the favoured lactam-lactim tautomer to be -12.2 Kcal, via Tomasi's polarisable continuum model (PCM) [26]. A more extensive treatment was performed by Burton et al. [27], who conducted gas-phase calculations to MP4/6-31g(d,p)//6-31g(d,p), made comparisons with experiments based on the pKas of known methyl-derivative model compounds, and estimated solvent effects by three methods: molecular dynamics simulation (free-energy perturbation, FEP), self-consistent reaction field (SCRF) in an ellipsoidal cavity, and the PCM model.
The aim of the current study of hydroxy-pyridazine 1-oxides 2, 3, 4, 5 is to describe the complex tautomerism of these molecules by applying and building on the theoretical methods previously applied to similar systems. Additionally, the same methods are applied to their anionic conjugate bases, likely to be present at physiological pH. In order to lay theoretical ground for medicinal chemistry studies, aqueous solution-phase modelling is performed, including semi-empirical methods not previously applied to pyridazines, but described in detail for the tautomerism of 5-isoxazolols by Cramer and Truhlar [28]. For completeness and comparison, the non-acidic non-tautomeric pyridazine-1,2-dioxide 1 is also included in the study.
Following the description of heterocyclic tautomerism given by Elguero et al. [24], a set of all possible tautomers and their conjugate bases were generated. Firstly, non-substituted ring carbons were assumed to be at least singly protonated. This defines the four anions 2-, 3-, 4-, 5-. To each of these, the addition of a proton to the unsubstituted nitrogen, either of the two oxygens, or to a carbon alpha- to an oxygen substituent (forming a methylene), defines a complete set of structures; five neutral tautomers for 2, six for 3, five for 4, and four for 5. No tautomerism is possible, nor may an anion be formed of 1 under this system. Chart 1 provides a key to these structures.
Initially, all structures were built in planar symmetry. Rotation of 180° about the carbon-oxygen bond generates two distinct conformers in the case of hydroxy tautomers. Optimisations of both conformers of each hydroxy tautomer were performed at HF / 6-31G(d), and the results summarised in Table 1. The higher energy conformer of each hydroxy tautomer was discarded from the study set.
Optimisation and frequency calculations were performed at HF / 6-31G(d) on all tautomers, including the recommended scaling factor of 0.8929 for frequencies and thermochemistry calculated at 298.15K. Structures 3a, 3f, and 5c were found to have imaginary frequencies at this level, and were re-optimised in C1 symmetry. A high energy minimum for structure 5c was found in C1 symmetry at HF / 3-21G, but optimisation failed to find a minimum at HF / 6-31G(d). In addition, methylene-containing tautomers 2d, 2e, 3c and 4c were found to have particularly small lowest vibrational frequencies, e.g. 34 cm-1 for 2d, indicating very soft bending of the methylene out of the ring plane at this level of theory and basis set. Minima and energies are presented in Table 2.
Tautomers whose energies were greater than 25 Kcal above that of the lowest energy tautomer at HF / 6-31G(d), i.e. 3b, 3e, 3f, 5b and 5c, were then discarded from the set as being very unlikely to contribute significantly to tautomerism (cf [25]). A second form of 3a, the planar structure 3a', found to be a transition structure at this level, was added to the set for further study on account of the very small difference between the energies of the Cs and C1 structures. The remaining set of structures was optimised at MP2 / 6-31G(d,p), and this geometry was used in all subsequent ab initio calculations.
Ideally, high accuracy gas-phase energies would be obtained from a compound method such as G2 theory, or alternately a complete basis set method such as CBS-Q, but such calculations were too expensive to apply to each system. Instead, single point energies were calculated at MP2 theory with the large basis set 6-311++G(3df,3pd). The extra diffuse functions and high angular momentum functions are employed to treat as far as possible on equal footing anions and neutral species, and atoms of differing hybridisation, known to present problems for calculations of heterocyclic tautomerism [28]. Higher order correlation effects were considered by performing single point energy calculations at QCISD(T) / 6-31G(d,p). Compound energies were then calculated by correcting the MP2 / 6-31g(d,p) energy as follows, and are presented in Table 3.
E0 = ZPE{HF/6-31G(d)}*0.8929 + E{MP2/6-311++G(3df,3pd)} + E{QCISD(T)/6-31G(d,p)} - E{MP2/6-31G(d,p)}
The basis set for the base energy of this method is intermediate in level
between the G2(MP2,SVP) theory of Smith and Radom, and Curtiss et al. [29,30]
i.e. 6-31G(d), and G2(MP2) theory [31] i.e.
/6-311G(d,p). Both these methods include an empirical
EHLC
term in the calculation of E0, to correct for residual
basis set deficiencies. Since this term is a linear combination of the number
of alpha and beta electrons, which are identical for all of the isomeric
structures considered, this term cancels for relative energies, and has not
been included. Relative energies calculated by the above means are expected to
be more accurate than G2(MP2,SVP) theory, but probably not better than the more
expensive G2(MP2) theory, despite the larger basis set correction term
employed. In particular, poor or inconsistent treatment by MP2 at the base level
of calculation is a possible source of error.
From the scaled thermochemical properties calculated at HF / 6-31G(d), gas phase relative free energies of tautomerism and of proton loss at 298.15K were calculated, and these values are included in Table 3.
Since the compound E0 values and free energies calculated for 3a, 3a', and 3d were found to be very similar, a series of higher level calculations were performed on these three structures, summarised in Table 4. Modified and extended G2 and G2(MP2) energies are calculated, and relative free energies at 298.15K calculated from the compound energy derived from the highest level calculations available. The high cost of the MP4/6-311G(2df,p) calculation has prevented the direct calculation of G2 energies for the C1 symmetry structure 3a. For comparison, this energy is approximated as follows: the unknown effect of higher level (MP4 vs MP2) correlation on 3a at 6-311G(2df,p) is estimated as the known effect on the Cs symmetry conformer 3a', corrected by the known difference between these two effects calculated with the smaller basis set 6-311G(d,p)
Free energies of aqueous solvation were estimated by three different theories, and the results presented in Table 5. Firstly, the isodensity polarised continuum model (IPCM) [32], an extension to the PCM model used previously for solvation of pyridazinones [25,27], was used to calculate single point energies at MP2 / 6-31G(d,p), using a solvent dielectric constant of 78.4 and the recommended isodensity value of 0.001, from which were subtracted the gas-phase energies at this level. Secondly, both single point energies at the MP2 / 6-31G(d,p) geometries and optimisations were performed using the semi-empirical AM1, and AM1-SM2 aqueous solvation model of Truhlar and Cramer [33], and free energies of solvation determined as the difference between the two. Single point energies may be more reliable than optimised calculations, due to poor semi-empirical geometries, however both types of calculation are presented for comparison. For several optimisations, maximum gradients remained large, a known tendency of aromatic systems in this model [40]. Thirdly, the same calculations were performed using the PM3-SM3 model [34]. The PM3-based semi-empirical calculations are of interest in light of the known superior treatment of pyridazinediones by PM3 vs. AM1 [19,35]. The alternative AM1-SM1a model which requires explicit declaration of atom type, suggested by Cramer and Truhlar for the calculation of heterocyclic tautomerism [28], was not used because of the ambiguous hybridisation of certain atoms in the set, particularly with respect to anions.
Addition of the free-energies of solvation to the gas-phase free energies gives an estimate of the relative free energies of tautomers in aqueous solution at 298.15K. For two tautomers in equilibrium, if the free energies in solution are known to a sufficient degree of accuracy, the tautomeric ratio may be calculated:
HA1 <=> HA2 ; K2,1 = [HA2]/[HA1] = exp((G1-G2)/RT)
From the relative free energies, normalised tautomer ratios (neglecting tautomers of very low concentration) have been calculated by this method, and are presented in Table 6.
In principle, assuming one acid produces one base in solution, pKa or pKb can be calculated ab initio by modelling the following reactions in the gas phase, and correcting with free energy of solvation for each species
HA + H2O <=> H3O+
+ A- ; pKa = -log(exp(-
G/RT))
HA + OH- <=> A- + H2O
; pKb = -log(exp(-
G/RT))
The latter reaction has the attraction of equal numbers of charged species
on both sides, thus minimising systematic error. Unfortunately, the difficulty
of treating different anions equally in the gas phase, and in particular the
difficulty of modelling the hydroxide ion, and the difficulty of obtaining
sufficiently precise and consistent treatment of
G of aqueous solvation
across all species preclude the accurate ab initio calculation of pKa
by this approach at present.
By making the assumption that systematic errors are consistent across similar systems, we may write the following:
pKa = k - log(exp((GHA - GA-)/RT))
where k is assumed to be constant across a series of closely related acid-base pairs, and is the sum of systematic errors and free energy terms of neglected species. Given the tautomeric proportions ai calculated as above, if a single base gives rise to multiple acid tautomers, this may be extended:
pKa = k - a1*log(exp((GHA1 - GA-)/RT)) - a2*log(exp((GHA2 - GA-)/RT)) - ...
Estimations of relative pKas by this method are given in Table 7.
Some literature pKa values of closely related systems are given in Table 8 for comparison.
The anions and low energy tautomers optimised at MP2/6-31G(d,p) for each system are presented as spheres of the van der Waals radii. Gas phase free energies of tautomers relative to the anion are indicated in pink, calculated at 298.15K from compound energies (Table 3.) Electrostatic potential-derived charges from the Merz-Kollman-Singh scheme are overlayed on atoms in yellow. Dipoles are indicated by overlayed vectors, with units of 1 debye corresponding to distances of 1 angstrom.
A key to structures referred to in the text is given in Chart 1.
In order to evaluate these results, discussion of the major sources of error is in order. Smith and Radom [29], and Curtiss et al. [30] compare gas phase energies calculated by G2(MP2) and G2(MP2,SVP) with experiment for a large test suite of small acyclic species, including dissociation and ionisation energies, and electron and proton affinities. The relative gas phase energies calculated by the compound method described here should be intermediate in accuracy between G2(MP2) and G2(MP2,SVP). While we are actually considering proton affinities in this work, the test suite only concerned cations in this regard. Hence, electron affinities, incorporating inherently less reliable calculations on anions, are probably a better indication of accuracy. The results obtained were mean absolute deviations (MAD) from experiment of 0.6 Kcal and 0.7 Kcal for the two methods respectively for proton affinities, and 1.9 and 2.1 Kcal respectively for electron affinities. This gives an indication of the magnitudes of error in comparing the gas-phase acidities of different systems. The accuracy of the relative gas phase free energies of tautomerism, being between neutral tautomeric species within a single system, is likely to be handled somewhat better than this.
Energy calculations on these heterocycles are slow to converge with increasing basis set and correlation method; a glance at Table 4 shows variation in qualitative predictions to the highest levels trialed. Very high levels of theory may be necessary to evaluate relative energies to the degree required for accurate prediction of tautomer proportions; a difference of only 1.36 Kcal corresponds to a ten-fold difference in tautomeric ratio. None the less, other than for 3, we are confident that the required level of accuracy has been attained for qualitative prediction of gas phase tautomerism.
The highest level calculations performed have been unable to resolve the issue of which tautomer of 3 is preferred in the gas phase. Both the N-hydroxy 3a and 4-hydroxy 3d structures are low in energy. The qualitative predictions vary with basis set, correlation treatment, and compound method employed. The issue of whether the planar 3a' or non-planar 3a is preferred is also unresolved; larger basis sets seem to favour planar geometry, while better treatment of electron correlation favours non-planar. Optimisation with a larger basis set, i.e. at 6-311+G(d,p) seems to make a small difference, but between structures so close in energy, perhaps a significant one. The inference from these results is that any difference in energy between 3a and 3d is sufficiently small to indicate the likely presence of both structures in the gas phase. Entropy and thermal terms appear to favour 3a vs 3d slightly. Concerning the structure of 3a: if the planar 3a' actually represents a local maximum (saddle point), as indicated by the HF/6-31G(d) frequency calculation, then the energy barrier of inversion of conformation is exceedingly small. Conversely, if 3a' is lower in energy than 3a as suggested by compound theories, and is a truly a local minimum, then either the hydroxyl rotation out of the plane and ring puckering is an exceedingly soft vibrational mode, or the compound will display significant populations in both non-planar and planar minima.
By contrast, the energies of tautomers of 2, 4, and 5 are sufficiently separated to allow confident prediction of tautomerism. The 3-substituted system 2 shows a clear preference of 3.3 Kcal in the gas phase for the 3-hydroxy tautomer 2c over the NH oxo-oxo tautomer 2b, in accordance with the predictions of Itai [11], with all other tautomers being substantially less significant. Also predicted was the preference of 4 for the 5-hydroxy 1-oxide tautomer 4d over the NH oxo-oxo tautomer 4b, and this is indeed found to be the case by 6.5 Kcal at 0K. However, it is the N-hydroxy-oxo tautomer 4a which is found to be next lowest in energy, being 4.3 Kcal relative to 4d. The comparative stability of 4a was not suggested in Okusa's experimental study [8], and stands in contrast to the much less favoured analogous N-hydroxy-oxo tautomer 2a. 3.3 Kcal separates the N-hydroxy-oxo tautomer 5a from the 6-hydroxy 1-oxide 5d, in accordance with literature predictions, with no other significant contributing tautomers. Strong intramolecular hydrogen bonding is observed for this species, affecting the conformation, and reversing the normal preference for the HONN dihedral to be 0°. This bonding would also be likely to have implications for physical properties, e.g. solid state where hydrogen bonded dimers are recorded for similar systems, and acidity.
It is also noteworthy that while none of the methylene-containing tautomers appear to contribute significantly to tautomerism - Elguero [24] suggests these are more important in five-membered heterocyclic rings than six - none the less, some methylene tautomers are less unfavourable than expected. For example, 4(3H)-pyridazinone 1-oxide 3c is more stable than the NH isomer, 4(2H)-pyridazinone 1-oxide 3b.
The accurate prediction of solution-phase properties is a much more challenging task than the prediction of gas phase properties, particularly for an aqueous system such as this one, for which hydrogen bonding interactions in the first solvation sheath are likely to be crucial. Of the commonly used methods, molecular dynamics (FEP) treatments generally ignore the effects of solute polarisation, and to a large extent, solvent polarisation [28], while giving fair treatment of first solvation sheath effects. Conversely, continuum model (SCRF) treatments neglect first solvation sheath effects, but potentially consider polarisation effects accurately. The improvement of continuum models to include a cavity defined by an isodensity surface (IPCM) is recent in origin, and little data exists to indicate how well it handles these kinds of systems. Concern also exists about the suitability of IPCM for calculations on anions due to "charge leakage" which occurs with truncation of wavefunctions as the isodensity cavity is contracted. The semi-empirical SMx models hold the attraction of consideration of both polarisation terms, and through empirical atomic parametisation, first solvation sheath effects, and are therefore the methods of choice for this problem.
Across the different solvation methods trialed, the estimates of
G are generally in
fair agreement, with some anomalies (Table 5). As
expected, the IPCM method shows greater variability, while the semi-empirical
models are more concordant. Certain qualitative trends are reversed for IPCM,
e.g. the comparison of 3a and 3a', which alters the predicted
preferred structure in solution; IPCM predicts the non-planar structure, while
SMx more credibly predicts the stabilisation of the planar structure. The
solvation energies of the anions also appear to be under-estimated by IPCM
compared with semi-empirical methods, as is the solvation of 1, and
non-hydroxy tautomers in general. Despite the variablility of these figures,
the qualitative predictions of tautomer distribution (Table
6) remain substantially similar across all methods, with the exception of
3. The quoted [40] RMS errors for
neutral species and ions using AM1-SM2 and PM3-SM3, are 0.9 and 4.0, and 1.2 and
5.6 respectively. The indication is that qualitative comparisons from these
treatments between tautomers in the aqueous phase are probably valid, except
where tautomers lie very close in energy (3). Quantitative estimates,
taking into account both gas phase and solution phase errors, may easily be in
error by an order of magnitude or more. Comparisons of aqueous acidities between
systems, requiring consistent treatment of different anions vs. neutrals, are
unlikely to be accurate within 5 Kcal/mol at best. Given that a difference of
only 5 Kcal translates to a pKa difference of 3.7, the quantities in
Table 7 should be viewed with a great deal of
scepticism; current models are simply unable to make quantitative predictions by
this means. The literature only contains a single pKa value for a hydroxy
pyridazine N-oxide with which to judge these results - that of 2, which
has been determined to be 4.1 [6].
This would imply pKa values of well below zero for 3 and 4,
according to the semi-empirical solvation methods, which seems unlikely to be
correct. On the other hand, all methods trialed gave the same qualitative
predictions of the order of acidity, i.e. pKa 4 > 3 > 5
> 2. This prediction can be compared with experimental pKa values of
closely related hydroxy heterocycles (Table 8.)
Comparing 12 with 13, and 14 with 15 indicates
increased acidity of a hydroxyl in the 4- vs. the 3- position of pyridazine,
which is consistent with the prediction of 3 and 4 to be more
acidic than 2 and 5. N-oxidation of hydroxy-pyridines causes a
substantial drop in pKa, of between 2.3 and 5.7 units, and N-oxidation of 12
to 2 causes a drop of 5.8 units, thus the N-oxides of 13, i.e.
3 and 4, are expected to be strongly acidic, in accordance with
these qualitative predictions. The pKa of 10, having a hydroxyl meta-
to an N-oxide, is 0.5 units higher than that of the ortho- and para- isomers
9 and 11. This is in accordance with the prediction that 4
is more acidic than 3, and 5 is more acidic than 2.
Note also, that the predicted relative acidity of 5 increases from the
gas phase to solution, compared with other tautomers. This is to be expected if
the extra stability of strong intramolecular hydrogen bonding of the neutral
species in the gas-phase is diminished in aqueous solution, allowing easier
dissociation.
The literature experimentally-based predictions of the tautomerism of 3-, 5-, and 6- hydroxy-pyridazine 1-oxides 2, 4, 5 are confirmed by theory. The hydroxy N-oxide tautomers are predicted to be substantially favoured in the case of 2 and 4, both in solution and in the gas phase. The N-hydroxy-oxo tautomer of 5 is predicted to be substantially favoured in both phases, a result suggested by experiments on similar structures. No firm conclusion has been reached regarding 4-hydroxy-pyridazine 1-oxide 3, nor has sufficient experimental data been found with which to make comparisons. However, it seems likely that in solution both planar N-hydroxy and 4-hydroxy tautomers are present in significant quantities. Likewise, in the gas phase, both N-hydroxy and 4-hydroxy tautomers are likely to be significant, although the preferred structure of the N-hydroxy compound remains unresolved.
MP4 and QCISD(T) energies were calculated using Molpro 96 [39]. AM1-SM2 and PM3-SM3 were performed using AMSOL 5.4 [40]. All other calculations were performed using Gaussian 94 [41]. Platforms used were SGI R10000 Power Challenge and IBM RS6000.
We would like to thank the NSW Centre for Parallel Computing and the Australian National University Supercomputer Facility for support and donation of computing resources, and Malcolm Gillies for technical support.