
[Sponsors] 
January 15, 2016, 07:42 
Energy balance in solid using chtMRSF

#1 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
Hi all,
I simulate a heat exchanger using chtMultiRegionSimpleFoam and OF 2.4.0. The results are somewhat realistic but I think there might be a problem with thermal convergence in the SOLID because the heat balance doesn't sum up. heat balance solid: The heat transferred at water side (blue) fluctuates and tops the heat transferred from exhaust gas side (red) which makes no sense for me Other parameters as global temperature increase, pressure drop and the residuals point to a converging solution (become constant) although the maximum temperature in the fluid rises continously. I measure the transferred heat with a functionObject Code:
energyAbsorbedSOLID { type patchExpression; outputControlMode timeStep; outputInterval 1; accumulations ( sum ); patches (".*"); region SOLID; variables ("kvar=140;" ); expression "kvar*(snGrad(T))*area()"; verbose true; } 

January 15, 2016, 07:47 

#2 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
Another thing are the high time step continuity errors, which decrease but are still quite high
Could somebody point me to a good definition of those values so I can understand them better? 

January 26, 2016, 07:47 

#3 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
It seems to be a matter of convergence.
If I continue the run it becomes constant with a small difference. Increasing nNonOrthCorrectors in the Solid was also helpful, since more iterations are done for the solid then, which let the case converge faster. 

January 26, 2016, 10:44 

#4  
Senior Member
Join Date: Sep 2013
Posts: 353
Rep Power: 19 
Quote:
Check your residuals (v2.4.0 has a neat function object for this). If there is no change check the log file if the h equation does iterate. You might want to use minIter 1 to be safe. The residuals alone however are often misleading. Always use the wallHeatFlux utility to check if the heatflux from one region to the next is the same. Then check again in paraview. If the heat flux isn't smooth your mesh might just be to coarse (always check the y+ value. I'd recommend a value below 1). Problems arise if the mesh isn't identical on bothsides which can be significant. Using prism layers helps alot. Temperature on a few patches or min/max values should also be checked. For all of those there are function objects to use. In your case it is probably the relaxation factor. The relaxation factor for h (or e depending on your thermoPhysicalProperties) has a huge influence on convergence speed. I only use them to stop the simulation from crashing in the first few timesteps. I'd fully remove them if possible (in both regions). Another reason might be your mesh or discretisation. 

January 27, 2016, 02:58 

#5 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
Hi Stephan,
I am using OF 3.0 now, just to avoid confusion. fvSchemes for WATER Code:
ddtSchemes { default steadyState; } gradSchemes { default Gauss linear; grad(U) cellLimited Gauss linear 1; } divSchemes { default none; div(phi,U) bounded Gauss upwind; div(phi,h) bounded Gauss upwind; div(phi,k) bounded Gauss upwind; div(phi,K) bounded Gauss upwind; div(phi,epsilon) bounded Gauss upwind; div(phi,omega) bounded Gauss upwind; div(phi,R) bounded Gauss upwind; div(R) Gauss linear; div((muEff*dev2(T(grad(U))))) Gauss linear; div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; } laplacianSchemes { default Gauss linear limited 0.333; } interpolationSchemes { default linear; } snGradSchemes { default limited 0.333; } fluxRequired { default no; p_rgh; } wallDist { method meshWave; } Code:
solvers { rho { solver PCG preconditioner DIC; tolerance 1e7; relTol 0; } p_rgh { solver GAMG; tolerance 1e7; relTol 0.01; smoother DIC; cacheAgglomeration true; nCellsInCoarsestLevel 3126; agglomerator faceAreaPair; mergeLevels 1; maxIter 100; } "(Uhkepsilonomega)" { solver PBiCG; preconditioner DILU; tolerance 1e7; relTol 0.1; } } SIMPLE { momentumPredictor on; nNonOrthogonalCorrectors 2; pRefCell 0; pRefValue 100000; rhoMin rhoMin [1 3 0 0 0] 1000; rhoMax rhoMax [1 3 0 0 0] 1100; } relaxationFactors { fields { rho 1; p_rgh 0.5; } equations { U 0.3; h 0.7; nuTilda 0.7; k 0.7; epsilon 0.7; omega 0.7; "ILambda.*" 0.7; } } relaxationFactors { fields { rho 1; p_rgh 0.4; } equations { U 0.5; h 0.9; nuTilda 0.9; k 0.9; epsilon 0.9; omega 0.9; "ILambda.*" 0.9; } } fvSchemes for SOLID Code:
ddtSchemes { default steadyState; } gradSchemes { default Gauss linear; } divSchemes { default none; } laplacianSchemes { default none; laplacian(alpha,h) Gauss linear corrected; div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; } interpolationSchemes { default linear; } snGradSchemes { default uncorrected; } fluxRequired { default no; } wallDist { method meshWave; } Code:
solvers { h { solver PCG; preconditioner DIC; tolerance 1e06; relTol 0.1; } } SIMPLE { nNonOrthogonalCorrectors 10; } relaxationFactors { fields { } equations { h 0.7; } } 

