UltraScan III
List of all members | Classes | Static Public Member Functions
US_AstfemMath Class Reference

A group of static mathematical functions to support finite element calculations. More...

#include "us_astfem_math.h"

Classes

class  AstFemParameters
 Parameters for finite element solution. More...
 
class  ComponentRole
 Component Role. More...
 
class  MfemData
 A data set comprised of scans from one sample taken at constant speed. More...
 
class  MfemInitial
 Initial concentration conditions. More...
 
class  MfemScan
 A scan entry. More...
 
class  ReactionGroup
 Reaction Group. More...
 

Static Public Member Functions

static void interpolate_C0 (MfemInitial &, MfemInitial &)
 Interpolate first onto second. More...
 
static void interpolate_C0 (MfemInitial &, double *, QVector< double > &)
 Interpolate starting concentration QVector MfemInitial onto C0. More...
 
static void zero_2d (int, int, double **)
 Initialize a 2d matrix in memory to all zeros. More...
 
static void initialize_2d (int, int, double ***)
 Create a 2d matrix in memory and initilize to all zeros. More...
 
static void clear_2d (int, double **)
 Delete a 2d matrix in memory. More...
 
static double maxval (const QVector< double > &)
 Find the maximum value in a vector. More...
 
static double minval (const QVector< double > &)
 Find the minimum value in a vector. More...
 
static double maxval (const QVector< US_Model::SimulationComponent > &)
 Find the maximum s in a vector of SimulationComponent entries. More...
 
static double minval (const QVector< US_Model::SimulationComponent > &)
 Find the minimum s in a vector of SimulationComponent entries. More...
 
static void initialize_3d (int, int, int, double ****)
 Create a 3d matrix in memory and initilize to all zeros. More...
 
static void clear_3d (int, int, double ***)
 Delete a 3d matrix in memory. More...
 
static void tridiag (double *, double *, double *, double *, double *, int)
 Solve a Ax = b where A is tridiagonal. More...
 
static double cube_root (double, double, double)
 Find the positive cubic-root of a cubic polynomial
with a0 <= 0 and
a1, a2 >= 0. More...
 
static int GaussElim (int, double **, double *)
 Solve Ax = b using Gaussian Elimination. More...
 
static double find_C1_mono_Nmer (int, double, double)
 Solve f(x) = x + K * x^n - C using Newton's method.
This function needs to be renamed! More...
 
static int interpolate (MfemData &, MfemData &, bool)
 Interpolate one dataset onto another using time or omega^2t. More...
 
static int interpolate (MfemData &, MfemData &, bool, int, int)
 Interpolate one dataset onto another using time or omega^2t. More...
 
static void QuadSolver (double *, double *, double *, double *, double *, double *, int)
 Solve Quad-diagonal system. More...
 
static void IntQT1 (double *, double, double, double **, double)
 Integration on test function seperately on left Q, right T. More...
 
static void IntQTm (double *, double, double, double **, double)
 Integration on test function. More...
 
static void IntQTn2 (double *, double, double, double **, double)
 Integration on test function. More...
 
static void IntQTn1 (double *, double, double, double **, double)
 Integration on test function. More...
 
static void DefineFkp (int, double **)
 Define Lamm equation values. More...
 
static double AreaT (double *, double *)
 Compute the area of a triangle (v1, v2, v3) More...
 
static void BasisTS (double, double, double *, double *, double *)
 Computer basis on standard element (TS) More...
 
static void BasisQS (double, double, double *, double *, double *)
 Computer basis on standard element (QS) More...
 
static void BasisTR (double *, double *, double, double, double *, double *, double *)
 Computer basis on real element T at given(xs,ts) point. More...
 
static void BasisQR (double *, double, double, double *, double *, double *, double)
 Computer basis on real element Q at given(xs,ts) point. More...
 
static double Integrand (double, double, double, double, double, double, double, double)
 Integrand for Lamm equation. More...
 
static void DefineGaussian (int, double **)
 Define Gaussian. More...
 
static void initSimData (US_DataIO::RawData &, US_DataIO::EditedData &, double)
 Initialize simulation data from experimental. More...
 
static double variance (US_DataIO::RawData &, US_DataIO::EditedData &)
 Calculate variance for Simulation-Experimental difference. More...
 
