research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Added hooks in for filters for New step directions and for

Browse source 
new solutions.
This commit is contained in:
Harry Moffat 2010-11-10 05:01:22 +00:00
parent 813a98ddf1
commit 54e1280e08
4 changed files with 73 additions and 21 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

45
Cantera/src/numerics/NonlinearSolver.cpp
View file

@ -108,7 +108,6 @@ namespace Cantera {
m_numTotalLinearSolves(0),
m_numTotalNewtIts(0),
m_min_newt_its(0),
filterNewstep(0),
maxNewtIts_(100),
m_jacFormMethod(NSOLN_JAC_NUM),
m_nJacEval(0),
@ -178,7 +177,6 @@ namespace Cantera {
m_numTotalLinearSolves(0),
m_numTotalNewtIts(0),
m_min_newt_its(0),
filterNewstep(0),
maxNewtIts_(100),
m_jacFormMethod(NSOLN_JAC_NUM),
m_nJacEval(0),
@ -237,7 +235,6 @@ namespace Cantera {
m_numTotalLinearSolves = right.m_numTotalLinearSolves;
m_numTotalNewtIts = right.m_numTotalNewtIts;
m_min_newt_its = right.m_min_newt_its;
filterNewstep = right.filterNewstep;
maxNewtIts_ = right.maxNewtIts_;
m_jacFormMethod = right.m_jacFormMethod;
m_nJacEval = right.m_nJacEval;
@ -1261,6 +1258,12 @@ namespace Cantera {
m_normSolnFRaw = solnErrorNorm(DATA_PTR(stp), "Initial Undamped Step of the iteration", 0);
}
calcSolnToResNormVector();
/*
* Filter out bad directions
*/
doublereal normFilter = filterNewStep(time_curr, DATA_PTR(m_y_n), DATA_PTR(stp));
// Damp the Newton step
/*
@ -1328,7 +1331,7 @@ namespace Cantera {
bool m_filterIntermediate = false;
if (m_filterIntermediate) {
if (m == 0) {
(void) filterNewStep(time_n, DATA_PTR(y_new), DATA_PTR(ydot_new));
(void) filterNewSolution(time_n, DATA_PTR(y_new), DATA_PTR(ydot_new));
}
}
@ -1592,7 +1595,7 @@ namespace Cantera {
ydot[j] += dy * CJ;
}
/*
* Call the functon
* Call the function
*/
@ -1618,17 +1621,17 @@ namespace Cantera {
if (m_print_flag >= 7 && s_print_NumJac) {
if (neq_ < 30) {
printf("\t\tCurrent Matrix:\n");
printf("\t\tCurrent Matrix and Residual:\n");
printf("\t\t I,J | ");
for (j = 0; j < neq_; j++) {
printf(" %5d ", j);
}
printf("\n");
printf("| Residual \n");
printf("\t\t --");
for (j = 0; j < neq_; j++) {
printf("------------");
}
printf("\n");
printf("| -----------\n");
for (i = 0; i < neq_; i++) {
@ -1636,14 +1639,14 @@ namespace Cantera {
for (j = 0; j < neq_; j++) {
printf(" % 11.4E", J(i,j) );
}
printf("\n");
printf(" | % 11.4E\n", f[i]);
}
printf("\t\t --");
for (j = 0; j < neq_; j++) {
printf("------------");
}
printf("\n");
printf("--------------\n");
}
}
@ -1681,7 +1684,7 @@ namespace Cantera {
}
}
//====================================================================================================================
// Apply a filtering step
// Apply a filtering process to the new step
/*
* @param timeCurrent Current value of the time
* @param y_current current value of the solution
@ -1689,8 +1692,24 @@ namespace Cantera {
*
* @return Returns the norm of the value of the amount filtered
*/
doublereal NonlinearSolver::filterNewStep(const doublereal timeCurrent, doublereal * const y_current, doublereal *const ydot_current) {
return 0.0;
doublereal NonlinearSolver::filterNewStep(const doublereal timeCurrent,
const doublereal * const ybase, doublereal * const step0) {
doublereal tmp = m_func->filterNewStep(timeCurrent, ybase, step0);
return tmp;
}
//====================================================================================================================
// Apply a filtering process to the new solution
/*
* @param timeCurrent Current value of the time
* @param y_current current value of the solution
* @param ydot_current Current value of the solution derivative.
*
* @return Returns the norm of the value of the amount filtered
*/
doublereal NonlinearSolver::filterNewSolution(const doublereal timeCurrent,
doublereal * const y_current, doublereal *const ydot_current) {
doublereal tmp = m_func->filterSolnPrediction(timeCurrent, y_current);
return tmp;
}
//====================================================================================================================
// Compute the Residual Weights

18
Cantera/src/numerics/NonlinearSolver.h
View file

@ -247,7 +247,16 @@ namespace Cantera {
doublereal * const ydot, int num_newt_its);
//! Apply a filtering step
//! Apply a filtering process to the step
/*!
* @param timeCurrent Current value of the time
* @param ybase current value of the solution
*
* @return Returns the norm of the value of the amount filtered
*/
doublereal filterNewStep(const doublereal timeCurrent, const doublereal * const ybase, doublereal * const step0);
//! Apply a filter to the solution
/*!
* @param timeCurrent Current value of the time
* @param y_current current value of the solution
@ -255,8 +264,8 @@ namespace Cantera {
*
* @return Returns the norm of the value of the amount filtered
*/
doublereal filterNewStep(const doublereal timeCurrent, doublereal * const y_current, doublereal * const ydot_current);
doublereal filterNewSolution(const doublereal timeCurrent, doublereal * const y_current, doublereal * const ydot_current);
//! Return the factor by which the undamped Newton step 'step0'
//! must be multiplied in order to keep the update within the bounds of an accurate jacobian.
@ -543,9 +552,6 @@ namespace Cantera {
//! Minimum number of newton iterations to use
int m_min_newt_its;
//! Boolean that turns on solution filtering
int filterNewstep;
//! Maximum number of newton iterations
int maxNewtIts_;

17
Cantera/src/numerics/ResidJacEval.cpp
View file

@ -210,6 +210,20 @@ namespace Cantera {
}
//====================================================================================================================
// Filter the solution predictions
/*
* Codes might provide a predicted step change. This routine filters the predicted
* solution vector eliminating illegal directions.
*
* @param t Time (input)
* @param y Solution vector (input, output)
* @param step Proposed step in the solution that will be cropped
*/
doublereal ResidJacEval::filterNewStep(doublereal t, const doublereal * const ybase, doublereal * const step)
{
return 0.0;
}
//====================================================================================================================
// Filter the solution predictions
/*
* Codes might provide a predicted solution vector. This routine filters the predicted
* solution vector.
@ -217,8 +231,9 @@ namespace Cantera {
* @param t Time (input)
* @param y Solution vector (input, output)
*/
void ResidJacEval::filterSolnPrediction(doublereal t, doublereal * const y)
doublereal ResidJacEval::filterSolnPrediction(doublereal t, doublereal * const y)
{
return 0.0;
}
//====================================================================================================================
// Evalulate any stopping criteria other than a final time limit

14
Cantera/src/numerics/ResidJacEval.h
View file

@ -132,6 +132,18 @@ namespace Cantera {
*/
virtual void getInitialConditions(const doublereal t0, doublereal * const y, doublereal * const ydot);
//! Filter the solution predictions
/*!
* Codes might provide a predicted step change. This routine filters the predicted
* solution vector eliminating illegal directions.
*
* @param t Time (input)
* @param y Solution vector (input, output)
* @param step Proposed step in the solution that will be cropped
*/
virtual doublereal filterNewStep(const doublereal t, const doublereal * const ybase,
doublereal * const step);
//! Filter the solution predictions
/*!
* Codes might provide a predicted solution vector. This routine filters the predicted
@ -140,7 +152,7 @@ namespace Cantera {
* @param t Time (input)
* @param y Solution vector (input, output)
*/
virtual void filterSolnPrediction(const doublereal t, doublereal * const y);
virtual doublereal filterSolnPrediction(const doublereal t, doublereal * const y);
//! Set a global value of the absolute tolerance
/*!