Skip to content
Brand Logo
  1. Home
  2. CRYSTAL

CRYSTAL

Discuss features, updates, and general use of the CRYSTAL module

56 Topics 263 Posts

Subcategories

  • SCF, Convergence, Thresholds, Density Functionals, Spin

    10 Topics
    44 Posts
    Aleksundefined

    Dear Rams,

    From your input/output files, it doesn't seem you are doing anything particularly wrong. It is more that the system is perhaps a bit more challenging. Here are a few ideas that might be of use regarding speed up and convergence:

    i) You are having quite large lattice parameters (15 Angstrom in each direction). In this case having a 4x4x4 mesh of k-points is perhaps even more than needed to converge the total energy, so you can try reducing it to a 3x3x3 or 2x2x2 mesh to reduce the computational cost and get some speed up (without compromising quality). You could probably go to 1 k-point just as well, but then you might not get a correct Fermi level sampled. Which brings me to the second point...

    ii) Do you expect your system to be metallic? If so, you can test a few values of the smearing parameter to see whether it makes a difference to the convergence behaviour. If you know you need to get a semiconductor, then you might want to try shifting the occupied states down in energy (keyword LEVSHIFT to prevent the system from getting stuck in a conducting state and causing the SCF to get stuck...

    iii) You have the DIISALLK keyword in your input file, have you tested that as necessary for the convergence of the undoped cell? It does increase the memory/storage requirements, so could be worth trying with just the default (gamma-point) DIIS...

    iv) One alternative is to turn off the DIIS convergence accelerator and revert to "classical" Fock mixing. Albeit slower in convergence than DIIS, it can stabilize the SCF procedure avoiding strong, sudden oscillation in the solution (while DIIS is turned on, the FMIXING is used only in the very first SCF step)...

    v) Coming back to the charge oscillations you observed, they do indeed start appearing very quickly in your second run.
    Here are the corresponding total energies from your OUTPUT_Run1.d12:

    CYC 0 ETOT(AU) -1.128999989618E+05 DETOT -1.13E+05 tst 0.00E+00 PX 1.00E+00 CYC 1 ETOT(AU) -1.129553711002E+05 DETOT -5.54E+01 tst 0.00E+00 PX 1.00E+00 CYC 2 ETOT(AU) -1.129699205883E+05 DETOT -1.45E+01 tst 7.60E-05 PX 9.78E-02 CYC 3 ETOT(AU) -1.123277877060E+05 DETOT 6.42E+02 tst 4.05E+00 PX 8.33E+00 CYC 4 ETOT(AU) -1.123934069238E+05 DETOT -6.56E+01 tst 1.68E-01 PX 6.66E+00 CYC 5 ETOT(AU) -1.128505608381E+05 DETOT -4.57E+02 tst 2.53E-01 PX 1.17E+00 CYC 6 ETOT(AU) -1.126570396187E+05 DETOT 1.94E+02 tst 1.41E-01 PX 1.13E+00 CYC 7 ETOT(AU) -1.129397986733E+05 DETOT -2.83E+02 tst 1.36E-01 PX 1.16E+00 CYC 8 ETOT(AU) -1.128184205493E+05 DETOT 1.21E+02 tst 9.81E-02 PX 6.95E-01 CYC 9 ETOT(AU) -1.129501102455E+05 DETOT -1.32E+02 tst 9.49E-02 PX 6.74E-01 CYC 10 ETOT(AU) -1.129223250516E+05 DETOT 2.78E+01 tst 5.80E-01 PX 6.62E-01

    And then from OUTPUT_Run2.d12

    CYC 0 ETOT(AU) -1.129223250516E+05 DETOT -1.13E+05 tst 0.00E+00 PX 1.00E+00 CYC 1 ETOT(AU) -8.145313982597E+04 DETOT 3.15E+04 tst 0.00E+00 PX 1.71E+02 CYC 2 ETOT(AU) -1.058417907360E+05 DETOT -2.44E+04 tst 2.27E+00 PX 1.71E+02 CYC 3 ETOT(AU) -7.201760490303E+04 DETOT 3.38E+04 tst 9.37E+00 PX 3.62E+02 CYC 4 ETOT(AU) -3.830339502397E+04 DETOT 3.37E+04 tst 1.36E+02 PX 7.83E+02 CYC 5 ETOT(AU) -5.612117026766E+04 DETOT -1.78E+04 tst 2.27E+01 PX 5.54E+02 CYC 6 ETOT(AU) -1.006357459536E+05 DETOT -4.45E+04 tst 3.23E+01 PX 5.67E+02 CYC 7 ETOT(AU) -1.131945182713E+05 DETOT -1.26E+04 tst 9.68E+00 PX 2.13E+02 CYC 8 ETOT(AU) -1.203304661632E+05 DETOT -7.14E+03 tst 8.51E+00 PX 1.06E+02 CYC 9 ETOT(AU) -1.148704842219E+05 DETOT 5.46E+03 tst 1.77E+01 PX 6.53E+01 CYC 10 ETOT(AU) -1.235795679474E+05 DETOT -8.71E+03 tst 2.63E+00 PX 4.40E+01

    In the first run the total energy does not oscillate, hence why the atomic charges are still somewhat as expected. In the second run, the SCF solutions starts fluctuating much more (see difference in total energies) and as a result the charges get destabilized. The SCF might not recover from that at all (even though it seems it could), so I would suggest first trying speeding up the simulation to try and get the job done in one run. Alternatively you might have to go with slower mixing and more restarts, but that is not a guarantee for a solution either. Further options include increasing the DFT grid (XXLGRID, XXXLGRID) or further increasing the last number of the integral tolerance (even though yours is pretty high already).

    Hope any of this is helpful in the end!

    Cheers,
    Aleks

  • Input Format, Pseudopotentials

    2 Topics
    15 Posts
    job314undefined

    Alessandro, the setback is temporary. What I really appreciate is that you gave CRYSTAL more life via user support than it ever had. I am very grateful for it

  • Internal Coordinates, Constraints, Convergence

    9 Topics
    39 Posts
    fbodoundefined

    Dear esmuigors
    esmuigors said in Unexpectedly cannot get the geometry after optimization: KEYWORD EXTPRT NOT ALLOWED:

    Should I only take the atoms labelled by "T" if the lattice is tetragonal, and which ones are to be taken in case of an orthorombic one? Or does "T"'and "F" stand just for "TRUE"'and "FALSE"?

    You are correct T and F stand for 'TRUE' and 'FALSE'.

    esmuigors said in Unexpectedly cannot get the geometry after optimization: KEYWORD EXTPRT NOT ALLOWED:

    I am actually talking about the standard runPcry23 script. Perhaps it should not be like that? Or should I modify this script?

    That script already contain the few line of codes I suggested you, if not feel free to modify it and add those lines.

    esmuigors said in Unexpectedly cannot get the geometry after optimization: KEYWORD EXTPRT NOT ALLOWED:

    Actually, I have now tried using only the "T"-labeled atoms and got this error:

    ERROR **** geometry **** FORMAT ERROR IN INPUT DECK

    I checked the input on top of the CS2_B1WC.pob_tzvp_rev2_gamma_onlyT.out file it seems like you are defining six atoms in the asymmetric unit, but only two atomic positions are specified in the input.

    CRYSTAL 0 0 0 64 6.34350113 5.54788046 9.71085545 6. ! Definitions of the number of atoms in the asymmetric unit 6 0.000000000000e+00 0.000000000000e+00 -5.000000000000e-01 16 2.303037287581e-17 3.119080311720e-01 1.207617820468e-01 ATOMSYMM ...

    The number of atoms in the asymmetric unit and the position specified should always match.

    I hope this helps you,

    Best,

  • Hessian, Phonons, Quasi-Harmonic Approximation, Anharmonic Force Constants

    5 Topics
    20 Posts
    spiekermannundefined

    Dear Alessandro,
    thank you for this explanation!
    Kind regards,
    Georg

  • Harmonic and Anharmonic Vibrational Spectra, Born Tensor, Raman Activities, Phonon Density-of-States

    12 Topics
    57 Posts
    PeterRemotoundefined

    aerba I've noticed for my systems, these SCF convergence issues appear when I try to compute for intensities in a restart job. Then when removing intensities, they converge fine and progress forward.

    Is it possible to finish the frequency calculations first and then compute for the Raman intensities afterwards? Would that be a work around? Maybe with RAMANREA - because I have the TENS_RAMAN.DAT, but it just gets stuck at some point during HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH
    FORCE CONSTANT MATRIX - NUMERICAL ESTIMATE
    HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH
    and doesn't progress with restarts

  • Strain, Elastic Tensor, Seismic Waves Velocities

    2 Topics
    5 Posts
    aerbaundefined

    Hi,

    Direct piezoelectric constants of a 3D lattice in CRYSTAL are defined and computed as:
    $$
    e_{ci}^{3D} = \left( \frac{\partial P_c}{\partial \eta_i}\right) = \frac{1}{V}\left( \frac{\partial^2 E}{\partial E_c\partial \eta_i}\right)
    $$
    that is as first derivatives of Cartesian components of the polarization (c=x,y,z) with respect to strain components, or, equivalently as second derivatives of the energy density (V is the volume of the 3D lattice cell) with respect to Cartesian components of an electric field \(E_c\) and strain components, where the strain \( \eta \) is dimensionless and thus the direct piezoelectric constants have units of \( \textup{charge/length}^2 \).

    For 1D and 2D periodic lattices, as the volume (V) is not uniquely defined (or not defined at all in some cases), one may divide by the length \(l \) and area \( A\) of the lattice cell instead:
    $$
    e_{ci}^{1D} = \frac{1}{l}\left( \frac{\partial^2 E}{\partial E_c\partial \eta_i}\right) \quad \textup{and} \quad e_{ci}^{2D} = \frac{1}{A}\left( \frac{\partial^2 E}{\partial E_c\partial \eta_i}\right)
    $$
    that would thus be expressed in units of \( \textup{charge} \) or \( \textup{charge/length} \) for 1D and 2D lattices, respectively.

    However, in CRYSTAL for 1D and 2D lattices we do not divide by \(l \) or \( A\) , and just define and compute the piezoelectric constants as:
    $$
    e_{ci}^\textup{1D and 2D} = \left( \frac{\partial^2 E}{\partial E_c\partial \eta_i}\right)
    $$
    with units of \( \textup{charge}\cdot\textup{length} \).

    Yes, these constants are physically meaningful for 1D and 2D systems. For a 2D monolayer system, for instance, depending on what you need to compare with, you can do one of two things:

    keep them as they are printed in the CRYSTAL output (units of \( \textup{charge}\cdot\textup{length} \))

    divide the values you get in the CRYSTAL output by the area of the 2D cell (and thus express them in units of \( \textup{charge/length} \))

    I would not divide by a volume because I would not know the physical meaning of the volume of a 2D monolayer system.

    Hope this helps,

  • Stress Tensor

    2 Topics
    24 Posts
    aerbaundefined

    For the purpose of finding the minimum energy structure to then do Raman calculations, it is.

    EOS gives you much more than that of course: the p(V) or, equivalently, V(p) relation (i.e. structure as a function of pressure), the bulk modulus K(p), and allows to compute the enthalpy H(p).

  • Two-Component Spinors, Non-Collinear Magnetisation, Spin and Particle Currents

    4 Topics
    18 Posts
    LUXiuundefined

    Thanks Jacques! These are very informative.
    The information for specific basis sets is helpful, I will try them.

  • Input Format, Symmetry, Manipulation, Slabs, Nanotubes, Fullerenes, Helices

    5 Topics
    20 Posts
    job314undefined

    Hi Giu, it worked!

  • Electric Field, Polarizability, Dielectric Tensor, Hyper-Polarizabilities

    4 Topics
    18 Posts
    esmuigorsundefined

    To clarify better, does the symmetrization correspond to the frequency symmetry (which determines the equal components, e.g., xxyy = xyxy = xyyx, etc. for the static case?
    Or is this about the geometrical symmetry (and possible reorientation of the molecule to have the dipole moment along an axis)?

  • Questions that do not fit in other categories

    1 Topics
    3 Posts
    job314undefined

    Thank you very much, that worked (but not restart)