## Physical Chemistry: The Simple Hückel Method (Part II)

In the previous article, we have learn how to set up the Hückel determinant for an aromatic molecule based on the topology of the pi-bonds. In this second part, we are going to learn how to calculate from the determinantal equation both the eigenvalues and the eigenvectors, corresponding to the orbital energy and orbital functions of the molecular system.

EIGENVALUES FROM THE HÜCKEL DETERMINANT

The calculation of the determinat gives a so-called characteristic equation. Namely a polynomial equation whose roots ($x_i$) are the eigenvalues of the system. As described in the previous article, the eigenvalues are related to the energy of the system by the relation $x=\frac{\alpha-E}{\beta}.$ Let work out the determinant for the allyl molecule. The determinant is given by

that is solved as

The polynomial $(x^3-2x)=0$ has three root that can be easily found by rearranging it in $x(x^2-2)=0$ giving $x_1=0$ and $x_{2,3}=\pm\sqrt{2}.$

We can now use these value to calculate the energy of the Hückel orbitals using the relation $x=\frac{\alpha-E}{\beta}.$ Therefore for

$x_1=0=\frac{\alpha-E_1}{\beta}$, we obtain $E_1=\alpha$, and

$x_{2,3}=\pm\sqrt{2}=\frac{\alpha-E_{2,3}}{\beta}$, we obtain $E_{2,3}=\alpha \mp\sqrt{2}\beta$.

## Calculation of Molecular Orbitals Coefficients

We are going to calculate the wavefunctions of the Hückel’s orbitals using the eigenvalues obtained from the characteristic equation. The coefficients of the molecular orbitals are calculated by solving the system of equations and the normalization condition given in the following slides.

For $x_1 =0$, we get two coefficients equal and different from zero, and using the normalization condition $2c^2=1$, the value of $c_1=c_3==1/\sqrt{2}$ is obtained.

For $x_{2,3} =\pm \sqrt{2}$, the following values are obtained.

Finally, the eigenfunctions and their graphical representation are reported in the following slide.

For linear alkenes the eigenvalues and eigenvectors can be calculated using simple formulas reported in the following slide [1].

Here, the value for the allyl molecule are recalculated again using the formulas in the previous slide. A simple program, scritten in awk language is also listed in the appendix.

The numerical representation obtained by the formula are given as follows

 Number of atoms                    : 3
Number of electrons                : 2
Number of double occupied orbitals : 0
Number of single occupied orbitals : 0
Orbital Energies:
E1=alpha   -1.414beta
E2=alpha   -0.000beta
E3=alpha   +1.414beta
Total Electronic Energy:
Epi=   2.000alpha   +2.828beta
Wavefunctions:
Psi(1)=   +0.500 X(1)  +0.707 X(2)  +0.500 X(3)
Psi(2)=   +0.707 X(1)  +0.000 X(2)  -0.707 X(3)
Psi(3)=   +0.500 X(1)  -0.707 X(2)  +0.500 X(3)

The linear chain next to the allyl molecule is the 1,3-butadiene, composed of 4 carbon atoms and two conjugated double bonds—the Hückel topological matrix is shown in the next slide. The determinant gives a 4th-degree characteristic equation whose roots can easily be found by variable substitution, as shown in the slide.

By running the program in the appendix, we can find the value of the energy and of the coefficients of the wavefunction

Energy:
E1=alpha   -1.618beta
Number of atoms                    : 4
Number of electrons                : 4
Number of double occupied orbitals : 0
Number of single occupied orbitals : 0
Orbital Energies:
E1=alpha   -1.618beta
E2=alpha   -0.618beta
E3=alpha   +0.618beta
E4=alpha   +1.618beta
Total Electronic Energy:
Epi=   4.000alpha   +4.472beta
Wavefunctions:
Psi(1)=   +0.372 X(1)  +0.602 X(2)  +0.602 X(3)  +0.372 X(4)
Psi(2)=   +0.602 X(1)  +0.372 X(2)  -0.372 X(3)  -0.602 X(4)
Psi(3)=   +0.602 X(1)  -0.372 X(2)  -0.372 X(3)  +0.602 X(4)
Psi(4)=   +0.372 X(1)  -0.602 X(2)  +0.602 X(3)  -0.372 X(4)

A rotational barrier separates the s-cis- and s-trans conformers of the 1,3-Butadiene. We can try to make an estimation of the transition barrier using the HMO (see Ref. [2]). We can assume that the two double bonds are not conjugated during the rotation as for the allyl molecules. Therefore, the energy of the transition state is given by the sum of the two isolated double bonds’ energy. Being a total of four electrons in the $\pi$ orbitals, we can calculate the total energy as $E_{\pi}=4\alpha+4\beta$. In the cis and trans butadiene, the HMO energy is the same as no difference in the molecule’s topology and is equal to $E_{\pi}=4\alpha+4.472\beta$. Therefore, the transition state’s energy difference is equal to $0.472\beta,$ that correspond approximately to $8~kcal/mol$. The experimental value of the transition barrier has a value of $4~kcal/mol.$ This result shows that the limit of the HMO approximation largely overestimated the true value of 100%.

Even though the HMO method is simple and approximated, it provided crucial results to understand the organic reaction mechanism. In the next article, we will calculate the HMO for cyclic conjugated molecules and discuss other results of this method.

If you have found interesting and useful my article, do not forget to press “Like it” and subscribe for updates on new ones!

## REFERENCES

1. J. P. Lowe. Quantum Chemistry. 1993, Academic Press.
2. F.A. Carroll. Structure and Mechanism,1998, BROOKS/COLE Publishing Company.

## APPENDIX

Program in awk language to calculate the energy and coefficients of the Hückle orbitals. You may change the variable for the number of atoms (n) and the number of electron (ne) to your convenience before running the program using the command:

gawk -f HMOLinearMol.awk

#======================================================================
#
#  NAME:  HMOLinearMol.awk
#
#======================================================================
#  DESCRIPTION: Calculate the energy level and coefficients of the
#               wavefunctions
#======================================================================
# Copyright (C): 2021 Danilo Roccatano
#======================================================================
#

BEGIN {
n=4 # number of atoms in the chain
ne = 5  # number of electrons in the chain
pi=4*atan2(1,1)
#
printf "Number of atoms                    : %d \n",n
printf "Number of electrons                : %d \n",ne
printf "Number of double occupied orbitals : %d \n",npe
printf "Number of single occupied orbitals : %d \n",nue

# Calculate the eigenvalues
tEa=0
tEb=0
printf "Orbital Energies:\n"
for (k=1;k<=n;k++) {
x=-2*cos(k*pi/(n+1))
printf " E%d=alpha %+8.3fbeta\n",k,x
Eb[k]=x
}
#
# Calculate the total energy of the electronic system
#
npe=int(ne/2)
nue=ne-2*npe
for (i=1;i<=npe+nue;i++) {
fa=1
if (i<=npe)  fa=2
tEa+=fa
tEb+=fa*Eb[n-(i-1)]
}
printf "Total Electronic Energy:\n"
printf " Epi=%8.3falpha %+8.3fbeta\n",tEa,tEb
#
# Calculate the wavefunction coefficients
#
printf "Wavefunctions:\n"
for (i=1;i<=n;i++) {
printf " Psi(%d)= ",i
for (k=1;k<=n;k++) {
c[k]=sqrt(2/(n+1))*sin(k*i*pi/(n+1))
printf "%+8.3f X(%d)",c[k],k
}
printf"\n"
}
}


## A Practical Introduction to the C Language for Computational Chemistry. Part 3

Sphere. From Space, from Space, Sir: whence else?

Square. Pardon me, my Lord, but is not your Lordship already in Space, your Lordship and his humble servant, even at this moment?

Sphere. Pooh! what do you know of Space? Define Space.

Square. Space, my Lord, is height and breadth indefinitely prolonged.

Sphere. Exactly: you see you do not even know what Space is. You think it is of Two Dimensions only; but I have come to announce to you a Third — height, breadth, and length.

Square. Your Lordship is pleased to be merry. We also speak of length and height, or breadth and thickness, thus denoting Two Dimensions by four names.

Sphere. But I mean not only three names, but Three Dimensions.

Flatland: A romance of many dimensions by Edwin A. Abbott

In the part 2 of this tutorial, we have learnt how to use arrays and how to read atomic coordinates from a file. In th eappendix, youc an find an example of solution of the exercises given in the previous tutorial.

In this third part, we are going to learn how to generate three-dimensional coordinated of atoms in a cubic crystal lattice and how to calculate non-bonded molecular potential and the force acting among them.

For this purpose, we need to modify the subroutines for the calculation of one dimnsional Lennard-Jones potentials and their forces developed in the previous tutorial to account to use directly the atomic coordinate of particles located in a tridimensional space.

