Subject: summary: problems with G94 volume runs

From: Christoph Maerker <[email-masked]>
Date: Wed, 3 Feb 1999 19:54:46 +0100 (MET)

Universite Louis Pasteur
Postal-Address: 4, rue Blaise Pascal, F- 67000 Strasbourg FRANCE
Phone: +33/3-88-41-53-19
Fax:   +33/3-88-60-63-83
X-Mailer: ELM [version 2.4ME+ PL37 (25)]
MIME-Version: 1.0
Content-Type: text/plain; charset=US-ASCII
Content-Transfer-Encoding: quoted-printable
Content-Length: 8963      

Dear CClers,

here is the summary for the problems with volume runs in G9x. I have receiv=
ed
three very helpful answers listed below. It was true that the problems
occurred on SGI's. However, the failure seems to come from Gaussian,
not the OS. Larry Cuffe provided a nice patch.

Doug Stack gave the hint to Gaussian test input file test232.com.
I ran it both with G94 and G98 on the same SGI; G94 (RevE.2) gave the zero=
=20
volume as reported before, but G98 (RevA.6) gave non-zero volumes, although
different from those listed in the test232.log file.

I tried other platforms as well (HP-UX and Linux), and got always different
results. I do not list the numbers, because everybody can do it on his own=
=20
hardware within seconds. There seems to be considerable numerical noise in=
=20
the MC procedure applied, rendering the results machine-dependent.

The recommended radii are less affected, whereas the volumes are. Thus, tho=
se
who are interested in the cavity radii and subsequent scrf=3Ddipole runs sh=
ould do
all computations on the VERY same machine !

Regards,
Christoph Maerker

PS The answers are listed according to the authors' last names.

=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D
From: "Dale A. Braden" <[email-masked]>
Subject: Re: CCL:problem with G94 volume runs

Dear Christoph,

Are you running G94 on an SGI?  I had this problem also, and the folks at
Gaussian told me that it was due to a bug in the random number generator
in the SGI operating system.  It appears in some versions of IRIX but not
others.  Gaussian did not give me any fix for it, either, so I was forced
to calculate volumes on another platform.

Best wishes,

Dale

Dale Braden
Department of Chemistry
University of Oregon
Eugene, OR 97403-1253
[email-masked]

=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D=3D=3D
From: [email-masked] (Larry Cuffe)
Subject: Re: SGI problems

Here is the code with a reference for the random number generator at the end
good luck and if there are any problems dont hesitate to contact me.
Larry Cuffe

  The version running on the SGI failed to compute molecular volumes correc=
tly. This computation is made in L604 which performs a Monte Carlo integrat=
ion over a box containing the molecule in this process points within the vo=
lume of the box are sampled a=
t a number of points and the electron density is calculated. If the electro=
n density is above a threshold value it is deemed to be part of the molecul=
e  At the end of this process the volume of the box is multiplied by the ra=
tio of occupied to empty poin=
ts and this is taken to be the molecular volume.  On the SGI machine comput=
ed molecular volumes returned by the program corresponded to either total o=
r zero occupancy.  This problem was traced back to problems with the random=
 number generator which is ca=
lled as a machine dependant function Gran() from the module of machine depe=
ndent utilities mdutil.f  This source file was edited and the call to the m=
achine dependent routine was replaced by a random number generator based on=
 the!
 proposals
    of S. D. Stearns.1.
The code for the relevent changes is detailed below.

Calling Routine in l604.f

*Deck GnMCPt
      Subroutine GnMCPt(IOut,IPrint,MaxPt,NPt,NDone,XO,YO,ZO,X1,Y1,Z1,
     $  XYZ)
      Implicit Real*8 (A-H,O-Z)
C
C     Generate the list of points for Monte-Carlo calculation
C
      Dimension XYZ(MaxPt,3)
      Save One
      Data One/1.0d0/
 1000 Format(' Points from I=3D',I10,' through I=3D',I10)
 1010 Format(' IPt=3D',I10,' XYZ=3D',3F12.6)
C
      IFirst =3D NDone + 1
      ILast =3D NDOne + NPt
      Lim =3D NPt / 2
      Do 10 I =3D 1, Lim
   10   XYZ(I,1) =3D XO + GRan(0)*X1
      Do 20 I =3D 1, Lim
   20   XYZ(I,2) =3D YO + GRan(0)*Y1
      Do 30 I =3D 1, Lim
   30   XYZ(I,3) =3D ZO + GRan(0)*Z1
      Do 40 I =3D (Lim+1), NPt
   40   XYZ(I,1) =3D XO + (One-GRan(0))*X1
      Do 50 I =3D (Lim+1), NPt
   50   XYZ(I,2) =3D YO + (One-GRan(0))*Y1
      Do 60 I =3D (Lim+1), NPt
   60   XYZ(I,3) =3D ZO + (One-GRan(0))*Z1
C
      If(IPrint.ge.3) then
        Write(IOut,1000) IFirst,ILast               =20
        Write(IOut,1010) (I,(XYZ(I,J),J=3D1,3),I=3D1,NPt)
        endif
      Return
      End

Original Code in mdutil.f

*Deck GRan
      Function GRan(IDummy)
      implicit real*8(a-h,o-z)
C
C     Working-Precision random number.
C
      Logical Called
      Integer ITime
#ifdef IBM_RS6K
#define USE_DRAND48
#endif
#ifdef _I386_
#define USE_DRAND48
#endif
#ifdef _SGI_
#define USE_RAND
      External Time
