InjectionModel.C 20.3 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
andy's avatar
andy committed
5
    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
6
7
8
9
10
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

11
12
13
14
    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
15
16
17
18
19
20
21

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
22
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
23
24
25
26

\*---------------------------------------------------------------------------*/

#include "InjectionModel.H"
27
#include "mathematicalConstants.H"
28
#include "meshTools.H"
29
#include "volFields.H"
30

31
using namespace Foam::constant::mathematical;
Andrew Heather's avatar
Andrew Heather committed
32

33
// * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * * //
34

35
36
37
38
39
40
41
42
43
44
45
46
template<class CloudType>
bool Foam::InjectionModel<CloudType>::validInjection(const label parcelI)
{
    notImplemented
    (
        "bool Foam::InjectionModel<CloudType>::validInjection(const label)"
    );

    return false;
}


47
template<class CloudType>
48
bool Foam::InjectionModel<CloudType>::prepareForNextTimeStep
49
(
50
    const scalar time,
Andrew Heather's avatar
Andrew Heather committed
51
    label& newParcels,
andy's avatar
andy committed
52
    scalar& newVolumeFraction
53
54
55
)
{
    // Initialise values
Andrew Heather's avatar
Andrew Heather committed
56
    newParcels = 0;
andy's avatar
andy committed
57
    newVolumeFraction = 0.0;
58
    bool validInjection = false;
59
60

    // Return if not started injection event
61
    if (time < SOI_)
62
    {
63
        timeStep0_ = time;
64
        return validInjection;
65
66
67
68
    }

    // Make times relative to SOI
    scalar t0 = timeStep0_ - SOI_;
69
    scalar t1 = time - SOI_;
70
71

    // Number of parcels to inject
72
    newParcels = this->parcelsToInject(t0, t1);
73
74

    // Volume of parcels to inject
75
76
77
    newVolumeFraction =
        this->volumeToInject(t0, t1)
       /(volumeTotal_ + ROOTVSMALL);
78

andy's avatar
andy committed
79
    if (newVolumeFraction > 0)
80
    {
81
82
83
84
85
86
87
88
89
90
91
        if (newParcels > 0)
        {
            timeStep0_ = time;
            validInjection = true;
        }
        else
        {
            // injection should have started, but not sufficient volume to
            // produce (at least) 1 parcel - hold value of timeStep0_
            validInjection = false;
        }
92
93
94
    }
    else
    {
95
        timeStep0_ = time;
96
        validInjection = false;
97
    }
98
99

    return validInjection;
100
101
102
}


103
template<class CloudType>
104
bool Foam::InjectionModel<CloudType>::findCellAtPosition
105
106
(
    label& cellI,
107
108
    label& tetFaceI,
    label& tetPtI,
109
110
    vector& position,
    bool errorOnNotFound
111
112
)
{
113
    const volVectorField& cellCentres = this->owner().mesh().C();
114

115
    const vector p0 = position;
116

andy's avatar
andy committed
117
    this->owner().mesh().findCellFacePt
118
119
120
121
122
123
    (
        position,
        cellI,
        tetFaceI,
        tetPtI
    );
124

125
126
    label procI = -1;

127
128
    if (cellI >= 0)
    {
129
130
        procI = Pstream::myProcNo();
    }
131

132
    reduce(procI, maxOp<label>());
133

graham's avatar
graham committed
134
135
    // Ensure that only one processor attempts to insert this Parcel

136
137
138
    if (procI != Pstream::myProcNo())
    {
        cellI = -1;
139
140
        tetFaceI = -1;
        tetPtI = -1;
141
142
    }

graham's avatar
graham committed
143
144
    // Last chance - find nearest cell and try that one - the point is
    // probably on an edge
145
    if (procI == -1)
146
    {
147
        cellI = this->owner().mesh().findNearestCell(position);
graham's avatar
graham committed
148

149
150
        if (cellI >= 0)
        {
151
            position += SMALL*(cellCentres[cellI] - position);
152

153
            if (this->owner().mesh().pointInCell(position, cellI))
154
155
156
157
            {
                procI = Pstream::myProcNo();
            }
        }
158

159
        reduce(procI, maxOp<label>());
160

161
162
163
        if (procI != Pstream::myProcNo())
        {
            cellI = -1;
164
165
            tetFaceI = -1;
            tetPtI = -1;
166
167
168
        }
    }

169
    if (procI == -1)
170
    {
171
172
173
174
175
176
177
        if (errorOnNotFound)
        {
            FatalErrorIn
            (
                "Foam::InjectionModel<CloudType>::findCellAtPosition"
                "("
                    "label&, "
178
179
180
181
                    "label&, "
                    "label&, "
                    "vector&, "
                    "bool"
182
183
184
185
186
187
188
189
190
                ")"
            )   << "Cannot find parcel injection cell. "
                << "Parcel position = " << p0 << nl
                << abort(FatalError);
        }
        else
        {
            return false;
        }
191
    }
192
193

    return true;
194
195
196
}


