post.c

OpenFVM-v1.1 open source cfd code

  fprintf (fp, "\n$EndView\n");

  free (sf);

}

void
WriteBinaryScalarFace (FILE * fp, char *label, int iphi, double curtime)
{

  int i;

  int face;

  double vs;

  double *sf;

  vs = 0.0;
  
  fprintf (fp, "$View\n");

  sf = calloc (nbfaces, sizeof (double));

  for (i = 0; i < nbfaces; i++)
    {

      face = i;

      if (iphi == 0)
	vs = V_GetCmp (&xuf, face);
      if (iphi == 1)
	vs = V_GetCmp (&xvf, face);
      if (iphi == 2)
	vs = V_GetCmp (&xwf, face);
      if (iphi == 3)
	vs = V_GetCmp (&xpf, face);
      if (iphi == 4)
	vs = V_GetCmp (&xTf, face);
      if (iphi == 5)
	vs = V_GetCmp (&xsf, face);

      sf[face] = vs;

    }

  PrintBinaryHeaderScalarFace (fp, label);
  PrintBinaryTime (fp, curtime);
  PrintBinaryScalarFace (fp, sf);

  fprintf (fp, "\n$EndView\n");

  free (sf);

}

void
WriteBinaryVectorMagnitudeFace (FILE * fp, char *label,
				int iphix, int iphiy, int iphiz,
				double curtime)
{

  int i;

  int face;

  double vx, vy, vz;

  double *sf;
  
  vx = 0.0;
  vy = 0.0;
  vz = 0.0;

  sf = calloc (nbfaces, sizeof (double));
  
  for (i = 0; i < nbfaces; i++)
    {

      face = i;

      if (iphix == 0)
	vx = V_GetCmp (&xuf, face);
      if (iphix == 1)
	vx = V_GetCmp (&xvf, face);
      if (iphix == 2)
	vx = V_GetCmp (&xwf, face);
      if (iphix == 3)
	vx = V_GetCmp (&xpf, face);
      if (iphix == 4)
	vx = V_GetCmp (&xTf, face);
      if (iphix == 5)
	vx = V_GetCmp (&xsf, face);

      if (iphiy == 0)
	vy = V_GetCmp (&xuf, face);
      if (iphiy == 1)
	vy = V_GetCmp (&xvf, face);
      if (iphiy == 2)
	vy = V_GetCmp (&xwf, face);
      if (iphiy == 3)
	vy = V_GetCmp (&xpf, face);
      if (iphiy == 4)
	vy = V_GetCmp (&xTf, face);
      if (iphiy == 5)
	vy = V_GetCmp (&xsf, face);

      if (iphiz == 0)
	vz = V_GetCmp (&xuf, face);
      if (iphiz == 1)
	vz = V_GetCmp (&xvf, face);
      if (iphiz == 2)
	vz = V_GetCmp (&xwf, face);
      if (iphiz == 3)
	vz = V_GetCmp (&xpf, face);
      if (iphiz == 4)
	vz = V_GetCmp (&xTf, face);
      if (iphiz == 5)
	vz = V_GetCmp (&xsf, face);

      sf[face] = sqrt (vx * vx + vy * vy + vz * vz);

    }

  fprintf (fp, "$View\n");

  PrintBinaryHeaderScalarFace (fp, label);
  PrintBinaryTime (fp, curtime);
  PrintBinaryScalarFace (fp, sf);

  fprintf (fp, "\n$EndView\n");

  free (sf);

}

void
WriteBinaryScalarElement (FILE * fp, char *label, int iphi, double curtime)
{

  int i;

  int element;

  double vs;

  double *se;

  vs = 0.0;
  
  fprintf (fp, "$View\n");

  se = calloc (nbelements, sizeof (double));

  for (i = 0; i < nbelements; i++)
    {

      element = i;

      if (iphi == 0)
	vs = V_GetCmp (&xul, element);
      if (iphi == 1)
	vs = V_GetCmp (&xvl, element);
      if (iphi == 2)
	vs = V_GetCmp (&xwl, element);
      if (iphi == 3)
	vs = V_GetCmp (&xpl, element);
      if (iphi == 4)
	vs = V_GetCmp (&xTl, element);
      if (iphi == 5)
	vs = V_GetCmp (&xsl, element);

      se[element] = vs;

    }

  PrintBinaryHeaderScalarElement (fp, label);
  PrintBinaryTime (fp, curtime);
  PrintBinaryScalarElement (fp, se);

  fprintf (fp, "\n$EndView\n");

  free (se);

}

