The picture on the front cover is the orthorhombic perovskite structure projected down the Z axis, copied from the article on Computer simulation of mantle minerals
July 10 saw the culmination of many weeks of preparation as RAL opened it doors to invited guests, schoolchildren and the general public. RAL Open Days are held once every ten years; this year they coincided with SERC's twenty-fifth anniversary.
It was estimated that 4,400 people visited the laboratory site during the five days (10-14 July), many of whom toured the Atlas Centre. Traffic was a bit slow in the first days: at Atlas, we suffered somewhat from being in a remote location. However, by the final day, the trickle had turned into a flood - Atlas Centre staff who were demonstrating found themselves a bit hoarse by the end of it.
Ten years has seen many changes at RAL and Atlas is no exception. The tour of the computing facilities included a few exhibits that gave a sense of perspective to our progress. Visitors could examine the consoles of the Ferranti Atlas and the IBM 360/195 before continuing to see the IBM 3090 and the Cray X-MP. Another vivid display of the speed at which technology changes was found in the operations area. On one hand, visitors could watch operators manually mounting round tapes and square cartridges; on the other, a robot eye's view of the Automated Cartridge System whizzing around mounting tape.
Even if no one could sit on the Cray, one could spin the hard disk of a Pc. After touring the network nerve centre in the Telecoms area, one could then send an electronic mail message. The tour included several such demonstrations, which tried to help people understand how computers worked, and show the link between the hardware, the software and the results.
Besides showing a dismantled personal computer, the PC Support Group had two other demonstrations. One PC ran a demonstration of a package called Mathematica 386. This is a mathematical modelling and equation solving package that allows a scientist to solve complex equations more quickly than by normal means. A second demonstration showed a PC running with two output screens allowing interactive debugging on one screen and program output on another.
The Office Automation demonstration featured IBM's PRofessional OFfice System (PROFS). Visitors were invited to try the various PROFS utilities, including sending and receiving mail, working with a diary, scheduling meetings, and creating and editing documents.
A demonstration of administrative databases showed several projects, both existing and under development. Of particular interest was the EXchange of Information of Research ProjecTS, or EXIRPTS. This will contain information and abstracts of current (that is, not yet published) scientific research underway in eight countries.
Graphics Group had a variety of demonstrations. Using the IRIS Silicon Graphics Workstation, one could, among other things, use a mouse to move across the surface of the Indian Ocean. A constantly running video demonstrated the most recent output of the Atlas Video Facility. The sequences included a flight across the Indian Ocean, a sequence from the Fine Resolution Antarctic Model (FRAM) video, an atmospheric model showing how chlorofluorocarbons are dispersed, simulated action of a crank shaft, a simulation of fluid flow dynamics, and excerpts from William Latham's videos, The Conquest of Form and The Evolution of Form.
Stills from William Latham's latest works were also on display on an IBM 5080 graphics terminal. A VAXstation displayed samples from UNIRAS and a Mitsubishi colour printer busily generated souvenirs from all the graphics demonstrations.
Students from the Wantage Sixth Form were presenting their project on global warming, the greenhouse effect and ozone layer depletion on a Commodore Amiga. They are using data gathered at RAL for the project.
The Joint Network Team were giving demonstrations of electronic mail, file transfer and the newly implemented X500 Pilot Directory Service. This service allows one to look up address or other information about another user. Eventually it will provide worldwide access to electronic directories.
An Apple Macintosh was used to display some of the work of the Advanced Research Computing Unit. Two simulations were running: one simulating the effects of a laser on a metal strip, the other showing the turbulence caused by a fast moving corridor of liquid inside a stationary liquid.
In the library, visitors could use JANET to find books in various university libraries, or find the best way to get from A to B using the library's AUTOROUTE software.
The effort that went into preparing for Open Days was rewarded by the interest shown by most of the visitors to the Atlas Centre. It was an opportunity to describe and exhibit the many computing services offered. It was also an opportunity for some of our remote users to see what it is they have been using.
Since the last issue of FLAGSHIP went to press the UNICOS service on the Cray X-MP/416 has settled down, but still has some niggling problems for users of tapes and very large datasets. The major problems experienced during the migration were largely caused by trying to run COS and UNICOS services side by side with inadequate disk resources. Since COS was closed down the UNICOS service has been very much better.
The remaining problems centre around two areas, user tapes and automatic migration of datasets to the Automated Cartridge System (ACS). User jobs requesting tape mounts frequently "hang" during execution, usually with no error message. The problem is not properly understood but is probably due to the curious mix of UNICOS 5.0 releases that we are currently running. We will get better support from Cray once we have moved to UNICOS 5.1, probably in late August. In the meantime, if a job requesting tape mounts has started execution but fails to complete within an hour or two, contact the Service Line; it may be necessary for us to rerun the job.
A more serious problem in the short term is that the "data migration" software that swaps user datasets to and from the ACS is not reliable and has poor functionality. It insists on copying unchanged datasets back to the ACS and has multiple copies of many datasets on the ACS, with the result that the pool of allocated tapes keeps filling up. When the pool of tapes runs out it may take an hour or more to load the ACS with new tapes. During this time all user jobs are suspended. The unreliability of data migration is in part related to the problem of handling tapes and is compounded by the poor and frequently misleading error messages that are becoming a feature of UNICOS.
Cray promise a much better data migration product and more comprehensible error messages with UNICOS 6.0. The problem of multiple copies of datasets will be removed and the data migration criteria and controls will be more flexible. Don't hold your breath on this - it is likely to be the end of 1990 before it is available.
The timing of the move to UNICOS was influenced by the upgrade to a /416 provided by Cray to run their potential bureau service. So far we have been reluctant to increase the maximum job size beyond 7 Million words due to the problems of disk space management. As we become more confident that we can run a reliable service the relaxation of this memory limit will move up the priority list. In the meantime the UNICOS service generally shows less than 2% idle time, except when data migration is playing up.
In my last note I caused great offence by thanking some of the RAL staff but not others; there are many "forgotten heroes" and the list below should be more complete:
Thanks
An Exabyte is a low-speed, high-density tape drive that uses standard, domestic 8mm video cartridges. These are double spool cartridges containing 8mm tape, somewhat similar to an audio tape cassette but a little big bigger. Each cartridge can hold up to 2.3 Gigabytes of data. Exabyte drives and controllers for attaching them to VAX computers have been available for several years at relatively low cost. As a result of these two factors, the Exabyte has become a very popular device for taking dumps of magnetic disk filestores and for transferring large quantities of data between sites.
CCD has recently purchased a controller and two Exabyte drives for attachment to the IBM 3090. This controller is manufactured only by a company in the USA and was designed primarily for smaller IBM machines. Consequently there have been some teething problems in getting it to work correctly from the VM/CMS system on the 3090. We hope these problems have now been solved and we are undertaking the small amount of software development needed to support this device and provide a service.
The initial service will provide a facility for copying from 3480 cartridge tapes to an Exabyte. This will be under user control from a job running in SLAC Batch. The intention is to provide an EXEC which will automate the copying of files from one or many 3480s. Initially this will be a straight tape to tape copy but a later phase will stage the data via disk for efficiency. Facilities for copying from disk to Exabyte and from Exabyte to disk or tape are planned future developments. Further details of this service will be announced via CMS NEWS when it is ready for release and will be the subject of a future FLAGSHIP article.
Computer modelling techniques can now be successfully employed to investigate many important mineral properties including their phase relations, equations of state and diffusion mechanisms. The ability to model minerals is particularly useful because of the severe technical difficulties encountered in many types of experimental investigation and the long time periods (thousands, even millions of years) over which many mineral processes occur. Such practical problems are accentuated in the study of the high pressure - high temperature phases thought to occur in the Earth's mantle. (The Earth's mantle extends from about 100 km to 2900 km and has a pressure range of up to 136 GPa (1,342,216 atmospheres) and a temperature range of up to 3000 K. A cross section of the Earth is shown in figure 1.)
Ideally, a simulation of a solid would be based on a direct solution of the Schrodinger equation to describe the interactions of electrons and nuclei within a system. However, for all but the most simple minerals, such quantum mechanical studies are not currently feasible and, therefore, we have chosen to employ an atomistic model.
The atomistic simulation technique is described in Parker and Price (1989). It is based on choosing an algorithm or interatomic potential to describe the energy surface within a solid. To model silicate minerals, the interatomic potential is usually considered to be composed of a long-range Coulombic component (Uc )and a term to model the short range forces (Ur ). Additionally, a shell model may be included to model ionic polarisability (U s )and a bond-bending term may be used to describe, in some measure, the directional component of covalent bonding. The lattice energy (Ut )can be evaluated from the summation of these terms over pairs or groups of atoms:
Ut = Uc + Ur + Ug + UbThe Coulombic potential term has the following form for a pair of atoms i and j considered to have point charges of qi and qj which are separated by a distance rij
Uc = qiqj/rij
In our simulations of silicates, we use a Buckingham potential for the short-range term:
Ur = Aij. exp( -rij/ Bij) - Cij rij -6
where Aij , Bij and Cij are parameters which need to be derived for each pair of atoms i and j considered. This analytical expression is chosen as it can model the repulsion between the overlapping electron charge cloud and also the attractive van der Waals interaction.
The shell model, sketched in figure 2, provides a simple way of including the effects of ionic polarisability in a simulation. The ion is considered to consist of two components, a core, on which all the mass is centred, surrounded by a massless shell. The core and shell are each assigned a charge, the sum of which equals the total ionic charge. They are coupled by a harmonic spring constant Kg . Only the shells are included in the short range force calculations. For a core and shell with separation w , the core-shell interaction, Us , is given by:
Us = Ks.w2
The bond-bending term for a bond between atoms j - i- k has the form:
Ub = Kb (θ ijk - θ 0 ) 2
where Kb is the spring constant, θ ijk is the simulated bond angle and θ 0 is the equilibrium bond angle.
The reliability of the atomistic model depends on the quality of the interatomic potential parameters. Aij , Bij , and Cij can be derived using quantum mechanical methods or empirically by fitting them so as to reproduce experimental data. The spring constants for the shell model and bond-bending term can only be derived empirically. It has proved possible to transfer potential terms fitted to a particular interaction in one solid to model a similar interaction in another solid. This "transfer" approach has proved highly successful in the simulation of high pressure magnesium silicate phases.
To start a simulation, it is necessary to have a suitable set of interatomic potentials and to supply the atomic positions. Then the simulation proceeds by an iterative process in which the atomic positions are adjusted until the minimum energy configuration is achieved. The minimisation method used is that of Newton-Raphson which requires the evaluation of both the first and second derivatives of the lattice energy with respect to the atomic coordinates. It is the calculation and subsequent manipulation of the second derivative matrix that is the most computationally expensive part of the simulation, both in terms of cpu time used and memory required. For example, a simulation of Mg2SiO4 olivine (using a unit cell of 44 species) takes a few hundred cpu seconds and around 2 Megawords of memory whereas a simulation of the zeolite ZSM5 with 480 species in the unit cell requires more memory than is currently available on the CRAY at RAL.
Once the first and second derivatives are calculated, properties including the elastic and dielectric constants and the lattice dynamics can be evaluated.
PARAPOCS (Phonon Assisted Relaxation Applied to the Prediction of Crystal Structures) is one of a suite of atomistic simulation codes designed for modelling inorganic solids. The related programs (THBREL, THBPHON, THBFIT and CASCADE) are supported by Maurice Leslie at the Daresbury Laboratory as part of CCP5.
PARAPOCS was written and is being supported and developed at the University of Bath. This program provides a valuable tool for modelling minerals because it includes the effects of both hydrostatic pressure and temperature in the simulation. I will outline a little more of how PARAPOCS works before discussing some of the ways in which it is currently being used to investigate the properties of high pressure-high temperature minerals phases.
PARAPOCS is able to evaluate the internal pressure in a crystal (defined as the first derivative of the lattice energy with respect to strain) and also the elastic constants (the second derivative of the lattice energy with respect to strain). The effect of temperature is simulated by calculating the kinetic pressure where the kinetic pressure is defined as the derivative of the vibrational free energy with respect to volume. The vibrational free energy is a function of the phonon frequencies calculated at a number of points in the irreducible Brillouin zone. The kinetic, internal and applied hydrostatic pressure are summed to determine the total pressure (P). The strain on the crystal (e) is then given by:
e = P.C-1
where C is the elastic constant tensor. The calculation then proceeds iteratively to adjust the atomic positions to remove all the strain acting on the crystal and thus to determine the equilibrium structure at the specified pressure and temperature. In other words, the crystal structure is adjusted until the internal pressure and the applied pressure holding the crystal together exactly balances the kinetic pressure which is trying to force the crystal apart.
Besides the vibrational free energy, other important thermodynamic properties such as the heat capacity, entropy, volume coefficient of thermal expansion and Gruneisen parameters are also functions of the phonon frequencies and are, therefore, easily evaluated by PARAPOCS.
PARAPOCS is being used to investigate the properties of a wide range of crystalline solids including superconducting oxides and zeolites. Here, however, I will describe its application to the study of high pressure minerals predicted to be important phases in the Earth's mantle.
(Mg,Fe)SiO3 with the perovskite structure is thought to make up 70% - 80% of the Earth's lower mantle (and thus 40% of the whole Earth!). It has been synthesised in small quantities in specialist high pressure apparatus, but it is not yet experimentally possible to determine many of its geophysically important properties such as its structure, volume expansion coefficient and rheological properties at lower mantle pressures and temperatures.
PARAPOCS was used to investigate the structure of the end-member MgSiO3 perovskite at high pressure and high temperature. There are many variations of the ideal cubic perovskite structure. At low pressures MgSiO3 has been observed to have an orthorhombically distorted perovskite structure (figure 3) caused by rotation of the SiO6 octahedra. It is suggested that this distortion may change with increasing pressure and temperature.
The accuracy of any simulation will be critically dependent on the interatomic potentials used. A number of different empirically derived potentials were, therefore, tested. The success of each interatomic potential in reproducing the observed properties of MgSiO3 was used to evaluate the usefulness of that potential model. Our most successful potential was transferred from MgO and SiO2 quartz and then the A term for the Mg - O interaction was modified by fitting to take into account the change in coordination of the Mg ion from 6 in MgO to 8 (or 12) in MgSiO3 perovskite. Simulations using PARAPOCS show that the distortion of the MgSiO3 perovskite structure changes little with increasing pressure and temperature.
A knowledge of the effect of pressure on the thermal expansion is required to calculate the density of mantle mineral phases at mantle conditions. Recently, it has been shown experimentally that the thermal expansion coefficient a can be expressed as:
(α/α 0 ) = (ρ/ρ 0 ) -δT
where δT is known as the Anderson-Gruneisen parameter and ρ is the density. The Anderson-Gruneisen parameter has been determined experimentally for some mantle magnesium silicate phases (MgO and forsterite) to have a value of between 2 and 3. However, it has also been estimated from seismic wave velocities to have a value of between 6 and 7. As PARAPOCS has previously been shown to model well the Mg2SiO4 - MgSiO3 phases (Price et al. 1989), Reynard and Price (1990) have used PARAPOCS to calculate the Anderson-Gruneisen parameter for the mantle phases Mg2SiO4 forsterite, γ - Mg2SiO4 , MgO and MgSiO3 perovskite with the aim of investigating the discrepancy between the experimental and seismic estimates, whether the Anderson-Gruneisen is independent of pressure over the range of pressures characteristic of the Earth's mantle and also to calculate δT for phases for which no experimental values are yet available. Calculations were carried out at temperatures of 300 to 2000 K and pressures of 5 to 120 GPa. The Anderson-Gruneisen parameters calculated for all four phases are between 2 and 3 and are, therefore, in good agreement with the experimental results. δT is also confirmed to be constant with pressure up until the core-mantle boundary. The seismic determination of the Anderson-Gruneisen parameter now needs to be re-examined.
The previous two examples demonstrate how PARAPOCS is being used to extend and question experimental data. One of the major outstanding geophysical problems concerns the nature of mantle convection. It is, for example, not known whether the mantle convects as a whole or whether it is divided into two (or more) separately convecting layers. Currently work is in progress to use PARAPOCS to help investigate the rheological properties of mantle materials. Specifically, we are working on the calculation of the activation volume of diffusion (V*) for lower mantle phases. The activation volume is defined as:
V* = dG/dP
Although the activation volume is an important component in rheological models of the mantle, the activation volumes currently used in mantle rheological models are often based on estimates derived from data on materials such as metals and alkali halides and rarely include any direct input from the properties of the lower mantle minerals themselves.
The activation volume is being calculated using a supercell method in PARAPOCS in which a (charge neutral) defect such as an MgO Schottky pair is embedded in a large cell of atoms on their perfect lattice positions. This involves running simulations using many hundreds of atoms and requiring large amounts of computer memory. Most of our work is run on the CRAY X-MP, but some of these calculations will be run on the IBM 3090 to take advantage of the larger memory available. Preliminary results for MgO predict the activation volume of formation of an Schottky defect to be 5cm3mol-1 at 100 GPa pressure and a temperature of 1500 K. Such results should provide valuable input to mantle rheological models.
In the future, we hope to extend PARAPOCS to calculate the intensities of Raman and infrared active modes. Also we aim to include anharmonic effects to improve the accuracy of our simulations at high temperature. With these developments, PARAPOCS will continue to provide a valuable tool in the study of mineral properties.
Parker, S.C. and Price, G.D. (1989) Computer modelling of phase transitions in minerals. Advances in Solid-State Chemistry 1, 295-327
Price, G.D., Wall, A. and Parker, S.C. (1989) The properties and behaviour of mantle minerals: a computer-simulation approach. Phil. Tran. R. Soc. Lond. A 328, 391-407
Reynard, B. and Price. G.D. (1990) Thermal expansion of mantle minerals at high pressures - a theoretical study. Geophys. Research Lett. In press.