static double variance (US_DataIO::RawData &, US_DataIO::EditedData &, QList< int >)
 Calculate variance for Simulation-Experimental difference. More...
 
static double calc_bottom (double, double, double *)
 Calculate bottom radius from channel bottom and rotor coefficients. More...
 

Detailed Description

A group of static mathematical functions to support finite element calculations.

Definition at line 13 of file us_astfem_math.h.

Member Function Documentation

double US_AstfemMath::AreaT ( double *  xv,
double *  yv 
)
static

Compute the area of a triangle (v1, v2, v3)

Parameters
xvThe XV vector
yvThe YV vector

Definition at line 1684 of file us_astfem_math.cpp.

void US_AstfemMath::BasisQR ( double *  vx,
double  xs,
double  ts,
double *  phi,
double *  phix,
double *  phiy,
double  dt 
)
static

Computer basis on real element Q at given(xs,ts) point.

Parameters
vxThe VX vector
xsThe start X coordinate
tsThe start Y coordinate
phiThe Phi vector to fill
phixThe PhiX vector to fill
phiyThe PhiY vector to fill
dtThe d-t constant

Definition at line 1790 of file us_astfem_math.cpp.

void US_AstfemMath::BasisQS ( double  xi,
double  et,
double *  phi,
double *  phi1,
double *  phi2 
)
static

Computer basis on standard element (QS)

Parameters
xiThe XI constant
etThe ET constant
phiThe Phi vector to fill
phi1The Phi1 vector to fill
phi2The Phi2 vector to fill

Definition at line 1712 of file us_astfem_math.cpp.

void US_AstfemMath::BasisTR ( double *  vx,
double *  vy,
double  xs,
double  ys,
double *  phi,
double *  phix,
double *  phiy 
)
static

Computer basis on real element T at given(xs,ts) point.

Parameters
vxThe Vx vector
vyThe Vy vector
xsThe start X coordinate
ysThe start Y coordinate
phiThe Phi vector to fill
phixThe PhiX vector to fill
phiyThe PhiY vector to fill

Definition at line 1737 of file us_astfem_math.cpp.

void US_AstfemMath::BasisTS ( double  xi,
double  et,
double *  phi,
double *  phi1,
double *  phi2 
)
static

Computer basis on standard element (TS)

Parameters
xiThe XI constant
etThe ET constant
phiThe Phi vector to fill
phi1The Phi1 vector to fill
phi2The Phi2 vector to fill

Definition at line 1692 of file us_astfem_math.cpp.

double US_AstfemMath::calc_bottom ( double  rpm,
double  bottom_chan,
double *  rotorcoefs 
)
static

Calculate bottom radius from channel bottom and rotor coefficients.

Parameters
rpmRotor revolutions per minute.
bottom_chanInitial bottom for centerpiece channel
rotorcoefsArray of 2 rotor coefficients
Returns
The calculated bottom radius value.

Definition at line 2015 of file us_astfem_math.cpp.

void US_AstfemMath::clear_2d ( int  val1,
double **  matrix 
)
static

Delete a 2d matrix in memory.

Parameters
val1First dimension
matrixval1 x n matrix to be deleted

Definition at line 109 of file us_astfem_math.cpp.

void US_AstfemMath::clear_3d ( int  val1,
int  val2,
double ***  matrix 
)
static

Delete a 3d matrix in memory.

Parameters
val1First dimension
val2Second dimension
matrixDeleted val1 x val2 x n matrix

Definition at line 179 of file us_astfem_math.cpp.

double US_AstfemMath::cube_root ( double  a0,
double  a1,
double  a2 
)
static

Find the positive cubic-root of a cubic polynomial
with a0 <= 0 and
a1, a2 >= 0.

Parameters
a0The a0 value
a1The a0 value
a2The a0 value
Returns
Cube root

Definition at line 240 of file us_astfem_math.cpp.

void US_AstfemMath::DefineFkp ( int  npts,
double **  Lam 
)
static

Define Lamm equation values.

Parameters
nptsOrder of the equation
LamMatrix of Lamm values to fill

Definition at line 1407 of file us_astfem_math.cpp.

void US_AstfemMath::DefineGaussian ( int  nGauss,
double **  Gs2 
)
static

Define Gaussian.