The Cartesian distance of two points in a Euclidean space is given by:

$r_{ij}=\sqrt{ (x_i-x_j)^2+(y_i-y_j)^2+(z_i-z_i)^2}$

where $(x_i,y_i,z_i)$ and $(x_j,y_j,z_j)$ are the coordinates of the two point particles. This expression will be used to modify the LJ energy potential and the components of the atomic forces.  The latter are derived taking in account the vectorial nature of the forces in the following components:

$f(x_i)=(x_i-x_j)\frac{24\epsilon}{r_{ij}^2} \left[ 2\left(\frac{\sigma}{r_{ij}}\right)^{12} - \left(\frac{\sigma}{r_{ij}}\right)^{6} \right]$

$f(y_i)=(y_i-y_j)\frac{24\epsilon}{r_{ij}^2} \left[ 2\left(\frac{\sigma}{r_{ij}}\right)^{12} - \left(\frac{\sigma}{r_{ij}}\right)^{6} \right]$

$f(z_i)=(z_i-z_j)\frac{24\epsilon}{r_{ij}^2} \left[ 2\left(\frac{\sigma}{r_{ij}}\right)^{12} - \left(\frac{\sigma}{r_{ij}}\right)^{6} \right]$

these expressions give the components of forces acting on the atom i. The force on the atom is given by $\vec{f_j}=-\vec{f_i}$, by the Netwon’s third law:

Lex. III. Actioni contrariam semper et equalem esse reactionem: sive corporum duorum actiones in se mutuo semper esse aequales et in partes contrarias dirigi.

Isaac Newton. Philosophiae Naturalis Principia Mathematica. London, 1686

Finally, the square root function in C is calculated by the sqrt().

EXERCISE 2: modify the program in the part 2 of this tutorial to include the possibility to read a file containing the coordinated in xyz format of two LJ particles. Using those coordinates, calculate the potential and forces acting between the two particles. Change the value of the coordinates and report the same value for different distances.

Atoms in molecules are surrounded by an electron density that depend on by their tendency to delocalize their electron (these properties is usually referred to as electronegativity). The electron distribution produces and electric field around the molecule that can be matched by assigning partial electric point charges to each atom (I have explained in another article how this is accomplished). Therefore, we can add to the expression of the LJ potential a Coulomb potential term that describe the electrostatic pair interactions of charged atoms.

$V(r_{ij})_{coul} = \frac{q_i q_j}{4\epsilon_0\epsilon_r r_{ij}}$

The force components for the Coulomb term on the particle i interacting with the particle j are given by the expressions

$f(x_i)_{coul}=\frac{q_i q_j}{4\pi\epsilon_0 \epsilon_r r_{ij}^2} \frac{(x_i-x_j)}{r_{ij}}$

$f(y_i)_{coul}=\frac{q_i q_j}{4\pi\epsilon_0\epsilon_r r_{ij}^2} \frac{(y_i-y_j)}{r_{ij}}$

$f(z_i)_{coul}=\frac{q_i q_j}{4\pi\epsilon_0\epsilon_r r_{ij}^2} \frac{(z_i-z_j)}{r_{ij}}$

The $q_i and q_j$ are the partial charges and $\epsilon_0 and \epsilon_r$ the dielectrioc constant of the vacuum and of the medium. The force on the atom is given by the Netwon’s third law. To extend the calculation to many-particles systems, we need to perform a double sum is used to calculate all the pair-wise interactions as follows

\begin{aligned}V(r_1,r_2,\cdots,r_n) &= \\ &+ \sum_i \sum_{j>i} 4\epsilon_{ij}\left[\left(\frac{\sigma}{r_{ij}}\right)^{12}-\left(\frac{\sigma}{r_{ij}}\right)^{6}\right]\\ &+ \sum_i \sum_{j>i} \left(\frac{q_i q_j}{4\epsilon_0\epsilon_r r_{ij}} \right) \end{aligned}

The index i>j means that in the second sum we have to exclude all the pair already accounted by the first index. If we have four particles for example we need to calculate the following 4*(4-1)/2=6 pairs

i:   1    1    1   2    2     3

j:   2     3   4   3    4     4

when the index i=2 the index j start from 3 since 2-1 has been calculated and 2-2 does not make sense.

To program this double summation, we need to use a two nested loops with two indexes (i and j). The following pseudo-code shows how to make the double loops. The potential and forces are updated using the unary operators += and they need to be initialized to zero. The results will be printed at the end of the calculation using a loop on the number of atoms (na).

