[Vyssion]

Hello,

I am trying to simulate a fan rotating within a fixed duct just floating in space. As it is currently set up, the fan is 1m in diameter and sits in the centre of a 4x4x4m cube domain.

I am using snappyHexMesh and pimpleDyMFoam to try and simulate this rotation, and seem to have it working for the most part... until about 50 ish iterations in, when the time step begins to approach zero (i.e. 1e-13) and the maximum courant number increases above the defined maximum I have set. (I have said "weird" in the title because the average courant number that I have is on the order of 1e-05 give or take)... At this point, the simulation hits a floating point error and terminates. I am experienced in commercial packages for CFD, but have never used OpenFOAM prior to a fortnight ago, and so am looking for some help to get it running smoothly.

I should add that I also have tried a simple MRF setup (which I have managed to converge) however, what I dont like about this option is that I have to state the Cp and Ct values within the MRFProperties script when I am looking to measure the "lift" generated just by rotating some .stl geometry.


Here is the error I receive:

Courant Number mean: 9.258629e-06 max: 1.998997
deltaT = 1.771031e-06
Time = 5.24984e-05

solidBodyMotionFunctions::rotatingMotion::transformation(): Time = 5.249842e-05 transformation: ((0 0 0) (0.9999862 (0 -0.005249818 0)))
PIMPLE: iteration 1
smoothSolver:  Solving for Ux, Initial residual = 0.0005662521, Final residual = 6.517957e-13, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 5.593663e-05, Final residual = 3.318731e-14, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 0.001630861, Final residual = 1.829808e-12, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.8991979, Final residual = 0.008066495, No Iterations 1
[5] #0  Foam::error::printStack(Foam::Ostream&) at ??:?
[5] #1  Foam::sigFpe::sigHandler(int) at ??:?
[5] #2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
[5] #3  Foam::GAMGSolver::scale(Foam::Field<double>&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, Foam::Field<double> const&, unsigned char) const at ??:?
[5] #4  Foam::GAMGSolver::Vcycle(Foam::PtrList<Foam::lduMatrix::smoother> const&, Foam::Field<double>&, Foam::Field<double> const&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::PtrList<Foam::Field<double> >&, Foam::PtrList<Foam::Field<double> >&, unsigned char) const at ??:?
[5] #5  Foam::GAMGSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
[5] #6  Foam::fvMatrix<double>::solveSegregated(Foam::dictionary const&) at ??:?
[5] #7  Foam::fvMatrix<double>::solve(Foam::dictionary const&) at ??:?
[5] #8  ? at ??:?
[5] #9  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
[5] #10  ? at ??:?
[vmware-vm:07270] *** Process received signal ***
[vmware-vm:07270] Signal: Floating point exception (8)
[vmware-vm:07270] Signal code:  (-6)
[vmware-vm:07270] Failing at address: 0x3e800001c66
[vmware-vm:07270] [ 0] /lib/x86_64-linux-gnu/libc.so.6(+0x352f0) [0x7fa0de9382f0]
[vmware-vm:07270] [ 1] /lib/x86_64-linux-gnu/libc.so.6(gsignal+0x37) [0x7fa0de938267]
[vmware-vm:07270] [ 2] /lib/x86_64-linux-gnu/libc.so.6(+0x352f0) [0x7fa0de9382f0]
[vmware-vm:07270] [ 3] /opt/openfoam30/platforms/linux64GccDPInt32Opt/lib/libOpenFOAM.so(_ZNK4Foam10GAMGSolver5scaleERNS_5FieldIdEES3_RKNS_9lduMatrixERKNS_10FieldFieldIS1_dEERKNS_8UPtrListIKNS_17lduInterfaceFieldEEERKS2_h+0xa2) [0x7fa0dfb901b2]
[vmware-vm:07270] [ 4] /opt/openfoam30/platforms/linux64GccDPInt32Opt/lib/libOpenFOAM.so(_ZNK4Foam10GAMGSolver6VcycleERKNS_7PtrListINS_9lduMatrix8smootherEEERNS_5FieldIdEERKS8_S9_S9_S9_S9_S9_RNS1_IS8_EESD_h+0x783) [0x7fa0dfb93073]
[vmware-vm:07270] [ 5] /opt/openfoam30/platforms/linux64GccDPInt32Opt/lib/libOpenFOAM.so(_ZNK4Foam10GAMGSolver5solveERNS_5FieldIdEERKS2_h+0x671) [0x7fa0dfb95d21]
[vmware-vm:07270] [ 6] /opt/openfoam30/platforms/linux64GccDPInt32Opt/lib/libfiniteVolume.so(_ZN4Foam8fvMatrixIdE15solveSegregatedERKNS_10dictionaryE+0x16c) [0x7fa0e257e1ec]
[vmware-vm:07270] [ 7] pimpleDyMFoam(_ZN4Foam8fvMatrixIdE5solveERKNS_10dictionaryE+0x199) [0x4575c9]
[vmware-vm:07270] [ 8] pimpleDyMFoam() [0x41ee71]
[vmware-vm:07270] [ 9] /lib/x86_64-linux-gnu/libc.so.6(__libc_start_main+0xf0) [0x7fa0de923a40]
[vmware-vm:07270] [10] pimpleDyMFoam() [0x422529]
[vmware-vm:07270] *** End of error message ***
--------------------------------------------------------------------------
mpirun noticed that process rank 5 with PID 7270 on node vmware-vm exited on signal 8 (Floating point exception).

Here is my snappyHexMesh.dict

// Which of the steps to run
castellatedMesh		true;
snap			true;
addLayers		true;

// GEOMETRY - DEFINITION OF ALL SURFACES
geometry
{
	fan_blades.stl			// Filename
	{
		type	triSurfaceMesh;	// Type
		name	fan_blades;	// Name to be referred to as
		regions
		{
			fan_blades	// Name of Solid which is inside the .stl
			{
				name	fan_blades;
			}
		}
	}

	duct.stl			// Filename
	{
		type	triSurfaceMesh;	// Type
		name	duct;		// Name to be referred to as
		regions
		{
			duct		// Name of Solid which is inside the .stl
			{
				name	duct;
			}
		}
	}

	refinementCylOuter		// User Defined Region
	{
		type	searchableCylinder;
		point1	(0 -0.15 0);	// Base
		point2	(0 0.3 0);	// Top
		radius	0.65;
	}

	refinementCylinder		// User Defined Region
	{
		type	searchableCylinder;
		point1	(0 -0.1 0);    	// Base
		point2	(0 0.25 0);    	// Top
		radius	0.5;       
	}

	rotatingCylinder			// User Defined Region
	{
		type	searchableCylinder;
		point1	(0 -0.05 0);
		point2	(0 0.15 0);
		radius	0.5;
		name	rotatingCylinder;
	}
};



// SETTINGS FOR castellatedMeshControls GENERATION

castellatedMeshControls
{
// REFINEMENT PARAMETERS
// If local number of cells is >= maxLocalCells on any processor switches from
// refinement followed by balancing (current method) to (weighted) balancing
// before refinement.
	maxLocalCells		200000;

// Overall cell limit (approximately). Refinement will stop immediately
// upon reaching this number.
// NOTE: that this is the number of cells before removing the part which
// 	 is not 'visible' from the keepPoint. The final number of cells might
// 	 actually be a lot less.
	maxGlobalCells		15000000;

// The surface refinement loop might spend lots of iterations refining just a
// few cells. This setting will cause refinement to stop if <= minimumRefine
// are selected for refinement. Note: it will at least do one iteration
// (unless the number of cells to refine is 0)
	minRefinementCells	5;

// Allow a certain level of imbalance during refining
// (since balancing is quite expensive)
// Expressed as fraction of perfect balance (= overall number of cells /
// nProcs). 0=balance always.
	maxLoadUnbalance	0.10;

// Number of buffer layers between different levels.
// 1 means normal 2:1 refinement restriction, larger means slower refinement.
	nCellsBetweenLevels	3;


// EXPLICIT FEATURE EDGE REFINEMENT
// Specifies a level for any cell intersected by its edges.
// This is a featureEdgeMesh, read from constant/triSurface for now.
	features
	(
		{
			file	"fan_blades.eMesh";	// File Containing Mesh
			level	9;			// Level of Refinement
		}
		{
			file	"duct.eMesh";		// File Containing Mesh
			level	8;			// Level of Refinement
		}
	);


// SURFACE BASED REFINEMENT
// Specifies two levels for every surface. The first is the minimum level,
// every cell intersecting a surface gets refined up to the minimum level.
// The second level is the maximum level. Cells that 'see' multiple
// intersections where the intersections make an
// angle > resolveFeatureAngle get refined up to the maximum level.
	// DO NOT PUT THE FILE NAME HERE!!!!!! 
	// Inside the .stl file there will be a solid defined with a given name
	// THAT is the string that you type here...
	refinementSurfaces
	{
		fan_blades
		{
			level	(6 8);	// (<min> <max>) Surface Refinement
			patchInfo	// Optional specification of patch type
			{		// No constraint types (e.g. symm) allowed.
				type	wall;
			}
		}

		duct
        	{
			level	(5 8);
			patchInfo
			{
				type	wall;
			}
		}

		rotatingCylinder
		{
			level	(6 8);
			faceZone 	rotatingCylinder;
			cellZone 	rotatingCylinder;
			cellZoneInside	inside;
		}
	}


// RESOLVE SHARP ANGLES
// Cells above angle use the max ref level, cells below use min
	resolveFeatureAngle	30;


// REGION-WISE REFINEMENT
// Specifies refinement level for cells in relation to a surface. One of
// three modes
// 1.	distance	'levels' specifies per distance to the surface the
//			wanted refinement level. The distances need to be 
//			specified in descending order.
// 2.	inside		'levels' is only one entry and only the level is used.
//   			All cells inside the surface get refined up to the level. 
//   			The surface needs to be closed for this to be possible.
// 3.	outside		Same as inside, but cells outside.

	refinementRegions		// See "geometry" sub-Dictionary for
	{				// User defined refinement regions
		refinementCylOuter
		{
			mode	inside;		// inside, outside, distance?
			levels	((1e15 6));	// FORMAT (<distance> <level>)
				// NOTE: for "inside" and "outside", <distance>
				//	 is ignored, but MUST be specified!!
		}

		refinementCylinder
		{
			mode	inside;	
			levels	((1e15 7));
		} 
	
		fan_blades
		{
			mode	distance;	// Refinement level <l> 
			levels	((0.05 8));	// within <d> of .stl
		}

		duct
		{
			mode	distance;
			levels	((0.02 7)); 
		}

		rotatingCylinder
		{
			mode	inside;
			levels ((1e15 7));
		}
	}


// MESH SELECTION
// After refinement patches get added for all refinementSurfaces and
// all cells intersecting the surfaces get put into these patches. The
// section reachable from the locationInMesh is kept.
	// NOTE: This point should never be on a face, always inside a cell,
	// even after refinement.
	locationInMesh			(3.9999 3.9999 3.9999);

// Whether any faceZones (as specified in the refinementSurfaces)
// are only on the boundary of corresponding cellZones or also allow
// free-standing zone faces. Not used if there are no faceZones.
	allowFreeStandingZoneFaces	true;
}





// SETTINGS FOR THE SNAPPING PROCESS
snapControls
{
// Number of patch smoothing iterations before finding correspondence to surface
	nSmoothPatch		3;

// Relative distance for points to be attracted by surface feature point
// or edge. True distance is this factor times local
// maximum edge length.
	// Tolerance*cellsize on surface = distance to look for closest cell
	// distance relative to local cell size to look for adjacent cells
	tolerance		4.0;

// Number of mesh displacement relaxation iterations.
// Deforms entire mesh to make good quality iterations
	nSolveIter		30;

// Maximum number of snapping relaxation iterations. Should stop before upon
// reaching a correct mesh. Relaxing sharp intersection makes smooth.
	nRelaxIter		5;

// FEATURE SNAPPING
// Number of feature edge snapping iterations. Leave out altogether to disable.
	// TRY WITH THIS DISABLED?!
	// Do not use here since mesh resolution too low and baffles present
        nFeatureSnapIter	10;

        // Detect (geometric only) features by sampling surface (default=false)
        implicitFeatureSnap	true;

        // Use castellatedMeshControls::features (default = true)
        explicitFeatureSnap	true;

        // Detect points on multiple surfaces (only for explicitFeatureSnap)
        multiRegionFeatureSnap	true;
}

// SETTINGS FOR THE LAYER ADDITION PROCESS
addLayersControls
{
// Are the thickness parameters below relative to the undistorted
// size of the refined cell outside layer (true) or absolute sizes (false).
// Relative to local cell size at surface

	relativeSizes			true;

	layers
	{
		fan_blades
		{
			nSurfaceLayers	5;
		}

		duct
		{
 			nSurfaceLayers	5;
		}
	}

// Expansion factor for layer mesh
	expansionRatio			1.2;

// Wanted thickness of final added cell layer. If multiple layers
// is the thickness of the layer furthest away from the wall.
// Relative to undistorted size of cell outside layer.
// See relativeSizes parameter.
// If "relativeSizes" parameter is false you need to specify in meters the 
// finallayerthickness, with it on we just use a factor of the local cell size
// Per final patch (so not geometry!) the layer information
	finalLayerThickness		0.2;

// Minimum thickness of cell layer. If for any reason layer cannot be above
// minThickness do not add layer. Relative to undistorted size of cell outside 
// layer. Not enforced but seeks to meet this. 
	minThickness			0.01;

// If points get not extruded do nGrow layers of connected faces that are
// also not grown. This helps convergence of the layer addition process
// close to features.
// Note: changed(corrected) w.r.t 17x! (didn't do anything in 17x)
	nGrow				0;	
// Sometimes cells don't extrude away from surface due to cells meeting at a
// sharp angle this forces the cells to grow a layer away from the surface.


// ADVANCED SETTINGS
// When not to extrude surface. 0 is flat surface, 90 is when two faces
// are perpendicular
	featureAngle			180;

// At non-patched sides allow mesh to slip if extrusion direction makes
// angle larger than slipFeatureAngle.
	slipFeatureAngle		359;

// Maximum number of snapping relaxation iterations. Should stop
// before upon reaching a correct mesh.
	nRelaxIter			5;

// Number of smoothing iterations of surface normals
	nSmoothSurfaceNormals		3;

// Number of smoothing iterations of interior mesh movement direction
	nSmoothNormals			3;

// Smooth layer thickness over surface patches
	nSmoothThickness		10;

// Stop layer growth on highly warped cells
	maxFaceThicknessRatio		0.3;

// Reduce layer growth where ratio thickness to medial distance is large
	maxThicknessToMedialRatio	0.5;

// Angle used to pick up medial axis points
// NOTE: changed(corrected) w.r.t 17x! 90 degrees corresponds to 130 in 17x.
	minMedianAxisAngle		80;

// Create buffer region for new layer terminations
	nBufferCellsNoExtrude		0;

// Overall max number of layer addition iterations. The mesher will exit
// if it reaches this number of iterations; possibly with an illegal mesh
	nLayerIter			30;
}





// GENERIC MESH QUALITY SETTINGS
// At any undoable phase these determine where to undo.
meshQualityControls
{
	#include	"meshQualityDict"


// ADVANCED QUALITY CONTROLS
// Number of error distribution iterations
	nSmoothScale	4;
// Amount to scale back displacement at error points
	errorReduction	0.75;
}





// ADVANCED snappyHexMesh PARAMETERS

// Flags for optional output
// 0	Only write final meshes
// 1	Write intermediate meshes
// 2	Write volScalarField with cellLevel for postprocessing
// 4	Write current intersections as .obj files
debug	0;

writeFlags
(
	scalarLevels
	layerSets
	layerFields	// Write volScalarField for layer coverage
);



// MERGE TOLERANCE 
// Is a fraction of overall bounding box of initial mesh.
// NOTE: the write tolerance needs to be higher than this.
mergeTolerance	1.001e-6;
//

Here is my fvSchemes script:

ddtSchemes
{
	default		Euler;
}

gradSchemes
{
	default		Gauss linear;
	grad(p)		Gauss linear;
	grad(U)		cellLimited Gauss linear 1;
}

divSchemes
{
	default				none;
	div(phi,U)      		Gauss linearUpwind grad(U);
	div(phi,k)      		Gauss upwind;
	div(phi,omega)  		Gauss upwind;
	div((nuEff*dev2(T(grad(U)))))	Gauss linear;
}

laplacianSchemes
{
	default		Gauss linear limited corrected 0.33; // ????
}

interpolationSchemes
{
	default		linear;
}

snGradSchemes
{
	default		limited corrected 0.33; // ????
}

wallDist
{
	method		meshWave;
}

Here is my fvSolution script:

solvers
{
	pcorr
	{
		solver			GAMG;
		tolerance		1e-03;
		relTol			0.0;
		smoother		DICGaussSeidel;
		nPreSweeps		0;
		nPostSweeps		2;
		cacheAgglomeration	false;		// ????
		nCellsInCoarsestLevel	10;
		agglomerator		faceAreaPair;
		mergeLevels		1;
	}
	
	p
	{
		$pcorr;
		tolerance		1e-06;
		relTol			0.1;
	}

	pFinal
	{
		$p;
		tolerance		1e-08;
		relTol			0.1;
	}

	U
	{
		solver			smoothSolver;
		smoother		symGaussSeidel;
		nSweeps			2;
		tolerance		1e-08;
		relTol			0.1;
	}

	UFinal
	{
		$U;
	}

	k
	{
		solver			smoothSolver;
		smoother		symGaussSeidel;
		nSweeps			2;
		tolerance		1e-08;
		relTol			0.1;
	}

	kFinal
	{
		$k;
	}

	omega
	{
		solver			smoothSolver;
		smoother		symGaussSeidel;
		nSweeps			2;
		tolerance		1e-08;
		relTol			0.1;
	}
	omegaFinal
	{
		$omega;
	}
}

PIMPLE
{
	correctPhi			no;
	nOuterCorrectors		2;
	nCorrectors			1;
	nNonOrthogonalCorrectors	1;
}

relaxationFactors
{
	fields
	{
		p	0.5;
	}
	equations
	{
		U	0.5;
		k	0.5;
		omega	0.5;
	}
}

cache
{
	grad(U);
}

Here is my controlDict script:

application		pimpleDyMFoam;

startFrom		startTime;

startTime		0;

stopAt			endTime;

endTime			0.1;

deltaT			1e-5;

writeControl		adjustableRunTime;
adjustTimeStep		true;
writeInterval		0.001;

purgeWrite		0;

////- For testing with moveDynamicMesh
//deltaT          0.01;
//writeControl    timeStep;
//writeInterval   1;

writeFormat		ascii;

writePrecision		7;

writeCompression	compressed;

timeFormat		general;

timePrecision		6;

runTimeModifiable	true;

maxCo			2;

functions
{
	#include	"readFields"
	#include	"Q"
//	#include	"surfaces"
	#include	"forces"
}

Now it seems that the issue is arrising for the Pressure script within the "0" directory and so I have also included my "p" file and "omega" since I am using the k-omega SST turbulence model:

Initial Conditions (referenced in the "p" and "omega" scripts (yes I know that my initial omega was set to 500 -- I was getting immediate divergence with it set to 0.1 and read that setting it much higher would then fix itself as it iterated):

turbulentKE          0.1;	// kInlet Value
turbulentOmega       500.0;	//

"p" file:

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
	top			// Face 5
	{
		type 	zeroGradient;
	}
	bottom			// Face 4
	{
		type 	fixedValue;
		value	uniform 0;	
	}
	fan_blades
	{
		type	zeroGradient;
	}
	duct
	{
		type	zeroGradient;
	}

	#include "caseDicts/setConstraintTypes"
	#include "include/frontBackLeftRightPatches"
}

And "omega" file:

#include        "include/initialConditions"

dimensions      [0 0 -1 0 0 0 0];

internalField   uniform $turbulentOmega;

boundaryField
{
    //- Set patchGroups for constraint patches
//    #includeEtc "caseDicts/setConstraintTypes"

	top
	{
   		type  	fixedValue;
    		value 	uniform $turbulentOmega;
	}


    bottom
    {
        type            inletOutlet;
        inletValue      uniform $turbulentOmega;
        value           uniform $turbulentOmega;
    }

    fan_blades
    {
        type            omegaWallFunction;
        value           $internalField;
    }
    duct
    {
        type            omegaWallFunction;
        value           uniform $turbulentOmega;
    }
	AMI1
	{
		type	cyclicAMI;
	}
	AMI2
	{
		type	cyclicAMI;
	}

    #include "include/frontBackLeftRightPatches"
}

Any help at all would be welcome!!

Thanks

[Vyssion]

I looked into the code more closely and found the issue.

In my snappyHexMesh script, I forgot to define the faceType of the rotatingCylinder at all. Once I set that to "boundary", the rest of it executed fine.

However...

I am now at the point of attempting to "moveDynamicMesh -checkAMI" and am finding that I am getting aspect ratios on the order of >1e+94, in addition to a lot of "bad" (collectively") cells including negative cell volumes. I have heard that this issue arises due to the way the meshing traverses the AMI boundary and that one possible way to solve this issue is to mesh inside the AMI first, then mesh the rest of the domain and use a "combineMeshes" function (which makes sense seeing as the meshes should end perfectly at the boundaries). This may also explain why I was able to get about 50 iterations out of the pimpleDyMFoam before things blew up. I will try meshing separately and combining and let you know how it goes.

One final thing though; before I defined the faceType of the AMI as a boundary, I was getting fatal errors stating that "cannot find patchField entry for <name>" for p, U, etc. Once I defined it as a boundary, these errors disappeared... mostly. I now receive this error:

Create time

Create mesh for time = 0

Selecting dynamicFvMesh solidBodyMotionFvMesh
Selecting solid-body motion function rotatingMotion
Applying solid body motion to cellZone rotatingCylinder

PIMPLE: no residual control data found. Calculations will employ 2 corrector loops

Reading field p

Reading field U

Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting turbulence model type RAS
Selecting RAS turbulence model kOmegaSST
Selecting patchDistMethod meshWave


--> FOAM FATAL IO ERROR: 
Cannot find patchField entry for AMI1

file: /home/osboxes/OpenFOAM/osboxes-3.0.1/run/base_model_bestMesh_AMI/0/omega.boundaryField from line 30 to line 27.

    From function GeometricField<Type, PatchField, GeoMesh>::GeometricBoundaryField::readField(const DimensionedField<Type, GeoMesh>&, const dictionary&)
    in file /home/openfoam/OpenFOAM/OpenFOAM-3.0.1/src/OpenFOAM/lnInclude/GeometricBoundaryField.C at line 209.

FOAM exiting

Any suggestions as to why this is? I will keep poking around.

[Vyssion]

Figured I would keep updating this so that anyone in the future who has the same issue can have some ideas to try.

I have found that when I performed the "moveDynamicMesh -checkAMI" operation, I was getting meshing errors build up after about 4 time steps (of 1e-05).

Poking around, it seems to be the way that snappyHexMesh works that is the issue. Because I defined the AMI as a boundary, it will mesh through it without any consideration that it will be rotating. Hence, when it does eventually rotate, the cells are dragged along with it.

Playing around with some of the settings in snappyHexMesh, (such as keeping the min-max cell sizes at the interface the same, along with any refinements which overlap the AMI zone being limited in the same way), I have managed to get it to get it to run until about 9 time steps before the -checkAMI blows up and I can simulate through this until 0.00018 seconds.

Considering I need about 0.032 seconds for a full revolution, I am a ways off!!

I need to look into meshing separately and then merging them... Problem is, I have no clue of the nitty-gritty which I would need to do in order to set this up... Has anyone ever done this before or knows of a tutorial that would help??

P.S. I am wondering if part of the issue is that my fan blades touch the duct itself (Im not interested in capturing tip clearance shear interactions at this point) and so right at the AMI boundary, there are zero cells at some points... Could this be an issue??

[EricLEE]

I ran the similar case(Re=1.4*10^5) with the standard k-epsilon turbulence model and I got a higher St=0.25~0.26. I'm confused about y+. I think that the min y+ should be higher than 30 when using wallfunction to make the first point near the cylinder wall locate in the logarithmic layer.

And I realize that y+ is different with that when calculated in initial situation. Which y+ should be based on to decide the distance of the first point to the wall.

Can you give me some opinions? Thanks a lot!