mirror of
https://github.com/NGSolve/netgen.git
synced 2024-09-22 17:52:50 +05:00
405 lines
7.2 KiB
C++
405 lines
7.2 KiB
C++
#ifndef FILE_GEOMOPS
|
|
#define FILE_GEOMOPS
|
|
|
|
/* *************************************************************************/
|
|
/* File: geomops.hpp */
|
|
/* Author: Joachim Schoeberl */
|
|
/* Date: 20. Jul. 02 */
|
|
/* *************************************************************************/
|
|
|
|
|
|
namespace netgen
|
|
{
|
|
|
|
/*
|
|
|
|
Point - Vector operations
|
|
|
|
*/
|
|
|
|
|
|
template <int D, typename T>
|
|
inline Vec<D,T> operator+ (Vec<D,T> a, Vec<D,T> b)
|
|
{
|
|
Vec<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = a(i) + b(i);
|
|
return res;
|
|
}
|
|
|
|
|
|
|
|
template <int D, typename T>
|
|
inline Point<D,T> operator+ (Point<D,T> a, Vec<D,T> b)
|
|
{
|
|
Point<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = a(i) + b(i);
|
|
return res;
|
|
}
|
|
|
|
|
|
|
|
template <int D, typename T>
|
|
inline Vec<D,T> operator- (Point<D,T> a, Point<D,T> b)
|
|
{
|
|
Vec<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = a(i) - b(i);
|
|
return res;
|
|
}
|
|
|
|
template <int D, typename T>
|
|
inline Point<D,T> operator- (Point<D,T> a, Vec<D,T> b)
|
|
{
|
|
Point<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = a(i) - b(i);
|
|
return res;
|
|
}
|
|
|
|
template <int D, typename T>
|
|
inline Vec<D,T> operator- (Vec<D,T> a, Vec<D,T> b)
|
|
{
|
|
Vec<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = a(i) - b(i);
|
|
return res;
|
|
}
|
|
|
|
|
|
|
|
template <int D, typename T>
|
|
inline Vec<D,T> operator* (T s, Vec<D,T> b)
|
|
{
|
|
Vec<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = s * b(i);
|
|
return res;
|
|
}
|
|
|
|
template <int D, int SW>
|
|
inline Vec<D, SIMD<double,SW>> operator* (SIMD<double,SW> s, Vec<D,double> b)
|
|
{
|
|
Vec<D,SIMD<double,SW>> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = s * b(i);
|
|
return res;
|
|
}
|
|
|
|
|
|
template <int D>
|
|
inline double operator* (Vec<D> a, Vec<D> b)
|
|
{
|
|
double sum = 0;
|
|
for (int i = 0; i < D; i++)
|
|
sum += a(i) * b(i);
|
|
return sum;
|
|
}
|
|
|
|
|
|
|
|
template <int D, typename T>
|
|
inline Vec<D,T> operator- (Vec<D,T> b)
|
|
{
|
|
Vec<D,T> res;
|
|
for (int i = 0; i < D; i++)
|
|
res(i) = -b(i);
|
|
return res;
|
|
}
|
|
|
|
|
|
template <int D, typename T>
|
|
inline Point<D,T> & operator+= (Point<D,T> & a, Vec<D,T> b)
|
|
{
|
|
for (int i = 0; i < D; i++)
|
|
a(i) += b(i);
|
|
return a;
|
|
}
|
|
|
|
template <int D, typename T>
|
|
inline Vec<D,T> & operator+= (Vec<D,T> & a, Vec<D> b)
|
|
{
|
|
for (int i = 0; i < D; i++)
|
|
a(i) += b(i);
|
|
return a;
|
|
}
|
|
|
|
|
|
template <int D>
|
|
inline Point<D> & operator-= (Point<D> & a, const Vec<D> & b)
|
|
{
|
|
for (int i = 0; i < D; i++)
|
|
a(i) -= b(i);
|
|
return a;
|
|
}
|
|
|
|
template <int D>
|
|
inline Vec<D> & operator-= (Vec<D> & a, const Vec<D> & b)
|
|
{
|
|
for (int i = 0; i < D; i++)
|
|
a(i) -= b(i);
|
|
return a;
|
|
}
|
|
|
|
|
|
|
|
template <int D, typename T1, typename T2>
|
|
inline Vec<D,T1> & operator*= (Vec<D,T1> & a, T2 s)
|
|
{
|
|
for (int i = 0; i < D; i++)
|
|
a(i) *= s;
|
|
return a;
|
|
}
|
|
|
|
|
|
template <int D>
|
|
inline Vec<D> & operator/= (Vec<D> & a, double s)
|
|
{
|
|
for (int i = 0; i < D; i++)
|
|
a(i) /= s;
|
|
return a;
|
|
}
|
|
|
|
|
|
|
|
|
|
// Matrix - Vector operations
|
|
|
|
/*
|
|
template <int H, int W>
|
|
inline Vec<H> operator* (const Mat<H,W> & m, const Vec<W> & v)
|
|
{
|
|
Vec<H> res;
|
|
for (int i = 0; i < H; i++)
|
|
{
|
|
res(i) = 0;
|
|
for (int j = 0; j < W; j++)
|
|
res(i) += m(i,j) * v(j);
|
|
}
|
|
return res;
|
|
}
|
|
*/
|
|
|
|
// thanks to VC60 partial template specialization features !!!
|
|
|
|
inline Vec<2> operator* (const Mat<2,2> & m, const Vec<2> & v)
|
|
{
|
|
Vec<2> res;
|
|
for (int i = 0; i < 2; i++)
|
|
{
|
|
res(i) = 0;
|
|
for (int j = 0; j < 2; j++)
|
|
res(i) += m(i,j) * v(j);
|
|
}
|
|
return res;
|
|
}
|
|
|
|
inline Vec<2> operator* (const Mat<2,3> & m, const Vec<3> & v)
|
|
{
|
|
Vec<2> res;
|
|
for (int i = 0; i < 2; i++)
|
|
{
|
|
res(i) = 0;
|
|
for (int j = 0; j < 3; j++)
|
|
res(i) += m(i,j) * v(j);
|
|
}
|
|
return res;
|
|
}
|
|
|
|
|
|
inline Vec<3> operator* (const Mat<3,2> & m, const Vec<2> & v)
|
|
{
|
|
Vec<3> res;
|
|
for (int i = 0; i < 3; i++)
|
|
{
|
|
res(i) = 0;
|
|
for (int j = 0; j < 2; j++)
|
|
res(i) += m(i,j) * v(j);
|
|
}
|
|
return res;
|
|
}
|
|
|
|
|
|
inline Vec<3> operator* (const Mat<3,3> & m, const Vec<3> & v)
|
|
{
|
|
Vec<3> res;
|
|
for (int i = 0; i < 3; i++)
|
|
{
|
|
res(i) = 0;
|
|
for (int j = 0; j < 3; j++)
|
|
res(i) += m(i,j) * v(j);
|
|
}
|
|
return res;
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
template <int H1, int W1, int H2, int W2>
|
|
inline Mat<H1,W2> operator* (const Mat<H1,W1> & a, const Mat<H2,W2> & b)
|
|
{
|
|
Mat<H1,W2> m;
|
|
for (int i = 0; i < H1; i++)
|
|
for (int j = 0; j < W2; j++)
|
|
{
|
|
double sum = 0;
|
|
for (int k = 0; k < W1; k++)
|
|
sum += a(i,k) * b(k, j);
|
|
m(i,j) = sum;
|
|
}
|
|
return m;
|
|
}
|
|
*/
|
|
|
|
template <typename T>
|
|
inline Mat<2,2,T> operator* (const Mat<2,2,T> & a, const Mat<2,2,T> & b)
|
|
{
|
|
Mat<2,2,T> m;
|
|
for (int i = 0; i < 2; i++)
|
|
for (int j = 0; j < 2; j++)
|
|
{
|
|
T sum(0);
|
|
for (int k = 0; k < 2; k++)
|
|
sum += a(i,k) * b(k, j);
|
|
m(i,j) = sum;
|
|
}
|
|
return m;
|
|
}
|
|
|
|
inline Mat<2,2> operator* (const Mat<2,3> & a, const Mat<3,2> & b)
|
|
{
|
|
Mat<2,2> m;
|
|
for (int i = 0; i < 2; i++)
|
|
for (int j = 0; j < 2; j++)
|
|
{
|
|
double sum = 0;
|
|
for (int k = 0; k < 3; k++)
|
|
sum += a(i,k) * b(k, j);
|
|
m(i,j) = sum;
|
|
}
|
|
return m;
|
|
}
|
|
|
|
template <typename T>
|
|
inline Mat<3,2,T> operator* (const Mat<3,2,T> & a, const Mat<2,2,T> & b)
|
|
{
|
|
Mat<3,2,T> m;
|
|
for (int i = 0; i < 3; i++)
|
|
for (int j = 0; j < 2; j++)
|
|
{
|
|
T sum(0.0);
|
|
for (int k = 0; k < 2; k++)
|
|
sum += a(i,k) * b(k, j);
|
|
m(i,j) = sum;
|
|
}
|
|
return m;
|
|
}
|
|
|
|
|
|
|
|
inline Mat<2,3> operator* (const Mat<2,2> & a, const Mat<2,3> & b)
|
|
{
|
|
Mat<2,3> m;
|
|
for (int i = 0; i < 2; i++)
|
|
for (int j = 0; j < 3; j++)
|
|
{
|
|
double sum = 0;
|
|
for (int k = 0; k < 2; k++)
|
|
sum += a(i,k) * b(k, j);
|
|
m(i,j) = sum;
|
|
}
|
|
return m;
|
|
}
|
|
|
|
template <typename T>
|
|
inline Mat<3,3,T> operator* (const Mat<3,3,T> & a, const Mat<3,3,T> & b)
|
|
{
|
|
Mat<3,3,T> m;
|
|
for (int i = 0; i < 3; i++)
|
|
for (int j = 0; j < 3; j++)
|
|
{
|
|
T sum = T(0);
|
|
for (int k = 0; k < 3; k++)
|
|
sum += a(i,k) * b(k, j);
|
|
m(i,j) = sum;
|
|
}
|
|
return m;
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
template <int H, int W>
|
|
inline Mat<W,H> Trans (const Mat<H,W> & m)
|
|
{
|
|
Mat<W,H> res;
|
|
for (int i = 0; i < H; i++)
|
|
for (int j = 0; j < W; j++)
|
|
res(j,i) = m(i,j);
|
|
return res;
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
template <int D, typename T>
|
|
inline ostream & operator<< (ostream & ost, const Vec<D,T> & a)
|
|
{
|
|
ost << "(";
|
|
for (int i = 0; i < D-1; i++)
|
|
ost << a(i) << ", ";
|
|
ost << a(D-1) << ")";
|
|
return ost;
|
|
}
|
|
|
|
template <int D, typename T>
|
|
inline ostream & operator<< (ostream & ost, const Point<D,T> & a)
|
|
{
|
|
ost << "(";
|
|
for (int i = 0; i < D-1; i++)
|
|
ost << a(i) << ", ";
|
|
ost << a(D-1) << ")";
|
|
return ost;
|
|
}
|
|
|
|
template <int D>
|
|
inline ostream & operator<< (ostream & ost, const Box<D> & b)
|
|
{
|
|
ost << b.PMin() << " - " << b.PMax();
|
|
return ost;
|
|
}
|
|
|
|
template <int H, int W, typename T>
|
|
inline ostream & operator<< (ostream & ost, const Mat<H,W,T> & m)
|
|
{
|
|
ost << "(";
|
|
for (int i = 0; i < H; i++)
|
|
{
|
|
for (int j = 0; j < W; j++)
|
|
ost << m(i,j) << " ";
|
|
ost << endl;
|
|
}
|
|
return ost;
|
|
}
|
|
|
|
|
|
}
|
|
|
|
#endif
|