feram: MD simulator for bulk and thin-film ferroelectrics
The parameter file for feram is a text file consisting of comment lines and 'tag = value(s)' lines. Filenames of parameter files may be foo123 or foo123.feram. In these cases, filenames of output files will be foo123.*,e.g. foo123.log, foo123.0000050000.cood, etc.
How to determine the parameters is described in parameters/parameters.html.
Lines beginning with '#' are ignored. Blank lines are also ignored.
# This is a comment line. # Here are two more # comment lines.
You can write a shebang at the first line of a .feram input file and execute the file as a script.
You must put ' = ', space-equal-space, between tag and value(s) as:
tag = 1.0 tag = -2.0 -3.0 -4.0 tag = 5.0 6.0 7.0
Followings are currently available tags.
verbose tag determines how much messages will be written to the .log file. (Details of .log file is described below.)
verbose = 0
(NOT IMPLEMENTED YET) .log file will not be created.
verbose = 1
1 is the default value for verbose tag. Energies of the first and last iterations will be reported in the .log file.
verbose = 2
When 2 is given, energies will be reported in the .log file for EVERY iterations.
verbose = 3
When 3 is given, detailed messages will be written into the .log file.
verbose = 4
When 4 is given, more detailed messages will be written into the .log file.
method = 'md'
'md' is for a molecular-dynamics simulation in the canonical ensemble using the Nosé-Poincaré thermostat.
method = 'vs'
'vs' is for a molecular-dynamics simulation in the canonical ensemble using the velocity scaling thermostat.
method = 'lf'
'lf' is for a molecular-dynamics simulation in the microcanonical ensemble using the leap-frog method.
method = 'mc'
'mc' is for a Monte Carlo simulation, but it is NOT IMPLEMENTED YET.
method = 'hl'
This 'hl' feature will be obsolete. For simulations of hysteresis loops, use n_E_wave_period and E_wave_type instead.
Pressure in GPa unit.
GPa = -5.0
You can apply uniaxial normal stress by giving three values to GPa.
GPa = 0.0 0.0 -3.0
Temperature in Kelvin.
kelvin = 100
Effective mass of u(R), the length of dipole on the unit cell R. In atomic mass unit.
If acoustic_mass_amu is not set, acoustic displacements, i.e. inhomogeneous strains, are optimized according to the {u(R)}. If it is set to the effective mass of acoustic displacements
acoustic_mass_amu = 46.44
molecular dynamics is performed on acoustic displacements. Its unit is atomic mass unit. Generally, it is average of mass of 5 atoms. You cannot use method = 'md' together with acoustic_mass_amu, currently.
Q_Nose for Nose-Poincare thermostat. For large system and high temperature, use large value. For Small system and low temperature, use small value.
Q_Nose = 14.4
Specifies the structure of the system. 'bulk' for infinitely periodic bulk 'film' for free standing thin film. 'epit' for epitaxially strained thin film. See epi_strain tag.
bulk_or_film = 'epit'
System size, L_x, L_y and L_z. They must be equal or larger than 4.
L = 32 32 4
Padding for speeding up the calculations. You do not have to give this value, if you cannot understand it. Best padding_y is system-dependent, but for 64x64x64 and 96x96x96
padding_y = 5
for 128x128x128
padding_y = 1
may be good. Its default value is 3.
0: There is no dead layer. 1: Single dead layer is on the bottom electrode. 2: There are dead layers both on the bottom and top electrodes.
gap_id = 1
Thickness of dead layer is always 1 unit cell.
Polarization of dead layer(s) in Angstrom unit.
gap_dipole_u = 0.0 0.0 0.16
Lattice constant for the perovskite ferroelectrics in Angstrom unit.
a0 = 3.99
Epitaxial strain.
epi_strain = -0.01
Note that epi_strain is effective only when bulk_or_film = 'epit'. See src/optimize-homo-strain.F for more details.
Time step in pico second.
dt = 0.002
The number of time steps of thermalizing, averaging, frequency of taking snapshots, and frequency of writing to .hl file.
n_thermalize = 40000 n_average = 10000 n_coord_freq = 50000 n_hl_freq = 5000
n_E_wave_period is the period of alternating external electric field which have external_E_field amplitude. Set 'triangular_sin' or 'triangular_cos' to E_wave_type. See example files in src/27example-BaTiO3-new-param-E_wave/ . n_E_wave_period must be a positive integer divisible by 4 for 'triangular_sin' or 'triangular_cos'. You can also use 'ramping_off' and 'ramping_on'. n_E_wave_period must no be divisible by 4 for 'ramping_on' and 'ramping_off' cases.
You can restart and continue a hysteresis-loop calculation from the n_hysteresis_loop_continue-th iteration.
n_hysteresis_loop_continue = 20000
If this tag is set as
coord_directory = 'coord_files'
.coord file will be written into coord_files/foo123.0000020000.coord . For example, use this tag to store huge .coord files into a local HDD.
If you do not want to get any .coord files, set it as
coord_directory = 'never'
If this tag is set as
distribution_directory = 'distribution_files'
.distribution and .distribution3d files will be written into distribution_files/foo123.distribution and distribution_files/foo123.distribution3d . For example, use this tag to store the huge .distribution3d file into a local HDD.
If you do not want to get any .distribution and .distribution3d files, set it as
distribution_directory = 'never'
This can also reduce computational time.
Vector of the external electric field in the unit of [V/Angstrom].
external_E_field = 0.00 0.00 -0.01
External electric field will be constant in this value when method is 'md', 'vs' or 'lf'. Or, if E_wave_type is specified, E will be altering in period n_E_wave_period.
Coefficients for the 4-th order polynomial.
P_kappa2 = 5.502 [eV/Angstrom^2] # P_4(u) = kappa2*u^2 + alpha*u^4 P_alpha = 110.4 [eV/Angstrom^4] # + gamma*(u_y*u_z+u_z*u_x+u_x*u_y), P_gamma = -163.1 [eV/Angstrom^4] # where u^2 = u_x^2 + u_y^2 + u_z^2
Coefficients for the 6-th order part
Coefficient for u^8
Short-range inter-site interaction coefficients
j = -2.648 3.894 0.898 -0.789 0.562 0.358 0.179 j(i) [eV/Angstrom^2]
Elastic Constants
B11 = 126. B12 = 44.9 B44 = 50.3 [eV]
Elastic Coupling
B1xx = -211. [eV/Angstrom^2] B1yy = -19.3 [eV/Angstrom^2] B4yz = -7.75 [eV/Angstrom^2]
Two integers for Marsaglia-Tsang 64-bit universal RNG (random number generator).
seed = 123384 8984847
Its default values are 123456789 987654321. It may be a good idea to generate the seed with jot(1) command.
$ jot -r 2 1 2147483648 639779840 1297280417
You can also use $((RANDOM)) of bash(1).
$ echo $((RANDOM*65536+RANDOM)) $((RANDOM*65536+RANDOM)) 2146723680 1863789493
Values for initial random setting of dipoles. The normal-distributed random number generator N(mu,sigma^2) is used.
For example, if you only have four dipoles and
init_dipo_avg = 0.10 0.00 0.00 [Angstrom] # Average (mu) of initial dipole displacements init_dipo_dev = 0.02 0.02 0.02 [Angstrom] # Deviation (sigma) of initial dipole displacements,
the dipoles may be
0.09 0.01 0.02 0.11 -0.02 0.01 0.08 0.02 -0.03 0.12 -0.01 0.00 .
If there is a .restart file in the same directory, there values are not used and the .restart file will be read.
Effective charge per site.
Z_star = 9.956
epsilon_inf = 5.24
This is for debug. You do not have to set this tag. Arbitrary kappa for Ewald summation in src/dipole-dipole-long-range.F.
kappa = 0.15
If
plot_dispersion = .true.
data files for plotting dispersion relation will be written. Default is .false..
For plot_dispersion = .true., L_x, L_y and L_z must be a same even number.
L = 32 32 32
may be a good idea.
See an example of src/08example-BaTiO3-dispersion/ or try
$ tar xf SOMEWHERE/feram-X.YY.ZZ.tar.xz
$ cd feram-X.YY.ZZ
$ mkdir builddirectory
$ cd builddirectory
$ ../configure && make -j4
$ make check TESTS='dispersion_check.sh dispersion_32x32x32.inhomo-K.gp'
$ less dispersion_32x32x32.feram # the input file
$ less dispersion_32x32x32.gp # the GNUPLOT script
$ gv src/dispersion_32x32x32.long+short.energy.eps &
$ gv src/dispersion_32x32x32.long+short.eps &
$ gv src/dispersion_32x32x32.long+short.interaction.eps &
$ gv src/dispersion_32x32x32.inhomo-K.eps &
# Previews of plots. You may also use evince(1) instead of gv(1).
Compare the plots with Fig.3 in [Takeshi Nishimatsu el al.: Phys. Rev. B 82, 134106 (2010), http://dx.doi.org/10.1103/PhysRevB.82.134106 ] and Fig.5 in parameters/parameters.html.
If feram is executed as
$ feram foo.feram bar.feram
and, in bar.feram,
continue = .true.
is declared, simulation of bar.feram uses results of foo.feram as initial configurations. The .restart file will not be used. Default value of continue is .false..
If you executed feram with a parameter file "foo123.feram" as
$ feram foo123.feram
feram searches for the "foo123.restart" file as an initial configurations. If there is no such file, feram starts from a random configurations.
Quadratic inter-atomic force constant (IFC) matrices can be read from .quadratic file. Set all j1-j7 to ZERO. Use P_kappa2 for *elevation*. (Lx/2+1)*Ly*Lz 3x3 matrices are required. See src/25example-BaTiO3-read-quadratic/ for more details.
If you executed feram with a parameter file "foo123.feram" as
$ feram foo123.feram
feram searches for the "foo123.defects" file as information of defects. If there is no such file, feram do not introduce defects.
The .defects input file gives positions and fixed dipole moments u for simulations of defects in ferroelectrics.
0 0 0 0.1 -0.1 0.7 1 0 0 0.1 0.1 0.7 2 0 0 0.1 0.2 0.7 : position_x position_y position_z fixed_u_x fixed_u_y fixed_u_z [Angstrom] :
Use src/feram_defects_maker.rb to make this .defects file.
If you executed feram with a parameter file "foo123.feram" as
$ feram foo123.feram
feram searches for the "foo123.localfield" file as information of local electric field. Format of .localfield file is:
4 6 13 0.00 0.00 0.10 : ix iy iz Ex Ey Ez :
The unit of local field is V/Angstrom. If feram cannot find both ionic.configuration and .localfield, it will set the local field to ZERO for all sites.
Use attached feram_random_field_generator for generating random fields.
$ ./feram_random_field_generator Lx Ly Lz mu_x mu_y mu_z sigma_x sigma_y sigma_z seed1 seed2 $ ./feram_random_field_generator 32 32 243 0.0 0.0 0.0 0.02 0.02 0.02 123456789 987654321
.modulation file gives acoustic modulation, a(R), of each site. Energy from acoustic modulation is
E_acoustic_modulation = modulation_constant Sigma_R_alpha eta_alpha(R) a(R) .
By preparing an FFTW wisdom file in your current directory, you can reduce the total computational time. Sometimes, however, this trick fails and the computational times increases. You should optimize the computational condition with trial-and-error.
You can generate an FFTW wisdom file with `feram_fftw_wisdom` command, which will be built in the src/ directory besides `feram` command. In the src/18example-benchmark/ directory, for example, you can generate an FFTW wisdom file of system size Lx*Ly*Lz=32*32*243 for the forward.feram input file,
$ w # Before executing feram_fftw_wisdom, confirm that load average is zero and there is no background job.
$ OMP_NUM_THREADS=6 ../feram_fftw_wisdom 1000 32 32 243 3
:
It takes a few minutes or a few hours depending on the system size.
'1000' in the command line arguments is the number of iterations of
FFT benchmarks. So, do not set a large number for it.
:
$ mv wisdom_new wisdom
$ OMP_NUM_THREADS=6 ../feram forward.feram
$ grep FFTW_WISDOM forward.log
feram_common.F: 47: FFTW_WISDOM: Successfully imported FFTW wisdom in current directory.
$ grep TIMING_REPORT forward.log
Command line arguments for feram_fftw_wisdom are
feram_fftw_wisdom <N_TIMES> <Lx> <Ly> <Lz> <padding_y>
Check http://www.fftw.org/doc/Other-Important-Topics.html for more details on FFTW wisdom. If you are using the MKL for FFTW3 wrappers, you do not have to prepare the wisdom file.
You can find input files and their simulated results in following directories in the software package.
feram-X.YY.ZZ/src/03example-BaTiO3-bulk-phase-transition/ feram-X.YY.ZZ/src/09example-BaTiO3-bulk-leap-frog/ feram-X.YY.ZZ/src/17example-PbTiO3-100-900K/ feram-X.YY.ZZ/src/20example-BaTiO3-new-param-pyro/ feram-X.YY.ZZ/src/21example-KNbO3/ feram-X.YY.ZZ/src/22example-BaTiO3-new-param-hl/ feram-X.YY.ZZ/src/23example-BaTiO3-new-param-hl-film-acoustic-mass/ feram-X.YY.ZZ/src/24example-BaTiO3-new-param-heating-cooling/ feram-X.YY.ZZ/src/25example-BaTiO3-read-quadratic/ feram-X.YY.ZZ/src/26example-BaTiO3-acoustic-MD/ feram-X.YY.ZZ/src/27example-BaTiO3-new-param-E_wave/ feram-X.YY.ZZ/src/28example-PbTiO3-elastocaloric-770K/ feram-X.YY.ZZ/src/29example-BaTiO3-new-param-defects/ feram-X.YY.ZZ/src/30example-BaTiO3-new-param-ramping-off/ feram-X.YY.ZZ/src/31example-BaTiO3-new-param-thermal-conductivity/ feram-X.YY.ZZ/src/32example-BaTiO3-old-64x64x1024-benchmark/ feram-X.YY.ZZ/src/34example-BST/
This project was partially supported by:
http://loto.sourceforge.net/feram/ is hosted by .
