xmd Code

/*

SOURCE CODE FOR READING AND WRITING RCV FILES
Uses particle structure, modelled after Pascal version of PLOTL

READRCV : read a single step from a rcv file.
WRITERCV: write a single step to a rcv file.

 6 Jun 1991
28 Oct 1991
 6 Nov 1991 - add code to leave types alone if selkeep flag on
15 Nov 1991 - set trans to 0 if boundrep is on
17 Nov 1991 - streamline file read
21 Nov 1991 - error with trans, was setting to zero when boundrep data
				  was present rather than vice versa.
15 May 1992 - error occured in WRITERCV with line maxtbl += 1.5;.  Computer
				  would hang when compiled with CLIPRCV.	Changed to maxtbl += 1;.
29 Jul 1992 - modified RCVIO0 (stepdata based) to particle based for PMD
 4 Feb 1993 - modify NSTR value (from 81 defined in PARTICLE.H to 256).
	Dec 1993 - change cur to double *, array of x y z x y z .. interleaved
 9 Dec 1993 - make default type 0
18 Mar 1996 - changed to allocate initial memory
*/

/*
************************************************************************
Compile Switches
************************************************************************
*/
/*  Use dummy RDF table (values increases linearly from 0 to 99)  */
#define USE_DUMMY_RDF 


/*
************************************************************************
File Includes
************************************************************************
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "strsub.h"
#include "cdalloc.h"
#include "rcvio.h"

/*
************************************************************************
Defines
************************************************************************
*/
#define NSTR 256
#define ANG_TO_CM 1e8

#define X 0
#define Y 1
#define Z 2

#define BOOLEAN int
#define FALSE 0
#define TRUE  1

#define TEST_MACRO \
	{ \
	int i,n,NumRead; \
	char *test; \
	test = "1 2 3 4"; \
	NumRead = sscanf (test, "%i %i", &i, &n); \
	printf ("File %s Line %i: NumRead,i,n = %i %i %i\n", \
		__FILE__, __LINE__, NumRead, i, n);  \
	}

#define WRITEMSG \
	printf ("In FILE %s Line %i\n", __FILE__, __LINE__);


