From 389b44d0e5ff2e4572fc5dee9ec4007891001aef Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 29 Jul 2022 13:46:33 +0200 Subject: [PATCH 01/23] [doc] constDiffSmoothing doc update --- doc/smoothingModel_constDiffSmoothing.txt | 20 +++++++++++++++----- 1 file changed, 15 insertions(+), 5 deletions(-) diff --git a/doc/smoothingModel_constDiffSmoothing.txt b/doc/smoothingModel_constDiffSmoothing.txt index 9ff82ecd..bfefed13 100644 --- a/doc/smoothingModel_constDiffSmoothing.txt +++ b/doc/smoothingModel_constDiffSmoothing.txt @@ -15,17 +15,22 @@ dictionary. smoothingModel constDiffSmoothing; constDiffSmoothingProps \{ - lowerLimit number1; - upperLimit number2; - smoothingLength lengthScale; - smoothingLengthReferenceField lengthScaleRefField; + lowerLimit number1; + upperLimit number2; + smoothingLength lengthScale; + smoothingLengthReferenceField lengthScaleRefField; + smoothingLengthFieldName fieldName1; + smoothingLengthReferenceFieldName fieldName2; verbose; \} :pre {number1} = scalar fields will be bound to this lower value :ulb,l {number2} = scalar fields will be bound to this upper value :l {lengthScale} = length scale over which the exchange fields will be smoothed out :l -{lengthScaleRefField} = length scale over which reference fields (e.g., the average particle velocity) will be smoothed out. Should be always larger than lengthScale. If not specified, will be equal to lengthScale. :l +{lengthScaleRefField} = (optional) length scale over which reference fields (e.g., the average particle velocity) will be smoothed out. Should be always larger than lengthScale. If not specified, will be equal to lengthScale. :l +{fieldName1} = (optional) name of scalar field to be used as local smoothing length. :l +{fieldName2} = (optional) name of scalar field to be used as local smoothing length for reference fields. :l + {verbose} = (optional, default false) flag for debugging output :l :ule @@ -51,6 +56,11 @@ which these reference fields are not specified. Values calculated in the cells e.g. the average particle velocity, which are not specified in all cells in case the flow is rather dilute. +Alternative to {smoothingLength} and {smoothingLengthReferenceField}, +{smoothingLengthFieldName} and/or {smoothingLengthReferenceFieldName} can be used +to define spatial variation of the smoothing lengths. Either the scalar or field +options must be used, giving both will result in errors. + [Restrictions:] This model is tested in a limited number of flow situations. From 3107a9ce10915caea4c13973c23240edab3ad81c Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 29 Jul 2022 13:57:58 +0200 Subject: [PATCH 02/23] Update smoothingModel_constDiffSmoothing.txt --- doc/smoothingModel_constDiffSmoothing.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/smoothingModel_constDiffSmoothing.txt b/doc/smoothingModel_constDiffSmoothing.txt index bfefed13..34b0c3e4 100644 --- a/doc/smoothingModel_constDiffSmoothing.txt +++ b/doc/smoothingModel_constDiffSmoothing.txt @@ -18,7 +18,7 @@ constDiffSmoothingProps lowerLimit number1; upperLimit number2; smoothingLength lengthScale; - smoothingLengthReferenceField lengthScaleRefField; + smoothingLengthReference lengthScaleRefField; smoothingLengthFieldName fieldName1; smoothingLengthReferenceFieldName fieldName2; verbose; From e95e1ec6adb184979edc8c3b67535ffcbb2c50a2 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 09:57:58 +0200 Subject: [PATCH 03/23] clean up ParmarBassetForce --- .../ParmarBassetForce/ParmarBassetForce.C | 38 ------------------- 1 file changed, 38 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index 79a72b32..b311f7f6 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -134,10 +134,6 @@ ParmarBassetForce::ParmarBassetForce particleCloud_.probeM().vectorFields_.append("ParmarBassetForce"); //first entry must the be the force particleCloud_.probeM().vectorFields_.append("Urel"); particleCloud_.probeM().vectorFields_.append("ddtUrel"); - // - particleCloud_.probeM().vectorFields_.append("UrelNoSmooth"); - particleCloud_.probeM().vectorFields_.append("ddtUrelNoSmooth"); - // particleCloud_.probeM().vectorFields_.append("Fshort"); particleCloud_.probeM().vectorFields_.append("Flong1"); particleCloud_.probeM().vectorFields_.append("Flong2"); @@ -220,22 +216,12 @@ void ParmarBassetForce::setForce() const Urel_ = Us_ - U_; ddtUrel_ = fvc::ddt(Us_) - fvc::ddt(U_) - (Us_ & fvc::grad(U_)); - // - volVectorField UrelNoSmooth_ = Urel_; - volVectorField ddtUrelNoSmooth_ = ddtUrel_; - // - smoothingM().smoothen(Urel_); smoothingM().smoothen(ddtUrel_); interpolationCellPoint UrelInterpolator_(Urel_); interpolationCellPoint ddtUrelInterpolator_(ddtUrel_); - // - interpolationCellPoint UrelNoSmoothInterpolator_(UrelNoSmooth_); - interpolationCellPoint ddtUrelNoSmoothInterpolator_(ddtUrelNoSmooth_); - // - for(int index = 0;index < particleCloud_.numberOfParticles(); index++) { vector ParmarBassetForce(0,0,0); @@ -412,34 +398,10 @@ void ParmarBassetForce::setForce() const Flong2[i] = FHist_[1][0][index][i]; } - // - // relative velocity (m/s) - vector UrelNoSmooth; - vector ddtUrelNoSmooth; - - if(forceSubM(0).interpolation()) - UrelNoSmooth = UrelNoSmoothInterpolator_.interpolate(position,cellI); - else - UrelNoSmooth = UrelNoSmooth_[cellI]; - - // acceleration (m/s2) - if(forceSubM(0).interpolation()) - ddtUrelNoSmooth = ddtUrelNoSmoothInterpolator_.interpolate(position,cellI); - else - ddtUrelNoSmooth = ddtUrelNoSmooth_[cellI]; - - UrelNoSmooth /= mps; - ddtUrelNoSmooth /= mpss; - // - #include "setupProbeModelfields.H" vValues.append(ParmarBassetForce); //first entry must the be the force vValues.append(Urel); vValues.append(ddtUrel); - // - vValues.append(UrelNoSmooth); - vValues.append(ddtUrelNoSmooth); - // vValues.append(Fshort); vValues.append(Flong1); vValues.append(Flong2); From d35249aa3c661440190c42b7415bb202a89ce415 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 10:03:11 +0200 Subject: [PATCH 04/23] clean up ParmarBassetForce --- .../forceModel/ParmarBassetForce/ParmarBassetForce.C | 10 ---------- 1 file changed, 10 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index b311f7f6..1d8ff351 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -184,8 +184,6 @@ void ParmarBassetForce::setForce() const vector Urel(0,0,0); vector ddtUrel(0,0,0); - scalar t0min = 0.00; - const volScalarField& nufField = forceSubM(0).nuField(); const volScalarField& rhoField = forceSubM(0).rhoField(); @@ -294,14 +292,6 @@ void ParmarBassetForce::setForce() const update_ddtUrelHist(ddtUrel,index); // add current dim.less ddtUrel to history update_rHist(r,index); // add current r to history - // warning and reset for too small t0 - if (t0(FHistSize_,NULL)), // FHist_[k={1,2}-1][ndt in past][particle ID][dim] - gH0_(NULL), - tRef_(NULL), - mRef_(NULL), - lRef_(NULL), + ddtUrelHistRegName_(typeName + "ddtUrelHist"), // indexed as: ddtUrelHist_[particle ID][iDim*nHist_+iHist] + rHistRegName_(typeName + "rHist"), // indexed as: rHist_[particle ID][iHist] + FHistRegName_(typeName + "FHist"), // indexed as: ddtUrelHist_[particle ID][iDim*FHistSize_*2+iHist*2+k] + gH0RegName_(typeName + "gH0"), + tRefRegName_(typeName + "tRef"), + mRefRegName_(typeName + "mRef"), + lRefRegName_(typeName + "lRef"), Urel_ ( IOobject ( @@ -116,6 +116,15 @@ ParmarBassetForce::ParmarBassetForce ) { + // allocate particle properties + particleCloud_.registerParticleProperty(ddtUrelHistRegName_, 3*nHist_,NOTONCPU,false); + particleCloud_.registerParticleProperty( rHistRegName_, nHist_,NOTONCPU,false); + particleCloud_.registerParticleProperty( FHistRegName_,6*FHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( gH0RegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( tRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( mRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( lRefRegName_, 1 ,NOTONCPU,false); + // init force sub model setForceSubModels(propsDict_); // define switches which can be read from dict @@ -154,20 +163,6 @@ ParmarBassetForce::ParmarBassetForce ParmarBassetForce::~ParmarBassetForce() { - particleCloud_.dataExchangeM().destroy(gH0_, 1); - particleCloud_.dataExchangeM().destroy(tRef_, 1); - particleCloud_.dataExchangeM().destroy(mRef_, 1); - particleCloud_.dataExchangeM().destroy(lRef_, 1); - - for (int i=0; i(ddtUrelHistRegName_); + double**& rHist_ = particleCloud_.getParticlePropertyRef( rHistRegName_); + double**& FHist_ = particleCloud_.getParticlePropertyRef( FHistRegName_); + double**& gH0_ = particleCloud_.getParticlePropertyRef( gH0RegName_); + double**& tRef_ = particleCloud_.getParticlePropertyRef( tRefRegName_); + double**& mRef_ = particleCloud_.getParticlePropertyRef( mRefRegName_); + double**& lRef_ = particleCloud_.getParticlePropertyRef( lRefRegName_); vector position(0,0,0); vector Urel(0,0,0); @@ -256,8 +255,8 @@ void ParmarBassetForce::setForce() const if (r<0.25 || r>2.0) { - gH0_[index][0] = NOTONCPU; //reset reference - ddtUrelHist_[0][index][0] = NOTONCPU; //reset ddtU history (only component used for checking nKnown) + gH0_[index][0] = NOTONCPU; //reset reference + ddtUrelHist_[index][0] = NOTONCPU; //reset ddtU history (only component used for checking nKnown) } } @@ -269,7 +268,7 @@ void ParmarBassetForce::setForce() const scalar Vs = rs*rs*rs*M_PI*4/3; scalar mRef = Vs*rhoField[cellI] * gH * 5.2863; // 9/(2*sqrt(pi))*(256/pi)^(1/6) = 5.2863 (Eq. 3.2) - gH0_[index][0] = gH; + gH0_[index][0] = gH; tRef_[index][0] = tRef; mRef_[index][0] = mRef; lRef_[index][0] = rs; @@ -289,14 +288,14 @@ void ParmarBassetForce::setForce() const ddtUrel /= mpss; // update ddtUrel and r history - update_ddtUrelHist(ddtUrel,index); // add current dim.less ddtUrel to history - update_rHist(r,index); // add current r to history + update_ddtUrelHist(ddtUrelHist_,ddtUrel,index); // add current dim.less ddtUrel to history + update_rHist(rHist_,r,index); // add current r to history // check length of known history int nKnown = 0; for (int j=0; j0) { for (int i=0; i<3; i++) // loop over dimensions - Fshort[i] = -calculateK0(r,dt) * ddtUrelHist_[0][index][i]; + Fshort[i] = -calculateK0(r,dt) * ddtUrelHist_[index][i*nHist_]; } // int_dt^t0 K(r,xi) * ddtU(t-xi) dxi (trapezoid rule) @@ -320,32 +319,32 @@ void ParmarBassetForce::setForce() const for (int j=1; j0; j--) // loop over past times - ddtUrelHist_[j][index][i] = ddtUrelHist_[j-1][index][i]; + ddtUrelHist_[index][i*nHist_+j] = ddtUrelHist_[index][i*nHist_+j-1]; - ddtUrelHist_[0][index][i] = ddtUrel[i]; + ddtUrelHist_[index][i*nHist_] = ddtUrel[i]; } } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -void Foam::ParmarBassetForce::update_rHist(scalar r, int index) const +void Foam::ParmarBassetForce::update_rHist(double**& rHist_, scalar r, int index) const { for (int j=nHist_-1; j>0; j--) // loop over past times - rHist_[j][index][0] = rHist_[j-1][index][0]; + rHist_[index][j] = rHist_[index][j-1]; - rHist_[0][index][0] = r; + rHist_[index][0] = r; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -void Foam::ParmarBassetForce::update_FHist(const vector& F1, const vector& F2, int index) const +void Foam::ParmarBassetForce::update_FHist(double**& FHist_, const vector& F1, const vector& F2, int index) const { for (int i=0; i<3; i++) // loop over dimensions { for (int k=0; k<2; k++) // loop over F1, F2 for (int j=FHistSize_-1; j>0; j--) // loop over past times - FHist_[k][j][index][i] = FHist_[k][j-1][index][i]; + FHist_[index][i*FHistSize_*2+j*2+k] = FHist_[index][i*FHistSize_*2+(j-1)*2+k]; - FHist_[0][0][index][i] = F1[i]; - FHist_[1][0][index][i] = F2[i]; + FHist_[index][i*FHistSize_*2 ] = F1[i]; + FHist_[index][i*FHistSize_*2+1] = F2[i]; } } @@ -516,18 +494,18 @@ void Foam::ParmarBassetForce::calculateCoeffs(int k, scalar t0, scalar r, double } } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -void Foam::ParmarBassetForce::solveFlongODE(int k, double C[4], scalar dt, int index) const +void Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUrelHist_, int k, double C[4], scalar dt, int index) const { if (discOrder_==1) { for (int i=0; i<3; i++) // loop over dimensions { - FHist_[k][0][index][i] = + FHist_[index][i*FHistSize_*2+k] = ( - ( C[1]/dt+2/(dt*dt)) * FHist_[k][1][index][i] - - ( 1/(dt*dt)) * FHist_[k][2][index][i] - - (C[2]+C[3]/dt ) * ddtUrelHist_[nInt_ ][index][i] - + ( C[3]/dt ) * ddtUrelHist_[nInt_+1][index][i] + ( C[1]/dt+2/(dt*dt)) * FHist_[index][i*FHistSize_*2+1*2+k] + - ( 1/(dt*dt)) * FHist_[index][i*FHistSize_*2+2*2+k] + - (C[2]+C[3]/dt ) * ddtUrelHist_[index][i*nHist_+nInt_ ] + + ( C[3]/dt ) * ddtUrelHist_[index][i*nHist_+nInt_+1] ) / (C[0]+C[1]/dt+1/(dt*dt)); // Eq. 3.20 using first order temporal discretisation } } @@ -535,41 +513,21 @@ void Foam::ParmarBassetForce::solveFlongODE(int k, double C[4], scalar dt, int i { for (int i=0; i<3; i++) // loop over dimensions { - FHist_[k][0][index][i] = + FHist_[index][i*FHistSize_*2+k] = ( - ( 4*C[1]/(2*dt) + 24/(4*dt*dt)) * FHist_[k][1][index][i] - - ( C[1]/(2*dt) + 22/(4*dt*dt)) * FHist_[k][2][index][i] - + ( 8/(4*dt*dt)) * FHist_[k][3][index][i] - - ( 1/(4*dt*dt)) * FHist_[k][4][index][i] + ( 4*C[1]/(2*dt) + 24/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+1*2+k] + - ( C[1]/(2*dt) + 22/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+2*2+k] + + ( 8/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+3*2+k] + - ( 1/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+4*2+k] - - (C[2] + 3*C[3]/(2*dt) ) * ddtUrelHist_[nInt_ ][index][i] - + ( 4*C[3]/(2*dt) ) * ddtUrelHist_[nInt_+1][index][i] - - ( C[3]/(2*dt) ) * ddtUrelHist_[nInt_+2][index][i] + - (C[2] + 3*C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_ ] + + ( 4*C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_+1] + - ( C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_+2] ) / (C[0] + 3*C[1]/(2*dt) + 9/(4*dt*dt)); // Eq. 3.20 using second order temporal discretisation } } } -void Foam::ParmarBassetForce::rescaleHist(scalar tScale, scalar mScale, scalar lScale, scalar rScale, int index) const -{ - for (int i=0; i<3; i++) // loop over dimensions - { - // rescale ddtU history - for (int j=0; j ddtUrelHist_; + const word ddtUrelHistRegName_; - mutable List rHist_; + const word rHistRegName_; - mutable List> FHist_; + const word FHistRegName_; - mutable double** gH0_; + const word gH0RegName_; - mutable double** tRef_; + const word tRefRegName_; - mutable double** mRef_; + const word mRefRegName_; - mutable double** lRef_; + const word lRefRegName_; mutable volVectorField Urel_; @@ -101,17 +101,15 @@ private: scalar trapWeight(int i, int n) const; - void update_ddtUrelHist(const vector& ddtUrel, int index) const; + void update_ddtUrelHist(double**& ddtUrelHist_, const vector& ddtUrel, int index) const; - void update_rHist(scalar r, int index) const; + void update_rHist(double**& rHist_, scalar r, int index) const; - void update_FHist(const vector& F1, const vector& F2, int index) const; + void update_FHist(double**& FHist_, const vector& F1, const vector& F2, int index) const; void calculateCoeffs(int k, scalar t0, scalar r, double c[2][4][3], double chi[2][4][2], double C[4]) const; - void solveFlongODE(int k, double C[4], scalar dt, int index) const; - - void rescaleHist(scalar tScale, scalar mScale, scalar lScale, scalar rScale, int index) const; + void solveFlongODE(double**& FHist_, double**& ddtUrelHist_, int k, double C[4], scalar dt, int index) const; public: @@ -136,9 +134,6 @@ public: // Member Functions void setForce() const; - - void reAllocArrays() const; - inline const smoothingModel& smoothingM() const { return smoothingModel_; From c3f609b0e9180955cd1f5c430a5cf568570daab3 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 15:05:48 +0200 Subject: [PATCH 06/23] allocate uRelOld only when needed, save memory --- .../forceModel/virtualMassForce/virtualMassForce.C | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/virtualMassForce/virtualMassForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/virtualMassForce/virtualMassForce.C index 4ef36cf9..70789880 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/virtualMassForce/virtualMassForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/virtualMassForce/virtualMassForce.C @@ -99,7 +99,12 @@ virtualMassForce::virtualMassForce ) ) { - particleCloud_.registerParticleProperty(UrelOldRegName_,3,NOTONCPU,false); + // allocate UrelOld only if needed + int UrelOldSize = 0; + if(!splitUrelCalculation_ || !useUs_ ) + UrelOldSize = 3; + + particleCloud_.registerParticleProperty(UrelOldRegName_,UrelOldSize,NOTONCPU,false); // init force sub model setForceSubModels(propsDict_); From 7ce49bf21ce3b46f81fd5d53169314e47e6dfb18 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 15:20:29 +0200 Subject: [PATCH 07/23] minor cleanup --- .../subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C | 1 - 1 file changed, 1 deletion(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index d8603504..7263f28f 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -34,7 +34,6 @@ Description #include "ParmarBassetForce.H" #include "addToRunTimeSelectionTable.H" #include "smoothingModel.H" -#include "constDiffSmoothing.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // From 19e780f93bebacdf29b17865272c637d04080bd6 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 15:48:50 +0200 Subject: [PATCH 08/23] added verbose mode for ParmarBassetForce --- .../ParmarBassetForce/ParmarBassetForce.C | 20 +++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index 7263f28f..604cf816 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -129,6 +129,7 @@ ParmarBassetForce::ParmarBassetForce // define switches which can be read from dict forceSubM(0).setSwitchesList(SW_TREAT_FORCE_EXPLICIT,true); // activate treatExplicit switch forceSubM(0).setSwitchesList(SW_INTERPOLATION,true); // activate search for interpolate switch + forceSubM(0).setSwitchesList(SW_VERBOSE,true); // activate search for verbose switch forceSubM(0).readSwitches(); //Extra switches/settings @@ -372,6 +373,25 @@ void ParmarBassetForce::setForce() const } ParmarBassetForce *= newton; + if (forceSubM(0).verbose() && index >= 0 && index < 2) + { + vector Flong1(0,0,0); + vector Flong2(0,0,0); + + for (int i=0; i<3; i++) // loop over dimensions + { + Flong1[i] = FHist_[index][i*FHistSize_*2 ]*newton; + Flong2[i] = FHist_[index][i*FHistSize_*2+1]*newton; + } + + Pout << "cellI = " << cellI << endl; + Pout << "index = " << index << endl; + Pout << "Fshort = " << Fshort*newton << endl; + Pout << "Flong1 = " << Flong1 << endl; + Pout << "Flong2 = " << Flong2 << endl; + Pout << "Ftotal = " << ParmarBassetForce << endl; + } + // Set value fields and write the probe if(probeIt_) { From a10c773e3125f44fda58adc60b2dbc290c664c3f Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 15:49:16 +0200 Subject: [PATCH 09/23] [DOC] added verbose mode for ParmarBassetForce --- doc/forceModel_ParmarBassetForce.txt | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/doc/forceModel_ParmarBassetForce.txt b/doc/forceModel_ParmarBassetForce.txt index f18206ca..bc022f53 100644 --- a/doc/forceModel_ParmarBassetForce.txt +++ b/doc/forceModel_ParmarBassetForce.txt @@ -23,7 +23,8 @@ ParmarBassetForceProps nIntegral scalar1; discretisationOrder scalar2; treatForceExplicit switch1; - interpolation switch2; + verbose switch2; + interpolation switch3; smoothingModel "smoothingModel"; \} :pre @@ -33,6 +34,7 @@ ParmarBassetForceProps {scalar2} = (1 or 2) discretisation order of the far history differential equations :l {switch1} = (optional, default true) sub model switch, see "forceSubModel"_forceSubModel.html for details :l {switch2} = (optional, default false) sub model switch, see "forceSubModel"_forceSubModel.html for details :l +{switch3} = (optional, default false) sub model switch, see "forceSubModel"_forceSubModel.html for details :l {smoothingModel} = name of smoothing model for the dU/dt field :l :ule From dc2be65fe3cfc1c889fe8296b57d61ecd899a16b Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Fri, 21 Oct 2022 16:22:00 +0200 Subject: [PATCH 10/23] reorder FHist to save memory --- .../ParmarBassetForce/ParmarBassetForce.C | 67 +++++++------------ .../ParmarBassetForce/ParmarBassetForce.H | 2 +- 2 files changed, 26 insertions(+), 43 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index 604cf816..d841f3d1 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -73,7 +73,7 @@ ParmarBassetForce::ParmarBassetForce nInt_(readLabel(propsDict_.lookup("nIntegral"))), discOrder_(readLabel(propsDict_.lookup("discretisationOrder"))), nHist_(nInt_+discOrder_+1), - FHistSize_(2*discOrder_+1), + FHistSize_(2*discOrder_), ddtUrelHistRegName_(typeName + "ddtUrelHist"), // indexed as: ddtUrelHist_[particle ID][iDim*nHist_+iHist] rHistRegName_(typeName + "rHist"), // indexed as: rHist_[particle ID][iHist] FHistRegName_(typeName + "FHist"), // indexed as: ddtUrelHist_[particle ID][iDim*FHistSize_*2+iHist*2+k] @@ -223,6 +223,7 @@ void ParmarBassetForce::setForce() const { vector ParmarBassetForce(0,0,0); vector Fshort(0,0,0); + vector Flong[2]={vector::zero, vector::zero}; label cellI = particleCloud_.cellIDs()[index][0]; if (cellI > -1) // particle Found @@ -328,9 +329,6 @@ void ParmarBassetForce::setForce() const //********* long term force computing (differential form) *********// - // update F1, F2 history - update_FHist(FHist_,vector::zero,vector::zero,index); - // initialise ddtUrel(t0) and Flong(:) as 0 and r(t0) as 1 if (nKnown == nInt_) { @@ -358,37 +356,28 @@ void ParmarBassetForce::setForce() const calculateCoeffs(k,t0,rHist_[index][nInt_],c,chi,C); // solve Eq. 3.20 - solveFlongODE(FHist_,ddtUrelHist_,k,C,dt,index); + Flong[k] = solveFlongODE(FHist_,ddtUrelHist_,k,C,dt,index); } } + // update F1, F2 history + update_FHist(FHist_,Flong[0],Flong[1],index); + //********* total force *********// // sum and convert to N - for (int i=0; i<3; i++) // loop over dimensions - { - ParmarBassetForce[i] = Fshort[i]; - for (int k=0; k<2; k++) // loop over F1, F2 - ParmarBassetForce[i] += FHist_[index][i*FHistSize_*2+k]; - } + ParmarBassetForce = Fshort; + for (int k=0; k<2; k++) // loop over F1, F2 + ParmarBassetForce += Flong[k]; ParmarBassetForce *= newton; if (forceSubM(0).verbose() && index >= 0 && index < 2) { - vector Flong1(0,0,0); - vector Flong2(0,0,0); - - for (int i=0; i<3; i++) // loop over dimensions - { - Flong1[i] = FHist_[index][i*FHistSize_*2 ]*newton; - Flong2[i] = FHist_[index][i*FHistSize_*2+1]*newton; - } - Pout << "cellI = " << cellI << endl; Pout << "index = " << index << endl; Pout << "Fshort = " << Fshort*newton << endl; - Pout << "Flong1 = " << Flong1 << endl; - Pout << "Flong2 = " << Flong2 << endl; + Pout << "Flong1 = " << Flong[0]*newton << endl; + Pout << "Flong2 = " << Flong[1]*newton << endl; Pout << "Ftotal = " << ParmarBassetForce << endl; } @@ -397,22 +386,13 @@ void ParmarBassetForce::setForce() const { scalar ReRef = 0.75/(gH0_[index][0]-0.105); - vector Flong1(0,0,0); - vector Flong2(0,0,0); - - for (int i=0; i<3; i++) // loop over dimensions - { - Flong1[i] = FHist_[index][i*FHistSize_*2 ]; - Flong2[i] = FHist_[index][i*FHistSize_*2+1]; - } - #include "setupProbeModelfields.H" vValues.append(ParmarBassetForce); //first entry must the be the force vValues.append(Urel); vValues.append(ddtUrel); vValues.append(Fshort); - vValues.append(Flong1); - vValues.append(Flong2); + vValues.append(Flong[0]); + vValues.append(Flong[1]); sValues.append(ReRef); sValues.append(tRef_[index][0]); sValues.append(mRef_[index][0]); @@ -513,16 +493,18 @@ void Foam::ParmarBassetForce::calculateCoeffs(int k, scalar t0, scalar r, double } } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -void Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUrelHist_, int k, double C[4], scalar dt, int index) const +vector Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUrelHist_, int k, double C[4], scalar dt, int index) const { + vector Flong = vector::zero; + if (discOrder_==1) { for (int i=0; i<3; i++) // loop over dimensions { - FHist_[index][i*FHistSize_*2+k] = + Flong[i] = ( - ( C[1]/dt+2/(dt*dt)) * FHist_[index][i*FHistSize_*2+1*2+k] - - ( 1/(dt*dt)) * FHist_[index][i*FHistSize_*2+2*2+k] + ( C[1]/dt+2/(dt*dt)) * FHist_[index][i*FHistSize_*2 +k] + - ( 1/(dt*dt)) * FHist_[index][i*FHistSize_*2+2+k] - (C[2]+C[3]/dt ) * ddtUrelHist_[index][i*nHist_+nInt_ ] + ( C[3]/dt ) * ddtUrelHist_[index][i*nHist_+nInt_+1] ) / (C[0]+C[1]/dt+1/(dt*dt)); // Eq. 3.20 using first order temporal discretisation @@ -532,12 +514,12 @@ void Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUrelH { for (int i=0; i<3; i++) // loop over dimensions { - FHist_[index][i*FHistSize_*2+k] = + Flong[i] = ( - ( 4*C[1]/(2*dt) + 24/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+1*2+k] - - ( C[1]/(2*dt) + 22/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+2*2+k] - + ( 8/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+3*2+k] - - ( 1/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+4*2+k] + ( 4*C[1]/(2*dt) + 24/(4*dt*dt)) * FHist_[index][i*FHistSize_*2 +k] + - ( C[1]/(2*dt) + 22/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+2+k] + + ( 8/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+4+k] + - ( 1/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+6+k] - (C[2] + 3*C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_ ] + ( 4*C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_+1] @@ -545,6 +527,7 @@ void Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUrelH ) / (C[0] + 3*C[1]/(2*dt) + 9/(4*dt*dt)); // Eq. 3.20 using second order temporal discretisation } } + return Flong; } } // End namespace Foam diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H index 03051678..df6e20a8 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H @@ -109,7 +109,7 @@ private: void calculateCoeffs(int k, scalar t0, scalar r, double c[2][4][3], double chi[2][4][2], double C[4]) const; - void solveFlongODE(double**& FHist_, double**& ddtUrelHist_, int k, double C[4], scalar dt, int index) const; + vector solveFlongODE(double**& FHist_, double**& ddtUrelHist_, int k, double C[4], scalar dt, int index) const; public: From ab73bf86ac1ef834243f1b4653c680ea07182c23 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Mon, 24 Oct 2022 09:42:25 +0200 Subject: [PATCH 11/23] reorder ddtUrelHist and rHist to save memory --- .../ParmarBassetForce/ParmarBassetForce.C | 48 ++++++++----------- 1 file changed, 21 insertions(+), 27 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index d841f3d1..17052c57 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -72,7 +72,7 @@ ParmarBassetForce::ParmarBassetForce Us_(sm.mesh().lookupObject (UsFieldName_)), nInt_(readLabel(propsDict_.lookup("nIntegral"))), discOrder_(readLabel(propsDict_.lookup("discretisationOrder"))), - nHist_(nInt_+discOrder_+1), + nHist_(nInt_+discOrder_), FHistSize_(2*discOrder_), ddtUrelHistRegName_(typeName + "ddtUrelHist"), // indexed as: ddtUrelHist_[particle ID][iDim*nHist_+iHist] rHistRegName_(typeName + "rHist"), // indexed as: rHist_[particle ID][iHist] @@ -282,18 +282,12 @@ void ParmarBassetForce::setForce() const scalar dt = U_.mesh().time().deltaT().value() / tRef_[index][0]; // dim.less scalar t0 = nInt_*dt; // dim.less - //********* update histories *********// - // non-dimensionlise Urel /= mps; ddtUrel /= mpss; - // update ddtUrel and r history - update_ddtUrelHist(ddtUrelHist_,ddtUrel,index); // add current dim.less ddtUrel to history - update_rHist(rHist_,r,index); // add current r to history - // check length of known history - int nKnown = 0; + int nKnown = 1; // we always know the current step for (int j=0; j0) - { - for (int i=0; i<3; i++) // loop over dimensions - Fshort[i] = -calculateK0(r,dt) * ddtUrelHist_[index][i*nHist_]; - } + for (int i=0; i<3; i++) // loop over dimensions + Fshort[i] = -calculateK0(r,dt) * ddtUrel[i]; // int_dt^t0 K(r,xi) * ddtU(t-xi) dxi (trapezoid rule) if (nShort>2) { - for (int j=1; j Date: Mon, 24 Oct 2022 11:47:51 +0200 Subject: [PATCH 12/23] remove unneeded component of rHist to save memory --- .../ParmarBassetForce/ParmarBassetForce.C | 53 ++++++++++--------- .../ParmarBassetForce/ParmarBassetForce.H | 6 ++- 2 files changed, 31 insertions(+), 28 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index 17052c57..2446796b 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -72,9 +72,10 @@ ParmarBassetForce::ParmarBassetForce Us_(sm.mesh().lookupObject (UsFieldName_)), nInt_(readLabel(propsDict_.lookup("nIntegral"))), discOrder_(readLabel(propsDict_.lookup("discretisationOrder"))), - nHist_(nInt_+discOrder_), + ddtUrelHistSize_(nInt_+discOrder_), + rHistSize_(nInt_), FHistSize_(2*discOrder_), - ddtUrelHistRegName_(typeName + "ddtUrelHist"), // indexed as: ddtUrelHist_[particle ID][iDim*nHist_+iHist] + ddtUrelHistRegName_(typeName + "ddtUrelHist"), // indexed as: ddtUrelHist_[particle ID][iDim*ddtUrelHistSize_+iHist] rHistRegName_(typeName + "rHist"), // indexed as: rHist_[particle ID][iHist] FHistRegName_(typeName + "FHist"), // indexed as: ddtUrelHist_[particle ID][iDim*FHistSize_*2+iHist*2+k] gH0RegName_(typeName + "gH0"), @@ -116,13 +117,13 @@ ParmarBassetForce::ParmarBassetForce { // allocate particle properties - particleCloud_.registerParticleProperty(ddtUrelHistRegName_, 3*nHist_,NOTONCPU,false); - particleCloud_.registerParticleProperty( rHistRegName_, nHist_,NOTONCPU,false); - particleCloud_.registerParticleProperty( FHistRegName_,6*FHistSize_,NOTONCPU,false); - particleCloud_.registerParticleProperty( gH0RegName_, 1 ,NOTONCPU,false); - particleCloud_.registerParticleProperty( tRefRegName_, 1 ,NOTONCPU,false); - particleCloud_.registerParticleProperty( mRefRegName_, 1 ,NOTONCPU,false); - particleCloud_.registerParticleProperty( lRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty(ddtUrelHistRegName_,3*ddtUrelHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( rHistRegName_, rHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( FHistRegName_, 6*FHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( gH0RegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( tRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( mRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( lRefRegName_, 1 ,NOTONCPU,false); // init force sub model setForceSubModels(propsDict_); @@ -288,7 +289,7 @@ void ParmarBassetForce::setForce() const // check length of known history int nKnown = 1; // we always know the current step - for (int j=0; j0; j--) // loop over past times - ddtUrelHist_[index][i*nHist_+j] = ddtUrelHist_[index][i*nHist_+j-1]; + for (int j=ddtUrelHistSize_-1; j>0; j--) // loop over past times + ddtUrelHist_[index][i*ddtUrelHistSize_+j] = ddtUrelHist_[index][i*ddtUrelHistSize_+j-1]; - ddtUrelHist_[index][i*nHist_] = ddtUrel[i]; + ddtUrelHist_[index][i*ddtUrelHistSize_] = ddtUrel[i]; } } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // void Foam::ParmarBassetForce::update_rHist(double**& rHist_, scalar r, int index) const { - for (int j=nHist_-1; j>0; j--) // loop over past times + for (int j=rHistSize_-1; j>0; j--) // loop over past times rHist_[index][j] = rHist_[index][j-1]; rHist_[index][0] = r; @@ -499,8 +500,8 @@ vector Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUre ( ( C[1]/dt+2/(dt*dt)) * FHist_[index][i*FHistSize_*2 +k] - ( 1/(dt*dt)) * FHist_[index][i*FHistSize_*2+2+k] - - (C[2]+C[3]/dt ) * ddtUrelHist_[index][i*nHist_+nInt_-1] - + ( C[3]/dt ) * ddtUrelHist_[index][i*nHist_+nInt_ ] + - (C[2]+C[3]/dt ) * ddtUrelHist_[index][i*ddtUrelHistSize_+nInt_-1] + + ( C[3]/dt ) * ddtUrelHist_[index][i*ddtUrelHistSize_+nInt_ ] ) / (C[0]+C[1]/dt+1/(dt*dt)); // Eq. 3.20 using first order temporal discretisation } } @@ -515,9 +516,9 @@ vector Foam::ParmarBassetForce::solveFlongODE(double**& FHist_, double**& ddtUre + ( 8/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+4+k] - ( 1/(4*dt*dt)) * FHist_[index][i*FHistSize_*2+6+k] - - (C[2] + 3*C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_-1] - + ( 4*C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_ ] - - ( C[3]/(2*dt) ) * ddtUrelHist_[index][i*nHist_+nInt_+1] + - (C[2] + 3*C[3]/(2*dt) ) * ddtUrelHist_[index][i*ddtUrelHistSize_+nInt_-1] + + ( 4*C[3]/(2*dt) ) * ddtUrelHist_[index][i*ddtUrelHistSize_+nInt_ ] + - ( C[3]/(2*dt) ) * ddtUrelHist_[index][i*ddtUrelHistSize_+nInt_+1] ) / (C[0] + 3*C[1]/(2*dt) + 9/(4*dt*dt)); // Eq. 3.20 using second order temporal discretisation } } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H index df6e20a8..1d89695d 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H @@ -71,9 +71,11 @@ private: label discOrder_; //ODE discretisation order - label nHist_; //no of timesteps to save ddtUrel history for + label ddtUrelHistSize_; //no of timesteps to save ddtUrel history for - label FHistSize_; + label rHistSize_; //no of timesteps to save r history for + + label FHistSize_; //no of timesteps to save Flong history for const word ddtUrelHistRegName_; From 2aa5c7880b26b329042951bc5674480244b4eb6c Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Mon, 24 Oct 2022 12:21:29 +0200 Subject: [PATCH 13/23] rework input --- doc/forceModel_ParmarBassetForce.txt | 6 ++--- .../ParmarBassetForce/ParmarBassetForce.C | 26 ++++++++++--------- .../ParmarBassetForce/ParmarBassetForce.H | 14 +++++----- 3 files changed, 24 insertions(+), 22 deletions(-) diff --git a/doc/forceModel_ParmarBassetForce.txt b/doc/forceModel_ParmarBassetForce.txt index bc022f53..a5e41730 100644 --- a/doc/forceModel_ParmarBassetForce.txt +++ b/doc/forceModel_ParmarBassetForce.txt @@ -30,12 +30,12 @@ ParmarBassetForceProps {U} = name of the finite volume fluid velocity field :ulb,l {Us} = name of the finite volume cell averaged particle velocity field :l -{scalar1} = number of timesteps considered in the near history :l -{scalar2} = (1 or 2) discretisation order of the far history differential equations :l +{scalar1} = number of timesteps considered in the near history (integer > 1):l +{scalar2} = (optional, default 1) discretisation order of the far history differential equations (1 or 2) :l {switch1} = (optional, default true) sub model switch, see "forceSubModel"_forceSubModel.html for details :l {switch2} = (optional, default false) sub model switch, see "forceSubModel"_forceSubModel.html for details :l {switch3} = (optional, default false) sub model switch, see "forceSubModel"_forceSubModel.html for details :l -{smoothingModel} = name of smoothing model for the dU/dt field :l +{smoothingModel} = name of smoothing model for the dUrel/dt field :l :ule [Examples:] diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index 2446796b..6ccc631b 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -70,8 +70,8 @@ ParmarBassetForce::ParmarBassetForce U_(sm.mesh().lookupObject (velFieldName_)), UsFieldName_(propsDict_.lookup("granVelFieldName")), Us_(sm.mesh().lookupObject (UsFieldName_)), - nInt_(readLabel(propsDict_.lookup("nIntegral"))), - discOrder_(readLabel(propsDict_.lookup("discretisationOrder"))), + nInt_(propsDict_.lookupOrDefault("nIntegral", -1)), + discOrder_(propsDict_.lookupOrDefault("discretisationOrder", 1)), ddtUrelHistSize_(nInt_+discOrder_), rHistSize_(nInt_), FHistSize_(2*discOrder_), @@ -115,16 +115,6 @@ ParmarBassetForce::ParmarBassetForce ) ) { - - // allocate particle properties - particleCloud_.registerParticleProperty(ddtUrelHistRegName_,3*ddtUrelHistSize_,NOTONCPU,false); - particleCloud_.registerParticleProperty( rHistRegName_, rHistSize_,NOTONCPU,false); - particleCloud_.registerParticleProperty( FHistRegName_, 6*FHistSize_,NOTONCPU,false); - particleCloud_.registerParticleProperty( gH0RegName_, 1 ,NOTONCPU,false); - particleCloud_.registerParticleProperty( tRefRegName_, 1 ,NOTONCPU,false); - particleCloud_.registerParticleProperty( mRefRegName_, 1 ,NOTONCPU,false); - particleCloud_.registerParticleProperty( lRefRegName_, 1 ,NOTONCPU,false); - // init force sub model setForceSubModels(propsDict_); // define switches which can be read from dict @@ -139,6 +129,18 @@ ParmarBassetForce::ParmarBassetForce if (discOrder_ < 1 || discOrder_ > 2) FatalError << "Parmar Basset Force: Discretisation order > 2 not implemented!" << abort(FatalError); + if (nInt_ < 1) + FatalError << "Parmar Basset Force: nIntegral missing or invalid, must be > 1" << abort(FatalError); + + // allocate particle properties + particleCloud_.registerParticleProperty(ddtUrelHistRegName_,3*ddtUrelHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( rHistRegName_, rHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( FHistRegName_, 6*FHistSize_,NOTONCPU,false); + particleCloud_.registerParticleProperty( gH0RegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( tRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( mRefRegName_, 1 ,NOTONCPU,false); + particleCloud_.registerParticleProperty( lRefRegName_, 1 ,NOTONCPU,false); + //Append the field names to be probed particleCloud_.probeM().initialize(typeName, typeName+".logDat"); particleCloud_.probeM().vectorFields_.append("ParmarBassetForce"); //first entry must the be the force diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H index 1d89695d..122e9172 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.H @@ -59,23 +59,23 @@ class ParmarBassetForce private: dictionary propsDict_; - word velFieldName_; + const word velFieldName_; const volVectorField& U_; - word UsFieldName_; + const word UsFieldName_; const volVectorField& Us_; - label nInt_; //no of timesteps to solve full integral + const int nInt_; //no of timesteps to solve full integral - label discOrder_; //ODE discretisation order + const int discOrder_; //ODE discretisation order - label ddtUrelHistSize_; //no of timesteps to save ddtUrel history for + const int ddtUrelHistSize_; //no of timesteps to save ddtUrel history for - label rHistSize_; //no of timesteps to save r history for + const int rHistSize_; //no of timesteps to save r history for - label FHistSize_; //no of timesteps to save Flong history for + const int FHistSize_; //no of timesteps to save Flong history for const word ddtUrelHistRegName_; From 20e75cf64fc44d125d8a247c91f11521a4da7a5c Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Mon, 24 Oct 2022 12:31:19 +0200 Subject: [PATCH 14/23] eliminate pow() for performance --- .../forceModel/ParmarBassetForce/ParmarBassetForce.C | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C index 6ccc631b..5dcb7e39 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/ParmarBassetForce/ParmarBassetForce.C @@ -255,7 +255,8 @@ void ParmarBassetForce::setForce() const if (gH0_[index][0]!=NOTONCPU) { - r = pow(gH0_[index][0]/gH,1.5); // Eq. 3.4 + scalar gHratio = gH0_[index][0]/gH; + r = gHratio*sqrt(gHratio); // gHratio^1.5, Eq. 3.4 if (r<0.25 || r>2.0) { @@ -314,7 +315,8 @@ void ParmarBassetForce::setForce() const for (int j=0; j<(nShort-1); j++) // we don't use the current step here, hence nShort-1 { scalar xi = (j+1)*dt; - scalar K = pow((pow(xi,.25) + rHist_[index][j]*xi),-2.); // Eq. 3.4 + scalar invsqrtK = sqrt(sqrt(xi)) + rHist_[index][j]*xi; // K^-0.5 + scalar K = 1./(invsqrtK*invsqrtK); // Eq. 3.4 for (int i=0; i<3; i++) // loop over dimensions Fshort[i] -= trapWeight(j,nShort-1) * K * ddtUrelHist_[index][i*ddtUrelHistSize_+j] * dt; @@ -416,7 +418,7 @@ void ParmarBassetForce::setForce() const scalar Foam::ParmarBassetForce::calculateK0(scalar r, scalar dt) const { scalar cbrtr = cbrt(r); // cube root of r - scalar gamma = cbrtr*pow(dt,0.25); + scalar gamma = cbrtr*sqrt(sqrt(dt)); /* scalar K0 = 2./(9.*pow(r,0.666)) * From 4ddf691936fe0bd7173ea2ba6415bea203184b33 Mon Sep 17 00:00:00 2001 From: tmjnijssen Date: Thu, 14 Apr 2022 15:29:10 +0200 Subject: [PATCH 15/23] add transfer of fluid properties to LIGGGHTS --- src/lagrangian/cfdemParticle/Make/files | 1 + .../transferFluidProperties.C | 108 ++++++++++++++++++ .../transferFluidProperties.H | 80 +++++++++++++ src/lagrangian/cfdemParticleComp/Make/files | 1 + 4 files changed, 190 insertions(+) create mode 100644 src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C create mode 100644 src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H diff --git a/src/lagrangian/cfdemParticle/Make/files b/src/lagrangian/cfdemParticle/Make/files index 00d33dbd..b0f4449f 100644 --- a/src/lagrangian/cfdemParticle/Make/files +++ b/src/lagrangian/cfdemParticle/Make/files @@ -96,6 +96,7 @@ $(forceModels)/potentialRelaxation/potentialRelaxation.C $(forceModels)/BeetstraDrag/BeetstraDrag.C $(forceModels)/BeetstraDragPoly/BeetstraDragPoly.C $(forceModels)/dSauter/dSauter.C +$(forceModels)/transferFluidProperties/transferFluidProperties.C $(forceModels)/Fines/Fines.C $(forceModels)/Fines/FinesFields.C $(forceModels)/Fines/FanningDynFines.C diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C new file mode 100644 index 00000000..5faf375d --- /dev/null +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C @@ -0,0 +1,108 @@ +/*---------------------------------------------------------------------------*\ +License + This 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. + This code 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 + along with this code. If not, see . + + Copyright (C) 2015- Thomas Lichtenegger, JKU Linz, Austria + +Description + transfer fluid properties to LIGGGHTS + +SourceFiles + transferFluidProperties.C +\*---------------------------------------------------------------------------*/ + +#include "error.H" + +#include "transferFluidProperties.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +defineTypeNameAndDebug(transferFluidProperties, 0); + +addToRunTimeSelectionTable +( + forceModel, + transferFluidProperties, + dictionary +); + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +// Construct from components +transferFluidProperties::transferFluidProperties +( + const dictionary& dict, + cfdemCloud& sm +) +: + forceModel(dict,sm), + propsDict_(dict.subDict(typeName + "Props")), + verbose_(propsDict_.lookupOrDefault("verbose",false)) +{ + particleCloud_.registerParticleProperty("fluidDensity",1); + particleCloud_.registerParticleProperty("fluidViscosity",1); + + // init force sub model + setForceSubModels(propsDict_); +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +transferFluidProperties::~transferFluidProperties() +{ +} + +// * * * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * public Member Functions * * * * * * * * * * * * * // + +void transferFluidProperties::setForce() const +{ + double**& fluidDensity_ = particleCloud_.getParticlePropertyRef("fluidDensity"); + double**& fluidViscosity_ = particleCloud_.getParticlePropertyRef("fluidViscosity"); + + const volScalarField& rhoField = forceSubM(0).rhoField(); + const volScalarField& nufField = forceSubM(0).nuField(); + + label cellI = 0; + + for(int index = 0; index < particleCloud_.numberOfParticles(); ++index) + { + cellI = particleCloud_.cellIDs()[index][0]; + if (cellI >= 0) + { + fluidDensity_[index][0] = rhoField[cellI]; + fluidViscosity_[index][0] = nufField[cellI] * rhoField[cellI]; + } + } + + particleCloud_.dataExchangeM().giveData("fluidDensity","scalar-atom",fluidDensity_); + particleCloud_.dataExchangeM().giveData("fluidViscosity","scalar-atom",fluidViscosity_); + + if (verbose_) Info << "give data done" << endl; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H new file mode 100644 index 00000000..9a2853ae --- /dev/null +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H @@ -0,0 +1,80 @@ +/*---------------------------------------------------------------------------*\ +License + This 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. + This code 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 + along with this code. If not, see . + + Copyright (C) 2015- Thomas Lichtenegger, JKU Linz, Austria + +Description + transfer fluid properties to LIGGGHTS + +SourceFiles + transferFluidProperties.C +\*---------------------------------------------------------------------------*/ + +#ifndef transferFluidProperties_H +#define transferFluidProperties_H + +#include "forceModel.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class transferFluidProperties Declaration +\*---------------------------------------------------------------------------*/ + +class transferFluidProperties +: + public forceModel +{ +private: + + dictionary propsDict_; + + bool verbose_; + +public: + + //- Runtime type information + TypeName("transferFluidProperties"); + + // Constructors + + //- Construct from components + transferFluidProperties + ( + const dictionary& dict, + cfdemCloud& sm + ); + + // Destructor + + ~transferFluidProperties(); + + + // Member Functions + void setForce() const; + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/lagrangian/cfdemParticleComp/Make/files b/src/lagrangian/cfdemParticleComp/Make/files index 2907a3ea..dfccc6e6 100644 --- a/src/lagrangian/cfdemParticleComp/Make/files +++ b/src/lagrangian/cfdemParticleComp/Make/files @@ -90,6 +90,7 @@ $(forceModels)/directedDiffusiveRelaxation/directedDiffusiveRelaxation.C $(forceModels)/BeetstraDrag/BeetstraDrag.C $(forceModels)/BeetstraDragPoly/BeetstraDragPoly.C $(forceModels)/dSauter/dSauter.C +$(forceModels)/transferFluidProperties/transferFluidProperties.C $(forceModels)/Fines/Fines.C $(forceModels)/Fines/FinesFields.C $(forceModels)/Fines/FanningDynFines.C From d6dab59bfde879c7e3f75db4de1e6d06a8416322 Mon Sep 17 00:00:00 2001 From: tmjnijssen Date: Thu, 14 Apr 2022 15:50:22 +0200 Subject: [PATCH 16/23] add interpolation --- .../transferFluidProperties.C | 31 ++++++++++++++++--- .../transferFluidProperties.H | 3 +- 2 files changed, 27 insertions(+), 7 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C index 5faf375d..5c2695d5 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.C @@ -52,14 +52,17 @@ transferFluidProperties::transferFluidProperties ) : forceModel(dict,sm), - propsDict_(dict.subDict(typeName + "Props")), - verbose_(propsDict_.lookupOrDefault("verbose",false)) + propsDict_(dict.subDict(typeName + "Props")) { particleCloud_.registerParticleProperty("fluidDensity",1); particleCloud_.registerParticleProperty("fluidViscosity",1); // init force sub model setForceSubModels(propsDict_); + // define switches which can be read from dict + forceSubM(0).setSwitchesList(SW_VERBOSE,true); // activate search for verbose switch + forceSubM(0).setSwitchesList(SW_INTERPOLATION,true); // activate search for interpolate switch + forceSubM(0).readSwitches(); } @@ -82,22 +85,40 @@ void transferFluidProperties::setForce() const const volScalarField& rhoField = forceSubM(0).rhoField(); const volScalarField& nufField = forceSubM(0).nuField(); + interpolationCellPoint rhoInterpolator_(rhoField); + interpolationCellPoint nufInterpolator_(nufField); + label cellI = 0; + double rho = 0.; + double nuf = 0.; + vector position(0,0,0); for(int index = 0; index < particleCloud_.numberOfParticles(); ++index) { cellI = particleCloud_.cellIDs()[index][0]; if (cellI >= 0) { - fluidDensity_[index][0] = rhoField[cellI]; - fluidViscosity_[index][0] = nufField[cellI] * rhoField[cellI]; + if(forceSubM(0).interpolation()) + { + position = particleCloud_.position(index); + rho = rhoInterpolator_.interpolate(position,cellI); + nuf = nufInterpolator_.interpolate(position,cellI); + } + else + { + rho = rhoField[cellI]; + nuf = nufField[cellI]; + } + + fluidDensity_[index][0] = rho; + fluidViscosity_[index][0] = nuf*rho; } } particleCloud_.dataExchangeM().giveData("fluidDensity","scalar-atom",fluidDensity_); particleCloud_.dataExchangeM().giveData("fluidViscosity","scalar-atom",fluidViscosity_); - if (verbose_) Info << "give data done" << endl; + if (forceSubM(0).verbose()) Info << "give data done" << endl; } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H index 9a2853ae..1dc52f85 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/transferFluidProperties/transferFluidProperties.H @@ -24,6 +24,7 @@ SourceFiles #define transferFluidProperties_H #include "forceModel.H" +#include "interpolationCellPoint.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -42,8 +43,6 @@ private: dictionary propsDict_; - bool verbose_; - public: //- Runtime type information From 57c8c1c7623ea061801c0c83e21ac774840424a8 Mon Sep 17 00:00:00 2001 From: tmjnijssen Date: Thu, 14 Apr 2022 15:50:39 +0200 Subject: [PATCH 17/23] add doc --- doc/CFDEMcoupling_models.txt | 1 + doc/forceModel_transferFluidProperties.txt | 42 ++++++++++++++++++++++ 2 files changed, 43 insertions(+) create mode 100644 doc/forceModel_transferFluidProperties.txt diff --git a/doc/CFDEMcoupling_models.txt b/doc/CFDEMcoupling_models.txt index 00ce8e28..073082dd 100644 --- a/doc/CFDEMcoupling_models.txt +++ b/doc/CFDEMcoupling_models.txt @@ -124,6 +124,7 @@ particleDeformation, potentialRelaxation, "surfaceTensionForce"_forceModel_surfaceTensionForce.html, terminalVelocity, +"transferFluidProperties"_forceModel_transferFluidProperties.html, turbulentDispersion, turbulentVelocityFluctuations, "virtualMassForce"_forceModel_virtualMassForce.html, diff --git a/doc/forceModel_transferFluidProperties.txt b/doc/forceModel_transferFluidProperties.txt new file mode 100644 index 00000000..8606da78 --- /dev/null +++ b/doc/forceModel_transferFluidProperties.txt @@ -0,0 +1,42 @@ +"CFDEMproject Website"_lws - "Main Page"_main :c + +:link(lws,http://www.cfdem.com) +:link(main,CFDEMcoupling_Manual.html) + +:line + +forceModel transferFluidProperties command :h3 + +[Syntax:] + +Defined in "couplingProperties"_CFDEMcoupling_dicts.html#couplingProperties +dictionary. + +forceModels +( + transferFluidProperties +); +transferFluidPropertiesProps +\{ + verbose switch1; + interpolation switch2; +\} :pre + +{switch1} = (optional, default false) sub model switch, see "forceSubModel"_forceSubModel.html for details :ulb,l +{switch2} = (optional, default false) sub model switch, see "forceSubModel"_forceSubModel.html for details :l +:ule + +[Description:] + +This "force model" does not influence the particles or the flow - it transfer to fluid density and (dynamic) +viscosity from OpenFOAM to LIGGGHTS. + + +[Restrictions:] + +This model requires {fix cfd/coupling/fluidproperties} to work. + +[Related commands:] + +"forceModel"_forceModel.html + From 8ba87ddedb181f2a47ebbddb6ca93e0f6de1d554 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Wed, 26 Oct 2022 16:01:58 +0200 Subject: [PATCH 18/23] probe bugfixes --- .../subModels/forceModel/BeetstraDrag/BeetstraDrag.C | 1 - .../subModels/probeModel/particleProbe/particleProbe.C | 2 +- 2 files changed, 1 insertion(+), 2 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C index b52f2a7b..345b27c7 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C @@ -77,7 +77,6 @@ BeetstraDrag::BeetstraDrag particleCloud_.probeM().vectorFields_.append("dragForce"); //first entry must be the force particleCloud_.probeM().vectorFields_.append("Urel"); particleCloud_.probeM().scalarFields_.append("Rep"); - particleCloud_.probeM().scalarFields_.append("betaP"); particleCloud_.probeM().scalarFields_.append("voidfraction"); particleCloud_.probeM().writeHeader(); diff --git a/src/lagrangian/cfdemParticle/subModels/probeModel/particleProbe/particleProbe.C b/src/lagrangian/cfdemParticle/subModels/probeModel/particleProbe/particleProbe.C index 5a0441d9..6060b38b 100644 --- a/src/lagrangian/cfdemParticle/subModels/probeModel/particleProbe/particleProbe.C +++ b/src/lagrangian/cfdemParticle/subModels/probeModel/particleProbe/particleProbe.C @@ -219,7 +219,7 @@ void particleProbe::writeProbe(int index, Field sValues, Field v *sPtr << setprecision(writePrecision_); forAll(vValues, iter) { - // if(!probeDebug_ && iter>0) break; + if(!probeDebug_ && iter>0) break; *sPtr << vValues[iter][0] << " "; *sPtr << vValues[iter][1] << " "; *sPtr << vValues[iter][2] << " "; From 3b7724003a1da69e9c40f2f2c638f41602f58b05 Mon Sep 17 00:00:00 2001 From: Tim MJ Nijssen <37873965+tmjnijssen@users.noreply.github.com> Date: Wed, 2 Nov 2022 16:32:00 +0100 Subject: [PATCH 19/23] Update smoothingModel_constDiffSmoothing.txt --- doc/smoothingModel_constDiffSmoothing.txt | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/doc/smoothingModel_constDiffSmoothing.txt b/doc/smoothingModel_constDiffSmoothing.txt index 34b0c3e4..99c06c58 100644 --- a/doc/smoothingModel_constDiffSmoothing.txt +++ b/doc/smoothingModel_constDiffSmoothing.txt @@ -18,7 +18,7 @@ constDiffSmoothingProps lowerLimit number1; upperLimit number2; smoothingLength lengthScale; - smoothingLengthReference lengthScaleRefField; + smoothingLengthReference lengthScaleRef; smoothingLengthFieldName fieldName1; smoothingLengthReferenceFieldName fieldName2; verbose; @@ -27,7 +27,7 @@ constDiffSmoothingProps {number1} = scalar fields will be bound to this lower value :ulb,l {number2} = scalar fields will be bound to this upper value :l {lengthScale} = length scale over which the exchange fields will be smoothed out :l -{lengthScaleRefField} = (optional) length scale over which reference fields (e.g., the average particle velocity) will be smoothed out. Should be always larger than lengthScale. If not specified, will be equal to lengthScale. :l +{lengthScaleRef} = (optional) length scale over which reference fields (e.g., the average particle velocity) will be smoothed out. Should be always larger than lengthScale. If not specified, will be equal to lengthScale. :l {fieldName1} = (optional) name of scalar field to be used as local smoothing length. :l {fieldName2} = (optional) name of scalar field to be used as local smoothing length for reference fields. :l @@ -56,7 +56,7 @@ which these reference fields are not specified. Values calculated in the cells e.g. the average particle velocity, which are not specified in all cells in case the flow is rather dilute. -Alternative to {smoothingLength} and {smoothingLengthReferenceField}, +Alternative to {smoothingLength} and {smoothingLengthReference}, {smoothingLengthFieldName} and/or {smoothingLengthReferenceFieldName} can be used to define spatial variation of the smoothing lengths. Either the scalar or field options must be used, giving both will result in errors. From 41b516c1e8cadd61a3e8b379333ec12d73abf40a Mon Sep 17 00:00:00 2001 From: Thomas Lichtenegger Date: Tue, 8 Nov 2022 17:41:17 +0100 Subject: [PATCH 20/23] Check out fields manually checked into registry. --- .../standardRecModel/standardRecModel.C | 19 ++++++++++++++++++- 1 file changed, 18 insertions(+), 1 deletion(-) diff --git a/src/recurrence/recModel/standardRecModel/standardRecModel.C b/src/recurrence/recModel/standardRecModel/standardRecModel.C index 36407588..d8b5cf00 100644 --- a/src/recurrence/recModel/standardRecModel/standardRecModel.C +++ b/src/recurrence/recModel/standardRecModel/standardRecModel.C @@ -159,7 +159,24 @@ standardRecModel::standardRecModel // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // standardRecModel::~standardRecModel() -{} +{ + const objectRegistry& objReg = base_.mesh().thisDb(); + + for(int j=0; j Date: Thu, 17 Nov 2022 14:35:17 +0100 Subject: [PATCH 21/23] rename cfdemSolverMultiphaseScalar/multiphaseMixture to cfdemSolverMultiphaseScalar/multiphaseMixtureScalar --- .../cfdemSolverMultiphaseScalar/Allwclean | 2 +- .../cfdemSolverMultiphaseScalar/Allwmake | 2 +- .../cfdemSolverMultiphaseScalar/Make/options | 2 +- .../cfdemSolverMultiphaseScalar.C | 2 +- .../createFields.H | 2 +- .../Make/files | 2 +- .../Make/options | 0 .../alphaContactAngleFvPatchScalarField.C | 0 .../alphaContactAngleFvPatchScalarField.H | 8 +-- .../multiphaseMixtureScalar.C} | 58 +++++++++---------- .../multiphaseMixtureScalar.H} | 18 +++--- .../phase/phase.C | 0 .../phase/phase.H | 2 +- etc/library-list.txt | 2 +- 14 files changed, 50 insertions(+), 50 deletions(-) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture => multiphaseMixtureScalar}/Make/files (83%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture => multiphaseMixtureScalar}/Make/options (100%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture => multiphaseMixtureScalar}/alphaContactAngle/alphaContactAngleFvPatchScalarField.C (100%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture => multiphaseMixtureScalar}/alphaContactAngle/alphaContactAngleFvPatchScalarField.H (96%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture/multiphaseMixture.C => multiphaseMixtureScalar/multiphaseMixtureScalar.C} (93%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture/multiphaseMixture.H => multiphaseMixtureScalar/multiphaseMixtureScalar.H} (96%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture => multiphaseMixtureScalar}/phase/phase.C (100%) rename applications/solvers/cfdemSolverMultiphaseScalar/{multiphaseMixture => multiphaseMixtureScalar}/phase/phase.H (98%) diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/Allwclean b/applications/solvers/cfdemSolverMultiphaseScalar/Allwclean index 1c8a5a68..8f8e3fe7 100755 --- a/applications/solvers/cfdemSolverMultiphaseScalar/Allwclean +++ b/applications/solvers/cfdemSolverMultiphaseScalar/Allwclean @@ -2,7 +2,7 @@ cd ${0%/*} || exit 1 # Run from this directory set -x -wclean libso multiphaseMixture +wclean libso multiphaseMixtureScalar wclean #------------------------------------------------------------------------------ diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/Allwmake b/applications/solvers/cfdemSolverMultiphaseScalar/Allwmake index fbe71a59..2e60481b 100755 --- a/applications/solvers/cfdemSolverMultiphaseScalar/Allwmake +++ b/applications/solvers/cfdemSolverMultiphaseScalar/Allwmake @@ -6,7 +6,7 @@ targetType=libso . $WM_PROJECT_DIR/wmake/scripts/AllwmakeParseArguments set -x -wmake $targetType multiphaseMixture +wmake $targetType multiphaseMixtureScalar wmake #------------------------------------------------------------------------------ diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/Make/options b/applications/solvers/cfdemSolverMultiphaseScalar/Make/options index 72a134ad..ce909f3b 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/Make/options +++ b/applications/solvers/cfdemSolverMultiphaseScalar/Make/options @@ -6,7 +6,7 @@ include $(CFDEM_ADD_LIBS_DIR)/additionalLibs EXE_INC = \ $(PFLAGS) \ -I$(CFDEM_OFVERSION_DIR) \ - -ImultiphaseMixture/lnInclude \ + -ImultiphaseMixtureScalar/lnInclude \ -I$(LIB_SRC)/transportModels \ -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \ diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/cfdemSolverMultiphaseScalar.C b/applications/solvers/cfdemSolverMultiphaseScalar/cfdemSolverMultiphaseScalar.C index 504450fc..0ce0ae16 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/cfdemSolverMultiphaseScalar.C +++ b/applications/solvers/cfdemSolverMultiphaseScalar/cfdemSolverMultiphaseScalar.C @@ -30,7 +30,7 @@ Description \*---------------------------------------------------------------------------*/ #include "fvCFD.H" -#include "multiphaseMixture.H" +#include "multiphaseMixtureScalar.H" #include "turbulentTransportModel.H" #include "pimpleControl.H" #include "fvOptions.H" diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/createFields.H b/applications/solvers/cfdemSolverMultiphaseScalar/createFields.H index 009bd16d..d3b65ec2 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/createFields.H +++ b/applications/solvers/cfdemSolverMultiphaseScalar/createFields.H @@ -88,7 +88,7 @@ surfaceScalarField phi linearInterpolate(U*voidfraction) & mesh.Sf() ); -multiphaseMixture mixture(U, phi, voidfraction); +multiphaseMixtureScalar mixture(U, phi, voidfraction); // Need to store rho for ddt(rho, U) volScalarField rho diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/Make/files b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/Make/files similarity index 83% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/Make/files rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/Make/files index 998255e4..91bbae6b 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/Make/files +++ b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/Make/files @@ -1,5 +1,5 @@ phase/phase.C alphaContactAngle/alphaContactAngleFvPatchScalarField.C -multiphaseMixture.C +multiphaseMixtureScalar.C LIB = $(CFDEM_LIB_DIR)/libcfdemMultiphaseInterFoamScalar diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/Make/options b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/Make/options similarity index 100% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/Make/options rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/Make/options diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.C b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/alphaContactAngle/alphaContactAngleFvPatchScalarField.C similarity index 100% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.C rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/alphaContactAngle/alphaContactAngleFvPatchScalarField.C diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.H b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/alphaContactAngle/alphaContactAngleFvPatchScalarField.H similarity index 96% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.H rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/alphaContactAngle/alphaContactAngleFvPatchScalarField.H index 09249c72..bbad25b3 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.H +++ b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/alphaContactAngle/alphaContactAngleFvPatchScalarField.H @@ -26,7 +26,7 @@ Class Description Contact-angle boundary condition for multi-phase interface-capturing - simulations. Used in conjuction with multiphaseMixture. + simulations. Used in conjuction with multiphaseMixtureScalar. SourceFiles alphaContactAngleFvPatchScalarField.C @@ -37,7 +37,7 @@ SourceFiles #define alphaContactAngleFvPatchScalarField_H #include "zeroGradientFvPatchFields.H" -#include "multiphaseMixture.H" +#include "multiphaseMixtureScalar.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -117,8 +117,8 @@ public: typedef HashTable < interfaceThetaProps, - multiphaseMixture::interfacePair, - multiphaseMixture::interfacePair::hash + multiphaseMixtureScalar::interfacePair, + multiphaseMixtureScalar::interfacePair::hash > thetaPropsTable; diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/multiphaseMixture.C b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/multiphaseMixtureScalar.C similarity index 93% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/multiphaseMixture.C rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/multiphaseMixtureScalar.C index 9d6cc130..0445368f 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/multiphaseMixture.C +++ b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/multiphaseMixtureScalar.C @@ -18,7 +18,7 @@ License \*---------------------------------------------------------------------------*/ -#include "multiphaseMixture.H" +#include "multiphaseMixtureScalar.H" #include "alphaContactAngleFvPatchScalarField.H" #include "Time.H" #include "subCycle.H" @@ -31,13 +31,13 @@ License // * * * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * // -const Foam::scalar Foam::multiphaseMixture::convertToRad = +const Foam::scalar Foam::multiphaseMixtureScalar::convertToRad = Foam::constant::mathematical::pi/180.0; // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -void Foam::multiphaseMixture::calcAlphas() +void Foam::multiphaseMixtureScalar::calcAlphas() { scalar level = 0.0; alphas_ == 0.0; @@ -51,7 +51,7 @@ void Foam::multiphaseMixture::calcAlphas() Foam::tmp -Foam::multiphaseMixture::calcNu() const +Foam::multiphaseMixtureScalar::calcNu() const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -74,7 +74,7 @@ Foam::multiphaseMixture::calcNu() const } Foam::tmp -Foam::multiphaseMixture::calcStf() const +Foam::multiphaseMixtureScalar::calcStf() const { tmp tstf ( @@ -134,7 +134,7 @@ Foam::multiphaseMixture::calcStf() const // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // -Foam::multiphaseMixture::multiphaseMixture +Foam::multiphaseMixtureScalar::multiphaseMixtureScalar ( const volVectorField& U, const surfaceScalarField& phi, @@ -230,7 +230,7 @@ Foam::multiphaseMixture::multiphaseMixture // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // Foam::tmp -Foam::multiphaseMixture::rho() const +Foam::multiphaseMixtureScalar::rho() const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -247,7 +247,7 @@ Foam::multiphaseMixture::rho() const Foam::tmp -Foam::multiphaseMixture::rho(const label patchi) const +Foam::multiphaseMixtureScalar::rho(const label patchi) const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -264,9 +264,9 @@ Foam::multiphaseMixture::rho(const label patchi) const Foam::tmp -Foam::multiphaseMixture::mu() const +Foam::multiphaseMixtureScalar::mu() const { - Info << "In multiphasemixture mu()" << endl; + Info << "In multiphaseMixtureScalar mu()" << endl; return rho()*nu(); // PtrDictionary::const_iterator iter = phases_.begin(); @@ -283,7 +283,7 @@ Foam::multiphaseMixture::mu() const Foam::tmp -Foam::multiphaseMixture::mu(const label patchi) const +Foam::multiphaseMixtureScalar::mu(const label patchi) const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -306,7 +306,7 @@ Foam::multiphaseMixture::mu(const label patchi) const Foam::tmp -Foam::multiphaseMixture::muf() const +Foam::multiphaseMixtureScalar::muf() const { return nuf()*fvc::interpolate(rho()); @@ -327,13 +327,13 @@ Foam::multiphaseMixture::muf() const Foam::tmp -Foam::multiphaseMixture::nu() const +Foam::multiphaseMixtureScalar::nu() const { return nu_; } Foam::tmp -Foam::multiphaseMixture::nu(const label patchi) const +Foam::multiphaseMixtureScalar::nu(const label patchi) const { //return nu_.boundaryField()[patchi]; PtrDictionary::const_iterator iter = phases_.begin(); @@ -355,7 +355,7 @@ Foam::multiphaseMixture::nu(const label patchi) const Foam::tmp -Foam::multiphaseMixture::nuf() const +Foam::multiphaseMixtureScalar::nuf() const { //return muf()/fvc::interpolate(rho()); PtrDictionary::const_iterator iter = phases_.begin(); @@ -374,7 +374,7 @@ Foam::multiphaseMixture::nuf() const } Foam::tmp -Foam::multiphaseMixture::Cp() const +Foam::multiphaseMixtureScalar::Cp() const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -395,7 +395,7 @@ Foam::multiphaseMixture::Cp() const } Foam::tmp -Foam::multiphaseMixture::kf() const +Foam::multiphaseMixtureScalar::kf() const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -417,7 +417,7 @@ Foam::multiphaseMixture::kf() const } Foam::tmp -Foam::multiphaseMixture::D() const +Foam::multiphaseMixtureScalar::D() const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -439,7 +439,7 @@ Foam::multiphaseMixture::D() const } Foam::tmp -Foam::multiphaseMixture::Cs() const +Foam::multiphaseMixtureScalar::Cs() const { PtrDictionary::const_iterator iter = phases_.begin(); @@ -456,7 +456,7 @@ Foam::multiphaseMixture::Cs() const } Foam::tmp -Foam::multiphaseMixture::diffusionCorrection() const +Foam::multiphaseMixtureScalar::diffusionCorrection() const { surfaceScalarField numerator @@ -517,7 +517,7 @@ Foam::multiphaseMixture::diffusionCorrection() const return correction; } -void Foam::multiphaseMixture::solve() +void Foam::multiphaseMixtureScalar::solve() { correct(); @@ -570,7 +570,7 @@ void Foam::multiphaseMixture::solve() } -void Foam::multiphaseMixture::correct() +void Foam::multiphaseMixtureScalar::correct() { forAllIter(PtrDictionary, phases_, iter) { @@ -579,7 +579,7 @@ void Foam::multiphaseMixture::correct() } -Foam::tmp Foam::multiphaseMixture::nHatfv +Foam::tmp Foam::multiphaseMixtureScalar::nHatfv ( const volScalarField& alpha1, const volScalarField& alpha2 @@ -605,7 +605,7 @@ Foam::tmp Foam::multiphaseMixture::nHatfv } -Foam::tmp Foam::multiphaseMixture::nHatf +Foam::tmp Foam::multiphaseMixtureScalar::nHatf ( const volScalarField& alpha1, const volScalarField& alpha2 @@ -622,7 +622,7 @@ Foam::tmp Foam::multiphaseMixture::nHatf // The dynamic contact angle is calculated from the component of the // velocity on the direction of the interface, parallel to the wall. -void Foam::multiphaseMixture::correctContactAngle +void Foam::multiphaseMixtureScalar::correctContactAngle ( const phase& alpha1, const phase& alpha2, @@ -726,7 +726,7 @@ void Foam::multiphaseMixture::correctContactAngle } -Foam::tmp Foam::multiphaseMixture::K +Foam::tmp Foam::multiphaseMixtureScalar::K ( const phase& alpha1, const phase& alpha2 @@ -742,7 +742,7 @@ Foam::tmp Foam::multiphaseMixture::K Foam::tmp -Foam::multiphaseMixture::nearInterface() const +Foam::multiphaseMixtureScalar::nearInterface() const { tmp tnearInt ( @@ -768,7 +768,7 @@ Foam::multiphaseMixture::nearInterface() const } -void Foam::multiphaseMixture::solveAlphas +void Foam::multiphaseMixtureScalar::solveAlphas ( const scalar cAlpha ) @@ -901,7 +901,7 @@ void Foam::multiphaseMixture::solveAlphas } -bool Foam::multiphaseMixture::read() +bool Foam::multiphaseMixtureScalar::read() { if (transportModel::read()) { diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/multiphaseMixture.H b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/multiphaseMixtureScalar.H similarity index 96% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/multiphaseMixture.H rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/multiphaseMixtureScalar.H index 5fe5a939..f1fbd153 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/multiphaseMixture.H +++ b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/multiphaseMixtureScalar.H @@ -17,10 +17,10 @@ License Copyright (C) 2018- Mathias Vångö, JKU Linz, Austria Class - multiphaseMixture + multiphaseMixtureScalar Description - This class is based on the OpenFOAM(R) Foam::multiphaseMixture class, + This class is based on the OpenFOAM(R) Foam::multiphaseMixtureScalar class, which is an incompressible multi-phase mixture with built in solution for the phase fractions with interface compression for interface-capturing. It has been extended to include the void fraction in the volume fraction @@ -33,11 +33,11 @@ Description between each phase-pair. SourceFiles - multiphaseMixture.C + multiphaseMixtureScalar.C \*---------------------------------------------------------------------------*/ -#ifndef multiphaseMixture_H -#define multiphaseMixture_H +#ifndef multiphaseMixtureScalar_H +#define multiphaseMixtureScalar_H #include "incompressible/transportModel/transportModel.H" #include "IOdictionary.H" @@ -52,10 +52,10 @@ namespace Foam { /*---------------------------------------------------------------------------*\ - Class multiphaseMixture Declaration + Class multiphaseMixtureScalar Declaration \*---------------------------------------------------------------------------*/ -class multiphaseMixture +class multiphaseMixtureScalar : public IOdictionary, public transportModel @@ -191,7 +191,7 @@ public: // Constructors //- Construct from components - multiphaseMixture + multiphaseMixtureScalar ( const volVectorField& U, const surfaceScalarField& phi, @@ -200,7 +200,7 @@ public: //- Destructor - virtual ~multiphaseMixture() + virtual ~multiphaseMixtureScalar() {} diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/phase/phase.C b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/phase/phase.C similarity index 100% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/phase/phase.C rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/phase/phase.C diff --git a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/phase/phase.H b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/phase/phase.H similarity index 98% rename from applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/phase/phase.H rename to applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/phase/phase.H index 8237f374..0bab2d57 100644 --- a/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/phase/phase.H +++ b/applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/phase/phase.H @@ -26,7 +26,7 @@ Class Description Single incompressible phase derived from the phase-fraction. - Used as part of the multiPhaseMixture for interface-capturing multi-phase + Used as part of the multiphaseMixtureScalar for interface-capturing multi-phase simulations. SourceFiles diff --git a/etc/library-list.txt b/etc/library-list.txt index aa1be375..d26022e4 100644 --- a/etc/library-list.txt +++ b/etc/library-list.txt @@ -3,4 +3,4 @@ lagrangian/cfdemParticleComp/dir recurrence/dir finiteVolume/dir ../applications/solvers/cfdemSolverMultiphase/multiphaseMixture/dir -../applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixture/dir +../applications/solvers/cfdemSolverMultiphaseScalar/multiphaseMixtureScalar/dir From 6af30fe0e5cc1a41103dda4d5f0da1875b82ea6e Mon Sep 17 00:00:00 2001 From: Thomas Lichtenegger Date: Tue, 6 Dec 2022 08:15:50 +0100 Subject: [PATCH 22/23] Allow displacement field computation for polydisperse systems. --- .../displacementField/displacementField.C | 171 ++++++++++++------ .../CFD/constant/displacementProperties | 8 + 2 files changed, 127 insertions(+), 52 deletions(-) diff --git a/applications/utilities/displacementField/displacementField.C b/applications/utilities/displacementField/displacementField.C index 6f7be11e..ad73aaed 100644 --- a/applications/utilities/displacementField/displacementField.C +++ b/applications/utilities/displacementField/displacementField.C @@ -37,11 +37,19 @@ Application // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // void findPairs(labelList &, labelList &, labelPairList &); void findPairsUnordered(labelList &, labelList &, labelPairList &); -void fillEmptyCells(fvMesh &, label, label, labelList &, volVectorField &, volVectorField &, scalarList &, volVectorField &, volVectorField &, bool, scalar); -void nearestNeighborCells(fvMesh &, label, label, label, labelList &, labelList &); -void normalizeFields(labelList &, volVectorField &, volVectorField &); -void readDump(std::string, labelList &, vectorList &); +void fillEmptyCells(fvMesh &, label, label, scalarList &, volVectorField &, volVectorField &, scalarList &, volVectorField &, volVectorField &, bool, scalar); +void nearestNeighborCells(fvMesh &, label, label, label, scalarList &, labelList &); +void normalizeFields(scalarList &, volVectorField &, volVectorField &); +void readDump(std::string, labelList &, scalarList &, vectorList &); scalar weightFun(scalar); +label maxNumParticles = 1000000; +scalar minVol = 1e-12; +scalar Pi43 = 4.1888; +label posIndex = -1; +label posRad = -1; +label posX = -1; +label posY = -1; +label posZ = -1; int main(int argc, char *argv[]) { @@ -88,6 +96,11 @@ int main(int argc, char *argv[]) label dumpIndexDisplacementIncrement(readLabel(displacementProperties.lookup("dumpIndexDisplacementIncrement"))); label nNeighMin(readLabel(displacementProperties.lookup("nNeighMin"))); label maxSearchLayers(displacementProperties.lookupOrDefault