Parameters
nGaussOrder of the matrix
Gs2Matrix to be filled

Definition at line 1840 of file us_astfem_math.cpp.

double US_AstfemMath::find_C1_mono_Nmer ( int  n,
double  K,
double  CT 
)
static

Solve f(x) = x + K * x^n - C using Newton's method.
This function needs to be renamed!

Parameters
nThe power exponent
KThe scalar constant
CTThe test C value

Definition at line 285 of file us_astfem_math.cpp.

int US_AstfemMath::GaussElim ( int  n,
double **  a,
double *  b 
)
static

Solve Ax = b using Gaussian Elimination.

Parameters
nThe order of the square matrix
aThe n x n matrix
bThe b vector
Returns
Success flag: -1 -> singular, no solution; 1 -> success

Definition at line 326 of file us_astfem_math.cpp.

void US_AstfemMath::initialize_2d ( int  val1,
int  val2,
double ***  matrix 
)
static

Create a 2d matrix in memory and initilize to all zeros.

Parameters
val1First dimension
val2Second dimension
matrixInitialized val1 x val2 matrix

Definition at line 95 of file us_astfem_math.cpp.

void US_AstfemMath::initialize_3d ( int  val1,
int  val2,
int  val3,
double ****  matrix 
)
static

Create a 3d matrix in memory and initilize to all zeros.

Parameters
val1First dimension
val2Second dimension
val3Third dimension
matrixInitialized val1 x val2 x val3 matrix

Definition at line 158 of file us_astfem_math.cpp.

void US_AstfemMath::initSimData ( US_DataIO::RawData simdata,
US_DataIO::EditedData editdata,
double  concval1 = 0.0 
)
static

Initialize simulation data from experimental.

Parameters
simdataReference to simulation Raw Data to initialize.
editdataReference to experimental Edited Data to mirror.
concval1Optional constant concentration value for first scan.

Definition at line 1913 of file us_astfem_math.cpp.

double US_AstfemMath::Integrand ( double  x_gauss,
double  D,
double  sw2,
double  u,
double  ux,
double  ut,
double  v,
double  vx 
)
static

Integrand for Lamm equation.

Parameters
x_gaussThe X-gaussian constant
DThe D constant
sw2The SW2 constant
uThe U constant
uxThe Ux constant
utThe Ut constant
vThe V constant
vxThe Vx constant

Definition at line 1831 of file us_astfem_math.cpp.

int US_AstfemMath::interpolate ( MfemData expdata,
MfemData simdata,
bool  use_time 
)
static

Interpolate one dataset onto another using time or omega^2t.

Parameters
expdataExperimental data to create, sized on input
simdataSimulation from which to create modeled experiment
use_timeFlag of whether to use time interpolation
Returns
Success flag: 0 -> success

Definition at line 390 of file us_astfem_math.cpp.

int US_AstfemMath::interpolate ( MfemData expdata,
MfemData simdata,
bool  use_time,
int  fscan,
int  lscan 
)
static

Interpolate one dataset onto another using time or omega^2t.

Parameters
expdataExperimental data to create, sized on input
simdataSimulation from which to create modeled experiment
use_timeFlag of whether to use time correction
fscanFirst update expdata scan index
lscanLast update expdata scan index plus one
Returns
Success flag: 0 -> success

Definition at line 466 of file us_astfem_math.cpp.

void US_AstfemMath::interpolate_C0 ( MfemInitial C0,
MfemInitial C1 
)
static

Interpolate first onto second.

Parameters
C0Input MfemInitial
C1MfemInitial with interpolated concentrations

Definition at line 48 of file us_astfem_math.cpp.

void US_AstfemMath::interpolate_C0 ( MfemInitial C0,
double *  C1,
QVector< double > &  xvec 
)
static

Interpolate starting concentration QVector MfemInitial onto C0.

Parameters
C0Input MfemInitial
C1First scan with interpolated concentrations
xvecStart x (radius) values for each scan

Definition at line 6 of file us_astfem_math.cpp.

void US_AstfemMath::IntQT1 ( double *  vx,
double  D,
double  sw2,
double **  Stif,
double  dt 
)
static

Integration on test function seperately on left Q, right T.

