 Okay, thanks. So it's my pleasure to start the second day of this open force field meeting I will talk in the next 45 minutes about electrostatics So just briefly introduce myself. So I'm Mikhail Schalpel. I'm here in the ghost lab in San Diego and I'm here on the fellowship from the Austrian Science Fund and in that fellowship I proposed to work on polarizable force fields and that's why I'm working on polarization and I did that for the first couple of months and then we found out that actually a lot of the things I'm doing for the polarizable force fields are also useful for fixed charge models That's why I also involved in that part of the open force field effort And that's the reason why I'm here today and present to you the ideas for our charge models So I want to start off with first slide you already saw yesterday And it's the general overview slide just to give you an idea where we are So we are now Defining our force field. So setting what kind of parameters we have and what kind of parameters we want to optimize in our optimizer So now a bit more detailed what I will cover in this talk So I will first talk about fixed charge models in Smirnov and our plans for implementing new fixed charge models There is like our Short-term goal a new version of REST, which we call REST2 Afterwards, Sue from Leaping's lab will briefly introduce their REST pipe program and then I will continue on with how we want to come up with a new generation of AM1 BCC's and In the second part then I will talk about polarizable charge models where we also have our short-term solution based on REST and one solution based on AM1 BCC so as I said, I want to start off with the fixed charge models and To do that, I would just want to give you a brief overview. What's the current status of fixed charge models? So as Smirnov is somehow a sibling of Gaff you can run them with the same charge models as Gaff so REST and AM1 BCC's and You usually use for our calculations mainly AM1 BCC up to now With the new or generation of Smirnov We also want to have new and better charge models and one idea is based on REST So that should be a new implementation of REST and the other thing is a new generation of AM1 BCC's So I want to start off with the REST based approach To do that. I want to briefly Summarize the idea of REST So the idea is we are using a gas phase QM calculation at a quite low level. So hard refog 631 g star And then we calculate the electrostatic potential around A molecule and then we fit charges which are reproducing the electrostatic potential And the nice thing is that actually hard refog with that low basis set is somehow over polarizing the gas phase so that we can use those charges for condensed phase simulations um That's actually just a lucky error compensation and not really Very well defined. So probably we can come up with a scheme which covers polarization better than this approach And While we are not using higher just higher level QM calculations The problem is if we are doing that and do that in the gas phase We end up with charges which are too low. So we are not accounting for the polarization in condensed phase And if we are running an implicit solvent calculation So the same molecule just with implicit solvent and use the dielectric constant of order we end up with charges Which are too high Because we do not account for something like polarization costs or so Um, so what can we do when we have one charge sets? That's too low too low and the other one is too high We can mix them together And our idea for rest two is uh that we are using two high level QM calculations So one is based on uh is a gas phase calculation one is an implicit solvent calculation Then we are doing two individual rest fits and then we Mix those two charge sets together Based on a mixing parameter So we are taking a fraction delta from the gas phase charges and a fraction one minus delta of the implicit solvent calculations And all that should then be implemented in within respite in the near future I want to talk about two main things now. First, I want to talk what kind of QM Calculation we want to use and the other thing is I want to talk about the mixing parameter delta and how we can fit this mixing parameter I will start off with What would be nice or what kind of QM level we want to use? So it would be nice if we could just use a very Good method like cc std Unfortunately, that's pretty expensive and we don't want to do that um But luckily a lot of people studied cc std calculation for small molecules and they uh compared other functional forms or other QM methods to the cc std result and they found out that for at least for electronic properties Double hybrids are also quite accurate double hybrids are still uh quite expensive and we don't want to use them for Uh direct like molecules or parameterize a lot of molecules So we are just using this double hybrid as our gold standard and compare what kind of cheaper functional is still good and is somehow Agrees with the double hybrid calculations and to do that We tested different functionals and different basis sets On a test set of 71 molecules and then compared the arrows to the double hybrid Calculations and you see here on the left side just what kind of molecules we included in our test set um We then did all those calculations and first calculate the molecular dipole moment And what we saw is okay if we compare to double hybrids Hardware foc is uh doing another really good job. So the arrows are quite high and for mp2 and um The free tested functions we are doing more or less Similar so it seems that uh, they're the dft functions are similar in the currency than mp2 Yes um That's the first property we calculated the other thing is we actually want to use the electrostatics potentials Uh for fitting charges. So we compared how the asp's are doing again. Hardware foc is doing not a great job mp2 is actually quite good as long as we're using our high enough basis set and Actually the function is doing only slightly worse than mp2, but are Much cheaper and it seems that when we compare the result of both Uh properties and also account how we want our how the speed of the method is we Think that with something like a bw6 e 95 with our augmented double Uh, see the basis that we can get quite good performances Uh, the double hybrids esp Was rough Yeah, uh, so the reference was uh double hybrid Uh This pbp 86 with quadruple set basis set Holy that the people from zoom can also hear what you asked Are you doing a norm over all esp grid points? Uh, yes. Okay. Okay. I understand. Thank you And i'm using uh emergency Coleman seeing uh kind surface Um, yeah, that was was what I want to tell you about the the cure method The next thing I want to talk about is the mixing parameter delta So we want to use as a starting value 0.5 There is some physical Justification to use 0.5 and it's also the ipal cure approach But we won't just use that as a starting value and then uh co-optimize this Mixing parameter together with the letter jones parameters using leaping force balance And all the experimental data we can access through the property calculator and the nice thing with that Parameter is that we are only introducing one additional parameter in the optimization But we can use that parameter to tune the polarity of our force fields And that's quite a nice idea we think Um to sum up this Rest two approach I just want to highlight the pros and cons about it So we think that we should be able to describe these be better Due to the higher level of qm We have that nice parameter to fit the the polarity against condensed phase experimental data And uh, we should be able to capture the polarizabilities of molecule more accurate because we are not Relying on the error compensation of hundred fox six three one g Uh, the downside is if we are using that new charge model with the current letter jones parameters We probably don't get better results. So we have to refit them Uh, but we should be able to do that within the open force field effort And the other thing is uh, the higher level qm calculations are more expensive And I will talk about how we can Tackle this problem that they are computed computational a bit more demanding in a few minutes um Yes john This is really cool stuff Does the delta have to be a constant for all molecules? Could you learn some sort of representation that maps individual molecules to individual deltas based upon functional groups or something like that? That's a good question. I want to say you this is something you already mentioned Okay, um actually At least at the beginning. I think it is nice to tune just one to get us in a right starting point and The problem is if you have different deltas for different functional groups when you have then a molecule with two Different functional groups. How do you combine them? so I think we can talk about a bit later when I explained what our plans for the bcc's are Then uh, you see where we can probably individual tune The for the electrostatics if I could just follow up on that. I think again like we talked about that if you Have their points of different functional groups. You won't maintain uh the charge of the molecule without imposing another constraint But yeah, why not it? so, uh It should think so is it yes, okay? so, uh, I I understand how you figured out your uh, your vacuum side of the calculation But I think I missed how you were doing the condensed phase component like you're mixing the condensed phase as well Did you mention how you were doing the condensed phase? calculation the condensed phase should be just an implicit solvent calculation with Our dielectric constant of 80. Okay, and which uh, which dielectric method are you using? Um, we still have to test what fits our needs best. Okay, maybe you have strong opinions on that invariably Okay, um, I just want to briefly mention our estimated timeline for this part of the project So we are currently Deciding what qm methods we want to use so what should be our reference and also what kind of implicit solvent model we want to use then we Want to test our to get a baseline how good our charge model does without refitting the mixing parameter and the lennard jones parameters And within the end of the year, we hope that we can refit lennard jones and mixing parameters And hopefully implement that within a year and then we want to improve that charge model until The end of the open force field effort and Yeah, now it's time for sue who is presenting leaping's respite approach because it fits better into the resp Is this working? Yeah, I think it works. Yeah so So the rest method and the am 1 bcc Are considered to be the standard method for atomic Atomic partial charge calculations But the existing tools are not interoperable with other software tools or not open enough So our main concepts are implementation of open source version of resp method In python and adding new features to improve our current charge models So i'm initially focused on reproducing the previous results And the package we encoded is called respite And these are major features of respite package user can set the restraints For the fittings and they can fit charges for multiple molecules and multiple conformers at the same times And also they can select which kind of which Greed selection schemes they they are going to use and also they can set the flexible constraint on charges for example, they can force symmetry On residues or across the residues So this is data data flow chart for of the package They have two parts one is ESP generators and yet the another one is resp optimizer So if we input the coordinates of molecules and the conformers with a setting for a molecule in With the input files containing molecule informations and grid settings ESP generator can take the inputs and generate the molecule object and they generate the ESP grid points around the molecules And they run cipher calculation to calculate the electrostatic potential and electric field on each grid And then they didn't they generate the ESP data Which has an appropriate directory structures for the resp optimizer and the resp optimizer takes the ESP data as an input and With the input file containing the models and constraints user going to use They can generate the fitted charges for each molecules So using the respite package we did we did some comparison between Comparison of performance between different charge models and these are preliminary results we got so For the comparison of the performance between the different charge models There are several ways we can try one is the one is the rendering of residuals Between the QM ESP and the MM ESP based on the fitted charges and to uh check if we can Reproduce the figures in the Christopher's 1993 paper. We We re-rendered the residuals around the Methanol molecules From the different charge fitting method So if you see these figures are the Colors of the its spheres are the values of the values of the ESP Point on each grid points and the radius means the magnitude of the residuals So if you see the like white big spheres on the Upper right figures that means there's a large residuals on the non-polar areas So by rendering these figures we can compare between different charge models and The if you see the upper left UNFR is the unrestrained fitting with with No force symmetries and the UNAP means the unrestrained fitting with Averaging charges after fittings and UNEQ means the unrestrained fitting with force symmetries And wk.fr st eq is the standard two-stage fit Method we um people usually use so you can see the two-stage fit can um improve the improve the charge fitting Uh results Especially around the polar areas when you compare the two-stage fitting figures with the UN dot eq And also we can calculate the relative salvation energy differences between two different set of charges Um to compare the performance between different charge models and uh, we Just for the simplicity of the notation. I just used LASP two RASP one RASP one and RASP two And RASP two means the two-stage fit to electrostatic potentials and RASP one means single stage fitting And RASP one means single stage fitting to electric field So when we try to calculate the relative salvation energy differences on methanol molecule You can see the table and when we switch from the two-stage fit restrained electrostatic potential to the single stage electric field electric fit electric field fitting the relative salvation energy differences was smaller than the one from one The one the red one uh from the single stage RASP fitting Also RASP two means the two-stage restrained electrostatic potential fit. I just um, yeah Oh RASP two. Oh, sorry. RASP two means the two-stage fit to electric field So from this calculation, we could um We could check that this could read to the read to Simplification in how we derive RASP charges And also we conducted some relative mm energies between conformers and compare them with different charge models And if you see the left figures, that's the that's the charge to charge comparison and the Black lines is the charges from the two-stage RASP fitting and the red dots are from the single stage RASP fitting Fitting method and the blue dots are from the single stage electric field fitting And you can see the you can these are from the 26 type peptides with five conformers itch and from these figures we could see the Electric field fitting lower the non-polar charges with Only with single-stage fitting compared to the single-stage fit Uh a single-stage electrostatic potential fit So these are so for these are our preliminary results for the calculations and We are going to do we are going to like focus on the another Larger sets of molecules to see to compare the performance between charges Thanks so I'll Now continue again With uh, or I told you before to talked about respite that I want to talk about That RASP is a bit more computational demanding in how we can come up with a solution for that and actually Uh a few years ago also hard to fuck was quite expensive and then chris bailey came up with The great idea of am1 bcc's And we want to do like a new version of am1 bcc's to reduce the costs of RASP too And to do that. I just want to briefly start with Explaining what are am1 bcc or what is am1 bcc and the idea is Uh Thanks, so yeah, so I'd be interested to know The how much of the cost difference between using am1 bcc and using Uh quantum approaches what what about cost difference actually is and how much that really impacts people because I know in the amber community A lot of people use am1 bcc purely because it's automated Not because it's cheaper, right? So you can use anti chamber You can tell it use am1 bcc. It's a single command line. It just runs That if you want to use qm you have to put in the qm package and it's a manual process so people don't do it, right? If that I feel like if that was automated The the cost of a qm part wouldn't I mean compared to the cost of actually doing the simulations later It still seems trivial I mean that's True now, but probably I don't know am1 bcc was developed 15 years ago and 15 years ago. It was still an issue Yeah, that's when then when it was I mean no the reasons it was done my sense of the time I'm just wondering are all those reasons still valid today Uh, probably not for hard to fox 631 g but if you're going to a higher level qm calculation Which is more expensive than probably those uh effects Player role. Yeah, I just think we get to quantify that that's what Within a set I think this is true though only if you're talking about simulations not If we're talking about something like conformational searches where a small molecule needs to be minimized in seconds Yeah, I was going to say something along the same lines like you know, if you if you're Main application is free energy calculations or something and sure you can just do a full qm, but maybe Maybe you're going to be doing docking or conformational searches And so you need something that's a fraction of a second per molecule So so we are going to need a fast method Anyway, even though for free energy calculations, you might not So I guess very quickly. So so we've always had the option of doing qm or am1 bcc, right? uh Is it is now not the time to standardize completely and then fit all of the parameters based on that single choice of Electro in that there will be error cancellation. I'd rather than we have kind of a fast method and a Yeah, I'm wondering I think that's a really interesting question I'm wondering if that would be I think we'll discuss electrostatics this afternoon and that's a broad enough question I'm wondering if we should raise it then okay Another possibility is to have we've always talked about using the same input data to train multiple like expense Different expenses of force fields So one could have a very expensive force field where you have to do some sort of high quality REST like electrostatics with high quality qm That would be maybe more useful for individual free energy calculations And then there could also be a very fast version for a1 bcc Trained with the same data and maybe with slight accuracy differences So there could be two different forms of the same force field that are supposed to be consistent Hi, this is thomas just a question Um, we talk about rest versus pcc. Is there any Any good study that shows us a priority of one of them over the other or it's just you get numbers with with a rest Rest charge to get another number with pcc charges But at the end does it really make a difference? I think I think I can partly answer that and I think the answer is that Because of um a wizard who's in this room M1 bcc is in some cases slightly better than rest and at least not usually much worse Um, but that's because of the wizardry in part So yeah, I was gonna say keep in mind know that the dihedral parameters Are a function of the charges Yeah, right So if if the dihedral parameters will fit to a1 bcc and then you use a high level qm method without refitting the dihedral's you you may not get Yeah, I think the am in amber You know family we sort of treated partial charges as swappable Which we are not allowed to do but we do it anyway And so this is part of our opportunity to kind of fix that in the way that john just described and the Sorry, yeah, I was gonna Mike you want to go really quick. We should probably get back to the talk and and say some of this for the breakout but yeah, uh well, I guess I don't myself or whatever, but I guess another thing is that the bcc's give you a set of handles to tune against empirical data That the rest doesn't immediately do aside from the delta that michael the suture I will tell you how we want to deal with or where we see advantages using bcc's Uh for the parameterization um yeah Uh, something else. I I just looked a lot of different charges in the last two months and it seems that Sometimes the rest charges are a bit odd and bcc's are a bit more consistent. So that's just my impression um But I will Follow up. So what the basic idea of a m1 bcc is so you use a very cheap method like a m1 In that case and then you get an initial set of rough charges And then you somehow correct those charges with one charge corrections And one charge corrections are nothing else than you allow to transfer electrons from one atom forming a bond to the other atom So i'm reducing The charge on one atom forming a bond and adding the same amount of charge on the other atom And uh, how much electrons i'm moving is stored in the force field file Which where we have that bcc's parameterized Ah, sorry. Yeah, this um, this kind of just a this is a pretty dumb question Can't give a can't give an order of magnitude estimate on how large the bcc's are are they around 0.1 or 0.3 or as large as 1 Probably chris bailey knows that by heart So um, and the answer is some are big and some are small but they're uh, generally I would say on the order of if you looked at a distribution. I'd say a lot of them are Uh, less than 0.1 But then there are a few really important ones that are probably between 0.1 and 0.2 and I'm trying to remember there could be No, I'd say most of them fall within that range Hi, actually i'm online can someone repeat that question Leaping ask uh, what the magnitude of the bcc's was And chris bailey answered that they were mainly in the range of 0.1 with a few outliers where the bcc's were higher To sum up the discussion So, uh, why do we now think that we can train a new version of bcc's? So the first thing is this am1 bcc was trained To reproduce resp Calculations sort of hardware fucker calculations. So we are now using new potentials for resp. So we can also train new bcc models on that new esp potentials The other thing is that we think that the smirner of representation of force fields give us a great framework to define the bcc's and That we can more easily see where we can combine bcc's and where we have to add bcc's The other thing is that we can improve the training set that so the qm training set and we want to do that by looking at Databases and pull out the frequently occurring substituents and there are also from industries the Request to include a lot of hetero aromatics and we know also where am1 bcc's Has its deficiencies. So we would like to include boron data and phosphorous and sulfurs more extensively and Chris Bailey did it to a small extent for the am1 bcc to include experimental reference data And we want to include much more experimental data And that should be one of the main advantages of the new generation of am1 bcc's And you see I put the am1 In parentheses because actually we want to probably replace am1 with another method just to reduce the magnitude of the bcc's so if we are already with our fast method quite close to the Target charges we can really use the bcc's to fine tune our charge model and that would be an advantage in our opinion and Uh, what goes in hand with that is the charge population methods So we would probably one don't want to use molecule charges anymore So probably want to reduce that with something else like for example hirschfold charges As those calculations are quite cheap. We can probably test a lot of different combinations here So yes We have another question Yeah, you mentioned a different qm method. Which ones did you have in mind? um so Of course, there's am1 probably with a lovedy and charges that are ever Anything else others mm prilicon methods or even a cheap hard to fuck calculation would be something we could do And probably we also want to do a resp fit for that cheap calculations. We will see so I'll just mention one of our infrastructure issues is that there's not really good Fully open easily kind of installable am1 implementations around So if we were able to find something that was in sci for for example Yeah, it was a cheap hard to fuck that could be from an infrastructure perspective much better Because then we could easily build on it and it would be fast So as I said, we will test something and probably if we really think that something which is not already implemented in sci for Probably we can get it implemented if we really think it's much better than anything else Otherwise, we will stick with the things we have implemented in the open source packages so Actually for to parameterize the molecule not a lot of things should change So you still run a low-level calculation, which is Probably not am1 anymore Any other cheap Qm method and then you look up again the bcc's and our first field file And this is a new set of bcc so with different values And then you should end up with quite decent force load charges And how do we want to train those bcc's? So the starting point is fit them to the electrostatic potential Of resp 2 so taking a fraction of the gas phase esp and the fraction of the implicit solvent esp and then solve that equation so that we reduce The error in the esp reproduction from the charges And then we should get an initial set of bcc's and then we want to use those bcc's and Include a lot of experimental training data. So fit those bcc's or optimize those bcc's together with lennard jones parameters and There the question from john codara from previous comes in play So here we think we can tune individual groups So we can tune individual bcc's and that should give us more flexibility for Tuning to experimental training data and again that should be implemented with property with the property calculation calculator and force balance So to also give you an estimated timeline for this part of our project So we somehow rely on the the rescue potential. So we have to Define them first and then we want to parameterize different flavors of am1 bcc's so am1 It's probably replaced with some other qm method and then we see which gives Which method is closest to our target value and probably we Start with those values or death that method and then we are going and fitting only the mixing parameter first and the lennard jones parameters to get us in a Better starting position for then a fine tuning of bcc's together with lennard jones parameters And we hope that we then implement that in smirnoff within the end of next year That would be the plan and then as more and more experimental data comes available We should be able to improve the bcc's further and further so Something else. I want to talk about fixed charges Including off-center charges. So there was or there is the idea that we can improve much when we include off-center charges and most of the infrastructure is already in place and we are doing are including off-center charges that Should go parallel to the other effort and only we have to implement implement that it's still in the toolkit But then it should be This part of the project also starting That was all I want to talk about fixed charge models and I want to now talk a bit about polarizable charge models and to start off I want to Just briefly give an idea why we need explicit polarizations are One reason so there's probably a few other reasons is that when we have a molecule And once it's in water once it's a membrane and even when it's going to a protein binding bucket It's getting even more complicated one set of charges. It's probably not representative for All electrostatic potentials of those molecules. So probably we can improve the description when we include polarization explicitly The problem of polarizable force fields is that they are more expensive and We are trying to making a polarizable force field are and also a parameter system scheme Which is as cheap as possible. So that's our idea and that's what I want to explain to you in the next Few minutes. Yes You you may cover it in the next few slides, but one thing that always comes up is The hardware these days is very different to what it used to look like 20 years ago, right? So so everything has to be naturally parallel So something can be Mathematically very expensive, but computationally very cheap Or it can be mathematically cheap, but computationally very expensive or Uh, they don't map directly to each other. So I would say When we you know one of the things in the vote right now is it's not it the old days used to be able to say Well, what's the best polarization model to use within how much computing I can do And that's what you choose. But these days we need to pick polarization scheme that is inherently parallel Which means no kind of scf cycles if not if possible, right or No, essentially avoid any kind of iteration inside and in a loop I'm just wondering if you're keeping that in mind when we yeah I will talk about it in a minute what we have in in mind for that Probably you're the expert in really implementing it and you can yeah, I wish I had time to implement stuff Sorry, I just have a quick question to the fixed charge model. You said You plan on finishing like end of next year. So 2020 You don't mean end of this year for so the idea is we have REST ready until the end of this year and am 1 bcc is probably implemented within the 2020 release of smirnoff Whenever that release will be so that would be the idea or the plan up to now So again, I want to talk about do two different schemes to get parameters for the electrostatics one is again based on rest and Extend the rest with polarization and the other one is an extension of am 1 bcc with polarization um Again, want to start with what should it look like or how do we want to parameterize one molecule? So To parameterize one molecule, we would have to run one gas phase qm calculation uh, and then we solving something like a rest like equation and The equation is very similar to the rest equation. It the only difference is the last term so we are including the intramolecular polarization already and We to do that We need force field parameters for the polarize abilities and our idea is that we pre parameterize those polarize abilities And two things I want to mention here is first all the data I show you in the following slides are based on mp2 calculations But probably we switch to the same method as we are using for rest too And the other thing I want to mention how is we or how do we obtain those polarizability parameters? There are different ways one would be use experimental values and the other way is qm derived Polarize abilities and we are focusing on a qm derived method So just to give you uh, an idea what I talk about in the next few slides So we are fitting polarizabilities based on mp2 calculations with a triple cedar basis set We used a relatively small test set of 22 molecules which could cover a subset of the chemical space and uh for every molecule we used multiple confirmations and uh for every confirmation Then we used seven single point calculations where we have one have the the gas phase calculation and six polarized version of the same confirmation And then we propose a sequential fit of polarizabilities and charges So we do in the first step only fit uh fit of the polarizabilities So we fit the polarizabilities to the difference of a polarized and unpolarized version of the same confirmation and the same molecule And then we should add enough just with the effect of the external electric field and uh, we then assign polarizability parameters or Based on smurfs patterns and we find out that we actually only need very few parameters So if we do uh different Different creation between elements and hybridization states we already quite Good in agreement with the qm calculations and then in the second step We do the fit for the charges and the the charge fit already include the preset polarizabilities which were fitted in the first step And we only fit them to the unpolarized molecule with the rest like equation. I showed you One slide before or if we are fitting a bcc that's an am1 bcc like equation So and we did that sequential fit and we also did a Co-optimization of charges and polarizabilities and then we compared the results So if we are co-optimizing charges and polarizabilities Uh, we end up with an error of 4.11 which doesn't mean a lot. So we Uh, normalize that to 100 because that's the smallest error we can get with that representation And then we looked how good the sequential Approach does and actually you see our error our additional error is almost nothing So we only have 1 percent additional error in reproducing dsp Um, we did a similar thing for bcc So just to give you the idea what would be the idea or how you are parameterizing a molecule with a polarizable version of bcc You again using one low-level qm calculation and then you look up the bcc and the polarizabilities In a force field file and again you could fit the bcc and the polarizabilities at the same time or fit them individually and When you do that you see when you're going from individual charges to rough charges with bcc You introduce an additional error of 60 percent in dsp, so that's similar to the error between resp and am1 bcc, but it's just dsp so it doesn't really mean that it's worse for other properties and you see if we are going to sequential polarizabilities with bcc's the additional error is 10 percent So it's much less than introducing or using bcc's So it seems to be valid that we can use those sequential optimized polarizabilities for both charge models And the nice thing about that is that the polarizabilities are fitted before the charge model So they have never seen the charge model So this is exactly the same polarizability parameters for both charge models and we probably can come up with Any other charge model and stick the polarizabilities on top of that as long as we include the intramolecular polarization during the fitting process of the charges So you're looking at the error with respect to the electrostatic potential, right? But that condensed phase property errors might be totally different And we're missing things like nuclear quantum effects that that may actually be relevant here and might need to be incorporated in an average way So is there any real drawback to co-optimizing things like with the rest of the force field and liner jones? The problem is when we are co-optimizing Okay, if we're co-optimizing charges and Polarizabilities, we need a lot of molecules to get the polarizabilities, right? Otherwise you have somehow ill-defined polarizability, at least when you're only fitting the resb so you end up with Polarizability shifted for example from a carbon in the cx3 group to the hydrogen So you cannot end up with a negative polarizability on the carbon and You could have positivity constraints, but you're only using elemental typings, right? Yes, more or less Every individual atom of every individual molecule, right? Just you have like 12 classes of Polarizabilities and you're fitting to hundreds of different condensed phase properties as well Yeah, the problem is when we're fitting charges and liner jones and polarizabilities for charges we have a lot of parameters So you have a lot of molecules and a lot of parameters because every charge is individual If I may I think I think I'm feeling you're actually going to address john's point because I think john is really asking about fitting to him condensed phase data Oh, yeah, I will talk about that in a minute. Okay. Sorry So but it makes our life when we are fitting only to the sp much easier because in the first step we use a lot of molecules To get the polarizabilities, right? And then we can split it up for the charges and when for the charges We have a lot of parameters, but every molecule is separate And co-optimization here means co-optimization of charges and polarizabilities Just to the okay. Sorry. Yeah, so who gets the Yeah um So as I told you we are looking forward to making a polarizable force field as efficient as possible and That comes to the question from ross So there is the way of using a self-consistent field approach for calculating Diples at the idea of that is the charges are inducing an electric field and then you use that electric field To calculate dipoles and the dipoles again are inducing an electric field which are inducing dipoles And so you have to solve that in the self-consistent field approach Then there is the direct approach which we favor and in the direct approach You only consider the electric field from the point charges and neglect the dipole-dipole interactions and We think that At least the faster method it is less accurate But we want to know where we are less accurate and we That looked at the few molecules and we saw that we only do worse for molecules with aromatic rings And As we now know where this direct approach has is loss or deficiencies. We can probably correct them During the the fitting process yes Just go back one slider. So yes so To those two approaches In the limit that the direct approach you you essentially repeat multiple times Uh, like if you did multiple time steps, right? Does the direct approach Approximate the self-consistent approach and if Yeah, it's zero kelvin when you like it like i'm trying to think if you optimize Structure. No, you don't but you don't go to the same thing. No because uh, you probably Even if the structure is always the same You always end up with the same dipole which probably induce the exact same dipole on the other or and that You never account for that. Okay, and that and that and it can't cancel completely, right? No I think I understand what you're alluding to for example Let's say in in time step one everything gets polarized by just the charges now You have those polarizations in time step two you could then say well, everybody feels All those induced polarizations as well as a charge polarization. So so it sort of build up over time I You'd have a time lag. I mean, I think I think it would I think it would be It would be some non-equilibrium thing, right? Maybe it's interesting idea. I never thought about it But yeah, there's something we could think of a feeling it would mess up the distribution Okay, so so we had this problem with amoeba where the essentially If you didn't converge dipoles all the way or if you were using the last set of dipoles You would get this weird thing where the forces are no longer the gradient of the potential And that causes all sorts of stat net problems. So we don't we don't want that But we want something that it's a one-to-one mapping from the current configuration to a polarizability But there's an intermediate as well You've been following the work of bernie brooks for example where they from a single evaluation They can extrapolate essentially 99 of all of the polarization and still have it be completely deterministic So you don't need any iterative stuff at all So there there might be dial like a dial between the two that you can make it just a bit more expensive But still get almost all of that polarization um as At least our initial plan is following the direct approach because it is much as at least a bit faster than the self consistent field approach And we know where the problems are so we think we can correct them And it's not only that we can run them the simulations faster it's also using us or helping us for the parameterization as We can run the simulations faster. We can also Generate more converged results for the benchmarking and the parameterization of our results And the direct approach is also very much in the same spirit as a m1 bcc Where you give up a bit of the currency for speed up and that's a nice, uh, analogone um With that I just want to briefly mention what other things we are planning on doing So we are aware that we have to scale the polarizabilities from gas phase to condensed phase and we want to do that First with one parameter Two experimental data, so against experimental data and together with lennard jones parameters And then in the second step, we want to tune individual polarizabilities together with the lennard jones parameters But the first step is probably necessary to get us to our more suitable starting point and then as soon as we improve the fixed charge model like when we include offset The charges we also want to do that for the polarizable force fields And with that I also want to briefly mention the timeline for this project so again The first point is the same thing and then we want to parameterize this independent set of polarizabilities And then we start off with fitting One scaling parameters and the current lineage with the current lennard jones parameters Then we are doing the co-optimization of the scaling parameters with the Lennard jones parameters and then we Assume we are running in a lot of troubles and have to do a lot of testing And probably a realistic goal is that we implement that 2021 or that would be the ideal case And with that I just want to thank the the whole electrostatics group who is always discussing All those issues and the gilsson lab and the people involved in the gilsson lab in the open force field effort and the infrastructure group for implementing the offset the charges soon and All other members for the of the open force fields and the funding agencies for funding my stay here at ucst Yeah, so uh I just want to make a a quick comment that um With resp and a1 bcc having been around for such a long time I'm just so happy to see these models being actively worked on in the laboratories of Lee ping and mike gilsson, and I'm really super excited about this work having this direct polarization model, which is Having charge model independent polarizabilities is just I think super exciting and very much in the spirit of this of the open force field But just thinking about that talk one flavor at the end of it may may have been a bit of Where we're doing all the charge models, which maybe sort of is true But I think part of the spirit of this initiative is There's going to be a lot of science that allows us to explore Which is some of what michael talked about it and then some of that will end up becoming Different families of force fields down the line And so you're getting sort of a flavor of that where it might be headed And so some of that will be sort of the core fixed charge force field that A lot of us are are most immediately focused on but then there'll be these these other directions as well Any other last last questions before we go on? Yeah, I have a question online Yeah, hi, this is dory from Bristol Myers So thanks for that talk I just kind of wanted to make the point that I I know that Obviously, this is a really difficult problem Um, but I do Have some concerns when people talk about co-optimizing the charge parameters alongside the Leonard Jones Because I feel like one thing that would be great to be able to keep is to be able to physically Understand what different energetic contributions are coming from whether it's electrostatics Or you know dispersion Or polarization and when we fit the things together I'm afraid that that those distinctions Might not be accurate anymore I want to briefly comment on that. So that's the idea of Having only one parameter, which is fitting only the polarity of the the force fields So if we only have one parameter, we are not Messing completely the the whole physical sense of it so that's At least one reason I think it is Or if I agree that if we're fitting all the charges and all the Leonard Jones parameters We probably overfit the whole problem and end up with something which is not Has nothing to do with the physics anymore is that the reality be The radius part of Leonard Jones does affect the electrostatics massively But it determines how this process can get to each other One of the things that I've sort of come to tentatively is an issue with all this space that I Okay, so since she's not saying anything I think she was happy with the answer Okay Yeah, okay. Let me try that again. That's why I shouldn't say anything I was just saying is that one of the things that we've come to the basically the The radius part this the repulsive part of the Leonard Jones potential Does affect electrostatics because it determines how close charges can get to each other And one of the impressions I've developed hacking around with force fields using force balance in our lab This is work of sam canton and is that um the density is largely determined by the radii by the essentially sigma But then you're all those are also determining the electrostatics And if we want to be able to independently tune densities and cohesion energies We probably need to have a dial not only on sigma, but also on the level of polarity in some way Yeah, and I think also I would add that In a way They already are sort of partly co-optimized. I mean so okay, we didn't fit at least we didn't we didn't fit the electrostatics To condense phase properties in the amber family of force fields, but um the Leonard Jones have been optimized a lot based on things like densities and so on and And the bcc is in a m1 bcc one reason why it's sometimes better than resp is because those take into account some aspects of condensed phase properties to a limited extent And in in fact these initial optimizations in many cases were under constrained, which is what we found when we went back to look at dielectric constants that you can make Dielectric constants better by tweaking Leonard Jones and electrostatics parameters without impacting densities and heat of vaporization Because the original optimization was under constrained so so they are coupled and I don't yeah, and We can make force fields better by carefully Working through that and I don't think when we talk about co-optimizing them That means we want to throw out all the physics or anything. I think we want to Just do it carefully in a way that adds some additional data to constrain the optimization