Wish list for future versions of MULTEM

[thomasaarholt]

This is a draft. For when I remember them, I thought I'd write a little wish list for future versions of MULTEM.
Hopefully these are things that I can help with in the coming weeks/months:

• Loading a .cif (standardised crystal format) file instead of a text file, and specifying the Debye-Waller factor as a constant for all elements AND/OR writing a conversion plugin for cif to txt/.mat
• In MATLAB, having a progress bar indicating the progress and average time remaining. Sometimes one runs the code with no real idea of how long things are gonna take.
• Some documentation indicating exactly how the results are saved in output_multislice, and how best to save them.

[Quantumstud]

Is there a package already developed for conversion of the cif file into the text file that can be used for the MATLAB based simulation?

[thomasaarholt]

I can provide a script that will do this, but it is written in python.

[Quantumstud]

Thanks a lot! That will be really helpful.

[thomasaarholt]

I will upload it as soon as I'm home from work.

…

On Wed, 10 Oct 2018, 18:50 Tara Prasad Mishra, ***@***.***> wrote: Thanks a lot! That will be really helpful. — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub <#13 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/ACmGj1TZ4cqvSvFjuQ_hZg9NgDNNf7brks5ujiVPgaJpZM4O8ubB> .

[thomasaarholt]

Use at your own risk. It will convert cif and xyz files in whatever orientation they are in (c-axis faces the viewer after simulation).

def save_cell_multem_txt(cif_cell, rms3d=0.085):
    '''Convert cif file to the multem txt format
    
    Can also take the xyz format as input (and really anything that ASE can take)
    
    Requires the ASE Atomic Simulation Environment Python package
    
    I recommend using Vesta (windows software) for manipulating 
    the cif file to enlarge it in the dimensions one wants
    - use Edit -> Edit Data -> Unit Cell -> Transform.
    
    The rms3d value is related to the debye-waller factor. It can
    be given as a constant for all 
    atoms, or as a dict of {atomic_number:rms3d_value}.
    It is a function of temperature, and can either be looked up in tables
    or calculated by DFT
    
    ie: (these are not real values) for carbon and zinc
    rms3d = {
    8:0.1,
    30:0.09
    }
    
    Warning: There may be issues if the cif file has been rotated (by 90 degs). Check carefully.
    '''
    import numbers
    import os
    from pathlib import Path
    import numpy as np
    from ase.io import read

    from scipy import io as out

    filepath = Path(cif_cell)
    name = filepath.stem
    path = filepath.parent
    
    cif_cell = read(cif_cell)
    
    sample_x = cif_cell.cell[0,0]
    sample_y = cif_cell.cell[1,1]
    sample_z = cif_cell.cell[2,2]

    lx = [sample_x]
    ly = [sample_y]
    dz = [sample_z/2] # Preferably some fraction of the unit cell, or just a small distance
    
    header = [lx + ly + dz + 5*[0]]
    
    occ = [1]
    label = [0]
    charge = [0]
    if isinstance(rms3d, numbers.Number):
        rms3d_list = [rms3d for Z in cif_cell.numbers]
    else:
        rms3d_list = [rms3d[Z] for Z in cif_cell.numbers]
    
    data = [[n] + list(xyz) + [rms3d] + occ + label + charge for n, xyz, rms3d in zip(cif_cell.numbers, cif_cell.positions, rms3d_list)]
    
    total = np.array(header + data)
    
    np.savetxt(path / (name + ".txt"), total, fmt='%.8f', newline="\n", delimiter=" ")

[Quantumstud]

Thank you. This works wonderfully well to generate an .txt file, which can be used in the MULTEM GUI. How do we use this script to generate the crystal by layers which can be used to generate the crystals as shown in the crystalline_materials so that we can use it in the MULTEM Matlab version?

[thomasaarholt]

My apologies, when I read "text file" in your first post, my brain assumed that you meant for the gui version, despite your clarification at the end. I will upload the version for matlab as well in a bit.

In terms of building the crystal, I recommend doing that in either ASE or Vesta, and exporting it as a .xyz file and then converting it using my scripts.

I explained some of this in the docstring in the above function.

[thomasaarholt]

As I'm now remembering things (a little while since I wrote the scripts):

