From dce684df5ed27edabf867d13f254cacb16b9684c Mon Sep 17 00:00:00 2001
From: Laura Pfeiffer <laura.pfeiffer@kit.edu>
Date: Mon, 5 Aug 2024 17:49:56 +0200
Subject: [PATCH] hopeffully fixed problems with adaptive
---
mpp | 2 +-
src/coupled/solvers/SegregatedSolver.cpp | 22 +++++++++++---------
src/elasticity/solvers/ElasticityNewmark.cpp | 1 -
src/elasticity/solvers/ElasticityNewmark.hpp | 2 +-
tools/solving/CoupledHoreka.py | 2 +-
5 files changed, 15 insertions(+), 14 deletions(-)
diff --git a/mpp b/mpp
index e04e54352..20f96e0f4 160000
--- a/mpp
+++ b/mpp
@@ -1 +1 @@
-Subproject commit e04e54352c9cfe78ab47c31d4cb360e6c91f4d4d
+Subproject commit 20f96e0f440ee73a40da6cacdcc95c6ac77a4494
diff --git a/src/coupled/solvers/SegregatedSolver.cpp b/src/coupled/solvers/SegregatedSolver.cpp
index 92f6ad19d..274c7c3fc 100644
--- a/src/coupled/solvers/SegregatedSolver.cpp
+++ b/src/coupled/solvers/SegregatedSolver.cpp
@@ -272,11 +272,10 @@ bool adaptiveSegregatedSolver::Step(IElphyAssemble &elphyAssemble, IElasticity &
elphyAssemble.ResetTime(mechAssemble.Time(), mechAssemble.LastTStep(), elphyDeltaTime);
double gammaInfNorm=elphyAssemble.InftyNormOnVertices(v_stretch(elphyValues));
- vout(2)<<"gammaInfNorm "<<gammaInfNorm<<endl;
int elphyStepCounter=0;
- double minGammaInfNorm =0.00001;
- if (gammaInfNorm<minGammaInfNorm) {
- while (elphyStepCounter != elphySteps && gammaInfNorm < minGammaInfNorm) {
+
+ if (gammaInfNorm<=0.0) {
+ while (elphyStepCounter != elphySteps && gammaInfNorm <=0.0) {
//while (!elphyAssemble.IsFinished()) {
elphySolver->printSolverNameInStep(elphyAssemble, ongoingElphyStep);
elphySolver->Step(elphyAssemble, elphyValues);
@@ -285,12 +284,12 @@ bool adaptiveSegregatedSolver::Step(IElphyAssemble &elphyAssemble, IElasticity &
ongoingElphyStep += 1;
elphyStepCounter += 1;
}
- }else if(minGammaInfNorm<=gammaInfNorm && gammaInfNorm<factorMaxGammaInfnorm*maxGammaInfNorm) {
+ }else if(0.0<gammaInfNorm && gammaInfNorm<factorMaxGammaInfnorm*maxGammaInfNorm) {
int numberOfInbetweenTimSteps=int(log2(elphySteps))-1;
vout(2)<<" numberOfInbetweenTimSteps "<<numberOfInbetweenTimSteps;
for (int i=0;i<numberOfInbetweenTimSteps;i++){
- double leftboundary= (double(i)/double(numberOfInbetweenTimSteps))*factorMaxGammaInfnorm*maxGammaInfNorm+minGammaInfNorm;
+ double leftboundary= (double(i)/double(numberOfInbetweenTimSteps))*factorMaxGammaInfnorm*maxGammaInfNorm;
double rightboundary= (double(i+1)/double(numberOfInbetweenTimSteps))*factorMaxGammaInfnorm*maxGammaInfNorm;
vout(2)<<" leftboundary "<<leftboundary<<" rightboundary "<<rightboundary;
if ( leftboundary< gammaInfNorm && gammaInfNorm<= rightboundary){
@@ -310,15 +309,18 @@ bool adaptiveSegregatedSolver::Step(IElphyAssemble &elphyAssemble, IElasticity &
ongoingElphyStep += 1;
elphyStepCounter=1;
}
- // Projects stretch from elphy to mech
- cardiacTransfer->Project(u_stretch(mechValues), v_stretch(elphyValues));
- mechAssemble.UpdateStretch(u_stretch(mechValues));
-
if(elphyStepCounter != oldEplpyStepCounter){
mechAssemble.ResetTime(mechAssemble.Time(), mechAssemble.LastTStep(), getNewMechDeltaTime(elphyStepCounter));
mechSolver->nonlinearSolver->InitializeMatrices(mechAssemble,u_displacement(mechValues));
+ mechSolver->nonlinearSolver->updateRHSValues(getNewMechDeltaTime(elphyStepCounter));
}
oldEplpyStepCounter=elphyStepCounter;
+
+ // Projects stretch from elphy to mech
+ cardiacTransfer->Project(u_stretch(mechValues), v_stretch(elphyValues));
+ mechAssemble.UpdateStretch(u_stretch(mechValues));
+
+
vout(2) << "Solving MechStep from " << mechAssemble.Time() << " to "
<< mechAssemble.NextTimeStep(false) << endl;
diff --git a/src/elasticity/solvers/ElasticityNewmark.cpp b/src/elasticity/solvers/ElasticityNewmark.cpp
index 637d7bfbf..8a84a297c 100644
--- a/src/elasticity/solvers/ElasticityNewmark.cpp
+++ b/src/elasticity/solvers/ElasticityNewmark.cpp
@@ -46,7 +46,6 @@ void NewmarkSolver::InitializeMatrices(const IElasticity &assemble, Vector &u){
*residualMatrix += *displacementMatrix;
}
void NewmarkSolver::operator()(const IElasticity &assemble, Vector &u) {
- vout(2)<<"NewMarkSolver operator dt "<<assemble.StepSize()<<endl;
double t = assemble.Time();
NewtonT<IElasticity>::operator()(assemble, u);
updateOldValues(assemble.StepSize(), u);
diff --git a/src/elasticity/solvers/ElasticityNewmark.hpp b/src/elasticity/solvers/ElasticityNewmark.hpp
index 5d933374b..19e24df7d 100644
--- a/src/elasticity/solvers/ElasticityNewmark.hpp
+++ b/src/elasticity/solvers/ElasticityNewmark.hpp
@@ -57,7 +57,7 @@ class NewmarkSolver : public NewtonT<IElasticity> {
NewmarkValues oldValues{};
- void updateRHSValues(double dt);
+ void updateRHSValues(double dt) override;
void updateOldValues(double dt, const Vector &u);
diff --git a/tools/solving/CoupledHoreka.py b/tools/solving/CoupledHoreka.py
index f0aa90389..e9b6c5d96 100644
--- a/tools/solving/CoupledHoreka.py
+++ b/tools/solving/CoupledHoreka.py
@@ -211,7 +211,7 @@ def startAll():
def startAdaptiveTimeStep():
lL=[2]
- tL=[0,2]
+ tL=[0]
nL=['SI']
sdt=0.0004
startJobs(lL,tL,nL,sdt,'adaptive/',{'Mesh':'TestBiventricleSimpleExcitation_smooth','ElphyProblem':'DeformedBiventricleCoarse','EndTime':0.15,'MechDeltaTime':0.0016,'AdaptiveTimeStepForMech':'true','MaximalGammaInf':0.3,'FactorMaximalGammaInf':0.8})
--
GitLab