From afe6dc088e2a9018ee4c94af02de2589b03675fe Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Fri, 19 Jan 2024 16:53:40 +0100 Subject: [PATCH 01/16] clean writing when file is provided --- calphy/input.py | 21 +++++++++++++-------- 1 file changed, 13 insertions(+), 8 deletions(-) diff --git a/calphy/input.py b/calphy/input.py index 6c7a29f..2c7d08f 100644 --- a/calphy/input.py +++ b/calphy/input.py @@ -329,6 +329,7 @@ def _validate_all(self) -> 'Input': self._element_dict[element]['mass'] = self.mass[count] self._element_dict[element]['count'] = 0 self._element_dict[element]['composition'] = 0.0 + self._element_dict[element]['atomic_number'] = mendeleev.element(element).atomic_number #generate temporary filename if needed write_structure_file = False @@ -336,7 +337,7 @@ def _validate_all(self) -> 'Input': if self.lattice == "": #fetch from dict if len(self.element) > 1: - raise ValueError("Cannot create lattice for more than one element") + raise ValueError("Cannot create lattice for more than one element, provide a lammps-data file explicitly") if self.element[0] in element_dict.keys(): self.lattice = element_dict[self.element[0]]['structure'] self.lattice_constant = element_dict[self.element[0]]['lattice_constant'] @@ -364,6 +365,9 @@ def _validate_all(self) -> 'Input': write_structure_file = True elif self.lattice.lower() in structure_dict.keys(): + if len(self.element) > 1: + raise ValueError("Cannot create lattice for more than one element, provide a lammps-data file explicitly") + #this is a valid structure if self.lattice_constant == 0: #we try try to get lattice_constant @@ -398,20 +402,21 @@ def _validate_all(self) -> 'Input': if not os.path.exists(self.lattice): raise ValueError(f'File {self.lattice} could not be found') if self.file_format == 'lammps-data': - aseobj = read(self.lattice, format='lammps-data', style='atomic') - structure = System(aseobj, format='ase') + #create atomic numbers for proper reading + Z_of_type = dict([(count+1, self._element_dict[element]['atomic_number']) for count, element in enumerate(self.element)]) + structure = read(self.lattice, format='lammps-data', style='atomic', Z_of_type=Z_of_type) + #structure = System(aseobj, format='ase') else: raise TypeError('Only lammps-data files are supported!') #extract composition - typelist = structure.atoms.types - #convert to species - typelist = [self.element[x-1] for x in typelist] - types, typecounts = np.unique(typelist, return_counts=True) + #this is the types read in from the file + types, typecounts = np.unique(structure.get_chemical_symbols(), return_counts=True) for c, t in enumerate(types): self._element_dict[t]['count'] = typecounts[c] self._element_dict[t]['composition'] = typecounts[c]/np.sum(typecounts) - self._natoms = structure.natoms + + self._natoms = len(structure) self._original_lattice = os.path.basename(self.lattice) self.lattice = os.path.abspath(self.lattice) From 08adaab1e44cc978d19feeda7038a0477190c287 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Fri, 19 Jan 2024 16:55:57 +0100 Subject: [PATCH 02/16] clean writing when lattice is created --- calphy/input.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/calphy/input.py b/calphy/input.py index 2c7d08f..25fbf80 100644 --- a/calphy/input.py +++ b/calphy/input.py @@ -380,10 +380,10 @@ def _validate_all(self) -> 'Input': lattice_constant=self.lattice_constant, repetitions=self.repeat, element=self.element) + structure = structure.write.ase() #extract composition - typelist = structure.atoms.species - types, typecounts = np.unique(typelist, return_counts=True) + types, typecounts = np.unique(structure.get_chemical_symbols(), return_counts=True) for c, t in enumerate(types): self._element_dict[t]['count'] = typecounts[c] @@ -393,7 +393,7 @@ def _validate_all(self) -> 'Input': #concdict_frac = {str(t): typecounts[c]/np.sum(typecounts) for c, t in enumerate(types)} #self._composition = concdict_frac #self._composition_counts = concdict_counts - self._natoms = structure.natoms + self._natoms = len(structure) self._original_lattice = self.lattice.lower() write_structure_file = True From 8b29d43aa5e86ef8b30c641294827de40c28ea65 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Fri, 19 Jan 2024 16:56:47 +0100 Subject: [PATCH 03/16] clean writing when no lattice is provided --- calphy/input.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/calphy/input.py b/calphy/input.py index 25fbf80..0cce0c4 100644 --- a/calphy/input.py +++ b/calphy/input.py @@ -351,16 +351,16 @@ def _validate_all(self) -> 'Input': lattice_constant=self.lattice_constant, repetitions=self.repeat, element=self.element) + structure = structure.write.ase() #extract composition - typelist = structure.atoms.species - types, typecounts = np.unique(typelist, return_counts=True) + types, typecounts = np.unique(structure.get_chemical_symbols(), return_counts=True) for c, t in enumerate(types): self._element_dict[t]['count'] = typecounts[c] self._element_dict[t]['composition'] = typecounts[c]/np.sum(typecounts) - self._natoms = structure.natoms + self._natoms = len(structure) self._original_lattice = self.lattice.lower() write_structure_file = True From b23a0a7d498653ea381e7ac4b55c638714f257af Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Fri, 19 Jan 2024 16:58:49 +0100 Subject: [PATCH 04/16] fix writing --- calphy/input.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/calphy/input.py b/calphy/input.py index 0cce0c4..cff326c 100644 --- a/calphy/input.py +++ b/calphy/input.py @@ -424,7 +424,7 @@ def _validate_all(self) -> 'Input': if write_structure_file: structure_filename = ".".join([self.create_identifier(), str(self.kernel), "data"]) structure_filename = os.path.join(os.getcwd(), structure_filename) - structure.write.file(structure_filename, format='lammps-data') + write(structure_filename, structure, format='lammps-data') self.lattice = structure_filename if self.mode == 'composition_scaling': From 8595d551a4d7f0987eae7985bc620eab32d8573b Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Fri, 19 Jan 2024 17:01:34 +0100 Subject: [PATCH 05/16] fix comp scale structure read --- calphy/input.py | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/calphy/input.py b/calphy/input.py index cff326c..3910b53 100644 --- a/calphy/input.py +++ b/calphy/input.py @@ -428,9 +428,6 @@ def _validate_all(self) -> 'Input': self.lattice = structure_filename if self.mode == 'composition_scaling': - aseobj = read(self.lattice, format='lammps-data', style='atomic') - structure = System(aseobj, format='ase') - #we also should check if actual contents are present input_chem_comp = {} for key, val in self._element_dict.items(): @@ -444,7 +441,7 @@ def _validate_all(self) -> 'Input': natoms1 = np.sum([val for key, val in self.composition_scaling._input_chemical_composition.items()]) natoms2 = np.sum([val for key, val in self.composition_scaling.output_chemical_composition.items()]) - if not (natoms1==natoms2==structure.natoms): + if not (natoms1==natoms2): raise ValueError(f"Input and output number of atoms are not conserved! Input {self.dict_to_string(self.input_chemical_composition)}, output {self.dict_to_string(self.output_chemical_composition)}, total atoms in structure {structure.natoms}") return self From 7367b836b7ef3fdcbec71ac333fa22856ff527d6 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 14:54:23 +0100 Subject: [PATCH 06/16] avoid nan spring constants --- calphy/solid.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/calphy/solid.py b/calphy/solid.py index 012c2c9..0369795 100644 --- a/calphy/solid.py +++ b/calphy/solid.py @@ -109,9 +109,13 @@ def analyse_spring_constants(self): k_std = [] for i in range(self.calc.n_elements): quant = np.loadtxt(file, usecols=(i+1, ), unpack=True)[-ncount+1:] - quant = 3*kb*self.calc._temperature/quant - k_mean.append(np.round(np.mean(quant), decimals=2)) - k_std.append(np.round(np.std(quant), decimals=2)) + mean_quant = np.round(np.mean(quant), decimals=2) + std_quant = np.round(np.std(quant), decimals=2) + if mean_quant == 0: + mean_quant = 1.00 + mean_quant = 3*kb*self.calc._temperature/mean_quant + k_mean.append(mean_quant) + k_std.append(std_quant) return k_mean, k_std From fba7320ab4bf20d9e61adf9b1fdf139671010777 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 15:02:42 +0100 Subject: [PATCH 07/16] dont divide by natoms in lammps script --- calphy/solid.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/calphy/solid.py b/calphy/solid.py index 0369795..8785566 100644 --- a/calphy/solid.py +++ b/calphy/solid.py @@ -423,7 +423,7 @@ def run_integration(self, iteration=1): lmp.command("variable step equal step") lmp.command("variable dU1 equal pe/atoms") for i in range(self.calc.n_elements): - lmp.command("variable dU%d equal f_ff%d/v_count%d"%(i+2, i+1, i+1)) + lmp.command("variable dU%d equal f_ff%d"%(i+2, i+1)) lmp.command("variable lambda equal f_ff1[1]") @@ -447,6 +447,9 @@ def run_integration(self, iteration=1): str2 = [] for i in range(self.calc.n_elements): str2.append("${dU%d}"%(i+2)) + for i in range(self.calc.n_elements): + str2.append("${count%d}"%(i+1)) + str2.append("${lambda}\"") str2 = " ".join(str2) str3 = " screen no file forward_%d.dat"%iteration @@ -472,6 +475,9 @@ def run_integration(self, iteration=1): str2 = [] for i in range(self.calc.n_elements): str2.append("${dU%d}"%(i+2)) + for i in range(self.calc.n_elements): + str2.append("${count%d}"%(i+1)) + str2.append("${lambda}\"") str2 = " ".join(str2) str3 = " screen no file backward_%d.dat"%iteration From ca61abc4aed0236a8deaf12c3589e1d07eb26112 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 16:40:01 +0100 Subject: [PATCH 08/16] modernise integrators --- calphy/integrators.py | 75 +++++++++++++++++++++++++++++-------------- calphy/solid.py | 18 +++++------ 2 files changed, 59 insertions(+), 34 deletions(-) diff --git a/calphy/integrators.py b/calphy/integrators.py index 41c5f79..bc8dbab 100644 --- a/calphy/integrators.py +++ b/calphy/integrators.py @@ -46,10 +46,12 @@ # TI PATH INTEGRATION ROUTINES #-------------------------------------------------------------------- -def integrate_path(fwdfilename, bkdfilename, - nelements=1, concentration=[1,], - usecols=(0, 1, 2), solid=True, - alchemy=False, composition_integration=False): +def integrate_path(calc, + fwdfilename, + bkdfilename, + solid=True, + alchemy=False, + composition_integration=False): """ Get a filename with columns du and dlambda and integrate @@ -72,23 +74,35 @@ def integrate_path(fwdfilename, bkdfilename, q : float heat dissipation during switching of system """ + natoms = np.array([calc._element_dict[x]['count'] for x in calc.element]) + concentration = np.array([calc._element_dict[x]['composition'] for x in calc.element]) + fdata = np.loadtxt(fwdfilename, unpack=True, comments="#") + bdata = np.loadtxt(bkdfilename, unpack=True, comments="#") + if solid: - fdui = np.loadtxt(fwdfilename, unpack=True, comments="#", usecols=(0,)) - bdui = np.loadtxt(bkdfilename, unpack=True, comments="#", usecols=(0,)) + fdui = fdata[0] + bdui = bdata[0] fdur = np.zeros(len(fdui)) bdur = np.zeros(len(bdui)) - for i in range(nelements): - fdur += concentration[i]*np.loadtxt(fwdfilename, unpack=True, comments="#", usecols=(i+1,)) - bdur += concentration[i]*np.loadtxt(bkdfilename, unpack=True, comments="#", usecols=(i+1,)) + for i in range(calc.n_elements): + fdur += concentration[i]*fdata[i+1]/natoms[i] + bdur += concentration[i]*bdata[i+1]/natoms[i] - flambda = np.loadtxt(fwdfilename, unpack=True, comments="#", usecols=(nelements+1,)) - blambda = np.loadtxt(bkdfilename, unpack=True, comments="#", usecols=(nelements+1,)) + flambda = fdata[calc.n_elements+1] + blambda = bdata[calc.n_elements+1] else: - fdui, fdur, flambda = np.loadtxt(fwdfilename, unpack=True, comments="#", usecols=usecols) - bdui, bdur, blambda = np.loadtxt(bkdfilename, unpack=True, comments="#", usecols=usecols) + fdui = fdata[0] + bdui = bdata[0] + + fdur = fdata[1] + bdur = bdata[1] + + flambda = fdata[2] + blambda = bdata[2] + fdu = fdui - fdur bdu = bdui - bdur @@ -106,11 +120,12 @@ def integrate_path(fwdfilename, bkdfilename, return w, q, flambda -def find_w(mainfolder, nelements=1, - concentration=[1,], nsims=5, - full=False, usecols=(0,1,2), +def find_w(mainfolder, + calc, + full=False, solid=True, - alchemy=False, composition_integration=False): + alchemy=False, + composition_integration=False): """ Integrate the irreversible work and dissipation for independent simulations @@ -142,15 +157,20 @@ def find_w(mainfolder, nelements=1, ws = [] qs = [] - for i in range(nsims): + for i in range(calc.n_iterations): fwdfilestring = 'forward_%d.dat' % (i+1) fwdfilename = os.path.join(mainfolder,fwdfilestring) + bkdfilestring = 'backward_%d.dat' % (i+1) bkdfilename = os.path.join(mainfolder,bkdfilestring) - w, q, flambda = integrate_path(fwdfilename, bkdfilename, - nelements=nelements, concentration=concentration, - usecols=usecols, solid=solid, - alchemy=alchemy, composition_integration=composition_integration) + + w, q, flambda = integrate_path(calc, + fwdfilename, + bkdfilename, + solid=solid, + alchemy=alchemy, + composition_integration=composition_integration) + ws.append(w) qs.append(q) @@ -444,7 +464,11 @@ def integrate_dcc(folder1, folder2, nsims=1, scale_energy=True, # REF. STATE ROUTINES: SOLID #-------------------------------------------------------------------- -def get_einstein_crystal_fe(temp, natoms, mass, vol, k, concentration, cm_correction=True): +def get_einstein_crystal_fe( + calc, + vol, + k, + cm_correction=True): """ Get the free energy of einstein crystal @@ -475,7 +499,10 @@ def get_einstein_crystal_fe(temp, natoms, mass, vol, k, concentration, cm_correc """ #convert mass first for single particle in kg - mass = (np.array(mass)/Na)*1E-3 + mass = np.array([calc._element_dict[x]['mass'] for x in calc.element]) + mass = (mass/Na)*1E-3 + natoms = np.sum([calc._element_dict[x]['count'] for x in calc.element]) + concentration = np.array([calc._element_dict[x]['composition'] for x in calc.element]) #convert k from ev/A2 to J/m2 k = np.array(k)*(eV2J/1E-20) diff --git a/calphy/solid.py b/calphy/solid.py index 8785566..4783b8e 100644 --- a/calphy/solid.py +++ b/calphy/solid.py @@ -447,8 +447,6 @@ def run_integration(self, iteration=1): str2 = [] for i in range(self.calc.n_elements): str2.append("${dU%d}"%(i+2)) - for i in range(self.calc.n_elements): - str2.append("${count%d}"%(i+1)) str2.append("${lambda}\"") str2 = " ".join(str2) @@ -475,8 +473,6 @@ def run_integration(self, iteration=1): str2 = [] for i in range(self.calc.n_elements): str2.append("${dU%d}"%(i+2)) - for i in range(self.calc.n_elements): - str2.append("${count%d}"%(i+1)) str2.append("${lambda}\"") str2 = " ".join(str2) @@ -520,13 +516,15 @@ def thermodynamic_integration(self): Calculates the final work, energy dissipation and free energy by matching with Einstein crystal """ - f1 = get_einstein_crystal_fe(self.calc._temperature, - self.natoms, [val['mass'] for key, val in self.calc._element_dict.items()], - self.vol, self.k, [val['composition'] for key, val in self.calc._element_dict.items()]) + f1 = get_einstein_crystal_fe( + self.calc, + self.vol, + self.k) + w, q, qerr = find_w(self.simfolder, - nelements=self.calc.n_elements, - concentration=[val['composition'] for key, val in self.calc._element_dict.items()], nsims=self.calc.n_iterations, - full=True, solid=True) + self.calc, + full=True, + solid=True) self.fref = f1 self.w = w From 6ae7875bd1c2fb65051c459fa618e9dbf4512514 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 16:42:29 +0100 Subject: [PATCH 09/16] fix temperature variable in integrators --- calphy/integrators.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/calphy/integrators.py b/calphy/integrators.py index bc8dbab..7ed59b3 100644 --- a/calphy/integrators.py +++ b/calphy/integrators.py @@ -515,12 +515,12 @@ def get_einstein_crystal_fe( F_cm = 0 for count, om in enumerate(omega): - F_harm += concentration[count]*np.log((hbar*om)/(kb*temp)) + F_harm += concentration[count]*np.log((hbar*om)/(kb*calc._temperature)) if cm_correction: - F_cm += np.log((natoms*concentration[count]/vol)*(2*np.pi*kbJ*temp/(natoms*concentration[count]*k[count]))**1.5) + F_cm += np.log((natoms*concentration[count]/vol)*(2*np.pi*kbJ*calc._temperature/(natoms*concentration[count]*k[count]))**1.5) #F_cm = 0 - F_harm = 3*kb*temp*F_harm - F_cm = (kb*temp/natoms)*F_cm + F_harm = 3*kb*calc._temperature*F_harm + F_cm = (kb*calc._temperature/natoms)*F_cm F_harm = F_harm + F_cm From d751893f654df1152d55e5992363e9230e78cac2 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 16:44:55 +0100 Subject: [PATCH 10/16] set guards to prevent log(0) call --- calphy/integrators.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/calphy/integrators.py b/calphy/integrators.py index 7ed59b3..eab48a8 100644 --- a/calphy/integrators.py +++ b/calphy/integrators.py @@ -515,9 +515,10 @@ def get_einstein_crystal_fe( F_cm = 0 for count, om in enumerate(omega): - F_harm += concentration[count]*np.log((hbar*om)/(kb*calc._temperature)) - if cm_correction: - F_cm += np.log((natoms*concentration[count]/vol)*(2*np.pi*kbJ*calc._temperature/(natoms*concentration[count]*k[count]))**1.5) + if concentration[count] > 0: + F_harm += concentration[count]*np.log((hbar*om)/(kb*calc._temperature)) + if cm_correction: + F_cm += np.log((natoms*concentration[count]/vol)*(2*np.pi*kbJ*calc._temperature/(natoms*concentration[count]*k[count]))**1.5) #F_cm = 0 F_harm = 3*kb*calc._temperature*F_harm F_cm = (kb*calc._temperature/natoms)*F_cm From 7f03691f8271cab13b025ac0186d548f3447f3fe Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 16:46:57 +0100 Subject: [PATCH 11/16] set guards to prevent div by 0 in integration --- calphy/integrators.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/calphy/integrators.py b/calphy/integrators.py index eab48a8..8fe0898 100644 --- a/calphy/integrators.py +++ b/calphy/integrators.py @@ -87,8 +87,9 @@ def integrate_path(calc, bdur = np.zeros(len(bdui)) for i in range(calc.n_elements): - fdur += concentration[i]*fdata[i+1]/natoms[i] - bdur += concentration[i]*bdata[i+1]/natoms[i] + if natoms[i] > 0: + fdur += concentration[i]*fdata[i+1]/natoms[i] + bdur += concentration[i]*bdata[i+1]/natoms[i] flambda = fdata[calc.n_elements+1] blambda = bdata[calc.n_elements+1] From 076624a0f93629c9bbebb3d0fc69c3889f56ff1a Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 17:08:40 +0100 Subject: [PATCH 12/16] fix ts mode --- calphy/phase.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/calphy/phase.py b/calphy/phase.py index e607398..fd4aeb4 100644 --- a/calphy/phase.py +++ b/calphy/phase.py @@ -56,11 +56,6 @@ def __init__(self, calculation=None, simfolder=None, log_to_screen=False): logfile = os.path.join(self.simfolder, "calphy.log") self.logger = ph.prepare_log(logfile, screen=log_to_screen) - self.logger.info("---------------input file----------------") - self.logger.info("commented out as causes crash when we're expanding the T range after a fail run") - # self.logger.info(yaml.safe_dump(self.calc.to_dict())) - self.logger.info("------------end of input file------------") - if self.calc._pressure is None: pressure_string = "None" else: From 9233e38d3e1888ae84f446d2026cfd3dc771c9d1 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 17:14:20 +0100 Subject: [PATCH 13/16] remove ghost elements from set mass function --- calphy/alchemy.py | 4 ++-- calphy/helpers.py | 12 +++--------- calphy/liquid.py | 4 ++-- calphy/phase.py | 8 ++++---- calphy/solid.py | 6 +++--- 5 files changed, 14 insertions(+), 20 deletions(-) diff --git a/calphy/alchemy.py b/calphy/alchemy.py index e78fbea..d114854 100644 --- a/calphy/alchemy.py +++ b/calphy/alchemy.py @@ -84,7 +84,7 @@ def run_averaging(self): #set up potential lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #add some computes lmp.command("variable mvol equal vol") @@ -157,7 +157,7 @@ def run_integration(self, iteration=1): #set up hybrid potential #here we only need to set one potential lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #lmp = ph.set_double_hybrid_potential(lmp, self.options, self.calc._pressureair_style, self.calc._pressureair_coeff) diff --git a/calphy/helpers.py b/calphy/helpers.py index afbabd2..e6295ab 100644 --- a/calphy/helpers.py +++ b/calphy/helpers.py @@ -118,19 +118,13 @@ def create_structure(lmp, calc): return lmp -def set_mass(lmp, options, ghost_elements=0): - count = 1 +def set_mass(lmp, options): for i in range(options.n_elements): lmp.command(f'mass {i+1} {options.mass[i]}') - count += 1 - - for i in range(ghost_elements): - lmp.command(f'mass {count+i} 1.00') - return lmp -def set_potential(lmp, options, ghost_elements=0): +def set_potential(lmp, options): """ Set the interatomic potential @@ -149,7 +143,7 @@ def set_potential(lmp, options, ghost_elements=0): lmp.command(f'pair_style {options._pair_style_with_options[0]}') lmp.command(f'pair_coeff {options.pair_coeff[0]}') - lmp = set_mass(lmp, options, ghost_elements=ghost_elements) + lmp = set_mass(lmp, options) return lmp diff --git a/calphy/liquid.py b/calphy/liquid.py index 8289bc0..fe7db4b 100644 --- a/calphy/liquid.py +++ b/calphy/liquid.py @@ -120,7 +120,7 @@ def run_averaging(self): #set up potential lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #Melt regime for the liquid lmp.velocity("all create", self.calc._temperature_high, np.random.randint(1, 10000)) @@ -189,7 +189,7 @@ def run_integration(self, iteration=1): #set hybrid ufm and normal potential #lmp = ph.set_hybrid_potential(lmp, self.options, self.eps) lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #remap the box to get the correct pressure lmp = ph.remap_box(lmp, self.lx, self.ly, self.lz) diff --git a/calphy/phase.py b/calphy/phase.py index fd4aeb4..3fd84c3 100644 --- a/calphy/phase.py +++ b/calphy/phase.py @@ -732,7 +732,7 @@ def reversible_scaling(self, iteration=1): #set up potential lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #remap the box to get the correct pressure lmp = ph.remap_box(lmp, self.lx, self.ly, self.lz) @@ -821,7 +821,7 @@ def reversible_scaling(self, iteration=1): else: self.check_if_solidfied(lmp, "traj.temp.dat") - lmp = ph.set_potential(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_potential(lmp, self.calc) #reverse scaling lmp.command("variable flambda equal ramp(${li},${lf})") @@ -923,7 +923,7 @@ def temperature_scaling(self, iteration=1): #set up potential lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #remap the box to get the correct pressure lmp = ph.remap_box(lmp, self.lx, self.ly, self.lz) @@ -1015,7 +1015,7 @@ def pressure_scaling(self, iteration=1): #set up potential lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #remap the box to get the correct pressure lmp = ph.remap_box(lmp, self.lx, self.ly, self.lz) diff --git a/calphy/solid.py b/calphy/solid.py index 4783b8e..053a143 100644 --- a/calphy/solid.py +++ b/calphy/solid.py @@ -226,7 +226,7 @@ def run_interactive_averaging(self): lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') else: lmp.command("include %s"%self.calc.potential_file) - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #add some computes lmp.command("variable mvol equal vol") @@ -305,7 +305,7 @@ def run_minimal_averaging(self): lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') else: lmp.command("include %s"%self.calc.potential_file) - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #add some computes @@ -388,7 +388,7 @@ def run_integration(self, iteration=1): lmp.command(f'pair_coeff {self.calc.pair_coeff[0]}') else: lmp.command("include %s"%self.calc.potential_file) - lmp = ph.set_mass(lmp, self.calc, ghost_elements=self.calc._ghost_element_count) + lmp = ph.set_mass(lmp, self.calc) #remap the box to get the correct pressure lmp = ph.remap_box(lmp, self.lx, self.ly, self.lz) From 2d020929bce3904ec9155b3046498c8c71372f21 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 17:15:35 +0100 Subject: [PATCH 14/16] formatting in liquid --- calphy/liquid.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/calphy/liquid.py b/calphy/liquid.py index fe7db4b..a04a611 100644 --- a/calphy/liquid.py +++ b/calphy/liquid.py @@ -331,7 +331,8 @@ def thermodynamic_integration(self): self.rho, 50, 1.5) #Get ideal gas fe - f2 = get_ideal_gas_fe(self.calc._temperature, self.rho, + f2 = get_ideal_gas_fe(self.calc._temperature, + self.rho, self.natoms, [val['mass'] for key, val in self.calc._element_dict.items()], [val['composition'] for key, val in self.calc._element_dict.items()]) From 062a68c55bd3cee7254dc1d2df5b4590541a28fa Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 17:17:04 +0100 Subject: [PATCH 15/16] remove unused ghost atom variable from input --- calphy/input.py | 1 - 1 file changed, 1 deletion(-) diff --git a/calphy/input.py b/calphy/input.py index 3910b53..2cee92e 100644 --- a/calphy/input.py +++ b/calphy/input.py @@ -202,7 +202,6 @@ class Calculation(BaseModel, title='Main input class'): #add second level options; for example spring constants spring_constants: Annotated[List[float], Field(default = None)] - _ghost_element_count: int = PrivateAttr(default=0) #structure items _structure: Any = PrivateAttr(default=None) From fcc9526a53767c60fc0fe9718220f9d49377c4a9 Mon Sep 17 00:00:00 2001 From: Sarath Menon Date: Tue, 20 Feb 2024 17:19:17 +0100 Subject: [PATCH 16/16] modify tests --- tests/test_integrators.py | 4 ---- 1 file changed, 4 deletions(-) diff --git a/tests/test_integrators.py b/tests/test_integrators.py index b7cb5af..802e2ad 100644 --- a/tests/test_integrators.py +++ b/tests/test_integrators.py @@ -8,7 +8,3 @@ def test_ideal_gas(): def test_uf(): a = get_uhlenbeck_ford_fe(1000, 0.07, 50, 2) assert np.abs(a-5.37158083028874) < 1E-5 - -def test_solid_ref(): - a = get_einstein_crystal_fe(1000, 1000, [26], 200, [1.2], [1]) - assert np.abs(a+0.4727067261423942) < 1E-5