197
198
199
template<class CloudType>
Foam::scalar Foam::InjectionModel<CloudType>::setNumberOfParticles
(
Andrew Heather's avatar
Andrew Heather committed
200
    const label parcels,
andy's avatar
andy committed
201
    const scalar volumeFraction,
Andrew Heather's avatar
updates    
Andrew Heather committed
202
    const scalar diameter,
Andrew Heather's avatar
updates    
Andrew Heather committed
203
    const scalar rho
204
205
206
207
208
209
210
)
{
    scalar nP = 0.0;
    switch (parcelBasis_)
    {
        case pbMass:
        {
andy's avatar
andy committed
211
212
213
214
            scalar volumep = pi/6.0*pow3(diameter);
            scalar volumeTot = massTotal_/rho;

            nP = (volumeFraction*volumeTot + delayedVolume_)/(parcels*volumep);
215
216
217
218
            break;
        }
        case pbNumber:
        {
219
220
221
222
223
            nP = massTotal_/(rho*volumeTotal_);
            break;
        }
        case pbFixed:
        {
224
            nP = nParticleFixed_;
225
226
227
228
229
230
231
            break;
        }
        default:
        {
            nP = 0.0;
            FatalErrorIn
            (
Andrew Heather's avatar
Andrew Heather committed
232
                "Foam::scalar "
Andrew Heather's avatar
Andrew Heather committed
233
234
                "Foam::InjectionModel<CloudType>::setNumberOfParticles"
                "("
Andrew Heather's avatar
Andrew Heather committed
235
236
237
238
                    "const label, "
                    "const scalar, "
                    "const scalar, "
                    "const scalar"
Andrew Heather's avatar
Andrew Heather committed
239
                ")"
240
241
242
243
244
245
246
247
248
249
            )<< "Unknown parcelBasis type" << nl
             << exit(FatalError);
        }
    }

    return nP;
}


template<class CloudType>
250
251
252
253
254
void Foam::InjectionModel<CloudType>::postInjectCheck
(
    const label parcelsAdded,
    const scalar massAdded
)
255
{
256
257
258
    const label allParcelsAdded = returnReduce(parcelsAdded, sumOp<label>());

    if (allParcelsAdded > 0)
259
    {
260
        Info<< nl
261
262
            << "--> Cloud: " << this->owner().name()
            << " injector: " << this->modelName() << nl
263
            << "    Added " << allParcelsAdded << " new parcels" << nl << endl;
264
265
266
    }

    // Increment total number of parcels added
267
268
269
270
    parcelsAddedTotal_ += allParcelsAdded;

    // Increment total mass injected
    massInjected_ += returnReduce(massAdded, sumOp<scalar>());
271
272

    // Update time for start of next injection
273
    time0_ = this->owner().db().time().value();
274
275
276
277
278
279
280
281

    // Increment number of injections
    nInjections_++;
}


// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

282
283
284
template<class CloudType>
Foam::InjectionModel<CloudType>::InjectionModel(CloudType& owner)
:
285
    SubModelBase<CloudType>(owner),
