NavierStokes3DJacobian.c File Reference

Contains the functions compute flux Jacobians for the 3D Navier-Stokes system. More...

#include <arrayfunctions.h>
#include <mathfunctions.h>
#include <matmult_native.h>
#include <physicalmodels/navierstokes3d.h>

Go to the source code of this file.

Functions
int	NavierStokes3DJacobian (double *Jac, double *u, void *p, int dir, int nvars, int upw)
int	NavierStokes3DStiffJacobian (double *Jac, double *u, void *p, int dir, int nvars, int upw)

Detailed Description

Contains the functions compute flux Jacobians for the 3D Navier-Stokes system.

Author

Debojyoti Ghosh

Definition in file NavierStokes3DJacobian.c.

Function Documentation

NavierStokes3DJacobian()

int NavierStokes3DJacobian	(	double *	Jac,
		double *	u,
		void *	p,
		int	dir,
		int	nvars,
		int	upw
	)

Function to compute the flux Jacobian of the 3D Navier-Stokes equations, given the solution at a grid point. The Jacobian is square matrix of size nvar=5, and is returned as a 1D array (double) of 25 elements in row-major format.

Parameters

Jac	Jacobian matrix: 1D array of size nvar^2 = 25
u	solution at a grid point (array of size nvar = 5)
p	object containing the physics-related parameters
dir	dimension (0 -> x, 1 -> y, 2 -> z)
nvars	number of vector components
upw	0 -> send back complete Jacobian, 1 -> send back Jacobian of right(+)-moving flux, -1 -> send back Jacobian of left(-)-moving flux

Definition at line 15 of file NavierStokes3DJacobian.c.

{
  NavierStokes3D *param = (NavierStokes3D*) p;
  static double R[_MODEL_NVARS_*_MODEL_NVARS_], D[_MODEL_NVARS_*_MODEL_NVARS_], 
                L[_MODEL_NVARS_*_MODEL_NVARS_], DL[_MODEL_NVARS_*_MODEL_NVARS_];

  /* get the eigenvalues and left,right eigenvectors */
  _NavierStokes3DEigenvalues_      (u,_NavierStokes3D_stride_,D,param->gamma,dir);
  _NavierStokes3DLeftEigenvectors_ (u,_NavierStokes3D_stride_,L,param->gamma,dir);
  _NavierStokes3DRightEigenvectors_(u,_NavierStokes3D_stride_,R,param->gamma,dir);

  int aupw = absolute(upw), k;
  k = 0;  D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  k = 6;  D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  k = 12; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  k = 18; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  k = 24; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );

  MatMult5(_MODEL_NVARS_,DL,D,L);
  MatMult5(_MODEL_NVARS_,Jac,R,DL);

  return(0);
}

NavierStokes3DStiffJacobian()

int NavierStokes3DStiffJacobian	(	double *	Jac,
		double *	u,
		void *	p,
		int	dir,
		int	nvars,
		int	upw
	)

Function to compute the Jacobian of the fast flux (representing the acoustic waves) of the 3D Navier-Stokes equations, given the solution at a grid point. The Jacobian is square matrix of size nvar=5, and is returned as a 1D array (double) of 25 elements in row-major format.

Parameters

Jac	Jacobian matrix: 1D array of size nvar^2 = 25
u	solution at a grid point (array of size nvar = 5)
p	object containing the physics-related parameters
dir	dimension (0 -> x, 1 -> y, 2 -> z)
nvars	number of vector components
upw	0 -> send back complete Jacobian, 1 -> send back Jacobian of right(+)-moving flux, -1 -> send back Jacobian of left(-)-moving flux

Definition at line 53 of file NavierStokes3DJacobian.c.

{
  NavierStokes3D *param = (NavierStokes3D*) p;
  static double R[_MODEL_NVARS_*_MODEL_NVARS_], D[_MODEL_NVARS_*_MODEL_NVARS_], 
                L[_MODEL_NVARS_*_MODEL_NVARS_], DL[_MODEL_NVARS_*_MODEL_NVARS_];

  /* get the eigenvalues and left,right eigenvectors */
  _NavierStokes3DEigenvalues_      (u,_NavierStokes3D_stride_,D,param->gamma,dir);
  _NavierStokes3DLeftEigenvectors_ (u,_NavierStokes3D_stride_,L,param->gamma,dir);
  _NavierStokes3DRightEigenvectors_(u,_NavierStokes3D_stride_,R,param->gamma,dir);

  int aupw = absolute(upw), k;
  if (dir == _XDIR_) {
    k = 0;  D[k] = 0.0;
    k = 6;  D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
    k = 12; D[k] = 0.0;
    k = 18; D[k] = 0.0;
    k = 24; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  } else if (dir == _YDIR_) {
    k = 0;  D[k] = 0.0;
    k = 6;  D[k] = 0.0;
    k = 12; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
    k = 18; D[k] = 0.0;
    k = 24; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  } else if (dir == _ZDIR_) {
    k = 0;  D[k] = 0.0;
    k = 6;  D[k] = 0.0;
    k = 12; D[k] = 0.0;
    k = 18; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
    k = 24; D[k] = absolute( (1-aupw)*D[k] + 0.5*aupw*(1+upw)*max(0,D[k]) + 0.5*aupw*(1-upw)*min(0,D[k]) );
  }

  MatMult5(_MODEL_NVARS_,DL,D,L);
  MatMult5(_MODEL_NVARS_,Jac,R,DL);

  return(0);
}