/usr/include/gmsh/STensor33.h is in libgmsh-dev 2.8.5+dfsg-1.1+b1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 | // Gmsh - Copyright (C) 1997-2014 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
//
// Contributor(s):
// Eric Bechet
//
#ifndef _STENSOR33_H_
#define _STENSOR33_H_
#include "STensor3.h"
#include "fullMatrix.h"
#include "Numeric.h"
// concrete class for general 3rd-order tensor in three-dimensional space
class STensor33 {
protected:
// 000 001 002 010 ... 211 212 220 221 222
double _val[27];
public:
inline int getIndex(int i, int j, int k) const
{
static int _index[3][3][3] = {{{0,1,2},{3,4,5},{6,7,8}},{{9,10,11},{12,13,14},{15,16,17}},{{18,19,20},{21,22,23},{24,25,26}}};
return _index[i][j][k];
}
STensor33(const STensor33& other)
{
for (int i = 0; i < 27; i++) _val[i] = other._val[i];
}
// default constructor, null tensor
STensor33(const double v = 0.0)
{
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
_val[getIndex(i, j, k)]=v;
}
inline double &operator()(int i, int j,int k)
{
return _val[getIndex(i, j, k)];
}
inline double operator()(int i, int j, int k) const
{
return _val[getIndex(i, j, k)];
}
STensor33 operator + (const STensor33 &other) const
{
STensor33 res(*this);
for (int i = 0; i < 27; i++) res._val[i] += other._val[i];
return res;
}
STensor33& operator += (const STensor33 &other)
{
for (int i = 0; i < 27; i++) _val[i] += other._val[i];
return *this;
}
STensor33& operator -= (const STensor33 &other)
{
for (int i = 0; i < 27; i++) _val[i] -= other._val[i];
return *this;
}
STensor33& operator *= (const double &other)
{
for (int i = 0; i < 27; i++) _val[i] *= other;
return *this;
}
STensor33 transpose (int n, int m) const
{
STensor33 ithis;
if ((n==0 && m==1) || (n==1 && m==0))
{
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
ithis(i,j,k) = (*this)(j,i,k);
return ithis;
}
if ((n==0 && m==2) || (n==2 && m==0))
{
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
ithis(i,j,k) = (*this)(k,j,i);
return ithis;
}
if ((n==1 && m==2) || (n==2 && m==1))
{
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
ithis(i,j,k) = (*this)(i,k,j);
return ithis;
}
return ithis+=(*this);
}
/* STensor33& operator *= (const STensor33 &other)
{
// to be implemented
return *this;
}*/
void print(const char *) const;
};
// tensor product
inline void tensprod(const SVector3 &a, const STensor3 &b, STensor33 &c)
{
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
c(i,j,k)=a(i)*b(j,k);
}
inline void tensprod(const STensor3 &a, const SVector3 &b, STensor33 &c)
{
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
c(i,j,k)=a(i,j)*b(k);
}
inline double dot(const STensor33 &a, const STensor33 &b)
{
double prod=0;
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
prod+=a(i,j,k)*b(i,j,k);
return prod;
}
// full contracted product
inline STensor33 operator*(const STensor33 &t, double m)
{
STensor33 val(t);
val *= m;
return val;
}
inline STensor33 operator*(double m,const STensor33 &t)
{
STensor33 val(t);
val *= m;
return val;
}
inline STensor3 operator*(const STensor33 &t, const SVector3 &m)
{
STensor3 val(0.);
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
val(i,j)+=t(i,j,k)*m(k);
return val;
}
inline STensor3 operator*( const SVector3 &m , const STensor33 &t)
{
STensor3 val(0.);
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
val(j,k)+=m(i)*t(i,j,k);
return val;
}
inline SVector3 operator*(const STensor33 &t, const STensor3 &m)
{
SVector3 val(0.);
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
val(i)+=t(i,j,k)*m(k,j);
return val;
}
inline SVector3 operator*( const STensor3 &m , const STensor33 &t)
{
SVector3 val(0.);
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
val(k)+=m(j,i)*t(i,j,k);
return val;
}
inline double operator*(const STensor33 &m, const STensor33 &t)
{
double val(0.);
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
for (int k = 0; k < 3; k++)
val+=m(i,j,k)*t(k,j,i);
return val;
}
#endif
|