math_util.h

Go to the documentation of this file.

/*
 * Copyright (c) 2009-2015:  G-CSC, Goethe University Frankfurt
 * Author: Sebastian Reiter
 * 
 * This file is part of UG4.
 * 
 * UG4 is free software: you can redistribute it and/or modify it under the
 * terms of the GNU Lesser General Public License version 3 (as published by the
 * Free Software Foundation) with the following additional attribution
 * requirements (according to LGPL/GPL v3 §7):
 * 
 * (1) The following notice must be displayed in the Appropriate Legal Notices
 * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
 * 
 * (2) The following notice must be displayed at a prominent place in the
 * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
 * 
 * (3) The following bibliography is recommended for citation and must be
 * preserved in all covered files:
 * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
 *   parallel geometric multigrid solver on hierarchically distributed grids.
 *   Computing and visualization in science 16, 4 (2013), 151-164"
 * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
 *   flexible software system for simulating pde based models on high performance
 *   computers. Computing and visualization in science 16, 4 (2013), 165-179"
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
 * GNU Lesser General Public License for more details.
 */
 
#ifndef __H__UGMATH__MATH_UTIL__
#define __H__UGMATH__MATH_UTIL__
 
#include <cstdlib>
#include "math_constants.h"
#include "../math_vector_matrix/math_vector.h"
#include "../math_vector_matrix/math_vector_functions.h"
#include "../ugmath_types.h"
#include "../../ug_config.h"
#include "eigenvalues.h"
 
// this is quite hacky...
// some environments (on windows) seem to define a small type or variable.
// this causes errors with parameter names etc. It is thus undefined below.
#ifdef small
  #undef small
#endif
 