Multem (Matlab) creates a crystalline sample based on some unit cell. That wasn't enough for me, since I wanted to build some larger and more complicated structures (think defected cells, grain boundaries etc).

Therefore I bypass the loading mechanism in the matlab version, and load a matlab .m file, replacing the equivalent of the Au001Crystal().

[thomasaarholt]

If you just want to create the equivalent of the Au001Crystal(), but for your material, I suggest you just duplicate a file in crystalling_materials, and then edit it to produce your own crystal.

For the cases of cubic, tetragonal and orthorhombic crystals, this is not hard. I am not sure how one would approach it for a non-orthogonal crystal.

[thomasaarholt]

def save_cell_multem(cif_cell, rms3d=0.085):
    '''Convert cif file to the multem matlab format
    
    Requires a `load(cif_cell.mat)` code in matlab, and should overwrite any a b c, spec_atoms etc.
    
    Can also take the xyz format as input (and really anything that ASE can take)
    
    Requires the ASE Atomic Simulation Environment Python package
    
    I recommend using Vesta (windows software) for manipulating 
    the cif file to enlarge it in the dimensions one wants
    - use Edit -> Edit Data -> Unit Cell -> Transform.
    
    The rms3d value is related to the debye-waller factor. It can
    be given as a constant for all 
    atoms, or as a dict of {atomic_number:rms3d_value}.
    It is a function of temperature, and can either be looked up in tables
    or calculated by DFT
    
    ie: (these are not real values) for carbon and zinc
    rms3d = {
    8:0.1,
    30:0.09
    }
    
    Warning: There may be issues if the cif file has been rotated (by 90 degs). Check carefully.
    '''
    import numbers
    import os
    from pathlib import Path
    import numpy as np
    from scipy import io as out
    from ase.io import read

    filepath = Path(cif_cell)
    name = filepath.stem
    path = filepath.parent
    
    cif_cell = read(cif_cell)

    occ = [1]
    label = [0]
    charge = [0]
    
    if isinstance(rms3d, numbers.Number):
        rms3d_list = [rms3d for Z in cif_cell.numbers]
    else:
        rms3d_list = [rms3d[Z] for Z in cif_cell.numbers]
    
    data = [[Z] + list(xyz) + [rms3d_list[Z]] +
            occ + label + charge for Z, xyz in zip(cif_cell.numbers, cif_cell.positions)]

    sample_x = cif_cell.cell[0, 0]
    sample_y = cif_cell.cell[1, 1]
    sample_z = cif_cell.cell[2, 2]

    matlabdict = {
        "spec_atoms": data,
        "spec_lx": sample_x,
        "spec_ly": sample_y,
        "spec_lz": sample_z,
        "spec_dz": sample_z/4,  # Preferably some fraction of the unit cell, or just a small distance
        "a": sample_x,
        "b": sample_y,
        "c": sample_z,
    }
    out.savemat(os.path.join(path, name), matlabdict)

Below is an example multem matlab script that loads the above output correctly. It also supports having vacuum around the edges, if that should be of interest.
STEM_example_aarholt.zip

[Quantumstud]

Wonderful! Thank you so much. The python script and the stem example scripts were beautiful. Previously all my attempts had resulted in MATLAB crashing altogether. The only changes I had to make was material specific. I sincerely hope the above example and the python code could be included in the original MULTEM package. This would greatly help new learners.

[Quantumstud]

I didn't notice this before, whenever I am converting a .cif file to the corresponding .mat file the number of atoms that are generated in the spec_atoms variables is 112. The actual atoms were about 64. (both the .cif and .mat file obtained after the conversion are attached).
I tried for a lot of other structures and always found out that the always the maximum number of atoms that are generated is 112.
On the other hand this convertor works wonderfully well for a .xyz crystal file but it doesn't take into account the lattice paramters for a .xyz file.
GaN.zip

[thomasaarholt]

How odd. Let me take a look.

[thomasaarholt]

First off - thanks for good problem solving and giving good feedback! Makes it much easier for me!

Solution below - I'm leaving my own discussion up so you can see what I've tried.