#define WRITEVAR(VAR_NAME,VAR_TYPE) \
   printf ("FILE %s LINE %i :", __FILE__, __LINE__); \
   printf ("%s = ", #VAR_NAME); \
   printf (#VAR_TYPE, (VAR_NAME) ); \
   printf ("\n");



/*
************************************************************************
Module Variables
************************************************************************
*/
/* TRANSLATE: CONVERTS TWO Code byte VALUES INTO integer FROM 0 TO 9999 */
/* THIS TABLE TRANSLATES BETWEEN charACTERS SENT FROM PURDUE
 CDC 205 AND THE CRAY AT LLL TO THE  IBM PC.  THE SEND PROGRAM
 TRANSLATES THE FOLLOWING charACTERS
 CHARACTER		 REASON
 ---------		 -------
 27 -> 17		 CRAY
 28 -> 14
 30 -> 15
 31 -> 25		 CYBER LINE FEED CONTROL CHARACTER
 92 -> 16
(94 -> 20)		 ASCII TO EBCDIC CONVERSION  )
(124-> 21)		 ASCII TO EBCDIC CONVERSION  )
 126-> 20		 VAX MAIL UTILITY
*/

int off	=	 27;
int base =	100; /* NOTE:	Base is fixed in values of table2  */
int off2 = 2727; /* off2 = off * (base + 1)	*/
int table[128] =
	{	0,   1,	 2,	3,   4,	 5,	6,   7,	 8,	9,  10,
			 11,	12,  13,  28,	30,  92,  27,	18,  19, 126,
			 21,	22,  23,  24,	31,  26,  27,	28,  29,  30,
			 31,	32,  33,  34,	35,  36,  37,	38,  39,  40,
			 41,	42,  43,  44,	45,  46,  47,	48,  49,  50,
			 51,	52,  53,  54,	55,  56,  57,	58,  59,  60,
			 61,	62,  63,  64,	65,  66,  67,	68,  69,  70,
			 71,	72,  73,  74,	75,  76,  77,	78,  79,  80,
			 81,	82,  83,  84,	85,  86,  87,	88,  89,  90,
			 91,	92,  93,  94,	95,  96,  97,	98,  99, 100,
			101, 102, 103, 104, 105, 106, 107, 108, 109, 110,
			111, 112, 113, 114, 115, 116, 117, 118, 119, 120,
			121, 122, 123, 124, 125, 126, 127					 };

int table2[128] =
 {  000,   100,	200,	 300,   400,	500,	 600,   700,	800,	 900,  1000,
			 1100,  1200,	1300,  2800,  3000,	9200,  2700,  1800,	1900, 12600,
			 2100,  2200,	2300,  2400,  3100,	2600,  2700,  2800,	2900,  3000,
			 3100,  3200,	3300,  3400,  3500,	3600,  3700,  3800,	3900,  4000,
			 4100,  4200,	4300,  4400,  4500,	4600,  4700,  4800,	4900,  5000,
			 5100,  5200,	5300,  5400,  5500,	5600,  5700,  5800,	5900,  6000,
			 6100,  6200,	6300,  6400,  6500,	6600,  6700,  6800,	6900,  7000,
			 7100,  7200,	7300,  7400,  7500,	7600,  7700,  7800,	7900,  8000,
			 8100,  8200,	8300,  8400,  8500,	8600,  8700,  8800,	8900,  9000,
			 9100,  9200,	9300,  9400,  9500,	9600,  9700,  9800,	9900, 10000,
			10100, 10200, 10300, 10400, 10500, 10600, 10700, 10800, 10900, 11000,
			11100, 11200, 11300, 11400, 11500, 11600, 11700, 11800, 11900, 12000,
			12100, 12200, 12300, 12400, 12500, 12600, 12700 				 };


/*
************************************************************************
Local function prototypes
************************************************************************
*/
void *ReAllocate
(
void *Pointer,
int Number,
int Size,
char *FunctionName,
char *VariableName
);


/*
************************************************************************
Exported Functions
************************************************************************
*/
int readrcv (FILE *fin, Particle_t *a)
	{
	int res	= 10000;
	double	  scale[3], tablefactor;
	unsigned  np,i, n,iline, ntable;
	long		 maxtable;
	int		 teof,xlo,xhi,ylo,yhi,zlo,zhi;
	int done;
	char instr[NSTR], *tstr;

	/* READ FILE HEADER								*/
	xlo = getc (fin); /*  GOBBLE UP FIRST CHARACTER  */
	if (xlo==EOF)
		return(EOF);
	fgets (instr, NSTR, fin);
	teof = sscanf (instr, "%li %li %lf %lf %lf %i %i\n",
		&a->run, &a->step, &a->time, &a->tempc, &a->temp,
		&np, &a->ndim);
	if (teof==EOF)
		return(EOF);

	/*  AJUST NUMBER OF POINTS  */
	np--;
	for (i=0;i<8;i++) {
		fgets (instr, NSTR, fin);
		/*  COPY FIRST NTITLESTR CHARACTERS TO TITLE  */
		n = 0;
		while (instr[n] && instr[n]!='\n' && n<NTITLESTR-1)
			{
			a->title[i][n] = instr[n];
			n++;
			}
		a->title[i][n] = 0;
	}
						  /*	READ MISC DATA LINE	*/
	fgets  (instr, NSTR, fin);
	tstr = instr;
	/*  Skip mass (obsolete)  */
	dblstr (&tstr);
	a->vibp	 = dblstr (&tstr);
	a->dtime  = dblstr (&tstr);
	a->cutoff = dblstr (&tstr);
	/*  Don't read in strain or stress (obsolete)  */

	/*  READ BOX  */
	fgets  (instr, NSTR, fin);
	tstr = instr;
	for (i=0;i<3;i++)
		a->bcur[i] = 1E-8*dblstr(&tstr);

	/*  READ Trans IF PRESENT	*/
	if (*tstr!=0)
		for (i=0;i<3;i++)
			a->trans[i] = 1E-8*dblstr(&tstr);
	else
		for (i=0;i<3;i++)
			 a->trans[i] = 0;

	/*  BOUNDARY IF PRESENT  */
	if (*tstr!=0) {
		for (i=0;i<3;i++) {
			a->surf[i] = intstr(&tstr);
			a->surf[i] = !a->surf[i];
		}
	}
	else
		for (i=0;i<3;i++) {
			a->surf [i] = 1;
			a->trans[i] = 0.0;
		}

	/*  Calculate scale for converting from code to coordinates  */
	for (i=0;i<3;i++)
		scale[i] = a->bcur[i] / res;

	/*  Reallocate particles if new number different from previous number  */
	if (np != a->np)
		ReallocateParticle (a, np);


	/*  Read COORDINATES IN Coded FORM	 */
	n = 0;
   iline = 0;
	do {	 /* UNTIL END-OF-FILE  OR	UNTIL ASCII 0 Code FOUND */
		i = 0;
	   fgets (instr, NSTR, fin);
		tstr = instr;
      do {		 /* UNTIL 13 POSITIONS Read IN	OR  END-OF-FILE  OR	0 FOUND*/
			xlo = *(tstr++);
			if (xlo!='\n' && xlo!=0)   {
				xhi = *(tstr++);
				ylo = *(tstr++);
				yhi = *(tstr++);
				zlo = *(tstr++);
				zhi = *(tstr++);
				if (feof(fin)) {
					fprintf (stderr, "ERROR (readrcv):  Premature end of file in step %i.\n", a->step);
					exit(1);
				}
				/*  TEST FOR N TOO BIG	*/
				if (n==np) {
					fprintf (stderr, "ERROR (readrcv): Data is corrupted.\n");
					fprintf (stderr,
						" Number of particles in header (%i) and body (%i)\n", np, n);
					fprintf (stderr,
						" for run (%i) step (%i) do not agree.\n", a->run, a->step);
					exit(1);
				}
				a->cur[X + 3*n] =
					scale[X] * (table[xlo]+table2[xhi]-off2) - a->trans[X];
				a->cur[Y + 3*n] =
					scale[Y] * (table[ylo]+table2[yhi]-off2) - a->trans[Y];
				a->cur[Z + 3*n] =
					scale[Z] * (table[zlo]+table2[zhi]-off2) - a->trans[Z];
/*				if (!a->selkeep)
					a->type[n] = 0;
*/
				n++;
				i++;
			}
		} while (i!=13 && xlo !=0 && xlo!='\n' && n!=np);
		iline++;
	} while (xlo !=0 && xlo!='\n' && n!=np);

	if (n<np) {
		fprintf (stderr, "ERROR (readrcv): Only %i of %i points read for step %i.\n",
								n, np, a->step);
		exit(1);
	}

	/*  Read RDF table  */
	fscanf (fin, "%li\n", &maxtable);
	tablefactor = maxtable/10000.0;
	ntable = 0;
	do {	  /* UNTIL END-OF-FILE	OR  UNTIL ASCII 0 Code FOUND */
		iline = 0;
		do {	  /*		UNTIL 13 POSITIONS Read IN
					 OR  END-OF-FILE	OR  0 FOUND 	*/
			xlo = getc(fin);
			done = (xlo==13 || xlo==0 || xlo==EOF);
			if (!done) {
				xhi = getc(fin);
				iline++;
				if (xhi!=0 && xhi!=13 && xhi!=EOF)
					a->rdftable[ntable++] =
							  tablefactor * (table[xlo]+table2[xhi]-off2);
			}
			done = (done || ntable==100);
		} while (iline<39 && !done);
		while (getc(fin)!='\n');
	} while (!done);

/*  EVERYTHING OK:  RETURN 0	*/
	return(0);
}
/*  END READRCV  */


void readtype (FILE *fin, Particle_t *a)
	{
	int itype;
	int NumType;
	int InputType;
	char InputString[NSTR];
	char *TempString;

	/*  Read header  TYPELIST n	*/
	fgets (InputString, NSTR, fin);
	TempString = InputString;

	/* Strip leading token	*/
	strhed (&TempString);

	/*  Read number  */
	NumType = intstr (&TempString);

	/*  READ TYPES  */
	itype = 0;
	while
	(
	EOF != fscanf(fin,"%i", &InputType)  &&
	itype<NumType
	)
		{
		a->type[itype] = InputType;
		itype++;
		}

	/*  Check for corret number of types read  */
	if (itype<NumType)
		{
		printf ("ERROR (readtype):  premature end of file.\n");
		exit(1);
		}
	}

void readilist(FILE *fin, Particle_t *a)
	{
	unsigned i,n,is,nalloc,*tind;
	fscanf (fin, "ILIST%i", &n);

	/*  ALLOCATE TEMPORARY INDEX STORAGE  */
	tind = (unsigned *)	calloc (n, sizeof(unsigned));
	if (tind==NULL)
		{
		printf ("ERROR (readilist):  ");
		printf ("Insufficient memory for allocating temporary indices.\n");
		exit(1);
		}

	/*  READ ILIST  */
	i = 0;
	while (EOF!=fscanf(fin,"%i", &is)  &&  i<n)
		{
		tind[i] = is;
		i++;
		}
	if (i<n)
		{
		printf ("ERROR (readilist):  premature end of file.\n");
		exit(1);
		}

	/*  FIND MAXIMUM INDEX	*/
	nalloc = tind[0];
	for (i=1;i<n;i++)
		if (tind[i]>nalloc)
			nalloc = tind[i];

	/*  Check for maximum index value too large	*/
	if (nalloc > a->np)
		{
		printf ("ERROR (readilist):  ");
		printf ("Value of largest index exceeds current number of particles.\n");
		exit(1);
		}

	/*  FREE STORAGE	*/
	free(tind);
	}

/*
************************************************************************
Local functions
************************************************************************
*/
/*
Re-allocate or Newly Allocate particle coordinates and
all arrays whose size depends on number of particles.
(if they are already in use, determined if pointer not NULL)
Arrays such as neighbor list arrays are not affected.
*/
/*
Allocate or reallocate variable storage.
If no room available then print error message and exit.
*/
void *ReAllocate
(
void *Pointer,
int Number,
int Size,
char *FunctionName,
char *VariableName
)
	{
	void *NewPointer;

	/*  If request is for 0 bytes, then just free allocation  */
	if (Size==0 || Number==0)
		{
		if (Pointer!=NULL)
			free (Pointer);
		return NULL;
		}

	/*  If pointer current null, just do allocation  */
	if (Pointer==NULL)
		NewPointer = calloc (Number, Size);

	else
		NewPointer = realloc (Pointer, Number*Size);

	if (NewPointer==NULL)
		{
		printf ("ERROR (%s):  Cannot allocate variable %s.\n",
			FunctionName, VariableName);
		exit(1);
		}

	return NewPointer;
	}