#endif
#ifdef _PMAX_
#define USE_RAND
      Real*4 Rand
#endif
#ifdef __CONVEX__
#define GRAN0_DONE
      Real*4 Secnds, R4Zero, Rand
      Save R4Zero
      Data R4Zero/0.0/
#endif
#ifdef __HPPA
#define GRAN0_DONE
$NOSTANDARD SYSTEM
      Real*4 R4Zero
      Integer ISeed
      Data R4Zero/0.0/
#endif
#ifdef IBM_PC
#define GRAN0_DONE
      Integer*2 I2Time(2), I2, J2
      Real*4 URand
      Equivalence (I2Time(1),ITime)
#endif
#ifndef GRAN0_DONE
      Integer Time
#endif
      Save Called, ITime
      Data Called/.False.
C
      If(Called) then
#ifdef __CONVEX__
#define GRAN1_DONE
        ITIME=3DIrand(ITime)
        GRan =3D Rand(ITime)
#endif
#ifdef __HPPA
#define GRAN1_DONE
        GRan =3D RAN(ISeed)
#endif
#ifdef IBM_PC
#define GRAN1_DONE
        GRan =3D URand(ITime)
#endif
#ifdef ALPHA_VMS
#define GRAN1_DONE
        Gran =3D getran()
#endif
#ifdef _SGI_
        ITime =3D Time()
        Call SRand(ITime)
#endif
#ifdef USE_RAND
#define GRAN1_DONE
        Gran =3D Rand()
#endif
#ifdef USE_DRAND48
#define GRAN1_DONE
        Gran =3D DRand48()
#endif
#ifndef GRAN1_DONE
        GRan =3D DRand(0)
#endif
      else
        Called =3D .True.
#ifdef __CONVEX__
#define GRAN2_DONE
        ITime =3D Int(20.0*Secnds(R4Zero))
        ITime =3D IRand(ITime)
#endif
#ifdef __HPPA
#define GRAN2_DONE
        ISeed =3D Int(20.0*Secnds(R4Zero))
#endif
#ifdef IBM_PC
#define GRAN2_DONE
        Call GetTim(I2,J2,I2Time(2),I2Time(1))
#endif
#ifdef ALPHA_OSF1
#define GRAN2_DONE
        ITime =3D Time()
        Call SRand(ITime)
#endif
#ifdef ALPHA_VMS
#define GRAN2_DONE
        GRan =3D getrans()
#endif
#ifndef GRAN2_DONE
        ITime =3D Time()
#endif
#ifdef __CONVEX__
#define GRAN3_DONE
        GRan =3D Rand(ITime)
#endif
#ifdef __HPPA
#define GRAN3_DONE
        GRan =3D Ran(ISeed)
#endif
#ifdef IBM_PC
#define GRAN3_DONE
        Gran =3D URand(ITime)
#endif
#ifdef ALPHA_VMS
#define GRAN3_DONE
#endif
#ifdef USE_DRAND48
#define GRAN3_DONE
        Call GSR48(ITime)
        Gran =3D DRand48()
#endif
#ifdef USE_RAND
#define GRAN3_DONE
        Call SRand(ITime)
        Gran =3D Rand()
#endif
#ifndef GRAN3_DONE
        GRan =3D DRand(ITime)
#endif
        endIf
      Return
      End

modified code=20

*Deck GRan
      Function GRan(IDummy)
      implicit real*8(a-h,o-z)
C
C     Working-Precision random number.
C     modified LPC 8 sept 98
      Logical Called
      Integer ITime
      Integer SEED
      Save SEED
      REAL*8 RANDX
      SEED=3D2045*SEED + 1
      SEED =3D SEED - (SEED/1048576)*1048576
      RANDX =3D REAL(SEED + 1)/1048577.0
      GRan=3DRANDX
      return
      end

1 S. D. Stearns, "A Portable Random Number Generator for Use in Signal Proc=
essing", Sandia Laboratory Technical Report (1981).

=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D
From: Douglas_Stack/CAS/UNO/[email-masked]

I've had this same problem w/ G94 running on an SGI Onynx sever (multiple
processors).  I ran a couple of test jobs that applied to the type of
calculations I was doing.  The test for molecular volume is test232.com, wh=
en I
performed this test, I did not get the same output as test232.log!  The
test232.log from Gaussian gives:

 Molar volume =3D  456.410 bohr**3/mol ( 40.729 cm**3/mol)
 Recommended a0 for SCRF calculation =3D  3.28 angstrom (  6.20 bohr)

My test232.log output was:

 Molar volume =3D 2389.438 bohr**3/mol (213.231 cm**3/mol)
 Recommended a0 for SCRF calculation =3D  5.33 angstrom ( 10.07 bohr)

This the only case where I've seen the output from a local job I've run
differing from the output of the Gaussian test file.

I obtained my a0 for SCRF calculations by importing the molecules into Spar=
tan,
creating an density isosurface at the 0.01 level, and having Spartan tell m=
e the
volume of that surface.  With the volume, I worked back to the radius, adde=
d 0.5
A and used that number in SCRF.

If Doug Fox see this tread, maybe he can give some insight to this problem.

Regards,

Douglas E. Stack
Assistant Professor
Department of Chemistry
University of Nebraska at Omaha
Omaha, NE 68182-0109
(402) 554-3647
(402) 544-3647 (fax)
[email-masked]