07 June 2013

441. Links to two good articles on windows, piracy and linux.

I don't like the idea of simply posting links to other blogs, but in this case I'll make an exception due to the quality of the articles.

These two articles from 2008 discuss the lack of adoption of linux based on the notion that windows is de facto free (gratis) to most people:
http://autotelic.com/windows_is_free
http://autotelic.com/no_really_-_windows_is_free

I find them well-written and well-reasoned.

And they address the basis of one of the more disparaging (although to some extent true) remarks by people who don't see what the fuss over linux is all about: That we can't get people to use Linux even by giving it away for free.

Anyway.From a more personal POV:

* Having lived in China I definitely agree with the idea that piracy of windows is one biggest problems when it comes to the adoption of free software in the developing world. Not only windows of course -- I could even buy SPSS, Origin and Matlab in my local computer software street stall or my local DVD store (all pirated, of course).

* Having lived in the developing world I also agree that Windows is free in the sense that you can hardly buy a computer without getting a copy of a recent windows version included (whether you want it or not).

* Working at a university -- and having worked at five in total -- I also agree that it's easy enough to get free access to most pieces of proprietary software and since the distinction between Home and Work is a bit fluid in academia, for all intents and purposes I have free access to Office, Windows, Photoshop, SPSS etc.

* And finally, having bought my first computer in as a teenager in 1993 (a 1.8 MHz 386SX, 2 Mb RAM, 28 Mb HDD -- second hand) I also grew up swapping floppies with windows (3.11 FTW!), DOS (the box had 5.0, but got DOS 6 from a class mate), and various pieces of free/shareware that we ordered via mailorder...(or bbs -- but anything over 100 kb took forever). I don't think teens of today look at things much differently from how we did back then.


06 June 2013

440. Briefly: Upgrading Arch -- consolidation of /sbin, /bin and /usr/sbin

Anyone running Arch and attempting to upgrade recently will have run into this:

(83/83) checking for file conflicts  [##########] 100%
error: failed to commit transaction (conflicting files)
filesystem: /bin exists in filesystem
filesystem: /sbin exists in filesystem
filesystem: /usr/sbin exists in filesystem
Errors occurred, no packages were upgraded.


There's a solution:
http://news.softpedia.com/news/Read-This-Article-Before-Updating-Your-Arch-Linux-358104.shtml

What I'm doing in this post is simply testing that solution. The intellectual contribution from me is very small.

1. Non-official packages
pacman -Qqo /bin /sbin /usr/sbin | pacman -Qm -
bootchart2-git 20130314-1
Not too bad.
sudo pacman -R bootchart2-git

2. Unofficial repositories
cat /etc/pacman.conf | grep '^\['
[options] [core] [extra] [community] [multilib

OK. Let's see whether there are any packages around:

paclist options | awk ' { print $1 } ' | pacman -Ql - | grep ' /s\?bin/\| /usr/sbin/'
apache /usr/sbin/ apache /usr/sbin/ab apache /usr/sbin/apachectl apache /usr/sbin/apxs [..] wpa_supplicant /usr/sbin/wpa_supplicant zvbi /usr/sbin/ zvbi /usr/sbin/zvbid
paclist core | awk ' { print $1 } ' | pacman -Ql - | grep ' /s\?bin/\| /usr/sbin/' paclist extra | awk ' { print $1 } ' | pacman -Ql - | grep ' /s\?bin/\| /usr/sbin/' paclist community | awk ' { print $1 } ' | pacman -Ql - | grep ' /s\?bin/\| /usr/sbin/' paclist multilib | awk ' { print $1 } ' | pacman -Ql - | grep ' /s\?bin/\| /usr/sbin/'

I think 'options' would qualify as an Official repository, so we are fine.
3. Orphan packages
find /bin /sbin /usr/sbin -exec pacman -Qo -- {} + >/dev/null

comes up blank. Good.

4. Upgrade
sudo pacman -Syu --ignore filesystem,bash
[..] (27/81) upgrading libgdm [################] 100% (28/81) upgrading gdm [################] 100% warning: directory permissions differ on /var/log/gdm/ filesystem: 711 package: 1770 (29/81) upgrading libx11 [..] (67/81) upgrading pacman-mirrorlist [################] 100% warning: /etc/pacman.d/mirrorlist installed as /etc/pacman.d/mirrorlist.pacnew (68/81) upgrading pm-utils [..] (80/81) upgrading wicd-gtk [################] 100% (81/81) upgrading zvbi [################] 100%

So far, so good.

sudo pacman -S bash
resolving dependencies... looking for inter-conflicts... Packages (1): bash-4.2.045-4 Total Installed Size: 3.51 MiB Net Upgrade Size: -0.10 MiB :: Proceed with installation? [Y/n] :: Proceed with installation? [Y/n] (1/1) checking keys in keyring [###################] 100% (1/1) checking package integrity [###################] 100% (1/1) loading package files [###################] 100% (1/1) checking for file conflicts [###################] 100% (1/1) checking available disk space [###################] 100% (1/1) upgrading bash
sudo pacman -Su
:: Starting full system upgrade... resolving dependencies... looking for inter-conflicts... Packages (1): filesystem-2013.05-2 Total Installed Size: 0.01 MiB Net Upgrade Size: -0.30 MiB :: Proceed with installation? [Y/n] (1/1) checking keys in keyring [###########] 100% (1/1) checking package integrity [###########] 100% (1/1) loading package files [###########] 100% (1/1) checking for file conflicts [###########] 100% (1/1) checking available disk space [###########] 100% (1/1) upgrading filesystem [###########] 100%

So far so done!

05 June 2013

439. Calculate frequencies from a hessian file from NWChem: example in Octave (matlab)

I wanted to calculate normal modes (frequencies) for specific atoms in a calculation, and so I had to write my own code.

This Octave code calculates frequencies for the first N atoms, where N is given in the input.mass file.

Background
The format that NWChem uses for the Hessian is that of a flat, triangular matrix i.e. a triangular matrix such as
1  
2 3 
4 5 6
is represented as
1
2
3
4
5
6

The Hessian is symmetric around the diagonal, so the full Hessian matrix is
1 2 4
2 3 5
4 5 6

The Hessian is independent of the masses of the atom pairs, while the frequencies are heavily dependent on the masses (isotope effects are quite visible for light elements).

To get the mass-weighted matrix we divide by the square root of the product of the masses (H * /(sqrt(m1*m2))). Note that the matrix reported in the nwchem output ("MASS-WEIGHTED NUCLEAR HESSIAN (Hartree/Bohr/Bohr/Kamu)") is multiplied by 1,000.

Once you have the mass-weighted hessian you need to calculate the eigenvalues, sort them and convert them to cm-1 using a scaling factor.

That's it.

The code:
See below for example input.mass and input.hess

%% prepare
clear;
format long

%%Calculate conversion factor from H/B/B/amu to cm-1

%% csi=299792458; %speed of light, m/s 
%% t2au=2.418884326505E-17; % seconds per a.u.
%% Better to do it by hand to avoid rounding errors:
cau=(2.99792458 * 2.418884326505)*1E-9; %c in metres per t(a.u.)

%% 1 electron (au)=9.10938291E-31 kg
%% 1 amu = 1.66053892E-27 kg
%% Better to do by hand to avoid rounding errors:
amu2au=(1.66053892/9.10938291)*1E4;% 1 amu in a.u. (via kgs)
%% For clarity
cmtom=1/100; %m per cm
%% And finally we get our scaling factor:
scaling=cmtom*(1/(2*pi*cau*sqrt(amu2au))); %( m/cm * 1/((m/au) * au) = m/cm * 1/m = 1/cm)


%%read masses
% The mass file contains the masses of the atoms
% The first line is the number of atoms in the file
% The remaining lines are the atom masses in the same order
% as the atoms are given in the nwchem input
protomasses=fopen("input.mass");
natoms=str2num(fgetl(protomasses));
for i = 1:natoms
 mass(end+1)=str2num(fgetl(protomasses));
end
fclose(protomasses);

%% Read and construct hessian from flat hessian in .hess file
%% The .hess file provided by nwchem is flat (i.e. one
%% dimensional) and is the triangular form (i.e half) of 
%% the full hessian. We use fgetl/str2num so that we can deal 
%% with instances of scientific notation in the hessian file.
%% While we"re at it we construct the mass-weighted force matrix too.
protohessian=fopen("input.hess"); 
hessian=zeros(3*natoms);
massweighted=zeros(3*natoms);

for i = 1:3*natoms
 for j=1:i
  hessian(i,j)=str2num(fgetl(protohessian));
  massweighted(i,j)=hessian(i,j)/sqrt( mass(ceil(i/3))*mass(ceil(j/3)));
 end
end

for i=1:3*natoms
 for j=1:i
  hessian(j,i)=hessian(i,j);
  massweighted(j,i)=massweighted(i,j);
 end
end

%% Diagonalize and compute frequencies in cm^{-1}
eigen=sort(eig(massweighted));
freqs=sqrt(eigen).*scaling;

%% Make imaginary frequencies negative and store them 
%% in a new array
for n=1:size(freqs,1)
 if imag(freqs(n))==0
  frequencies(end+1,1)=real(freqs(n));
 else
  frequencies(end+1,1)=-imag(freqs(n));
 end
end

%% Echo frequencies to stdout
printf("%10.4f \n",frequencies)
%% Save frequencies as well to modes.outs
outfile=fopen("normal.out","w");
fprintf(outfile,"%i \n",natoms);
fprintf(outfile,"%10.10f \n",frequencies);
fclose(outfile);
%save 'modes.out' -ascii  frequencies

input.mass (for water):
3
1.5994910D+01
1.0078250D+00
1.0078250D+00

input.hess (this one has imaginary frequencies as well):
     6.6177469151D-01
    -5.8658669668D-12
    -1.0013075598D-05
     1.0754299967D-09
     4.5060920407D-10
     3.6644723357D-01
    -3.3088202114D-01
     2.1099357839D-10
     1.6617441386D-01
     3.6163164885D-01
     2.5270659061D-12
     4.0920019206D-06
     3.2209366184D-11
     1.6382988861D-11
     8.3427731090D-07
     2.3904755566D-01
    -2.2311539742D-10
    -1.8322029567D-01
    -2.0261099118D-01
     1.1292349908D-10
     1.7796238990D-01
    -3.3088202212D-01
    -2.4469194991D-10
    -1.6617441477D-01
    -3.0749615389D-02
     1.2368245322D-10
    -3.6436594678D-02
     3.6163164980D-01
     2.5272503844D-12
     4.0920029550D-06
     3.2022391582D-11
     1.6289326095D-11
    -4.9229359909D-06
    -1.5407535297D-11
    -1.8816632580D-11
     8.3427660670D-07
    -2.3904755666D-01
    -2.2750774181D-10
    -1.8322029575D-01
     3.6436523006D-02
     1.1512005611D-10
     5.2580053385D-03
     2.0261099171D-01
     1.1238793371D-10
     1.7796238961D-01

Output:
 
  -11.0036 
   -1.6327 
    3.1676 
    3.9298 
    7.5811 
   12.2862 
 1619.0207 
 3616.0904 
 3781.1341

c.f.
 ----------------------------------------------------------------------------
 Normal Eigenvalue ||                 Infra Red Intensities
  Mode   [cm**-1]  || [atomic units] [(debye/angs)**2] [(KM/mol)] [arbitrary]
 ------ ---------- || -------------- ----------------- ---------- -----------
    1      -11.004 ||    0.426523           9.840       415.796      59.477
    2       -1.633 ||    0.000029           0.001         0.028       0.004
    3        3.168 ||    0.000003           0.000         0.003       0.000
    4        3.930 ||    0.000700           0.016         0.682       0.098
    5        7.581 ||    0.134394           3.101       131.014      18.741
    6       12.286 ||    0.000000           0.000         0.000       0.000
    7     1619.021 ||    0.070174           1.619        68.409       9.786
    8     3616.091 ||    0.004517           0.104         4.404       0.630
    9     3781.135 ||    0.009065           0.209         8.837       1.264
 ----------------------------------------------------------------------------