[hrishikeshsivanandan]

Dear all,
I am faced with the problem of having bad math being solved upon trying to set up Rayleigh damping layers in the domain. The error is as follows:
#1 Foam::sigFpe::sigHandler(int) at ??:?
#2 ? in "/lib/x86_64-linux-gnu/libc.so.6"
#3 Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) at ??:?
#4 Foam::symGaussSeidelSmoother::smooth(Foam::Field<d ouble>&, Foam::Field<double> const&, unsigned char, int) const at ??:?
#5 Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#6 Foam::fvMatrix<double>::solveSegregated(Foam::dict ionary const&) at ??:?
#7 Foam::fvMatrix<double>::solve(Foam::dictionary const&) in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"
#8 Foam::fvMatrix<double>::solve() in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"
#9 5? in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"
#10 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#11 ?0 in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"


Can someone help?

[mAlletto]

Can you provide a more detailed description. At least how was your simulation set up and how did you implement your damping.

I think it is quite obvious that with a description you provide nobody can infer the cause of your problem

[HPE]

agree with mAlletto.

please provide sufficient detail for us to understand your problem.

[hrishikeshsivanandan]

I am trying to simulate buoyancy waves in large wind farms using a custom made modified k-\epsilon CFD solver based on OpenFOAM and SOWFA environments.

With regard to the boundary conditions,
I am trying to implement timeVaryingMappedFixedValues at the inlet face for U,epsilon,k,T. Further, the top wall is modelled as a slip patch for all variables except T, which is set as a fixed Gradient matching the atmospheric lapse rate.
At the exit, all the variables are set to zeroGradient.
The lateral faces are modelled as cyclic boundaries for all variables.

The initial field for these variables are set to match the atmospheric flow called based on variation of these flow variables with height.