Vt=0; /* total potential */
for (i=0;i<na;i++) {
F[i]=0.0;
V[i] =0.0; } /* initialize the forces */
for (i=0;i<na;i++) {
for (j=i+1; j< na; j++) {

calculate the distance using xx[i],xx[j],yy[i],yy[j],zz[i],zz[j]
calculate the LJ potential and forces

V+= V[i] + ….
Vt+= V[i] + ….
F[i]+= F[i] + …
}
}
/* Add here the printout of the potential for each atom and the total potential */



EXERCISE 3: implement the double loop to calculate the potential energy and forces for the coordinates generated by the program from EXERCISE 1.

### CUBIC CRYSTAL LATTICE COORDINATES GENERATION

Atomic coordinates for simple regular systems, such as a cubic crystal lattice, can be easily generated using nested for-loops. Let start with the example of a simple one-dimensional lattice, where the atoms are lined up along the x-axis. This means that coordinates of these atoms along the other two directions are null. The coordinate of a given atom along the x direction are generated by multiply the lattice index for the lattice size.  The following piece of code show how to implement the algorithm using for-loops. The variable xaxis is the lattice space in Angstrom units and nx is the number of lattice points along the x direction.

xaxis = 3.0; /* this value can be given in input */
nx = 5;
xx[0] = 0;
for (i=1;i<=nx;i++) {
xx[i] = i*xaxis;
}


We can extend the method for a 2 D lattice using a nested for-loop as in this code snippet:

xaxis = 3.;
yaxis = 3.;
nx=5;
ny=5;
xx[0] = 0.;
yy[0] = 0.;
for (i=0;i<nx;i++) {
xx[i] = i*xaxis;
for (j=0;j<ny;j++) {
yy[j] = j*yaxis;
}
}


Now putting together what we have learnt from the last tutorial, we can write a code that generate the coordinates of a squared two-dimensional lattice and save them in xyz format. The code is reported in the appendix. The coordinated in xyz format (file name out.xyz) can be visualized them using the program rasmol. If you run from the comand line type rasmol –xyz  2dlattice.xyz.

With nx=10, and ny=1 you should get the the 1D-dimensional lattice that visualized with rasmol (background white, spacefill) will look like the following picture

With nx=10, and ny=10 you should get the the 2D-dimensional lattice that visualized with rasmol (background white, ball&stick) will look like the following picture

EXERCISE 1. Extend the 2D lattice program to calculate coordinates of atoms in a 3D lattice for a given number of lattice point along the 3 axis (nx, ny, nz) and for given lattice sizes (xaxis, yaxis, zaxis). You may want also to add the possibility to input these parameters from the keyboard to make the program more user-friendly.

## APPENDIX

TABLE I: LENNARD-JONES PARAMETERS FOR NOBLE GASES

### PROGRAM I: Two dimensional cubic lattice generator

/*
PROGRAM NAME : 2D reticulum generator
DESCRIPTION  : the program generate one or two dimensional
lattices having nx,ny alternating atoms in the x
and y direction, respectively. The atoms are separated
ati equal distances xsize,ysize
number of atoms in each dimension to be set
AUTHOR :  Danilo Roccatano
DATE   : BREMEN, 2007

*/
#include <stdio.h>
#include <math.h>
#include <string.h>
#include <unistd.h>

int  main() {
char sym1[]="Na";
char sym2[]="Cl";
char csym[3];
float xaxis,yaxis;
int nx,ny;
int nn,n,nsel;
float xx[10];
float yy[10];
float zz;

/*
* Initialize variables
*/

xaxis=2.82;
yaxis=2.82;
nx=5;
ny=5;
zz=0.0;

/*Opening output file*/
FILE *fout;
if (!(fout=fopen("2dlattice.xyz","w"))){
printf ("Error opening file %s\n", "2dlattice.xyz");
}
n=nx*ny;
fprintf(fout,"%i\n",n);
fprintf(fout,"%i %s%s-atoms in two dimensional rectangular lattice\n",n,sym1,sym2);
nn=0;
/*Computing and writing the lattice coordinates*/
for(int i=0;i<nx;i++){
xx[i]=i*xaxis;
for(int j=0;j<ny;j++){
yy[j]=j*yaxis;
nn+=1;
if (nn%2) strcpy(csym, sym1);
else strcpy(csym,sym2);
fprintf(fout,"%2s %8.3f %8.3f %8.3f %5d\n", csym, xx[i], yy[j], zz,nn);
}
}
printf("Generated a rectangular lattice of %i %s%s-atoms. \n",n,sym1,sym2);
fclose(fout);
}



