Okay, so welcome to the practical part of the tutorial. So here in the first part, we will go through the basic system setup with only quantum mechanics using a newly developed interface between Gromax and CP2K. So basically the idea of the interface is that most of the people know how to use Gromax. Yeah, most of you probably know how to use Gromax, but there is a very limited number of people who can actually perform simulations also using CP2K. So our idea of developing the following interface is to provide the people with knowledge of Gromax to also be able to use quantum mechanics, especially in periodic systems. So yeah, let's start the practical part. So first I will do a quick lecture recap just for most important parts from it. Then I will describe how the Gromax-CP2K interface works. So then we will do a setup, then we'll do exercise on setting up this for now pure QM calculation using that interface. And I will explain what in CP2K input is kind of important to know. And yeah, next you can do exercises again. So yeah, and we finish with some more advanced usage, how to use the Gromax, for example, pulling a code to generate to do umbrella sampling for QM system. Okay, so let's start with the lecture recap. It's really short, it will be really short. So as might you know, that Gromax is molecular mechanics program. So it uses a classical so-called molecular force field to perform simulation of the system. What means that molecular force field, classical force field is that your total energy of your system and that's the gradient, that's the movement of the system is described by a number of quasi-classical parameters. So they are split usually into bonded and non-bonded. So bonded interactions or parameters is represented by bonds, angles, torsions and several actually kind of torsions in reality, which have a form of usually harmonic oscillator function. And non-bonded parameters are usually two types of them. So first of all, slender-johns or wonder-wise parameters or repulsion in atoms. And the second one is columbic interaction between your atoms, which also have a partial charge. So columbic interaction, so-called. And this is how you do it in Gromax. So the difference with CP2K now is that instead of a number of these parameters, instead of large number of parameters, you now describe a system on a fully quantum level. So starting from the initial, from the wave function basically. And what you do here is that you first have some quantum system. So for each atom in that quantum system, you provide the basis set and you also choose a method, for example, PFT, yeah? And then you calculate the forces and energies on the basis of just to that part, to the parameters. So how in CP2K it's working, it's working quite differently from normal QM software like Gaussian, for example. So the CP2K uses quick step algorithm, which is kind of extension of the, which kind of extension it is. It uses the so-called mixed Gaussian play-wave basis set, what it is. So it was explained in the lecture, but in the short.  So you have plain waves in the box. Yeah, like a box, particle in the box wave function like as a basis. And what you want to do, so you basically construct your initial density matrix for the DFT using Gaussian basis set. And then you map that basis onto the real space, multi-grid. So grids and plain waves are basically the same determination. You have the same definition. Usually grids is a real space grids and for the reciprocal space, yeah, you can use plain waves, so basically it's the same thing. You map basically that on the grid in the real space. Then you do this, you pass this into the reciprocal space, you do fast Fourier transform. And then you minimize energy, why you need to do it in the reciprocal space because then you have the periodicity. And then you can obtain a real density matrix in the reciprocal space. Then you can solve for a constant equation in the reciprocal space, obtain new density matrix, check if you converge or not in a density. So if your density changes just like the previous iteration is small enough, I'm saying, okay, I'm converged. And then from that converged density, you calculate energy forces, any other properties you want. If not, then you go onto the second cycle, you can expand your density in Gaussian basis set, you construct again constant matrix, map it back to the real space and so on. And the cycle continues. So this was a short recap of what you should pick up from the test lecture. And now let's move to how it actually works when you go to the, so our newly developed Gromach CP2K interface. So what Gromach CP2K essentially is doing, it uses a machinery of Gromach's for classical MD and methods from Gromach's on classical MD to perform also simulation with quantum forces, which can be derived from the CP2K. So here's more general slides. So here you have all three terms. So we have part from the molecular mechanics prior to the quantum mechanics and coupling parameters. It will be covered in the next lecture. So now just forget about this. But what is important here that the interface could provide you a full KMM description of the region. Normally as you do, but using the same topologies, force fields, parameters, methods, which you can normally use in Gromach's. So these are features which interface provides you in general. So it take a normal Gromach's topology, it automatically converts it from classical MD to QMM. So it made some modifications to charges, bonds. It automatically set up link atoms. It will be covered in the next lecture also. So it also provide some set of validated parameters, QM parameters for starting your simulation. It's not 100% correct for all cases. So we need to sometimes to check use of cells. But in most cases, it will be more or less suitable for biological systems. So it's also compatible with most of the simulations available in Gromach's. Maybe for now with except of the free energy perturbation. Yes. And also fully compatible with any Gromach's tools, which is in Gromach's and it supports also it can be highly parallelizable in some terms like if you're using ensemble methods. Okay. So how you should set up your QMM question in general? So first you need to make your system like you normally do it in Gromach's. So you just do a normal classical MD system for MD simulation. So you can equilibrate it. You can maybe this classical MD, you can do whatever you want. And then you just switch several parameters in the MDP file. And the system will automatically go into the QMM mode. And then you can do this QMM calculations or pure QMM calculations in the Gromach's. There's some examples of what you can do now. Okay. So let's start with our practical. So please open the QMM tutorial page and open the eight from the episodes choose the practical Gromach CP2K part one. There is a guide I chose it. So, and please do the following thing. So please do the set it up to your own environment part. And yeah, and that's all for now. Yes. I think we will continue in five minutes. By the way, if someone is interested, you can just copy paste some comments into directly from the site from the episode into the terminal window. Yes. It also may be much faster than printing it each time with the hands. Okay. I see that most of the people managed to do so. So let's continue. So the first exercise would be setting up a very simple QM system in a very small box like shown on the right side of the slide. So it is a very small molecule called N-Meteala-Cetamide with only 12 atoms. There will be no M-M system for now in this practical. So it have a charge zero of that part and multiplicity is one of other fundamental properties being multiplicity to be one. We will use a DFTB functional for that exercise. So the box is really small. It's only one molecule. In reality, it's periodic box. So you have images of the molecule here and there. And yeah, let's do the starting with the first exercise. Do the steps one to five. And once you have done also raise your hand. Good. If you have any problems with setting up the simulation, especially doing this JMFs, PDB to JMFs, but I suppose most of you are already familiar with gromaxel, there shouldn't be any problems. Yes. Do you want to say anything about those of these parameters in the MDP file? Yes, I will. Yes, yes, yes. And if you don't specify them, there are defaults and things like that. Yes. So the idea is that, yes, on the next slide, I will show multiplicity. Yeah. What does QM multiplicity one mean? So basically your QM system have another fundamental property which is called SPIN, yeah? Which is not the case in classical system, but in quantum mechanical system, you have that property as a SPIN. And when it comes to the quantum mechanics, you always need to specify SPIN state of your system. One means a singlet state. Basically you are informing that you want to calculate your system in a singlet state. Singlet state means that all electrons is paired up, basically. Yes, there is another states like triplet, for example, states mostly. Okay. So now we probably can go through the, okay, let's wait another minute. Okay. Now, for the basic set, for now you don't need to specify anything. Functional PB. So PB is one of the DFT functionals, which is kind of standard for the, one of the standard DFT functionals, which is used a lot in the computations. Usually it can be different. For example, a sample of you probably know about the bitrelip functionals, yeah. I mean, tomorrow I think Holly will also show you how you can use that in CP2K, whereas it is a bit much slower than the PB function. Computational mean, but yeah. This is basically the method how you calculate your system. How to explore it more simply. It's just a method. Just a pre-parameterized method to use. Okay, spin restricted or unrestricted. So yes, for now the interface automatically set up if it's restricted or unrestricted. So on the basis of the multiplicity you put it in, you put in and if it's multiplicity one, it will be restricted. If it's more than one, it will be unrestricted. But I mean tomorrow I will show you how to use non-standard parameters. For example, you can use unrestricted, if you want to use unrestricted coefficients with spin multiplicity one. So with single state it's also possible, but I will show you tomorrow how you can do that in the interface. Okay, let's continue. So please see the MDP parameters. Please open the EM.MDP file and you will see the following picture. So this is a very simple setup, which we use just to minimize the energy. So we use the steepest descent optimization algorithm here. This is sub-tolerance parameters, 100 maximum steps. So there is output frequency, normal parameters, so grow box, like we want to coordinate each one step, we want to put into log file each step, we want to calculate energy each step, and we want to output energy each step. There is cut-offs, they are set it up to very small value because we have a very small box and we don't have any mechanical parts, so they're fairly not so important, these parameters. And what is important is that part of the parameters. This is how you activate the QMM calculation in CP2K or Engromus, a two-CP2K interface. So, yes, you can also use another integrators like CG, yes, right, correct. But yeah, steepest just, yeah, there is a specific reason of putting parameters to 0.2 because I mean, for pure QMM calculation, that cut-offs is not needed. It's needed for MMM calculation, which is not case here. If you put standard cut-offs, what will happen? The Gromach still complains that it cannot, that the system is too small, so it cannot use the number of processors you want to ask him. It cannot just decompose the system. So, I've just put a very little small parameters because it's very small system. It's like a box of one nanometer by one nanometer by one nanometer. And it is too small for the standard parameters. So, I've just reduced them. But they are not important because it's MMM parameters. But there is no MMM system. Okay, as if you do a standard parameters, you then need to adjust the number of a Bayer ranks of open and bitreds and so on, which is just a bit painful. Okay, so how you activate the QMM parameters, you just put QMM active to true, by default it is false. So, it activates a module. As that you need to specify which group of atoms will be your QM system. So, in that case, whole system is QM. So, we just use here system, the standard group in Gromach. Then there is three parameters defining which method we want to use. So, here we use the QMM method PBE. It's just been that Gromach still automatically set up the CP2K parameters for the PBE function with a double set up basis set. And also you need to put charge and multiplicity of your system as I said before. It's a pure quantum related parameters so you need to know beforehand what is them. So, next, I will show if you open your CP2K input file, if you do a plus NMA.info, then what you see, it will be a CP2K parameters which have been automatically generated by the Gromachs. So, the third global section, of course it's just a standard head-heading parameters. We want to go with energy and force because we integrate that in the Gromachs. The sixth section is the most important for several evals section. So, here is the several subsections. For that, we use a method QMM code. So, it means that it's quick step with some external charges is there present. So, next section is DFT section which contains the parameters of your DFT then QMM section contains parameters for your QMM set up. MMS section which defines your main charges and so on is there present and subsets system coordinates, basis sets and so on. Okay, I will shortly go through that sections, what you see. So, method we're using QMM. Here is a charge and multiplicity. So, I pass it just direct to CP2K. This is basis set file name and potential file name. This is just a basis set database of CP2K. Usually they are provided with CP2K installation that kind of files, but you also can use your own if you want. There is a set up of the multi grids which Emiliana talked about in his lecture. So, basically the CP2K is multi grids approach. So, you not use one grid but use several grids with the different cut-offs. The largest one here, for example, using five grids with the largest grid cut-off 450 reburts. And relative cut-off is, for now it's not used but it will be used for defining the projection of your Gaussian basis functions onto the grid. So, it basically decides on which grid your Gaussian will be projected. Yeah, and this aligns the grids. So, a CFGF restart, it means that CP2K first will try to search anyway function file existing in your directory for a start. Or if it's not presented, it will use just a standard atomic case. And this is a accuracy of your convergence in density matrix space. Okay. So, next section is set up of actual potential or the potential of a DFT functional. So, first parameter is XC. It's actually exchange correlation code. Yeah, so it's a DFT setup. So, here we use the DFT functional with cut-off, it's precision parameters which automatically set up by the interface. In most cases, they will be okay. It's quite high actually, even. So, the method is GPV for quick steps. There is quick step parameters. We use the method GPV, which Emiliano talked about. Again, accuracy parameter, extrapolation method, it means how your wave function will be reused from previous step. So, for example, if you calculate step-by-step, yeah? So, what CP2K will do, it will use a previously converged wave function for the previous step onto the next step to accelerate convergence, I can say, of the new, at the new step. And what is also in the next section is subsystem section, so what you need to look at. So, this subsection specifies a system which we're using. So, here it's a full system cell. Full system cell, in that case, is one by one by one nanometer because they're in angstrom, so it's 10 by 10 by 10 angstroms. And you have a fully periodic cell, like usually in gromax. This, yeah, actually this cell should be the same as a gromax cell. It's automatically generated like that. And next, go in the basis set setup for each atom kind in your system. So, for example, here is kind hydrogen, yeah? So, it's element hydrogen and we have basis set for that and we have, we set up the potential which will be also at the later stage, will be explained what it is. So, there is kind of carbon, there should be kind of nitrogen and for oxygen, which are set it up in the same way. And now you can do the, oops, sorry, if you do the point seven, part seven, step seven of this exercise, you will find something like that. So, this is the result of your energy minimization. It's basically potential energy with respect to the step, optimization step. You see that it emerges as quite fast, actually in 20 steps, it's already almost at the minimum and then it almost has no change. So, yeah, this is kind of the result of the energy minimization. So, next, what we will do instead of energy minimization, we will do molecular dynamics with NVT ensemble. So, please do the remaining steps of exercise one. So, for the sake of the recording, I'll just read out the question in the answer. So, we have that. So, the question was Matthew asked, for the interface, is it required that Gromax is compiled with double precision? And Dimitri, you said no, could also be single precision, but in general for QMM, you want to have full double precision in the coordinates? Yeah, so, I can explain briefly what I meant. So, precision does not matter for Gromax. You can have both single and double precision. CP2 cable anyway work in double precision as most of the quantum codes, if not all of them. So, yes, single precision generate less trajectories. Trajectors will be smaller, for example, but in QMM, it is not the limiting step. I mean, usually your trajectories will be not so big as classical MD, because you cannot have so many steps usually. And thus, yeah, both will work, but double precision coordinates will be better, in my opinion. Can anyone share the links pertaining to the integrated Gromax CP2 cable so we can install our local clusters? Arna, did you get the link into the document? I guess somewhere, or on the site? It was on the registration page, but I did not put it with the course materials. I can put it into the chat. Yes, I guess, yeah. I guess into the chat or on the course page later, we'll put it, yeah, because I remember that Arna put it at some point, which I don't remember where. Yes, it is available into versions with plant and without also plant, if you want to use plant. Yeah, but it is available, but I warned you that it is just, I can say alpha version. So it is not included any support from the Gromax directly. For now, once it will be integrated into the Gromax, it will go into better beta version probably. I hope before the July, it could happen. And then it will be, yes. Yeah, this version, which I use in the larger, is also includes plant support. But in general, it is 2021.1 Gromax with KMM. Also, I think this week or next week will be Gromax 2021.2 release, and I will also provide that version. Okay, so if someone finished with steps eight to 11, please raise your hand so we can continue. Yes, we will continue anyway in four minutes. So, okay, so anyway, I will continue in a moment. You will always have a time to finish your exercise if you want, as I've already said, it will be arch accounts will be active. So if you want, you can go, and also the web page will be online, so you can go through the tutorial once again. If something is not clear, probably I will continue because most of the people I see manage to do this. So results of this will be like that. So basically, we are trying to use NVT ensemble, but, and here, like your dynamics will look like if you download trajectory, and for example, render it in PyMol or VMD or hotel program you're using. And the temperature plot will look like that. So there is a really high oscillations. It's oscillating again around 300 gallons, but it's really high oscillations. Why it is so? Because of course, the system was too small. You don't have much degrees of freedom. I will, later, we'll see an example that in the bigger system in the full QMM system, it will be much smaller oscillations, but basically, this is the dynamics still works. Okay, so I think this is the end of the XSS one. So this is how you set up the really small QMM system with the interface. And I think if someone is already done, let's move to the exercise two. Exercise two will be a bit more involving and I really want you to kind of make attention, pay attention to what you're doing because here we want to use your collective kind of force to do something meaningful. And what we will do here, we will again use purely quantum system without a M-sub system. But what we want to, our objective will be to make a free energy profile of the steel beam isomerization reaction using the QM, using the same functional PV. So it's already bigger system as you see. And what is important that it could isomerize around this double bond here from trans to cis configuration. And we want to produce the energy profile, how we will do the steel beam roll assembling. I will explain later what we will do or how it works and what will be the result. But for now, I want you to start exercise two and do the steps one to seven. Yes, for the step two, you need to do your grow file which you need to pick up. It can be picked up either in the document which Arnold linked or I can share a bit in the video here for a moment, Arnold can cut up. So here is your username on the archer and here is a structure you need to pick up on the step two and here is the hydro angle which you need to pick up. I'm just trying to explain it. So basically, so first you need to type zero then type enter and then the option seven will also appear. You need to copy basically the zero group into seven. So just everything which is typed you just need to type literally as it is. So there is zero enter, then name seven QM atoms enter then Q enter. Yes, in reality, there is no need, but yeah, just, yeah. You just to show how you need to specify the QM atoms, how you can specify the QM atoms using the GMX make index option. Because of course for bigger system like proteins you want to just pick up not the whole system, but several atoms and this is how you should do that. Okay, let's continue this. You'll have a time to finish your tutorial later. I mean, the whole week will be yours, that's steps. So what we're gonna do and what is now simulating is so-called umbrella sampling. So what is umbrella sampling? It's a short introduction. So basically, if your system have some reactivity so you want to study the action and you have two states, a state A and state B and there is some barrier between them and you want to know the energy profile, free energy profile in reality of that system.  You need some method to do that and umbrella sampling is one of that methods. So basically how it works, you simulate your system in the number of points along the pathway. So you should know your reaction coordinate beforehand, but it's in most cases it is like that. And you simulate your system in number of states. So in how you do that, you basically put a harmonic potential in each of that points, which are here indicated by arrows, you put a harmonic potential which will attract your system or your reaction coordinate towards that point, basically put into that point and you simulate it with that umbrella with that harmonic potential for some time. And then you take a distribution in each of that individual frames, in distribution of your coordinate, reaction coordinate initial frames and from that and from knowing the force constant of your harmonic potential and distribution, you can convolve your free energy profile. So for that, for example, here is the distribution which frame, yeah. It should look in reality in a best case, it should look like Gaussian, but the most important thing is they need to overlap. You should see the overlap between the individual frames. That's what we need to check. It need to be sufficient, sufficient, yeah. How sufficient it is very arbitrarily sync in reality. I will show you later why. And the energy profile, you can integrate from the coordinate using the GMX-HUN tool. So you just basically need gather all cool x.xvg files and tpr files for each simulation in one directory. And then with the GMX-HUN tool, you can integrate the whole profile and that's what we're going to do now. But for further information, if you're interested in how to use umbrella sampling in Gromax, there is a tutorial on that here. So, yeah, probably the best way for now, how you can do the reaction, how you can study the reactions in the Gromax. Yeah, so it was worth to try and worth to learn how to use it. Let's move to the next slide. Yes. So how it's looking in the MDP parameters, you can open that less with less MDP file and go to the section called pool. And in pooling code, it's basically pooling code. It's basically the molecular dynamics with applied additional force. And in pooling code, you just say set up that you want to pool. We want to pool here one coordinate, dihedral angle, that dihedral angle in Gromax it consists of four coordinates. It's defined by four coordinates. So here we define four groups with the coordinates. So first and fourth group is the center of masses of that ring coordinate. And two and three, it's position of these two atoms. And then we want to pool using umbrella. So it's basically harmonic potential and pool geometry is dihedral and pool dimensions is all dimensions. It's four. Yes, yes, yes. So we pool in all dimensions. And this is how you define your dihedral angle. So basically the dihedral angle will be spawned between the two planes defined by pairs of vectors. So basically between pairs between vector one to two and two to three, it's one plane. And the second plane is two to three and three to four vectors spawned on that vectors. So basically this dihedral angle between these two kind of planes defined by these vectors. So here it's some initial value but in reality you need to set up your own here. And you need to scan over that coordinate. In reality, you need to launch like many simulations starting with a different pooling coordinates. So the pooling rate is zero because we just want to restrain it in place. This is K value, which also need for to convolve the profile in kilojoules per multiliter square nanometer. And I think it's nanometer here, it's radian. Oh no, yeah, it's nanometer, yeah. And this is, you informed Gromov that you want to output force and coordinates of that pooling, of that reaction coordinate each step. Okay, so here is your QMM parameters. So it's basically the same as it was in previous case but here you set it up your QM atoms group first time. And this will be important in the further simulations with real QMM systems. So be careful usually with that because with that group, with that kind of groups you can set up your QM system, especially in the bigger protein systems. And again, if you have used the same PV functional zero charge and singlet spin state. Yes, so beforehand I've done, so this how your coordinate looks like in dynamics. And beforehand I've done this simulation with just normal MD, classical MD setup, use number 14 force field. And I've simulated for one and a second for each frame. It's very important, in my opinion. And you will get the following can most profile. So basically, this is your trans structure, this is your C structure. It's the hydro angle value and this free energy. And you see that it is really high barrier. It's 160 kilojoules per mole is an organization barrier. Yeah, let's see how it will change if we, instead of MM force field, we'll use QM system. So please do the part nine, copy your output result. You can finish the exercise to part nine. And after the lunch, I will convolve the profile and show you the result. So just do up to the part nine and then go for the for lunch. Yeah, that's the idea. Please do not forget to copy also pull x.xg files because I see that some people forgetting to do so because there is 10 deeper files and only eight xg files, we need both of them. So what is x axis and y axis? Yeah, so x axis is a time step and the y axis is coordinate. If you use pull x, it will be coordinate. If you look into the pull f dot xg file, it will be a force. That's why it's in kilojoules per mole. It should be per mole per nanometer, I guess. But I am now looking at the pull x files and they're all time with respect to the position. So we are probably, we are looking into the wrong file. Yes, I suppose it's about this slide and about this group one, group two, group three, group four. How they defined, they defined in the index file as always in Gromax. So if you go to the, I will show you now. If you go to your directory with a steel bitmarked directory, there is index dot ndx file and in that file, I will open it. Yes, so that file will look like that. So oops, sorry, yes. Okay, so this is the file index. So, and here is a group defined. So here they are. I just provided them beforehand. Yeah, in reality, the group one, as you see the many of atoms, it just ring one. As a group four is the second ring, a center of masses of these things, of course, defines. And these two is just basically of one atom, of one carbon atom here. So they defined in that file, where you also need to define, for example, QM atoms. Yeah, but this is a general Gromax thing. Always when you see that kind of stuff, you usually need to look them up in the index file. And another question is what is the pulling dimension means, I guess? So that pull, yes, second question, yeah. Ah, what has to be considered to set the force constant and pull rate in this case? So the force constant should be pretty arbitrary, I think. The only thing should be that your force constant should be stronger now that you can overcome the barrier. It's need to be guessed by you, so it's up to you. If you see that your system cannot cross the barrier near the barrier, so you always stick to the sides, then you need to increase the force constant. Pulling rate zero is that because you don't want to really pull the system. So what pulling rate will does, it will basically each peak a second will decrease or increase this value. Depending on the rate by some number. And so it means that your system will start pulling. There is another option, how you can do the free energy barrier using, for example, the Erzien scheme equality. And for that, you will need the pulling rate, but for umbrella sampling, pulling rate should be zero. So this is defined by umbrella sampling. So I guess the second question, the third question, in the GDP file dimension y, y, y corresponds to the group connections for the hero. No, y, y, y here means that you are not restricting pulling in any Cartesian dimensions. So for example, you can pull in Gromax, you can pull atom only x dimension, yeah, not in y and z. And you mean that if your pulling force will be in Cartesian space, not in x direction, it will just project out the components in y and z and pull into x direction. And here we want, of course, to pull in all three dimensions. So your forces will be not restricted to one dimension. Forces which you apply to pull, will be not restricted to one dimension, they will be pulled in all three dimensions, fully freely. Yeah, it's all, usually it's the case what we need to use in numeral sampling. Yeah, if you just don't want to pull in one dimension in Cartesian space. For example, in membrane systems, you usually want to pull along the z direction than you set only z direction, yeah, to y, to z. So about the results of the first part. And now you can go to the directory where you copied your PULIX and TPR files. It is tutorial shared tutorial umbrella. And I've generated for you the following command. So I'm just pull here is a command as you can do. So it's basically doing ham. So this is that file which contains names of all TPR files. This is that file which contains names of all PULIX files. And here is some options. So you want to also build histogram, I will show you later why. You need skill of JOLS per mole, minimum minus 180, maximum 20, ham beans, it's how it samples, and starting from zero time. Okay, and there is a two files now. One is profile.xug file. The second is histo.xug file. So the profile.xug file looks like that. I can show you. So this is your profile. So something happened apparently on the middle, but I will show you what happened later. So you can do it yourself. So this is what it looks like. And yeah, in reality, what's happened here? And that's why you always need to be careful when you're doing any simulations and check always, do a sanity check of your simulations. So of course you know that it cannot be that any other barrier it just drops to zero. And I can show you what's happened. So in reality, if you check the umbrella sampling simulation tutorial, it will be said that this heavily under-sampled region. It means that you have basically in none of your simulations, none of the trajectories ever reached that point. It's basically zero sample region. And if you built a histo, from the histo.xug, you can build that kind of histograms. This is basically for each of your umbrellas, for each of your frames, you can build a distribution in that frame. And ideally it should be Gaussian. Also it's not the case here. So they're not looking like a normal Gaussian. But the important idea is that exactly at that point you have really zero samples. So you basically did not. And that's why Ham says that, okay, if you don't have anything, it's probably zero here. Yeah, but what is also important here is that if you check the barrier, it's already not so big as it was in the MMFORCE field. It's much, much lower. And that means that your MMFORCE field problem not well described in this isomerization also. But what will happen, how to solve that problem? So to solve that problem, you need to take not 100 femtosecond time steps, but 100 femtosecond, so 100 femtosecond, but you need a picosecond, so you need much, much more simulations to be done. And that's what's happened when you do 10,000 steps. Now you see that histograms already looking like Gaussians really. So they're really Gaussians, they're living Gaussians. Still, that region is not very well sampled, but that's what's also connected to the question of Oliver, that probably the... I would say that if I do it again, I will adjust my force constant near that point. So I will increase that region force constant so we can pull it much further away. And maybe put some additional frames in that region. So I will put several other frames with a bit different angles, the actual angle values. Okay, but now you see that the profile again changed and now it's much more converted, it's smooth. So here is, yeah, also the idea is that if you do this for the protein system, I can warn you that you need to do a lot, the simulation should be done in much, much further. Your simulation should be extended even further. Typically for the protein system we are doing in our practical work, we are doing no less than 30 picoseconds simulation for each window. And this is only after the several milliseconds of MMP equilibration. So this is not really so fast, yeah, but this is the best way you can obtain your reaction free energy to estimate your reaction kinetics and thermodynamics. Okay, so this on one slide, all three approaches. So what we use, so this black one is a MMP force field. The orange one is the one with a very low sampling and the blue one is came with a high sampling. So always should be, you always should check if you convert your profile if you're doing an umbrella sampling. So how it's usually done, you, for example, do 30 picoseconds, okay, let's do 10 picosecond simulation. Then you extend this by another like five picoseconds and check if your profile changing. So you only stop doing simulation after your profile stay on the place. So it's not changing when you increase your trajectory further. Okay.