4075885785

Build Type

push

github

Committed by Shyue Ping Ong

Commit Message

Merge branch 'master' of github.com:materialsproject/pymatgen

Run Details

96 of 96 new or added lines in 27 files covered. (100.0%)

81013 of 102710 relevant lines covered (78.88%)

0.79 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

13.07

/pymatgen/command_line/vampire_caller.py

# Copyright (c) Pymatgen Development Team.
# Distributed under the terms of the MIT License.

"""
This module implements an interface to the VAMPIRE code for atomistic
simulations of magnetic materials.

This module depends on a compiled vampire executable available in the path.
Please download at https://vampire.york.ac.uk/download/ and
follow the instructions to compile the executable.

If you use this module, please cite the following:

"Atomistic spin model simulations of magnetic nanomaterials."
R. F. L. Evans, W. J. Fan, P. Chureemart, T. A. Ostler, M. O. A. Ellis
and R. W. Chantrell. J. Phys.: Condens. Matter 26, 103202 (2014)
"""

from __future__ import annotations

import logging
import subprocess
from shutil import which

import pandas as pd
from monty.dev import requires
from monty.json import MSONable

from pymatgen.analysis.magnetism.heisenberg import HeisenbergMapper

__author__ = "ncfrey"
__version__ = "0.1"
__maintainer__ = "Nathan C. Frey"
__email__ = "ncfrey@lbl.gov"
__status__ = "Development"
__date__ = "June 2019"

VAMPEXE = which("vampire-serial")