Parameters
vxThe vx vector
DThe D value
sw2The sw2 value
StifThe Stif matrix
dtThe dt value

Definition at line 846 of file us_astfem_math.cpp.

void US_AstfemMath::IntQTm ( double *  vx,
double  D,
double  sw2,
double **  Stif,
double  dt 
)
static

Integration on test function.

Parameters
vxThe vx vector
DThe D value
sw2The sw2 value
StifThe Stif matrix
dtThe dt value

Definition at line 990 of file us_astfem_math.cpp.

void US_AstfemMath::IntQTn1 ( double *  vx,
double  D,
double  sw2,
double **  Stif,
double  dt 
)
static

Integration on test function.

Parameters
vxThe vx vector
DThe D value
sw2The sw2 value
StifThe Stif matrix
dtThe dt value

Definition at line 1332 of file us_astfem_math.cpp.

void US_AstfemMath::IntQTn2 ( double *  vx,
double  D,
double  sw2,
double **  Stif,
double  dt 
)
static

Integration on test function.

Parameters
vxThe vx vector
DThe D value
sw2The sw2 value
StifThe Stif matrix
dtThe dt value

Definition at line 1161 of file us_astfem_math.cpp.

double US_AstfemMath::maxval ( const QVector< double > &  value)
static

Find the maximum value in a vector.

Parameters
valueVector whose maximum is found
Returns
Vector maximum value

Definition at line 137 of file us_astfem_math.cpp.

double US_AstfemMath::maxval ( const QVector< US_Model::SimulationComponent > &  value)
static

Find the maximum s in a vector of SimulationComponent entries.

Parameters
valueVector of components whose maximum s value is found
Returns
Vector maximum value

Definition at line 148 of file us_astfem_math.cpp.

double US_AstfemMath::minval ( const QVector< double > &  value)
static

Find the minimum value in a vector.

Parameters
valueVector whose minimum is found
Returns
Vector minimum value

Definition at line 116 of file us_astfem_math.cpp.

double US_AstfemMath::minval ( const QVector< US_Model::SimulationComponent > &  value)
static

Find the minimum s in a vector of SimulationComponent entries.

Parameters
valueVector of components whose minimum s value is found
Returns
Vector minimum value

Definition at line 127 of file us_astfem_math.cpp.

void US_AstfemMath::QuadSolver ( double *  ai,
double *  bi,
double *  ci,
double *  di,
double *  cr,
double *  solu,
int  N 
)
static

Solve Quad-diagonal system.

Parameters
aiThe initial a vector
biThe initial b vector
ciThe initial c vector
diThe initial d vector
crThe Cr vector
soluThe calculated solution vector
NThe length of vectors

Definition at line 780 of file us_astfem_math.cpp.

void US_AstfemMath::tridiag ( double *  a,
double *  b,
double *  c,
double *  r,
double *  u,
int  N 
)
static

Solve a Ax = b where A is tridiagonal.

Parameters
aArray of a values
bArray of b values
cArray of c values
rArray of r values
uArray of u values
NLength of vectors

Definition at line 194 of file us_astfem_math.cpp.

double US_AstfemMath::variance ( US_DataIO::RawData simdata,
US_DataIO::EditedData editdata 
)
static

Calculate variance for Simulation-Experimental difference.

Parameters
simdataReference to simulation Raw Data.
editdataReference to experimental Edited Data.
Returns
The variance (average of differences squared) between the simulated and experimental readings values.

Definition at line 1962 of file us_astfem_math.cpp.

double US_AstfemMath::variance ( US_DataIO::RawData simdata,
US_DataIO::EditedData editdata,
QList< int >  exclScans 
)
static

Calculate variance for Simulation-Experimental difference.

Parameters
simdataReference to simulation Raw Data.
editdataReference to experimental Edited Data.
exclScansList of excluded scans.
Returns
The variance (average of differences squared) between the simulated and experimental readings values.

Definition at line 1972 of file us_astfem_math.cpp.

void US_AstfemMath::zero_2d ( int  val1,
int  val2,
double **  matrix 
)
static

Initialize a 2d matrix in memory to all zeros.

Parameters
val1First dimension
val2Second dimension
matrixInitialized val1 x val2 matrix

Definition at line 88 of file us_astfem_math.cpp.


The documentation for this class was generated from the following files: