#include <DGMSolver.hpp>

Static Public Member Functions
static unsigned int	nchoosek (unsigned int n, unsigned int k)
	compute combinational number, n choose k, maximum output is std::numeric_limits<unsigned int>::max(); More...

static unsigned int	compute_numcols_vander_multivar (unsigned int kvars, unsigned int degree)
	compute the number of columns for a multivariate vandermonde matrix, given certen degree More...

static void	gen_multivar_monomial_basis (const int kvars, const double *vars, const int degree, std::vector< double > &basis)
	compute the monomial basis of mutiple variables, up to input degree, lexicographically ordered More...

static void	gen_vander_multivar (const int mrows, const int kvars, const double *us, const int degree, std::vector< double > &V)
	compute multivariate vandermonde matrix, monomial basis ordered in the same way as gen_multivar_monomial_basis More...

static void	rescale_matrix (int mrows, int ncols, double V, double ts)

static void	compute_qtransposeB (int mrows, int ncols, const double Q, int bncols, double bs)

static void	qr_polyfit_safeguarded (const int mrows, const int ncols, double V, double D, int *rank)

static void	backsolve (int mrows, int ncols, double R, int bncols, double bs, double *ws)

static void	backsolve_polyfit_safeguarded (int dim, int degree, const bool interp, int mrows, int ncols, double R, int bncols, double bs, const double ws, int degree_out)

static void	vec_dotprod (const int len, const double a, const double b, double *c)

static void	vec_scalarprod (const int len, const double a, const double c, double b)

static void	vec_crossprod (const double a[3], const double b[3], double(&c)[3])

static double	vec_innerprod (const int len, const double a, const double b)

static double	vec_2norm (const int len, const double *a)

static double	vec_normalize (const int len, const double a, double b)

static double	vec_distance (const int len, const double a, const double b)

static void	vec_projoff (const int len, const double a, const double b, double *c)

static void	vec_linear_operation (const int len, const double mu, const double a, const double psi, const double b, double *c)

static void	get_tri_natural_coords (const int dim, const double cornercoords, const int npts, const double currcoords, double *naturalcoords)

Private Member Functions
	DGMSolver ()

	~DGMSolver ()

Detailed Description

Definition at line 8 of file DGMSolver.hpp.

Constructor & Destructor Documentation

◆ DGMSolver()

moab::DGMSolver::DGMSolver ( )

inlineprivate

Definition at line 10 of file DGMSolver.hpp.

10 {};

◆ ~DGMSolver()

moab::DGMSolver::~DGMSolver ( )

inlineprivate

Definition at line 11 of file DGMSolver.hpp.

11 {};

Member Function Documentation

◆ backsolve()

void moab::DGMSolver::backsolve	(	int	mrows,
		int	ncols,
		double *	R,
		int	bncols,
		double *	bs,
		double *	ws
	)

static

Definition at line 256 of file DGMSolver.cpp.

 {
 #if 0
     std::cout.precision(16);
     std::cout<<"Before backsolve  "<<std::endl;
     std::cout<<" V = "<<std::endl;
     for (int k=0; k< ncols; k++){
         for (int j=0; j<mrows; ++j){
             std::cout<<R[mrows*k+j]<<std::endl;
           }
         std::cout<<std::endl;
       }
     std::cout<<std::endl;
  
     //std::cout<<"#pnts = "<<mrows<<std::endl;
     std::cout<<"bs = "<<std::endl;
     std::cout<<std::endl;
     for (int k=0; k< bncols; k++){
         for (int j=0; j<mrows; ++j){
             std::cout<<"  "<<bs[mrows*k+j]<<std::endl;
           }
       }
     std::cout<<std::endl;
 #endif
  
     for( int k = 0; k < bncols; k++ )
     {
         for( int j = ncols - 1; j >= 0; j-- )
         {
             for( int i = j + 1; i < ncols; ++i )
                 bs[mrows * k + j] = bs[mrows * k + j] - R[mrows * i + j] * bs[mrows * k + i];
  
             assert( R[mrows * j + j] != 0 );
             bs[mrows * k + j] = bs[mrows * k + j] / R[mrows * j + j];
         }
     }
  
     for( int k = 0; k < bncols; k++ )
     {
         for( int j = 0; j < ncols; ++j )
             bs[mrows * k + j] = bs[mrows * k + j] / ws[j];
     }
 }

References moab::R.

Referenced by moab::HiReconstruction::eval_vander_bivar_cmf(), and moab::HiReconstruction::eval_vander_univar_cmf().

◆ backsolve_polyfit_safeguarded()

void moab::DGMSolver::backsolve_polyfit_safeguarded	(	int	dim,
		int	degree,
		const bool	interp,
		int	mrows,
		int	ncols,
		double *	R,
		int	bncols,
		double *	bs,
		const double *	ws,
		int *	degree_out
	)

static

Definition at line 300 of file DGMSolver.cpp.

 {
     if( ncols < 1 ) std::cout << "ERROR: Invalid input to safeguarded polyfit backsolve routine.\n";
 #if 0
     std::cout.precision(12);
     std::cout<<"Before backsolve  "<<std::endl;
     std::cout<<" V = "<<std::endl;
     for (int k=0; k< ncols; k++){
         for (int j=0; j<mrows; ++j){
             std::cout<<R[mrows*k+j]<<std::endl;
           }
         std::cout<<std::endl;
       }
     std::cout<<std::endl;
  
     //std::cout<<"#pnts = "<<mrows<<std::endl;
     std::cout<<"bs = "<<std::endl;
     std::cout<<std::endl;
     for (int k=0; k< bncols; k++){
         for (int j=0; j<mrows; ++j){
             std::cout<<"  "<<bs[mrows*k+j]<<std::endl;
         }
     }
     std::cout<<std::endl;
  
     //std::cout<<" ] "<<std::endl;
  
     std::cout<<"Input ws = [ ";
     for (int k=0; k< ncols; k++){
             std::cout<<ws[k]<<", ";
           }
     std::cout<<" ] "<<std::endl;
  
     std::cout << "R: " << R << "size: [" << mrows << "," << ncols << "]" << std::endl;
     std::cout << "bs: " << bs << "size: [" << mrows << "," << bncols << "]" << std::endl;
     std::cout << "ws: " << ws << "size: [" << ncols << "," << 1 << "]" << std::endl;
     std::cout << "degree_out: " << degree_out << std::endl;
 #endif
  
     int deg, numcols = 0;
  
     for( int k = 0; k < bncols; k++ )
     {
         deg = degree;
         /* I think we should consider interp = true/false -Xinglin*/
         if( dim == 1 )
             numcols = deg + 1 - interp;
         else if( dim == 2 )
             numcols = ( deg + 2 ) * ( deg + 1 ) / 2 - interp;
  
         assert( numcols <= ncols );
  
         std::vector< double > bs_bak( numcols );
  
         if( deg >= 2 )
         {
             for( int i = 0; i < numcols; i++ )
             {
                 assert( mrows * k + i < mrows * bncols );
                 bs_bak.at( i ) = bs[mrows * k + i];
             }
         }
  
         while( deg >= 1 )
         {
             int cend       = numcols - 1;
             bool downgrade = false;
  
             // The reconstruction can be applied only on edges (2-d) or faces (3-d)
             assert( cend >= 0 );
             assert( dim > 0 && dim < 3 );
  
             for( int d = deg; d >= 0; d-- )
             {
                 int cstart = 0;
                 if( dim == 1 )
                 {
                     cstart = d;
                 }
                 else if( dim == 2 )
                 {
                     cstart = ( ( d + 1 ) * d ) / 2;
                     // cstart = ((d*(d+1))>>1)-interp;
                 }
  
                 // Solve for  bs
                 for( int j = cend; j >= cstart; j-- )
                 {
                     assert( mrows * k + j < mrows * bncols );
                     for( int i = j + 1; i < numcols; ++i )
                     {
                         assert( mrows * k + i < mrows * bncols );
                         assert( mrows * i + j < mrows * ncols );  // check R
                         bs[mrows * k + j] = bs[mrows * k + j] - R[mrows * i + j] * bs[mrows * k + i];
                     }
                     assert( mrows * j + j < mrows * ncols );  // check R
                     bs[mrows * k + j] = bs[mrows * k + j] / R[mrows * j + j];
                 }
  
                 // Checking for change in the coefficient
                 if( d >= 2 && d < deg )
                 {
                     double tol;
  
                     if( dim == 1 )
                     {
                         tol = 1e-06;
                         assert( mrows * cstart + cstart < mrows * ncols );  // check R
                         double tb = bs_bak.at( cstart ) / R[mrows * cstart + cstart];
                         assert( mrows * k + cstart < mrows * bncols );
                         if( fabs( bs[mrows * k + cstart] - tb ) > ( 1 + tol ) * fabs( tb ) )
                         {
                             downgrade = true;
                             break;
                         }
                     }
  
                     else if( dim == 2 )
                     {
                         tol = 0.05;
  
                         std::vector< double > tb( cend - cstart + 1 );
                         for( int j = 0; j <= ( cend - cstart ); j++ )
                         {
                             tb.at( j ) = bs_bak.at( cstart + j );
                         }
  
                         for( int j = cend; j >= cstart; j-- )
                         {
                             int jind = j - cstart;
  
                             for( int i = j + 1; i <= cend; ++i )
                             {
                                 assert( mrows * i + j < mrows * ncols );  // check R
                                 tb.at( jind ) = tb.at( jind ) - R[mrows * i + j] * tb.at( i - cstart );
                             }
                             assert( mrows * j + j < mrows * ncols );  // check R
                             tb.at( jind ) = tb.at( jind ) / R[mrows * j + j];
                             assert( mrows * k + j < mrows * bncols );
                             double err = fabs( bs[mrows * k + j] - tb.at( jind ) );
  
                             if( ( err > tol ) && ( err >= ( 1 + tol ) * fabs( tb.at( jind ) ) ) )
                             {
                                 downgrade = true;
                                 break;
                             }
                         }
  
                         if( downgrade ) break;
                     }
                 }
  
                 cend = cstart - 1;
             }
  
             if( !downgrade )
                 break;
             else
             {
                 deg = deg - 1;
                 if( dim == 1 )
                     numcols = deg + 1;
                 else if( dim == 2 )
                     numcols = ( deg + 2 ) * ( deg + 1 ) / 2;
  
                 for( int i = 0; i < numcols; i++ )
                 {
                     assert( mrows * k + i < mrows * bncols );
                     bs[mrows * k + i] = bs_bak.at( i );
                 }
             }
         }
         assert( k < bncols );
         degree_out[k] = deg;
  
         for( int i = 0; i < numcols; i++ )
         {
             // assert(mrows*k+i < mrows*bncols);
             // assert(i < ncols);
             bs[mrows * k + i] = bs[mrows * k + i] / ws[i];
         }
  
         for( int i = numcols; i < mrows; i++ )
         {
             // assert(mrows*k+i < mrows*bncols);
             bs[mrows * k + i] = 0;
         }
     }
 }

