BUGFIX: Seg-fault in multiphaseInterFoam. Author: Henrik Rusche. Merge: Hrvoje Jasak
[foam-extend-3.2.git] / src / coupledMatrix / coupledFvMatrices / coupledFvScalarMatrix / coupledFvScalarMatrix.C
blob6ea52dff8995c8045d43e08dac5a55908ecdcb51
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | foam-extend: Open Source CFD
4    \\    /   O peration     | Version:     3.2
5     \\  /    A nd           | Web:         http://www.foam-extend.org
6      \\/     M anipulation  | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
8 License
9     This file is part of foam-extend.
11     foam-extend 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.
16     foam-extend is distributed in the hope that it will be useful, but
17     WITHOUT ANY WARRANTY; without even the implied warranty of
18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19     General Public License for more details.
21     You should have received a copy of the GNU General Public License
22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
24 Description
25      Finite-Volume scalar matrix member functions and operators
27 \*---------------------------------------------------------------------------*/
29 #include "coupledFvScalarMatrix.H"
30 #include "zeroGradientFvPatchFields.H"
31 #include "lduInterfaceFieldPtrsList.H"
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 namespace Foam
38 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
40 template<>
41 coupledSolverPerformance coupledFvMatrix<scalar>::solve
43     const dictionary& solverControls
46     if (debug)
47     {
48         InfoIn("coupledFvMatrix<Type>::solve(const dictionary)")
49             << "solving coupledFvMatrix<Type>" << endl;
50     }
52     typedef FieldField<Field, scalar> scalarFieldField;
54     PtrList<lduMatrix>& matrices = *this;
56     // Make a copy of the diagonal and complete the source
57     scalarFieldField saveDiag(this->size());
58     scalarFieldField psi(this->size());
59     FieldField<Field, scalar> source(this->size());
60     lduInterfaceFieldPtrsListList interfaces(this->size());
62     PtrList<FieldField<Field, scalar> > bouCoeffs(this->size());
63     PtrList<FieldField<Field, scalar> > intCoeffs(this->size());
65     // Prepare block solution
66     forAll (matrices, rowI)
67     {
68         fvScalarMatrix& curMatrix =
69             static_cast<fvScalarMatrix&>(matrices[rowI]);
71         saveDiag.set(rowI, new scalarField(curMatrix.diag()));
72         // HR 17/Feb/2013
73         // Need to be able to compare references to support hacks such as in jumpCyclic
74         // psi.set(rowI, new scalarField(curMatrix.psi()));
75         psi.set(rowI, &curMatrix.psi());
76         source.set(rowI, new scalarField(curMatrix.source()));
78         curMatrix.addBoundarySource(source[rowI], 0);
80         interfaces[rowI] = curMatrix.psi().boundaryField().interfaces();
82         curMatrix.addBoundaryDiag(curMatrix.diag(), 0);
84         bouCoeffs.set
85         (
86             rowI,
87             new FieldField<Field, scalar>
88             (
89                 curMatrix.boundaryCoeffs()
90             )
91         );
93         intCoeffs.set
94         (
95             rowI,
96             new FieldField<Field, scalar>
97             (
98                 curMatrix.internalCoeffs().component(0)
99             )
100         );
101     }
103     // Solver call
104     coupledSolverPerformance solverPerf = coupledLduSolver::New
105     (
106         this->coupledPsiName(),
107         *this,
108         bouCoeffs,
109         intCoeffs,
110         interfaces,
111         solverControls
112     )->solve(psi, source, 0);
114     solverPerf.print();
116     // HR 17/Feb/2013
117     // Not needed since reference is used
118     // Update solution
119     //forAll (matrices, rowI)
120     //{
121     //    fvScalarMatrix& curMatrix =
122     //        static_cast<fvScalarMatrix&>(matrices[rowI]);
123     //
124     //    curMatrix.psi().internalField() = psi[rowI];
125     //}
127     // Correct boundary conditions
128     forAll (matrices, rowI)
129     {
130         fvScalarMatrix& curMatrix =
131             static_cast<fvScalarMatrix&>(matrices[rowI]);
133         curMatrix.psi().correctBoundaryConditions();
134     }
136     //HR 17.2.2013: Clear references to internal field without deleting the objects
137     forAll (matrices, rowI)
138     {
139         psi.set(rowI, NULL).ptr();
140     }
142     return solverPerf;
146 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
148 } // End namespace Foam
150 // ************************************************************************* //