## AN EXAMPLE OF SOLUTION FOR THE PART 2

This is an example of solution of the exercises for the last practical. The program read the parameters from the file “ljparam.dat” containing the following data.

LJ parameters
3
Atom_name   sigma(A)       epsilon(K)
Ar           3.336000      141.200000
Kr           3.575000      191.400000
Xe           3.924000      257.400000 

The program ask to choose two atoms types, then it calculates the combined parameters and ask for the maximum value of distance to calculate the value of the potential and the force. The values are stored in the file “Results.out“, and they can be plotted with xmagrace (see the graph at the end).

#include <stdio.h>
#include <math.h>
#include <string.h>
int main() {
char line[160],title[160];
char aty[10];
char yn[2];
int k=1;
int ia1,ia2;
int np; // number of parametes in the library
float V,s,s1,s2,e,e1,e2,F;
float x,rmin,xmin,xmax,xinc;
float C6,C12;
float s6,ix6,ix;
float sig[10];
float eps[10];
float kb=8.3145e-3; /*Boltzmann constant*Avogadro Number  */
FILE *fd;
FILE *fout;

/*
* READ LJ PARAMETERS FROM A LIBRARY FILE
*/

printf("\nREADING THE LJ PARAMETER FILE \n");
if (!(fd=fopen("ljparam.dat","r"))){
printf("\nError opening file %s\n", "ljparam.dat");
}

fgets (title, 80, fd);
fgets (line, 80, fd);
sscanf (line,"%d",&np);
fgets (line, 80, fd);
printf (" %s  \n",line);
while (fgets (line, 80, fd) != NULL){
if (line[0] != '\n' && strlen(line) > 10) {
sscanf (line,"%s%f%f",&aty[2*k],&sig[k],&eps[k]);
printf ("%i %s  %8.3f   %8.3f\n",k, &aty[2*k],sig[k],eps[k]);
k++;
}
}
fclose(fd);
/*
*   OPEN THE OUTPUT FILE
*/

if (!(fout=fopen("Results.out","w"))){
printf ("Error opening file %s\n", "Results.out");
}

/*Input for the atom types*/
printf("\n");
printf("Enter the atomtype for particle 1: ");
scanf("%i",&ia1);
printf("Enter the atomtype for particle 2: ");
scanf("%i",&ia2);
printf("\n");

/*Berthelot-Lorenz mixing, etc.*/

s=0.5*(sig[ia1]+sig[ia2])/10.;
e=kb*sqrt(eps[ia1]*eps[ia2]);
s6=pow(s,6);
C6=4*e*s6;
C12=4*e*s6*s6;
rmin=pow(2.,1./6.)*s;

/*Output the results on the screen*/
printf("Mixed Lennard-Jones parameters:\n");
printf(" Epsilon:  %f kJ/mol \n",e);
printf(" R min  :  %f nm     \n",rmin);
printf(" Sigma  :  %f nm\n",s);
printf(" C6     :  %f kJ*nm^6/mol \n",C6);
printf(" C12    :  %f kJ*nm^12/mol \n\n",C12);
/* Estimate the minimum distance to calculate the LJ function */
xmin=s-0.025;
ix=1./xmin;
ix6=1./pow(xmin,6.);
V=C12*ix6*ix6-C6*ix6;
printf(" Value of Potential V(xmin=%f nm):  %f kJ/mol \n",xmin,V);

/*Output the results in the file  */
fprintf(fout,"#Lennard-Jones parameters:\n\n");
fprintf(fout,"# Atom Type 1  :  %c%c \n",aty[2*ia1],aty[2*ia1 + 1]);
fprintf(fout,"#     Sigma    :  %f nm/10\n",sig[ia1]);
fprintf(fout,"#   Epsilon    :  %f K\n\n",eps[ia1]);
fprintf(fout,"# Atom Type 2  :  %c%c \n",aty[2*ia2],aty[2*ia2 + 1]);
fprintf(fout,"#     Sigma    :  %f nm/10\n",sig[ia2]);
fprintf(fout,"#   Epsilon    :  %f K\n\n",eps[ia2]);
fprintf(fout,"# Mixed Sigma  :  %f nm\n",s);
fprintf(fout,"# Mixed Epsilon:  %f kJ/mol \n",e);
fprintf(fout,"#    C6        :  %f kJ*nm^6/mol \n",C6);
fprintf(fout,"#   C12        :  %f kJ*nm^12/mol \n\n",C12);
fprintf(fout,"# Distance [nm] Potential [kJ/mol] Force [nN]\n");

/*Input for the value of the maximum distance range*/
do {
printf("\nEnter the value of the maximum distance distance (in nm): ");
scanf("%f",&xmax);
printf("\nEnter the distance increment (in nm): ");
scanf("%f",&xinc);
printf("Finish?(y/n)\n");
scanf("%s",yn);
} while (strcmp(yn,"y"));

for (x=xmin;x<=xmax;x+=xinc){
/*Calculate the LJ potential and the force between xmin and xmax  x*/
ix=1./x;
ix6=1./pow(x,6.);
V=C12*ix6*ix6-C6*ix6;
F=(12.*ix*C12*ix6*ix6-6.*ix*C6*ix6)/602.2; /*to obtain [nN]*/
fprintf(fout," %f    %f    %f\n", x,V,F);
}
fclose(fout);
return 0;
}