References dim, and moab::R.

Referenced by moab::HiReconstruction::eval_vander_bivar_cmf(), and moab::HiReconstruction::eval_vander_univar_cmf().

◆ compute_numcols_vander_multivar()

unsigned int moab::DGMSolver::compute_numcols_vander_multivar	(	unsigned int	kvars,
		unsigned int	degree
	)

static

compute the number of columns for a multivariate vandermonde matrix, given certen degree

If degree = 0, out put is 1; If kvars = 1, degree = k, output is k+1; If kvars = 2, degree = k, output is (k+2)*(k+1)/2;

Definition at line 42 of file DGMSolver.cpp.

 {
     unsigned int mcols = 0;
     for( unsigned int i = 0; i <= degree; ++i )
     {
         unsigned int temp = nchoosek( kvars - 1 + i, kvars - 1 );
         if( !temp )
         {
             std::cout << "overflow to compute nchoosek n= " << kvars - 1 + i << " k= " << kvars - 1 << std::endl;
             return 0;
         }
         mcols += temp;
     }
     return mcols;
 }

References nchoosek().

Referenced by gen_multivar_monomial_basis(), and gen_vander_multivar().

◆ compute_qtransposeB()

void moab::DGMSolver::compute_qtransposeB	(	int	mrows,
		int	ncols,
		const double *	Q,
		int	bncols,
		double *	bs
	)