void
WriteBinaryVectorElement (FILE * fp, char *label, int iphix, int iphiy,
			  int iphiz, double curtime)
{

  fprintf (fp, "$View\n");

  PrintBinaryHeaderVectorElement (fp, label);
  PrintBinaryTime (fp, curtime);
  PrintBinaryVectorElement (fp, iphix, iphiy, iphiz);

  fprintf (fp, "\n$EndView\n");

}

void
WriteBinaryStreamFunction (FILE * fp, char *label, double curtime)
{

  double *sn;

  fprintf (fp, "$View\n");

  sn = calloc (nbnodes, sizeof (double));

  CalculateStreamFunction (sn);

  PrintBinaryHeaderScalarElement (fp, label);
  PrintBinaryTime (fp, curtime);
  PrintBinaryScalarNode (fp, sn);

  fprintf (fp, "\n$EndView\n");

  free (sn);
}

void
WriteBinaryVorticity (FILE * fp, char *label, int axis, double curtime)
{

  int i;

  int element;

  double omega = 0.0; 

  double *se;

  msh_vector gradu, gradv, gradw;

  fprintf (fp, "$View\n");

  se = calloc (nbelements, sizeof (double));

  for (i = 0; i < nbelements; i++)
    {

      element = i;

      gradu = Gradient (&xul, &xuf, LOGICAL_FALSE, element);
      gradv = Gradient (&xvl, &xvf, LOGICAL_FALSE, element);
      gradw = Gradient (&xwl, &xwf, LOGICAL_FALSE, element);

      if (axis == 0)
	{
	  omega = gradw.y - gradv.z;
	}

      if (axis == 1)
	{
	  omega = gradu.z - gradw.x;
	}

      if (axis == 2)
	{
	  omega = gradv.x - gradu.y;
	}

      se[element] = omega;

    }

  PrintBinaryHeaderScalarElement (fp, label);
  PrintBinaryTime (fp, curtime);
  PrintBinaryScalarElement (fp, se);

  fprintf (fp, "\n$EndView\n");

  free (se);

}

void
WriteResults (FILE * fp, int onFace, int inElement, double curtime)
{

  VecGhostGetLocalForm (xu, &xul);
  VecGhostGetLocalForm (xv, &xvl);
  VecGhostGetLocalForm (xw, &xwl);
  VecGhostGetLocalForm (xp, &xpl);
  VecGhostGetLocalForm (xT, &xTl);
  VecGhostGetLocalForm (xs, &xsl);

  if (parameter.savflux == LOGICAL_TRUE)
    {
      if (parameter.wbinary == LOGICAL_TRUE)
	{
	  sprintf (viewname, "Flux_[%s/%s]-%03d", parameter.ulength,
		   parameter.utime, processor);
	  WriteBinaryFlux (fp, viewname, curtime);
	}
      else
	{
	  sprintf (viewname, "Flux_[%s/%s]-%03d", parameter.ulength,
		   parameter.utime, processor);
	  WriteAsciiFlux (fp, viewname, curtime);
	}
    }

  if (onFace == LOGICAL_TRUE)
    {

      if (parameter.wbinary == LOGICAL_TRUE)
	{

	  // Write velocity-U to file 
	  sprintf (viewname, "VelocityU(Face)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.fsav[iu] == LOGICAL_TRUE)
	    WriteBinaryScalarFace (fp, viewname, iu, curtime);

	  // Write velocity-V to file 
	  sprintf (viewname, "VelocityV(Face)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.fsav[iv] == LOGICAL_TRUE)
	    WriteBinaryScalarFace (fp, viewname, iv, curtime);

	  // Write velocity-W to file 
	  sprintf (viewname, "VelocityW(Face)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.fsav[iw] == LOGICAL_TRUE)
	    WriteBinaryScalarFace (fp, viewname, iw, curtime);

	  // Write pressure to file
	  sprintf (viewname, "Pressure(Face)_[%s/(%s.%s.%s)]-%03d",
		   parameter.umass, parameter.ulength, parameter.utime,
		   parameter.utime, processor);
	  if (parameter.fsav[ip] == LOGICAL_TRUE)
	    WriteBinaryScalarFace (fp, viewname, ip, curtime);

	  // Write temperature to file
	  sprintf (viewname, "Temperature(Face)_[%s]-%03d",
		   parameter.utemperature, processor);
	  if (parameter.fsav[iT] == LOGICAL_TRUE)
	    WriteBinaryScalarFace (fp, viewname, iT, curtime);

	  // Write indicator function to file
	  sprintf (viewname, "Gamma(Face)_[-]-%03d", processor);
	  if (parameter.fsav[is] == LOGICAL_TRUE)
	    WriteBinaryScalarFace (fp, viewname, is, curtime);

	  // Write velocity to file 
	  sprintf (viewname, "Velocity(Face)_[%s/%s]-%03d", parameter.ulength,
		   parameter.utime, processor);
	  if (parameter.fvec == LOGICAL_TRUE)
	    WriteBinaryVectorMagnitudeFace (fp, viewname, iu, iv, iw,
					    curtime);

	}
      else
	{

	  // Write velocity-U to file 
	  sprintf (viewname, "VelocityU(Face)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.fsav[iu] == LOGICAL_TRUE)
	    WriteAsciiScalarFace (fp, viewname, iu, curtime);

	  // Write velocity-V to file 
	  sprintf (viewname, "VelocityV(Face)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.fsav[iv] == LOGICAL_TRUE)
	    WriteAsciiScalarFace (fp, viewname, iv, curtime);

	  // Write velocity-W to file 
	  sprintf (viewname, "VelocityW(Face)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.fsav[iw] == LOGICAL_TRUE)
	    WriteAsciiScalarFace (fp, viewname, iw, curtime);

	  // Write pressure to file
	  sprintf (viewname, "Pressure(Face)_[%s/(%s.%s.%s)]-%03d",
		   parameter.umass, parameter.ulength, parameter.utime,
		   parameter.utime, processor);
	  if (parameter.fsav[ip] == LOGICAL_TRUE)
	    WriteAsciiScalarFace (fp, viewname, ip, curtime);

	  // Write temperature to file
	  sprintf (viewname, "Temperature(Face)_[%s]-%03d",
		   parameter.utemperature, processor);
	  if (parameter.fsav[iT] == LOGICAL_TRUE)
	    WriteAsciiScalarFace (fp, viewname, iT, curtime);

	  // Write indicator function to file
	  sprintf (viewname, "Gamma(Face)_[-]-%03d", processor);
	  if (parameter.fsav[is] == LOGICAL_TRUE)
	    WriteAsciiScalarFace (fp, viewname, is, curtime);

	  // Write velocity to file 
	  sprintf (viewname, "Velocity(Face)_[%s/%s]-%03d", parameter.ulength,
		   parameter.utime, processor);
	  if (parameter.fvec == LOGICAL_TRUE)
	    WriteBinaryVectorMagnitudeFace (fp, viewname, iu, iv, iw,
					    curtime);

	}
    }

  if (inElement == LOGICAL_TRUE)
    {

      if (parameter.wbinary == LOGICAL_TRUE)
	{

	  // Write velocity-U to file 
	  sprintf (viewname, "VelocityU(Cell)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.csav[iu] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, iu, curtime);

	  // Write velocity-V to file 
	  sprintf (viewname, "VelocityV(Cell)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.csav[iv] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, iv, curtime);

	  // Write velocity-W to file 
	  sprintf (viewname, "VelocityW(Cell)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.csav[iw] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, iw, curtime);

	  // Write pressure to file
	  sprintf (viewname, "Pressure(Cell)_[%s/(%s.%s.%s)]-%03d",
		   parameter.umass, parameter.ulength, parameter.utime,
		   parameter.utime, processor);
	  if (parameter.csav[ip] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, ip, curtime);

	  // Write indicator function to file
	  sprintf (viewname, "Temperature(Cell)_[%s]-%03d",
		   parameter.utemperature, processor);
	  if (parameter.csav[iT] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, iT, curtime);

	  // Write indicator function to file
	  sprintf (viewname, "Gamma(Cell)_[-]-%03d", processor);
	  if (parameter.csav[is] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, is, curtime);

	  // Write velocity to file 
	  sprintf (viewname, "Velocity(Cell)_[%s/%s]-%03d", parameter.ulength,
		   parameter.utime, processor);
	  if (parameter.cvec == LOGICAL_TRUE)
	    WriteBinaryVectorElement (fp, viewname, iu, iv, iw, curtime);

	  // Write vorticity to file                      
	  sprintf (viewname, "VorticityU(Cell)_[1/%s]-%03d", parameter.utime,
		   processor);
	  if (parameter.vortex[0] == LOGICAL_TRUE)
	    WriteBinaryVorticity (fp, viewname, 0, curtime);

	  sprintf (viewname, "VorticityV(Cell)_[1/%s]-%03d", parameter.utime,
		   processor);
	  if (parameter.vortex[1] == LOGICAL_TRUE)
	    WriteBinaryVorticity (fp, viewname, 1, curtime);

	  sprintf (viewname, "VorticityW(Cell)_[1/%s]-%03d", parameter.utime,
		   processor);
	  if (parameter.vortex[2] == LOGICAL_TRUE)
	    WriteBinaryVorticity (fp, viewname, 2, curtime);

	  sprintf (viewname, "StreamFunction(Cell)_[-]-%03d", processor);
	  if (parameter.streamf == LOGICAL_TRUE)
	    WriteBinaryStreamFunction (fp, viewname, curtime);

	}
      else
	{

	  // Write velocity-U to file 
	  sprintf (viewname, "VelocityU(Cell)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.csav[iu] == LOGICAL_TRUE)
	    WriteAsciiScalarElement (fp, viewname, iu, curtime);

	  // Write velocity-V to file 
	  sprintf (viewname, "VelocityV(Cell)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.csav[iv] == LOGICAL_TRUE)
	    WriteAsciiScalarElement (fp, viewname, iv, curtime);

	  // Write velocity-W to file                     
	  sprintf (viewname, "VelocityW(Cell)_[%s/%s]-%03d",
		   parameter.ulength, parameter.utime, processor);
	  if (parameter.csav[iw] == LOGICAL_TRUE)
	    WriteAsciiScalarElement (fp, viewname, iw, curtime);

	  // Write pressure to file
	  sprintf (viewname, "Pressure(Cell)_[%s/(%s.%s.%s)]-%03d",
		   parameter.umass, parameter.ulength, parameter.utime,
		   parameter.utime, processor);
	  if (parameter.csav[ip] == LOGICAL_TRUE)
	    WriteAsciiScalarElement (fp, viewname, ip, curtime);

	  // Write indicator function to file
	  sprintf (viewname, "Temperature(Cell)_[%s]-%03d",
		   parameter.utemperature, processor);
	  if (parameter.csav[iT] == LOGICAL_TRUE)
	    WriteAsciiScalarElement (fp, viewname, iT, curtime);

	  // Write indicator function to file
	  sprintf (viewname, "Gamma(Cell)_[-]-%03d", processor);
	  if (parameter.csav[is] == LOGICAL_TRUE)
	    WriteAsciiScalarElement (fp, viewname, is, curtime);

	  // Write velocity to file 
	  sprintf (viewname, "Velocity(Cell)_[%s/%s]-%03d", parameter.ulength,
		   parameter.utime, processor);
	  if (parameter.cvec == LOGICAL_TRUE)
	    WriteAsciiVectorElement (fp, viewname, iu, iv, iw, curtime);

	  // Write vorticity to file                      
	  sprintf (viewname, "VorticityU(Cell)_[1/s]-%03d", processor);
	  if (parameter.vortex[0] == LOGICAL_TRUE)
	    WriteAsciiVorticity (fp, viewname, 0, curtime);

	  sprintf (viewname, "VorticityV(Cell)_[1/s]-%03d", processor);
	  if (parameter.vortex[1] == LOGICAL_TRUE)
	    WriteAsciiVorticity (fp, viewname, 1, curtime);

	  sprintf (viewname, "VorticityW(Cell)_[1/s]-%03d", processor);
	  if (parameter.vortex[2] == LOGICAL_TRUE)
	    WriteAsciiVorticity (fp, viewname, 2, curtime);

	  sprintf (viewname, "StreamFunction(Cell)_[-]-%03d", processor);
	  if (parameter.streamf == LOGICAL_TRUE)
	    WriteAsciiStreamFunction (fp, viewname, curtime);

	}

    }

  VecGhostRestoreLocalForm (xu, &xul);
  VecGhostRestoreLocalForm (xv, &xvl);
  VecGhostRestoreLocalForm (xw, &xwl);
  VecGhostRestoreLocalForm (xp, &xpl);
  VecGhostRestoreLocalForm (xT, &xTl);
  VecGhostRestoreLocalForm (xs, &xsl);

}

void
WriteResidual (FILE * fp, int iter, double *res)
{

  int i;

  PetscFPrintf (PETSC_COMM_WORLD, fp, "%d", iter);

  for (i = 0; i < nphi; i++)
    {
      PetscFPrintf (PETSC_COMM_WORLD, fp, " \t %E", res[i]);
    }

  PetscFPrintf (PETSC_COMM_WORLD, fp, "\n");

}