View Issue Details
ID | Project | Category | View Status | Date Submitted | Last Update |
---|---|---|---|---|---|
0001425 | OpenFOAM | Bug | public | 2014-10-27 22:47 | 2014-12-12 16:27 |
Reporter | Assigned To | henry | |||
Priority | normal | Severity | major | Reproducibility | always |
Status | closed | Resolution | no change required | ||
Platform | Supermicro | OS | Centos | OS Version | 6.5 |
Summary | 0001425: faceToPointInterpolate does not interpolate across a cyclic patch (and perhaps also a processor patch) | ||||
Description | faceToPointInterpolate uses the nearest face value for points at the boundaries of a region. This is incorrect for cyclic boundaries; the points from the opposite face should be used as well. I have tested a number of scalar and vector fields with the same result. | ||||
Steps To Reproduce | Take a sinusoidal boundary surface, say with normal in the x direction (on average) and cyclic boundaries in the other direction (y). The pointNormals to the sinusoidal boundary on both sides of the cyclic patch should point in the same direction, but since they are taken from the interior face values they are at a slight angle. If the values from both sides of the cyclic patch were included in the interpolation, the correct result would be obtained. I think this is an error because the solution is not strictly periodic. It causes serious problems when the (sinusoidal) boundary is moving which we are currently patching within our solver. But I think this is a generic problem with the interpolation function and a general solution is needed. | ||||
Additional Information | This also seems to affect the processor boundary patch in the same way, presumably for the same reason. Results in a parallel simulation do not agree exactly with results obtained in serial near a cyclic patch. | ||||
Tags | No tags attached. | ||||
|
faceToPointInterpolate operates on a PrimitivePatch level - it does not necessary relate to a (decomposed) mesh - it could e.g. operate on a triSurface. The problem you are seeing is from the normal being calculated using the local patch only. As a work around maybe you can specify a fixed normal (e.g. fixedNormalSlip). |
|
Its not just the normals but fields as well. We use gradC (C is a scalar field) interpolated to points to determine the motion of the (initially) sinusoidal boundary. We are working around this in the cyclic case - just doing the additional interpolation by averaging. Harder with the parallel case. Did you see the previous bug report (1424)? This may be related. In general it seems to me that interpolation for any coupled patch should involve both regions (at least as an option if not by default), even if there is a jump condition at the boundary. Are we missing a function that interpolates correctly at coupled patches? |
|
Get yourself a pointVectorField with the correct bcs, e.g. pf, and const pointConstraints& pcs = pointConstraints::New(pf.mesh()); pcs.constrainDisplacement(pf); |
2014-11-11 02:23
|
moveMesh.H (2,057 bytes)
// Include code for dynamic mesh motion { // Create point vector fields pointVectorField& intPointDisp = const_cast<pointVectorField&> ( mesh.objectRegistry::lookupObject<pointVectorField> ( "pointDisplacement" ) ); // Get patch ID for int2ext label intPatch = mesh.boundaryMesh().findPatchID("int2ext"); // Correct point normals on cyclic boundaries (1,0,0) vectorField pointN = mesh.boundaryMesh()[intPatch].pointNormals(); vector N(1.0, 0.0, 0.0); for(std::map<int, std::pair<int,int> >::iterator iter = cList.begin(); iter != cList.end(); ++iter){ pointN[iter->second.first] = N; pointN[iter->second.second] = N; } // Calculate surface motion - interpolate faceGradC to points scalarField faceGradC = gradC.boundaryField()[intPatch]; primitivePatchInterpolation patchInterpolator(mesh.boundaryMesh()[intPatch]); scalarField pointGradC = patchInterpolator.faceToPointInterpolate(faceGradC); vectorField pointU = -pointGradC*pointN/Pe.value(); // Scale by Pe vectorField intDisp = pointU*runTime.deltaTValue(); // Correct points for interpolation error across cyclic boundary for(std::map<int, std::pair<int,int> >::iterator iter = cList.begin(); iter != cList.end(); ++iter){ vector avgX(0.5*(intDisp[iter->second.first].x() + intDisp[iter->second.second].x()), 0,0); intDisp[iter->second.first] = avgX; intDisp[iter->second.second] = avgX; } vectorField &intPatchDisp = refCast<vectorField>(intPointDisp.boundaryField()[intPatch]); intPointDisp.boundaryField()[intPatch] == intPatchDisp + intDisp; // Mesh relaxation for(int cycle=1; cycle <= meshCycles; cycle++){ Info << nl << "Mesh update cycle: " << cycle << endl; mesh.update(); } Info<< "Mesh update: ExecutionTime = " << runTime.elapsedCpuTime() << " s" << " ClockTime = " << runTime.elapsedClockTime() << nl << endl; } |
|
Mattijs I was not sure how to implement your suggestion since the fields seem to be of different type (vectorField vs pointVectorField). Attached is my code to move a boundary in response to a concentration gradient (dissolution). As you can see there are two "fixes": one for the normals and one for the displacements. Both are needed because faceToPointInterpolate does not handle coupled boundary conditions correctly. I recently noticed the same problem with the symmetryPlane boundary condition. The point normals have a component normal to the symmetryPlane just as they do with cyclic. Do you see this as a bug? Missing feature? User error? It would be helpful to know. Thanks. |
|
The displacement of a point is determined by all the patches it is on. One of the neighbouring patches might be a constraint patch (e.g. symmetryPlane). Thus you cannot just use a single patch to determine the motion - you have to take all the constraints into account. This (and the parallel synchronisation) is handled by pointConstraints. |
|
mattijs I think you can close this out. We agree this is not a bug. Two points were confusing to us. 1) pointNormal being just a vector field does not know about boundary conditions 2) more confusing was that the pointMotionU field did not correct for velocities normal to a cyclic patch. We later realized that this was because we specified a fixed value on an adjacent patch. OpenFOAM does not correct the edge points for the cyclic boundary condition but allows the fixedValue patch to be the decider. This seems to us to be backwards but maybe there is a reason. In any case we can accept that if we modify a fixed value patch we should make sure it satisfies the boundary conditions on the edges with adjacent cyclic patches. pointConstraint does not fix this, but it seems that perhaps the syncTools class is what we need. Thanks for your help. Sorry for the mistaken report. |
Date Modified | Username | Field | Change |
---|---|---|---|
2014-10-27 22:47 |
|
New Issue | |
2014-10-28 10:57 |
|
Note Added: 0003268 | |
2014-10-28 12:28 |
|
Note Added: 0003271 | |
2014-10-30 09:42 |
|
Note Added: 0003274 | |
2014-11-11 02:23 |
|
File Added: moveMesh.H | |
2014-11-11 02:31 |
|
Note Added: 0003279 | |
2014-11-12 10:08 |
|
Note Added: 0003285 | |
2014-11-26 03:12 |
|
Note Added: 0003297 | |
2014-12-12 16:27 | henry | Status | new => closed |
2014-12-12 16:27 | henry | Assigned To | => henry |
2014-12-12 16:27 | henry | Resolution | open => no change required |