I had a similar problem for my borat glass systems. I needed to create and break bonds between my Borat atoms and my oxygen atoms according to the distance. I ended up writing a new unit for Lammps which checks the bond distances during runtime.
What I actually did was I started the program by defining an additioanl fake set of potential parameters which give always zero energy. In the input file every atom is given the maximum number of bonds I expected to occur during the simulation, with every bond using that fake set of potential parameters. During runtime Lammps now checks in a new routine every n steps the distances between my borate and oxygen atoms, and according to the coordination of the atoms this routine overrites all internal bond parameters (which atoms are belonging to a bond, what potential is to be used for that bond and so on).
This way of doing it by just changing the parameters of existing bonds during runtime seemed to me easier than actually creating new bonds.
It works quiet smooth for me.
P.S. In my special problem the potential I wanted to use is an angle depended potential. Depending on the coordination of the borate (3 or 4 oxygens nearby) it needs to use different parameters. So i mainly have to change angle lists. But I imagine that one could program something similar for all other type of bonds as well.
-------- Original-Nachricht --------