Restrain the distance between atoms when refining xyzpars

[Nils Prinz]

Hello,

I am working on some refinements, that use a tetragonal Maghemite cell (P 43 3 2 symmetry). When refining the xyzpars it happens, that atoms are placed too close to each other. Is there a way to tell that the minimum distance between atoms needs to be for example >1.9A ?

If this is possible, it opens up many possibilities to refine disordered structures but stay in the physical boundaries :)

All the best,

Nils

[Simon Billinge]

There is in general infrastructure to introduce soft restraints (costs for exceeding some limit, e.g. too large of a displacement away from an average value) and also hard constraints, where regression updates are disallowed if they exceed some limit, so I believe that this should be possible.  Has anyone done it already and have a recipe?

--
You received this message because you are subscribed to the Google Groups "diffpy-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to diffpy-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/diffpy-users/ccbee070-c4d7-4d42-901c-2eebdbec82afn%40googlegroups.com.

--

Simon Billinge

Professor, Columbia University

Physicist, Brookhaven National Laboratory

[Mikkel Juelsholt]

Hi Nils 

There are a couple of way to do it, but only in Diffpy-CMI, not in PDFgui. Usually when I restrain atomic positions I do something like this: 

#Import the phase from the PDFgenerator

phase_X = pdfgenerator_X.phase

sgpars_X = phase_X.sgpars

#Global limits for the fractional coordinates.

pos_restrainO = 0.1

pos_restrainFe= 0.2

    