## A Practical Introduction to the C Language for Computational Chemistry. Part 2

In the first part of this introduction to C language, we have learnt the basic of the C language by writing simple programs for the calculation of the non-bonded interaction between two particles at variable distances. Some solutions to the first part exercises are reported in the appendix of this article.

In this second tutorial, we will learn how to use arrays data types and how to load them with a set of data read from a file. We will also use these data to perform numerical calculations and write results in output files.

Arrays and Pointers Datatypes

The program that calculate the energy of interaction between two particle doe not take in account the actual position in space of the two particle but only their distance. If we want to study the dynamics of a system composed by multiple atoms in a tridimensional space, it is way more convenient to represent the and calculate their interactions by using the coordinates directly to evaluate the distances.

## The Dandelion (Taraxacum Officinalis) and OpenCV

The dandelion’s pallid tube
Astonishes the grass,
And winter instantly becomes
An Infinite Alas —

The tube uplifts a signal Bud
And then a shouting Flower, —
The Proclamation of the Suns
That septulture is o’er.

– Emily Dickinson

The yellow flowers and the delicate and beautiful florescence of Dandelion catch the attention of both romantic and curious souls. The aerial consistency of the fine silk decorated seeds that glance to the sunlight as crystalline material became the favourite subject of inspired photographer and the toy of amused children. Beside the grace of its forms, other interesting and curious secret is hidden in its phloem fluids. In fact, if you cut one of the stems of the plant a milky sticky liquid will flow out of the wound resection. This latex is going to polymerize at 30-35 oC in few minutes in a yellow-brown quite solid mass. Around the year 1982, I have annotated this observation but I could not find in my later notes further follow-ups study on the topics. It was a casual observation but I didn’t know at that time that this latex is indeed very useful. A variety of the Taraxacum (Taraxacum koksaghyz, Russian Dandelion) was used in Russian and American to produce a replacement of the natural rubber from Brazil during the WWII that was in shortage because of the war. Now days, many researches are in progress to exploit the lattice of Taraxacum and also Taraxacum brevicorniculatumas a convenient replacement of the rubber plant lattice. A recent study has shown the presence of rubber particles in the lattice of these plant in 32% proportion composed prevalently by poly(cis-1,4-isoprene) at >95% of purity (www.biomedcentral.com/1471-2091/11/11). The brownish lattice condensate that, as I reported in my note, forms after exposing for several minutes the latex to the air, is caused by the presence of the polyphenol oxidase (PPO) enzyme that produces the fast coagulation of the latex by catalysing the oxidation of polyphenols. Genetic engineer the plant, it is possible to reduce the amount of PPO in the latex making economically competitive the use of this resource for the production of latex.

## Chimica Fisica: La Termodinamica, La Meravigliosa Cattedrale Della Scienza. Parte I.

Indipendentemente dai motivi del culto, le antiche cattedrali invitano ad un’ammirata contemplazione, ispirano rispetto e quiete. Anche il visitatore piu’ disinvolto non si esimere dal moderare la voce, non insiste in argomenti futili: delle navate, l’eco delle sue stesse parole sembra destare insolite suggestioni. L’impegno di generazioni di architetti e di artigiani e’ stato dimenticato, le loro impalcature sono state rimosse ormai da lungo tempo, i loro stessi errori sono stati cancellati dai secoli. Il monumento che essi crearono, ora compiuto e perfetto, ci appare come la testimonianza di un disegno sopraordinario. Se evochiamo in noi il ricordo di un cantiere in attività, con il rumore ritmato dei martelli, le voci ed i gesti degli operai, l’odore stantio del legno e di tabacco, alle splendide structure che ora ammiriamo non possiamo attribuire altro significato che quello di essere il frutto di un ordine imposto alla mera fatica umana.