static

Definition at line 174 of file DGMSolver.cpp.

 {
     for( int k = 0; k < ncols; k++ )
     {
         for( int j = 0; j < bncols; j++ )
         {
             double t2 = 0;
             for( int i = k; i < mrows; i++ )
                 t2 += Q[mrows * k + i] * bs[mrows * j + i];
             t2 = t2 + t2;
  
             for( int i = k; i < mrows; i++ )
                 bs[mrows * j + i] -= t2 * Q[mrows * k + i];
         }
     }
 }

Referenced by moab::HiReconstruction::eval_vander_bivar_cmf(), and moab::HiReconstruction::eval_vander_univar_cmf().

◆ gen_multivar_monomial_basis()

void moab::DGMSolver::gen_multivar_monomial_basis	(	const int	kvars,
		const double *	vars,
		const int	degree,
		std::vector< double > &	basis
	)

static

compute the monomial basis of mutiple variables, up to input degree, lexicographically ordered

if degree = 0, output basis = {1} If kvars = 1, vars = {u}, degree = k, basis = {1,u,...,u^k} If kvars = 2, vars = {u,v}, degree = k, basis = {1,u,v,u^2,uv,v^2,u^3,u^2*v,uv^2,v^3,...,u^k,u^k-1*v,...,uv^k-1,v^k} If kvars = 3, vars = {u,v,w}, degree = k, basis = {1,u,v,w,u^2,uv,uw,v^2,v*w,w^2,...,u^k,u^k-1v,u^k-1w,...,v^k,v^k-1w,...,vw^k-1,w^k}

Parameters

kvars	Integer, number of variables
vars	Pointer to array of doubles, size = kvars, variable values
degree	Integer, maximum degree
basis	Reference to vector, user input container to hold output which is appended to existing data; users don't have to preallocate for basis, this function will allocate interally

Definition at line 58 of file DGMSolver.cpp.

 {
     unsigned int len = compute_numcols_vander_multivar( kvars, degree );
     basis.reserve( len - basis.capacity() + basis.size() );
     size_t iend = basis.size();
 #ifndef NDEBUG
     size_t istr = basis.size();
 #endif
     basis.push_back( 1 );
     ++iend;
     if( !degree )
     {
         return;
     }
     std::vector< size_t > varspos( kvars );
     // degree 1
     for( int ivar = 0; ivar < kvars; ++ivar )
     {
         basis.push_back( vars[ivar] );
         varspos[ivar] = iend++;
     }
     // degree 2 to degree
     for( int ideg = 2; ideg <= degree; ++ideg )
     {
         size_t preend = iend;
         for( int ivar = 0; ivar < kvars; ++ivar )
         {
             size_t varpreend = iend;
             for( size_t ilast = varspos[ivar]; ilast < preend; ++ilast )
             {
                 basis.push_back( vars[ivar] * basis[ilast] );
                 ++iend;
             }
             varspos[ivar] = varpreend;
         }
     }
     assert( len == iend - istr );
 }

References compute_numcols_vander_multivar().

◆ gen_vander_multivar()

void moab::DGMSolver::gen_vander_multivar	(	const int	mrows,
		const int	kvars,
		const double *	us,
		const int	degree,
		std::vector< double > &	V
	)

static

compute multivariate vandermonde matrix, monomial basis ordered in the same way as gen_multivar_monomial_basis

if degree = 0, V = {1,...,1}'; If kvars = 1, us = {u1;u2;..,;um}, degree = k, V = {1 u1 u1^2 ... u1^k;1 u2 u2^2 ... u2^k;...;1 um um^2 ... um^k}; *If kvars = 2, us = {u1 v1;u2 v2;...;um vm}, degree = k, V = {1 u1 v1 u1^2 u1v1 v1^2;...;1 um vm um^2 umvm vm^2};

Parameters

mrows	Integer, number of points to evaluate Vandermonde matrix
kvars	Integer, number of variables
us	Pointer to array of doubles, size = mrow*kvars, variable values for all points. Stored in row-wise, like {u1 v1 u2 v2 ...}
degree	Integer, maximum degree
basis	Reference to vector, user input container to hold Vandermonde matrix which is appended to existing data; users don't have to preallocate for basis, this function will allocate interally; the Vandermonde matrix is stored in an array, columnwise, like {1 ... 1 u1 ...um u1^2 ... um^2 ...}

Definition at line 100 of file DGMSolver.cpp.

 {
     unsigned int ncols = compute_numcols_vander_multivar( kvars, degree );
     V.reserve( mrows * ncols - V.capacity() + V.size() );
     size_t istr = V.size(), icol = 0;
     // add ones, V is stored in an single array, elements placed in columnwise order
     for( int irow = 0; irow < mrows; ++irow )
     {
         V.push_back( 1 );
     }
     ++icol;
     if( !degree )
     {
         return;
     }
     std::vector< size_t > varspos( kvars );
     // degree 1
     for( int ivar = 0; ivar < kvars; ++ivar )
     {
         for( int irow = 0; irow < mrows; ++irow )
         {
             V.push_back( us[irow * kvars + ivar] );  // us stored in row-wise
         }
         varspos[ivar] = icol++;
     }
     // from 2 to degree
     for( int ideg = 2; ideg <= degree; ++ideg )
     {
         size_t preendcol = icol;
         for( int ivar = 0; ivar < kvars; ++ivar )
         {
             size_t varpreend = icol;
             for( size_t ilast = varspos[ivar]; ilast < preendcol; ++ilast )
             {
                 for( int irow = 0; irow < mrows; ++irow )
                 {
                     V.push_back( us[irow * kvars + ivar] * V[istr + irow + ilast * mrows] );
                 }
                 ++icol;
             }
             varspos[ivar] = varpreend;
         }
     }
     assert( icol == ncols );
 }

References compute_numcols_vander_multivar().

Referenced by moab::HiReconstruction::eval_vander_bivar_cmf(), and moab::HiReconstruction::eval_vander_univar_cmf().

◆ get_tri_natural_coords()

void moab::DGMSolver::get_tri_natural_coords	(	const int	dim,
		const double *	cornercoords,
		const int	npts,
		const double *	currcoords,
		double *	naturalcoords
	)

static

Definition at line 667 of file DGMSolver.cpp.

 {
     assert( dim == 2 || dim == 3 );
     double a = 0, b = 0, d = 0, tol = 1e-12;
     for( int i = 0; i < dim; ++i )
     {
         a += ( cornercoords[dim + i] - cornercoords[i] ) * ( cornercoords[dim + i] - cornercoords[i] );
         b += ( cornercoords[dim + i] - cornercoords[i] ) * ( cornercoords[2 * dim + i] - cornercoords[i] );
         d += ( cornercoords[2 * dim + i] - cornercoords[i] ) * ( cornercoords[2 * dim + i] - cornercoords[i] );
     }
     double det = a * d - b * b;
     assert( det > 0 );
     for( int ipt = 0; ipt < npts; ++ipt )
     {
         double e = 0, f = 0;
         for( int i = 0; i < dim; ++i )
         {
             e += ( cornercoords[dim + i] - cornercoords[i] ) * ( currcoords[ipt * dim + i] - cornercoords[i] );
             f += ( cornercoords[2 * dim + i] - cornercoords[i] ) * ( currcoords[ipt * dim + i] - cornercoords[i] );
         }
         naturalcoords[ipt * 3 + 1] = ( e * d - b * f ) / det;
         naturalcoords[ipt * 3 + 2] = ( a * f - b * e ) / det;
         naturalcoords[ipt * 3]     = 1 - naturalcoords[ipt * 3 + 1] - naturalcoords[ipt * 3 + 2];
         if( naturalcoords[ipt * 3] < -tol || naturalcoords[ipt * 3 + 1] < -tol || naturalcoords[ipt * 3 + 2] < -tol )
         {
             std::cout << "Corners: \n";
             std::cout << cornercoords[0] << "\t" << cornercoords[1] << "\t" << cornercoords[3] << std::endl;
             std::cout << cornercoords[3] << "\t" << cornercoords[4] << "\t" << cornercoords[5] << std::endl;
             std::cout << cornercoords[6] << "\t" << cornercoords[7] << "\t" << cornercoords[8] << std::endl;
             std::cout << "Candidate: \n";
             std::cout << currcoords[ipt * dim] << "\t" << currcoords[ipt * dim + 1] << "\t" << currcoords[ipt * dim + 2]
                       << std::endl;
             exit( 0 );
         }
         assert( fabs( naturalcoords[ipt * 3] + naturalcoords[ipt * 3 + 1] + naturalcoords[ipt * 3 + 2] - 1 ) < tol );
         for( int i = 0; i < dim; ++i )
         {
             assert( fabs( naturalcoords[ipt * 3] * cornercoords[i] +
                           naturalcoords[ipt * 3 + 1] * cornercoords[dim + i] +
                           naturalcoords[ipt * 3 + 2] * cornercoords[2 * dim + i] - currcoords[ipt * dim + i] ) < tol );
         }
     }
 }