286
287
288
    SOI_(0.0),
    volumeTotal_(0.0),
    massTotal_(0.0),
289
    massFlowRate_(owner.db().time(), "massFlowRate"),
290
291
    massInjected_(this->template getModelProperty<scalar>("massInjected")),
    nInjections_(this->template getModelProperty<label>("nInjections")),
292
293
    parcelsAddedTotal_
    (
294
        this->template getModelProperty<scalar>("parcelsAddedTotal")
295
    ),
296
    parcelBasis_(pbNumber),
297
    nParticleFixed_(0.0),
298
    time0_(0.0),
299
300
    timeStep0_(this->template getModelProperty<scalar>("timeStep0")),
    delayedVolume_(0.0)
301
{}
302
303


304
305
306
307
308
template<class CloudType>
Foam::InjectionModel<CloudType>::InjectionModel
(
    const dictionary& dict,
    CloudType& owner,
309
310
    const word& modelName,
    const word& modelType
311
)
312
:
313
    SubModelBase<CloudType>(modelName, owner, dict, typeName, modelType),
314
    SOI_(0.0),
315
    volumeTotal_(0.0),
316
    massTotal_(0.0),
317
    massFlowRate_(owner.db().time(), "massFlowRate"),
318
319
    massInjected_(this->template getModelProperty<scalar>("massInjected")),
    nInjections_(this->template getModelProperty<scalar>("nInjections")),
320
321
    parcelsAddedTotal_
    (
322
        this->template getModelProperty<scalar>("parcelsAddedTotal")
323
    ),
324
    parcelBasis_(pbNumber),
325
    nParticleFixed_(0.0),
326
    time0_(owner.db().time().value()),
327
328
    timeStep0_(this->template getModelProperty<scalar>("timeStep0")),
    delayedVolume_(0.0)
329
{
Andrew Heather's avatar
Andrew Heather committed
330
331
332
333
334
335
    // Provide some info
    // - also serves to initialise mesh dimensions - needed for parallel runs
    //   due to lazy evaluation of valid mesh dimensions
    Info<< "    Constructing " << owner.mesh().nGeometricD() << "-D injection"
        << endl;

336
337
338
339
    if (owner.solution().transient())
    {
        this->coeffDict().lookup("massTotal") >> massTotal_;
        this->coeffDict().lookup("SOI") >> SOI_;
340
        SOI_ = owner.db().time().userTimeToTime(SOI_);
341
342
343
    }
    else
    {
344
345
        massFlowRate_.reset(this->coeffDict());
        massTotal_ = massFlowRate_.value(owner.db().time().value());
346
347
    }

348
    const word parcelBasisType = this->coeffDict().lookup("parcelBasisType");
349

350
    if (parcelBasisType == "mass")
351
352
353
    {
        parcelBasis_ = pbMass;
    }
354
    else if (parcelBasisType == "number")
355
356
357
    {
        parcelBasis_ = pbNumber;
    }
358
359
360
361
    else if (parcelBasisType == "fixed")
    {
        parcelBasis_ = pbFixed;

362
        Info<< "    Choosing nParticle to be a fixed value, massTotal "
363
364
365
            << "variable now does not determine anything."
            << endl;

366
        nParticleFixed_ = readScalar(this->coeffDict().lookup("nParticle"));
367
    }
368
369
370
371
    else
    {
        FatalErrorIn
        (
Andrew Heather's avatar
Andrew Heather committed
372
373
374
375
376
            "Foam::InjectionModel<CloudType>::InjectionModel"
            "("
                "const dictionary&, "
                "CloudType&, "
                "const word&"
Andrew Heather's avatar
Andrew Heather committed
377
            ")"
378
        )<< "parcelBasisType must be either 'number', 'mass' or 'fixed'" << nl
379
380
381
382
383
         << exit(FatalError);
    }
}