namespace ug
{
 
//  deg_to_rad
template <class TNumber>
inline TNumber
deg_to_rad(TNumber deg);
 
//  rad_to_deg
template <class TNumber>
inline TNumber
rad_to_deg(TNumber rad);
 
//  urand
template <class TNumber>
TNumber
urand(TNumber lowerBound, TNumber upperBound);
 
//  clip
template <class TNumber>
TNumber
clip(TNumber val, TNumber lowerBound, TNumber upperBound);
 
template <class TNumber>
inline TNumber sq(TNumber val);
 
UG_API
bool FindNormal(vector3& normOut, const vector3& v);
 
 
UG_API
bool ConstructOrthonormalSystem(matrix33& matOut, const vector3& v,
                size_t vColInd);
 
template <class vector_t>
void CalculateCenter(vector_t& centerOut, const vector_t* pointSet,
           size_t numPoints);
 
template <class vector_t>
vector_t
TriangleBarycenter(const vector_t& p1, const vector_t& p2, const vector_t& p3);
 
 
UG_API
bool TriangleCircumcenter(vector2& centerOut, const vector2& p1,
              const vector2& p2, const vector2& p3);
 
UG_API
bool TriangleCircumcenter(vector3& centerOut, const vector3& p1,
              const vector3& p2, const vector3& p3);
 
UG_API
void CalculateCovarianceMatrix(matrix33& matOut, const vector3* pointSet,
                const vector3& center, size_t numPoints);
 
 
UG_API
bool FindClosestPlane(vector3& centerOut, vector3& normalOut,
            const vector3* pointSet, size_t numPoints);
 
 
UG_API
bool TransformPointSetTo2D(vector2* pointSetOut, const vector3* pointSet,
              size_t numPoints);
 
 
template <class vector_t>
number DropAPerpendicular(vector_t& vOut, const vector_t& v,
              const vector_t& v0, const vector_t& v1);
 
 
template <class vector_t>
vector_t PointOnRay(const vector_t& from, const vector_t& dir, number s);
 
 
template <class vector_t>
number ProjectPointToRay(vector_t& vOut, const vector_t& v,
              const vector_t& from, const vector_t& dir);
 
 
 
template <class vector_t>
number ProjectPointToLine(vector_t& vOut, const vector_t& v,
              const vector_t& from, const vector_t& to);
 
 
/*
 * the coordinates of the point are given through v.
 * the line-segment is defined by its endpoints v1 and v2.
 *
 * The method returns the distance.
 */
template <class vector_t>
inline
number DistancePointToLine(const vector_t& v, const vector_t& v1,
              const vector_t& v2);
 
template <class vector_t>
number DistancePointToLine(number& tOut, const vector_t& v,
               const vector_t& v1, const vector_t& v2);
              
/*
 * the coordinates of the point are given through v.
 * the ray is defined by from (an arbitrary point on the ray) and
 * dir (the direction of the ray).
 *
 * The method returns the distance.
 */
template <class vector_t>
inline
number DistancePointToRay(const vector_t& v, const vector_t& from,
              const vector_t& dir);
 
template <class vector_t>
inline
number DistancePointToRay(number& tOut, const vector_t& v,
              const vector_t& from, const vector_t& dir);
 
 
template <class vector_t>
number DistancePointToTriangle(vector_t& vOut, number& bc1Out, number& bc2Out,
              const vector_t& p, const vector_t& v1, const vector_t& v2,
              const vector_t& v3, const vector_t& n);
 
template <class vector_t>
number DistancePointToPlane(const vector_t& v, const vector_t& p,
              const vector_t& n);
 
 
template <class vector_t>
void ProjectPointToPlane(vector_t& vOut, const vector_t& v,
            const vector_t& p, const vector_t& n);
 
//  RayPlaneIntersection
template <class vector_t>
bool RayPlaneIntersection(vector_t& vOut, number& tOut,
              const vector_t& rayFrom, const vector_t& rayDir,
              const vector_t& p, const vector_t& n);
 
//  RayRayIntersection2d
 
template <class vector_t>
bool RayRayIntersection2d(vector_t &vOut, number& t0Out, number& t1Out,
               const vector_t &p0, const vector_t &dir0,
               const vector_t &p1, const vector_t &dir1);
 
//  RayLineIntersection2d
 
template <class vector_t>
bool RayLineIntersection2d(vector_t &vOut, number& bcOut, number& tOut,
               const vector_t &p0, const vector_t &p1,
               const vector_t &vFrom, const vector_t &vDir,
               number sml = 0);
 
 
//  LineLineIntersection2d
 
template <class vector_t>
bool LineLineIntersection2d(vector_t &vOut, number& t0Out, number& t1Out,
               const vector_t &from0, const vector_t &to0,
               const vector_t &from1, const vector_t &to1,
               const number threshold = 0);
 
 
 
UG_API
bool RayRayIntersection3d(vector3& aOut, vector3& bOut,
              const vector3& p0, const vector3& dir0,
              const vector3& p1, const vector3& dir1);
 
 
UG_API
bool LineLineIntersection3d(vector3& aOut, vector3& bOut,
              const vector3& a1, const vector3& a2,
                const vector3& b1, const vector3& b2);
 
 
template <class vector_t>
bool RayRayProjection(number& t1Out, number& t2Out,
            const vector_t& from1, const vector_t& dir1,
            const vector_t& from2, const vector_t& dir2);
 
 
template <class vector_t>
bool LineLineProjection(number& t1Out, number& t2Out,
              const vector_t& a1, const vector_t& a2,
              const vector_t& b1, const vector_t& b2);
 
 
UG_API
number DistanceLineToLine(const vector3& a1, const vector3& a2,
              const vector3& b1, const vector3& b2);
 
//  RayTriangleIntersection
 
template <class vector_t>
bool RayTriangleIntersection(vector_t &vOut, number& bc1Out, number& bc2Out, number& tOut,
               const vector_t &p0, const vector_t &p1, const vector_t &p2, 
               const vector_t &vFrom, const vector_t &vDir,
               const number small = SMALL);
 
//  RayTriangleIntersection
 
template <class vector_t> inline
bool RayTriangleIntersection(vector_t &vOut, const vector_t &p0,
               const vector_t &p1, const vector_t &p2, 
               const vector_t &vFrom, const vector_t &vDir);
 
//  RayBoxIntersection
 
template <class vector_t>
bool RayBoxIntersection(const vector_t& rayFrom, const vector_t& rayDir,
            const vector_t& boxMin, const vector_t& boxMax,
            number* tNearOut = NULL, number* tFarOut = NULL);
 
 
template <class vector_t>
bool LineBoxIntersection(const vector_t& v1, const vector_t& v2,
            const vector_t& boxMin, const vector_t& boxMax);
 
 
template <class vector_t>
int RaySphereIntersection(number& s1Out, number& s2Out,
              const vector_t& v, const vector_t& dir,
              const vector_t& center, number radius);
 
template <class vector_t>
int LineSphereIntersection(number& s1Out, number& s2Out,
              const vector_t& v1, const vector_t& v2,
              const vector_t& center, number radius);
 
 
UG_API
bool RayCylinderIntersection(number& tMinOut, number& tMaxOut, const vector3& rayFrom,
               const vector3& rayDir, const vector3& cylCenter,
               const vector3& cylAxis, number cylRadius);
 
//  TriangleTriangleIntersection
 
UG_API
bool TriangleTriangleIntersection(const MathVector<3>& p0, const MathVector<3>& p1,
                  const MathVector<3>& p2, const MathVector<3>& q0,
                  const MathVector<3>& q1, const MathVector<3>& q2,
                  MathVector<3>* ip1Out = NULL,
                  MathVector<3>* ip2Out = NULL,
                  number snapThreshold = SMALL);
 
//  TriBoxIntersection
//  opposed to most of the other methods here, TriBoxIntersection is not
//  a template method. Instead it directly operates on the built in
//  vector3 type. This is due to its implementation, which is taken
//  from 'Game Programming Gems'. Please have a look at the copyright notice
//  in 'tri_box.cpp'.
 
UG_API
bool TriangleBoxIntersection(const MathVector<3>& p0, const MathVector<3>& p1,
               const MathVector<3>& p2,
               const MathVector<3>& boxMin, const MathVector<3>& boxMax);
            
 
template <class vector_t>
bool BoxBoxIntersection(const vector_t& box1Min, const vector_t& box1Max,
            const vector_t& box2Min, const vector_t& box2Max);
 
//  BoxBoundProbe
 
template <class vector_t>
bool BoxBoundProbe(const vector_t& v, const vector_t& boxMin,
          const vector_t& boxMax);
 
template <class vector_t>
number TriangleArea(const vector_t& p1, const vector_t& p2,
          const vector_t& p3);
 
template <class vector_t>
number QuadrilateralArea(const vector_t& p1, const vector_t& p2,
             const vector_t& p3, const vector_t& p4);
 
 
template <class vector_t>
number GeometricApproximationDegree(vector_t& n1, vector_t& n2, vector_t& n3,
                  vector_t& tn);
 
 
template <class vector_t>
number TriangleQuality_Area(const vector_t& p1, const vector_t& p2,
              const vector_t& p3);
 
 
template <class vector_t>
bool PointIsInsideTriangle(const vector_t& v, const vector_t& v0,
              const vector_t& v1, const vector_t& v2);
 
 
template <class vector_t>
bool PointIsInsideTriangle_HighAcc(const vector_t& v, const vector_t& v0,
                   const vector_t& v1, const vector_t& v2);
 
 
template <class vector_t>
bool PointIsInsideQuadrilateral(const vector_t& v, const vector_t& v0,
                const vector_t& v1, const vector_t& v2,
                const vector_t& v3);
 
//  PointIsInsideTetrahedron
 
template <class vector_t>
bool PointIsInsideTetrahedron(const vector_t& v, const vector_t& v0, const vector_t& v1,
                const vector_t& v2, const vector_t& v3);
 
//  CalculateTetrahedronVolume - mstepnie
UG_API
number CalculateTetrahedronVolume(const vector3& a, const vector3& b,
                  const vector3& c, const vector3& d);
 
//  CalculatePyramidVolume - mscherer
number CalculatePyramidVolume(const vector3& v0, const vector3& v1,
    const vector3& v2, const vector3& v3, const vector3& v4);
 
//  CalculatePrismVolume - mscherer
number CalculatePrismVolume(const vector3& v0, const vector3& v1,
    const vector3& v2, const vector3& v3, const vector3& v4,
    const vector3& v5);
 
//  CalculateHexahedronVolume - mscherer
number CalculateHexahedronVolume(const vector3& v0, const vector3& v1,
    const vector3& v2, const vector3& v3, const vector3& v4,
    const vector3& v5, const vector3& v6, const vector3& v8);
 
//  CalculateHexahedronVolume - mstepnie
number CalculateOctahedronVolume(const vector3& a, const vector3& b,
    const vector3& c, const vector3& d, const vector3& e,
    const vector3& f);
 
//  BinomialCoefficient
 
UG_API
int BinomCoeff(int n, int k);
 
//  ReflectVectorAtPlane
 
template <class vector_t>
void ReflectVectorAtPlane(vector_t& vReflectedOut, const vector_t& v,
                          const vector_t& n, const vector_t& r0);
 
// end group ugbase_math
 
}// end of namespace
 
//  include implementation
#include "math_util_impl.hpp"
 
#endif