class VampireCaller:
    """
    Run Vampire on a material with magnetic ordering and exchange parameter information to compute the critical
    temperature with classical Monte Carlo.
    """

    @requires(
        VAMPEXE,
        "VampireCaller requires vampire-serial to be in the path."
        "Please follow the instructions at https://vampire.york.ac.uk/download/.",
    )
    def __init__(
        self,
        ordered_structures=None,
        energies=None,
        mc_box_size=4.0,
        equil_timesteps=2000,
        mc_timesteps=4000,
        save_inputs=False,
        hm=None,
        avg=True,
        user_input_settings=None,
    ):
        """
        user_input_settings is a dictionary that can contain:
        * start_t (int): Start MC sim at this temp, defaults to 0 K.
        * end_t (int): End MC sim at this temp, defaults to 1500 K.
        * temp_increment (int): Temp step size, defaults to 25 K.

        Args:
            ordered_structures (list): Structure objects with magmoms.
            energies (list): Energies of each relaxed magnetic structure.
            mc_box_size (float): x=y=z dimensions (nm) of MC simulation box
            equil_timesteps (int): number of MC steps for equilibrating
            mc_timesteps (int): number of MC steps for averaging
            save_inputs (bool): if True, save scratch dir of vampire input files
            hm (HeisenbergModel): object already fit to low energy
                magnetic orderings.
            avg (bool): If True, simply use <J> exchange parameter estimate.
                If False, attempt to use NN, NNN, etc. interactions.
            user_input_settings (dict): optional commands for VAMPIRE Monte Carlo

        Parameters:
            sgraph (StructureGraph): Ground state graph.
            unique_site_ids (dict): Maps each site to its unique identifier
            nn_interactions (dict): {i: j} pairs of NN interactions
                between unique sites.
            ex_params (dict): Exchange parameter values (meV/atom)
            mft_t (float): Mean field theory estimate of critical T
            mat_name (str): Formula unit label for input files
            mat_id_dict (dict): Maps sites to material id # for vampire
                indexing.

        TODO:
            * Create input files in a temp folder that gets cleaned up after run terminates
        """
        self.mc_box_size = mc_box_size
        self.equil_timesteps = equil_timesteps
        self.mc_timesteps = mc_timesteps
        self.save_inputs = save_inputs
        self.avg = avg

        if not user_input_settings:  # set to empty dict
            self.user_input_settings = {}
        else:
            self.user_input_settings = user_input_settings

        # Get exchange parameters and set instance variables
        if not hm:
            hmapper = HeisenbergMapper(ordered_structures, energies, cutoff=3.0, tol=0.02)

            hm = hmapper.get_heisenberg_model()

        # Attributes from HeisenbergModel
        self.hm = hm
        self.structure = hm.structures[0]  # ground state
        self.sgraph = hm.sgraphs[0]  # ground state graph
        self.unique_site_ids = hm.unique_site_ids
        self.nn_interactions = hm.nn_interactions
        self.dists = hm.dists
        self.tol = hm.tol
        self.ex_params = hm.ex_params
        self.javg = hm.javg

        # Full structure name before reducing to only magnetic ions
        self.mat_name = hm.formula

        # Switch to scratch dir which automatically cleans up vampire inputs files unless user specifies to save them
        # with ScratchDir('/scratch', copy_from_current_on_enter=self.save_inputs,
        #                 copy_to_current_on_exit=self.save_inputs) as temp_dir:
        #     os.chdir(temp_dir)

        # Create input files
        self._create_mat()
        self._create_input()
        self._create_ucf()

        # Call Vampire
        with subprocess.Popen(["vampire-serial"], stdout=subprocess.PIPE, stderr=subprocess.PIPE) as process:
            stdout, stderr = process.communicate()
            stdout = stdout.decode()

        if stderr:
            vanhelsing = stderr.decode()
            if len(vanhelsing) > 27:  # Suppress blank warning msg
                logging.warning(vanhelsing)

        if process.returncode != 0:
            raise RuntimeError(f"Vampire exited with return code {process.returncode}.")

        self._stdout = stdout
        self._stderr = stderr

        # Process output
        nmats = max(self.mat_id_dict.values())
        parsed_out, critical_temp = VampireCaller.parse_stdout("output", nmats)
        self.output = VampireOutput(parsed_out, nmats, critical_temp)

    def _create_mat(self):
        structure = self.structure
        mat_name = self.mat_name
        magmoms = structure.site_properties["magmom"]

        # Maps sites to material id for vampire inputs
        mat_id_dict = {}

        nmats = 0
        for key in self.unique_site_ids:
            spin_up, spin_down = False, False
            nmats += 1  # at least 1 mat for each unique site

            # Check which spin sublattices exist for this site id
            for site in key:
                m = magmoms[site]
                if m > 0:
                    spin_up = True
                if m < 0:
                    spin_down = True

            # Assign material id for each site
            for site in key:
                m = magmoms[site]
                if spin_up and not spin_down:
                    mat_id_dict[site] = nmats
                if spin_down and not spin_up:
                    mat_id_dict[site] = nmats
                if spin_up and spin_down:
                    # Check if spin up or down shows up first
                    m0 = magmoms[key[0]]
                    if m > 0 and m0 > 0:
                        mat_id_dict[site] = nmats
                    if m < 0 and m0 < 0:
                        mat_id_dict[site] = nmats
                    if m > 0 > m0:
                        mat_id_dict[site] = nmats + 1
                    if m < 0 < m0:
                        mat_id_dict[site] = nmats + 1

            # Increment index if two sublattices
            if spin_up and spin_down:
                nmats += 1

        mat_file = [f"material:num-materials={nmats}"]

        for key in self.unique_site_ids:
            i = self.unique_site_ids[key]  # unique site id

            for site in key:
                mat_id = mat_id_dict[site]

                # Only positive magmoms allowed
                m_magnitude = abs(magmoms[site])

                if magmoms[site] > 0:
                    spin = 1
                if magmoms[site] < 0:
                    spin = -1

                atom = structure[i].species.reduced_formula

                mat_file += [f"material[{mat_id}]:material-element={atom}"]
                mat_file += [
                    f"material[{mat_id}]:damping-constant=1.0",
                    f"material[{mat_id}]:uniaxial-anisotropy-constant=1.0e-24",  # xx - do we need this?
                    f"material[{mat_id}]:atomic-spin-moment={m_magnitude:.2f} !muB",
                    f"material[{mat_id}]:initial-spin-direction=0,0,{spin}",
                ]

        mat_file = "\n".join(mat_file)
        mat_file_name = mat_name + ".mat"

        self.mat_id_dict = mat_id_dict

        with open(mat_file_name, "w") as f:
            f.write(mat_file)

    def _create_input(self):
        structure = self.structure
        mcbs = self.mc_box_size
        equil_timesteps = self.equil_timesteps
        mc_timesteps = self.mc_timesteps
        mat_name = self.mat_name

        input_script = [f"material:unit-cell-file={mat_name}.ucf"]
        input_script += [f"material:file={mat_name}.mat"]

        # Specify periodic boundary conditions
        input_script += [
            "create:periodic-boundaries-x",
            "create:periodic-boundaries-y",
            "create:periodic-boundaries-z",
        ]

        # Unit cell size in Angstrom
        abc = structure.lattice.abc
        ucx, ucy, ucz = abc[0], abc[1], abc[2]

        input_script += [f"dimensions:unit-cell-size-x = {ucx:.10f} !A"]
        input_script += [f"dimensions:unit-cell-size-y = {ucy:.10f} !A"]
        input_script += [f"dimensions:unit-cell-size-z = {ucz:.10f} !A"]

        # System size in nm
        input_script += [
            f"dimensions:system-size-x = {mcbs:.1f} !nm",
            f"dimensions:system-size-y = {mcbs:.1f} !nm",
            f"dimensions:system-size-z = {mcbs:.1f} !nm",
        ]

        # Critical temperature Monte Carlo calculation
        input_script += [
            "sim:integrator = monte-carlo",
            "sim:program = curie-temperature",
        ]

        # Default Monte Carlo params
        input_script += [
            f"sim:equilibration-time-steps = {equil_timesteps}",
            f"sim:loop-time-steps = {mc_timesteps}",
            "sim:time-steps-increment = 1",
        ]

        # Set temperature range and step size of simulation
        if "start_t" in self.user_input_settings:
            start_t = self.user_input_settings["start_t"]
        else:
            start_t = 0

        if "end_t" in self.user_input_settings:
            end_t = self.user_input_settings["end_t"]
        else:
            end_t = 1500

        if "temp_increment" in self.user_input_settings:
            temp_increment = self.user_input_settings["temp_increment"]
        else:
            temp_increment = 25

        input_script += [
            f"sim:minimum-temperature = {start_t}",
            f"sim:maximum-temperature = {end_t}",
            f"sim:temperature-increment = {temp_increment}",
        ]

        # Output to save
        input_script += [
            "output:temperature",
            "output:mean-magnetisation-length",
            "output:material-mean-magnetisation-length",
            "output:mean-susceptibility",
        ]

        input_script = "\n".join(input_script)

        with open("input", "w") as f:
            f.write(input_script)

    def _create_ucf(self):
        structure = self.structure
        mat_name = self.mat_name

        abc = structure.lattice.abc
        ucx, ucy, ucz = abc[0], abc[1], abc[2]

        ucf = ["# Unit cell size:"]
        ucf += [f"{ucx:.10f} {ucy:.10f} {ucz:.10f}"]

        ucf += ["# Unit cell lattice vectors:"]
        a1 = list(structure.lattice.matrix[0])
        ucf += [f"{a1[0]:.10f} {a1[1]:.10f} {a1[2]:.10f}"]
        a2 = list(structure.lattice.matrix[1])
        ucf += [f"{a2[0]:.10f} {a2[1]:.10f} {a2[2]:.10f}"]
        a3 = list(structure.lattice.matrix[2])
        ucf += [f"{a3[0]:.10f} {a3[1]:.10f} {a3[2]:.10f}"]

        nmats = max(self.mat_id_dict.values())

        ucf += ["# Atoms num_materials; id cx cy cz mat cat hcat"]
        ucf += [f"{len(structure)} {nmats}"]

        # Fractional coordinates of atoms
        for site, r in enumerate(structure.frac_coords):
            # Back to 0 indexing for some reason...
            mat_id = self.mat_id_dict[site] - 1
            ucf += [f"{site} {r[0]:.10f} {r[1]:.10f} {r[2]:.10f} {mat_id} 0 0"]

        # J_ij exchange interaction matrix
        sgraph = self.sgraph
        ninter = 0
        for idx in range(len(sgraph.graph.nodes)):
            ninter += sgraph.get_coordination_of_site(idx)

        ucf += ["# Interactions"]
        ucf += [f"{ninter} isotropic"]

        iid = 0  # counts number of interaction
        for idx in range(len(sgraph.graph.nodes)):
            connections = sgraph.get_connected_sites(idx)
            for c in connections:
                jimage = c[1]  # relative integer coordinates of atom j
                dx = jimage[0]
                dy = jimage[1]
                dz = jimage[2]
                j = c[2]  # index of neighbor
                dist = round(c[-1], 2)

                # Look up J_ij between the sites
                if self.avg is True:  # Just use <J> estimate
                    j_exc = self.hm.javg
                else:
                    j_exc = self.hm._get_j_exc(idx, j, dist)

                # Convert J_ij from meV to Joules
                j_exc *= 1.6021766e-22

                j_exc = str(j_exc)  # otherwise this rounds to 0

                ucf += [f"{iid} {idx} {j} {dx} {dy} {dz} {j_exc}"]
                iid += 1

        ucf = "\n".join(ucf)
        ucf_file_name = mat_name + ".ucf"

        with open(ucf_file_name, "w") as f:
            f.write(ucf)

    @staticmethod
    def parse_stdout(vamp_stdout, nmats):
        """Parse stdout from Vampire.

        Args:
            vamp_stdout (txt file): Vampire 'output' file.
            nmats (int): Num of materials in Vampire sim.

        Returns:
            parsed_out (DataFrame): MSONable vampire output.
            critical_temp (float): Calculated critical temp.
        """
        names = ["T", "m_total"] + ["m_" + str(i) for i in range(1, nmats + 1)] + ["X_x", "X_y", "X_z", "X_m", "nan"]

        # Parsing vampire MC output
        df = pd.read_csv(vamp_stdout, sep="\t", skiprows=9, header=None, names=names)
        df.drop("nan", axis=1, inplace=True)

        parsed_out = df.to_json()

        # Max of susceptibility <-> critical temp
        critical_temp = df.iloc[df.X_m.idxmax()]["T"]

        return parsed_out, critical_temp


