CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM (http://www.cfd-online.com/Forums/openfoam/)
-   -   How the boundary conditions are disposed into the matrix coefficient? (http://www.cfd-online.com/Forums/openfoam/113490-how-boundary-conditions-disposed-into-matrix-coefficient.html)

mxylondon February 20, 2013 12:39

How the boundary conditions are disposed into the matrix coefficient?
 
Hi guys, does any one know?
How the boundary conditions are disposed into the matrix A of Ax=B?
In other word, which section of source code I need to find about it?

Cheers

santiagomarquezd March 1, 2013 17:03

Hi the BC contributions (internalCoffs and boundaryCoeffs) to the matrix and RHS are calculated where the discrete version of the differential operators is assembled. These contributions are stored to be inserted into the matrix prior to the solution. In the case of a case of scalar transport and the div operator the lines are 98 in gaussConvectionScheme.C

Code:

00066 template<class Type>
00067 tmp<fvMatrix<Type> >
00068 gaussConvectionScheme<Type>::fvmDiv
00069 (
00070    const surfaceScalarField& faceFlux,
00071    const GeometricField<Type, fvPatchField, volMesh>& vf
00072 ) const
00073 {
00074    tmp<surfaceScalarField> tweights = tinterpScheme_().weights(vf);
00075    const surfaceScalarField& weights = tweights();
00076
00077    tmp<fvMatrix<Type> > tfvm
00078    (
00079        new fvMatrix<Type>
00080        (
00081            vf,
00082            faceFlux.dimensions()*vf.dimensions()
00083        )
00084    );
00085    fvMatrix<Type>& fvm = tfvm();
00086
00087    fvm.lower() = -weights.internalField()*faceFlux.internalField();
00088    fvm.upper() = fvm.lower() + faceFlux.internalField();
00089    fvm.negSumDiag();
00090
00091    forAll(vf.boundaryField(), patchI)
00092    {
00093        const fvPatchField<Type>& psf = vf.boundaryField()[patchI];
00094        const fvsPatchScalarField& patchFlux = faceFlux.boundaryField()[patchI];
00095        const fvsPatchScalarField& pw = weights.boundaryField()[patchI];
00096
00097        fvm.internalCoeffs()[patchI] = patchFlux*psf.valueInternalCoeffs(pw);
00098        fvm.boundaryCoeffs()[patchI] = -patchFlux*psf.valueBoundaryCoeffs(pw);
00099    }
00100
00101    if (tinterpScheme_().corrected())
00102    {
00103        fvm += fvc::surfaceIntegrate(faceFlux*tinterpScheme_().correction(vf));
00104    }
00105
00106    return tfvm;
00107 }

and 106 and 109 in fvScalarMatrix.C

Code:

00095 template<>
00096 Foam::lduMatrix::solverPerformance Foam::fvMatrix<Foam::scalar>::fvSolver::solve
00097 (
00098    const dictionary& solverControls
00099 )
00100 {
00101    GeometricField<scalar, fvPatchField, volMesh>& psi =
00102        const_cast<GeometricField<scalar, fvPatchField, volMesh>&>
00103        (fvMat_.psi());
00104
00105    scalarField saveDiag(fvMat_.diag());
00106    fvMat_.addBoundaryDiag(fvMat_.diag(), 0);
00107
00108    scalarField totalSource(fvMat_.source());
00109    fvMat_.addBoundarySource(totalSource, false);
00110
00111    // assign new solver controls
00112    solver_->read(solverControls);
00113
00114    lduMatrix::solverPerformance solverPerf = solver_->solve
00115    (
00116        psi.internalField(),
00117        totalSource
00118    );
00119
00120    solverPerf.print();
00121
00122    fvMat_.diag() = saveDiag;
00123
00124    psi.correctBoundaryConditions();
00125
00126    psi.mesh().setSolverPerformance(psi.name(), solverPerf);
00127
00128    return solverPerf;
00129 }

Regards.


All times are GMT -4. The time now is 01:00.