384
385
386
387
388
389
390
391
392
393
template<class CloudType>
Foam::InjectionModel<CloudType>::InjectionModel
(
    const InjectionModel<CloudType>& im
)
:
    SubModelBase<CloudType>(im),
    SOI_(im.SOI_),
    volumeTotal_(im.volumeTotal_),
    massTotal_(im.massTotal_),
394
    massFlowRate_(im.massFlowRate_),
395
396
397
398
399
400
    massInjected_(im.massInjected_),
    nInjections_(im.nInjections_),
    parcelsAddedTotal_(im.parcelsAddedTotal_),
    parcelBasis_(im.parcelBasis_),
    nParticleFixed_(im.nParticleFixed_),
    time0_(im.time0_),
401
402
    timeStep0_(im.timeStep0_),
    delayedVolume_(im.delayedVolume_)
403
404
405
{}


406
407
408
409
410
411
412
413
414
// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

template<class CloudType>
Foam::InjectionModel<CloudType>::~InjectionModel()
{}


// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

415
416
417
418
419
420
421
template<class CloudType>
void Foam::InjectionModel<CloudType>::updateMesh()
{
    // do nothing
}


422
423
424
425
426
427
428
429
430
431
432
433
template<class CloudType>
Foam::scalar Foam::InjectionModel<CloudType>::timeEnd() const
{
    notImplemented
    (
        "Foam::scalar Foam::InjectionModel<CloudType>::timeEnd() const"
    );

    return 0.0;
}


434
435
436
437
438
439
440
441
442
443
444
445
446
template<class CloudType>
Foam::label Foam::InjectionModel<CloudType>::parcelsToInject
(
    const scalar time0,
    const scalar time1
)
{
    notImplemented
    (
        "Foam::label Foam::InjectionModel<CloudType>::parcelsToInject"
        "("
            "const scalar, "
            "const scalar"
andy's avatar
andy committed
447
        ")"
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
    );

    return 0;
}


template<class CloudType>
Foam::scalar Foam::InjectionModel<CloudType>::volumeToInject
(
    const scalar time0,
    const scalar time1
)
{
    notImplemented
    (
        "Foam::scalar Foam::InjectionModel<CloudType>::volumeToInject"
        "("
            "const scalar, "
            "const scalar"
andy's avatar
andy committed
467
        ")"
468
469
470
471
472
473
474
475
476
    );

    return 0.0;
}


template<class CloudType>
Foam::scalar Foam::InjectionModel<CloudType>::averageParcelMass()
{
477
478
479
480
481
482
483
484
485
486
    label nTotal = 0.0;
    if (this->owner().solution().transient())
    {
        nTotal = parcelsToInject(0.0, timeEnd() - timeStart());
    }
    else
    {
        nTotal = parcelsToInject(0.0, 1.0);
    }

487
488
489
490
    return massTotal_/nTotal;
}


