Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 82
Текст из файла (страница 82)
14.3).Element DElement AFig. 14.3 Schematic of a boundary mesh with stretched elementsSolutionSolution starts by computing the equivalent E factor for the applied relaxationfactorE¼11¼¼51 k 1 0:8thus the pseudo time step for each of the elements isDt ¼ E / Dtq VC¼5 CaCTherefore the relative pseudo time steps for elements A and D can becomputed asDtA¼DtB VCAVCDVCAa CD1¼ 0:05=¼ 0:12a CAa CDVCDa CAimplying that the pseudo time step for element D is nearly 20 times that ofelement A.54414.2.3.314Discretization of the Source Term, Relaxation, and Other DetailsFalse Transient RelaxationThe false transient relaxation method [3] is a modification of the Euler first orderimplicit transient method, wherein the previous iteration values are used instead ofthe old time step values.
As in the Euler method the diagonal dominance of thealgebraic equation is increased through the addition of the pseudo-transient term,aoC /C , to the diagonal coefficient and the pseudo old time step term, aoC /C , to theright hand side. With these modifications the equation becomesXaC þ aoC /C þaF /F ¼ bC þ aoC /Cð14:17ÞFNBðCÞwhere the coefficient aoC is computed asaoC ¼qC VCDtð14:18Þand is basically equal to the transient coefficient obtained from the first orderimplicit Euler discretization of the transient term, qC is the density, VC the elementvolume, and Dt a user-defined false time step.
For large values of Dt, the added termis negligible, and under-relaxation effects are negligible. Therefore the solution ofthe equation is the same as the original unrelaxed one. For very small values of Dt,the value of aoC becomes large dominating other terms. The solution is heavilyunder-relaxed leading to minute change in the value of /C (i.e., /C /C ).In addition to allowing the solution to advance consistently over the entirecomputational domain, the false transient method ensures the addition of a non-zerocontribution to the diagonal coefficient even in extreme cases when the diagonalcoefficient is zero.There is no general rule for assigning optimum under-relaxation factors as valuesused in one case may not work properly for another case.
In addition, differentequations may be assigned different under-relaxation factors. Further, it is notnecessary to use the same under-relaxation value throughout the computationaldomain. Furthermore, under-relaxation values may vary from iteration to iteration.For the SIMPLE algorithm to be described in Chaps. 15 and 16, Raithby andSchneider [4] derived an optimum relation between the under-relaxation factors forthe velocity and pressure fields, which will be presented in the next chapter.14.3Residual Form of the EquationThe discretized algebraic equation has been written so far in its “direct” or “standard” form as14.3Residual Form of the Equation545aC / C þXaF /F ¼ bC ;ð14:19ÞFNBðC Þwhich is the form used in OpenFOAM®.Equation (14.19) can also be written in “correction” or “residual” form byrearranging its terms so as to solve for the correction needed to satisfy the equation.Thus if /C and /0C are the previous iteration value of /C and the correction neededto satisfy Eq.
(14.19), respectively, then the solution is given by/C ¼ /C þ /0Cð14:20Þand Eq. (14.19) can be rewritten asXaC /C þ /0C þaF /F þ /0F ¼ bCð14:21ÞFNBðC ÞoraC /0C þX0aF /0F ¼ bC @aC /C þFNBðC ÞX1aF /F Að14:22ÞFNBðC ÞNote that the right hand side of Eq. (14.22) represents the residual error of theequation for the field /C . Denoting this residual over element C by Res/C ,Eq. (14.22) becomesaC /0C þXaF /0F ¼ Res/Cð14:23ÞFNBðC ÞIt is worth noting that for the exact field Res/C would be zero.While mathematically equivalent to Eq. (14.19), Eq.
(14.23) has one numericaladvantage. In this form, numerical errors during the solution of the equation areslightly less than those associated with the standard form for cases when smallvariations are expected for large values of /.14.3.1 Residual Form of Patankar’s Under-RelaxationThe residual form of the implicit Patankar relaxation equation is derived byrewriting Eq. (14.19) in the following form:54614Discretization of the Source Term, Relaxation, and Other Details0aC /C þ /0C ¼ k/ @bC 1 aF /F þ /0F A þ 1 k/ aC /CXð14:24ÞFNBðC Þwhich can be simplified toaC /0C þ k/X20XaF /0F ¼ k/ 4bC @aC /C þFNBðC Þ13aF /F A5ð14:25ÞFNBðC ÞNoticing that the right hand side of the above equation represents the residual ofthe original equation, Eq. (14.25) can be written asaCk//0C þXaF /0F ¼ Res/Cð14:26ÞFNBðC Þimplying that under-relaxing the equation in residual form necessitates modifyingonly the diagonal coefficient.14.4Residuals and Solution ConvergenceIn any iterative solution process it is important to be able to determine when thesolution can be considered good enough, or when the error can be estimated to bebelow a certain tolerance, or even to what precision the conservation equations havebeen satisfied.
Having the tools to answer any of the above questions is animportant ingredient for any CFD code. This can be rephrased as how to evaluatethe degree of convergence of the solution field without knowing the final solution.To this end a number of indicators were proposed over the years, from the simplemonitoring of a point at a location over a number of iterations, to the monitoring ofan integrated value such as the drag coefficient, the total mass flow, the wall shearstress, and so on, or more commonly the monitoring of some type of equationresidual.
The challenge being that the method has to be used for a wide range offlow parameters and for a variety of geometries and boundary conditions.14.4.1 ResidualsAs a solution to the discretized system of equations represented by Eq. (14.1) issought, the error in the balance equation is quantified by defining the elementresidual error as14.4Residuals and Solution Convergence0Res/C ¼ bC @aC /C þ547X1aF /F A:ð14:27ÞFNBðCÞIt is clear that when the solution is reached and the equation satisfied, the valueof Res/C will be zero. Using Res/C a number of residual indicators for the wholedomain can be derived as detailed next.14.4.2 Absolute ResidualAs defined by Eq.
(14.27), the residual may be a positive or a negative quantity.Since the sign is immaterial, the absolute value of the Res/C , denoted by R/C , is usedto decide whether a solution has converged or not. If R/C decreases with iterationsthen the solution will be converging otherwise it will be diverging. The value of R/Cat point C is defined as01XR/C ¼ bC @aC /C þaF /F AFNBðC Þð14:28Þ14.4.3 Maximum ResidualThe solution is assumed to have converged when the maximum value of theabsolute residuals, defined as,R/C; max01 X@A¼ max bC aC /C þaF /F ¼ max R/Call cells all cellsFNBðC Þð14:29Þover the domain drops below a vanishing quantity e, i.e.,R/C; max e ) solution has converged:ð14:30Þ14.4.4 Root-Mean Square ResidualAnother parameter used as a convergence indicator, is the average of the square rootof the sum of the squares of the absolute residuals R/C;rms , mathematically given by54814R/C; rms ¼Discretization of the Source Term, Relaxation, and Other Detailsvffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi!!2ffiuuPPubC aC /C þaF / FutCall elementsFNBðC Þnumber of elementsvffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiuP / 2uRut Call cells C¼:number of elementsð14:31ÞIn this case the convergence criteria is written asR/C; rms e ) solution has converged:ð14:32Þ14.4.5 Normalization of the ResidualThe level of the absolute residual is a strong function of the variable /.
Thereforedifferent variables result in different levels of R/C . This makes it difficult to discernwhether a solution has converged or not. In such cases a better insight can be gainedthrough scaling the different residuals by dividing them by their respective maximum fluxes. Recalling that aC represents the sum of the fluxes over the element, theresiduals are scaled relative to the local value of the property / to obtain a relativeerror by dividing them by the maximum value of aC /C over the domain such thatPaF /F bC aC /C þFNBðC ÞR/C; scaled ¼:ð14:33Þmax jaC /C jall cellsThe solution is assumed to have converged when the maximum value of thescaled absolute residuals has dropped below a vanishing quantity e, i.e.,max R/C;scaled e ) solution has convergedall cellsð14:34ÞIt is common to require e for scaled residuals to be of the order of 103 to 105or less for convergence.In addition to using either the absolute or scaled residuals, it is always insightfulto monitor integrated quantities, as mentioned above, before concluding that thesolution has converged.
Always ensure proper convergence before declaring that asolution has converged as non-converged solutions can be misleading.14.5Computational Pointers14.5549Computational PointersIn this section, the source term linearization and relaxation techniques used inuFVM and OpenFOAM® are discussed.14.5.1 uFVM14.5.1.1Source Term LinearizationThe linearization and assembly of the source term in uFVM is implemented in thefunction cfdAssembleSourceTerm displayed in Listing 14.1. The function relies onthe linearization set by the user, which has to supply the constant part of the source,Sb, and the linearized part, Sc. These terms are then added to FLUXV and FLUXC,as in Eq. (14.4).It is worth noting that in uFVM the equations are solved in residual form, whichexplains the implementation of the total source in FLUXTE, rather than its constantpart only.theEquationField = cfdGetMeshField(theEquationName);phi = theEquationField.phi(iElements);%Sb = cfdComputeFormulaAtLocale(theTerm.Sb,'Interior Elements')';Sc = cfdComputeFormulaAtLocale(theTerm.Sc,'Interior Elements')';%volume = [theMesh.elements.volume];%% Assemble Source Term%pos = zeros(1,size(phi));pos(Sc<0) = 1;theFluxes.FLUXCE = -pos .* Sc .* volume;theFluxes.FLUXTE = -(Sb +Sc .*phi) .* volume;Listing 14.1 Linearization and implementation of the source term14.5.1.2Under-RelaxationThe implicit under-relaxation method of Patankar is implemented in ufvm.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
