Menu

A few general quaries.

SOUMEN
2023-09-22
2023-10-09
  • SOUMEN

    SOUMEN - 2023-09-22

    Query 1: A multi-well reaction that we are currently studying has all its transition states below the isolated reactant. The highest barrier lies below 17 kcal/mol from the isolated reactant and the rest are around ~ 50 kcal/mol below the isolated reactant energy. In that reaction, while dealing with intermediate uni-molecular step we are using Eckart tunneling with RRKM technique.
    We face difficulties like Chemically significant eigenvalues (CSE) not well separated from internal energy relaxation eigenvalues (IEREs) and also at low temperatures for a few steps Bartis-Widom Phenomenological Rate Coefficients become zero.
    That problem vanishes when tunneling is removed from the uni-molecular step.
    (So far we have understood no flux is settled down in the stationery points and all the flux reaches one of the final products while tunneling is considered; that may cause the problem).

    So is there an option instead of using RRKM we supply temperature-dependent TST rate constant (CVT/SCT) externally for those intermediate steps so that the mesmer can process them??
    Although RRKM is the best we can have to make use of the energy release in each step to propagate further.

    Query 2: Here the situation is kind of the opposite. In this reaction one barrier is above the isolated reactant energy (~ 6 kcal/mol) and the rest of them are below the isolated reactant energy. This reaction is an atmospheric one. how do choose the energy above the top hill here?
    If say I allow an excess of ~10KT above the top hill (i.e. it is above the highest energy TS which is above the isolated reactant energy); how to justify that excess energy in case of ordinary atmospheric reaction case ??
    (is it better to take 0 KT ??)
    Or is there any option to include temperature-dependent CVT/SCT rate directly here externally ??
    That may solve the problem of excess energy.

    Thanks in advance. If necessary I will attach input files on the upcoming conversation.

     
    • struanr

      struanr - 2023-09-23

      Hello,

      Thank-you for using MESMER.

      Query 1: If the CVT/SCT rate coefficient can be well fitted to a modified Arrhenius form, then it can be incorporated into MESMER using the inverse Laplace transform (ILT) method. This method will take the canonical rate coefficient and obtain from it an effective microcanonical rate coefficient function that can be used in a master equation calculation.

      Query 2: The maximum energy cut-off that should be used in a MESMER calculation can be set automatically using the keyword me:automaticallySetMaxEne , however, it is wise to check this by doing a number of different calculations with different maximum energies until the lower eigenvalues have converged. Incorporating a CVT/SCT canonical rate coefficient could be done by the same procedure as described for Query 1.

      For what it is worth, atmospheric systems are usually run at low temperatures and it maybe necessary to go to higher precision.

      I would be happy to take a look at these systems and try to give further advice.

      With regards, Struan

       
  • SOUMEN

    SOUMEN - 2023-09-26

    Thank you very much for your reply.
    Here I am attaching two files. Both are usually the same. One of them has Eckart tunneling in on (File1_withT) and in the other one Eckart is off (File2_withoutT).
    The system is like the Query 1. All the barriers are significantly below the isolated reactant energy. Also, there is a branching of products.
    From the output, it can be observed when tunneling was considered one of the branches getting zero Bartis-Widom phenomenological rate constants. Whereas when tunneling was not considered thats the same branch got a very small Bartis-Widom phenomenological rate constants (as that one has a relatively large barrier height i.e. consistent with PES).

    So I want your valuable suggestion about the consideration of tunneling. is it really necessary to consider the tunneling correction when this PES is considered? If so; then how to justify the zero B-W PRC.

    I will send the PES regarding Query 2 in near future. Also, thanks in advance.

     
    • struanr

      struanr - 2023-09-28

      Hello,

      "From the output, it can be observed when tunneling was considered one of the branches getting zero Bartis-Widom phenomenological rate constants. Whereas when tunneling was not considered thats the same branch got a very small Bartis-Widom phenomenological rate constants (as that one has a relatively large barrier height i.e. consistent with PES).

      So I want your valuable suggestion about the consideration of tunneling. is it really necessary to consider the tunneling correction when this PES is considered? If so; then how to justify the zero B-W PRC."

      Thank-you for sending your files. I have looked at them and I believe I have located the issue and it concerns the imaginary frequency associated with reaction transition state TS-5 which is entered as 78.05 cm-1. This value is very small and because of this the parameter "c" (see W.H. Miller, JACS, 101(23), 1979) is very large, too large in fact to be represented by a standard double precision variable and as a consequence the tunnelling coefficients, and so therefore the rate coefficients, are set to zero. This small value of the imaginary frequency should be checked because it is a measure of the curvature of the barrier and it suggests that this barrier has a small curvature which does not seem to be consistent with the height of the barrier. I attach a .zip file containing the results I obtained when I re-run your examples but for the tunnelling case I increased the imaginary frequency to 780.05 cm-1 as a test, you will see this gives more reasonable results.

      Incidentally, I made some minor changes to the definition of reaction R1 expressing the MesmerILT section in the current format.

      Finally, I note there are a number of negative rate coefifcients. I suspect that these occur because the barriers between the species Add-I, Add-II and Add-III are very low. Indeed the barrier heights are less than barrier to internal rotation in ethane (about 10 kJ/mol). For such low barriers, even at the low temperature of 240 K, equilibrium between these species will be established very quickly and under such circumstances it is hard to define meaningful rate coefficients. I would suggest that you consider treating these species as a single species and combining their states.

      I hope this helpful.

      With regards, Struan

       
  • SOUMEN

    SOUMEN - 2023-09-29

    Thank you very much for your reply.
    So, you were correct. I somehow put a wrong value of imaginary frequency for TS-5. The actual value of it is around 1006.17 cm-1. That solves the issue with tunneling.
    Now to solve the issue with negative rate coefficients your suggestion is to treat intermediate species having low barriers as a single one and combine their states. I want your guide here.
    If I understood correctly one possible option is to remove Add-II and Add-III and directly connect Add-I with PC-2 via TS-5. That solves the problem of negative rate coefficients. Although initially, we tried to remove Add-II and connect Add-I with Add-III via TS-3 but still some negative coefficients were still there. Is there any other way?
    Also, how wrong my result is when a few of the rate coefficients are negative?
    I am attaching the corrected Mesmer input file where the wrong imaginary frequency is corrected (namely File1_withT) for your reference.
    Also, I am attaching another two files. Both are modified accordingly to remove the negative rate coefficients.
    The file ‘File1_withT_mod1’ is modified from ‘File1_withT’ by removing the small barrier step as mentioned above. Also, this one is for your reference only.
    The second file namely ‘File2_withT_mod1’ has barrierless steps connecting R1 and R2 with PC-1 without any TS and MESMER ILT is applied here (i.e. a new channel is added). So far, all the necessary corrections and modifications have been made to this file. But still, this has the problem of negative rate coefficients (only two negative values). Is the problem arising due to the wrong declaration of the MESMER ILT parameter or so ? (Negative rate coefficients arise between the interconversion of two final product complexes (PC-1 and PC-2)). I set the final products to the sinks (that may cause inconvenient results). Here my goal is to know the rate coefficients of individual branches or channels.
    When I remove ‘sinks’ from all the final PC and set them to ‘modelled’ the problem of negative rate coefficients vanishes but I lose the information about the rate coefficients of two different channels. I want to report rate coefficients in two different channels in from the Mesmer results.

     
    • struanr

      struanr - 2023-10-01

      Hello,
      “So, you were correct. I somehow put a wrong value of imaginary frequency for TS-5. The actual value of it is around 1006.17 cm-1. That solves the issue with tunneling.”

      I am pleased to read this.

      “Now to solve the issue with negative rate coefficients your suggestion is to treat intermediate species having low barriers as a single one and combine their states. I want your guide here. If I understood correctly one possible option is to remove Add-II and Add-III and directly connect Add-I with PC-2 via TS-5. That solves the problem of negative rate coefficients. Although initially, we tried to remove Add-II and connect Add-I with Add-III via TS-3 but still some negative coefficients were still there. Is there any other way?”

      It depends on what the species Add-I, Add-II and Add-III are. For instance, if they are conformers connected by an internal rotation barrier then it is possible to model this motion as hindered internal rotation.

      “Also, how wrong my result is when a few of the rate coefficients are negative?”

      This is an interesting question and depends on how rate coefficients are interpreted. The prevailing interpretation (the one you find in textbooks) is that individual rate coefficients describe the rate of change of one specific species to another and so all rate coefficients must be positive. There is some debate about this and some discussion that it is the eigenvalues describing the relaxation that are important and the details of the rate coefficients are less important. However, I would suggest that one should endeavour to achieve the former.

      “The file ‘File1_withT_mod1’ is modified from ‘File1_withT’ by removing the small barrier step as mentioned above.”

      This seems reasonable to me, providing we account for the other isomers/confers as discussed above.

      “The second file namely ‘File2_withT_mod1’ has barrierless steps connecting R1 and R2 with PC-1 without any TS and MESMER ILT is applied here (i.e. a new channel is added). So far, all the necessary corrections and modifications have been made to this file. But still, this has the problem of negative rate coefficients (only two negative values). Is the problem arising due to the wrong declaration of the MESMER ILT parameter or so ? (Negative rate coefficients arise between the interconversion of two final product complexes (PC-1 and PC-2)). I set the final products to the sinks (that may cause inconvenient results). Here my goal is to know the rate coefficients of individual branches or channels.
      When I remove ‘sinks’ from all the final PC and set them to ‘modelled’ the problem of negative rate coefficients vanishes but I lose the information about the rate coefficients of two different channels. I want to report rate coefficients in two different channels in from the Mesmer results.”

      I think the problem here is that MESMER is getting confused by the different channels to the PC-1 product. There might be a bug, but I am not sure and will need to investigate this some more. In the meantime, the workaround is to define sink species with different names but are in fact are the same which is safe to do for this example if both reactions to form PC-1 are irreversible. I attach the results that I have obtained by adding a new sink PC-1b. Not surprisingly, most of the input flux appears in one form or other of PC-1 and of these PC-1b is the larger component which is also expected. I ran the calculation at double-double and quad-double precision and did not see any substantial difference. There are still some negative rate coefficients but these are very small in magnitude and can be considered numerical noise.

      With regards, Struan

       
  • SOUMEN

    SOUMEN - 2023-10-03

    Thank you for your time and effort to clarify all my doubts. I believe we are at the end of this conversation trail.
    So far I have learned a lot from our conversation.
    I have my last queries regarding the PES we have been talking about so far.
    The file "File1_withT" contains the full PES for your reference. In the PES Add-I, Add-II and Add-III are rotational isomers. A COOH group rotates about an angle slightly larger than 180° ( dihedrally) via TS-2 and TS-3. we already discussed output of this file has quite a few negative rate constants and that is due to those small rotational barriers (TS-2 and TS-3).
    As a way out you suggested the use of the hindered internal rotation model. Here my question is even if I model the Add-I, Add-II, Add-III TS-2, and TS-3 with a hindered internal rotation model will the problem of negative rate coefficients be solved??

    Or is it better to remove Add-II, Add-III, TS-2, and TS-5 and connect directly TS-5 with Add-I as shown in file 'File2_withT_modified'? So far this approach seems best to me because only two negative rate coefficients were there and they were very small and can be ignored.
    Although you told me about it previously I am confirming about it.

    For now, considering we did all the calculations by removing all those low-barrier rotamers then how should I justify this to the reviewer or anyone else?
    Since we can see that all the important rate constants are almost the same for both the output and also seems reasonable to us. Then the shortcut process as shown in file 'File2_withT_modified' be okay to use.

     
    • struanr

      struanr - 2023-10-04

      Hello,

      “I have my last queries regarding the PES we have been talking about so far.The file "File1_withT" contains the full PES for your reference. In the PES Add-I, Add-II and Add-III are rotational isomers. A COOH group rotates about an angle slightly larger than 180° ( dihedrally) via TS-2 and TS-3. we already discussed output of this file has quite a few negative rate constants and that is due to those small rotational barriers (TS-2 and TS-3). As a way out you suggested the use of the hindered internal rotation model. Here my question is even if I model the Add-I, Add-II, Add-III TS-2, and TS-3 with a hindered internal rotation model will the problem of negative rate coefficients be solved??”

      Without trying one cannot know for sure, but I strongly suspect that it would eliminate the negative eigenvalues, because the internal rotational wells would be accounted for within the density of states. There is a facility within MESMER to treat hindered rotors, this approach will need information about the hindered rotor potential, which you appear to have already, and the internal rotor moment of inertia which can be obtained from the geometry of the adduct.

      “Or is it better to remove Add-II, Add-III, TS-2, and TS-5 and connect directly TS-5 with Add-I as shown in file 'File2_withT_modified'? So far this approach seems best to me because only two negative rate coefficients were there and they were very small and can be ignored. Although you told me about it previously I am confirming about it.”

      If the motion is an internal motion it would be best to treat it as such.

      “For now, considering we did all the calculations by removing all those low-barrier rotamers then how should I justify this to the reviewer or anyone else? Since we can see that all the important rate constants are almost the same for both the output and also seems reasonable to us. Then the shortcut process as shown in file 'File2_withT_modified' be okay to use."

      If you use the hindered rotor approach above, this include the rotamers correctly and thus easy to justify as this is how virtually all calculations are done for systems with internal rotations are done. I strongly recommend you consider this.

      I hope this helps.

      Regards, Struan

       
  • SOUMEN

    SOUMEN - 2023-10-09
     

    Last edit: SOUMEN 2023-10-09
  • SOUMEN

    SOUMEN - 2023-10-09

    Dear Struanr;
    Thank you very much for your reply and guidance. I am glad to say after considering hindered rotors treatment for R1, R2, Add-I, Add-II, Add-III, TS-2, and TS-3 we have been able to get rid of most of the negative frequencies. There are only two (or we can say one) negative rate coefficients so far. One is R1 (deficient reactant) to ADD-I (association product) , and the other one is the reverse of it i.e., ADD-I to R1. We can very easily neglect the second one as it is of the order of 10^-40. I will be glad if you look into the file and provide more insight about it.

    I am not attaching the file here as we are considering publishing the result. As you are aware the hindered rotor treatment required all the optimized geometry so it exposes all the details of PES to the public. Therefore I request you provide your email ID or any other communication channel so that I can send the input file.

    With Regards
    Soumen

     
    • struanr

      struanr - 2023-10-09

      Dear Soumen,

      You can send the input files to struanhrobertson@outlook.com.

      With regards, Struan

       

Log in to post a comment.