References dim.

◆ nchoosek()

unsigned int moab::DGMSolver::nchoosek	(	unsigned int	n,
		unsigned int	k
	)

static

compute combinational number, n choose k, maximum output is std::numeric_limits<unsigned int>::max();

Definition at line 19 of file DGMSolver.cpp.

 {
     if( k > n )
     {
         return 0;
     }
     unsigned long long ans = 1;
     if( k > ( n >> 1 ) )
     {
         k = n - k;
     }
     for( unsigned int i = 1; i <= k; ++i )
     {
         ans *= n--;
         ans /= i;
         if( ans > std::numeric_limits< unsigned int >::max() )
         {
             return 0;
         }
     }
     return ans;
 }

Referenced by compute_numcols_vander_multivar().

◆ qr_polyfit_safeguarded()

void moab::DGMSolver::qr_polyfit_safeguarded	(	const int	mrows,
		const int	ncols,
		double *	V,
		double *	D,
		int *	rank
	)

static

Definition at line 191 of file DGMSolver.cpp.

 {
     double tol = 1e-8;
     *rank      = ncols;
     double* v  = new double[mrows];
  
     for( int k = 0; k < ncols; k++ )
     {
         int nv = mrows - k;
  
         for( int j = 0; j < nv; j++ )
             v[j] = V[mrows * k + ( j + k )];
  
         double t2 = 0;
  
         for( int j = 0; j < nv; j++ )
             t2 = t2 + v[j] * v[j];
  
         double t    = sqrt( t2 );
         double vnrm = 0;
  
         if( v[0] >= 0 )
         {
             vnrm = sqrt( 2 * ( t2 + v[0] * t ) );
             v[0] = v[0] + t;
         }
         else
         {
             vnrm = sqrt( 2 * ( t2 - v[0] * t ) );
             v[0] = v[0] - t;
         }
  
         if( vnrm > 0 )
         {
             for( int j = 0; j < nv; j++ )
                 v[j] = v[j] / vnrm;
         }
  
         for( int j = k; j < ncols; j++ )
         {
             t2 = 0;
             for( int i = 0; i < nv; i++ )
                 t2 = t2 + v[i] * V[mrows * j + ( i + k )];
  
             t2 = t2 + t2;
  
             for( int i = 0; i < nv; i++ )
                 V[mrows * j + ( i + k )] = V[mrows * j + ( i + k )] - t2 * v[i];
         }
  
         D[k] = V[mrows * k + k];
  
         for( int i = 0; i < nv; i++ )
             V[mrows * k + ( i + k )] = v[i];
  
         if( ( fabs( D[k] ) ) < tol && ( *rank == ncols ) )
         {
             *rank = k;
             break;
         }
     }
  
     delete[] v;
 }

Referenced by moab::HiReconstruction::eval_vander_bivar_cmf(), and moab::HiReconstruction::eval_vander_univar_cmf().

◆ rescale_matrix()

void moab::DGMSolver::rescale_matrix	(	int	mrows,
		int	ncols,
		double *	V,
		double *	ts
	)

static