491
492
493
494
template<class CloudType>
template<class TrackData>
void Foam::InjectionModel<CloudType>::inject(TrackData& td)
{
495
    if (!this->active())
Andrew Heather's avatar
Andrew Heather committed
496
497
498
499
    {
        return;
    }

500
    const scalar time = this->owner().db().time().value();
501
502

    // Prepare for next time step
503
504
    label parcelsAdded = 0;
    scalar massAdded = 0.0;
Andrew Heather's avatar
Andrew Heather committed
505
    label newParcels = 0;
andy's avatar
andy committed
506
    scalar newVolumeFraction = 0.0;
507

andy's avatar
andy committed
508
    if (prepareForNextTimeStep(time, newParcels, newVolumeFraction))
509
    {
510
511
        scalar delayedVolume = 0;

andy's avatar
andy committed
512
513
514
515
        const scalar trackTime = this->owner().solution().trackTime();
        const polyMesh& mesh = this->owner().mesh();
        typename TrackData::cloudType& cloud = td.cloud();

516
517
518
        // Duration of injection period during this timestep
        const scalar deltaT =
            max(0.0, min(trackTime, min(time - SOI_, timeEnd() - time0_)));
519

520
521
        // Pad injection time if injection starts during this timestep
        const scalar padTime = max(0.0, SOI_ - time0_);
522

523
524
525
526
        // Introduce new parcels linearly across carrier phase timestep
        for (label parcelI = 0; parcelI < newParcels; parcelI++)
        {
            if (validInjection(parcelI))
527
            {
528
529
530
531
532
533
534
535
                // Calculate the pseudo time of injection for parcel 'parcelI'
                scalar timeInj = time0_ + padTime + deltaT*parcelI/newParcels;

                // Determine the injection position and owner cell,
                // tetFace and tetPt
                label cellI = -1;
                label tetFaceI = -1;
                label tetPtI = -1;
Andrew Heather's avatar
Andrew Heather committed
536

537
                vector pos = vector::zero;
538

539
                setPositionAndCell
540
                (
541
542
543
                    parcelI,
                    newParcels,
                    timeInj,
544
545
546
547
548
                    pos,
                    cellI,
                    tetFaceI,
                    tetPtI
                );
549

550
551
552
                if (cellI > -1)
                {
                    // Lagrangian timestep
553
                    const scalar dt = time - timeInj;
554

555
556
                    // Apply corrections to position for 2-D cases
                    meshTools::constrainToMeshCentre(mesh, pos);
557

558
                    // Create a new parcel
andy's avatar
andy committed
559
560
                    parcelType* pPtr =
                        new parcelType(mesh, pos, cellI, tetFaceI, tetPtI);
561

562
563
564
565
566
567
568
569
                    // Check/set new parcel thermo properties
                    cloud.setParcelThermoProperties(*pPtr, dt);

                    // Assign new parcel properties in injection model
                    setProperties(parcelI, newParcels, timeInj, *pPtr);

                    // Check/set new parcel injection properties
                    cloud.checkParcelProperties(*pPtr, dt, fullyDescribed());
570

571
572
                    // Apply correction to velocity for 2-D cases
                    meshTools::constrainDirection
573
                    (
574
575
576
                        mesh,
                        mesh.solutionD(),
                        pPtr->U()
577
                    );
Andrew Heather's avatar
Andrew Heather committed
578

579
580
581
582
583
                    // Number of particles per parcel
                    pPtr->nParticle() =
                        setNumberOfParticles
                        (
                            newParcels,
andy's avatar
andy committed
584
                            newVolumeFraction,
585
586
587
588
                            pPtr->d(),
                            pPtr->rho()
                        );

589
                    if (pPtr->nParticle() >= 1.0)
590
                    {
591
592
593
594
                        parcelsAdded++;
                        massAdded += pPtr->nParticle()*pPtr->mass();

                        if (pPtr->move(td, dt))
595
596
597
598
599
600
601
                        {
                            td.cloud().addParticle(pPtr);
                        }
                        else
                        {
                            delete pPtr;
                        }
602
                    }
603
604
605
606
607
                    else
                    {
                        delayedVolume += pPtr->nParticle()*pPtr->volume();
                        delete pPtr;
                    }
608
                }
Andrew Heather's avatar
Andrew Heather committed
609
            }
610
        }
611
612

        delayedVolume_ = delayedVolume;
613
614
    }

615
    postInjectCheck(parcelsAdded, massAdded);
616
617
618
}


