Go to the documentation of this file.
44 for (
label i = 0; i < m; ++i)
47 scalar largestCoeff =
mag(tmpMatrix(iMax, i));
50 for (
label j = i + 1; j < m; ++j)
52 if (
mag(tmpMatrix(j, i)) > largestCoeff)
55 largestCoeff =
mag(tmpMatrix(iMax, i));
63 Swap(tmpMatrix(i,
k), tmpMatrix(iMax,
k));
65 Swap(sourceSol[i], sourceSol[iMax]);
69 if (
mag(tmpMatrix(i, i)) < 1
e-20)
77 for (
label j = i + 1; j < m; ++j)
79 sourceSol[j] -= sourceSol[i]*(tmpMatrix(j, i)/tmpMatrix(i, i));
84 tmpMatrix(i,
k)*tmpMatrix(j, i)/tmpMatrix(i, i);
90 for (
label j = m - 1; j >= 0; --j)
96 ntempvec += tmpMatrix(j,
k)*sourceSol[
k];
99 sourceSol[j] = (sourceSol[j] - ntempvec)/tmpMatrix(j, j);
130 for (
label i = 0; i < m; ++i)
132 label ip = pivotIndices[i];
133 Type
sum = sourceSol[ip];
134 sourceSol[ip] = sourceSol[i];
135 const scalar* __restrict__ luMatrixi = luMatrix[i];
139 for (
label j = ii - 1; j < i; ++j)
141 sum -= luMatrixi[j]*sourceSol[j];
152 for (
label i = m - 1; i >= 0; --i)
154 Type
sum = sourceSol[i];
155 const scalar* __restrict__ luMatrixi = luMatrix[i];
157 for (
label j = i + 1; j < m; ++j)
159 sum -= luMatrixi[j]*sourceSol[j];
162 sourceSol[i] =
sum/luMatrixi[i];
178 for (
label i = 0; i < m; ++i)
180 Type
sum = sourceSol[i];
181 const scalar* __restrict__ luMatrixi = luMatrix[i];
185 for (
label j = ii - 1; j < i; ++j)
187 sum -= luMatrixi[j]*sourceSol[j];
195 sourceSol[i] =
sum/luMatrixi[i];
198 for (
label i = m - 1; i >= 0; --i)
200 Type
sum = sourceSol[i];
201 const scalar* __restrict__ luMatrixi = luMatrix[i];
203 for (
label j = i + 1; j < m; ++j)
205 sum -= luMatrixi[j]*sourceSol[j];
208 sourceSol[i] =
sum/luMatrixi[i];
238 template<
class Form,
class Type>
249 <<
"A and B must have identical inner dimensions but A.n = "
250 <<
A.n() <<
" and B.m = " <<
B.m()
256 for (
label i = 0; i <
A.m(); ++i)
258 for (
label j = 0; j <
B.n(); ++j)
260 for (
label l = 0; l <
B.m(); ++l)
262 ans(i, j) +=
A(i, l)*
B(l, j);
void LUsolve(scalarSquareMatrix &matrix, List< Type > &source)
static constexpr const zero Zero
Global zero.
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
void Swap(DynamicList< T, SizeMin1 > &a, DynamicList< T, SizeMin2 > &b)
A templated (m x n) matrix of objects of <T>. The layout is (mRows x nCols) - row-major order:
Swap arguments as per std::swap, but in Foam namespace.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
errorManip< error > abort(error &err)
A templated 2D square symmetric matrix of objects of <T>, where the n x n matrix dimension is known a...
void LUBacksubstitute(const scalarSquareMatrix &luMmatrix, const labelList &pivotIndices, List< Type > &source)
errorManipArg< error, int > exit(error &err, const int errNo=1)
void LUDecompose(scalarSquareMatrix &matrix, labelList &pivotIndices)
LU decompose the matrix with pivoting.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label m() const noexcept
The number of rows.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Various functions to operate on Lists.
void multiply(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
const volScalarField & psi