Now, your cif file definitely contains 64 atoms (try opening it in notepad++ or some other text editor.

I believe there are some interesting things going on with the specific cif file or the format, since it's supposed to be a repeatable cell. The ASE software, which is what I'm using for reading in the crystal, sees your cif file as containing 112 atoms. It also raises some issues about how certain atom positions are equivalent.

from ase.io import read
cr = read("NewGaN.cif")
len(cr.numbers)
112

Now, if I immediately export your cif file from ase as a new cif file and also as a xyz file, we see what is wrong, although I'm not sure I can explain it.

If I import your cif in Vesta, and export as .xyz, then import that file in ASE, it reads as having 105 atoms!

from ase.io import read
cr = read("NewGaN.xyz")
len(cr.numbers)
105

The 105 number should be correct, since now the atoms that are exactly on the border count twice rather than just once. 64 + 15 + 13 +13 = 105. (The numbers I am adding are the numbers on the three cell faces)

I observed that there are some numbers that are either very small and are just above or below zero in your original cif file (i.e. -6.16298e-33). I tried changing these to 0, but it did not have any effect on importing the cif file again. I have previously encountered this as an issue - for cif files especially, I recommend manually changing any value that is either extremely close to 1 or 0, equal to 0.

Solution
Right, having just tried opening the file in Vesta and exporting back into cif and importing this new cif file into ASE, it now reads having 64 atoms. That is what you expected, but probably not what you want, since the edge atoms are now gone.

This leads me to conclude that whatever method you are generating the cif file with is using either an old format or is exporting it incorrectly. Or perhaps ASE is importing it incorrectly, but I think that's unlikely seeing as it's under active development.

For what you should do: Just stick with converting to xyz and reading that into my script. That should do it. It makes more sense to be using xyz files anyway, since our model in its current form shouldn't be repeated in 3D-space in multem, and hence you probably want it to appear "complete". Therefore we want to use the non-repeating "xyz" format rather than the repeating "cif" format.

[thomasaarholt]

So with that giant post out of the way - how are you generating your cif files?

[Quantumstud]

So with that giant post out of the way - how are you generating your cif files?

I am generating my .cif files using virtual nano lab.

[Quantumstud]

First off - thanks for good problem solving and giving good feedback! Makes it much easier for me!

Solution below - I'm leaving my own discussion up so you can see what I've tried.

Now, your cif file definitely contains 64 atoms (try opening it in notepad++ or some other text editor.

I believe there are some interesting things going on with the specific cif file or the format, since it's supposed to be a repeatable cell. The ASE software, which is what I'm using for reading in the crystal, sees your cif file as containing 112 atoms. It also raises some issues about how certain atom positions are equivalent.

from ase.io import read
cr = read("NewGaN.cif")
len(cr.numbers)
112

Now, if I immediately export your cif file from ase as a new cif file and also as a xyz file, we see what is wrong, although I'm not sure I can explain it.

If I import your cif in Vesta, and export as .xyz, then import that file in ASE, it reads as having 105 atoms!

from ase.io import read
cr = read("NewGaN.xyz")
len(cr.numbers)
105

The 105 number should be correct, since now the atoms that are exactly on the border count twice rather than just once. 64 + 15 + 13 +13 = 105. (The numbers I am adding are the numbers on the three cell faces)

I observed that there are some numbers that are either very small and are just above or below zero in your original cif file (i.e. -6.16298e-33). I tried changing these to 0, but it did not have any effect on importing the cif file again. I have previously encountered this as an issue - for cif files especially, I recommend manually changing any value that is either extremely close to 1 or 0, equal to 0.

Solution
Right, having just tried opening the file in Vesta and exporting back into cif and importing this new cif file into ASE, it now reads having 64 atoms. That is what you expected, but probably not what you want, since the edge atoms are now gone.

This leads me to conclude that whatever method you are generating the cif file with is using either an old format or is exporting it incorrectly. Or perhaps ASE is importing it incorrectly, but I think that's unlikely seeing as it's under active development.

For what you should do: Just stick with converting to xyz and reading that into my script. That should do it. It makes more sense to be using xyz files anyway, since our model in its current form shouldn't be repeated in 3D-space in multem, and hence you probably want it to appear "complete". Therefore we want to use the non-repeating "xyz" format rather than the repeating "cif" format.

Thank you for your prompt reply. I think the .cif created using the virtual nanolab has some kind of weird problem in Vesta. It creates some additional atoms which are read in the ase. As you had advised I exported the .cif back from the Vesta and read it again. This shows up no error.

Thank you once again, I am now able to run the simulations. It is quite wonderful.