619
620
621
622
623
624
625
626
627
628
629
630
631
632
template<class CloudType>
template<class TrackData>
void Foam::InjectionModel<CloudType>::injectSteadyState
(
    TrackData& td,
    const scalar trackTime
)
{
    if (!this->active())
    {
        return;
    }

    const polyMesh& mesh = this->owner().mesh();
633
    typename TrackData::cloudType& cloud = td.cloud();
634

635
636
    massTotal_ = massFlowRate_.value(mesh.time().value());

637
638
639
640
641
642
643
644
645
646
647
648
    // Reset counters
    time0_ = 0.0;
    label parcelsAdded = 0;
    scalar massAdded = 0.0;

    // Set number of new parcels to inject based on first second of injection
    label newParcels = parcelsToInject(0.0, 1.0);

    // Inject new parcels
    for (label parcelI = 0; parcelI < newParcels; parcelI++)
    {
        // Volume to inject is split equally amongst all parcel streams
andy's avatar
andy committed
649
        scalar newVolumeFraction = 1.0/scalar(newParcels);
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675

        // Determine the injection position and owner cell,
        // tetFace and tetPt
        label cellI = -1;
        label tetFaceI = -1;
        label tetPtI = -1;

        vector pos = vector::zero;

        setPositionAndCell
        (
            parcelI,
            newParcels,
            0.0,
            pos,
            cellI,
            tetFaceI,
            tetPtI
        );

        if (cellI > -1)
        {
            // Apply corrections to position for 2-D cases
            meshTools::constrainToMeshCentre(mesh, pos);

            // Create a new parcel
andy's avatar
andy committed
676
677
            parcelType* pPtr =
                new parcelType(mesh, pos, cellI, tetFaceI, tetPtI);
678

andy's avatar
andy committed
679
            // Check/set new parcel thermo properties
680
            cloud.setParcelThermoProperties(*pPtr, 0.0);
681

682
683
684
            // Assign new parcel properties in injection model
            setProperties(parcelI, newParcels, 0.0, *pPtr);

andy's avatar
andy committed
685
            // Check/set new parcel injection properties
686
687
            cloud.checkParcelProperties(*pPtr, 0.0, fullyDescribed());

688
            // Apply correction to velocity for 2-D cases
andy's avatar
andy committed
689
            meshTools::constrainDirection(mesh, mesh.solutionD(), pPtr->U());
690
691
692
693
694
695

            // Number of particles per parcel
            pPtr->nParticle() =
                setNumberOfParticles
                (
                    1,
andy's avatar
andy committed
696
                    newVolumeFraction,
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
                    pPtr->d(),
                    pPtr->rho()
                );

            // Add the new parcel
            td.cloud().addParticle(pPtr);

            massAdded += pPtr->nParticle()*pPtr->mass();
            parcelsAdded++;
        }
    }

    postInjectCheck(parcelsAdded, massAdded);
}


713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
template<class CloudType>
void Foam::InjectionModel<CloudType>::setPositionAndCell
(
    const label parcelI,
    const label nParcels,
    const scalar time,
    vector& position,
    label& cellOwner,
    label& tetFaceI,
    label& tetPtI
)
{
    notImplemented
    (
        "void Foam::InjectionModel<CloudType>::setPositionAndCell"
        "("
            "const label, "
            "const label, "
            "const scalar, "
            "vector&, "
            "label&, "
            "label&, "
            "label&"
        ")"
    );
}


template<class CloudType>
void Foam::InjectionModel<CloudType>::setProperties
(
    const label parcelI,
    const label nParcels,
    const scalar time,
    typename CloudType::parcelType& parcel
)
{
    notImplemented
    (
        "void Foam::InjectionModel<CloudType>::setProperties"
        "("
            "const label, "
            "const label, "
            "const scalar, "
            "typename CloudType::parcelType&"
        ")"
    );
}


template<class CloudType>
bool Foam::InjectionModel<CloudType>::fullyDescribed() const
{
    notImplemented
    (
        "bool Foam::InjectionModel<CloudType>::fullyDescribed() const"
    );

    return false;
}


775
template<class CloudType>
776
void Foam::InjectionModel<CloudType>::info(Ostream& os)
777
{
778
    os  << "    " << this->modelName() << ":" << nl
andy's avatar
andy committed
779
780
        << "        number of parcels added     = " << parcelsAddedTotal_ << nl
        << "        mass introduced             = " << massInjected_ << nl;
781

782
    if (this->outputTime())
783
    {
784
785
786
787
        this->setModelProperty("massInjected", massInjected_);
        this->setModelProperty("nInjections", nInjections_);
        this->setModelProperty("parcelsAddedTotal", parcelsAddedTotal_);
        this->setModelProperty("timeStep0", timeStep0_);
788
    }
789
790
791
}


792
793
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

794
#include "InjectionModelNew.C"
795
796

// ************************************************************************* //