NPT relaxation issue

[matukumilli]

Dear Developers of SSCHA,

When I do a NPT relaxation using the script from SimpleVCRelax example, the cell volume increases monotonously.
I have to do this in this manual way as I am using MLFF generated in VASP.
Both the gradients (Str, FC) are staying far from zero but slowly coming down.
Volume simply increasing without any fluctuations. When I set target_pressure to a small finite value, there is a slight change in rate of volume increase but no trend of convergence.
Please help me to understand what's going wrong.
Please note that this MLFF worked well to minimise free energy at fixed volume.

&relax
type = "vc-relax"
start_pop = 1
max_pop_id = 1
generate_ensemble = .false.
fix_volume = .false.
target_pressure = 0 ! [GPa]
bulk_modulus = 8 ! [GPa]
n_configs = 100
&end

&inputscha
n_random = 100
n_random_eff = 50
data_dir = "pop1/data"
population = 1
fildyn_prefix = "Ddyn/dynq_after_"
nqirr = 8
supercell_size = 2 2 2
Tg = 200
T = 200
meaningful_factor = 1e-4
gradi_op = "all"
print_stress = .true.
eq_energy = -30.7753
lambda_a = 0.5
lambda_w = 0.5
root_representation = "normal"
preconditioning = .true.
/

I use the above script along with ensemble prepared via MLFF in loop.

Regards,
Prasad

[diegomartinez2]

Hi, is not your bulk_modulus too small?.
Set to a higher value to get a lower change in the cell shape, something like x10 times what you got.
Also I think you can use "recalc" to use finite differences to get the static bulk modulus.

[mesonepigreco]

Hi,
Unfortunately, doing one step at time the code is not able to automatically detect the best unit cell step, and uses the value of the bulk modulus provided.
In your case, it seems that the bulk modulus of your system is even smaller than the one you selected there. You can have an estimation by hand.
Compare the pressure between different calculation at different volumes, how much do they change? Then guess how much bigger the step in volume should have been to make that pressure go to zero, since the step in the volume is ~ 1/B (B is the Bulk Modulus), edit the value of B accordingly to get a bigger step and converge much faster.

In practice, this means doing by hand what the cell relaxation code does between different steps to reduce the number of calculations.

[mesonepigreco]

You should remove the target_pressure flag if you fix the volume to true. Sorry, probably is not so clear from that input

[matukumilli]

Please have a look at my input and output files. I already removed target_pressure line.

vcr-fixVol.tar.gz

[matukumilli]

Even after removing target_pressure line, it is still printing P = 0.0000 GPa.

Are those attached input output files are of any use to understand this behaviour?

Could you please help me in NPT relaxation?

Thanks,
Prasad

[mesonepigreco]

Hi, there is a problem in your calculation:

 ==== STRESS TENSOR [GPa] ==== 
      0.15670724      0.00000000      0.00000000                0.00118771      0.00000000      0.00000000
      0.00000000      0.15841902      0.00000000    +-          0.00000000      0.00133348      0.00000000
      0.00000000      0.00000000      0.15966600                0.00000000      0.00000000      0.00123602

 Ab initio average stress [GPa]:
      0.00000000      0.00000000      0.00000000
      0.00000000      0.00000000      0.00000000
      0.00000000      0.00000000      0.00000000


 ======================
 ENTHALPIC CONTRIBUTION
 ======================

 P = 0.0000 GPa   V = 414.9813 A^3

 P V = 0.00000000e+00 eV

 Helmoltz Free energy = -4.6072217522e+01 eV <--
 Gibbs Free energy = 1.6588394861e+02 eV 
 Zero energy = -2.1195616613e+02 eV

The code is actually performing a fixed volume minimization, as it is optimizing the Helmoltz free energy instead of the Gibbs free energy, so that part works correctly (the small arrow indicates which is the value that is important for the minimization).

However, you have an Ab initio average stress of exactly 0, which probably means that the code is not correctly loading the stress tensors from the ab'initio calculation.

In data_dir = "../pop1/data", did you store the pressures_population1_X.dat with the stress tensors in Ry/Bohr^3 from the ab-initio calculations?

[matukumilli]

yes, you are right, I missed the pressure_..dat files generation part , only extracted energies and forces. Thanks for pointing it out.

Now, even after having all required files in data_dir, it still printing-out P as zero.

 ==== STRESS TENSOR [GPa] ====
     -1.21298413      0.00000000      0.00000000                0.01117302      0.00000000      0.00000000
      0.00000000     -1.07123555      0.00000000    +-          0.00000000      0.01770616      0.00000000
      0.00000000      0.00000000     -0.99850654                0.00000000      0.00000000      0.01385735

 Ab initio average stress [GPa]:
     -1.29622476      0.00000000      0.00000000
      0.00000000     -1.15592997      0.00000000
      0.00000000      0.00000000     -1.08162879


 ======================
 ENTHALPIC CONTRIBUTION
 ======================

 P = 0.0000 GPa   V = 410.3506 A^3

 P V = 0.00000000e+00 eV

 Helmoltz Free energy = -4.4830219540e+01 eV <--
 Gibbs Free energy = 1.6712594659e+02 eV
 Zero energy = -2.1195616613e+02 eV

Ab initio average stresses are finite now.
Why pressure still showing as 0, could you please help me here? I tried with both ON/OFF fix_volume tag, but no change: P = 0.0000 GPa .

Warm regards,
Prasad

[mesonepigreco]

Hi, 0 is the target pressure, so it computes the enthalpy (or gibbs free energy) based on the target pressure.
In fact, the stress tensor is not 0 but approximately -1 GPa, which is your correct pressure.

Maybe is a bit confusing. I will add an issue to clarify that P is the target pressure in the output.
The reason why we use target pressure is that it reduces the error on the final estimation of the Free energy.

Your output now indicates that everything is correct and VC relax should work as intended (even with the full vc relax without fix_volume)

[matukumilli]

Got it. Thanks Dr. Lorenzo for your continued support in resolving this issue.

Regards,
Prasad

[chenyubi14]

Dear developers,

I am currently exploring NPT relaxation using VASP. In the provided tutorial, relax.vc_relax(target_press = 0) is used for this purpose.
However, my preference leans towards manual submission with VASP because relax.vc_relax executes VASP one job at a time, which is very slow. Instead, I prefer to submit multiple VASP jobs in the ensemble, and process SSCHA relaxation by

## generate ensemble
ensemble.generate(number_of_configurations)
ensemble.save(directory, population_id)

## load ensemble
ensemble.load(directory, population = population_id, N = number_of_configurations)
minim = sscha.SchaMinimizer.SSCHA_Minimizer(ensemble)
minim.init()

This process is specifically tailored for NVT relaxation. However, I am encountering difficulties in configuring a similar setup for NPT relaxation. In minim.init(), there is no setup like target_press for NPT.
From my understanding of the code, relax.vc_relax is the only method for NPT, but relax.vc_relax is not compatible with the manual submission that I prefer.

Additionally, I noticed the first user of this discussion is able to run NPT by VASP. His method uses a file containing the SSCHA inputs, rather than using current Python code tutorials.

&relax
...
&end

&inputscha
...
/

How can I use NPT for VASP manual submission? Specifically, could you provide me with an example code to adapt minim.init() for an NPT ensemble?

Thanks a lot,
Yubi