Confusion regarding ramp variable

Hi everyone

I am using fix adept/fep to adept different values of the lambda parameter in the lj/cut/soft pair-style. To make the transition from 0 to 1 or 1 to 0, I am using ramp function to evaluate the new value of the parameter after every certain number of timesteps.

The problem is, the final value of the ramp function seems off by a little, depending on the number of steps in the ramp.

This is the section of my input file where I define relevant variables and the fix. Below that is the output file written at the end of the run.

variable nsteps equal 5
variable simlength equal 2000
variable totalsimlength equal ${simlength}*(${nsteps}+1)
variable delta_lambda equal 1.000000000/${nsteps}
variable lambda equal ramp(0.000000000,1.000000000+${delta_lambda})

fix change_lambda all adapt/fep ${simlength} pair lj/cut/soft lambda * * v_lambda after yes

run ${totalsimlength} 

write_data data_final.lmp
data_final.lmp: 

Pair Coeffs # lj/cut/soft                                                                    
1 0.352644 2.1595 1.00008                                                        
2 0.012785 4.8305 1.00008                                                                   
3 0.155402 3.166 1.00008                                                                     
4 0 1 1.00008 

According to the rule for calculation of ramp variable the value of the variable lambda in the final leg of the run should be 1.00000 and not 1.00008 as printed in the data file.

Can someone help me why is it not 1. Moreover, the value 1.00008 becomes 1.00002 if I increase the number of steps from 5 to 20.

This is because you are using the fix adapt/fep keyword “after yes” that means that the variable “lambda” is not evaluated on steps 2000, 4000, 6000, 8000, 10000 as you seem to expect but on steps 2001, 4001, 6001, 8001, 10001 and there the value of lambda is higher by one increment.

If you change nsteps, you have more steps in total, so the increment of the ramped variable per timestep is smaller. If you increase the number of steps by a factor of 4 the increment is a quarter.

Thanks Axel for pointing that out.

It would be great if it is clearly mentioned in the documentation. For now, it only says, adoption of the new value will be one step later. It doesn’t say clearly that evaluation is delayed too.

Thanks again
Chaitanya

Hi @akohlmey

Another way I realized this is using the ramp function in a similar manner to enable the charges in my system. In this case, while increasing the charges from 0 to 1, the final charges as output in the dumpfile are 1.00005, which I understand is because it is incremented by an extra timestep.
But what I do not understand is why then the initial adapted value is 0? Shouldn’t it be 0.00005 as we would get if the adaptation is at the first step? Similar is the case even for the earlier example but I couldn’t exactly print what the first adapted value was hence I am presenting the charge example.

The input file commands and the dump file output (for the first configuration printed) is shown below:

variable nsteps equal 20
variable simlength equal 1000
variable dumpfrequency equal 100
variable totalsimlength equal ${simlength}*(${nsteps}+1)
variable delta_lambda equal 1.0/${nsteps}
variable qCl equal -1.0*ramp(0,1+${delta_lambda})

fix scale_charges solute adapt/fep ${simlength} atom charge 2 v_qCl after yes

thermo 10000
thermo_style custom step etotal ke pe temp

dump TRAJ all custom ${dumpfrequency} dump.lammpstrj id mol type element q xu yu zu vx vy vz ix iy iz
dump_modify TRAJ element Na Cl O H sort id

run_style respa 2 8

run ${totalsimlength}                    #### Production run
34 34 2 Cl -0

Thanks for your patience.

Regards

This is the same problem and thus you have the same issues.

This is basically due to your “abusing” of fix adapt. If I would want to do this kind of change in rather large steps, I would just do a loop and then change the pair_coeff settings or charge values (via the set command), run for a chunk and then do the next loop iteration. That bypasses a lot of the rather confusing hackery that you are doing to make fix adapt to the stepwise change you want.

You should also take note that changing the interactions so massively so suddenly, will require significant re-equilibration time. This is why fix adapt is usually used with a continuous ramp. In fact, depending on the geometry and properties of your system, it may even be needed to use fix adapt to gradually move the interaction settings from the value for one step to the next to avoid significant disruptions and resulting disturbances (and sudden significant change of interactions is a bit like a small explosion).

Dear Axel,

Yes, this is the same problem because the confusion still persists.

I am a little confused by your use of word ‘abuse’. Fix adapt/fep documentation says that it is especially designed to use with compute fep or compute ti computes. That’s exactly what I am trying to do.

Also, by continuous ramp do you mean ramping the variable at each timestep?

Regards,
Chaitanya

I have nothing to add to my previous responses.

1 Like

To answer directly, your first dump frame is dumped on step 0, and on step 0 the value of ramp(0,$something) is zero.

Hi @srtee ,

Thanks for the response.

I went ahead and tested some things. In conclusion, I found that the operation of fix adapt/fep is non consistent.

When setting the ‘after’ as yes, it calculates the variable values and adapts it at one timestep later, EXCEPT for the first time. It adapts the value calculated at 0th timestep at the 0th timestep itself and not at timestep 1.

I do not understand the intention behind this set up, because this renders using ramp function with this fix almost useless, unless you make the frequency of adaptation very small, which unnecessarily adds computational expense.

Regards,
Chaitanya

Well, I don’t understand the intention either, but you can take that up with the feature’s original authors:

In any case, it is not too difficult to either quickly patch a custom ramp of your choosing:

variable offset equal #something, you can calculate it pretty easily
variable lambda equal (step==0)?0:ramp(${offset},1+${delta_lambda}-${offset})

or rewrite the C++ code to something of your own choosing and recompile as needed. That’s the beauty of open-source software! You can even contribute a bugfix if you’d like! There are many options – venting on an internet forum is the least productive one.

1 Like

Obviously, the author of this feature is thinking differently and there are several examples in the examples/PACKAGES/fep folder for validation. Have you considered that it may be you that is not using this feature correctly and making some assumptions that are not correct?

Perhaps, you should stop with your system and try to reproduce the data of the examples in “your way” to see, if you are getting consistent results.

… or even contact the author. He’s been quite helpful in the past with people that were asking (respectful) questions. Search results for 'Agilio Padua order:latest' - Materials Science Community Discourse

Hi @srtee ,

I did not intend to vent it. I understand there are ways to bypass the issues. I even took Axel’s suggestion to use loop which proved to be easier. I wanted to only know if there are some alternative benefits to make the first adaptation at timestep 0 instead of at 1 when ‘after’ keyword is yes. Also, since it is not mentioned in the documentation, I thought of putting it on this forum thinking the conversation around it may enlighten me and help others if they come this way.

However, I apologise if my response came across as rude or arrogant.

Regards,
Chaitanya

1 Like