#include <stdio.h>
#include <stdlib.h>
#include "surferC.h"

void save_(double *u, double *v, double *w, double *c, double *p, int *nx, int *ny, int *nz, int *t)
{
  int i, j, k;
  float val;
  int valt;
  char s[256];
  FILE *fptr;
 
#ifdef HAVE_MPI
  sprintf(s, "%s.%d.%s", "dump",*PPRANK,"dat");
#else
  sprintf(s, "%s", "dump.dat");
#endif

  fptr = fopen(s, "wb");

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        val = u[i + (*nx)*j + (*nx)*(*ny)*k];
        fwrite(&val, sizeof(float), 1, fptr);
      }

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        val = v[i + (*nx)*j + (*nx)*(*ny)*k];
        fwrite(&val, sizeof(float), 1, fptr);
      }

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        val = w[i + (*nx)*j + (*nx)*(*ny)*k];
        fwrite(&val, sizeof(float), 1, fptr);
      }

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        val = c[i + (*nx)*j + (*nx)*(*ny)*k];
        fwrite(&val, sizeof(float), 1, fptr);
      }            

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        val = p[i + (*nx)*j + (*nx)*(*ny)*k];
        fwrite(&val, sizeof(float), 1, fptr);
      }

  valt = *t;
  fwrite(&valt, sizeof(int), 1, fptr);

  fclose(fptr);
}

void iniget_(double *u, double *v, double *w, 
			double *c, double *p, int *nx, 
			int *ny, int *nz, int *t0)
{
  int i, j, k;
  float val;
  int valt;
  char s[256];
  FILE *fptr;

#ifdef HAVE_MPI
  sprintf(s, "%s.%d.%s", "dump",*PPRANK,"dat");
#else
  sprintf(s, "%s", "dump.dat");
#endif

  fptr = fopen(s, "rt");

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        fread(&val, sizeof(float), 1, fptr);
	u[i + (*nx)*j + (*nx)*(*ny)*k] = val;
      }

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        fread(&val, sizeof(float), 1, fptr);
	v[i + (*nx)*j + (*nx)*(*ny)*k] = val;
      }

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        fread(&val, sizeof(float), 1, fptr);
	w[i + (*nx)*j + (*nx)*(*ny)*k] = val;
      }

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        fread(&val, sizeof(float), 1, fptr);
	c[i + (*nx)*j + (*nx)*(*ny)*k] = val;
      }            

  for (i = 0; i < *nx; i++)
    for (j = 0; j < *ny; j++)
      for (k = 0; k < *nz; k++) {
        fread(&val, sizeof(float), 1, fptr);
	p[i + (*nx)*j + (*nx)*(*ny)*k] = val;
      }

  fread(&valt, sizeof(int), 1, fptr);
  *t0 = valt;
  
  fclose(fptr);
}