Skip to content

Commit 02152f6

Browse files
committed
Add sensitivity parType dependence multiplex for parameters with BndCompTypeSec dependence (#419)
1 parent 3661655 commit 02152f6

File tree

5 files changed

+147
-49
lines changed

5 files changed

+147
-49
lines changed

src/libcadet/model/GeneralRateModel.cpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -2843,7 +2843,7 @@ bool GeneralRateModel<ConvDispOperator>::setParameter(const ParameterId& pId, do
28432843
return true;
28442844
if (multiplexCompTypeSecParameterValue(pId, hashString("PAR_DIFFUSION"), _parDiffusionMode, _parDiffusion, _disc.nParType, _disc.nComp, value, nullptr))
28452845
return true;
2846-
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, value, nullptr))
2846+
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, value, nullptr))
28472847
return true;
28482848
const int mpIc = multiplexInitialConditions(pId, value, false);
28492849
if (mpIc > 0)
@@ -2923,7 +2923,7 @@ void GeneralRateModel<ConvDispOperator>::setSensitiveParameterValue(const Parame
29232923
return;
29242924
if (multiplexCompTypeSecParameterValue(pId, hashString("PAR_DIFFUSION"), _parDiffusionMode, _parDiffusion, _disc.nParType, _disc.nComp, value, &_sensParams))
29252925
return;
2926-
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, value, &_sensParams))
2926+
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, value, &_sensParams))
29272927
return;
29282928
if (multiplexInitialConditions(pId, value, true) != 0)
29292929
return;
@@ -2989,7 +2989,7 @@ bool GeneralRateModel<ConvDispOperator>::setSensitiveParameter(const ParameterId
29892989
return true;
29902990
}
29912991

2992-
if (multiplexBndCompTypeSecParameterAD(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, adDirection, adValue, _sensParams))
2992+
if (multiplexBndCompTypeSecParameterAD(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, adDirection, adValue, _sensParams))
29932993
{
29942994
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
29952995
return true;

src/libcadet/model/GeneralRateModel2D.cpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -2341,7 +2341,7 @@ bool GeneralRateModel2D::setParameter(const ParameterId& pId, double value)
23412341
return true;
23422342
if (multiplexCompTypeSecParameterValue(pId, hashString("PAR_DIFFUSION"), _parDiffusionMode, _parDiffusion, _disc.nParType, _disc.nComp, value, nullptr))
23432343
return true;
2344-
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, value, nullptr))
2344+
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, value, nullptr))
23452345
return true;
23462346
const int mpIc = multiplexInitialConditions(pId, value, false);
23472347
if (mpIc > 0)
@@ -2381,7 +2381,7 @@ void GeneralRateModel2D::setSensitiveParameterValue(const ParameterId& pId, doub
23812381
return;
23822382
if (multiplexCompTypeSecParameterValue(pId, hashString("PAR_DIFFUSION"), _parDiffusionMode, _parDiffusion, _disc.nParType, _disc.nComp, value, &_sensParams))
23832383
return;
2384-
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, value, &_sensParams))
2384+
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, value, &_sensParams))
23852385
return;
23862386
if (multiplexInitialConditions(pId, value, true) != 0)
23872387
return;
@@ -2432,7 +2432,7 @@ bool GeneralRateModel2D::setSensitiveParameter(const ParameterId& pId, unsigned
24322432
return true;
24332433
}
24342434

2435-
if (multiplexBndCompTypeSecParameterAD(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, adDirection, adValue, _sensParams))
2435+
if (multiplexBndCompTypeSecParameterAD(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, adDirection, adValue, _sensParams))
24362436
{
24372437
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
24382438
return true;

src/libcadet/model/GeneralRateModelDG.cpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -2145,7 +2145,7 @@ bool GeneralRateModelDG::setParameter(const ParameterId& pId, double value)
21452145
return true;
21462146
if (multiplexCompTypeSecParameterValue(pId, hashString("PAR_DIFFUSION"), _parDiffusionMode, _parDiffusion, _disc.nParType, _disc.nComp, value, nullptr))
21472147
return true;
2148-
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, value, nullptr))
2148+
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, value, nullptr))
21492149
return true;
21502150
const int mpIc = multiplexInitialConditions(pId, value, false);
21512151
if (mpIc > 0)
@@ -2222,7 +2222,7 @@ void GeneralRateModelDG::setSensitiveParameterValue(const ParameterId& pId, doub
22222222
return;
22232223
if (multiplexCompTypeSecParameterValue(pId, hashString("PAR_DIFFUSION"), _parDiffusionMode, _parDiffusion, _disc.nParType, _disc.nComp, value, &_sensParams))
22242224
return;
2225-
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, value, &_sensParams))
2225+
if (multiplexBndCompTypeSecParameterValue(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, value, &_sensParams))
22262226
return;
22272227
if (multiplexInitialConditions(pId, value, true) != 0)
22282228
return;
@@ -2287,7 +2287,7 @@ bool GeneralRateModelDG::setSensitiveParameter(const ParameterId& pId, unsigned
22872287
return true;
22882288
}
22892289

2290-
if (multiplexBndCompTypeSecParameterAD(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, adDirection, adValue, _sensParams))
2290+
if (multiplexBndCompTypeSecParameterAD(pId, hashString("PAR_SURFDIFFUSION"), _parSurfDiffusionMode, _parSurfDiffusion, _disc.nParType, _disc.nComp, _disc.strideBound, _disc.nBound, _disc.boundOffset, _disc.nBoundBeforeType, adDirection, adValue, _sensParams))
22912291
{
22922292
LOG(Debug) << "Found parameter " << pId << ": Dir " << adDirection << " is set to " << adValue;
22932293
return true;

src/libcadet/model/ParameterMultiplexing.cpp

Lines changed: 134 additions & 38 deletions
Original file line numberDiff line numberDiff line change
@@ -143,7 +143,7 @@ MultiplexMode readAndRegisterMultiplexCompTypeSecParam(IParameterProvider& param
143143

144144
values = std::move(p);
145145

146-
for (std::size_t s = 0; s < values.size() / nComp; ++s)
146+
for (std::size_t s = 0; s < values.size() / nComp / nParType; ++s)
147147
{
148148
for (unsigned int i = 0; i < nComp; ++i)
149149
parameters[makeParamId(nameHash, uoi, i, ParTypeIndep, BoundStateIndep, ReactionIndep, s)] = &values[s * nParType * nComp + i];
@@ -206,7 +206,32 @@ bool multiplexCompTypeSecParameterValue(const ParameterId& pId, StringHash nameH
206206
return true;
207207
}
208208
case MultiplexMode::ComponentType:
209+
{
210+
if ((pId.component == CompIndep) || (pId.particleType == ParTypeIndep) || (pId.boundState != BoundStateIndep)
211+
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
212+
return false;
213+
214+
if (sensParams && !contains(*sensParams, &data[pId.particleType * nComp + pId.component]))
215+
return false;
216+
217+
for (unsigned int s = 0; s < data.size() / nComp / nParType; ++s)
218+
data[s * nComp * nParType + pId.particleType * nComp + pId.component].setValue(value);
219+
220+
return true;
221+
}
209222
case MultiplexMode::ComponentSectionType:
223+
{
224+
if ((pId.component == CompIndep) || (pId.particleType == ParTypeIndep) || (pId.boundState != BoundStateIndep)
225+
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
226+
return false;
227+
228+
if (sensParams && !contains(*sensParams, &data[pId.section * nComp * nParType + pId.particleType * nComp + pId.component]))
229+
return false;
230+
231+
data[pId.section * nComp * nParType + pId.particleType * nComp + pId.component].setValue(value);
232+
233+
return true;
234+
}
210235
case MultiplexMode::RadialSection:
211236
case MultiplexMode::Independent:
212237
case MultiplexMode::ComponentRadial:
@@ -231,33 +256,56 @@ bool multiplexCompTypeSecParameterAD(const ParameterId& pId, StringHash nameHash
231256
switch (mode)
232257
{
233258
case MultiplexMode::Component:
234-
{
235-
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState != BoundStateIndep)
236-
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
237-
return false;
259+
{
260+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState != BoundStateIndep)
261+
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
262+
return false;
238263

239-
sensParams.insert(&data[pId.component]);
264+
sensParams.insert(&data[pId.component]);
240265

241-
for (unsigned int i = 0; i < nParType; ++i)
242-
data[i * nComp + pId.component].setADValue(adDirection, adValue);
266+
for (unsigned int i = 0; i < nParType; ++i)
267+
data[i * nComp + pId.component].setADValue(adDirection, adValue);
243268

244-
return true;
245-
}
269+
return true;
270+
}
246271
case MultiplexMode::ComponentSection:
247-
{
248-
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState != BoundStateIndep)
249-
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
250-
return false;
272+
{
273+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState != BoundStateIndep)
274+
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
275+
return false;
251276

252-
sensParams.insert(&data[pId.section * nComp * nParType + pId.component]);
277+
sensParams.insert(&data[pId.section * nComp * nParType + pId.component]);
253278

254-
for (unsigned int i = 0; i < nParType; ++i)
255-
data[i * nComp + pId.section * nComp * nParType + pId.component].setADValue(adDirection, adValue);
279+
for (unsigned int i = 0; i < nParType; ++i)
280+
data[i * nComp + pId.section * nComp * nParType + pId.component].setADValue(adDirection, adValue);
256281

257-
return true;
258-
}
282+
return true;
283+
}
259284
case MultiplexMode::ComponentType:
285+
{
286+
if ((pId.component == CompIndep) || (pId.particleType == ParTypeIndep) || (pId.boundState != BoundStateIndep)
287+
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
288+
return false;
289+
290+
sensParams.insert(&data[pId.particleType * nComp + pId.component]);
291+
292+
for (unsigned int s = 0; s < data.size() / nComp / nParType; ++s)
293+
data[s * nComp * nParType + pId.particleType * nComp + pId.component].setADValue(adDirection, adValue);
294+
295+
return true;
296+
}
260297
case MultiplexMode::ComponentSectionType:
298+
{
299+
if ((pId.component == CompIndep) || (pId.particleType == ParTypeIndep) || (pId.boundState != BoundStateIndep)
300+
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
301+
return false;
302+
303+
sensParams.insert(&data[pId.section * nComp * nParType + pId.particleType * nComp + pId.component]);
304+
305+
data[pId.section * nComp * nParType + pId.particleType * nComp + pId.component].setADValue(adDirection, adValue);
306+
307+
return true;
308+
}
261309
case MultiplexMode::RadialSection:
262310
case MultiplexMode::Independent:
263311
case MultiplexMode::ComponentRadial:
@@ -404,7 +452,7 @@ MultiplexMode readAndRegisterMultiplexBndCompTypeSecParam(IParameterProvider& pa
404452
}
405453

406454
bool multiplexBndCompTypeSecParameterValue(const ParameterId& pId, StringHash nameHash, MultiplexMode mode, std::vector<active>& data,
407-
unsigned int nParType, unsigned int nComp, unsigned int const* strideBound, unsigned int const* nBound, unsigned int const* boundOffset, double value, std::unordered_set<active*> const* sensParams)
455+
unsigned int nParType, unsigned int nComp, unsigned int const* strideBound, unsigned int const* nBound, unsigned int const* boundOffset, unsigned int const* nBoundBeforeType, double value, std::unordered_set<active*> const* sensParams)
408456
{
409457
if (pId.name != nameHash)
410458
return false;
@@ -443,7 +491,32 @@ bool multiplexBndCompTypeSecParameterValue(const ParameterId& pId, StringHash na
443491
return true;
444492
}
445493
case MultiplexMode::ComponentType:
494+
{
495+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
496+
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
497+
return false;
498+
499+
if (sensParams && !contains(*sensParams, &data[nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState]))
500+
return false;
501+
502+
for (unsigned int s = 0; s < data.size() / nComp / nParType; ++s)
503+
data[s * strideBound[nParType] + nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState].setValue(value);
504+
505+
return true;
506+
}
446507
case MultiplexMode::ComponentSectionType:
508+
{
509+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
510+
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
511+
return false;
512+
513+
if (sensParams && !contains(*sensParams, &data[pId.section * strideBound[nParType] + nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState]))
514+
return false;
515+
516+
data[pId.section * strideBound[nParType] + nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState].setValue(value);
517+
518+
return true;
519+
}
447520
case MultiplexMode::RadialSection:
448521
case MultiplexMode::Independent:
449522
case MultiplexMode::ComponentRadial:
@@ -461,7 +534,7 @@ bool multiplexBndCompTypeSecParameterValue(const ParameterId& pId, StringHash na
461534
}
462535

463536
bool multiplexBndCompTypeSecParameterAD(const ParameterId& pId, StringHash nameHash, MultiplexMode mode, std::vector<active>& data,
464-
unsigned int nParType, unsigned int nComp, unsigned int const* strideBound, unsigned int const* nBound, unsigned int const* boundOffset, unsigned int adDirection, double adValue, std::unordered_set<active*>& sensParams)
537+
unsigned int nParType, unsigned int nComp, unsigned int const* strideBound, unsigned int const* nBound, unsigned int const* boundOffset, unsigned int const* nBoundBeforeType, unsigned int adDirection, double adValue, std::unordered_set<active*>& sensParams)
465538
{
466539
if (pId.name != nameHash)
467540
return false;
@@ -472,33 +545,56 @@ bool multiplexBndCompTypeSecParameterAD(const ParameterId& pId, StringHash nameH
472545
switch (mode)
473546
{
474547
case MultiplexMode::Component:
475-
{
476-
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
477-
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
478-
return false;
548+
{
549+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
550+
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
551+
return false;
479552

480-
sensParams.insert(&data[boundOffset[pId.component] + pId.boundState]);
553+
sensParams.insert(&data[boundOffset[pId.component] + pId.boundState]);
481554

482-
for (unsigned int i = 0; i < nParType; ++i)
483-
data[boundOffset[pId.component] + pId.boundState + i * strideBound[0]].setADValue(adDirection, adValue);
555+
for (unsigned int i = 0; i < nParType; ++i)
556+
data[boundOffset[pId.component] + pId.boundState + i * strideBound[0]].setADValue(adDirection, adValue);
484557

485-
return true;
486-
}
487-
case MultiplexMode::ComponentSection:
488-
{
489-
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
490-
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
491-
return false;
558+
return true;
559+
}
560+
case MultiplexMode::ComponentSection: // note that due to ParTypeIndep, boundOffsets and strideBounds are the same for all particle types
561+
{
562+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
563+
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
564+
return false;
492565

493566
sensParams.insert(&data[pId.section * strideBound[nParType] + boundOffset[pId.component] + pId.boundState]);
494567

495-
for (unsigned int i = 0; i < nParType; ++i)
568+
for (unsigned int i = 0; i < nParType; ++i)
496569
data[pId.section * strideBound[nParType] + boundOffset[pId.component] + pId.boundState + i * strideBound[0]].setADValue(adDirection, adValue);
497570

498-
return true;
499-
}
571+
return true;
572+
}
500573
case MultiplexMode::ComponentType:
574+
{
575+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
576+
|| (pId.reaction != ReactionIndep) || (pId.section != SectionIndep))
577+
return false;
578+
579+
sensParams.insert(&data[nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState]);
580+
581+
for (unsigned int s = 0; s < data.size() / nComp / nParType; ++s)
582+
data[s * strideBound[nParType] + nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState].setADValue(adDirection, adValue);
583+
584+
return true;
585+
}
501586
case MultiplexMode::ComponentSectionType:
587+
{
588+
if ((pId.component == CompIndep) || (pId.particleType != ParTypeIndep) || (pId.boundState == BoundStateIndep)
589+
|| (pId.reaction != ReactionIndep) || (pId.section == SectionIndep))
590+
return false;
591+
592+
sensParams.insert(&data[pId.section * strideBound[nParType] + nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState]);
593+
594+
data[pId.section * strideBound[nParType] + nBoundBeforeType[pId.particleType] + boundOffset[pId.particleType * nComp + pId.component] + pId.boundState].setADValue(adDirection, adValue);
595+
596+
return true;
597+
}
502598
case MultiplexMode::RadialSection:
503599
case MultiplexMode::Independent:
504600
case MultiplexMode::ComponentRadial:

0 commit comments

Comments
 (0)