src/engine/engineTopoChangerMesh/thoboisSliding/thoboisSlidingMove.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation, either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Class
  25     verticalValvesGambit
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "thoboisSliding.H"
  30 #include "slidingInterface.H"
  31 #include "layerAdditionRemoval.H"
  32 #include "attachDetach.H"
  33 #include "surfaceFields.H"
  34 #include "regionSplit.H"
  35 #include "componentMixedTetPolyPatchVectorField.H"
  36 #include "mapPolyMesh.H"
  37
  38 #include "tetPolyMesh.H"
  39 #include "tetPointFields.H"
  40 #include "elementFields.H"
  41 #include "fixedValueTetPolyPatchFields.H"
  42 #include "slipTetPolyPatchFields.H"
  43
  44 #include "tetFem.H"
  45
  46 #include "symmetryFvPatch.H"
  47 #include "wedgeFvPatch.H"
  48 #include "emptyFvPatch.H"
  49 #include "zeroGradientTetPolyPatchFields.H"
  50 #include "tetMotionSolver.H"
  51
  52 #include "fixedValueTetPolyPatchFields.H"
  53 #include "mixedTetPolyPatchFields.H"
  54 #include "slipTetPolyPatchFields.H"
  55 #include "zeroGradientTetPolyPatchFields.H"
  56
  57 #include "zeroGradientFvPatchFields.H"
  58 #include "fvCFD.H"
  59
  60 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  61
  62 void Foam::thoboisSliding::makeLayersLive()
  63 {
  64     // Enable layering
  65     forAll (topoChanger_, modI)
  66     {
  67         if (isA<layerAdditionRemoval>(topoChanger_[modI]))
  68         {
  69             topoChanger_[modI].enable();
  70         }
  71         else if (isA<slidingInterface>(topoChanger_[modI]))
  72         {
  73             topoChanger_[modI].disable();
  74         }
  75         else if (isA<attachDetach>(topoChanger_[modI]))
  76         {
  77             topoChanger_[modI].enable();
  78         }
  79         else
  80         {
  81             FatalErrorIn("void Foam::engineTopoFvMesh::makeLayersLive()")
  82                 << "Don't know what to do with mesh modifier "
  83                 << modI << " of type " << topoChanger_[modI].type()
  84                 << abort(FatalError);
  85         }
  86     }
  87 }
  88
  89 void Foam::thoboisSliding::makeSlidersLive()
  90 {
  91     // Enable sliding interface
  92     forAll (topoChanger_, modI)
  93     {
  94         if (isA<layerAdditionRemoval>(topoChanger_[modI]))
  95         {
  96             topoChanger_[modI].disable();
  97         }
  98         else if (isA<slidingInterface>(topoChanger_[modI]))
  99         {
 100             topoChanger_[modI].enable();
 101         }
 102         else if (isA<attachDetach>(topoChanger_[modI]))
 103         {
 104             topoChanger_[modI].enable();
 105         }
 106         else
 107         {
 108             FatalErrorIn("void Foam::engineTopoFvMesh::makeLayersLive()")
 109                 << "Don't know what to do with mesh modifier "
 110                 << modI << " of type " << topoChanger_[modI].type()
 111                 << abort(FatalError);
 112         }
 113     }
 114
 115 }
 116
 117 bool Foam::thoboisSliding::attached() const
 118 {
 119     const polyTopoChanger& morphs = topoChanger_;
 120
 121     bool result = false;
 122
 123     forAll (morphs, modI)
 124     {
 125         if (typeid(morphs[modI]) == typeid(slidingInterface))
 126         {
 127             result =
 128                 result
 129              || refCast<const slidingInterface>(morphs[modI]).attached();
 130         }
 131     }
 132
 133     // Check thal all sliders are in sync (debug only)
 134     forAll (morphs, modI)
 135     {
 136         if (typeid(morphs[modI]) == typeid(slidingInterface))
 137         {
 138             if
 139             (
 140                 result
 141              != refCast<const slidingInterface>(morphs[modI]).attached()
 142             )
 143             {
 144                 FatalErrorIn("bool movingSquaresTM::attached() const")
 145                     << "Slider " << modI << " named " << morphs[modI].name()
 146                     << " out of sync: Should be" << result
 147                     << abort(FatalError);
 148             }
 149         }
 150     }
 151
 152     return result;
 153 }
 154
 155
 156 void Foam::thoboisSliding::valveDetach()
 157 {
 158     // Enable sliding interface
 159     forAll (topoChanger_, modI)
 160     {
 161         if (isA<attachDetach>(topoChanger_[modI]))
 162         {
 163             const attachDetach& ad =
 164                 refCast<const attachDetach>(topoChanger_[modI]);
 165
 166             const word masterName = ad.masterPatchID().name();
 167
 168             // Find the valve with that name
 169             label valveIndex = -1;
 170
 171             forAll (valves_, valveI)
 172             {
 173                 if
 174                 (
 175                     valves_[valveI].detachInCylinderPatchID().name()
 176                  == masterName
 177                 )
 178                 {
 179                     valveIndex = valveI;
 180                     break;
 181                 }
 182             }
 183
 184             if (valveIndex < 0)
 185             {
 186                 FatalErrorIn
 187                 (
 188                     "void Foam::engineTopoFvMesh::prepareValveDetach()"
 189                 )   << "Cannot match patch for attach/detach " << modI
 190                     << abort(FatalError);
 191             }
 192
 193             if (debug)
 194             {
 195                 Info<< " valveI: " << valveIndex << " attached: "
 196                     << ad.attached()
 197                     << " valve open: " << valves_[valveIndex].isOpen()
 198                     << endl;
 199             }
 200
 201             ad.setDetach();
 202         }
 203     }
 204 }
 205
 206 void Foam::thoboisSliding::valveAttach()
 207 {
 208     // Enable sliding interface
 209     forAll (topoChanger_, modI)
 210     {
 211         if (isA<attachDetach>(topoChanger_[modI]))
 212         {
 213             const attachDetach& ad =
 214                 refCast<const attachDetach>(topoChanger_[modI]);
 215
 216             const word masterName = ad.masterPatchID().name();
 217
 218             // Find the valve with that name
 219             label valveIndex = -1;
 220
 221             forAll (valves_, valveI)
 222             {
 223                 if
 224                 (
 225                     valves_[valveI].detachInCylinderPatchID().name()
 226                  == masterName
 227                 )
 228                 {
 229                     valveIndex = valveI;
 230                     break;
 231                 }
 232             }
 233
 234             if (valveIndex < 0)
 235             {
 236                 FatalErrorIn
 237                 (
 238                     "void Foam::engineTopoFvMesh::prepareValveDetach()"
 239                 )   << "Cannot match patch for attach/detach " << modI
 240                     << abort(FatalError);
 241             }
 242
 243             if (debug)
 244             {
 245                 Info<< " valveI: " << valveIndex << " attached: "
 246                     << ad.attached()
 247                     << " valve open: " << valves_[valveIndex].isOpen()
 248                     << endl;
 249             }
 250
 251             ad.setAttach();
 252
 253         }
 254     }
 255 }
 256
 257 void Foam::thoboisSliding::prepareValveDetach()
 258 {
 259     // Enable sliding interface
 260     forAll (topoChanger_, modI)
 261     {
 262         if (isA<attachDetach>(topoChanger_[modI]))
 263         {
 264             const attachDetach& ad =
 265                 refCast<const attachDetach>(topoChanger_[modI]);
 266
 267             const word masterName = ad.masterPatchID().name();
 268
 269             // Find the valve with that name
 270             label valveIndex = -1;
 271
 272             forAll (valves_, valveI)
 273             {
 274                 if
 275                 (
 276                     valves_[valveI].detachInCylinderPatchID().name()
 277                  == masterName
 278                 )
 279                 {
 280                     valveIndex = valveI;
 281                     break;
 282                 }
 283             }
 284
 285             if (valveIndex < 0)
 286             {
 287                 FatalErrorIn
 288                 (
 289                     "void Foam::engineTopoFvMesh::prepareValveDetach()"
 290                 )   << "Cannot match patch for attach/detach " << modI
 291                     << abort(FatalError);
 292             }
 293
 294             if (debug)
 295             {
 296                 Info<< " valveI: " << valveIndex << " attached: "
 297                     << ad.attached()
 298                     << " valve open: " << valves_[valveIndex].isOpen()
 299                     << endl;
 300             }
 301
 302             if (valves_[valveIndex].isOpen())
 303             {
 304                 ad.setAttach();
 305             }
 306             else
 307             {
 308                 ad.setDetach();
 309             }
 310         }
 311     }
 312 }
 313
 314
 315 bool Foam::thoboisSliding::update()
 316 {
 317     tetMotionSolver& mSolver =
 318         refCast<tetMotionSolver>(msPtr_());
 319
 320     // Detaching the interfacethobois2DSlidingDeform
 321     if (attached())
 322     {
 323         Info << "Decoupling sliding interfaces" << endl;
 324         makeSlidersLive();
 325         valveDetach();
 326         autoPtr<mapPolyMesh> topoChangeMap1 = topoChanger_.changeMesh();
 327
 328         Info << "sliding interfaces successfully decoupled!!!" << endl;
 329
 330         if (topoChangeMap1->morphing())
 331         {
 332             mSolver.updateMesh(topoChangeMap1());
 333         }
 334     }
 335     else
 336     {
 337         Info << "Sliding interfaces decoupled" << endl;
 338         valveDetach();
 339     }
 340
 341     Info << "Executing layer action" << endl;
 342
 343     scalar deltaZ = engTime().pistonDisplacement().value();
 344
 345     // Piston Layering
 346     makeLayersLive();
 347
 348
 349     {
 350         // Activate bottom layer
 351         Info << "Valve layering mode" << endl;
 352         forAll(valves_, valveI)
 353         {
 354             if(valves_[valveI].bottomPatchID().active())
 355             {
 356                 // Find piston mesh modifier
 357                 const label valveLayerID2 =
 358                     topoChanger_.findModifierID
 359                     (
 360                         "valveBottomLayer"
 361                       + Foam::name(valveI + 1)
 362                     );
 363
 364                 topoChanger_[valveLayerID2].enable();
 365             }
 366         }
 367     }
 368
 369     virtualPistonPosition() += deltaZ;
 370
 371     {
 372         forAll(valves_, valveI)
 373         {
 374             if(!realDeformation())
 375             {
 376                 scalar valveDisplacement =
 377                     valves_[valveI].curVelocity()*
 378                     valves_[valveI].cs().axis().z()*engTime().deltaT().value();
 379
 380                 Info<< "valveDisplacement = " << valveDisplacement << nl
 381                     << "valve velocity =" << valves_[valveI].curVelocity()
 382                     << endl;
 383
 384                 valveTopPosition_[valveI] += valveDisplacement;
 385                 valveBottomPosition_[valveI] += valveDisplacement;
 386                 valvePistonPosition_[valveI] += deltaZ;
 387             }
 388         }
 389     }
 390
 391     forAll(valves_,valveI)
 392     {
 393         if(valves_[valveI].curLift() > valves_[valveI].deformationLift())
 394         {
 395             if(valves_[valveI].poppetPatchID().active())
 396             {
 397                 // Find valve layer mesh modifier
 398                 const label valveLayerID =
 399                     topoChanger_.findModifierID
 400                     (
 401                         "valvePoppetLayer"
 402                       + Foam::name(valveI + 1)
 403                     );
 404
 405                 if (valveLayerID < 0)
 406                 {
 407                         FatalErrorIn("void engineFvMesh::moveAndMorph()")
 408                            << "valve modifier not found."
 409                            << abort(FatalError);
 410                 }
 411
 412                 topoChanger_[valveLayerID].enable();
 413             }
 414         }
 415         else
 416         {
 417             if(valves_[valveI].poppetPatchID().active())
 418             {
 419                 // Find valve layer mesh modifier
 420                 const label valveLayerID =
 421                     topoChanger_.findModifierID
 422                     (
 423                         "valvePoppetLayer"
 424                       + Foam::name(valveI + 1)
 425                     );
 426
 427                 if (valveLayerID < 0)
 428                 {
 429                         FatalErrorIn("void engineFvMesh::moveAndMorph()")
 430                            << "valve modifier not found."
 431                            << abort(FatalError);
 432                 }
 433
 434                 topoChanger_[valveLayerID].disable();
 435             }
 436         }
 437     }
 438
 439     pointField newPoints = points();
 440
 441 //    pointField newPoints = allPoints();
 442
 443     autoPtr<mapPolyMesh> topoChangeMap2 = topoChanger_.changeMesh();
 444
 445     // Changing topology by hand
 446     if (topoChangeMap2->morphing())
 447     {
 448         mSolver.updateMesh(topoChangeMap2());
 449
 450         if (topoChangeMap2->hasMotionPoints())
 451         {
 452             movePoints(topoChangeMap2->preMotionPoints());
 453             newPoints = topoChangeMap2->preMotionPoints();
 454         }
 455         setV0();
 456         resetMotion();
 457     }
 458
 459     // Work array for new points position.
 460
 461 // NEW !
 462
 463     // Reset the position of layered interfaces
 464
 465     {
 466
 467 #       include "setValveMotionBoundaryConditionThoboisSliding.H"
 468
 469         dynamicLabelList constrainedPoints(mSolver.curPoints()().size()/100);
 470         DynamicList<vector> constrainedVelocity
 471         (
 472             mSolver.curPoints()().size()/100
 473         );
 474 #       include "setThoboisSlidingConstraint.H"
 475
 476         forAll (constrainedPoints, i)
 477         {
 478             mSolver.setConstraint
 479             (
 480                 constrainedPoints[i],
 481                 constrainedVelocity[i]
 482             );
 483         }
 484
 485         mSolver.solve();
 486
 487         newPoints = mSolver.curPoints();
 488         movePoints(newPoints);
 489
 490         setVirtualPositions();
 491         mSolver.clearConstraints();
 492
 493 //      layering
 494
 495 #       include "moveValvePointsThoboisSliding.H"
 496         movePoints(newPoints);
 497 //        newPoints = points();
 498         newPoints = allPoints();
 499         setVirtualPositions();
 500
 501     }
 502
 503     forAll(valves(), valveI)
 504     {
 505         Info << "Valve Top Position for valve n. " << valveI + 1 << " = " <<
 506         valveTopPosition_[valveI] << endl;
 507         Info << "Valve Bottom Position for valve n. " << valveI + 1 << " = " <<
 508         valveBottomPosition_[valveI] << endl;
 509     }
 510
 511     Info << "virtualPistonPosition = " << virtualPistonPosition()
 512     << ", deckHeight = " << deckHeight()
 513     << ", pistonPosition = " << pistonPosition()
 514     << endl;
 515
 516     pistonPosition() += deltaZ;
 517
 518     scalarField Vold = V();
 519
 520     volScalarField vMesh
 521     (
 522         IOobject
 523         (
 524             "vMesh",
 525             time().timeName(),
 526             *this,
 527             IOobject::NO_READ,
 528             IOobject::NO_WRITE
 529         ),
 530         *this,
 531         dimensionedScalar("zero", dimensionSet(0, 0, 0, 0, 0), 0.0),
 532         zeroGradientFvPatchScalarField::typeName
 533     );
 534
 535     vMesh.internalField() = V();
 536     vMesh.oldTime().internalField() = V0();
 537
 538
 539     // Coupling the interface Again
 540
 541     {
 542
 543         pointField newPointsNew = allPoints();
 544
 545         // Attach the interface
 546         Info << "Coupling sliding interfaces" << endl;
 547         makeSlidersLive();
 548         valveAttach();
 549         prepareValveDetach();
 550
 551         pointField oldPointsNew = oldAllPoints();
 552
 553         // Changing topology by hand
 554         autoPtr<mapPolyMesh> topoChangeMap3 = topoChanger_.changeMesh();
 555
 556         Info << "Sliding interfaces coupled: " << attached() << endl;
 557         if (topoChangeMap3->morphing())
 558         {
 559             Info << "mesh changed 3" << endl;
 560
 561             mSolver.updateMesh(topoChangeMap3());
 562
 563             if (topoChangeMap3->hasMotionPoints())
 564             {
 565                 movePoints(topoChangeMap3->preMotionPoints());
 566                 resetMotion();
 567                 setV0();
 568             }
 569
 570             {
 571                 // correct the motion after attaching the sliding interface
 572                 Info<< "oldPointsNew.size = " << oldPointsNew.size() << nl
 573                     << "allPointsNew.size = " << allPoints().size() << nl
 574                     << "pointsNew.size = " << points().size() << endl;
 575
 576                 pointField mappedOldPointsNew(allPoints().size());
 577
 578                 pointField newPoints = allPoints();
 579 //                pointField newPoints = points();
 580
 581                 mappedOldPointsNew.map
 582                 (
 583                     oldPointsNew,
 584                     topoChangeMap3->pointMap()
 585                 );
 586
 587                 movePoints(mappedOldPointsNew);
 588
 589                 resetMotion();
 590                 setV0();
 591
 592                 movePoints(mappedOldPointsNew);
 593                 Info<< "max mesh phi, 3 = " << max(phi()) << nl
 594                     << "min mesh phi, 3 = " << min(phi()) << nl
 595                     << "mappedOldPointsNew.size() = "
 596                     << mappedOldPointsNew.size() << nl
 597                     << "newPoints.size() = " << newPoints.size() << nl
 598                     << "topoChangeMap3->pointMap().size() = "
 599                     << topoChangeMap3->pointMap().size() << endl;
 600
 601                 movePoints(newPoints);
 602             }
 603         }
 604     }
 605
 606     return true;
 607 }
 608
 609