[ala]

Hi *,

I have read the following wiki article to lduMatrix and lduAddressing:
http://openfoamwiki.net/index.php/Matrices_in_OpenFOAM
This article is a good starting point for understanding the lduMatrix and the corresponding lduAddressing. Unfortunately this article does not cover the fvMatrix which has lduMatrix as a base class.
I also searched at this forum but did not find an answer.

And here are my questions concerning the fvMatrix:
- What is the meaning and the usage of the members internalCoeffs_ and boundaryCoeffs_?
- What does "pseudo-matrix coeff" mean and for what is it used?

fvMatrix.H, lines 129-135:

//- Boundary scalar field containing pseudo-matrix coeffs
//  for internal cells
FieldField<Field, Type> internalCoeffs_;

//- Boundary scalar field containing pseudo-matrix coeffs
//  for boundary cells
FieldField<Field, Type> boundaryCoeffs_;

Thanks in advance,
Alicja

[santiagomarquezd]

Hi, ala, internalCoeffs_ and boundaryCoeffs_ are the contribution of BC to matrix coefficients, they are calculated in diff operators and applied at solving time, for example, for divergence (gaussConvectionScheme.C):

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

Regards.

[ala]

Hi Santiago,

thanks for your reply. Ok, so that are matrix contributions from BC.
I took a look at fvMatrix.C, at the functions
- addToInternalField()
- addBoundaryDiag()
- addBoundarySource():

00061 template<class Type> 
00062 template<class Type2> 
00063 void Foam::fvMatrix<Type>::addToInternalField 
00064 ( 
00065     const unallocLabelList& addr, 
00066     const tmp<Field<Type2> >& tpf, 
00067     Field<Type2>& intf 
00068 ) const 
00069 { 
00070     addToInternalField(addr, tpf(), intf); 
00071     tpf.clear(); 
00072 } 
00073  

00114  
00115 template<class Type> 
00116 void Foam::fvMatrix<Type>::addBoundaryDiag 
00117 ( 
00118     scalarField& diag, 
00119     const direction solveCmpt 
00120 ) const 
00121 { 
00122     forAll(internalCoeffs_, patchI) 
00123     { 
00124         addToInternalField 
00125         ( 
00126             lduAddr().patchAddr(patchI), 
00127             internalCoeffs_[patchI].component(solveCmpt), 
00128             diag 
00129         ); 
00130     } 
00131 } 
00132  


00148  
00149 template<class Type> 
00150 void Foam::fvMatrix<Type>::addBoundarySource 
00151 ( 
00152     Field<Type>& source, 
00153     const bool couples 
00154 ) const 
00155 { 
00156     forAll(psi_.boundaryField(), patchI) 
00157     { 
00158         const fvPatchField<Type>& ptf = psi_.boundaryField()[patchI]; 
00159         const Field<Type>& pbc = boundaryCoeffs_[patchI]; 
00160  
00161         if (!ptf.coupled()) 
00162         { 
00163             addToInternalField(lduAddr().patchAddr(patchI), pbc, source); 
00164         } 
00165         else if (couples) 
00166         { 
00167             tmp<Field<Type> > tpnf = ptf.patchNeighbourField(); 
00168             const Field<Type>& pnf = tpnf(); 
00169  
00170             const unallocLabelList& addr = lduAddr().patchAddr(patchI); 
00171  
00172             forAll(addr, facei) 
00173             { 
00174                 source[addr[facei]] += cmptMultiply(pbc[facei], pnf[facei]); 
00175             } 
00176         } 
00177     } 
00178 } 
00179

I understand why I have to distinguish between the contribution of BC to matrix coefficients as boundaryCoeffs and internalCoeffs. Thank you for your hint : )

Now i don't know, why are this coefficients separately stored (as boundaryCoeffs and as internalCoeffs) and not immediately added to the 'real' matrix, and i have to (?) add this coefficients firstly during the solve( ) process.

Regards, Ala

[randolph]

Originally Posted by santiagomarquezd

Hi, ala, internalCoeffs_ and boundaryCoeffs_ are the contribution of BC to matrix coefficients, they are calculated in diff operators and applied at solving time, for example, for divergence (gaussConvectionScheme.C):

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

Regards.

Hi,

Is it possible to apply different boundary condition for the different operator? Any suggestions?

Thanks,
Rdf

[santiagomarquezd]

Hmm, the BC's are intended to be selected just by patch...