Anche la scienza ha i suoi templi, costruiti con gli sforzi di pochi architetti e di molto operai, e di fronte ad essi proviamo lo stesso sentimento. Anche in questi templi l’atmosfera e’ solenne, e forse lo e’ a tal punto da condizionare l’espressione stessa del pensiero scientifico, che una lunga tradizione vuole assai severo e formale.

G.N, Lewis, M. Randall -“Thermodinamica”, Leonardo Edizioni Scientifiche, Roma (1971).

## INTRODUZIONE

Gilbert Newton Lewis e Merle Randall nella introduzione alla prima edizione del loro autorevole testo di termodinamica chimica descrivono  la termodinamica come la cattedrale della scienza. La loro non è solo una concessione poetica di un’epoca ancora permeata dal romanticismo scientifico, ma una meravigliosa analogia per questa fondamentale disciplina della scienza che più di ogni altra contiene  le leggi arcane che goverano il nostro Universo e il suo destino.

In questa serie di articoli riporto alcuni appunti su argomenti vari di termodinamica chimica che possono essere utili come riferimento o come materiale didattico.

## Programming in Awk Language. LiStLA: Little Statistics Library in Awk. Part I.

In the following previous Awk programming articles

The AWK Programming Language

Awk Programming II: Life in a Shell

Awk Programming III: the One-Dimensional Cellular Automaton

I have given a short introduction to this very useful Unix program by also showing two example of elaborate applications. In this fourth article of the series, I am going to show a little library of functions that can be used for basic statistical analysis of data sets. I have written (and rewritten) many of these function but I have spent little time to collect them in a library that can be used by other user. So this article give me the motivation to achieve this target. I didn’t extensively test the library so I am realeasing it as alpha version. If you spot errors or you have improved it then please just send me your modified code!

We start with a function that can be used to read data from a text file (ascii format). A good data reader should be able to read common data format such as comma separated (cvs) or space separated data files. It should also be able to spik blank lines or lines starting with special characters. It would be also handy to select the columns that need to be read and also check and skip lines with inconsistent data sets (missing data or NaNs). This is what exacty work the function ReadData() given in the Appendix. But shall we see it more in details.

The function read the data from a file with name provided in the variable filename. The program skips all empty record, those starting with one of the characters contained in the regular expression skipchar. For example, a regular expressions such as skipchr=”@|#|;” skips the occurrence of the characters “at” or “hash” or semicolomn. The variable warn is used to check the behavior of the program if alphabetic characters or NaN or INF values are present in the data. If the variable is set to 0, the function gives a warning without stop the program, if set to 1 then the function terminate the program after the first warning.

The field separator is specified in fsep and it is used to set the awk internal variable FS and define the separator between data. The variable can be assigned with single character such as fsep=” “ or fsep=”,” or ESC codes such as fsep=FS=”\t” for tab-delimited.

The column in the data record can be read in two ways by set the element zero of the array range[]. For range[0]=0, a adjoint range of data is specified by setting the first element is at range[1] the last one in range[2]. For range[0]=1, the first element in range[1] is the number of data to read followed by the specific field in the record where the data is located.

## Molekulare Maschinen: Die Coronavirus SARS-CoV-2 Bedrohung, Teil I.

Was Freunde mit und für uns tun, ist auch ein Erlebtes; denn es stärkt und fördert unsere Persönlichkeit. Was Feinde gegen uns unternehmen, erleben wir nicht, wir erfahren’s nur, lehnen’s ab und schützen uns dagegen wie gegen Frost, Sturm, Regen und Schloßenwetter oder sonst äußere Übel, die zu erwarten sind.

Johann Wolfgang von Goethe (1749-1832), Maximen und Reflexionen. Aphorismen und Aufzeichnungen.

