TiAl potential (ACE1pack-julia)

In this tutorial we repeat what was done in [Fitting a TiAL potential][TiAl.md], but only using ACE1pack.

ACE1pack.jl has two purposes: (1) to import and re-export ACE1.jl, IPFitting.jl, JuLIP.jl with guaranteed version compatibility; and (2) to have several convenience wrappers for setting up the least-squares problem (ACE1.jl & JuLIP.jl) and solving it (IPFitting.jl). For full documentation see ACE1pak overview.

First import ACE1pack

using ACE1pack

First, we need to construct various parameters' dictionaries that define various aspects of fitting an ACE potential. We use various *params() functions that return these dictionaries and let us only specify mandatory and non-default parameter values.

data_param_dict = data_params(
    fname = joinpath(ACE1pack.artifact("TiAl_tutorial"), "TiAl_tutorial.xyz"),
    energy_key = "energy",
    force_key = "force",
    virial_key = "virial")
Dict{String, String} with 4 entries:
  "force_key"  => "force"
  "energy_key" => "energy"
  "fname"      => "/home/runner/.julia/artifacts/b437d7d5fac4424b8203d0afc31732…
  "virial_key" => "virial"

N.B. there is no way to sub-select data (should there be?) and here this tutorial diverges from [Fitting a TiAL potential][TiAl.md], but the the .xyz file is small enough.

r0 = 2.88
species = ["Ti", "Al"]     # symbols (:Ti, :Al) also work
2-element Vector{String}:
 "Ti"
 "Al"

basis_params of type="ace" can optionally have radial part defined.

ace_radial_params = basis_params(
    type = "radial",
    r0 = r0,
    rin = 0.6 * r0,
    rcut = 5.5)
Dict{String, Any} with 5 entries:
  "rcut" => 5.5
  "rin"  => 1.728
  "pin"  => 2
  "pcut" => 2
  "type" => "radial"

Construct ACE basis

ACE_basis_param_dict = basis_params(
    type = "ace",
    species = species,
    N = 3,
    maxdeg = 6,
    r0 = r0,
    radial = ace_radial_params)
Dict{String, Any} with 7 entries:
  "N"         => 3
  "maxdeg"    => 6
  "radial"    => Dict{String, Any}("rcut"=>5.5, "rin"=>1.728, "pin"=>2, "pcut"=…
  "degree"    => Dict{String, Any}("csp"=>1.0, "chc"=>0.0, "wL"=>1.5, "bhc"=>0.…
  "type"      => "ace"
  "transform" => Dict{String, Any}("r0"=>2.88, "type"=>"polynomial", "p"=>2)
  "species"   => ["Ti", "Al"]

and pair basis.

pair_basis_param_dict = basis_params(
    type = "pair",
    species = species,
    maxdeg = 6,
    r0 = r0,
    rcut = 7.0)
Dict{String, Any} with 8 entries:
  "rcut"      => 7.0
  "rin"       => 0.0
  "maxdeg"    => 6
  "pin"       => 0
  "pcut"      => 2
  "type"      => "pair"
  "transform" => Dict{String, Any}("r0"=>2.88, "type"=>"polynomial", "p"=>2)
  "species"   => ["Ti", "Al"]

The keys in the following dictionary are for reference, the basis function kind is defined by the type parameter. This way, it's possible to specify multiple "ACE" and/or "pair", etc basis.

basis_param_dicts = Dict(
    "pair" => pair_basis_param_dict,
    "ace" => ACE_basis_param_dict)
Dict{String, Dict{String, Any}} with 2 entries:
  "ace"  => Dict("N"=>3, "maxdeg"=>6, "radial"=>Dict{String, Any}("rcut"=>5.5, …
  "pair" => Dict("rcut"=>7.0, "rin"=>0.0, "maxdeg"=>6, "pin"=>0, "pcut"=>2, "ty…

We also need to give the "isolated atom" energies that will be subtracted from total energies before the fit.

e0 = Dict(
    "Ti" => -1586.0195,
    "Al" =>  -105.5954)
Dict{String, Float64} with 2 entries:
  "Ti" => -1586.02
  "Al" => -105.595

weights are given in a dictionary as before

weights = Dict(
        "FLD_TiAl" => Dict("E" => 30.0, "F" => 1.0 , "V" => 1.0 ),
        "TiAl_T5000" => Dict("E" => 5.0, "F" => 1.0 , "V" => 1.0 ))
Dict{String, Dict{String, Float64}} with 2 entries:
  "FLD_TiAl"   => Dict("V"=>1.0, "E"=>30.0, "F"=>1.0)
  "TiAl_T5000" => Dict("V"=>1.0, "E"=>5.0, "F"=>1.0)

The fit will be done using LSQR solver. lsqr_atol has a default value of 1e-6, so we can skip it here.

solver_param_dict = solver_params(
    type = "lsqr",
    lsqr_damp = 1e-2)

and define parameters for smoothness prior.

smoothness_prior_params = regularizer_params(
    type = "laplacian",
    rlap_scal = 3.0)        # default
Dict{String, Any} with 2 entries:
  "rlap_scal" => 3.0
  "type"      => "laplacian"

Finally, let's upt everything together. Note that exporting to pacemaker code isn't supported (should it be?).

ace_fit_params = fit_params(
    data = data_param_dict,
    basis = basis_param_dicts,
    solver = solver_param_dict,
    e0 = e0,
    weights = weights,
    P = smoothness_prior_params,
    ACE_fname = "ACE.json"  # change to `nothing` if you don't want to save the potential
)

fitted_potential, lsqinfo = ACE1pack.fit_ace(ace_fit_params)

fitted_potential is the fitted potential which can be evaluated in julia and lsqinfo contains information about the fit. The potential will also be saved to the file ACE.json which can be read in python or julia. If you want to export the potential to LAMMPS, use

ACE1pack.ExportMulti.export_ACE("TiAl_tutorial_pot.yace", fitted_potential)

This page was generated using Literate.jl.