src/mesh/cfMesh/utilities/helperClasses/matrices/matrix3DI.H

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | cfMesh: A library for mesh generation
   4    \\    /   O peration     |
   5     \\  /    A nd           | Author: Franjo Juretic (franjo.juretic@c-fields.com)
   6      \\/     M anipulation  | Copyright (C) Creative Fields, Ltd.
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of cfMesh.
  10
  11     cfMesh is free software; you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation; either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     cfMesh is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with cfMesh.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Class
  25     matrix3D
  26
  27 Description
  28     Implementation od 3 x 3 matrix
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32
  33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  34
  35 namespace Foam
  36 {
  37
  38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  39
  40 inline void matrix3D::calculateDeterminant()
  41 {
  42     if( calculatedDet_ )
  43         return;
  44
  45     const FixedList<FixedList<scalar, 3>, 3>& mat = *this;
  46
  47     det_ =
  48         mat[0][0] * (mat[1][1] * mat[2][2] - mat[1][2] * mat[2][1]) -
  49         mat[0][1] * (mat[1][0] * mat[2][2] - mat[1][2] * mat[2][0]) +
  50         mat[0][2] * (mat[1][0] * mat[2][1] - mat[1][1] * mat[2][0]);
  51
  52     calculatedDet_ = true;
  53 }
  54
  55 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  56
  57 inline matrix3D::matrix3D()
  58 :
  59     det_(),
  60     calculatedDet_(false)
  61 {}
  62
  63 inline matrix3D::matrix3D(const matrix3D& m)
  64 :
  65     FixedList<FixedList<scalar, 3>, 3>(m),
  66     det_(m.det_),
  67     calculatedDet_(m.calculatedDet_)
  68 {}
  69
  70 inline matrix3D::~matrix3D()
  71 {}
  72
  73 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  74
  75 inline scalar matrix3D::determinant()
  76 {
  77     calculateDeterminant();
  78
  79     return det_;
  80 }
  81
  82 inline matrix3D matrix3D::inverse()
  83 {
  84     calculateDeterminant();
  85
  86     const FixedList<FixedList<scalar, 3>, 3>& mat = *this;
  87
  88     matrix3D imat;
  89     imat[0][0] = (mat[1][1] * mat[2][2] - mat[1][2] * mat[2][1]) / det_;
  90     imat[0][1] = (mat[0][2] * mat[2][1] - mat[0][1] * mat[2][2]) / det_;
  91     imat[0][2] = (mat[0][1] * mat[1][2] - mat[0][2] * mat[1][1]) / det_;
  92     imat[1][0] = (mat[1][2] * mat[2][0] - mat[1][0] * mat[2][2]) / det_;
  93     imat[1][1] = (mat[0][0] * mat[2][2] - mat[0][2] * mat[2][0]) / det_;
  94     imat[1][2] = (mat[0][2] * mat[1][0] - mat[0][0] * mat[1][2]) / det_;
  95     imat[2][0] = (mat[1][0] * mat[2][1] - mat[1][1] * mat[2][0]) / det_;
  96     imat[2][1] = (mat[0][1] * mat[2][0] - mat[0][0] * mat[2][1]) / det_;
  97     imat[2][2] = (mat[0][0] * mat[1][1] - mat[0][1] * mat[1][0]) / det_;
  98
  99     return imat;
 100 }
 101
 102 inline FixedList<scalar, 3> matrix3D::solve
 103 (
 104     const FixedList<scalar, 3>& source
 105 )
 106 {
 107     FixedList<scalar, 3> result;
 108     result[0] = solveFirst(source);
 109     result[1] = solveSecond(source);
 110     result[2] = solveThird(source);
 111
 112     return result;
 113 }
 114
 115 inline scalar matrix3D::solveFirst(const FixedList<scalar, 3>& source)
 116 {
 117     calculateDeterminant();
 118
 119     const FixedList<FixedList<scalar, 3>, 3>& mat = *this;
 120
 121     return
 122     (
 123         (mat[1][1] * mat[2][2] - mat[1][2] * mat[2][1]) * source[0] +
 124         (mat[0][2] * mat[2][1] - mat[0][1] * mat[2][2]) * source[1] +
 125         (mat[0][1] * mat[1][2] - mat[0][2] * mat[1][1]) * source[2]
 126     ) / det_;
 127 }
 128
 129 inline scalar matrix3D::solveSecond(const FixedList<scalar, 3>& source)
 130 {
 131     calculateDeterminant();
 132
 133     const FixedList<FixedList<scalar, 3>, 3>& mat = *this;
 134
 135     return
 136     (
 137         (mat[1][2] * mat[2][0] - mat[1][0] * mat[2][2]) * source[0] +
 138         (mat[0][0] * mat[2][2] - mat[0][2] * mat[2][0]) * source[1] +
 139         (mat[0][2] * mat[1][0] - mat[0][0] * mat[1][2]) * source[2]
 140     ) / det_;
 141 }
 142
 143 inline scalar matrix3D::solveThird(const FixedList<scalar, 3>& source)
 144 {
 145     calculateDeterminant();
 146
 147     const FixedList<FixedList<scalar, 3>, 3>& mat = *this;
 148
 149     return
 150     (
 151         (mat[1][0] * mat[2][1] - mat[1][1] * mat[2][0]) * source[0] +
 152         (mat[0][1] * mat[2][0] - mat[0][0] * mat[2][1]) * source[1] +
 153         (mat[0][0] * mat[1][1] - mat[0][1] * mat[1][0]) * source[2]
 154     ) / det_;
 155 }
 156
 157 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 158
 159 } // End namespace Foam
 160
 161 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //