What’s computational physics about?

Is it a subfield of theory? What about  - uhmm - “numerical experiments” (if there’s something like that)?

One could consider numerical simulation a third paradigm of science besides theory and experiment (see e.g. this link).
The term “numerical experiment”, however, would suggest a black box concept, i.e. you feed some input into the code and get numbers out, but don’t necessarily understand how and why it works and therefore not even if the output means anything.
So I would prefer the term “simulation” and make computational physics a part of theory and apply it to problems which are too complex to be solved by hand. Along this line scientific computing can aim at providing new “tools of discovery” of complex phenomena.

What do you think “scientific computing” is about and do you think it is important at all to decide which field it belongs to?

Bandstructure code from 1972

Electronic bandstructure of cubic YZn calculated with an APW code from the early 70's.

Are there still source codes of other old bandstructure programs available?


See further details about how to compile and use the code on nowadays hardware below.

In 1972, Hoffstein, Ray, and Belakhovsky have published a symmetrized augmented planewave code http://dx.doi.org/10.1016/0010-4655(72)90098-7 (Comput. Phys. Commun. 4 361 (1972)), which I've used to calculate the data shown in the figure above. Neither source nor binary must be distributed. You need access to the Computer Physics Communications Program Library to get the source code (tip: search for the authors to find the archive).

The program was written in Fortran IV (≈ Fortran 66) and the file you download is a punch code conversion containing the command line for the compilation, the source code, what is fed to standard input...
To extract the Fortran code and the input data run:

zcat acmj_v1_0.gz | tail +10 | head -1224 >sapw.f
zcat acmj_v1_0.gz | tail +1308 | head -313 >>sapw.f
zcat acmj_v1_0.gz | tail +1622 | head -506 >inph

We've skipped a few lines of the source code, which back in the days served for the calculation of a float encoded in a four character string in a special (but simple) way. The replacement of this routine for other hardware was actually worth a second publication by de Hosson in 1975 (http://dx.doi.org/10.1016/0010-4655(75)90091-0) ...

We'll take advantage of something else that was developed in 1972 - C - to implement this simple conversion:

/* conv.c */

double conv_(char *x)
{
    double sign = 1., c = 0., d;
    int i;

    for(i=0, d=1000.; i<4; i++, d/=10.) {
       switch(x[i]) {
       case '+':
           sign = 1.;
           break;
       case '-':
           sign = -1.;
           break;
       case 'A':
           c += sign * d * 0.86602540378;
           break;
       case 'B':
           c += sign * d * 0.5;
           break;
       case 'C':
           c += sign * d * 1.5;
       case 32:
           break;
       default:
           c += ((double) (x[i]-48)) * sign * d;
       }
    }

    return c;
}

Phew, what a masterpiece ;)

Compile with:

gcc -O2 -c conv.c
gfortran -O2 -fdefault-real-8 -o sapw sapw.f conv.o

The following python script can be used to calculate the dispersion from Γ to R:

from os import system,popen

e = []
for i in range(9):
  e.append({})

for s in (-2.0,-1.8,-1.5,-1.2,-1.0,-0.7,-0.5,-0.2,0.4,0.6,0.8):
  for i in range(9):
    k = float(i)/10.+0.1

    system('cat inph >inp')
    f = open('inp','a')
    print >>f,'       %g       %g       %g' % (k,k,k)
    print >>f,' -4  3  1  1  1 R 15 BANDS\n      %g     0.05   50' %s
    f.close()

    p = popen("./sapw <inp | grep EIGENVA | awk '{print $3}'")
    a = p.readline().strip()
    p.close()
    if len(a)>0:
      e[i][a] = 1

for i in range(9):
  k = float(i)/10.+0.1
  q = e[i].keys()
  q.sort()
  for a in q:
    print k, float(a)*2.0

You spin my subspace right round... - direct minimization in DFT for metallic systems

The article

Direct minimization technique for metals in density functional theory
Freysoldt, Boeck, Neugebauer
Phys. Rev. B 79, 241103(R) (2009)
http://link.aps.org/doi/10.1103/PhysRevB.79.241103

describes a method for the direct minimization of the expansion coefficients of the Kohn-Sham orbitals in a density functional theory scheme with simultaneous and consistent updates of occupation numbers and subspace rotations.
This extension to previous approaches allowed the authors to perform fast direct-minimization-DFT calculations for metallic systems (the convergence of which depends sensitively on changes in the Kohn-Sham Fermi surface during the iterations).
Since the KS orbitals have to be orthonormalized, the asymptotic complexity, however, is O(N³), i.e. the same as for self-consistent, iterative diagonalization of the KS Hamiltonian or any other KS orbital-based approach.