Further, the Rayleigh damping layers are set at the top of the domain having a 5km thickness and Rayleigh damping (c_{max} of 0.00046
Mathematically, these are of the form:
F_damping = c(U_{ref}-U), -----1
where c = 0.5*c_{max}(1-f*cos(\pi(z-z_{start})/w}-----2,
where in 2, 'w' is the width of the inversion layer, z_start is the altitude at which the layer starts and z is the altitude it self. Further f is +1 ( if set as step up) and -1 (if set as step down). c_{max} is the maximum viscosity (set by the user).

After a few iterations, I am seeing the errors, I mentioned earlier.

[HPE]

could you please attach the full error log and checkMesh output as well?

[hrishikeshsivanandan]

The checkMesh output is as follows:
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 6
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
Build : 6-fa1285188035
Exec : checkMesh
Date : Feb 10 2021
Time : 22:10:32
Host : "HRISHIKESH-G3-3579"
PID : 24750
I/O : uncollated
Case : /media/hrishi0509/DATA/Thesis/Simulations/thesis_hrishi/Julia_Tutorials/FullFarm1/FullFarm_UniformT/onlytop_0.0046_top.in.out
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 10)
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 513744
faces: 1432475
internal faces: 1362900
cells: 459800
faces per cell: 6.07954545455
boundary patches: 6
point zones: 1
face zones: 3
cell zones: 25

Overall number of cells of each type:
hexahedra: 454575
prisms: 0
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 0
polyhedra: 5225
Breakdown of polyhedra by number of faces:
faces number of cells
13 5225

Checking topology...
Boundary definition OK.
Cell to face addressing OK.
Point usage OK.
Upper triangular ordering OK.
Face vertices OK.
Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces...
Patch Faces Points Surface topology
lower 41800 42616 ok (non-closed singly connected)
upper 5225 5376 ok (non-closed singly connected)
west 2915 3024 ok (non-closed singly connected)
east 2915 3024 ok (non-closed singly connected)
north 8360 9174 ok (non-closed singly connected)
south 8360 9174 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (0 0 0) (19000 5500 10000)
Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
Mesh has 3 solution (non-empty) directions (1 1 1)
Boundary openness (-1.09793055398e-15 1.80618621633e-16 0) OK.
Max cell openness = 1.0848454272e-16 OK.
Max aspect ratio = 4 OK.
Minimum face area = 2400. Maximum face area = 70908.956528. Face area magnitudes OK.
Min volume = 240000. Max volume = 7090895.6528. Total volume = 1.045e+12. Cell volumes OK.
Mesh non-orthogonality Max: 42.3741905455 average: 5.00567636831
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 1.6864796634 OK.
Coupled point location match (average 0) OK.

Mesh OK.

End

[1]+ Done paraview --mesa dum.foam
************************************************** *************************************************

Further, the last few iterations in the solver log file is as follows:
Time = 80 Time Step = 80
updating force
updating force
returning force
smoothSolver: Solving for Ux, Initial residual = 0.389550534462, Final residual = 0.0334176710723, No Iterations 1
smoothSolver: Solving for Uy, Initial residual = 0.342880967445, Final residual = 0.0247859399143, No Iterations 1
smoothSolver: Solving for Uz, Initial residual = 0.780889902598, Final residual = 0.0720472041904, No Iterations 1
DICPCG: Solving for p_rgh, Initial residual = 0.0618551876474, Final residual = 0.000616615296909, No Iterations 49
Buoyancy term B: min(kEpsilonABL:B) [0 2 -3 0 0 0 0] -4.03194312755e+12, max(kEpsilonABL:B) [0 2 -3 0 0 0 0] 0.000921954386475
smoothSolver: Solving for epsilon, Initial residual = 0.890995287383, Final residual = 0.0484526507249, No Iterations 3
bounding epsilon, min: -0.0103846954547 max: 0.514556510583 average: 0.00046377757542
smoothSolver: Solving for k, Initial residual = 0.884980518772, Final residual = 0.0774153753853, No Iterations 2
bounding k, min: -5.98731119167 max: 63.0293175404 average: 0.0733353330506
returning force
smoothSolver: Solving for T, Initial residual = 0.0620273374252, Final residual = 0.00293887167671, No Iterations 1
Continuity Report:
-Local Cell Continuity Error:
minimum: 1.86014128849e-12
maximum: 0.0868597619827
weighted mean: 4.50138167532e-06
-Boundary Flux:
lower - flux: -1.6115109247e-11 / area: 104500000
upper - flux: 4.36557456851e-10 / area: 104500000
west - flux: -548072873.397 / area: 55000000
east - flux: 548072873.397 / area: 55000000
total - flux: -4.76837158203e-07 / area: 319000000
ExecutionTime = 770.55 s ClockTime = 888 s

Time = 81 Time Step = 81
updating force
updating force
returning force
smoothSolver: Solving for Ux, Initial residual = 0.434682013836, Final residual = 0.0416011959379, No Iterations 1
smoothSolver: Solving for Uy, Initial residual = 0.306563102169, Final residual = 0.0227160269001, No Iterations 1
smoothSolver: Solving for Uz, Initial residual = 0.708224305329, Final residual = 0.0677720940285, No Iterations 1
DICPCG: Solving for p_rgh, Initial residual = 0.966261251138, Final residual = 0.0084277138978, No Iterations 21
Buoyancy term B: min(kEpsilonABL:B) [0 2 -3 0 0 0 0] -4.7258784106e+14, max(kEpsilonABL:B) [0 2 -3 0 0 0 0] 0.001154448976
smoothSolver: Solving for epsilon, Initial residual = 0.990726867968, Final residual = 0.0582640849199, No Iterations 2
bounding epsilon, min: -0.0210408648249 max: 68964.0665679 average: 21.2863323864
smoothSolver: Solving for k, Initial residual = 0.00180223419247, Final residual = 7.46755985941e-05, No Iterations 2
bounding k, min: -0.476577612178 max: 1508422025.42 average: 35385.5552345
returning force
smoothSolver: Solving for T, Initial residual = 0.0282149934703, Final residual = 0.00213460899956, No Iterations 1
Continuity Report:
-Local Cell Continuity Error:
minimum: 3.96641527362e-10
maximum: 439225.270708
weighted mean: 271.363935157
-Boundary Flux:
lower - flux: 1.32445165946e-11 / area: 104500000
upper - flux: 9.60426405072e-10 / area: 104500000
west - flux: -548072873.397 / area: 55000000
east - flux: 548072873.397 / area: 55000000
total - flux: 1.66893005371e-06 / area: 319000000
ExecutionTime = 779.31 s ClockTime = 898 s

Time = 82 Time Step = 82
updating force
updating force
returning force
smoothSolver: Solving for Ux, Initial residual = 0.363114559086, Final residual = 0.00963046907832, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 0.0414754899336, Final residual = 0.0031153025679, No Iterations 1
smoothSolver: Solving for Uz, Initial residual = 0.233946049796, Final residual = 0.00702479250304, No Iterations 2
DICPCG: Solving for p_rgh, Initial residual = 0.998633687394, Final residual = 0.00966257801876, No Iterations 221
Buoyancy term B: min(kEpsilonABL:B) [0 2 -3 0 0 0 0] -1.76131734081e+15, max(kEpsilonABL:B) [0 2 -3 0 0 0 0] 486.869602501
smoothSolver: Solving for epsilon, Initial residual = 0.25359769274, Final residual = 0.00782482435813, No Iterations 2
smoothSolver: Solving for k, Initial residual = 0.988772110512, Final residual = 0.0387214151603, No Iterations 2
bounding k, min: 6.66133814775e-17 max: 2.3427438967e+18 average: 1.74886814659e+14
returning force
smoothSolver: Solving for T, Initial residual = 0.000298984998007, Final residual = 2.9805860157e-05, No Iterations 1
Continuity Report:
-Local Cell Continuity Error:
minimum: 1.13880505431e-10
maximum: 178386.835691
weighted mean: 256.820513864
-Boundary Flux:
lower - flux: -1.63265993344e-11 / area: 104500000
upper - flux: -2.18278728426e-10 / area: 104500000
west - flux: -548072873.397 / area: 55000000
east - flux: 548072873.397 / area: 55000000
total - flux: -7.15255737305e-07 / area: 319000000
ExecutionTime = 791.06 s ClockTime = 911 s

Time = 83 Time Step = 83
updating force
updating force
returning force
smoothSolver: Solving for Ux, Initial residual = 0.589327392615, Final residual = 0.0241378702493, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 0.290985231987, Final residual = 0.0211399670241, No Iterations 1
smoothSolver: Solving for Uz, Initial residual = 0.38420797035, Final residual = 0.0176802020221, No Iterations 2
DICPCG: Solving for p_rgh, Initial residual = 0.999962424413, Final residual = 0.00716361177681, No Iterations 953
Buoyancy term B: min(kEpsilonABL:B) [0 2 -3 0 0 0 0] -4.40708302596e+18, max(kEpsilonABL:B) [0 2 -3 0 0 0 0] 1.37061354786e+12
smoothSolver: Solving for epsilon, Initial residual = 0.00721233532312, Final residual = 0.000235223942019, No Iterations 1
smoothSolver: Solving for k, Initial residual = 2.41697900931e-05, Final residual = 2.24859772384e-06, No Iterations 2
bounding k, min: 6.66133814775e-17 max: 6.39922765685e+20 average: 2.19251575933e+16
returning force
smoothSolver: Solving for T, Initial residual = 4.80405221942e-08, Final residual = 5.92868648753e-09, No Iterations 1
Continuity Report:
-Local Cell Continuity Error:
minimum: 0.00420573618511
maximum: 1.54957233121e+13
weighted mean: 330762106.742
-Boundary Flux:
lower - flux: 1.75349864263e-17 / area: 104500000
upper - flux: 2.50910403565e-14 / area: 104500000
west - flux: -548072873.397 / area: 55000000
east - flux: 547427840 / area: 55000000
total - flux: -645033.396691 / area: 319000000
ExecutionTime = 811.62 s ClockTime = 934 s

Time = 84 Time Step = 84
updating force
updating force
returning force
smoothSolver: Solving for Ux, Initial residual = 0.588876408857, Final residual = 0.0100423090454, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 0.569679092668, Final residual = 0.0547516746918, No Iterations 1
smoothSolver: Solving for Uz, Initial residual = 0.0894108285395, Final residual = 0.00145537039155, No Iterations 2
DICPCG: Solving for p_rgh, Initial residual = 0.920701158782, Final residual = 13.1231897405, No Iterations 1000
Buoyancy term B: min(kEpsilonABL:B) [0 2 -3 0 0 0 0] -6.22981657465e+25, max(kEpsilonABL:B) [0 2 -3 0 0 0 0] 2.15298589571e+23
smoothSolver: Solving for epsilon, Initial residual = 0.000134228206104, Final residual = 5.67280762781e-06, No Iterations 2
smoothSolver: Solving for k, Initial residual = 1.30443379764e-13, Final residual = 1.30443379764e-13, No Iterations 0
returning force
#0 Foam::error:rintStack(Foam::Ostream&) at ??:?
#1 Foam::sigFpe::sigHandler(int) at ??:?
#2 ? in "/lib/x86_64-linux-gnu/libc.so.6"
#3 Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) at ??:?
#4 Foam::symGaussSeidelSmoother::smooth(Foam::Field<d ouble>&, Foam::Field<double> const&, unsigned char, int) const at ??:?
#5 Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#6 Foam::fvMatrix<double>::solveSegregated(Foam::dict ionary const&) at ??:?
#7 Foam::fvMatrix<double>::solve(Foam::dictionary const&) in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"
#8 Foam::fvMatrix<double>::solve() in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"
#9 ? in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"
#10 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#11 ? in "/home/hrishi0509/OpenFOAM/hrishi0509-6/platforms/linux64GccDPInt32Opt/bin/steadySolverWithoutTurbine"

[mAlletto]

What do you mean with step up and step down?

Is the damping source explicit or implicit?

By the way by clicking on the sigma symbol you can write in math mode and clicking on the # symbol you can write in code mode. This leads to a much more readable post

[Regis_]

Hrishikesh,
I implemented the Rayleigh damping layers in SOWFA and can help you debug your case. I see that you have a custom solver called `steadySolverWithoutTurbine`, though. Were you able to run the standard SOWFA solver with the damping layers? In any case, explain better what you did, what you modified, and we can go from there.