for par in sgpars_X.xyzpars:

    lclabel = par.par.obj.GetName().lower()

    lcsymbol = lclabel.rstrip(string.digits)

    # name the variables

    name="{}_{}_X".format(par.par.name, lclabel)

    tags = ['xyz', 'xyz_' + lclabel, 'xyz_' + lcsymbol]

    recipe.addVar(par, name=name, tags=tags)

    if name.count('x_o') == 1:

        recipe.restrain(name, lb = par.par.value - 1, ub = par.par.value + 1, sig = 0.0005)

    if name.count('y_o') == 1:

        recipe.restrain(name, lb = par.par.value - 1, ub = par.par.value + 1, sig = 0.0005) 

    if name.count('z_o') == 1:

        recipe.restrain(name, lb = par.par.value - 0.1, ub = par.par.value + 0.1, sig = 0.0005)

    if name.count('x_Fe) == 1:

        recipe.restrain(name, lb = par.par.value - pos_restrainM, ub = par.par.value + pos_restrainM, sig = 0.0005)

    if name.count('y_Fe') == 1:

        recipe.restrain(name, lb = par.par.value - pos_restrainM, ub = par.par.value + pos_restrainM, sig = 0.0005) 

    if name.count('z_Fe') == 1:

        recipe.restrain(name, lb = par.par.value - 0.1, ub = par.par.value + 0.1, sig = 0.0005)

This should restrain all the Fe coordinates and O to the set values. Remember they are fractional coordinates, so 0.1 will be a different value in Å along the a and c axis in your cell. I think this should do the trick. 

If you want to restrain the the distance between two atoms you could do something like this: 

phase_crystal1 = pdfgenerator_crystal1.phase

sgpars_crystal1 = phase_crystal1.sgpars

for par in sgpars_crystal1.xyzpars:

    lclabel = par.par.obj.GetName().lower()

    lcsymbol = lclabel.rstrip(string.digits)

    # name the variables

    name="{}_{}_bronze".format(par.par.name, lclabel)

    tags = ['xyz', 'xyz_' + lclabel, 'xyz_' + lcsymbol]

    recipe.addVar(par, name=name, tags=tags)

    

recipe.newVar('Dummy_variablex')

recipe.newVar('Dummy_variabley')

recipe.newVar('Dummy_variablez')

recipe.constrain('(xFe1**2-xO1**2)**(1/2)', 'Dummy variablex')

recipe.constrain('(yFe1**2-yO1**2)**(1/2)', 'Dummy variabley')

recipe.constrain('(zFe1**2-zO1**2)**(1/2)', 'Dummy variablez')

recipe.constrain('Dummy_variablex', lb = 2.0, ub = 3.0, sig = 0.0005)

recipe.constrain('Dummy_variabley', lb = 2.0, ub = 3.0, sig = 0.0005)

recipe.constrain('Dummy_variablez', lb = 2.0, ub = 3.0, sig = 0.0005)

I think this works, but I have not tried it out so there are likely a mistake somewhere. The idea is to create a variable which is always equal to the the x (y or z) distance between  the atoms of interest and the restrain the dummy variable to be between some values, in this case 2 and 3. I am not 100 % it makes sense to do it like this, this is just off top of my head. I am quite sure I mix up some geometry, but it should work if you get the right equations in the constrain and the right number in the restrain. 

Hope this helps!

Mikkel

[Simon Billinge]

thank you so much Mikkel, this is super helpful!

To view this discussion on the web visit https://groups.google.com/d/msgid/diffpy-users/24544c6b-33d3-4197-a203-8ab59935be40n%40googlegroups.com.

[Sandra Skjærvø]

This is very helpful for me too! Thanks, Mikkel!

[Mikkel Juelsholt]

I glad you can all use it. But I think I made a small mistake. 

Diffpy usually does not use capital letters for naming atomic positions, so the variable names I have used in the example should be fe and o, instead of Fe and O. But you should check the names of your atom variables in your own script. 

Also a handy trick to export the full refined structure can be done with this function: 

def export_cifs():

    name = name + '.cif'

    with open(name, 'w') as fp:

        crystal_structure.CIFOutput(fp)

    return

# uncomment to really export

export_cifs()

[Nils Prinz]

Hi Mikkel,

this is indeed very helpful. I will try it out and will let you know how it worked out.

Thank you!!

[Nils Prinz]

Hi,

the first script is running fine!

I tried to get this running:

phase_crystal1 = pdfgenerator_crystal1.phase

sgpars_crystal1 = phase_crystal1.sgpars

for par in sgpars_crystal1.xyzpars:

    lclabel = par.par.obj.GetName().lower()

    lcsymbol = lclabel.rstrip(string.digits)

    # name the variables

    name="{}_{}".format(par.par.name, lclabel)

    tags = ['xyz', 'xyz_' + lclabel, 'xyz_' + lcsymbol]

    recipe.addVar(par, name=name, tags=tags)

    

recipe.newVar('Dummy_variablex')

recipe.newVar('Dummy_variabley')

recipe.newVar('Dummy_variablez')

recipe.constrain('(xFe1**2-xO1**2)**(1/2)', 'Dummy variablex')

recipe.constrain('(yFe1**2-yO1**2)**(1/2)', 'Dummy variabley')

recipe.constrain('(zFe1**2-zO1**2)**(1/2)', 'Dummy variablez')

recipe.constrain('Dummy_variablex', lb = 2.0, ub = 3.0, sig = 0.0005)

recipe.constrain('Dummy_variabley', lb = 2.0, ub = 3.0, sig = 0.0005)

recipe.constrain('Dummy_variablez', lb = 2.0, ub = 3.0, sig = 0.0005)

Unfortunately python responds with ValueError: The parameter '(xFe1**2-xO1**2)**(1/2)' cannot be found. I already tried different names like x_Fe1, x_fe1 (I think x_fe1 should be the one that is coming from the lclabel name.

Maybe the names inside the string cannot be recognized as values?

All help is welcome!

Nils

[Mikkel Juelsholt]

Hi Nils

So you get the error (I just found out myself) because the constrain syntax wants an already defined parameter first and then the math that it needs to be equal to. Essentially you need to reverse the order of how we give it the parameter names. So the code should hopefully work if you do something like below. Also I noticed an error in my original post. The recipe.constrain('Dummy_variablex', lb = 2.0, ub = 3.0, sig = 0.0005) should be a restrain command. I fixed it below.

phase_crystal1 = pdfgenerator_crystal1.phase

sgpars_crystal1 = phase_crystal1.sgpars

for par in sgpars_crystal1.xyzpars:

    lclabel = par.par.obj.GetName().lower()

    lcsymbol = lclabel.rstrip(string.digits)

    # name the variables

    name="{}_{}".format(par.par.name, lclabel)

    tags = ['xyz', 'xyz_' + lclabel, 'xyz_' + lcsymbol]

    recipe.addVar(par, name=name, tags=tags)

    

recipe.newVar('Dummy_variablex')

recipe.newVar('Dummy_variabley')

recipe.newVar('Dummy_variablez')

recipe.constrain( 'Dummy variablex', '(xFe1**2-xO1**2)**(1/2)')

recipe.constrain('Dummy variabley', '(yFe1**2-yO1**2)**(1/2)',)

recipe.constrain('Dummy variablez', '(zFe1**2-zO1**2)**(1/2)')

recipe.restrain('Dummy_variablex', lb = 2.0, ub = 3.0, sig = 0.0005)

recipe.restrain('Dummy_variabley', lb = 2.0, ub = 3.0, sig = 0.0005)

recipe.restrain('Dummy_variablez', lb = 2.0, ub = 3.0, sig = 0.0005

Excited to hear if it works. 

Mikkel

[Nils Prinz]

Hi Mikkel,

sorry to reach to you so late. The code with the distance restraint works (python-wise). But I checked your formula (x_fe1**2-x_o1**2)**1/2 and this only yields the difference of the fractional coordinates, if I am not wrong. So we need to take into account the distance of points in a xyz-cartesian system, which is ((x_fe1-x_o1)**2+(y_fe1-y_o1)**2+(z_fe1-z_o1)**)**1/2. But we also need to consider the fractional nature of these values (so we need to take into account the a,b,c parameters, to get the real distances in Angstroem). If I did no mistake it should finally look like this:

((x_fe1*a-x_o1*a)**2+(y_fe1*b-y_o1*b)**2+(z_fe1*c-z_o1*c)**)**1/2

We could still run into problems, when x,y,z get negative, so probably this should be squared?

I now try to get this code running, but I am also pleased to hear what you think about this. I heard that Topas can also restrain distances. It would be of great value, if we can implement this in diffpy-cmi, as well. I will feedback any progress here.

Best,

Nils

[Mikkel Juelsholt]

Hi Nils 

No worries. I think everything looks okay. But you are correct I worked with the fractional coordinates above, so if you want a very specific distance in Å you need to take the unit cell into account. I do not see why you should run into problems when xyz gets negative. I think Diffpy should handle that no problem.

Topas can restrain distances as well. As far as I am aware there is only the rigid body macro in topas academic to restrain distances. I know this is written for organic molecules, but it should work for inorganic structures as well. But it is the same math. 

Less writing of math yourself since Topas uses more macro based input-file writing, but less flexibility because you are not in the Diffpy/python environment. But in your case I think both should do the trick. 

Cheers Mikkel