Predict dipole during MD

[elindgren]

It would be very useful to be able to compute the dipole during MD. My idea was to use dump_observer, where one would supply two models, nep_main and nep_dipole. nep_main would be used to run the MD and nep_dipole would be used to compute the dipole as an observer model.

In order to implement the above scheme we would need a function for computing dipoles in, for instance, force/nep3.cu. Also, it would be beneficial if there was a flag in the nep.txt file that indicated what type of model it was (regular, polarizability, dipole), which I don't believe is done atm? @brucefan1983 Have you been planning to implement a function for calculating dipoles during MD?

[brucefan1983]

These extensions are in principle possible, but I have not found time yet. Some extensions would break existing code for reading the nep.txt file and invole a lot of refactoring within GPUMD. Have you tried the predicition mode of nep, which should be very fast as it predicts in a batched way. Most of the time is simply spent on reading the xyz files and it is also easy to automate it within calorine by calling the nep executable and creating and destroying some tempeary train.xyz files.

[elindgren]

I've been playing around a little bit with the prediction mode for dipoles, but I don't seem to get consistent results with the results from predicting the dipole via NEP_CPU. The CPU-predicted dipole matches the target DFT dipole well, but the ones from prediction are about one order of magnitude smaller. @brucefan1983 is the output predicted dipoles in dipole_train.out perhaps normalized?

Input parameters for reference

Input or default parameters:
    (input)   train_mode = dipole.
    (input)   calculation mode = predict.
    (input)   use NEP version 4.
    (input)   number of atom types = 4.
        (default) type 0 (Si with Z = 14) has force weight of 1.
        (default) type 1 (F with Z = 9) has force weight of 1.
        (default) type 2 (C with Z = 6) has force weight of 1.
        (default) type 3 (H with Z = 1) has force weight of 1.
    (default) will not add the ZBL potential.
    (input)   radial cutoff = 8 A.
    (input)   angular cutoff = 6 A.
    (input)   n_max_radial = 15.
    (input)   n_max_angular = 8.
    (input)   basis_size_radial = 8.
    (input)   basis_size_angular = 8.
    (input)   l_max_3body = 4.
    (input)   l_max_4body = 2.
    (input)   l_max_5body = 0.
    (input)   number of neurons = 80.

nep4_dipole.txt

[elindgren]

Now that I think about it, the outputs are probably dipole/atom right?

[brucefan1983]

Yes, like energy and virial, it is /atom.

[brucefan1983]

Because we wanted to have nice parity plots.

[brucefan1983]

In GPUMD-v3.9, we will output dipole/polarizability within GPUMD run (need both NEP potential model and NEP dipole/polarizability model), and also calcualte IR and Raman spectra on the fly.

[elindgren]

In GPUMD-v3.9, we will output dipole/polarizability within GPUMD run (need both NEP potential model and NEP dipole/polarizability model), and also calcualte IR and Raman spectra on the fly.

Nice! That will be very useful :rocket: