Skip to content

Commit d87f214

Browse files
jrgparkinsonjrgparkinson
authored
Merge pull request #33 from jrgparkinson/development
Fix BC calculation
2 parents c032446 + ae68d9f commit d87f214

File tree

8 files changed

+156
-85
lines changed

8 files changed

+156
-85
lines changed

‎.gitignore

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -46,4 +46,5 @@ visitlog.py
4646
**/matlab/MethodsPaper/*.png
4747
*.csv
4848
/zzindex.html
49-
*.ipynb_checkpoints*
49+
*.ipynb_checkpoints*
50+
/requirements (copy).txt

‎params/agu2020/.gitignore

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
/inputs.2d.1x1.1024
2+
/inputs.2d.1x1.128
3+
/inputs.2d.1x1.256
4+
/inputs.2d.1x1.32
5+
/inputs.2d.1x1.512
6+
/inputs.2d.1x1.64
7+
/new
8+
/new_inputs

‎src/AMRNonLinearMultiCompOp.cpp

Lines changed: 1 addition & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -101,7 +101,6 @@ void AMRNonLinearMultiCompOp::computeDiffusedVar(FArrayBox& a_diffusedVar,
101101
m_diffusedVarBC(a_diffusedVar, m_domain.domainBox(), m_domain, m_dx, a_homogeneous);
102102
}
103103

104-
// this->m_bc
105104
}
106105

107106
void AMRNonLinearMultiCompOp::residualI(LevelData<FArrayBox>& a_lhs,
@@ -122,7 +121,7 @@ void AMRNonLinearMultiCompOp::residualI(LevelData<FArrayBox>& a_lhs,
122121

123122
for (dit.begin(); dit.ok(); ++dit)
124123
{
125-
m_bc(phi[dit], dbl[dit()],m_domain, dx, a_homogeneous);
124+
m_bc(phi[dit], dbl[dit()], m_domain, dx, a_homogeneous);
126125
}
127126
}
128127

@@ -140,8 +139,6 @@ void AMRNonLinearMultiCompOp::residualI(LevelData<FArrayBox>& a_lhs,
140139
// Is this quicker than the level data version?
141140
computeDiffusedVar(derivedVar[dit], phi[dit], dit(), a_homogeneous);
142141

143-
144-
145142
#if CH_SPACEDIM == 1
146143
FORT_NLVCCOMPUTERES1D
147144
#elif CH_SPACEDIM == 2

‎src/CellToEdge2.cpp

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,8 @@
1919

2020
#include "NamespaceHeader.H"
2121

22+
23+
2224
/// -----------------------------------------------------------
2325
void CellToEdge2(const LevelData<FArrayBox>& a_cellData,
2426
LevelData<FluxBox>& a_edgeData,
@@ -81,6 +83,7 @@ void CellToEdge2(const FArrayBox& a_cellData,
8183

8284

8385
}
86+
8487
}
8588
}
8689
}

‎src/MushyLayerParams.cpp

Lines changed: 88 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,26 @@ string MushyLayerParams::s_ADVECTIVE_VELOCITY_SCALE = "Advective";
3232
/// Darcy velocity scale
3333
string MushyLayerParams::s_DARCY_VELOCITY_SCALE = "Darcy";
3434

35+
void MushyLayerParams::printBCs(string bcName, Vector<int> bcTypeLo, Vector<int> bcTypeHi, RealVect bcValLo, RealVect bcValHi)
36+
{
37+
char buffer [100];
38+
sprintf(buffer, "%-10s", bcName.c_str());
39+
pout() << buffer;
40+
41+
Vector<string> bcTypeNames;
42+
bcTypeNames.push_back("Fixed");
43+
bcTypeNames.push_back("No flux");
44+
45+
for (int dir=0; dir < SpaceDim; dir ++)
46+
{
47+
sprintf(buffer, "|%-10s= %02.2f|%-10s= %02.2f", m_scalarBCTypes[bcTypeLo[dir]].c_str(), bcValLo[dir],
48+
m_scalarBCTypes[bcTypeHi[dir]].c_str(), bcValHi[dir]);
49+
pout() << buffer;
50+
// pout() << "| " << bcTypeLo[dir] << "= " << bcValLo[dir] << " | " << bcTypeHi[dir] << "= " << bcValHi[dir];
51+
}
52+
pout() << endl;
53+
}
54+
3555
// Initialise in initialization list
3656
MushyLayerParams::MushyLayerParams() : physicalProblem(PhysicalProblems::m_mushyLayer),
3757
viscosity(-999),
@@ -132,37 +152,69 @@ MushyLayerParams::MushyLayerParams() : physicalProblem(PhysicalProblems::m_mushy
132152
m_porosityFunction(ParamsPorosityFunctions::m_porosityConstant),
133153
m_viscosityFunction(ViscosityFunction::uniformViscosity),
134154
max_viscosity(1.0)
135-
{}
155+
{
156+
m_scalarBCTypes.push_back("fixed");
157+
m_scalarBCTypes.push_back("noflux");
158+
m_scalarBCTypes.push_back("open");
159+
m_scalarBCTypes.push_back("inflow");
160+
161+
m_vectorBCTypes.push_back("noflow");
162+
m_vectorBCTypes.push_back("inflowVelocity");
163+
m_vectorBCTypes.push_back("open");
164+
m_vectorBCTypes.push_back("outflownormal");
165+
m_vectorBCTypes.push_back("inflowoutflow");
166+
m_vectorBCTypes.push_back("noshear");
167+
m_vectorBCTypes.push_back("symmetry");
168+
m_vectorBCTypes.push_back("inflowPlume");
169+
m_vectorBCTypes.push_back("outflowPressureGrad");
170+
m_vectorBCTypes.push_back("pressureHead");
171+
}
136172

137173
MushyLayerParams::~MushyLayerParams() {
138174
}
139175

140176
void
141177
MushyLayerParams::computeDerivedBCs ()
142178
{
179+
pout() << "MushyLayerParams::computeDerivedBCs ThetaInitial=" << ThetaInitial
180+
<< ", H top = " << bcValEnthalpyHi[SpaceDim-1] << endl;
143181
// Now do defaults for the other boundary values
144182
for (int dir = 0; dir < SpaceDim; dir++) {
183+
Real bulkCHi = bcValBulkConcentrationHi[dir];
184+
Real bulkCLo = bcValBulkConcentrationLo[dir];
185+
186+
if (bcTypeBulkConcentrationHi[dir] == PhysBCUtil::Neumann)
187+
{
188+
bulkCHi = ThetaInitial;
189+
}
190+
191+
if (bcTypeBulkConcentrationLo[dir] == PhysBCUtil::Neumann)
192+
{
193+
bulkCLo = ThetaInitial;
194+
}
195+
145196
::computeBoundingEnergy (bcValEnthalpyHi[dir],
146-
bcValBulkConcentrationHi[dir], bcValSolidusHi[dir],
197+
bulkCHi, bcValSolidusHi[dir],
147198
bcValLiquidusHi[dir], bcValEutecticHi[dir],
148199
specificHeatRatio, stefan, compositionRatio,
149200
waterDistributionCoeff, thetaEutectic,
150201
ThetaEutectic);
151-
::computeEnthalpyVars (bcValEnthalpyHi[dir], bcValBulkConcentrationHi[dir],
202+
::computeEnthalpyVars (bcValEnthalpyHi[dir], bulkCHi,
152203
bcValPorosityHi[dir], bcValTemperatureHi[dir],
153204
bcValLiquidConcentrationHi[dir],
154205
bcValSolidConcentrationHi[dir], bcValSolidusHi[dir],
155206
bcValLiquidusHi[dir], bcValEutecticHi[dir],
156207
specificHeatRatio, stefan, compositionRatio,
157208
waterDistributionCoeff, thetaEutectic,
158209
ThetaEutectic);
159-
::computeBoundingEnergy (bcValEnthalpyLo[dir],
160-
bcValBulkConcentrationLo[dir], bcValSolidusLo[dir],
210+
211+
::computeBoundingEnergy (bcValEnthalpyLo[dir], bulkCLo,
212+
bcValSolidusLo[dir],
161213
bcValLiquidusLo[dir], bcValEutecticLo[dir],
162214
specificHeatRatio, stefan, compositionRatio,
163215
waterDistributionCoeff, thetaEutectic,
164216
ThetaEutectic);
165-
::computeEnthalpyVars (bcValEnthalpyLo[dir], bcValBulkConcentrationLo[dir],
217+
::computeEnthalpyVars (bcValEnthalpyLo[dir], bulkCLo,
166218
bcValPorosityLo[dir], bcValTemperatureLo[dir],
167219
bcValLiquidConcentrationLo[dir],
168220
bcValSolidConcentrationLo[dir], bcValSolidusLo[dir],
@@ -183,10 +235,24 @@ MushyLayerParams::computeDerivedBCs ()
183235
else
184236
{
185237
Hinitial = bcValEnthalpyLo[dir];
186-
;
187238
}
188239
}
189240
}
241+
242+
pout() << "After computeDerivedBCs: " << endl << " ";
243+
Vector<string> direction; direction.push_back("x"); direction.push_back("y"); direction.push_back("z");
244+
for (int dir=0; dir < SpaceDim; dir ++)
245+
{
246+
pout() << "| " << direction[dir] << "(lo) | " << direction[dir] << "(hi) ";
247+
}
248+
pout () << endl;
249+
250+
// printBCs("H", bcTypeEnthalpyLo, bcTypeEnthalpyHi, bcValEnthalpyLo, bcValEnthalpyHi);
251+
// printBCs("C", bcTypeBulkConcentrationLo, bcTypeBulkConcentrationHi, bcValBulkConcentrationLo, bcValBulkConcentrationHi);
252+
printBCs("T", bcTypeTemperatureLo, bcTypeTemperatureHi, bcValTemperatureLo, bcValTemperatureHi);
253+
printBCs("Sl", bcTypeLiquidConcentrationLo, bcTypeLiquidConcentrationHi, bcValLiquidConcentrationLo, bcValLiquidConcentrationHi);
254+
printBCs("porosity", bcTypePorosityLo, bcTypePorosityHi, bcValPorosityLo, bcValPorosityHi);
255+
190256
}
191257

192258
void MushyLayerParams::getParameters()
@@ -628,7 +694,6 @@ void MushyLayerParams::getParameters()
628694

629695

630696
// Now BC values
631-
632697
parseBCVals("liquidConcentrationLoVal", bcValLiquidConcentrationLo);
633698
parseBCVals("liquidConcentrationHiVal", bcValLiquidConcentrationHi);
634699

@@ -794,6 +859,21 @@ void MushyLayerParams::getParameters()
794859

795860
// uncomment to print the parameters
796861
// printParameters();
862+
863+
pout() << "BCs used: " << endl << " ";
864+
Vector<string> direction; direction.push_back("x"); direction.push_back("y"); direction.push_back("z");
865+
for (int dir=0; dir < SpaceDim; dir ++)
866+
{
867+
pout() << "| " << direction[dir] << "(lo) | " << direction[dir] << "(hi) ";
868+
}
869+
pout () << endl;
870+
871+
printBCs("H", bcTypeEnthalpyLo, bcTypeEnthalpyHi, bcValEnthalpyLo, bcValEnthalpyHi);
872+
printBCs("C", bcTypeBulkConcentrationLo, bcTypeBulkConcentrationHi, bcValBulkConcentrationLo, bcValBulkConcentrationHi);
873+
printBCs("T", bcTypeTemperatureLo, bcTypeTemperatureHi, bcValTemperatureLo, bcValTemperatureHi);
874+
printBCs("Sl", bcTypeLiquidConcentrationLo, bcTypeLiquidConcentrationHi, bcValLiquidConcentrationLo, bcValLiquidConcentrationHi);
875+
printBCs("porosity", bcTypePorosityLo, bcTypePorosityHi, bcValPorosityLo, bcValPorosityHi);
876+
797877
}
798878

799879
void MushyLayerParams::parseBCs(string a_name, Vector<int>* a_bcHolder, bool required)

‎src/MushyLayerParams.h

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -722,6 +722,9 @@ class MushyLayerParams {
722722
*/
723723
Real max_viscosity;
724724

725+
Vector<string> m_scalarBCTypes;
726+
Vector<string> m_vectorBCTypes;
727+
725728
/// Get all parameters from the inputs file
726729
void getParameters();
727730

@@ -781,6 +784,9 @@ class MushyLayerParams {
781784

782785
/// Compute boundary conditions for fields which can be found from the enthalpy and salinity via the phase diagram
783786
void computeDerivedBCs ();
787+
788+
void printBCs(string bcName, Vector<int> bcTypeLo, Vector<int> bcTypeHi, RealVect bcValLo, RealVect bcValHi);
789+
784790
};
785791

786792
#endif /* SRC_MUSHYLAYERPARAMS_H_ */

‎srcSubcycle/AMRLevelMushyLayerDiagnostics.cpp

Lines changed: 0 additions & 31 deletions
Original file line numberDiff line numberDiff line change
@@ -561,8 +561,6 @@ void AMRLevelMushyLayer::computeDiagnostics()
561561
}
562562

563563
// Work out mushy layer depth
564-
565-
566564
Real depth = computeMushDepth();
567565

568566
if (calcDiagnostics
@@ -580,35 +578,6 @@ void AMRLevelMushyLayer::computeDiagnostics()
580578
// Real mushAvBulkC = 0.0;
581579
Real mushAvPorosity = averageOverMushyRegion(ScalarVars::m_porosity, mushVol);
582580

583-
// int numMushyCells = 0;
584-
585-
// int lo_j = m_problem_domain.domainBox().smallEnd()[1];
586-
587-
// for (DataIterator dit = m_grids.dataIterator(); dit.ok(); ++dit)
588-
// {
589-
// FArrayBox& porosity = (*m_scalarNew[ScalarVars::m_porosity])[dit];
590-
// FArrayBox& bulkConc = (*m_scalarNew[ScalarVars::m_bulkConcentration])[dit];
591-
//
592-
// for (BoxIterator bit = BoxIterator(m_grids[dit]); bit.ok(); ++bit)
593-
// {
594-
// IntVect iv = bit();
595-
// RealVect loc;
596-
// ::getLocation(iv, loc, m_dx);
597-
//
598-
// bool is_sea_ice = loc[1] > (m_domainHeight-depth);
599-
//
600-
// if (is_sea_ice)
601-
// {
602-
// numMushyCells++;
603-
// mushAvBulkC += bulkConc(iv);
604-
// mushAvPorosity += porosity(iv);
605-
// }
606-
// }
607-
// }
608-
// mushAvBulkC = mushAvBulkC / numMushyCells;
609-
// mushAvPorosity = mushAvPorosity / numMushyCells;
610-
// mushVol = numMushyCells*m_dx*m_dx;
611-
612581
// Convert mush volume in number of cells to physical volume (in dimensionless units)
613582
mushVol *= pow(m_dx, SpaceDim);
614583

0 commit comments

Comments
 (0)