diff --git a/dart/src/wBLRegion.cpp b/dart/src/wBLRegion.cpp
index 2caf52a4c3c78e502ee26ff1c27ed61d0a00a75b..3b394a79bea12ee2a20daa9e97b5e270543f670c 100644
--- a/dart/src/wBLRegion.cpp
+++ b/dart/src/wBLRegion.cpp
@@ -96,7 +96,6 @@ void BLRegion::SetStagBC(double Re)
std::vector<double> ClosParam = closSolver->LaminarClosures(U[0], U[1], Ret[0], Me[0], name);
- std::cout << ClosParam[3] << std::endl;
Hstar[0] = ClosParam[0];
Hstar2[0] = ClosParam[1];
Hk[0] = ClosParam[2];
@@ -168,7 +167,7 @@ double BLRegion::GetDeltaStarExt(size_t iPoint){return DeltaStarExt[iPoint];}
void BLRegion::PrintSolution(size_t iPoint)
{
- std::cout << "Point " << iPoint << std::scientific;
+ std::cout << "\nPoint " << iPoint << " \n" << std::scientific;
std::cout << "Theta " << U[iPoint * nVar + 0] << " \n" << std::scientific;
std::cout << "H " << U[iPoint * nVar + 1] << " \n" << std::scientific;
std::cout << "N " << U[iPoint * nVar + 2] << " \n" << std::scientific;
diff --git a/dart/src/wTimeSolver.cpp b/dart/src/wTimeSolver.cpp
index 66402b83734a2f42ed021d62ae7dc65e5eed8fe3..e314c7285521f7b3c0f05d72bff983d4c11e96d3 100644
--- a/dart/src/wTimeSolver.cpp
+++ b/dart/src/wTimeSolver.cpp
@@ -56,9 +56,9 @@ void TimeSolver::InitialCondition(size_t iPoint, BLRegion *bl)
if (bl->Regime[iPoint] == 1 && bl->U[iPoint * nVar + 4] == 0){
bl->U[iPoint * nVar + 4] = 0.03;
} // End if
- else if (bl->Regime[iPoint] == 0){
- bl->U[iPoint*nVar+4] = 0;
- } // End else if
+ //else if (bl->Regime[iPoint] == 0){
+ // bl->U[iPoint*nVar+4] = 0;
+ //} // End else if
} // End InitialCondition
int TimeSolver::Integration(size_t iPoint, BLRegion *bl)
@@ -97,6 +97,13 @@ int TimeSolver::Integration(size_t iPoint, BLRegion *bl)
while (innerIters < maxIt){
+ if (bl->name == "lower" && iPoint == 425)
+ {
+ std::cout << "iter " << innerIters << std::endl;
+ bl->PrintSolution(iPoint-1);
+ std::cout << (normSysRes/normSysRes0) << " CFL" << CFL << std::endl;
+ }
+
if(innerIters % itFrozenJac == 0)
{
if (nErrors >= 5)
@@ -104,6 +111,10 @@ int TimeSolver::Integration(size_t iPoint, BLRegion *bl)
ResetSolution(iPoint, bl, Sln0);
return -1; // Convergence failed
}
+ if (std::isnan(CFL))
+ {
+ CFL = 0.5;
+ }
CFL = ControlIntegration(iPoint, bl, Sln0, CFL);
itFrozenJac = itFrozenJac0;
JacMatrix = SysMatrix->ComputeJacobian(iPoint, bl, SysRes, 1e-8);
@@ -127,7 +138,7 @@ int TimeSolver::Integration(size_t iPoint, BLRegion *bl)
// CFL adaptation.
CFL = CFL0* pow(normSysRes0/normSysRes, 0.7);
- CFL = std::max(CFL, 1.0);
+ CFL = std::max(CFL, 0.3);
timeStep = SetTimeStep(CFL, dx, bl->GetVtExt(iPoint));
IMatrix = Matrix<double, 5, 5>::Identity() / timeStep;
@@ -148,6 +159,7 @@ int TimeSolver::Integration(size_t iPoint, BLRegion *bl)
}
else
{
+ ResetSolution(iPoint, bl, Sln0, true);
return innerIters;
}
@@ -181,15 +193,23 @@ double TimeSolver::ControlIntegration(size_t iPoint, BLRegion *bl, std::vector<d
/* Check for eventual errors */
- bl->U[iPoint*nVar + 0] = std::max(bl->U[iPoint*nVar + 0], 1e-7);
+ bl->U[iPoint*nVar + 0] = std::max(bl->U[iPoint*nVar + 0], 1e-6);
bl->U[iPoint*nVar + 1] = std::max(bl->U[iPoint*nVar + 1], 1.0005);
+ if (bl->U[iPoint*nVar + 3] <= 0)
+ {
+ std::cout << "Negative Ue detected at point " << iPoint << std::endl;
+ ResetSolution(iPoint, bl, Sln0, true);
+ nErrors+=1;
+ return 0.2;
+ }
+
for (size_t k=0; k<nVar; ++k)
{
if (std::isnan(bl->U[iPoint*nVar + k]))
{
ResetSolution(iPoint, bl, Sln0, true);
- itFrozenJac0 = 1; // Impose continuous Jacobian evaluation.
+ //itFrozenJac0 = 1; // Impose continuous Jacobian evaluation.
nErrors+=1;
return 0.1; // New CFL.
}
@@ -244,6 +264,12 @@ void TimeSolver::ResetSolution(size_t iPoint, BLRegion *bl, std::vector<double>
bl->Us[iPoint] = turbParam[7];
} // End else if
+ if (bl->name=="lower" && iPoint == 425)
+ {
+ std::cout << "end Reset" << std::endl;
+ bl->PrintSolution(iPoint);
+ }
+
} // End ResetSolution
diff --git a/dart/src/wViscFluxes.cpp b/dart/src/wViscFluxes.cpp
index 83cd1a78b97b418758552ba0c14c7cebf45e19d1..a22bc6e2c101b643dec131a578ede7e460821267 100644
--- a/dart/src/wViscFluxes.cpp
+++ b/dart/src/wViscFluxes.cpp
@@ -78,15 +78,15 @@ Vector<double, 5> ViscFluxes::BLlaws(size_t iPoint, BLRegion *bl, std::vector<do
if (bl->Regime[iPoint] == 0){
std::vector<double> lamParam = bl->closSolver->LaminarClosures(u[0], u[1], bl->Ret[iPoint], bl->GetMe(iPoint), bl->name);
- bl->Hstar[iPoint] = lamParam[0];
- bl->Hstar2[iPoint] = lamParam[1];
- bl->Hk[iPoint] = lamParam[2];
- bl->Cf[iPoint] = lamParam[3];
- bl->Cd[iPoint] = lamParam[4];
- bl->Delta[iPoint] = lamParam[5];
- bl->Cteq[iPoint] = lamParam[6];
- bl->Us[iPoint] = lamParam[7];
- } // End if
+ bl->Hstar[iPoint] = lamParam[0];
+ bl->Hstar2[iPoint] = lamParam[1];
+ bl->Hk[iPoint] = lamParam[2];
+ bl->Cf[iPoint] = lamParam[3];
+ bl->Cd[iPoint] = lamParam[4];
+ bl->Delta[iPoint] = lamParam[5];
+ bl->Cteq[iPoint] = lamParam[6];
+ bl->Us[iPoint] = lamParam[7];
+ } // End if
else if (bl->Regime[iPoint] == 1){
std::vector<double> turbParam = bl->closSolver->TurbulentClosures(u[0], u[1], u[4], bl->Ret[iPoint], bl->GetMe(iPoint), bl->name);
bl->Hstar[iPoint] = turbParam[0];
@@ -143,8 +143,8 @@ Vector<double, 5> ViscFluxes::BLlaws(size_t iPoint, BLRegion *bl, std::vector<do
double Cteqa = bl->Cteq[iPoint];
double Usa = bl->Us[iPoint];
- /*if (bl->name=="upper" && iPoint == 1){
- std::cout << "Reta" << Reta << std::endl;
+ /*if (bl->name=="lower" && iPoint == 2){
+ std::cout << "Reta" << bl->Ret[iPoint] << std::endl;
std::cout << "Hstara" << Hstara << std::endl;
std::cout << "Hstar2a" << Hstar2a << std::endl;
std::cout << "Hka" << Hka << std::endl;
@@ -156,7 +156,7 @@ Vector<double, 5> ViscFluxes::BLlaws(size_t iPoint, BLRegion *bl, std::vector<do
/* Space part */
- /*if(bl->name=="lower" && iPoint==1){
+ /*if(bl->name=="lower" && iPoint==2){
std::cout << "due_dx" << due_dx << std::endl;
std::cout << "c" << c << std::endl;
std::cout << "u[1]" << u[1] << std::endl;
diff --git a/dart/src/wViscSolver.cpp b/dart/src/wViscSolver.cpp
index f4b16a9061054fa71f437cbe7cfbb9e365f1518f..5ffef5261bb69d22b908e74f76a21194480fc778 100644
--- a/dart/src/wViscSolver.cpp
+++ b/dart/src/wViscSolver.cpp
@@ -169,7 +169,7 @@ int ViscSolver::Run(unsigned int couplIter){
/* Prepare time integration */
if ((couplIter == 0) || (nbElmUpper != bl[0]->GetnMarkers()) || (stagPtMvmt)){
- tSolver->SetCFL0(0.5);
+ tSolver->SetCFL0(1);
tSolver->SetitFrozenJac(1);
tSolver->SetinitSln(true);
if (couplIter != 0 && (nbElmUpper != bl[0]->GetnMarkers())){
diff --git a/dart/tests/blicpp.py b/dart/tests/blicpp.py
index 3437e134e9f0719faaacb757a7f18c888c7bd3fd..3b756977788a3a050226708880c302a04bb748dc 100644
--- a/dart/tests/blicpp.py
+++ b/dart/tests/blicpp.py
@@ -62,7 +62,7 @@ def main():
meshFactor = 2
CFL0 = 1
- plotVar = 'Hk'
+ plotVar = 'Cf'
# ========================================================================================
diff --git a/dart/vCoupler.py b/dart/vCoupler.py
index 4d218f2b211e425cfd4de04d3024b10681f6844e..a93f9ece5f9cf530d8ad0f3c90820a3cdac1db59 100644
--- a/dart/vCoupler.py
+++ b/dart/vCoupler.py
@@ -35,7 +35,7 @@ class Coupler:
self.iSolver=iSolver
self.vSolver=vSolver
- self.maxCouplIter = 100
+ self.maxCouplIter = 150
self.couplTol = 1e-4
self.tms=fwk.Timers()
@@ -81,7 +81,9 @@ class Coupler:
self.group[n].connectListNodes, self.group[n].connectListElems, dataIn[n] = self.group[n].connectList()
fMeshDict, cMeshDict, data = GroupConstructor().mergeVectors(dataIn, 0, 1)
- if ((len(fMeshDict['locCoord'][:fMeshDict['bNodes'][1]])) != (nElmUpperPrev or couplIter==0)):
+ if ((len(fMeshDict['locCoord'][:fMeshDict['bNodes'][1]]) != nElmUpperPrev) or couplIter==0):
+
+ print('STAG POINT MOVED')
# Initialize mesh
self.vSolver.InitMeshes(0, fMeshDict['locCoord'][:fMeshDict['bNodes'][1]], fMeshDict['globCoord'][:fMeshDict['bNodes'][1]])
@@ -113,11 +115,12 @@ class Coupler:
error = abs((self.vSolver.Cdt - cdPrev)/self.vSolver.Cdt) if self.vSolver.Cdt != 0 else 1
"""# Plot State.
- plotVar = 'H'
- blScalars, blVectors = viscU.ExtractVars(self.vSolver)
- xtr=self.vSolver.Getxtr()
- wData = viscU.Dictionary(blScalars, blVectors, xtr)
- viscU.PlotVars(wData, plotVar)"""
+ if couplIter >= 12:
+ plotVar = 'H'
+ blScalars, blVectors = viscU.ExtractVars(self.vSolver)
+ xtr=self.vSolver.Getxtr()
+ wData = viscU.Dictionary(blScalars, blVectors, xtr)
+ viscU.PlotVars(wData, plotVar)"""
if error <= self.couplTol and exitCode==0:
diff --git a/dart/viscous/Master/DCoupler.py b/dart/viscous/Master/DCoupler.py
index b64b9ee563b5a3fce91e4e06209f117c35c3ac1b..295cfc5be8c0fd36f36b82748101a09511ea06ce 100755
--- a/dart/viscous/Master/DCoupler.py
+++ b/dart/viscous/Master/DCoupler.py
@@ -62,7 +62,7 @@ class Coupler:
# Run inviscid solver
print('+----------------------------- Inviscid solver --------------------------------+')
self.tms['inv'].start()
- self.isolver.reset()
+ #self.isolver.reset()
self.isolver.run()
self.tms['inv'].stop()
@@ -97,6 +97,10 @@ class Coupler:
# Check convergence and output coupling iteration.
+ print('tms @ it', it, self.tms)
+ if it == 84:
+ converged = 1
+
if error <= self.tol:
print(self.tms)
print('')