address Atgeirs comments.

This commit is contained in:
Robert Kloefkorn 2014-10-20 12:32:11 +02:00
parent e552acc244
commit 42bcb5f0c8
3 changed files with 16 additions and 9 deletions

View File

@ -26,8 +26,9 @@ namespace Opm {
std::string control = param.getDefault("timestep.control", std::string("pid") ); std::string control = param.getDefault("timestep.control", std::string("pid") );
const double tol = param.getDefault("timestep.control.tol", double(1e-3) ); const double tol = param.getDefault("timestep.control.tol", double(1e-3) );
if( control == "pid" ) if( control == "pid" ) {
timeStepControl_ = TimeStepControlType( new PIDTimeStepControl( tol ) ); timeStepControl_ = TimeStepControlType( new PIDTimeStepControl( tol ) );
}
else if ( control == "pid+iteration" ) else if ( control == "pid+iteration" )
{ {
const int iterations = param.getDefault("timestep.control.targetiteration", int(25) ); const int iterations = param.getDefault("timestep.control.targetiteration", int(25) );
@ -44,8 +45,9 @@ namespace Opm {
const double time, const double timestep ) const double time, const double timestep )
{ {
// init last time step as a fraction of the given time step // init last time step as a fraction of the given time step
if( last_timestep_ < 0 ) if( last_timestep_ < 0 ) {
last_timestep_ = initial_fraction_ * timestep ; last_timestep_ = initial_fraction_ * timestep ;
}
// create adaptive step timer with previously used sub step size // create adaptive step timer with previously used sub step size
AdaptiveSimulatorTimer timer( time, time+timestep, last_timestep_ ); AdaptiveSimulatorTimer timer( time, time+timestep, last_timestep_ );
@ -85,7 +87,7 @@ namespace Opm {
// also catch linear solver not converged // also catch linear solver not converged
} }
// (linearIterations < 0 means on convergence in solver) // (linearIterations < 0 means no convergence in solver)
if( linearIterations >= 0 ) if( linearIterations >= 0 )
{ {
// advance by current dt // advance by current dt
@ -104,7 +106,7 @@ namespace Opm {
last_state = state ; last_state = state ;
last_well_state = well_state; last_well_state = well_state;
} }
else // in case of no convergence else // in case of no convergence (linearIterations < 0)
{ {
// increase restart counter // increase restart counter
if( restarts >= solver_restart_max_ ) { if( restarts >= solver_restart_max_ ) {
@ -135,8 +137,9 @@ namespace Opm {
std::cout << "Last suggested step size = " << unit::convert::to( last_timestep_, unit::day ) << " (days)" << std::endl; std::cout << "Last suggested step size = " << unit::convert::to( last_timestep_, unit::day ) << " (days)" << std::endl;
} }
if( ! std::isfinite( last_timestep_ ) ) // check for NaN if( ! std::isfinite( last_timestep_ ) ) { // check for NaN
last_timestep_ = timestep; last_timestep_ = timestep;
}
} }
} }

View File

@ -50,11 +50,13 @@ namespace Opm
assert( state.saturation().size() == satSize ); assert( state.saturation().size() == satSize );
// compute u^n - u^n+1 // compute u^n - u^n+1
for( std::size_t i=0; i<pSize; ++i ) for( std::size_t i=0; i<pSize; ++i ) {
p0_[ i ] -= state.pressure()[ i ]; p0_[ i ] -= state.pressure()[ i ];
}
for( std::size_t i=0; i<satSize; ++i ) for( std::size_t i=0; i<satSize; ++i ) {
sat0_[ i ] -= state.saturation()[ i ]; sat0_[ i ] -= state.saturation()[ i ];
}
// compute || u^n - u^n+1 || // compute || u^n - u^n+1 ||
const double stateOld = euclidianNormSquared( p0_.begin(), p0_.end() ) + const double stateOld = euclidianNormSquared( p0_.begin(), p0_.end() ) +
@ -65,8 +67,9 @@ namespace Opm
euclidianNormSquared( state.saturation().begin(), state.saturation().end() ); euclidianNormSquared( state.saturation().begin(), state.saturation().end() );
// shift errors // shift errors
for( int i=0; i<2; ++i ) for( int i=0; i<2; ++i ) {
errors_[ i ] = errors_[i+1]; errors_[ i ] = errors_[i+1];
}
// store new error // store new error
const double error = stateOld / stateNew; const double error = stateOld / stateNew;

View File

@ -60,8 +60,9 @@ namespace Opm
double euclidianNormSquared( Iterator it, const Iterator end ) const double euclidianNormSquared( Iterator it, const Iterator end ) const
{ {
double product = 0.0 ; double product = 0.0 ;
for( ; it != end; ++it ) for( ; it != end; ++it ) {
product += ( *it * *it ); product += ( *it * *it );
}
return product; return product;
} }