Definition at line 150 of file DGMSolver.cpp.

 {
     // This function rescales the input matrix using the norm of each column.
     double* v = new double[mrows];
     for( int i = 0; i < ncols; i++ )
     {
         for( int j = 0; j < mrows; j++ )
             v[j] = V[mrows * i + j];
  
         // Compute norm of the column vector
         double w = vec_2norm( mrows, v );
  
         if( fabs( w ) == 0 )
             ts[i] = 1;
         else
         {
             ts[i] = w;
             for( int j = 0; j < mrows; j++ )
                 V[mrows * i + j] = V[mrows * i + j] / ts[i];
         }
     }
     delete[] v;
 }

References vec_2norm().

Referenced by moab::HiReconstruction::eval_vander_bivar_cmf(), and moab::HiReconstruction::eval_vander_univar_cmf().

◆ vec_2norm()

double moab::DGMSolver::vec_2norm	(	const int	len,
		const double *	a
	)

static

Definition at line 544 of file DGMSolver.cpp.

 {
     if( !a )
     {
         MB_SET_ERR_RET_VAL( "NULL Pointer", 0.0 );
     }
     double w = 0, s = 0;
     for( int k = 0; k < len; k++ )
         w = std::max( w, fabs( a[k] ) );
  
     if( w == 0 )
     {
         return 0;
     }
     else
     {
         for( int k = 0; k < len; k++ )
         {
             s += ( a[k] / w ) * ( a[k] / w );
         }
         s = w * sqrt( s );
     }
     return s;
 }

References MB_SET_ERR_RET_VAL.

Referenced by rescale_matrix(), and vec_projoff().

◆ vec_crossprod()

void moab::DGMSolver::vec_crossprod	(	const double	a[3],
		const double	b[3],
		double(&)	c[3]
	)

static

Definition at line 523 of file DGMSolver.cpp.

 {
     c[0] = a[1] * b[2] - a[2] * b[1];
     c[1] = a[2] * b[0] - a[0] * b[2];
     c[2] = a[0] * b[1] - a[1] * b[0];
 }

Referenced by moab::HiReconstruction::average_vertex_normal(), and moab::HiReconstruction::polyfit3d_surf_get_coeff().

◆ vec_distance()

double moab::DGMSolver::vec_distance	(	const int	len,
		const double *	a,
		const double *	b
	)

static

Definition at line 604 of file DGMSolver.cpp.

 {
     double res = 0;
     for( int i = 0; i < len; ++i )
     {
         res += ( a[i] - b[i] ) * ( a[i] - b[i] );
     }
     return sqrt( res );
 }

◆ vec_dotprod()

void moab::DGMSolver::vec_dotprod	(	const int	len,
		const double *	a,
		const double *	b,
		double *	c
	)

static

Definition at line 499 of file DGMSolver.cpp.

 {
     if( !a || !b || !c )
     {
         MB_SET_ERR_RET( "NULL Pointer" );
     }
     for( int i = 0; i < len; ++i )
     {
         c[i] = a[i] * b[i];
     }
 }

References MB_SET_ERR_RET.

◆ vec_innerprod()

double moab::DGMSolver::vec_innerprod	(	const int	len,
		const double *	a,
		const double *	b
	)

static

Definition at line 530 of file DGMSolver.cpp.

 {
     if( !a || !b )
     {
         MB_SET_ERR_RET_VAL( "NULL Pointer", 0.0 );
     }
     double ans = 0;
     for( int i = 0; i < len; ++i )
     {
         ans += a[i] * b[i];
     }
     return ans;
 }

References MB_SET_ERR_RET_VAL.

Referenced by moab::HiReconstruction::compute_weights(), moab::HiReconstruction::polyfit3d_curve_get_coeff(), moab::HiReconstruction::polyfit3d_surf_get_coeff(), vec_projoff(), moab::HiReconstruction::walf3d_curve_vertex_eval(), and moab::HiReconstruction::walf3d_surf_vertex_eval().

◆ vec_linear_operation()

void moab::DGMSolver::vec_linear_operation	(	const int	len,
		const double	mu,
		const double *	a,
		const double	psi,
		const double *	b,
		double *	c
	)

static

Definition at line 650 of file DGMSolver.cpp.

 {
     if( !a || !b || !c )
     {
         MB_SET_ERR_RET( "NULL Pointer" );
     }
     for( int i = 0; i < len; ++i )
     {
         c[i] = mu * a[i] + psi * b[i];
     }
 }