January 27, 2016, 03:11 

#6 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
My residuals go down and become constant. For h at 1e6 but it takes sometime.
I also monitored pressure drop, temperature increase and maximum temperature in fluid and solid which become constant as well. Heat flux I checked with a functionObject for the SOLID region Code:
energyAbsorbedWATER { type patchExpression; outputControlMode timeStep; outputInterval 1; accumulations ( sum ); patches ("SOLID_TO_WATER"); region WATER; variables ("kvar=0.5;" ); expression "(snGrad(T))*area()"; verbose true; } The average yPlus is below one for all walls. 

January 27, 2016, 03:15 

#7 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
I also wanted to check the heat balance in the fluid.
And to obtain the enthalpy flow difference between IN and OUTLET I used Code:
energyOutlet { type patchExpression; outputControlMode timeStep; outputInterval 1; accumulations ( sum ); patches ("OUTLET"); region WATER; variables ( ); expression "(A+B*TC*sqr(T))*T*phi"; verbose true; } But there seems to be some difference between what is going into the water and what I get from subtracting enthalpy flow for inlet and outlet. 

January 27, 2016, 11:31 

#8 
Senior Member
Join Date: Sep 2013
Posts: 353
Rep Power: 19 
If your residuals become constant (in one region) make sure that this is not due to not solving anymore.
Code:
Solving for h, Initial residual = .., Final residual = .., No Iterations 0 Code:
solvers { h { solver PCG; preconditioner DIC; tolerance 1e08; relTol 0.1; minIter 1; } } SIMPLE { nNonOrthogonalCorrectors 3; } relaxationFactors { fields { } equations { h 1.0; } } Do the same problems occur whitout using polynomial transport? 

January 28, 2016, 03:10 

