00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #include "Particle.H"
00027 #include <lagrangian/Cloud.H>
00028 #include <OpenFOAM/wedgePolyPatch.H>
00029 #include <OpenFOAM/symmetryPolyPatch.H>
00030 #include <OpenFOAM/cyclicPolyPatch.H>
00031 #include <OpenFOAM/processorPolyPatch.H>
00032 #include <OpenFOAM/wallPolyPatch.H>
00033 #include <OpenFOAM/transform.H>
00034
00035
00036
00037 template<class ParticleType>
00038 void Foam::Particle<ParticleType>::findFaces
00039 (
00040 const vector& position,
00041 DynamicList<label>& faceList
00042 ) const
00043 {
00044 const polyMesh& mesh = cloud_.polyMesh_;
00045 const labelList& faces = mesh.cells()[celli_];
00046 const vector& C = mesh.cellCentres()[celli_];
00047
00048 faceList.clear();
00049 forAll(faces, i)
00050 {
00051 label facei = faces[i];
00052 scalar lam = lambda(C, position, facei);
00053
00054 if ((lam > 0) && (lam < 1.0))
00055 {
00056 faceList.append(facei);
00057 }
00058 }
00059 }
00060
00061
00062 template<class ParticleType>
00063 void Foam::Particle<ParticleType>::findFaces
00064 (
00065 const vector& position,
00066 const label celli,
00067 const scalar stepFraction,
00068 DynamicList<label>& faceList
00069 ) const
00070 {
00071 const polyMesh& mesh = cloud_.pMesh();
00072 const labelList& faces = mesh.cells()[celli];
00073 const vector& C = mesh.cellCentres()[celli];
00074
00075 faceList.clear();
00076 forAll(faces, i)
00077 {
00078 label facei = faces[i];
00079 scalar lam = lambda(C, position, facei, stepFraction);
00080
00081 if ((lam > 0) && (lam < 1.0))
00082 {
00083 faceList.append(facei);
00084 }
00085 }
00086 }
00087
00088
00089 template<class ParticleType>
00090 template<class TrackData>
00091 void Foam::Particle<ParticleType>::prepareForParallelTransfer
00092 (
00093 const label patchi,
00094 TrackData& td
00095 )
00096 {
00097
00098 facei_ = patchFace(patchi, facei_);
00099 }
00100
00101
00102 template<class ParticleType>
00103 template<class TrackData>
00104 void Foam::Particle<ParticleType>::correctAfterParallelTransfer
00105 (
00106 const label patchi,
00107 TrackData& td
00108 )
00109 {
00110 const processorPolyPatch& ppp =
00111 refCast<const processorPolyPatch>
00112 (cloud_.pMesh().boundaryMesh()[patchi]);
00113
00114 celli_ = ppp.faceCells()[facei_];
00115
00116 if (!ppp.parallel())
00117 {
00118 if (ppp.forwardT().size() == 1)
00119 {
00120 const tensor& T = ppp.forwardT()[0];
00121 transformPosition(T);
00122 static_cast<ParticleType&>(*this).transformProperties(T);
00123 }
00124 else
00125 {
00126 const tensor& T = ppp.forwardT()[facei_];
00127 transformPosition(T);
00128 static_cast<ParticleType&>(*this).transformProperties(T);
00129 }
00130 }
00131 else if (ppp.separated())
00132 {
00133 if (ppp.separation().size() == 1)
00134 {
00135 position_ -= ppp.separation()[0];
00136 static_cast<ParticleType&>(*this).transformProperties
00137 (
00138 -ppp.separation()[0]
00139 );
00140 }
00141 else
00142 {
00143 position_ -= ppp.separation()[facei_];
00144 static_cast<ParticleType&>(*this).transformProperties
00145 (
00146 -ppp.separation()[facei_]
00147 );
00148 }
00149 }
00150
00151
00152 if (stepFraction_ > (1.0 - SMALL))
00153 {
00154 stepFraction_ = 1.0;
00155 facei_ = -1;
00156 }
00157 else
00158 {
00159 facei_ += ppp.start();
00160 }
00161 }
00162
00163
00164
00165
00166 template<class ParticleType>
00167 Foam::Particle<ParticleType>::Particle
00168 (
00169 const Cloud<ParticleType>& cloud,
00170 const vector& position,
00171 const label celli
00172 )
00173 :
00174 cloud_(cloud),
00175 position_(position),
00176 celli_(celli),
00177 facei_(-1),
00178 stepFraction_(0.0),
00179 origProc_(Pstream::myProcNo()),
00180 origId_(cloud_.getNewParticleID())
00181 {}
00182
00183
00184 template<class ParticleType>
00185 Foam::Particle<ParticleType>::Particle(const Particle<ParticleType>& p)
00186 :
00187 cloud_(p.cloud_),
00188 position_(p.position_),
00189 celli_(p.celli_),
00190 facei_(p.facei_),
00191 stepFraction_(p.stepFraction_),
00192 origProc_(p.origProc_),
00193 origId_(p.origId_)
00194 {}
00195
00196
00197
00198
00199 template<class ParticleType>
00200 template<class TrackData>
00201 Foam::label Foam::Particle<ParticleType>::track
00202 (
00203 const vector& endPosition,
00204 TrackData& td
00205 )
00206 {
00207 facei_ = -1;
00208
00209
00210 while(!onBoundary() && stepFraction_ < 1.0 - SMALL)
00211 {
00212 stepFraction_ += trackToFace(endPosition, td)*(1.0 - stepFraction_);
00213 }
00214
00215 return facei_;
00216 }
00217
00218
00219
00220 template<class ParticleType>
00221 Foam::label Foam::Particle<ParticleType>::track(const vector& endPosition)
00222 {
00223 int dummyTd;
00224 return track(endPosition, dummyTd);
00225 }
00226
00227 template<class ParticleType>
00228 template<class TrackData>
00229 Foam::scalar Foam::Particle<ParticleType>::trackToFace
00230 (
00231 const vector& endPosition,
00232 TrackData& td
00233 )
00234 {
00235 const polyMesh& mesh = cloud_.polyMesh_;
00236
00237 DynamicList<label>& faces = cloud_.labels_;
00238 findFaces(endPosition, faces);
00239
00240 facei_ = -1;
00241 scalar trackFraction = 0.0;
00242
00243 if (faces.empty())
00244 {
00245 trackFraction = 1.0;
00246 position_ = endPosition;
00247 }
00248 else
00249 {
00250 scalar lambdaMin = GREAT;
00251
00252 if (faces.size() == 1)
00253 {
00254 lambdaMin = lambda(position_, endPosition, faces[0], stepFraction_);
00255 facei_ = faces[0];
00256 }
00257 else
00258 {
00259
00260
00261
00262
00263
00264 forAll(faces, i)
00265 {
00266 scalar lam =
00267 lambda(position_, endPosition, faces[i], stepFraction_);
00268
00269 if (lam < lambdaMin)
00270 {
00271 lambdaMin = lam;
00272 facei_ = faces[i];
00273 }
00274 }
00275 }
00276
00277 bool internalFace = cloud_.internalFace(facei_);
00278
00279
00280
00281
00282
00283
00284 if (lambdaMin > 0.0)
00285 {
00286 if (lambdaMin <= 1.0)
00287 {
00288 trackFraction = lambdaMin;
00289 position_ += trackFraction*(endPosition - position_);
00290 }
00291 else
00292 {
00293 trackFraction = 1.0;
00294 position_ = endPosition;
00295 }
00296 }
00297 else if (static_cast<ParticleType&>(*this).softImpact())
00298 {
00299
00300
00301
00302 trackFraction = 1.0;
00303 position_ = endPosition;
00304 }
00305
00306
00307 if (internalFace)
00308 {
00309 if (celli_ == mesh.faceOwner()[facei_])
00310 {
00311 celli_ = mesh.faceNeighbour()[facei_];
00312 }
00313 else if (celli_ == mesh.faceNeighbour()[facei_])
00314 {
00315 celli_ = mesh.faceOwner()[facei_];
00316 }
00317 else
00318 {
00319 FatalErrorIn
00320 (
00321 "Particle::trackToFace(const vector&, TrackData&)"
00322 )<< "addressing failure" << nl
00323 << abort(FatalError);
00324 }
00325 }
00326 else
00327 {
00328 ParticleType& p = static_cast<ParticleType&>(*this);
00329
00330
00331 if (p.softImpact())
00332 {
00333 trackFraction = 1.0;
00334 position_ = endPosition;
00335 }
00336
00337 label patchi = patch(facei_);
00338 const polyPatch& patch = mesh.boundaryMesh()[patchi];
00339
00340 if (!p.hitPatch(patch, td, patchi))
00341 {
00342 if (isA<wedgePolyPatch>(patch))
00343 {
00344 p.hitWedgePatch
00345 (
00346 static_cast<const wedgePolyPatch&>(patch), td
00347 );
00348 }
00349 else if (isA<symmetryPolyPatch>(patch))
00350 {
00351 p.hitSymmetryPatch
00352 (
00353 static_cast<const symmetryPolyPatch&>(patch), td
00354 );
00355 }
00356 else if (isA<cyclicPolyPatch>(patch))
00357 {
00358 p.hitCyclicPatch
00359 (
00360 static_cast<const cyclicPolyPatch&>(patch), td
00361 );
00362 }
00363 else if (isA<processorPolyPatch>(patch))
00364 {
00365 p.hitProcessorPatch
00366 (
00367 static_cast<const processorPolyPatch&>(patch), td
00368 );
00369 }
00370 else if (isA<wallPolyPatch>(patch))
00371 {
00372 p.hitWallPatch
00373 (
00374 static_cast<const wallPolyPatch&>(patch), td
00375 );
00376 }
00377 else
00378 {
00379 p.hitPatch(patch, td);
00380 }
00381 }
00382 }
00383 }
00384
00385
00386
00387
00388
00389
00390
00391 if (trackFraction < SMALL)
00392 {
00393 position_ += 1.0e-3*(mesh.cellCentres()[celli_] - position_);
00394 }
00395
00396 return trackFraction;
00397 }
00398
00399 template<class ParticleType>
00400 Foam::scalar Foam::Particle<ParticleType>::trackToFace
00401 (
00402 const vector& endPosition
00403 )
00404 {
00405 int dummyTd;
00406 return trackToFace(endPosition, dummyTd);
00407 }
00408
00409 template<class ParticleType>
00410 void Foam::Particle<ParticleType>::transformPosition(const tensor& T)
00411 {
00412 position_ = transform(T, position_);
00413 }
00414
00415
00416 template<class ParticleType>
00417 void Foam::Particle<ParticleType>::transformProperties(const tensor&)
00418 {}
00419
00420
00421 template<class ParticleType>
00422 void Foam::Particle<ParticleType>::transformProperties(const vector&)
00423 {}
00424
00425
00426 template<class ParticleType>
00427 template<class TrackData>
00428 bool Foam::Particle<ParticleType>::hitPatch
00429 (
00430 const polyPatch&,
00431 TrackData&,
00432 const label
00433 )
00434 {
00435 return false;
00436 }
00437
00438
00439 template<class ParticleType>
00440 template<class TrackData>
00441 void Foam::Particle<ParticleType>::hitWedgePatch
00442 (
00443 const wedgePolyPatch& wpp,
00444 TrackData&
00445 )
00446 {
00447 vector nf = wpp.faceAreas()[wpp.whichFace(facei_)];
00448 nf /= mag(nf);
00449
00450 static_cast<ParticleType&>(*this).transformProperties(I - 2.0*nf*nf);
00451 }
00452
00453
00454 template<class ParticleType>
00455 template<class TrackData>
00456 void Foam::Particle<ParticleType>::hitSymmetryPatch
00457 (
00458 const symmetryPolyPatch& spp,
00459 TrackData&
00460 )
00461 {
00462 vector nf = spp.faceAreas()[spp.whichFace(facei_)];
00463 nf /= mag(nf);
00464
00465 static_cast<ParticleType&>(*this).transformProperties(I - 2.0*nf*nf);
00466 }
00467
00468
00469 template<class ParticleType>
00470 template<class TrackData>
00471 void Foam::Particle<ParticleType>::hitCyclicPatch
00472 (
00473 const cyclicPolyPatch& cpp,
00474 TrackData&
00475 )
00476 {
00477 label patchFacei_ = cpp.whichFace(facei_);
00478
00479 facei_ = cpp.transformGlobalFace(facei_);
00480
00481 celli_ = cloud_.polyMesh_.faceOwner()[facei_];
00482
00483 if (!cpp.parallel())
00484 {
00485 const tensor& T = cpp.transformT(patchFacei_);
00486
00487 transformPosition(T);
00488 static_cast<ParticleType&>(*this).transformProperties(T);
00489 }
00490 else if (cpp.separated())
00491 {
00492 position_ += cpp.separation(patchFacei_);
00493 static_cast<ParticleType&>(*this).transformProperties
00494 (
00495 cpp.separation(patchFacei_)
00496 );
00497 }
00498 }
00499
00500
00501 template<class ParticleType>
00502 template<class TrackData>
00503 void Foam::Particle<ParticleType>::hitProcessorPatch
00504 (
00505 const processorPolyPatch& spp,
00506 TrackData& td
00507 )
00508 {}
00509
00510
00511 template<class ParticleType>
00512 template<class TrackData>
00513 void Foam::Particle<ParticleType>::hitWallPatch
00514 (
00515 const wallPolyPatch& spp,
00516 TrackData&
00517 )
00518 {}
00519
00520
00521 template<class ParticleType>
00522 template<class TrackData>
00523 void Foam::Particle<ParticleType>::hitPatch
00524 (
00525 const polyPatch& spp,
00526 TrackData&
00527 )
00528 {}
00529
00530
00531
00532
00533 #include "ParticleIO.C"
00534
00535