class VampireOutput(MSONable):
    """
    This class processes results from a Vampire Monte Carlo simulation
    and returns the critical temperature.
    """

    def __init__(self, parsed_out=None, nmats=None, critical_temp=None):
        """
        Args:
            parsed_out (json): json rep of parsed stdout DataFrame.
            nmats (int): Number of distinct materials (1 for each specie and up/down spin).
            critical_temp (float): Monte Carlo Tc result.
        """
        self.parsed_out = parsed_out
        self.nmats = nmats
        self.critical_temp = critical_temp

1	# Copyright (c) Pymatgen Development Team.
2	# Distributed under the terms of the MIT License.
3
4	"""	1✔
5	This module implements an interface to the VAMPIRE code for atomistic
6	simulations of magnetic materials.
7
8	This module depends on a compiled vampire executable available in the path.
9	Please download at https://vampire.york.ac.uk/download/ and
10	follow the instructions to compile the executable.
11
12	If you use this module, please cite the following:
13
14	"Atomistic spin model simulations of magnetic nanomaterials."
15	R. F. L. Evans, W. J. Fan, P. Chureemart, T. A. Ostler, M. O. A. Ellis
16	and R. W. Chantrell. J. Phys.: Condens. Matter 26, 103202 (2014)
17	"""
18
19	from __future__ import annotations	1✔
20
21	import logging	1✔
22	import subprocess	1✔
23	from shutil import which	1✔
24
25	import pandas as pd	1✔
26	from monty.dev import requires	1✔
27	from monty.json import MSONable	1✔
28
29	from pymatgen.analysis.magnetism.heisenberg import HeisenbergMapper	1✔
30
31	__author__ = "ncfrey"	1✔
32	__version__ = "0.1"	1✔
33	__maintainer__ = "Nathan C. Frey"	1✔
34	__email__ = "ncfrey@lbl.gov"	1✔
35	__status__ = "Development"	1✔
36	__date__ = "June 2019"	1✔
37
38	VAMPEXE = which("vampire-serial")	1✔
39
40
41	class VampireCaller:	1✔
42	"""
43	Run Vampire on a material with magnetic ordering and exchange parameter information to compute the critical
44	temperature with classical Monte Carlo.
45	"""
46
47	@requires(	1✔
48	VAMPEXE,
49	"VampireCaller requires vampire-serial to be in the path."
50	"Please follow the instructions at https://vampire.york.ac.uk/download/.",
51	)
52	def __init__(	1✔
53	self,
54	ordered_structures=None,
55	energies=None,
56	mc_box_size=4.0,
57	equil_timesteps=2000,
58	mc_timesteps=4000,
59	save_inputs=False,
60	hm=None,
61	avg=True,
62	user_input_settings=None,
63	):
64	"""
65	user_input_settings is a dictionary that can contain:
66	* start_t (int): Start MC sim at this temp, defaults to 0 K.
67	* end_t (int): End MC sim at this temp, defaults to 1500 K.
68	* temp_increment (int): Temp step size, defaults to 25 K.
69
70	Args:
71	ordered_structures (list): Structure objects with magmoms.
72	energies (list): Energies of each relaxed magnetic structure.
73	mc_box_size (float): x=y=z dimensions (nm) of MC simulation box
74	equil_timesteps (int): number of MC steps for equilibrating
75	mc_timesteps (int): number of MC steps for averaging
76	save_inputs (bool): if True, save scratch dir of vampire input files
77	hm (HeisenbergModel): object already fit to low energy
78	magnetic orderings.
79	avg (bool): If True, simply use <J> exchange parameter estimate.
80	If False, attempt to use NN, NNN, etc. interactions.
81	user_input_settings (dict): optional commands for VAMPIRE Monte Carlo
82
83	Parameters:
84	sgraph (StructureGraph): Ground state graph.
85	unique_site_ids (dict): Maps each site to its unique identifier
86	nn_interactions (dict): {i: j} pairs of NN interactions
87	between unique sites.
88	ex_params (dict): Exchange parameter values (meV/atom)
89	mft_t (float): Mean field theory estimate of critical T
90	mat_name (str): Formula unit label for input files
91	mat_id_dict (dict): Maps sites to material id # for vampire
92	indexing.
93
94	TODO:
95	* Create input files in a temp folder that gets cleaned up after run terminates
96	"""
97	self.mc_box_size = mc_box_size	×
98	self.equil_timesteps = equil_timesteps	×
99	self.mc_timesteps = mc_timesteps	×
100	self.save_inputs = save_inputs	×
101	self.avg = avg	×
102
103	if not user_input_settings: # set to empty dict	×
104	self.user_input_settings = {}	×
105	else:
106	self.user_input_settings = user_input_settings	×
107
108	# Get exchange parameters and set instance variables
109	if not hm:	×
110	hmapper = HeisenbergMapper(ordered_structures, energies, cutoff=3.0, tol=0.02)	×
111
112	hm = hmapper.get_heisenberg_model()	×
113
114	# Attributes from HeisenbergModel
115	self.hm = hm	×
116	self.structure = hm.structures[0] # ground state	×
117	self.sgraph = hm.sgraphs[0] # ground state graph	×
118	self.unique_site_ids = hm.unique_site_ids	×
119	self.nn_interactions = hm.nn_interactions	×
120	self.dists = hm.dists	×
121	self.tol = hm.tol	×
122	self.ex_params = hm.ex_params	×
123	self.javg = hm.javg	×
124
125	# Full structure name before reducing to only magnetic ions
126	self.mat_name = hm.formula	×
127
128	# Switch to scratch dir which automatically cleans up vampire inputs files unless user specifies to save them
129	# with ScratchDir('/scratch', copy_from_current_on_enter=self.save_inputs,
130	# copy_to_current_on_exit=self.save_inputs) as temp_dir:
131	# os.chdir(temp_dir)
132
133	# Create input files
134	self._create_mat()	×
135	self._create_input()	×
136	self._create_ucf()	×
137
138	# Call Vampire
139	with subprocess.Popen(["vampire-serial"], stdout=subprocess.PIPE, stderr=subprocess.PIPE) as process:	×
140	stdout, stderr = process.communicate()	×
141	stdout = stdout.decode()	×
142
143	if stderr:	×
144	vanhelsing = stderr.decode()	×
145	if len(vanhelsing) > 27: # Suppress blank warning msg	×
146	logging.warning(vanhelsing)	×
147
148	if process.returncode != 0:	×
149	raise RuntimeError(f"Vampire exited with return code {process.returncode}.")	×
150
151	self._stdout = stdout	×
152	self._stderr = stderr	×
153
154	# Process output
155	nmats = max(self.mat_id_dict.values())	×
156	parsed_out, critical_temp = VampireCaller.parse_stdout("output", nmats)	×
157	self.output = VampireOutput(parsed_out, nmats, critical_temp)	×
158
159	def _create_mat(self):	1✔
160	structure = self.structure	×
161	mat_name = self.mat_name	×
162	magmoms = structure.site_properties["magmom"]	×
163
164	# Maps sites to material id for vampire inputs
165	mat_id_dict = {}	×
166
167	nmats = 0	×
168	for key in self.unique_site_ids:	×
169	spin_up, spin_down = False, False	×
170	nmats += 1 # at least 1 mat for each unique site	×
171
172	# Check which spin sublattices exist for this site id
173	for site in key:	×
174	m = magmoms[site]	×
175	if m > 0:	×
176	spin_up = True	×
177	if m < 0:	×
178	spin_down = True	×
179
180	# Assign material id for each site
181	for site in key:	×
182	m = magmoms[site]	×
183	if spin_up and not spin_down:	×
184	mat_id_dict[site] = nmats	×
185	if spin_down and not spin_up:	×
186	mat_id_dict[site] = nmats	×
187	if spin_up and spin_down:	×
188	# Check if spin up or down shows up first
189	m0 = magmoms[key[0]]	×
190	if m > 0 and m0 > 0:	×
191	mat_id_dict[site] = nmats	×
192	if m < 0 and m0 < 0:	×
193	mat_id_dict[site] = nmats	×
194	if m > 0 > m0:	×
195	mat_id_dict[site] = nmats + 1	×
196	if m < 0 < m0:	×
197	mat_id_dict[site] = nmats + 1	×
198
199	# Increment index if two sublattices
200	if spin_up and spin_down:	×
201	nmats += 1	×
202
203	mat_file = [f"material:num-materials={nmats}"]	×
204
205	for key in self.unique_site_ids:	×
206	i = self.unique_site_ids[key] # unique site id	×
207
208	for site in key:	×
209	mat_id = mat_id_dict[site]	×
210
211	# Only positive magmoms allowed
212	m_magnitude = abs(magmoms[site])	×
213
214	if magmoms[site] > 0:	×
215	spin = 1	×
216	if magmoms[site] < 0:	×
217	spin = -1	×
218
219	atom = structure[i].species.reduced_formula	×
220
221	mat_file += [f"material[{mat_id}]:material-element={atom}"]	×
222	mat_file += [	×
223	f"material[{mat_id}]:damping-constant=1.0",
224	f"material[{mat_id}]:uniaxial-anisotropy-constant=1.0e-24", # xx - do we need this?
225	f"material[{mat_id}]:atomic-spin-moment={m_magnitude:.2f} !muB",
226	f"material[{mat_id}]:initial-spin-direction=0,0,{spin}",
227	]
228
229	mat_file = "\n".join(mat_file)	×
230	mat_file_name = mat_name + ".mat"	×
231
232	self.mat_id_dict = mat_id_dict	×
233
234	with open(mat_file_name, "w") as f:	×
235	f.write(mat_file)	×
236
237	def _create_input(self):	1✔
238	structure = self.structure	×
239	mcbs = self.mc_box_size	×
240	equil_timesteps = self.equil_timesteps	×
241	mc_timesteps = self.mc_timesteps	×
242	mat_name = self.mat_name	×
243
244	input_script = [f"material:unit-cell-file={mat_name}.ucf"]	×
245	input_script += [f"material:file={mat_name}.mat"]	×
246
247	# Specify periodic boundary conditions
248	input_script += [	×
249	"create:periodic-boundaries-x",
250	"create:periodic-boundaries-y",
251	"create:periodic-boundaries-z",
252	]
253
254	# Unit cell size in Angstrom
255	abc = structure.lattice.abc	×
256	ucx, ucy, ucz = abc[0], abc[1], abc[2]	×
257
258	input_script += [f"dimensions:unit-cell-size-x = {ucx:.10f} !A"]	×
259	input_script += [f"dimensions:unit-cell-size-y = {ucy:.10f} !A"]	×
260	input_script += [f"dimensions:unit-cell-size-z = {ucz:.10f} !A"]	×
261
262	# System size in nm
263	input_script += [	×
264	f"dimensions:system-size-x = {mcbs:.1f} !nm",
265	f"dimensions:system-size-y = {mcbs:.1f} !nm",
266	f"dimensions:system-size-z = {mcbs:.1f} !nm",
267	]
268
269	# Critical temperature Monte Carlo calculation
270	input_script += [	×
271	"sim:integrator = monte-carlo",
272	"sim:program = curie-temperature",
273	]
274
275	# Default Monte Carlo params
276	input_script += [	×
277	f"sim:equilibration-time-steps = {equil_timesteps}",
278	f"sim:loop-time-steps = {mc_timesteps}",
279	"sim:time-steps-increment = 1",
280	]
281
282	# Set temperature range and step size of simulation
283	if "start_t" in self.user_input_settings:	×
284	start_t = self.user_input_settings["start_t"]	×
285	else:
286	start_t = 0	×
287
288	if "end_t" in self.user_input_settings:	×
289	end_t = self.user_input_settings["end_t"]	×
290	else:
291	end_t = 1500	×
292
293	if "temp_increment" in self.user_input_settings:	×
294	temp_increment = self.user_input_settings["temp_increment"]	×
295	else:
296	temp_increment = 25	×
297
298	input_script += [	×
299	f"sim:minimum-temperature = {start_t}",
300	f"sim:maximum-temperature = {end_t}",
301	f"sim:temperature-increment = {temp_increment}",
302	]
303
304	# Output to save
305	input_script += [	×
306	"output:temperature",
307	"output:mean-magnetisation-length",
308	"output:material-mean-magnetisation-length",
309	"output:mean-susceptibility",
310	]
311
312	input_script = "\n".join(input_script)	×
313
314	with open("input", "w") as f:	×
315	f.write(input_script)	×
316
317	def _create_ucf(self):	1✔
318	structure = self.structure	×
319	mat_name = self.mat_name	×
320
321	abc = structure.lattice.abc	×
322	ucx, ucy, ucz = abc[0], abc[1], abc[2]	×
323
324	ucf = ["# Unit cell size:"]	×
325	ucf += [f"{ucx:.10f} {ucy:.10f} {ucz:.10f}"]	×
326
327	ucf += ["# Unit cell lattice vectors:"]	×
328	a1 = list(structure.lattice.matrix[0])	×
329	ucf += [f"{a1[0]:.10f} {a1[1]:.10f} {a1[2]:.10f}"]	×
330	a2 = list(structure.lattice.matrix[1])	×
331	ucf += [f"{a2[0]:.10f} {a2[1]:.10f} {a2[2]:.10f}"]	×
332	a3 = list(structure.lattice.matrix[2])	×
333	ucf += [f"{a3[0]:.10f} {a3[1]:.10f} {a3[2]:.10f}"]	×
334
335	nmats = max(self.mat_id_dict.values())	×
336
337	ucf += ["# Atoms num_materials; id cx cy cz mat cat hcat"]	×
338	ucf += [f"{len(structure)} {nmats}"]	×
339
340	# Fractional coordinates of atoms
341	for site, r in enumerate(structure.frac_coords):	×
342	# Back to 0 indexing for some reason...
343	mat_id = self.mat_id_dict[site] - 1	×
344	ucf += [f"{site} {r[0]:.10f} {r[1]:.10f} {r[2]:.10f} {mat_id} 0 0"]	×
345
346	# J_ij exchange interaction matrix
347	sgraph = self.sgraph	×
348	ninter = 0	×
349	for idx in range(len(sgraph.graph.nodes)):	×
350	ninter += sgraph.get_coordination_of_site(idx)	×
351
352	ucf += ["# Interactions"]	×
353	ucf += [f"{ninter} isotropic"]	×
354
355	iid = 0 # counts number of interaction	×
356	for idx in range(len(sgraph.graph.nodes)):	×
357	connections = sgraph.get_connected_sites(idx)	×
358	for c in connections:	×
359	jimage = c[1] # relative integer coordinates of atom j	×
360	dx = jimage[0]	×
361	dy = jimage[1]	×
362	dz = jimage[2]	×
363	j = c[2] # index of neighbor	×
364	dist = round(c[-1], 2)	×
365
366	# Look up J_ij between the sites
367	if self.avg is True: # Just use <J> estimate	×
368	j_exc = self.hm.javg	×
369	else:
370	j_exc = self.hm._get_j_exc(idx, j, dist)	×
371
372	# Convert J_ij from meV to Joules
373	j_exc *= 1.6021766e-22	×
374
375	j_exc = str(j_exc) # otherwise this rounds to 0	×
376
377	ucf += [f"{iid} {idx} {j} {dx} {dy} {dz} {j_exc}"]	×
378	iid += 1	×
379
380	ucf = "\n".join(ucf)	×
381	ucf_file_name = mat_name + ".ucf"	×
382
383	with open(ucf_file_name, "w") as f:	×
384	f.write(ucf)	×
385
386	@staticmethod	1✔
387	def parse_stdout(vamp_stdout, nmats):	1✔
388	"""Parse stdout from Vampire.
389
390	Args:
391	vamp_stdout (txt file): Vampire 'output' file.
392	nmats (int): Num of materials in Vampire sim.
393
394	Returns:
395	parsed_out (DataFrame): MSONable vampire output.
396	critical_temp (float): Calculated critical temp.
397	"""
398	names = ["T", "m_total"] + ["m_" + str(i) for i in range(1, nmats + 1)] + ["X_x", "X_y", "X_z", "X_m", "nan"]	×
399
400	# Parsing vampire MC output
401	df = pd.read_csv(vamp_stdout, sep="\t", skiprows=9, header=None, names=names)	×
402	df.drop("nan", axis=1, inplace=True)	×
403
404	parsed_out = df.to_json()	×
405
406	# Max of susceptibility <-> critical temp
407	critical_temp = df.iloc[df.X_m.idxmax()]["T"]	×
408
409	return parsed_out, critical_temp	×
410
411
412	class VampireOutput(MSONable):	1✔
413	"""
414	This class processes results from a Vampire Monte Carlo simulation
415	and returns the critical temperature.
416	"""
417
418	def __init__(self, parsed_out=None, nmats=None, critical_temp=None):	1✔
419	"""
420	Args:
421	parsed_out (json): json rep of parsed stdout DataFrame.
422	nmats (int): Number of distinct materials (1 for each specie and up/down spin).
423	critical_temp (float): Monte Carlo Tc result.
424	"""
425	self.parsed_out = parsed_out	×
426	self.nmats = nmats	×
427	self.critical_temp = critical_temp	×

materialsproject / pymatgen / 4075885785

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous