Hi All,
I want to add stresses value to the vtu file for the el2p problem. Based on
geomechanics/poroselastic/main.cc, I tried to implementing the same for
el2p problem. But I getting error segmentation fault error.
// function to evaluate the element stresses
template< class StressType,
class Problem,
class GridGeometry,
class GridVariables,
class SolutionVector,
class SigmaStorage>
void assembleElementStresses(SigmaStorage& sigmaStorage,
SigmaStorage& effSigmaStorage,
const Problem& problem,
const GridGeometry& fvGridGeometry,
const GridVariables& gridVariables,
const SolutionVector& x)
{
for (const auto& element : elements(fvGridGeometry.gridView()))
{
auto fvGeometry = localView(fvGridGeometry);
auto elemVolVars = localView(gridVariables.curGridVolVars());
fvGeometry.bind(element);
elemVolVars.bind(element, fvGeometry, x);
// evaluate flux variables cache at cell center
using FluxVarsCache = typename
GridVariables::GridFluxVariablesCache::FluxVariablesCache;
FluxVarsCache fluxVarsCache;
fluxVarsCache.update(problem, element, fvGeometry, elemVolVars,
element.geometry().center());
//get lame parameters, the pressure and compute stress tensor
const auto sigma = StressType::stressTensor(problem, element,
fvGeometry, elemVolVars, fluxVarsCache); // this gives error
const auto effSigma = StressType::effectiveStressTensor(problem,
element, fvGeometry, elemVolVars, fluxVarsCache);
//std::vector<int> v{1,2,3};
const auto eIdx = fvGridGeometry.elementMapper().index(element);
// pass values into storage container
using FVGridGeometry = typename GridVariables::GridGeometry;
for (int dir = 0; dir < FVGridGeometry::GridView::dimension; ++dir)
{
sigmaStorage[dir][eIdx] = sigma[dir];
//sigmaStorage[dir][eIdx] = v[dir]; //sigma[dir];
effSigmaStorage[dir][eIdx] = effSigma[dir];
//effSigmaStorage[dir][eIdx] = v[dir]; //effSigma[dir];
}
}
}
In main method,
//containers to store the sigma
std::array<std::vector<ForceVector>, dim> sigmaStorage;
std::array<std::vector<ForceVector>, dim> effSigmaStorage;
const auto numCells = poroMechFvGridGeometry->gridView().size(0);
std::for_each(sigmaStorage.begin(), sigmaStorage.end(),
[numCells](auto& sigma){sigma.resize(numCells);});
std::for_each(effSigmaStorage.begin(), effSigmaStorage.end(),
[numCells](auto& effSigma){effSigma.resize(numCells);});
for (int dir = 0; dir < dim; ++dir)
{
poroMechVtkWriter.addField(sigmaStorage[dir], "sigma_" +
std::to_string(dir), PoroMechVtkOutputModule::FieldType::element);
poroMechVtkWriter.addField(sigmaStorage[dir], "effSigma_" +
std::to_string(dir), PoroMechVtkOutputModule::FieldType::element);
}
using StressType = GetPropType<PoroMechTypeTag, Properties::StressType>;
assembleElementStresses<StressType>(sigmaStorage, effSigmaStorage,
*poroMechProblem, *poroMechFvGridGeometry,
*poroMechGridVariables, x[poroMechId]);
Any help is appreciated.
Thank You!,
Regards,
Ranjeet
_______________________________________________
Dumux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
