Skip to content

Commit bc3df2a

Browse files
committed
Add a new explanation for the VC Relax
====================================== What is the target pressure and what the current pressure.
1 parent dcaf70a commit bc3df2a

1 file changed

Lines changed: 14 additions & 6 deletions

File tree

Modules/Relax.py

Lines changed: 14 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -582,8 +582,6 @@ def vc_relax(self, target_press = 0, static_bulk_modulus = 100,
582582
if kind_minimizer == "RPSD":
583583
# Compute the static bulk modulus
584584
sbm = GetStaticBulkModulus(self.minim.dyn.structure, self.calc)
585-
print ("BM:")
586-
print (sbm)
587585
BFGS = sscha.Optimizer.SD_PREC_UC(self.minim.dyn.structure.unit_cell, sbm)
588586

589587
# Generate the ensemble
@@ -603,16 +601,13 @@ def vc_relax(self, target_press = 0, static_bulk_modulus = 100,
603601
self.calc, True, stress_numerical, cluster = self.cluster)
604602
#self.minim.ensemble.get_energy_forces(self.calc, True, stress_numerical = stress_numerical)
605603

606-
print("RELAX force length:", len(self.minim.ensemble.force_computed))
607604

608605
if ensemble_loc is not None and self.save_ensemble:
609606
self.minim.ensemble.save_bin(ensemble_loc, pop)
610-
print("RELAX force length:", len(self.minim.ensemble.force_computed))
611607

612608
self.minim.population = pop
613609
self.minim.init(delete_previous_data = False)
614610

615-
print("RELAX force length:", len(self.minim.ensemble.force_computed))
616611
self.minim.run(custom_function_pre = self.__cfpre__,
617612
custom_function_post = self.__cfpost__,
618613
custom_function_gradient = self.__cfg__)
@@ -629,6 +624,14 @@ def vc_relax(self, target_press = 0, static_bulk_modulus = 100,
629624
# Get the pressure
630625
Press = np.trace(stress_tensor) / 3
631626

627+
# Get the error correctly accounting for symmetries
628+
stress_diagonal_err = np.diag(stress_err)
629+
press_err = 0
630+
for i in range(3):
631+
if not stress_diagonal_err[i] in stress_diagonal_err[:i]:
632+
press_err += stress_diagonal_err[i]**2
633+
press_err = np.sqrt(press_err)
634+
632635
# Get the volume
633636
Vol = self.minim.dyn.structure.get_volume()
634637

@@ -656,6 +659,11 @@ def vc_relax(self, target_press = 0, static_bulk_modulus = 100,
656659
ENTHALPIC CONTRIBUTION
657660
======================
658661
662+
Current pressure P = {:.4f} +- {:.4f} GPa
663+
Target pressure P = {:.4f} GPa
664+
665+
For enthalpy we use the target pressure.
666+
659667
P = {:.4f} GPa V = {:.4f} A^3
660668
661669
P V = {:.8e} eV
@@ -664,7 +672,7 @@ def vc_relax(self, target_press = 0, static_bulk_modulus = 100,
664672
Gibbs Free energy = {:.10e} eV {}
665673
Zero energy = {:.10e} eV
666674
667-
""".format(target_press , Vol,target_press_evA3 * Vol, helmoltz, mark_helmoltz, gibbs, mark_gibbs, self.minim.eq_energy)
675+
""".format(Press / CC.Units.GPA_TO_EV_PER_A3, press_err/ CC.Units.GPA_TO_EV_PER_A3, target_press , Vol,target_press_evA3 * Vol, helmoltz, mark_helmoltz, gibbs, mark_gibbs, self.minim.eq_energy)
668676
print(message)
669677
# print " ====================== "
670678
# print " ENTHALPIC CONTRIBUTION "

0 commit comments

Comments
 (0)