Estudio de la hernandulcina: conexión estructura-propiedad

Facultad de Ciencias Trabajo Fin de Grado Grado en Química Estudio de la hernandulcina: Conexión estructura-propiedad Autor: Alba Arribas Sanz Tutores: Iker León y José Luis Alonso 2 ÍNDICE 1. RESUMEN/ ABSTRACT ......................................................................................................................... 1 2. INTRODUCCIÓN .................................................................................................................................... 3 2.1. Importancia del dulzor .................................................................................................................. 3 2.2. Edulcorantes en la actualidad ....................................................................................................... 4 2.2.1. Uso de edulcorantes en la dieta ............................................................................................ 4 2.2.2. Uso de edulcorantes en fármacos ......................................................................................... 4 2.3. Percepción del dulzor ................................................................................................................... 5 2.3.1. Percepción general del gusto ................................................................................................ 5 2.3.2. Mecanismos de percepción del dulzor ................................................................................. 6 2.3.3. Relación estructura-propiedad del dulzor ............................................................................ 7 3. JUSTIFICACIÓN ................................................................................................................................... 10 4. OBJETIVOS .......................................................................................................................................... 12 5. METODOLOGÍA .................................................................................................................................. 13 5.1. Mecánica molecular .................................................................................................................... 13 5.1.1. Superficies de energía potencial ......................................................................................... 13 5.1.2. Campos de fuerzas .............................................................................................................. 14 5.2. Métodos mecanocuánticos ......................................................................................................... 17 5.2.1. Métodos basados en la función de onda ............................................................................ 18 5.2.1.1. Métodos perturbativos moller-plesset ........................................................................ 18 5.2.2. Métodos semiempíricos ...................................................................................................... 19 5.2.2.1. Parametric model 6 ..................................................................................................... 19 5.2.3. Teoría del funcional de la densidad (DFT). .......................................................................... 20 5.2.3.1. Funcional B3LYP ........................................................................................................... 21 5.2.3.2. Funcional M062X ......................................................................................................... 21 5.2.3.3. Dispersiones de grimme y de becke johnson .............................................................. 21 5.2.4. Bases de cálculo .................................................................................................................. 23 5.3. Espectroscopía de rotación ......................................................................................................... 25 5.3.1. Fundamentos....................................................................................................................... 25 5.3.2. Expansión supersónica y ablación láser .............................................................................. 27 5.3.3. Rotación interna .................................................................................................................. 28 4 6. RESULTADOS ...................................................................................................................................... 29 6.1. Búsqueda conformacional .......................................................................................................... 29 6.1.1. Comparativa de bases y niveles de cálculo ......................................................................... 33 6.2. Predicción del espectro ............................................................................................................... 40 6.2.1. Evaluación de la existencia de rotación interna .................................................................. 45 6.3. Correspondencia con el triángulo del dulzor .............................................................................. 47 7. CONCLUSIONES .................................................................................................................................. 50 8. BIBLIOGRAFÍA ..................................................................................................................................... 51 9. ÍNDICE DE FIGURAS, TABLAS Y ECUACIONES. ................................................................................... 54 1 1. RESUMEN En este trabajo de fin de grado se estudia el panorama conformacional de la hernandulcina, una molécula de origen natural que posee un gran potencial como edulcorante. Para buscar las estructuras más relevantes, se hace uso de métodos computacionales. En primer lugar, se buscan los posibles confórmeros con mecánica molecular. Posteriormente, las estructuras se optimizan con métodos del Teorema del Funcional de la Densidad. Esto nos permite realizar una predicción de los espectros de rotación y evaluar la existencia de rotación interna en la molécula, que será de gran ayuda de cara a una posterior experimentación. Finalmente, se comprobará si las estructuras más relevantes de la molécula cumplen con los requisitos de la teoría del dulzor de Shallenberger-AcreeKier. Para ello, se evaluará la correspondencia entre los parámetros geométricos de los confórmeros más estables con la teoría del dulzor. Por otro lado, se comparan las 30 estructuras más estables empleando distintas funciones de base y niveles de cálculo, incluyendo métodos semiempíricos, DFT y post-Hartree-Fock, cuyos resultados se procesan para dilucidar que estructuras proporcionan resultados más fiables y el coste computacional que implica cada metodología. Palabras clave: (Dulzor, Química Computacional, Espectroscopia de Rotación, Confórmeros) ABSTRACT In this work the conformational landscape of the hernandulcin, a molecule of natural source with great potential as sweetener, is studied. In order to find the most relevant structures, computational methods are used. In first place, molecular mechanics is employed to find the possible conformers. Afterwards, the structures are optimized using Density Functional Theory methods. This allows us to predict the rotational spectra and to evaluate the existence of internal rotation, which will be of great help in a possible future experimentation. Finally, whether the most relevant structures of the molecule fulfil the requirements of the sweet taste theory proposed by Shallenberger-Acree-Kier will be evaluated. In order to this, the correspondence between the geometrical parameters of the most stable conformers will be evaluated according to the sweet taste theory. Additionally, the 30 most stable structures are compared employing different basis sets and levels of theory, including semiempirical, DFT, and post-Hartree-Fock, and the results are processed to elucidate their reliability and their computational cost. Keywords: (Sweetness, Computational Chemistry, Rotational Spectroscopy, Conformers) 2 3 2. INTRODUCCIÓN El dulzor es uno de los cinco sabores básicos que universalmente se asocia a una sensación agradable. Lo que percibimos como dulce es generalmente causado por el sabor de los azúcares naturales, tales como los que podemos encontrar en frutas y miel. Existen otros tipos de sustancias que activan en el gusto la respuesta de las células sensoriales responsables del dulzor. Algunos ejemplos de ello son ciertos amino ácidos, alcoholes, aldehídos o cetonas. 2.1. IMPORTANCIA DEL DULZOR La razón de la preferencia por el dulzor tiene un sentido evolutivo, tal y como se describe en numerosas investigaciones como las de Hladick et al. 1 Estos estudios sugieren que la capacidad de identificar moléculas dulces permite a los seres vivos detectar fuentes de azúcares simples. El agrado por el dulzor está relacionado con la presencia de vitaminas y minerales, la madurez de las frutas y la presencia de azúcares fácilmente metabolizables, que proporcionan una fuente de energía rápidamente accesible para el organismo. Por otro lado, se cree que la utilidad evolutiva reside no sólo en el dulzor, sino en una escala dulce-amargo (sabores considerados opuestos), de forma que un sabor dulce se relaciona con una fuente energética, y un sabor amargo se corresponde con una señal de toxicidad. Gracias a esta escala, un animal puede considerar si una planta es beneficiosa y por tanto comestible o si debe ser rechazada por la presencia de potenciales toxinas2–4. La percepción del dulzor se relaciona con el consumo de vegetales debido a que, tras realizar estudios en animales, se ha llegado a la conclusión de que la preferencia por los sabores dulces tiene lugar en especies cuya dieta está basada en plantas, mientras que los estrictamente carnívoros no poseen sensibilidad al dulzor. Estos animales bien han perdido alguna de las proteínas que actúan como receptor, o bien los genes que las codifican han pseudogenizado. Se ha comprobado que hay ejemplos de animales distintos con antepasados comunes en los que unos han perdido los mecanismos de recepción del dulzor al prescindir de las plantas en su dieta, mientras que no lo hacen los que siguen consumiendo plantas. Por ejemplo, la nutria asiática (Aonyx cinereus), estrictamente carnívora, no es capaz de detectar compuestos que los humanos percibimos como dulces, mientras que el oso andino (Tremarctos ornatus), cuya dieta es omnívora, muestra una fuerte preferencia por los azúcares y algunos edulcorantes.2,5 Otro indicio de la importancia del dulzor desde un punto de vista evolutivo puede observarse en el hecho de que, en los recién nacidos, el gusto es el sentido más desarrollado junto con el dolor. Concretamente el dulzor es uno de los primeros estímulos que se pueden percibir positivamente de forma innata, produciendo respuesta incluso a disoluciones muy diluidas de azúcar. Los estudios de comparación genética pueden aportar información que ayude a entender cómo nuestra historia 4 evolutiva ha modificado nuestro sentido del gusto y dieta. Este conocimiento es necesario para comprender por qué en la actualidad poseemos una predilección por este tipo de sabor, que nos ha llevado a utilizar sustancias dulces en nuestra alimentación, muchas veces añadidas artificialmente. 6,7 2.2. EDULCORANTES EN LA ACTUALIDAD Debido a la preferencia que poseemos por el sabor dulce, se han buscado y desarrollado sustancias comestibles que tengan este sabor y no sean necesariamente azúcares. Un ejemplo claro son los aditivos alimentarios, que son sustancias que normalmente no se consumen como alimentos en sí mismos, sino que se añaden intencionalmente a los alimentos con un fin determinado como la conservación del alimento o la modificación de sus propiedades. 8 De forma más específica, los edulcorantes son sustancias naturales o artificiales con la propiedad de producir un sabor dulce similar al del azúcar. Los edulcorantes artificiales se definen en el Código Alimentario Español como sustancias sápidas sintéticas que, sin tener cualidades nutritivas, poseen un poder edulcorante superior al de cualquier hidrato de carbono al que sustituyen o refuerzan. 9 A continuación, se detallan dos de los usos más importantes que se dan a los edulcorantes hoy en día. 2.2.1. USO DE EDULCORANTES EN LA DIETA Debido a la preferencia por los alimentos dulces, existe una tendencia por recurrir a alimentos que posean esta cualidad. El efecto en los seres humanos ha sido tal que la abundancia de azúcar en la dieta moderna ha llegado a niveles nocivos en algunos casos, relacionándose directamente con el aumento de enfermedades como obesidad, diabetes, hipertensión y caries. Por ello, se ha buscado sustituirlo con edulcorantes no calóricos que ayuden a disminuir su incidencia, siendo la principal razón del uso y desarrollo de edulcorantes en comida, bebida y otros productos como pasta de dientes. 10 2.2.2. USO DE EDULCORANTES EN FÁRMACOS El sabor es el factor que más influye en cuanto a la aceptación de los fármacos administrados por vía oral, y es un requisito importante cuando se trata de pastillas o comprimidos que poseen de por sí un sabor amargo, especialmente cuando se destinan a pacientes pediátricos o geriátricos. Es también un factor importante a considerar desde el punto de vista comercial, razón por la que se utilizan edulcorantes, ya sean artificiales o naturales, en la composición de los medicamentos. 11 5 Así mismo, se debe tener en cuenta que las interacciones entre el principio activo y los excipientes pueden afectar a la naturaleza química, estabilidad y la biodisponibilidad del medicamento, que darían lugar a cambios en la eficiencia y seguridad del fármaco. 12–14 Algunos de los edulcorantes que podemos encontrar en este campo son: • La sacarosa es uno de los más utilizados, con ella se hacen siropes que aportan viscosidad y consistencia a los fluidos. • La lactosa es otro azúcar ampliamente utilizado en farmacia como excipiente, especialmente en pastillas. • La trehalosa se emplea en la conservación de la sangre destinada a las transfusiones ya que ayuda a incrementar el tiempo de vida de las plaquetas. También se utiliza para preservar los embriones sometidos a liofilización aumentando la duración de su viabilidad. • El xilitol se utiliza en recubrimientos de vitaminas o como expectorantes por el efecto de vaporización que produce al tomarlo. Así mismo, puede encontrarse en la composición de suplementos alimentarios como aminoácidos, elementos traza o azúcares no reductores. 2.3. PERCEPCIÓN DEL DULZOR En los apartados anteriores se ha manifestado la importancia de los edulcorantes, por ello, es fundamental conocer cómo funcionan para saber diseñar o encontrar nuevas sustancias dulces. Es decir, deben conocerse los mecanismos por los cuales percibimos los sabores. 2.3.1. PERCEPCIÓN GENERAL DEL GUSTO El gusto, junto con el olfato, es uno de los sentidos que se percibe por medio de interacciones químicas, por lo que el sabor de una sustancia depende de su naturaleza y actividad química. Existe un amplio espectro de sabores que implican sabores secundarios como el metálico o el de la grasa, así como sensaciones que implican la textura o la temperatura, pero para simplificar, nos centraremos en los cinco sabores básicos que podemos percibir los humanos: ácido, amargo, dulce, salado y umami. 15,16 • El sabor ácido es desencadenado por un cambio en la concentración de protones en la membrana celular (Figura 1B). • Para el sabor salado, el estímulo principal viene dado por la permeación de iones Na+ a través de los canales en las papilas gustativas, que producen una despolarización de las células receptoras (Figura 1C). 12 4. OBJETIVOS El objetivo principal de este trabajo es realizar un estudio computacional de la molécula de hernandulcina con el propósito de dilucidar qué disposiciones adopta en un entorno libre de interacciones, y determinar cuáles serán sus confórmeros más estables. Se pretende dar una información inicial de la molécula que sirva de base para posibles estudios posteriores, como por ejemplo en medios acuosos o formando complejos con unidades clave en las proteínas que actúan de receptor. Obtener los confórmeros más estables, sus momentos dipolares y constantes de rotación con el fin de predecir el espectro rotacional de la molécula, que facilite la posterior interpretación del espectro obtenido experimentalmente. También se preverá en qué región de frecuencias sería más conveniente realizar el experimento. Comprobar si la hernandulcina es una sustancia cuyo dulzor se explique por la teoría de Shallenberger-Acree-Kier, lo que fundamenta el diseño de nuevos edulcorantes partiendo de estructuras moleculares semejantes. Realizar una comparación de los resultados proporcionados por distintos niveles de cálculo, aprovechando que la molécula a estudiar es relativamente grande y existe poca información para estos casos, centrando más la investigación en los métodos DFT y distintas bases de Pople. Esto permitirá dilucidar cuales de los niveles de cálculo disponibles son los más adecuados para el estudio de moléculas con un gran número de átomos e interacciones intramoleculares. 13 5. METODOLOGÍA Los estudios estructurales de este tipo se suelen realizar siguiendo dos líneas de trabajo. En primer lugar, se efectúa un análisis teórico de la molécula con el fin de obtener una visión del panorama conformacional, así como de las energías y parámetros espectroscópicos relevantes. Después, con la información obtenida, se comienza el trabajo experimental. Por último, se concluiría comparando los resultados de ambos procedimientos para dilucidar cuáles son los confórmeros que presenta la molécula realmente. Para este trabajo de fin de grado, nos centraremos en la primera parte teórica, cuyos resultados serán relevantes de cara a una posible experimentación posterior. 5.1. MECÁNICA MOLECULAR Con el fin de buscar los distintos confórmeros estructurales de una molécula sin hacer uso de métodos computacionales, sería necesario hacer uso de la intuición química a fin de encontrar las estructuras más relevantes. Esto conllevaría una gran cantidad de tiempo además de correr el riesgo de no considerar todos los confórmeros. Una metodología ampliamente utilizada para la búsqueda conformacional son los métodos basados en campos de fuerza. Estos emplean un tratamiento puramente clásico, por lo que ofrecen una aproximación muy simplificada en cuanto al cálculo de energías y parámetros estructurales, pero resultan muy efectivos para estudiar los posibles confórmeros de sistemas con una elevada cantidad de estructuras o sistemas que posean un gran número de átomos en un tiempo reducido. Para ello, se buscan los mínimos de su curva de energía potencial mediante métodos de mecánica molecular. 5.1.1. SUPERFICIES DE ENERGÍA POTENCIAL Como bien sabemos, una molécula puede adoptar infinitas disposiciones en el espacio cambiando las medidas de sus parámetros geométricos, y cada una de esas disposiciones tendrá una determinada energía. La representación de esa energía respecto de los parámetros geométricos se denomina superficie de energía potencial (SEP). En el caso de una molécula diatómica, donde se tiene un solo grado de libertad (la longitud del enlace), la representación de la energía viene dada por la curva de energía potencial en dos dimensiones. Éste sería el caso más sencillo, sin embargo, cuando se tratan moléculas poliatómicas la representación requiere de múltiples dimensiones, dando lugar a una hipersuperficie de energía potencial. En esta superficie se hallan mínimos correspondientes a las estructuras estables, denominadas confórmeros, y todos los caminos de interconversión entre ellos incluyendo estados de transición, suele representarse como se observa en la Figura 8. Una 14 estructura de equilibrio se encuentra en un punto de la SEP para la que son cero todas las primeras derivadas de la energía con respecto a las coordenadas geométricas individuales y cuya representación diagonal de la matriz de derivadas segundas de la energía, tiene todos los elementos positivos. En términos simples, una estructura de equilibrio corresponde al fondo de un pozo en la superficie de energía potencial global.33 5.1.2. CAMPOS DE FUERZAS Los métodos de mecánica molecular (MM, molecular mechanics) calculan la energía mecánica involucrada en la deformación de las moléculas. Los átomos en este modelo se consideran masas esféricas y los enlaces que los unen como tensores. Para describir las interacciones entre las masas se utilizan funciones de potencial derivadas de la mecánica clásica. La expresión de la energía potencial de cada geometría es la suma de la tensión de enlaces, torsión alrededor de los enlaces simples y ángulos de flexión, y las fuerzas intermoleculares como los enlaces de hidrógeno o de van der Waals. Es obvio que este tipo de método no proporcionará valores de la energía precisos, pero es perfectamente válido para una primera búsqueda conformacional, especialmente para moléculas de gran tamaño como en el caso de este trabajo. Además, el coste computacional que supone es menor que el de otros métodos. De esta manera es posible evaluar la energía de los puntos de la superficie de energía potencial en un tiempo reducido. 34 Para comenzar, la información necesaria que se requiere es la estructura de la molécula, o de forma más específica, las posiciones relativas de los átomos que la componen en el espacio, así como una Figura 8: Representación habitual de una superficie de energía potencial. 15 serie de datos tales como las constantes de fuerza k de los distintos enlaces, ángulos de giro y torsiones que no pueden producirse, como por ejemplo en los dobles enlaces. Las constantes se encuentran en la base de datos del programa utilizado y han sido obtenidas experimentalmente o mediante cálculos ab initio. Esta base de datos de compuestos empleada durante la parametrización es fundamental para una búsqueda conformacional correcta, así como para modelar otras propiedades. El conjunto de parámetros obtenidos empíricamente para describir las interacciones entre átomos y las funciones de potencial utilizadas forman el campo de fuerzas (FF, force field) en el que se basan este tipo de métodos. Existen varios campos de fuerzas que pueden utilizarse, cada uno con sus propias características. Por ejemplo, la familia MM, son ampliamente utilizados para cálculos computacionales de moléculas pequeñas, existiendo distintas variaciones dentro de este tipo que se designan con números (por ejemplo, MM2, MM3…). Los denominados AMBER (Assisted Model Building and Energy Refinement), son muy utilizados para el estudio de algunos sistemas orgánicos, proteínas y ácidos nucleicos. En el caso de los CHARMM (Chemistry at Harvard Macromolecular Mechanics), se trata de campos de fuerzas parametrizados con datos experimentales y utilizados tanto para moléculas pequeñas como para complejos solvatados de macromoléculas biológicas. Uno de los campos de fuerza más extendidos debido a su gran reproducibilidad en sistemas con varios grupos funcionales es el MMFFs (Merck Molecular Force Field). Este último es el que se ha elegido para realizar este trabajo, por lo que se describe con más detalle a continuación. 34 5.1.2.1. MERCK MOLECULAR FORCE FIELD El campo de fuerzas MMFF (Merck Molecular Force Field) está diseñado para su uso en compuestos farmacéuticos y considera de forma precisa las energías conformacionales y las interacciones no enlazantes. Es adecuado para estudiar moléculas tanto en fase gas como en fases condensadas. Posee un gran número de términos cruzados, lo que es la razón principal de su gran transferibilidad, es decir, las propiedades calculadas no diferirán en gran medida al cambiar las condiciones en las que se encuentre la molécula. 34,35 En este caso se calcula la energía como la suma de las siguientes contribuciones. • Tensión de enlaces: Se utiliza la ley de Hooke a la que se le añaden términos adicionales del desarrollo en serie de Taylor. 𝐸𝑒𝑛𝑙𝑎𝑐𝑒 =𝑘𝑒𝑛𝑙𝑎𝑐𝑒(𝑟𝑖𝑗 −𝑟𝑖𝑗,𝑒𝑞)2[ 1 +𝑐𝑠(𝑟𝑖𝑗 −𝑟𝑖𝑗,𝑒𝑞)+ 7 12 𝑐𝑠2(𝑟𝑖𝑗 −𝑟𝑖𝑗,𝑒𝑞)2] (1) donde Kenlace , es la constante de fuerza, rij , es la distancia entre átomos, rij, eq es la distancia de equilibrio entre los átomos, y cs es la constante de estiramiento cúbica. 16 • Deformación de ángulos de valencia: Se utiliza una función análoga a la del caso anterior en la que se sustituyen las distancias de enlace por ángulos entre enlaces. 𝐸á𝑛𝑔𝑢𝑙𝑜 =𝐾  (𝜃𝑖𝑗𝑘 − 𝜃𝑖𝑗𝑘,𝑒𝑞)2[1+𝑐𝑏(𝜃𝑖𝑗𝑘 − 𝜃𝑖𝑗𝑘,𝑒𝑞)] (2) donde K es la constate de fuerza, 𝜃𝑖𝑗𝑘 es el ángulo de enlace entre los átomos i, j , y k , 𝜃𝑖𝑗𝑘,𝑒𝑞 es el ángulo de equilibrio, y cb es la constante cúbica de doblado. Para ángulos próximos a la linealidad se tiene la siguiente aproximación: 𝐸á𝑛𝑔𝑢𝑙𝑜,𝑙𝑖𝑛𝑒𝑎𝑙 =𝐾𝑖𝑗𝑘,𝑙𝑖𝑛𝑒𝑎𝑙 (1 + 𝑐𝑜𝑠  𝑖𝑗𝑘) (3) • Tensión-Doblado. 𝐸𝑒𝑠𝑡𝑖𝑟𝑎𝑚𝑖𝑒𝑛𝑡𝑜−𝑑𝑜𝑏𝑙𝑎𝑑𝑜 =[𝑘𝑖𝑗𝑘 (𝑟𝑖𝑗 −𝑟𝑖𝑗,𝑒𝑞)+𝑘𝑘𝑗𝑖 (𝑟𝑖𝑗 −𝑟𝑖𝑗,𝑒𝑞) ] (𝜃𝑖𝑗𝑘 − 𝜃𝑖𝑗𝑘,𝑒𝑞) (4) donde 𝑘𝑖𝑗𝑘 y 𝑘𝑘𝑗𝑖 son las constantes de fuerza que unen los estiramientos ij y kj al ángulo ijk . • Doblado fuera del plano (OOP, Out of plane). 𝐸𝑂𝑂𝑃 =𝑘𝑜𝑜𝑝 (  𝑖𝑗𝑘;𝑙)2 (5) donde el término  𝑖𝑗𝑘;𝑙 es el ángulo entre el enlace jl y el plano ijk , en el que j es el átomo central. • Torsiones: Debido a que la torsión es un movimiento periódico (circular), la función que la describe también lo es. 𝐸𝑡𝑜𝑟𝑠𝑖ó𝑛 =1 2 [𝑉1(1+𝑐𝑜𝑠  )+ 𝑉2(1+𝑐𝑜𝑠2  )+𝑉3(1+𝑐𝑜𝑠3  )] (6) donde los términos V son constantes de fuerza de la serie de Fourier y  es el ángulo diedro. • Interacciones de Van der Waals: Surge de las interacciones entre las nubes electrónicas de dos átomos no enlazados. 𝐸𝑉𝑑𝑊 =ϵij(1.07 𝑅𝑖𝑗 ∗ 𝑅𝑖𝑗 +0.07𝑅𝑖𝑗 ∗)7( 1.12 𝑅𝑖𝑗 ∗7 𝑅𝑖𝑗 7+0.07𝑅𝑖𝑗 ∗7−2) (7) donde 𝑅𝑖𝑗 es la distancia entre los átomos i y j , 𝑅𝑖𝑗 ∗ es la distancia cuya energía de interacción entre átomos es mínima y ϵij es la profundidad del pozo de energía potencial entre átomos. • Interacciones electrostáticas: Describen la interacción culómbica entre átomos con cargas parciales. 17 𝐸𝑖𝑗 = 𝑞𝑖𝑞𝑗 𝐷(𝑅𝑖𝑗 +  ) 𝑛 (8) donde D es la constante dieléctrica,  es la constante dieléctrica de amortiguación, y los términos q las cargas de los respectivos átomos. • Términos cruzados y términos no enlazantes adicionales: Describen el acoplamiento entre enlaces, ángulos y torsiones por lo que los funcionales que los describen son combinaciones de los funcionales para los términos por separado. Sus expresiones son más complejas. • Estrategias de parametrización: Son funciones que miden de la desviación entre los valores de energía predichos y los valores experimentales. La elección de estas funciones es arbitraria. 5.2. MÉTODOS MECANOCUÁNTICOS Los métodos de la mecánica cuántica se basan en la resolución de la ecuación de Schrödinger, la cual sólo puede resolverse de manera exacta en el caso en el que se tenga un sistema monoelectrónico. Debido a que en sistemas de más de un electrón las partículas no se mueven de forma independiente, si no que existe correlación electrónica, deben considerarse varias aproximaciones para resolver la ecuación de Schrödinger. 34 Las aproximaciones consisten en: • Aproximación de Born-Oppenheimer: Se basa en que la masa de los núcleos es mucho mayor que la de los electrones, por eso, los electrones se mueven a mayores velocidades respondiendo instantáneamente a los cambios producidos en la molécula. Con esta aproximación es posible suprimir los términos nucleares considerándolos constantes al calcular las derivadas, y parametrizar el término de atracción electrónica entre núcleos y electrones. • Modelo Hartree-Fock: La mayor parte de los métodos mecanocuánticos utilizan el cálculo Hartree-Fock como base. La metodología empleada consiste en resolver el Hamiltoniano con un solo determinante de Slater. Estos cálculos se fundamentan en el método de Hartree, que aplica la aproximación de la partícula independiente para los electrones, de manera que permite dividir el operador Hamiltoniano en operadores individuales para cada electrón. Aplicando la aproximación de la partícula independiente, el término de un electrón se considera independiente del resto. De esta forma se puede considerar que un electrón se mueve en el campo de energía potencial promedio formado por el resto. Es una modificación en la que se describe la función de onda de varios electrones como un producto antisimetrizado (determinante de Slater) de funciones de onda de un electrón. 18 • El comportamiento de los orbitales moleculares se aproxima al resultado de la combinación lineal de los orbitales atómicos que lo forman. De esta forma, los electrones se pueden considerar independientes asumiendo que cada uno estará en un orbital diferente. Así, la función de onda del sistema se puede escribir como el producto de las funciones de onda de cada partícula. • Aproximación de la separación de orbitales : Consiste en la separación de las contribuciones  y  a un enlace, de forma que permite caracterizar la contribución . 5.2.1. MÉTODOS BASADOS EN LA FUNCIÓN DE ONDA Tras plantear la ecuación de Schrödinger y aplicar la aproximación de Born-Oppenheimer, se divide el Hamiltoniano electrónico en distintos términos, de los cuales, el correspondiente a la repulsión interelectrónica no tiene solución exacta. Se emplea entonces el modelo de Hartree-Fock. El concepto fundamental en el que se basa la teoría Hartree-Fock (HF) consiste en que cada electrón percibe al resto como un campo promedio, lo que permite un enorme progreso en la realización de los cálculos de orbitales moleculares. Según el teorema variacional, la energía calculada resulta siempre mayor que la auténtica energía que posee el sistema, por lo que debe procurarse minimizar el resultado para cometer el menor error posible. Para ello se plantea el determinante secular cuya resolución puede realizarse por distintos métodos. 33 Los métodos post-HF realizan el cálculo de la energía de correlación. De nuevo, se encuentran distintas metodologías como variacionales, perturbativos y de agregados acoplados. Para este trabajo se ha empleado el método MP2 que detallaremos a continuación. 36 5.2.1.1. MÉTODOS PERTURBATIVOS MOLLER-PLESSET La teoría de perturbaciones de varios cuerpos (MBPT, many-body perturbation theory) trata la correlación electrónica como una perturbación en la función de onda de HF. Se comienza considerando el operador Hamiltoniano 𝐻  no perturbado y al que se le añade posteriormente una pequeña perturbación 𝑉  que representa una pequeña modificación en el sistema. El operador Hamiltoniano en este caso vendría dado por la expresión (9). 𝐻  =𝐻 0+ 𝜆 𝑉  (9) donde 𝐻 0 es el Hamiltoniano anterior a la introducción de la perturbación que se resuelve de forma exacta o aproximada, 𝜆 es el parámetro de perturbación que mide la magnitud de la perturbación, y 𝑉  es el operador de perturbación. 19 Se asume que el factor de corrección es pequeño comparado con el Hamiltoniano inicial por lo que la función de onda modificada y energía pueden expresarse en forma de la expresión de Taylor en función del parámetro de perturbación. 34 𝐸𝑛=𝐸0 𝑛+ 𝜆 𝐸0 𝑛+𝜆2𝐸2 𝑛+⋯ (10) Ψ𝑛=Ψ0 𝑛+ 𝜆 Ψ0 𝑛+𝜆2Ψ2 𝑛+⋯ (11) 5.2.2. MÉTODOS SEMIEMPÍRICOS Los métodos semiempíricos modifican los cálculos Hartree-Fock introduciendo funciones con parámetros obtenidos experimentalmente, que mejoran la calidad del cálculo computacional. Se basan en tres aproximaciones principales. 34 • Despreciar la contribución de los electrones internos: Dado que no contribuyen a la actividad química, suelen reemplazarse sus contribuciones al Hamiltoniano por una función parametrizada. • Uso de un conjunto de base mínimo: Las funciones utilizadas para describir los electrones de valencia emplean un número mínimo de funciones base. • Supresión o uso de aproximaciones para las integrales de dos electrones: Los métodos actuales poseen la aproximación modified neglect of differential overlap (MNDO). En este método, los parámetros se asignan para los distintos átomos y se encajan para reproducir propiedades como los calores de formación, variables geométricas, momentos dipolares y primeras energías de ionización. 5.2.2.1. PARAMETRIC MODEL 6 Existen múltiples modelos como el AM1, PM1, PM3… que van mejorando sucesivamente la parte MNDO. El método semiempírico que se utilizará en este trabajo es el parametric model 6 (PM6), que emplea un conjunto de datos mucho mayor para ajustar los parámetros e introduce mejoras en los términos núcleo-núcleo. Las correcciones que aplica se basan en introducir parámetros de pares de núcleos en lugar de parámetros específicos de un elemento, así como diferentes potenciales de repulsión internuclear para N-H, O-H, C-C y Si-O. También añade orbitales d para algunos elementos de forma similar al método MNDO/d.37 20 5.2.3. TEORÍA DEL FUNCIONAL DE LA DENSIDAD (DFT). El método de la teoría del funcional de la densidad se basa en que la energía total del estado fundamental de un sistema con N electrones se puede describir mediante un funcional de densidad electrónica. Los métodos actuales están basados en dos teoremas de Hohenberg y Kohn. El primero de ellos enuncia que cada observable de un sistema mecano-cuántico estacionario, incluyendo la energía, puede ser calculado a partir de la densidad electrónica del estado fundamental únicamente. El segundo afirma que la densidad electrónica del estado fundamental es aquella que minimiza el funcional de la energía, es decir, al igual que con la teoría de orbitales moleculares (MO, molecular orbitals) la densidad electrónica cumple el principio variacional.34,38 Sin embargo, los teoremas de Hohenberg y Kohn no proporcionan una expresión para el funcional de la energía y la utilidad de los cálculos DFT depende del uso de aproximaciones adecuadas. Para ello, el funcional se describe como la energía total según la expresión de Hartree añadiendo un término adicional, el de intercambio-correlación. Ésta es la expresión de Kohn-Sham (12) que consta de varios términos, todos ellos funcionales de la densidad electrónica, que comprenden el sumatorio de las interacciones de las distintas partículas dentro de la molécula. 𝐸[𝜌]=𝑇𝑠 [𝜌]+𝐸𝐻[𝜌]+𝐸𝑒𝑛[𝜌]+𝐸𝑛𝑛 [𝜌]+𝐸𝑥𝑐[𝜌] (12) donde 𝑇𝑠 [𝜌] es la energía cinética de los electrones sin interacción entre ellos, 𝐸𝐻[𝜌] es la energía de Hartree, 𝐸𝑒𝑛[𝜌] la energía de interacción electrón-núcleo, 𝐸𝑛𝑛 [𝜌] la energía de interacción núcleo-núcleo y 𝐸𝑥𝑐[𝜌] el funcional de intercambio-correlación.39 Existen varios tipos de funcionales de intercambio-correlación que sirven como aproximaciones detallando en mejor o peor medida la forma de la densidad electrónica. Se pueden clasificar de la siguiente manera: • LDA (Local Density Aproximation): Es el más sencillo de todos, se calcula el valor del funcional de intercambio-correlación en un punto r como el valor de la densidad electrónica ρ en dicha posición. Sería equivalente a tener un gas de electrones libres en esa densidad. Para los sistemas que incluyan polarización de spin se utiliza el funcional LSDA (Local Spin Density Aproximation). Tiene una expresión como la siguiente, (13). 𝐸𝑥𝑐[𝜌]=∫𝑑3𝑟 𝜌(𝑟)𝜀𝑥𝑐(𝜌(𝑟)) (13) 21 • GGA (Generalized Gradient Approximation): Esta aproximación se basa en la anterior que sólo tiene en cuenta los valores locales de la densidad electrónica, y añade el cambio local de la misma, es decir, su gradiente. 𝐸𝑥 𝑐 ⁄ 𝐺𝐺𝐴[𝜌(𝑟)]=𝐸𝑥 𝑐 ⁄ 𝐿𝑆𝐷𝐴[𝜌(𝑟)]+ Δ 𝜖𝑥 𝑐 ⁄ [|∇ρ(r) | 𝜌43 ⁄(𝑟)] (14) • Meta-GGA: Una de las deficiencias de los funcionales GGA es su incapacidad para describir correctamente el comportamiento de la densidad de energía y el potencial de intercambio simultáneamente en regiones alejadas del núcleo para obtener una mejora. Estos funcionales añaden a los anteriores la segunda derivada de la densidad electrónica aplicando el operador Laplaciano a la expresión. Un ejemplo de estos funcionales son los de tipo M06. • DFT Híbridos: Dado que un funcional puro muestra poca contribución de spin, se añaden contribuciones del funcional de intercambio proporcionado por la función de onda de Hartree-Fock. Los funcionales que se van a utilizar son de tipo GGA hibrido B3LYP e híbrido meta-GGA M06-2X. 36,40 5.2.3.1. FUNCIONAL B3LYP Son métodos híbridos que combinan funcionales para la correlación y para el intercambio. En este caso concreto, se tiene el funcional de intercambio de Becke (B) que añade tres parámetros a la expresión y el de correlación no local de Lee, Yang y Parr (LYP), que a diferencia de otros, no es una corrección del funcional LSDA (Local Spin Density Approximation), sino que está diseñado para proporcionar la energía de correlación por sí mismo. 41 5.2.3.2. FUNCIONAL M062X Como se ha mencionado antes, el funcional M062X se puede clasificar como un funcional híbrido de tipo meta GGA (hybrid meta-GGA). Tienen muy en cuenta el cambio del valor de la densidad electrónica en la molécula y además duplica el intercambio no local (2X). Se ha parametrizado para no metales, por lo que se recomienda para el estudio de interacciones no covalentes. 5.2.3.3. DISPERSIONES DE GRIMME Y DE BECKE JOHNSON La dispersión de Grimme es una modificación realizada a la teoría del funcional de la densidad de Kohn-Sham, que proporciona mayor precisión, mayor rango de aplicación y menor empirismo. Las principales modificaciones respecto de los métodos anteriores son los coeficientes de dispersión específicos para pares de átomos y radios de corte, ambos computados por primeros principios. La principal desventaja es que la corrección no es dependiente o no afecta a la estructura electrónica. 28 5.3.3. ROTACIÓN INTERNA La rotación interna es una vibración en la que el movimiento realizado por la molécula conecta varias disposiciones (n) equivalentes entre sí, siendo el ángulo de rotación interna (α) la coordenada que describe dicho movimiento. La función del potencial es simétrica y periódica con el ángulo α, siendo el máximo de potencial la barrera que impide el movimiento. Sin embargo, en los casos en los que la barrera es baja, el sistema puede pasar de un mínimo a otro por medio del efecto túnel. Este suceso causa que los niveles cercanos a la barrera se dividan en n componentes. Un caso concreto que va a aplicarse en este trabajo es la rotación de los grupos metilo. En este caso da lugar a una función de energía potencial en tres partes, con tres mínimos equivalentes y tres máximos equivalentes. Si la barrera de potencial es suficientemente pequeña, dicha rotación puede dar lugar a desdoblamientos hiperfinos en el espectro de rotación, dependiendo de la magnitud de la barrera de rotación interna. En esos casos, los niveles rotacionales se dividirán en un nivel no degenerado de componente A y un nivel doblemente degenerado de componente E, como se observa en la Figura 14 50. Entonces, cada transición rotacional aparece desdoblada en dos componentes. Como resultado de la presencia de rotación interna, el espectro obtenido se complica haciendo más compleja su resolución. 51,52 Figura 13: Formación y zonas del chorro supersónico. Figura 14: Función del potencial proporcionado por la rotación interna de un grupo metilo al cambiar el ángulo  29 6. RESULTADOS 6.1. BÚSQUEDA CONFORMACIONAL En este trabajo se ha estudiado la molécula de hernandulcina, cuya fórmula molecular es C15H24O2, y que posee en su estructura 41 átomos en total. Tal y como se muestra en la Figura 15, se trata de una estructura grande y con numerosos enlaces simples, por lo que cabe esperar que adopte una gran cantidad de confórmeros. Por esta razón, primero se recurre a los métodos de mecánica molecular para dilucidar las estructuras que poseerá la molécula, para posteriormente optimizarlas mediante distintos métodos. Su tamaño y flexibilidad le dotan de una importante complejidad desde el punto de vista espectroscópico y computacional. Por ello, con este estudio se pretende obtener una información básica que después podrá utilizarse para valorar la viabilidad de la molécula, conocer sus interacciones intramoleculares y así tener un punto de partida en el caso de realizar su estudio experimental. También será de ayuda en el caso de hacerse estudios posteriores como por ejemplo, los agregados moleculares del compuesto con agua o moléculas claves en el receptor del dulzor. Se comienza construyendo la estructura de la molécula con el programa Maestro de Schrödinger 53, en el que se pueden visualizar todas las posibles torsiones de la molécula como se indica en la Figura 15. Los indicadores azules representan giros a través de un enlace y en amarillo se indica la rotura de un enlace con el fin de modificar la estructura del ciclo. Seguidamente se procede a lanzar un cálculo de mecánica molecular empleando el campo de fuerzas de Merck (MMFF) y se establece un umbral de energía de 25 KJ/mol por encima del cual no se consideren nuevas estructuras. Una vez terminada la búsqueda conformacional se obtienen 109 confórmeros, cada uno de ellos descrito en un fichero que proporciona las coordenadas de los átomos, parámetros energéticos y otras propiedades de cada confórmero. Posteriormente, se optimizan las estructuras de los confórmeros obtenidos con un nivel de cálculo DFT, utilizando el funcional B3LYP y una base extendida, 6-311++G(d,p), añadiendo las dispersiones de Grimme y Becke Johnson. Este cálculo se realiza con el software Gaussian 16 54. Al comprobar los Figura 15: Posibles torsiones que pueden efectuarse en la estructura de la hernandulcina. 30 parámetros de todas las estructuras, se observa que las energías y las constantes de rotación A, B y C coinciden en algunas estructuras. Este hecho indica que esos confórmeros convergen a la misma estructura, por lo que realmente se obtienen 81 estructuras diferentes de las 109 iniciales. Tabla 1: Parámetros espectroscópicos calculados y energías relativas a nivel B3LYP/6-311++G(d,p) GD3BJ, correspondientes a los 10 confórmeros más estables de la hernandulcina. Los dos confórmeros menos estables se muestran a la derecha. A, B, y C representan las constantes rotacionales; µa, µb, y µc son los valores de las componentes del momento dipolar eléctrico (en Debyes). Por último, la energía electrónica relativa, E, la energía relativa en el punto cero, EZPE, y energías de Gibbs relativa (G) calculadas a 298 K son relativas con respecto al mínimo global. Los valores han sido calculados mediante B3LYP/6-311++G(d,p) GD3BJ. En la Figura 16 se muestran las estructuras de los diez confórmeros más estables junto con los dos menos estables, mientras que en la Tabla 1 se recogen los datos de las constantes de rotación, momentos dipolares y energías relativas para dichas estructuras. Los datos del resto de los confórmeros se encuentran en el anexo. La numeración de las estructuras sigue el proporcionado según el orden energético por el método de mecánica molecular. Como puede verse en la Figura 16, las estructuras más estables son aquellas que poseen un enlace de hidrógeno entre el grupo hidroxilo y el grupo cetona. Aquellas estructuras que no cuentan con dicho enlace, sufren de una desestabilización considerable, como es el caso de los confórmeros 91 y 96. Confórmero 5 1 66 6 4 14 3 57 102 11 91 96 A(MHz) 501 688 886 536 931 553 945 1157 636 523 598 941 B(MHz) 429 246 209 392 192 358 189 206 325 445 328 200 C(MHz) 290 196 186 263 177 290 170 193 286 322 278 191  A(D) 1.56 0.81 2.35 -0.20 -2.80 2.92 2.17 3.89 2.06 -3.50 -0.04 2.69  B(D) 4.19 5.29 4.38 5.06 4.51 2.82 5.02 3.41 3.53 1.16 -2.41 2.56  C(D) -2.72 -1.66 -1.24 -2.20 1.00 -3.46 -1.02 -1.08 -3.19 3.05 2.53 1.37 E (cm-1) 0 531 551 459 861 807 894 936 1025 1053 2725 2816 EZPE (cm-1) 0 470 489 491 776 779 798 906 1030 1045 2512 2586 G (cm-1) 0 140 112 528 318 782 335 755 983 1068 2204 2002 31 Confórmero 5 Confórmero 1 Confórmero 66 Confórmero 6 Confórmero 4 Confórmero 14 Confórmero 3 Confórmero 57 Confórmero 102 Confórmero 11 Confórmero 91 Confórmero 96 Figura 16: Los diez confórmeros de la hernandulcina más estables y los dos menos estables hallados con el método B3LYP/6-311++G(d,p) GD3BJ. 32 Entre las estructuras más estables es posible diferenciar dos grandes grupos: las estructuras que tienen la cadena carbonada plegada sobre el ciclo (confórmeros 5 y 6, por ejemplo) y aquellos cuya estructura se dispone de una forma más extendida (confórmeros 1 y 66, por ejemplo). La diferencia entre ambos grupos resulta muy interesante y, como se mencionará más adelante, dará lugar a diferencias entre los distintos métodos de cálculo. El grupo de estructuras de conformación plegada adopta esta configuración gracias a una estabilización debida a fuerzas dispersivas entre la cadena carbonada y el anillo. Por otro lado, las estructuras de conformación extendida son más complejas de evaluar, pero parece que su disposición se ve favorecida por una interacción muy débil de los hidrógenos con los oxígenos. La estabilidad final vendrá dada por un compromiso entre la estabilización por la suma de las fuerzas estabilizantes y las fuerzas repulsivas ejercidas por el impedimento estérico. Para corroborar esta observación entre los distintos grupos, se realiza un NCI Plot (Non-Covalent Interactions) a las estructuras 5 y 66 por ser las más estables de cada grupo. El NCI Plot es un programa que permite visualizar e identificar las interacciones no covalentes. Los cálculos están basados en la densidad electrónica y sus derivadas (picos que aparecen en el gradiente de densidad reducido a bajas densidades). El análisis de las interacciones se muestra en la Figura 17, en la que puede verse cómo el confórmero 5 presenta unas fuerzas dispersivas entre la cadena carbonada y el anillo, mientras que el confórmero 66 no posee dichas fuerzas, lo que conlleva una estabilización adicional de la primera. De hecho, la estabilización es tal que la estructura 5 es 500 cm-1 más baja en energía que el siguiente confórmero más estable. Confórmero 5 Confórmero 66 Figura 17: Confórmeros 5 y 66 de la hernandulcina junto con el análisis NCI Plot a color. Las superficies en azul corresponden a fuerzas atractivas, las repulsiones se muestran en rojo y en tonos verdes se representan las interacciones débiles. 33 6.1.1. COMPARATIVA DE BASES Y NIVELES DE CÁLCULO El éxito de un modelo depende de su capacidad para reproducir consistentemente los datos experimentales, y no es probable que un único modelo teórico sea el ideal para describir todos los sistemas moleculares. Si bien todos estos modelos pueden ser aplicados habitualmente a moléculas de tamaño considerable, difieren en el tiempo de computación requerido y en la eficacia para reproducir un sistema. De forma que es de gran importancia saber qué modelos proporcionan mejores resultados en menores periodos de tiempo o en qué situaciones son necesarios los modelos de mayor coste computacional. Es difícil cuantificar el tiempo de computación total de un cálculo debido a que depende de varios factores. Para moléculas orgánicas de tamaño moderado (de unos 10 átomos diferentes de H), los modelos HF/6-31G(d), B3LYP/6-31G(d) y MP2/6-31G(d) es de esperar que exhiban tiempos de computación global en una proporción de aproximadamente 1:1.5:10. El tiempo de computación de Hartree-Fock y los modelos B3LYP y MP2 con bases mayores que 6-31G(d) aumentan aproximadamente con el cubo (HF y B3LYP) y la quinta potencia (MP2) del número total de funciones de base. 33 La hernandulcina se presenta como un caso muy didáctico de cara a ilustrar la necesidad de las distintas metodologías. Por un lado, es una molécula relativamente grande por lo que nos mostrará los problemas relacionados con el tiempo de cálculo y tamaño de archivo generado durante el cálculo. Por otro lado, dispone de dos grupos predominantes del cual uno de ellos presenta estructuras estabilizadas mediante fuerzas dispersivas de largo alcance difíciles de modelar. Por ello y con el fin de realizar una comparativa entre distintos métodos computacionales disponibles, se toman las 30 estructuras más estables de las totales optimizadas y se reoptimizan como sigue: 1) Con el nivel B3LYP aplicando distintas bases: o B3LYP/6-31G(d) o B3LYP/6-31G(d,p) o B3LYP/6-311G(d) o B3LYP/6-311G(d,p) o B3LYP/6-31+G(d) o B3LYP/6-31++G(d) o B3LYP/6-311+G(d) o B3LYP/6-31+G(d,p) o B3LYP/6-31++G(d,p) o B3LYP/6-311++G(d) o B3LYP/6-311+G(d,p) o B3LYP/6-311++G(d,p) GD3 o B3LYP/6-311++G(d,p) GD3BJ o B3LYP/Def2TZVP 34 2) Con distintos niveles empleando la base 6-311++G(d,p) como referencia: o M06-2X, correspondiente a metodología DFT. o PM6, de tipo semiempírico. o Método Hartree-Fock (HF). o MP2, que utiliza metodología ab initio post-HF. Una vez terminados los cálculos, con el fin de contrastar la información se crea un script con Wolfram Mathematica55 para que automatice la tarea de extraer el tiempo de optimización promedio (tiempo dedicado a la optimización por cada número de pasos necesarios para optimizar y por CPU). Esta representación del tiempo requerido frente al nivel o base empleado se muestra en la Figura 18. Así mismo, la Figura 19 muestra una representación del tiempo necesario para el cálculo de frecuencias frente al nivel y base utilizado. Estos datos nos darán una idea estimada del coste computacional de cada método. Los resultados numéricos se muestran en la Tabla 2. Tabla 2: Datos obtenidos para el tiempo de cálculo de los distintos niveles de cálculo/base. t. opt/ h:min:s t. frec/ h:min:s PM6/6-311G++(d,p) 0:00:01 0:00:03 B3LYP/6-31G(d) 0:01:49 0:25:29 B3LYP/6-31G(d,p) 0:01:46 0:23:44 B3LYP/6-311G(d) 0:02:10 0:28:02 B3LYP/6-311G(d,p) 0:02:56 0:40:17 B3LYP/6-31+G(d) 0:04:12 0:51:36 B3LYP/6-31++G(d) 0:04:23 0:48:00 B3LYP/6-311+G(d) 0:05:11 0:56:26 B3LYP/6-31+G(d,p) 0:05:45 1:16:44 B3LYP/6-31++G(d,p) 0:05:58 1:12:10 B3LYP/6-311++G(d) 0:06:39 1:10:56 B3LYP/6-311+G(d,p) 0:07:33 1:28:27 HF/6-311++G(d,p) 0:06:19 1:03:19 B3LYP/6-311++G(d,p) GD3 0:11:42 2:35:47 B3LYP/6-311++G(d,p) GD3BJ 0:22:45 4:18:04 M062X/6-311++G(d,p) 0:16:19 4:10:09 B3LYP/Def2TZVP 0:15:18 3:00:43 MP2/6-311G++(d,p) 4:48:02 10:35:07 35 Figura 19: Tiempo necesario para el cálculo de frecuencias por CPU para cada nivel de teoría y base. Figura 18: Tiempo necesario para la optimización de geometrías para cada nivel de teoría y base, por paso y por CPU. 36 Como puede observarse, cualquier pequeña mejora que se introduzca en las bases aumenta ligeramente el tiempo de cálculo. El tiempo empleado en la base Pople más completa aumenta cuatro veces respecto de la más sencilla, lo cual parece bastante razonable teniendo en cuenta que la mejora es sustancial. No es extraño que la base 6-311++G(d,p) haya sido y sea ampliamente utilizada. Otro detalle importante es el apenas aumento del tiempo de cálculo de un nivel HF a un B3LYP. Teniendo en cuenta que los B3LYP, así como otros métodos DFT, suelen proporcionar muy buenos resultados generalmente, no es de extrañar que su desarrollo recibiera el Premio Nobel en 1998. Finalmente existe un incremento de tiempo notable, aunque no exagerado, al aplicar las dispersiones de Grimme, que tienen en cuenta las interacciones de largo alcance o los funcionales que las incorporan (M06-2X). En breve se discutirá si este incremento está justificado. En relación al método semiempírico PM6, si bien es cierto que reproduce un aparente orden energético prácticamente igual a como lo hacen los métodos DFT más completos en un tiempo varios ordenes de magnitud menor, presenta varios problemas. En primer lugar, no proporciona un resultado coherente en cuanto a las energías relativas. Mientras que los métodos DFT que incluyen dispersiones sitúan al confórmero más estable 400-500 cm-1 más bajo en energía que el siguiente, el PM6 obtiene que hay hasta 14 confórmeros por debajo de ese límite, lo cual es un factor que daría lugar a conclusiones erroneas de cara a la cantidad de estructuras que se espera encontrar en el experimento. Es decir, a pesar de que el tiempo de cálculo es muy rápido (segundos), los resultados son muy “sospechosos” y casi todas las estructuras obtenidas son isoenergéticas. Además, tal y como se muestra en la Tabla 3 para la estructura 66, las geometrías son muy “pobres” y, mientras que todos los métodos proporcionan unas constantes rotacionales similares, el método PM6 proporciona un valor de la constante rotacional A un 10% fuera del valor predicho por el resto de los métodos. Esto afectaría enormemente de cara a una asignación conformacional, sobre todo para las búsqueda de transiciones tipo b (las más intensas) y c. Todo ello indica que el método no es muy fiable en el contexto que nos movemos, que requiere una precisión notable tanto en los valores energéticos como en los estructurales. Esto ocurre con varias estructuras. Curiosamente, la discrepancia entre las geometrías es debido a una sobreestimación de las interacciones, como por ejemplo la interacción C-H···O-H en la estructura 66 que la precide a 2.35 Å tal y como se puede apreciar en la Figura 20. La diferencia es de casi un Amstrong con respecto al resto de los métodos. 37 Tabla 3: Parámetros espectroscópicos calculados con distintos métodos de cálculo, empleando la misma base, 6-311++G(d,p), correspondientes al confórmero 66. A, B, y C representan las constantes rotacionales; µa, µb, y µc son los valores de las componentes del momento dipolar eléctrico (en Debyes). Finalmente, los cálculos MP2 presentan varios problemas: por un lado, el tiempo empleado para la optimización de la geometría es muy alto, siendo al menos un orden de magnitud superior a los DFT más exigentes; además, hay que tener en cuenta que el número de pasos requerido para optimizar fue superior al resto y que, en algunos casos, la convergencia no era posible sin un criterio de seguimiento de las constantes de fuerza para la minimización de la estructura; por otro lado, el método no paraleliza bien, por lo que el empleo del doble de procesadores no conlleva a la mitad del tiempo necesario, en contraposición con los métodos DFT donde sí que se cumple esta regla; por último, el scratch necesario además de requerir un tamaño considerable (cada cálculo requiere de 0.15 TB durante el cálculo de la optimización y de 0.5 TB durante el cálculo de frecuencias), requiere además de una gestión del sistema de colas especial: la escriturá en disco es uno de los mayores problemas en los clusters, y el modo óptimo no se arregla añadiendo un disco más grande. La solución es añadir un sistema de discos (cluster) que puedan escribir la información en varios de ellos a la vez, ya que la escritura en disco, incluso para los discos SSD, es muy pesada y si se tiene dos o más calculos escribiendo en un único disco, es muy normal que los cálculos “mueran”. Este hecho se comenta Confórmero B3LYP-GD3BJ M06-2X MP2 PM6 A(MHz) 885.9 912.5 887.6 1004 B(MHz) 208.7 214.5 214.2 212 C(MHz) 185.7 190.6 190.1 194  A(D) 2.3 2.3 2.3 2.8  B(D) 4.4 4.0 3.8 4.3  C(D) -1.2 -1.3 -1.2 -1.5 B3LYP/6-311++G(d,p) GD3BJ PM6/6-311++G(d,p) Figura 20: Distancias entre el hidrógeno del grupo metilo ubicado al final de la cadena carbonada y el grupo hidroxilo, proporcionada por los niveles de cálculo B3LYP y PM6, usando la base 6-311++G(d,p), para la estructura 66. 44 Figura 24: Espectro completo predicho para el confórmero 5 de la hernandulcina. Figura 25: Transiciones de tipo a en la predicción del espectro 5 de la hernandulcina. 45 6.2.1. EVALUACIÓN DE LA EXISTENCIA DE ROTACIÓN INTERNA Si observamos la estructura de la molécula, vemos que existen cuatro grupos metilo que pueden dar lugar a rotación interna por el giro entorno al enlace sencillo que une el carbono metílico al resto de la molécula. Para comprobar si va a tener lugar la rotación interna a través de alguno de ellos, se realiza un escaneo de la superficie de energía potencial relajada (SEP) seleccionando un ángulo diedro que incluya el enlace de giro y los dos enlaces adyacentes como se muestra en la Figura 27. Las SEP proporcionan la energía de la molécula al variar el ángulo  de cada diedro, con lo que se obtienen gráficas de energía potencial cuya forma es periódica, y el máximo indica la barrera de energía potencial necesaria para la rotación. La altura de la barrera energética para cada caso se recoge en la Tabla 6, las gráficas se muestran en la Figura 28. Figura 27: Estructura de la hernandulcina donde se resaltan los ángulos diedros con una posible rotación interna. Figura 26: Transiciones de tipo c en la predicción del espectro de la hernandulcina, donde se muestra el valor de su separación. 46 Tabla 6: Datos de los diedros para los que se estudia la rotación interna. Como dato indicativo, la rotación interna se produce cuando la altura de la barrera es inferior a unos 300 cm-1. En los datos obtenidos para la hernandulcina, la energía máxima proporcionada por la rotación supera ampliamente en todos los casos el límite por el que puede producirse este efecto. Por ello, se puede decir que, afortunadamente, no existirá el fenómeno de rotación interna en esta molécula. Únicamente puede discutirse en el caso del diedro B, sin embargo, es poco probable que ocurra. Esta información es una buena noticia de cara a la resolución del espectro, ya que no se espera que las señales salgan desdobladas por la ruptura de la degeneración como se explica en el punto 5.3.3., lo cual complicaría la asignación además de que también tendría efectos en la intensidad de las señales, ya que disminuiría al repartirse entre los distintos estados y para cada metilo. Diedro A B C D Energía máxima 790 cm-1 380 cm-1 1090 cm-1 640 cm-1 Figura 28: Superficie de energía potencial relajada de las cuatro posibles torsiones que pueden dar lugar a rotación interna. 47 6.3. CORRESPONDENCIA CON EL TRIÁNGULO DEL DULZOR Sabiendo que la hernandulcina es dulce cabe preguntarse si la teoría de Shallenberger-Acree-Kier puede explicar esta relación estructura-propiedad. Así, una vez que se han determinado las estructuras más estables, y por lo tanto más relevantes, es posible comprobar si existen los puntos necesarios para formar el glucóforo en dichas estructuras. En principio, los grupos cetona e hidroxilo cumplen con las características de electronegatividad y distancias requeridas y corresponden a los puntos AH/B. Por ello, faltaría localizar dónde se encuentra el punto . Midiendo las distancias entre átomos en cada una de las estructuras tridimensionales más estables, y por lo tanto más relevantes, se comprueba qué átomo o grupo de átomos puede actuar como punto de anclaje hidrófobo con el receptor dando lugar al glucóforo propuesto por la teoría. Las distancias para el punto  de acuerdo a la teoría de Shallenberger-Acree-Kier se indican en la Figura 5 y se recuerdan en la Tabla 7. Esta medición se realiza para los confórmeros que aparecen como los más estables según todas las metodologías empleadas, y que además resultan ser los de menor energía relativa de Gibbs. Este último criterio es importante puesto que tiene en cuenta la estabilidad a una temperatura dada, en este caso a 298 K, y por ello la más relevante en el cuerpo humano. Los resultados se recogen en la Tabla 7, y el triángulo del dulzor para las estructuras más relevantes se muestra en Figura 29. Tabla 7: Distancias entre los tres puntos de interacción del glucóforo con el receptor medidos en Angstrom. Distancia A-B A- B- Teórica 2.60 3.50 5.50 Confórmero 1 2.71 3.58 5.71 Confórmero 3 2.71 3.51 5.98 Confórmero 4 2.71 3.51 5.95 Confórmero 5 2.71 3.40 5.14 Confórmero 66 2.71 3.57 5.71 48 Figura 29: Triángulo del dulzor indicado en los cinco confórmeros más estables. Confórmero 1 Confórmero 3 Confórmero 4 Confórmero 5 Confórmero 66 49 A la vista de los resultados, se llega a la conclusión de que el punto hidrófobo  se encuentra entre uno de los átomos de carbono que forma el doble enlace y el carbono alílico en la cadena carbonada. Para hallar las distancias aproximadas se toma como referencia un punto entre los carbonos indicados, representado en color rosa en las imágenes de la Figura 29. Como se puede observar, este punto no tiene por qué corresponder a un átomo o grupo funcional concreto, sino simplemente pertenecer a una zona de fuerzas dispersivas en la molécula cuyas distancias con los puntos AH/B sea adecuada. Es interesante ver como el NCI Plot de la Figura 17, refuerza esta observación y se observa que existen fuerzas dispersivas fuertes en esa región, lo que también lleva a pensar que será un punto de interacción favorable con el receptor responsable del dulzor en las papilas gustativas. Según los estudios de Compadre et al. 27, cuando se sintetizan derivados de esta molécula en los que no existe el doble enlace en la cadena carbonada, como se puede ver en la Figura 30, la sustancia se vuelve amarga, lo cual corrobora las conclusiones aquí obtenidas. Finalmente, se puede concluir que la molécula encaja con la teoría del sabor dulce de ShallenbergerAcree-Kier y posee un glucóforo que consta de tres puntos de anclaje con los receptores de las papilas gustativas. No solo es eso si no que, como se ha mencionado anteriormente, todas las estructuras relevantes en términos de energía libre de Gibbs cumplen dicha teoría. Estas son las estructuras que estarán interconvirtiéndose continuamente a la temperatura del cuerpo humano. Si el receptor resultase ser muy específico (modelo llave-cerradura), el hecho de tener diversas formas en las que pueden encontrarse los confórmeros (plegada, extendida…) que cumplen con los requisitos, haría que al menos una de ellas encaje. Es decir, la adaptación de la hernandulcina sería excelente. Si el receptor no fuese tan específico, entonces todas las estructuras relevantes cumplirían con los requisitos de la teoría. En cualquier caso, ambos hechos podrían explicar el elevado dulzor de la hernandulcina. Hernandulcina Derivado amargo 1 Derivado amargo 2 Figura 30: Estructura de la hernandulcina en comparación con la de dos de sus derivados amargos. 50 7. CONCLUSIONES Se ha realizado el estudio del panorama conformacional de la hernandulcina obteniendo 81 confórmeros estables por debajo de 25 KJ/mol. Entre los confórmeros obtenidos, se pueden diferenciar tres grupos. El primer grupo contiene aquellos confórmeros que presentan un enlace de hidrógeno y además se estabilizan adicionalmente por fuerzas dispersivas entre el anillo y la cadena carbonada, lo que las sitúa entre las estructuras más estables. El segundo grupo también contiene confórmeros estabilizados por enlaces de hidrógeno, pero disponen la cadena carbonada de una forma más extendida, que los desestabiliza ligeramente. Finalmente, el tercer grupo contiene aquellos confórmeros que no poseen un enlace de hidrógeno en su estructura, lo que les lleva a una alta desestabización. Gracias a esta búsqueda, se han determinado las estructuras más relevantes. Al ser una molécula grande, se ha usado la hernandulcina como modelo para comparar el coste computacional y fiabilidad de distintos niveles de cálculo. Los métodos B3LYP tienen un coste computacional similar al método Hartree-Fock. Además, al introducir funcionales que tienen en cuenta interacciones no-covalentes como la dispersión de Grimme, Becke-Johnson o al usar M06-2X que las tiene en cuenta, el coste computacional aumenta sólo ligeramente, pero con resultados más acordes. Los métodos semiempíricos son muy rápidos, pero no estiman correctamente las diferencias energéticas ni las geometrías, de manera que no sería correcto tomarlos como referencia a la hora de realizar una predicción. Por último, se demuestra que los cálculos MP2 requieren de mucho tiempo de cálculo, por lo que su uso no es viable para moléculas de gran tamaño. Se ha comprobado la viabilidad de la realización de un experimento de espectroscopía de rotación para la molécula. Por un lado, se ha demostrado que la molécula no presentará rotación interna. Por otro lado, se ha predicho el espectro que cabría esperar para dos de las estructuras más relevantes. Curiosamente, estas dos estructuras presentan distintos límites de simetría: una estructura es oblate y la otra prolate. En ambos casos, se encuentran una gran cantidad de transiciones debido a que las moléculas son activas en los tres ejes principales de inercia y se hace una predicción del espectro que servirá de cara a una posible experimentación futura. Las transiciones más intensas serán las de tipo b. Finalmente, se ha dilucidado cuáles serán los tres puntos necesarios para formar el glucóforo en la hernandulcina según la teoría del dulzor de Shallenberger-Acree-Kier. Los grupos cetona e hidroxilo cumplen con las características de electronegatividad y distancias requeridas y corresponden a los puntos AH/B. Por otro lado, el punto  se encuentra entre uno de los átomos de carbono que forma el doble enlace y el carbono alílico en la cadena carbonada. Todas las estructuras más estables poseen el triángulo del dulzor además de poseer distintas disposiciones espaciales, lo que explicaría el elevado dulzor de la hernandulcina. 51 8. BIBLIOGRAFÍA 1 C. M. Hladik, B. Simmen and P. Pasquet, Primatological and anthropological aspects of taste perception and the evolutionary interpretation of ”basic tastes”., Anthropology, 2003, 41, 67–74. 2 G. K. Beauchamp, Why do we like sweet taste: A bitter tale?, Physiol. Behav., 2016, 164, 432– 437. 3 P. Mergenthaler, U. Lindauer, G. A. Dienel and A. Meisel, Sugar for the brain: the role of glucose in physiological and pathological brain function, Trends Neurosci., 2013, 36, 587–597. 4 I. Ramirez, Why do sugars taste good?, Neurosci. Biobehav. Rev., 1990, 14, 125–134. 5 P. Jiang, J. Josue, X. Li, D. Glaser, W. Li, J. G. Brand, R. F. Margolskee, D. R. Reed and G. K. Beauchamp, Major taste loss in carnivorous mammals, Proc. Natl. Acad. Sci., 2012, 109, 4956– 4961. 6 The Science of Taste, https://foodinsight.org/the-science-of-taste/, (accessed 3 April 2019). 7 Preferencia por lo dulce, https://www.tellmegen.com/results/rasgos/preferencia-dulce/, (accessed 3 April 2019). 8 Reglamento (CE) no 1333/2008 sobre aditivos alimentarios, 2008. 9 Reglamento de la Unión Europea 1129/2011, 2011. 10 P. J. Rogers, The role of low-calorie sweeteners in the prevention and management of overweight and obesity: evidence v. conjecture, Proc. Nutr. Soc., 2018, 77, 230–238. 11 R. Chauhan, Taste Masking: A Unique Approach for Bitter Drugs, J. Stem Cell Biol. Transplant. , 2017, 1, 12. 12 S. S. Bharate, S. B. Bharate and A. N. Bajaj, Interactions and incompatibilities of pharmaceutical excipients with active pharmaceutical ingredients: a comprehensive review., J. Excipients Food Chem., 2010, 1, 3–26. 13 S. S. Bharate, S. B. Bharate and A. N. Bajaj, ChemInform Abstract: Interactions and Incompatibilities of Pharmaceutical Excipients with Active Pharmaceutical Ingredients, ChemInform, 2011, 42, 3–26. 14 K. Srikanth, K. Priya and V. R. Mohan-Gupta, Natural sweeteners: A complete review, J. Pharm. Res., 2011, 4, 2034–2039. 15 S. D. Roper and N. Chaudhari, The cell biology of taste, J. Cell Biol., 2010, 190, 285–296. 16 D. Liu, N. Archer, K. Duesing, G. Hannan and R. Keast, Mechanism of fat taste perception: Association with diet and obesity, Prog. Lipid Res., 2016, 63, 41–49. 17 T. Yamamoto, Brain mechanisms of sweetness and palatability of sugars., Nutr. Rev., 2003, 61, S5-9. 18 QIAGEN - Sweet Taste Signaling, https://www.qiagen.com/us/shop/genes-andpathways/pathway-details/?pwid=425, (accessed 9 April 2019). 19 R. S. Shallenberger, The AH,B glycophore and general taste chemistry, Food Chem., 1996, 56, 209–214. 52 20 R. S. Shallenberger, Taste chemistry, Springer-Sciende+Business Media, B.V., 1st edn., 1993. 21 F. Bruni, C. Di Mino, S. Imberti, S. E. McLain, N. H. Rhys and M. A. Ricci, Hydrogen Bond Length as a Key to Understanding Sweetness, J. Phys. Chem. Lett., 2018, 9, 3667–3672. 22 L. B. Kier, A Molecular Theory, J. Pharm. Sci., 1972, 61, 1394–1397. 23 C. Bermúdez, I. Peña, S. Mata and J. L. Alonso, Sweet Structural Signatures Unveiled in Ketohexoses, Chem. - A Eur. J., 2016, 22, 16829–16837. 24 E. R. Alonso, I. León, L. Kolesniková and J. L. Alonso, The Structural Signs of Sweetness in Artificial Sweeteners: A Rotational Study of Sorbitol and Dulcitol, ChemPhysChem, 2018, 19, 3334–3340. 25 E. R. Alonso, Biomolecules and interstellar molecules: structure, interactions and spectroscopic characterization, Tésis Dr. 26 Hernandulcin | C15H24O2 | ChemSpider, http://www.chemspider.com/ChemicalStructure.111731.html, (accessed 22 March 2019). 27 C. M. Compadre, R. A. Hussain, J. M. Pezzuto, R. L. de Compadre Lopez and A. D. Kinghorn, Analysis of structural features responsible for the sweetness of the sesquiterpene, hernandulcin, Experientia, 1988, 44, 447–449. 28 A. Arias Contreras, Medicinal Plants at Pura Vida Spa & Yoga Retreat: Sweet Aztec herb (Hierba dulce, orozus), http://wellnessplant.blogspot.com/2011/11/sweet-aztec-herb-hierba-dulceorozus.html, (accessed 1 May 2019). 29 F. A. Souto-Bachiller, M. De Jesus-Echevarría, O. E. Cárdenas-González, M. F. Acuña-Rodriguez, P. A. Meléndez and L. Romero-Ramsey, Terpenoid composition of Lippia dulcis, Phytochemistry, 1997, 44, 1077–1086. 30 A. D. Kinghorn and E. J. Kennelly, Discovery of Highly Sweet Compounds from Natural Sources, J. Chem. Educ., 2009, 72, 676–680. 31 C. M. Compadre, J. M. Pezzuto, A. D. Kinghorn and S. K. Kamath, Hernandulcin: An intensely sweet compound discovered by review of ancient literature, Science (80-. )., 1985, 227, 417–419. 32 C. M. Compadre, R. A. Hussain, R. L. Lopez de Compadre, J. M. Pezzuto and A. D. Kinghorn, The Intensely Sweet Sesquiterpene Hernandulcin: Isolation, Synthesis, Characterization, and Preliminary Safety Evaluation, J. Agric. Food Chem., 1987, 35, 273–279. 33 T. Engel, P. Reid and W. Hehre, in Química Física, Pearson Educación S.A., 2006, pp. 597–656. 34 K. I. Ramachandran, G. Deepa and K. Namboori, Computational Chemistry and Molecular Modeling, Springer-Verlag Berlin Heidelberg, 2008. 35 A. D. McNaught and A. Wilkinson, in IUPAC Compendium of Chemical Terminology, Blackwell Scientific Publications, Oxford, Research Triagle Park, NC, 2nd edn., 1997. 36 C. J. Cramer, Essentials of Computational Chemistry, John Wiley & Sons Ltd, 2004. 37 A. S. Christensen, T. Kubař, Q. Cui and M. Elstner, Semiempirical Quantum Mechanical Methods for Noncovalent Interactions for Chemical and Biochemical Applications, Chem. Rev., 2016, 116, 5301–5337. 38 P. Hohenberg and W. Kohn, Inhomogeneous Electron Gas, Phys. Rev., 1964, 46, 864–871. 53 39 M. Bockstedte, A. Kley, J. Neugebauer and M. Scheffler, Density-functional theory calculations for poly-atomic systems: electronic structure, static and elastic properties and ab initio molecular dynamics, Comput. Phys. Commun., 1997, 107, 187–222. 40 L. Singh, David J., Nordstrom, in Planewaves, Pseudopotentials and the LAPW Method, Springer US, 2nd edn., 1994, pp. 5–21. 41 Density Functional (DFT) Methods, http://gaussian.com/dft, (accessed 12 March 2019). 42 S. Grimme, J. Antony, S. Ehrlich and H. Krieg, A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu, J. Chem. Phys., 2010, 132, 1–20. 43 A. D. Becke and E. R. Johnson, A density-functional model of the dispersion interaction, J. Chem. Phys., , DOI:10.1063/1.2065267. 44 S. Grimme, S. Ehlich and L. Goerikg, Effect of the Damping Function in Dispersion Corrected Density Functional Theory, J. Comput. Chem., 2011, 32, 1456–1465. 45 J. Foresman and Æ. Frisch, Explor. Chem. with Electron. Struct. Methods, 1996, 302. 46 F. Weigend and R. Ahlrichs, Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: Design and assessment of accuracy, Phys. Chem. Chem. Phys., 2005, 7, 3297. 47 J. L. Alonso and J. C. López, in Encyclopedia of Signaling Molecules, 2015, pp. 2–4. 48 P. Atkins and J. de Paula, Química Física, Editorial Médica Panamericana, 8th edn., 2008. 49 Stefano Rampino, Schematic diagram of a supersonic expansion from a nozzle. The change... | Download Scientific Diagram, https://www.researchgate.net/figure/Schematic-diagram-of-asupersonic-expansion-from-a-nozzle-The-change-of-Mach-number-M_fig1_243480825, (accessed 4 July 2019). 50 P. Pinacho Morante, Microsolvation of biomolecular models by microwave spectroscopy: structure and cooperative effects., Tésis Dr. 51 C. C. Lin and J. D. Swalen, Internal Rotation and Microwave Spectroscopy, Rev. Mod. Phys., 1959, 31, 841–892. 52 V. D. G. Lister, J. N. Macdonald and N. L. Owen., Internal Rotation and Inversion. An Introduction to Large Amplitude Motions in Molecules. Academic Press, London 1978., Angew. Chemie, 1979, 91, 450–451. 53 Maestro | Schrödinger, https://www.schrodinger.com/maestro, (accessed 1 May 2019). 54 Gaussian.com | Expanding the limits of computational chemistry, https://gaussian.com/, (accessed 10 May 2019). 55 Wolfram Mathematica: Computación técnica moderna, http://www.wolfram.com/mathematica/, (accessed 1 July 2019). 56 JB95 Spectral fitting program | NIST, https://www.nist.gov/services-resources/software/jb95spectral-fitting-program, (accessed 4 July 2019).