Ein Virus ist Leben in der einfachsten Form. Es ist die minimalistische Reduktion eines Organismus auf seine wesentlichen Funktionselemente. Noch pragmatischer ist ein Virus ein Behälter mit genetischem Code mit einem effizienten molekularen Mechanismus, der es ihm ermöglicht, in eine Wirtszelle eines Organismus einzudringen, der sich selbstständig reproduzieren kann. Als molekulare Maschine kann ein Virus der Form und der zerstörerischen Kraft des Todessterns in der Star-Wars-Saga ähneln. Daher ist es eine Art molekulare Maschine, die wir absolut nicht in uns haben wollen!

Wie der große Goethe sagt, ist der Feind Teil unserer Erfahrung und wir müssen ihn jagen und uns tatsächlich vor anderen möglichen Feinden schützen. Dieser epische Naturkrieg veranlasste mich, diesen Blog zu starten, in dem ich mitteilen werde, was ich über diese gefährliche molekulare Maschine lerne.

## Le Macchine Molecolari: La minaccia del Coronavirus SARS-CoV-2. Parte I

Difficilmente è vinto colui che sa conoscere le forze sue e quelle del nemico.

Nicollò Machiavelli in Dell’arte della guerra (1519-1520)

Un virus è la vita nella forma più semplice. È la riduzione minimalista di un organismo ai suoi elementi essenziali di funzionalità. Più pragmaticamente, un virus è un contenitore di codice genetico dotato di un efficiente meccanismo molecolare che gli consente d’invadere una cellula ospite di un organismo capace di riprodursi autonomamente. Come macchina molecolare, un virus può assomigliare nella forma e potere distruttivo, alla Morte Nera della saga di Star Wars. Pertanto, è un tipo di macchina molecolare che non vogliamo assolutamente avere dentro di noi!

La diffusione del coronavirus SARS-CoV-2 (COVID-19) ha prodotto una nuova pandemia, ovvero una infezione causata da un agente patogeno che colpisce l’intera popolazione di una specie vivente, in questo caso quella umana. Questa situazione di emergenza globale è il risultato di una competizione naturale tra specie viventi che ci rammenta di essere ancora un tassello nell’ecosistema di Gaia. Tuttavia, anche se sia sempre arduo da credere visto lo stato in cui abbiamo ridotto il nostro pianeta, siamo la forma di vita più intelligente nell’universo conosciuto. Quindi sarebbe abbastanza imbarazzante essere sconfitti da un nemico invisibile.

## Retro programming nostalgia IV: L’Equilibrio e la Titolazione Acido/Base

La motivazione per questo articolo nasce dal mio interesse per il retro-computing connesso, da una parte, alla rivalutazione delle mie esplorazioni giovanili del calcolo scientifico in linguaggio BASIC e dall’altra, alla popolarità che, negli ultimi anni, stanno avendo nel settore amatoriale e della didattica i microcomputer su scheda singola  (single-board computer, quali, per esempio  il Raspberry Pi).  Questi piccoli computer hanno una potenza considerevolmente maggiore a un costo decisamente inferiore dei microcalcolatori degli anni 80. Questo ha reso possibile l’emulazione su questi calcolatori dei sistemi operativi di mitici modelli di home computer della Commodore e i modelli MSX.

Pertanto sta prendendo piede anche un rinnovato interesse nel linguaggio di programmazione BASIC. Questo interesse nel retro-computing riflette la nostalgia nelle grandi emozioni che negli anni 70-80 lo sviluppo della tecnologia informatica consumistica ha portato alla mia generazione. Ricordo che rimasi folgorato dalla creatività nell’uso e nella programmazione di questi microcomputer al punto che ha ridiretto i miei interessi scientifici e la mia carriera accademica.

Ho raccontato in altri articoli delle mie prime avventure di programmazione con  home computer della Commodore e i sistemi MSX alla fine degli anni ’80 e inizi degli anni ’90 e delle mie riscoperte di archeologia informatica. Tra i reperti ho rinvenuto un piccolo programma che ho usato per studiare le titolazioni acido/base sviluppato in MSX BASIC. Pertanto ho colto l’occasione per scrivere delle note sull’equilibrio acido base  e la titolazione e quindi fornire una versione restaurata e migliorate del mio programma, a gli studenti appassionati di programmazione  che sono  alle prese con questo importante concetto della chimica analitica.

$y=T(1+x)^{\frac{\lambda}{1+\lambda}}(1-x)^{\frac{1}{1+\lambda}}$,
where $T$ and $\lambda$ are two parameters, and rotate it around the central axis
\begin{aligned} x'&=&x\\ y' &=&y*cos(\theta) \\ z' &=& y*sin(\theta) \end{aligned}