#9  
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
1) Number of iterations
Quote:
You are right, thats what happens after some time. Code:
Solving for solid region SOLID DICPCG: Solving for h, Initial residual = 9.999899e07, Final residual = 9.999899e07, No Iterations 0 DICPCG: Solving for h, Initial residual = 9.999774e07, Final residual = 9.999774e07, No Iterations 0 DICPCG: Solving for h, Initial residual = 9.999774e07, Final residual = 9.999774e07, No Iterations 0 DICPCG: Solving for h, Initial residual = 9.999774e07, Final residual = 9.999774e07, No Iterations 0 DICPCG: Solving for h, Initial residual = 9.999774e07, Final residual = 9.999774e07, No Iterations 0 DICPCG: Solving for h, Initial residual = 9.999774e07, Final residual = 9.999774e07, No Iterations 0 DICPCG: Solving for h, Initial residual = 9.999774e07, Final residual = 9.999774e07, No Iterations 0 [...] Code:
residualControl { p_rgh 1e6; U 1e8; T 1e6; // possibly check turbulence fields "(kepsilonomega)" 1e7; } Or how I can add a minimum number of iterations? 

January 28, 2016, 12:24 

#10 
Senior Member
Join Date: Sep 2013
Posts: 353
Rep Power: 19 
This is done as mentioned above. You just need to lower your tolerance or add minIter (again highlighted below) in your fvSolution files. This is one of the reasons heat transfer simulations should not be judged by residual alone. Much more import are the temperature values or the heat fluxes. Residual Control only aborts a calculation once the inital residual reaches a certain value. The solver on the other side stops solving once the tolerance or the relative tolerance from one time step to the next is reached . It is however correct that residual control is not implemented in this solver, maybe in part for this reason.
Code:
solvers { h { .... tolerance 1e08; relTol 0.1; minIter 1; } } relaxationFactors { fields { } equations { h 1.0; } } 

February 8, 2016, 11:17 

#11 
Member
Join Date: May 2015
Posts: 68
Rep Power: 9 
Hi,
your suggestions have been pure gold as they worked out for my case pretty good. But now I am having another question: It seems like I am using a significant amount of energy (17 W) when coupling solid and fluid region with compressible::turbulentTemperatureCoupledBaffleMix ed. Some loss makes sense since I have Arbitrary Mesh Interfaces (AMI) which I couple using nearestPatchFaceAMI. But 17 W is too much and not confirmed by CFX calculations. So I wonder whether it is a matter of thermophysicalProperties I am using temperature dependent values for heat capacity (hPolynomial) and coupling is specified as Code:
{ type compressible::turbulentTemperatureCoupledBaffleMixed; value uniform 348.15; Tnbr T; kappa solidThermo; //fluidThermo kappaName none; } Is solidThermo temperature dependent or should I use the lookup entry? Description of coupling condition Code:
Description Mixed boundary condition for temperature, to be used for heattransfer on backtoback baffles. Specifies gradient and temperature such that the equations are the same on both sides:  refGradient = zero gradient  refValue = neighbour value  mixFraction = nbrKDelta / (nbrKDelta + myKDelta()) where KDelta is heattransfer coefficient K * deltaCoeffs Example usage: myInterfacePatchName { type compressible::turbulentTemperatureCoupledBaffleMixed; neighbourFieldName T; kappa lookup; kappaName kappa; value uniform 300; } Needs to be on underlying mapped(Wall)FvPatch. Note: kappa : heat conduction at patch. Gets supplied how to lookup calculate kappa:  'lookup' : lookup volScalarField (or volSymmTensorField) with name  'fluidThermo' : use fluidThermo and compressible::RASmodel to calculate kappa  'solidThermo' : use solidThermo kappa()  'directionalSolidThermo' directionalKappa() Note: runs in parallel with arbitrary decomposition. Uses mapped functionality to calculate exchange. 

February 8, 2016, 14:03 

#12 
Senior Member
Join Date: Sep 2013
Posts: 353
Rep Power: 19 
solidThermo is correct for solid regions and fluidThermo correct for fluid regions. This means that kappa is determined via the thermophysical model used. (In your case the polynomial ones) Since the "turbulentTemperatureCoupledBaffleMixed" boundary condition can be used with other solvers and solvers written by yourself it has the lookup functionality. lookup sets kappa to the value specified in the kappa file. (You need to create one and probably modify the chtmultiregion solvers to use it). The name of the file needs to be set by the kappaName entry. Which is why it right now is set to none. If no file is present lookup should however be the same as solidThermo.
One way to check where your looses occur is to use the mapFields utility to map one boundary field to the neighbour region and afterwards subtract them with foamCalc. With that you should see problematic areas. The syntax for this is a bit complicated but you could give it a try. Are the results using the nearestPatchFaceAMI better than with the other options available or why did you switch? And can you post a closeup of the boundary between the regions so I can get a feeling for the differences between the mesh points of the two regions? I do not expect your polynomial to influence your losses. 

Tags 
chtmrsf, chtmultiregionsimplefoam, convergence, multi regions 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
HELP needed: Energy balance in gaseous combustion  James Willie  FLUENT  11  August 10, 2005 04:58 
UDS for solid enegy and solid mass balance  Arvind Phate  FLUENT  2  July 20, 2005 22:53 
UDF Scalar Code: HT 1  Greg Perkins  FLUENT  8  October 20, 2000 12:40 
UDFs for Scalar Eqn  Fluid/Solid HT  Greg Perkins  FLUENT  0  October 13, 2000 23:03 
Why FVM for highRe flows?  Zhong Lei  Main CFD Forum  23  May 14, 1999 13:22 