References MB_SET_ERR_RET.

Referenced by moab::HiReconstruction::average_vertex_normal(), moab::HiReconstruction::average_vertex_tangent(), moab::HiReconstruction::polyfit3d_curve_get_coeff(), moab::HiReconstruction::polyfit3d_surf_get_coeff(), and moab::HiReconstruction::walf3d_curve_vertex_eval().

◆ vec_normalize()

double moab::DGMSolver::vec_normalize	(	const int	len,
		const double *	a,
		double *	b
	)

static

Definition at line 569 of file DGMSolver.cpp.

 {
     if( !a || !b )
     {
         MB_SET_ERR_RET_VAL( "NULL Pointer", 0.0 );
     }
     double nrm = 0, mx = 0;
     for( int i = 0; i < len; ++i )
     {
         mx = std::max( fabs( a[i] ), mx );
     }
     if( mx == 0 )
     {
         for( int i = 0; i < len; ++i )
         {
             b[i] = 0;
         }
         return 0;
     }
     for( int i = 0; i < len; ++i )
     {
         nrm += ( a[i] / mx ) * ( a[i] / mx );
     }
     nrm = mx * sqrt( nrm );
     if( nrm == 0 )
     {
         return nrm;
     }
     for( int i = 0; i < len; ++i )
     {
         b[i] = a[i] / nrm;
     }
     return nrm;
 }

References MB_SET_ERR_RET_VAL.

Referenced by moab::HiReconstruction::average_vertex_normal(), moab::HiReconstruction::average_vertex_tangent(), and moab::HiReconstruction::polyfit3d_surf_get_coeff().

◆ vec_projoff()

void moab::DGMSolver::vec_projoff	(	const int	len,
		const double *	a,
		const double *	b,
		double *	c
	)

static

Definition at line 614 of file DGMSolver.cpp.

 {
     if( !a || !b || !c )
     {
         MB_SET_ERR_RET( "NULL Pointer" );
     }
     // c = a-<a,b>b/<b,b>;
     double bnrm = vec_2norm( len, b );
     if( bnrm == 0 )
     {
         for( int i = 0; i < len; ++i )
         {
             c[i] = a[i];
         }
         return;
     }
     double innerp = vec_innerprod( len, a, b ) / bnrm;
  
     if( innerp == 0 )
     {
         if( c != a )
         {
             for( int i = 0; i < len; ++i )
             {
                 c[i] = a[i];
             }
         }
         return;
     }
  
     for( int i = 0; i < len; ++i )
     {
         c[i] = a[i] - innerp * b[i] / bnrm;
     }
 }

References MB_SET_ERR_RET, vec_2norm(), and vec_innerprod().

Referenced by moab::HiReconstruction::polyfit3d_surf_get_coeff().

◆ vec_scalarprod()

void moab::DGMSolver::vec_scalarprod	(	const int	len,
		const double *	a,
		const double	c,
		double *	b
	)

static

Definition at line 511 of file DGMSolver.cpp.

 {
     if( !a || !b )
     {
         MB_SET_ERR_RET( "NULL Pointer" );
     }
     for( int i = 0; i < len; ++i )
     {
         b[i] = c * a[i];
     }
 }

References MB_SET_ERR_RET.

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

Static Public Member Functions

Private Member Functions

Detailed Description

Constructor & Destructor Documentation

◆ DGMSolver()

◆ ~DGMSolver()

Member Function Documentation

◆ backsolve()

◆ backsolve_polyfit_safeguarded()

◆ compute_numcols_vander_multivar()

◆ compute_qtransposeB()

◆ gen_multivar_monomial_basis()

◆ gen_vander_multivar()

◆ get_tri_natural_coords()

◆ nchoosek()

◆ qr_polyfit_safeguarded()

◆ rescale_matrix()

◆ vec_2norm()

◆ vec_crossprod()

◆ vec_distance()

◆ vec_dotprod()

◆ vec_innerprod()

◆ vec_linear_operation()

◆ vec_normalize()

◆ vec_projoff()

◆ vec_scalarprod()