
[Sponsors] 
November 22, 2012, 10:58 
Some questions about Laplacian with nonorthogonal correction

#1 
New Member
Join Date: Jan 2012
Posts: 6
Rep Power: 7 
Hi,
I have taken a look at the nonorthogonal correction for the Laplacian for a scalar variable in Openfoam and compared the treatment of the corrections to the description of the diffusive term in Hrv's thesis. I will first write what i think is happening in the OF code and then have some questions, because not everthing is clear to me. The code lines in the following refer to the current version of the source guide. As far as i figured out the interesting things happen in gaussLaplacianSchemes.C after line 39. In line 56 with fvmLaplacianUncorrected(gammaMagSf, vf) the implicit part of the Laplacian is written into the matrix. I think this is the function defined in gaussLaplacianScheme.C in line 46. It writes into the ofdiagonal matrix coeficient gamma_f*s/d with nomenclature according to Hrv's thesis. Then the diagonal matrix entries are filled with the negative sum of the ofdiagonal enties (why the negative sum?). The boundaries are also taken care of. Then we jump back to gaussLaplacianSchemes.C. In line 66 the content of fvm.faceFluxCorrectionPtr() is filled with: gamma_f*s*correction correction is defined in correctedSnGrads.C in line 42 and just returns the values of fullGradCorrection(vf), which is defined in correctedSnGrad.C in line 44. fullGradCorrection(vf) a bit confused me, but i think it just returns the inner product of the interpolated gradient at the face with the nonOrthCorrectionVector: nonOrthCorrectionVector dot nabla(vf)_f The nonOrthCorrectionVectors and the geometric things needed for them are built in surfaceinterpolation.C as follows: deltaCoeffs = 1/d (line 232) nonOrthDeltaCoeffs= 1/(n dot d), (stabilised with the max...) (line 292) nonOrthCorrectionVectors= n  d/(n dot d) (line 340), which is just nonOrthCorrectionVectors= k/s= (s(ss/sd)d) /s with k according to Hrv's thesis k=s(ss/sd)d So gammaMagSf*this>tsnGradScheme_().correction(vf) in gaussLaplacianSchemes.C line 66 fills the faceFluxCorrectionPtr with: gamma_f*s*correction = gamma_f*s*k/s*nabla(vf)_f = gamma_f*k*nabla(vf)_f which is just the second part on the right side of formula 3.28 in Hrv's thesis multiplied by the interpolated diffusion coefficient gamma_f. The content of the faceFluxCorrectionPtr is summed up with fvc::div for the faces of every cell. It is multiplied with the cell volume (why?) and goes to the right hand side of the equation into the source term with a minus. If the flux is calculated from the matrix (line 938 in fvMatrix.C) the content of the faceFluxCorrectionPtr is just added. That's what i think is happening for the corrected laplacian of a scalar equation. If you think i got something wrong, i'll be pleased if you propose your correction. Now i have some questions. My first question is as follows: It has often been written in this forum that the approach taken in OF is the overrelaxed approach described in Hrv's PhD Thesis. But if what i wrote above is right, only the explicit part of the overrelaxed approach is used and the implicit part (the matrix coefficients) is taken from the orthogonal correction approach (they should be calculated according to the expression with (phi_Nphi_P) in formula (3.33) and (3.34) in Hrv's thesis, multiplied by gamma_f). But then in gaussLaplacianSchemes.C line 56, fvmLaplacianUncorrected would have to be called with gamma*sk, which is gamma*Delta with the Delta defined according to Hrv's thesis, formula (3.32). Does anybody know why the face area is taken instead of the magnitude of the Delta vector? My second question seems quite simple to answer (and i hope will turn out to be stupid : In gaussLaplacianSchemes.C (line 70) the divergence of the corrections is multiplied by the cell volume as far as i can see. But where does this cell volume come from? In my opinion this multiplication by the cell volume should not be there (i hope this question can easily be answered...). I found no way to formulate the discretization of a transport equation in such a way that the explicit flux contribution is multiplied by the volume. Can anybody give me a hint? Third question: In fvmLaplacianUncorrected there is a coefficient with a positive sign written into the ofdiagonal matrix entries and a contribution with negative sign into the diagonal entries. I would have expexted the signs of the diagonal and the ofdiagonal entries to be the other way around (ofdiagonal entry negative, diagonal positive, as described e.g. in the book by Ferziger and Peric). Can anybody explain? I hope that my question is understandable. If my notation is not clear, feel free to ask. Thanks in advance for your answers or comments. 

December 2, 2012, 22:18 

#2 
Senior Member
Fumiya Nozaki
Join Date: Jun 2010
Location: Yokohama, Japan
Posts: 205
Blog Entries: 1
Rep Power: 10 
Hi koderer,
Thank you for your clarifying this topic. I'm also trying to fully understand the discretization of laplacian term. Below are my answers to your questions: 1) This point is not clear to me. Can somebody help me with it? 2) fvc::div(*fvm.faceFluxCorrectionPtr())().internalF ield() is already divided by the cell volume, so multiplying this term by the cell volume equals Σ_f (gamma_f * k & nabla(vf)_f). 3) The orthogonal contribution is evaluated by gamma_f*s*((phi_Nphi_P))/d, so I think the offdiagonal coefficient is gamma_f*s/d and the diagonal coefficient is Σ_f (gamma_f*s/d) and the code is correct. Please correct the mistakes, if any. Hope that helps, fumiya 

December 12, 2012, 08:19 

#3 
New Member
Join Date: Jan 2012
Posts: 6
Rep Power: 7 
Thank you for your answer, fumiya.
1) I'm glad i'm not the only one wondering. 2) Can you tell me the code lines where this is done? I am still not sure, because in the implicit counter part fvm::div is always used in a way that makes only sense without a built in division by the cell volume. 3) What you write seems reasonable, but the code lines i refer to seem to say the opposite. Can somebody help? 

December 12, 2012, 18:30 

#4 
Senior Member
Fumiya Nozaki
Join Date: Jun 2010
Location: Yokohama, Japan
Posts: 205
Blog Entries: 1
Rep Power: 10 
Hi koderer,
1) In gaussLaplacianScheme.C Code:
00052 tmp<surfaceScalarField> tdeltaCoeffs = 00053 this>tsnGradScheme_().deltaCoeffs(vf); 00054 const surfaceScalarField& deltaCoeffs = tdeltaCoeffs(); Code:
00093 // Return the interpolation weighting factors for the given field 00094 virtual tmp<surfaceScalarField> deltaCoeffs 00095 ( 00096 const GeometricField<Type, fvPatchField, volMesh>& 00097 ) const 00098 { 00099 return this>mesh().nonOrthDeltaCoeffs(); 00100 } one to take into consideration the mesh nonorthogonality. Therefore, the coefficients calculated in the following line are gamma_f * Delta/d, where Delta is the vector defined in the Dr. Jasak's thesis(Eq. (3.32)). Code:
00066 fvm.upper() = deltaCoeffs.internalField()*gammaMagSf.internalField(); I think the fvc::div operation to the surfaceScalarField is caluculated using the surfaceIntegrate function defined in the fvcSurfaceIntegrate.C and in this file Code:
00075 ivf /= mesh.V(); When we solve the pressure poisson equation fvm::laplacian(rAU, p) == fvc::div(phi), the fvc::div(phi) term is calculated in the same way. For explicit source terms, they are incorporated into an equation simply as a field of values. In this example fvc::div(phi) is volScalarField(not multiplied by the cell volume). Last edited by fumiya; December 16, 2012 at 06:11. 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
About nonorthogonal correction in fvclaplacian  7islands  OpenFOAM Running, Solving & CFD  4  December 9, 2012 23:54 
Floating point exception error  Alan  OpenFOAM Running, Solving & CFD  10  April 6, 2012 14:02 
Nonorthogonal correction  maka  OpenFOAM Running, Solving & CFD  9  October 25, 2011 10:55 
using nonorthogonal correction for the laplacian term  stevenvanharen  OpenFOAM Running, Solving & CFD  4  January 10, 2011 08:58 
nonorthogonal correction in fvclaplacian  ssp  OpenFOAM  0  May 4, 2010 04:02 