post.c

OpenFVM-v1.1 open source cfd code

  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 + 1);
      if (iphix == 1)
	vx = V_GetCmp (&xvf, face + 1);
      if (iphix == 2)
	vx = V_GetCmp (&xwf, face + 1);
      if (iphix == 3)
	vx = V_GetCmp (&xpf, face + 1);
      if (iphix == 4)
	vx = V_GetCmp (&xTf, face + 1);
      if (iphix == 5)
	vx = V_GetCmp (&xsf, face + 1);

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

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

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

    }

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

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

  fprintf (fp, "$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 (&xu, element + 1);
      if (iphi == 1)
	vs = V_GetCmp (&xv, element + 1);
      if (iphi == 2)
	vs = V_GetCmp (&xw, element + 1);
      if (iphi == 3)
	vs = V_GetCmp (&xp, element + 1);
      if (iphi == 4)
	vs = V_GetCmp (&xT, element + 1);
      if (iphi == 5)
	vs = V_GetCmp (&xs, element + 1);

      se[element] = vs;

    }

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

  fprintf (fp, "$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, "$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, "$EndView\n");

  free (sn);
}

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

  int i;

  int element;

  double omega;

  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 (&xu, &xuf, LOGICAL_FALSE, element);
      gradv = Gradient (&xv, &xvf, LOGICAL_FALSE, element);
      gradw = Gradient (&xw, &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, "$EndView\n");

  free (se);

}

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

  int i, j;

  int element, face;

  double vs;

  double *se;

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

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

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

      element = i;

      vs = 0.0;

      for (j = 0; j < elements[element].nbfaces; j++)
	{

	  face = elements[element].face[j];

	  vs += V_GetCmp (&uf, face + 1) * faces[face].Aj;

	}

      se[element] = vs;

    }

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

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

  free (se);

}

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

  int i, j;

  int element, face;

  double vs;

  double *se;

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

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

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

      element = i;

      vs = 0.0;

      for (j = 0; j < elements[element].nbfaces; j++)
	{

	  face = elements[element].face[j];

	  vs += V_GetCmp (&uf, face + 1) * faces[face].Aj;

	}

      se[element] = vs;

    }

  PrintAsciiHeaderScalarElement (fp, label);
  PrintAsciiTime (fp, curtime);
  PrintAsciiScalarElement (fp, se);

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

  free (se);

}

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

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

  if (onFace == LOGICAL_TRUE)
    {

      if (parameter.wbinary == LOGICAL_TRUE)
	{

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

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

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

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

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

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

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

	}
      else
	{

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

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

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

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

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

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

	  // Write velocity to file 
	  sprintf (viewname, "Velocity(Face)_[%s/%s]", parameter.ulength,
		   parameter.utime);
	  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]", parameter.ulength,
		   parameter.utime);
	  if (parameter.csav[iu] == LOGICAL_TRUE)
	    WriteBinaryScalarElement (fp, viewname, iu, curtime);

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

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

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

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

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

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

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

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

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

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

	}
      else
	{

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

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

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

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

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

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

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

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

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

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

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

	}

    }



}

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

  int i;

  fprintf (fp, "%d", iter);

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

  fprintf (fp, "\n");

}