This file is indexed.

/usr/include/gmsh/STensor43.h is in libgmsh-dev 2.10.1+dfsg1-1ubuntu4.

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
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
// Gmsh - Copyright (C) 1997-2015 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 _STENSOR43_H_
#define _STENSOR43_H_

#include "STensor33.h"
#include "fullMatrix.h"
#include "Numeric.h"

// concrete class for general 3rd-order tensor in three-dimensional space

class STensor43 {
 protected:
  // 0000 0001 0002 0010 ... 2211 2212 2220 2221 2222
  double _val[81];
 public:
  inline int getIndex(int i, int j, int k, int l) const
  {
    static int _index[3][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}}},
                                     {{{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}}}};
    return _index[i][j][k][l];
  }
  STensor43(const STensor43& other)
  {
    for (int i = 0; i < 81; i++) _val[i] = other._val[i];
  }
  // default constructor, null tensor
  STensor43(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++)
          for (int l = 0; l < 3; l++)
            if ((i==k)&&(j==l))
              _val[getIndex(i, j, k, l)]=v;
            else
              _val[getIndex(i, j, k, l)]=0.0;
  }
  // Symmetric identity tensor
  STensor43(const double vik, const double vil)
  {
    for (int i = 0; i < 3; i++)
      for (int j = 0; j < 3; j++)
        for (int k = 0; k < 3; k++)
          for (int l = 0; l < 3; l++)
          {
            _val[getIndex(i, j, k, l)]= 0.;
            if ((i==k)&&(j==l))
              _val[getIndex(i, j, k, l)]+=0.5*vik;
            if ((i==l)&&(j==k))
              _val[getIndex(i, j, k, l)]+=0.5*vil;
          }
  }
  inline double &operator()(int i, int j,int k, int l)
  {
    return _val[getIndex(i, j, k, l)];
  }
  inline double operator()(int i, int j, int k, int l) const
  {
    return _val[getIndex(i, j, k ,l)];
  }
  STensor43 operator + (const STensor43 &other) const
  {
    STensor43 res(*this);
    for (int i = 0; i < 81; i++) res._val[i] += other._val[i];
    return res;
  }
  STensor43& operator += (const STensor43 &other)
  {
    for (int i = 0; i < 81; i++) _val[i] += other._val[i];
    return *this;
  }
  STensor43& operator -= (const STensor43 &other)
  {
    for (int i = 0; i < 81; i++) _val[i] -= other._val[i];
    return *this;
  }
  STensor43& operator *= (const double &other)
  {
    for (int i = 0; i < 81; i++) _val[i] *= other;
    return *this;
  }
  STensor43 transpose (int n, int m) const
  {
    STensor43 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++)
            for (int l = 0; l < 3; l++)
              ithis(i,j,k,l) = (*this)(j,i,k,l);
      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++)
            for (int l = 0; l < 3; l++)
              ithis(i,j,k,l) = (*this)(k,j,i,l);
      return ithis;
    }
    if ((n==0 && m==3) || (n==3 && m==0))
    {
      for (int i = 0; i < 3; i++)
        for (int j = 0; j < 3; j++)
          for (int k = 0; k < 3; k++)
            for (int l = 0; l < 3; l++)
              ithis(i,j,k,l) = (*this)(l,j,k,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++)
            for (int l = 0; l < 3; l++)
              ithis(i,j,k,l) = (*this)(i,k,j,l);
      return ithis;
    }
    if ((n==1 && m==3) || (n==3 && m==1))
    {
      for (int i = 0; i < 3; i++)
        for (int j = 0; j < 3; j++)
          for (int k = 0; k < 3; k++)
            for (int l = 0; l < 3; l++)
              ithis(i,j,k,l) = (*this)(i,l,k,j);
      return ithis;
    }
    if ((n==2 && m==3) || (n==3 && m==2))
    {
      for (int i = 0; i < 3; i++)
        for (int j = 0; j < 3; j++)
          for (int k = 0; k < 3; k++)
            for (int l = 0; l < 3; l++)
              ithis(i,j,k,l) = (*this)(i,j,l,k);
      return ithis;
    }
    return ithis+=(*this);
  }
/*  STensor43& operator *= (const STensor43 &other)
  {
// to be implemented
    return *this;
  }*/
  void print(const char *) const;
  const double* data() const {return _val;}
  double* data() {return _val;}
};

// tensor product
inline void tensprod(const STensor3 &a, const STensor3 &b, STensor43 &c)
{
    for (int i = 0; i < 3; i++)
      for (int j = 0; j < 3; j++)
        for (int k = 0; k < 3; k++)
          for (int l = 0; l < 3; l++)
            c(i,j,k,l)=a(i,j)*b(k,l);
}

inline void tensprod(const SVector3 &a, const STensor33 &b, STensor43 &c)
{
    for (int i = 0; i < 3; i++)
      for (int j = 0; j < 3; j++)
        for (int k = 0; k < 3; k++)
          for (int l = 0; l < 3; l++)
            c(i,j,k,l)=a(i)*b(j,k,l);
}
inline void tensprod(const STensor33 &a, const SVector3 &b, STensor43 &c)
{
    for (int i = 0; i < 3; i++)
      for (int j = 0; j < 3; j++)
        for (int k = 0; k < 3; k++)
          for (int l = 0; l < 3; l++)
            c(i,j,k,l)=a(i,j,k)*b(l);
}

inline double dot(const STensor43 &a, const STensor43 &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++)
        for (int l = 0; l < 3; l++)
      prod+=a(i,j,k,l)*b(i,j,k,l);
  return prod;
}

// full contracted product
inline STensor43 operator*(const STensor43 &t, double m)
{
  STensor43 val(t);
  val *= m;
  return val;
}
inline STensor43 operator*(double m,const STensor43 &t)
{
  STensor43 val(t);
  val *= m;
  return val;
}

inline STensor33 operator*(const STensor43 &t, const SVector3 &m)
{
  STensor33 val(0.);
  for (int i = 0; i < 3; i++)
    for (int j = 0; j < 3; j++)
      for (int k = 0; k < 3; k++)
        for (int l = 0; l < 3; l++)
          val(i,j,k)+=t(i,j,k,l)*m(l);
  return val;
}
inline STensor33 operator*( const SVector3 &m , const STensor43 &t)
{
  STensor33 val(0.);
  for (int i = 0; i < 3; i++)
    for (int j = 0; j < 3; j++)
      for (int k = 0; k < 3; k++)
        for (int l = 0; l < 3; l++)
          val(j,k,l)+=m(i)*t(i,j,k,l);
  return val;
}

inline STensor3 operator*(const STensor43 &t, const STensor3 &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++)
        for (int l = 0; l < 3; l++)
          val(i,j)+=t(i,j,k,l)*m(l,k);
  return val;
}
inline STensor3 operator*( const STensor3 &m , const STensor43 &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++)
        for (int l = 0; l < 3; l++)
          val(k,l)+=m(j,i)*t(i,j,k,l);
  return val;
}

inline SVector3 operator*(const STensor43 &t, const STensor33 &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++)
        for (int l = 0; l < 3; l++)
          val(i)+=t(i,j,k,l)*m(l,k,j);
  return val;
}
inline SVector3 operator*( const STensor33 &m , const STensor43 &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++)
        for (int l = 0; l < 3; l++)
          val(l)+=m(k,j,i)*t(i,j,k,l);
  return val;
}

inline double operator*( const STensor43 &m , const STensor43 &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++)
        for (int l = 0; l < 3; l++)
          val+=m(i,j,k,l)*t(l,k,j,i);
  return val;
}

#endif