|
@@ -24,7 +24,7 @@ Street, Fifth Floor, Boston, MA 02110-1301, USA. |
|
|
#include <time.h> |
|
|
#include <time.h> |
|
|
|
|
|
|
|
|
/*Cantera include files*/ |
|
|
/*Cantera include files*/ |
|
|
#include <cantera/IdealGasMix.h> |
|
|
|
|
|
|
|
|
#include "IdealGasMix.h" |
|
|
|
|
|
|
|
|
/*Sundials include files*/ |
|
|
/*Sundials include files*/ |
|
|
#include <cvode/cvode.h> /* prototypes for CVODE fcts., consts. */ |
|
|
#include <cvode/cvode.h> /* prototypes for CVODE fcts., consts. */ |
|
@@ -140,6 +140,7 @@ typedef struct { |
|
|
bool IVPSuccess; |
|
|
bool IVPSuccess; |
|
|
bool BVPSuccess; |
|
|
bool BVPSuccess; |
|
|
bool sensSuccess; |
|
|
bool sensSuccess; |
|
|
|
|
|
bool solveBVP; //The BVP does not HAVE to be solved: use this flag to skip solving it! |
|
|
int nsp,neq,npts,nvar; //key quantities required in for loops: |
|
|
int nsp,neq,npts,nvar; //key quantities required in for loops: |
|
|
//nsp-> no: of chemical species |
|
|
//nsp-> no: of chemical species |
|
|
//neq-> number of governing equations |
|
|
//neq-> number of governing equations |
|
@@ -522,21 +523,23 @@ int main(int argc, char *argv[]) |
|
|
ier=KINSetPrintLevel(mem,2); |
|
|
ier=KINSetPrintLevel(mem,2); |
|
|
|
|
|
|
|
|
/* Solve the BVP! */ |
|
|
/* Solve the BVP! */ |
|
|
ier = KINSol(mem, /* KINSol memory block */ |
|
|
|
|
|
yp, /* initial guess on input; solution vector */ |
|
|
|
|
|
KIN_NONE,//KIN_LINESEARCH,/* global strategy choice */ |
|
|
|
|
|
scale, /* scaling vector, for the variable cc */ |
|
|
|
|
|
scale); /* scaling vector for function values fval */ |
|
|
|
|
|
|
|
|
|
|
|
if (check_flag(&ier, "KINSol", 1)){ |
|
|
|
|
|
data->BVPSuccess=false; |
|
|
|
|
|
}else{ |
|
|
|
|
|
data->BVPSuccess=true; |
|
|
|
|
|
/* Get scaled norm of the system function */ |
|
|
|
|
|
ier = KINGetFuncNorm(mem, &data->ftol); |
|
|
|
|
|
//ier = check_flag(&ier, "KINGetfuncNorm", 1); |
|
|
|
|
|
printf("\nComputed solution (||F|| = %g):\n\n",data->ftol); |
|
|
|
|
|
printf("KinSOL Ignition Delay: %15.6es\n", taup(1)); |
|
|
|
|
|
|
|
|
if(data->solveBVP){ |
|
|
|
|
|
ier = KINSol(mem, /* KINSol memory block */ |
|
|
|
|
|
yp, /* initial guess on input; solution vector */ |
|
|
|
|
|
KIN_NONE,//KIN_LINESEARCH,/* global strategy choice */ |
|
|
|
|
|
scale, /* scaling vector, for the variable cc */ |
|
|
|
|
|
scale); /* scaling vector for function values fval */ |
|
|
|
|
|
|
|
|
|
|
|
if (check_flag(&ier, "KINSol", 1)){ |
|
|
|
|
|
data->BVPSuccess=false; |
|
|
|
|
|
}else{ |
|
|
|
|
|
data->BVPSuccess=true; |
|
|
|
|
|
/* Get scaled norm of the system function */ |
|
|
|
|
|
ier = KINGetFuncNorm(mem, &data->ftol); |
|
|
|
|
|
//ier = check_flag(&ier, "KINGetfuncNorm", 1); |
|
|
|
|
|
printf("\nComputed solution (||F|| = %g):\n\n",data->ftol); |
|
|
|
|
|
printf("KinSOL Ignition Delay: %15.6es\n", taup(1)); |
|
|
|
|
|
} |
|
|
} |
|
|
} |
|
|
//ier=check_flag(&ier, "KINSol", 1); |
|
|
//ier=check_flag(&ier, "KINSol", 1); |
|
|
|
|
|
|
|
@@ -549,7 +552,7 @@ int main(int argc, char *argv[]) |
|
|
|
|
|
|
|
|
/********************************************************/ |
|
|
/********************************************************/ |
|
|
/*Begin sensitivity analysis here:*/ |
|
|
/*Begin sensitivity analysis here:*/ |
|
|
if(data->BVPSuccess){ |
|
|
|
|
|
|
|
|
if(data->BVPSuccess || !data->solveBVP){ |
|
|
/* Create banded SUNMatrix for use in linear solves; bandwidths |
|
|
/* Create banded SUNMatrix for use in linear solves; bandwidths |
|
|
* are dictated by the dependence of the solution at a given |
|
|
* are dictated by the dependence of the solution at a given |
|
|
* time on one time-step ahead and one time-step behind:*/ |
|
|
* time on one time-step ahead and one time-step behind:*/ |
|
@@ -777,6 +780,8 @@ static int parseInput(UserData data, int argc, char *argv[]){ |
|
|
data->TIgnFac=0.20; |
|
|
data->TIgnFac=0.20; |
|
|
/*Sorting*/ |
|
|
/*Sorting*/ |
|
|
data->sort=false; |
|
|
data->sort=false; |
|
|
|
|
|
/* Solve two-point BVP: */ |
|
|
|
|
|
data->solveBVP=true; |
|
|
/*****************************************************/ |
|
|
/*****************************************************/ |
|
|
|
|
|
|
|
|
int ier; |
|
|
int ier; |
|
@@ -785,7 +790,7 @@ static int parseInput(UserData data, int argc, char *argv[]){ |
|
|
char comp[BUFSIZE+1]; |
|
|
char comp[BUFSIZE+1]; |
|
|
bool enteredT0, enteredP, enteredMech, enteredComp; |
|
|
bool enteredT0, enteredP, enteredMech, enteredComp; |
|
|
enteredT0=enteredP=enteredMech=enteredComp=false; |
|
|
enteredT0=enteredP=enteredMech=enteredComp=false; |
|
|
while((opt=getopt(argc,argv,"a:r:f:T:P:m:c:t:q:vsod")) != -1){ |
|
|
|
|
|
|
|
|
while((opt=getopt(argc,argv,"a:r:f:T:P:m:c:t:q:vsodn")) != -1){ |
|
|
switch(opt){ |
|
|
switch(opt){ |
|
|
case 'a': |
|
|
case 'a': |
|
|
data->atol=RCONST(atof(optarg)); |
|
|
data->atol=RCONST(atof(optarg)); |
|
@@ -834,6 +839,9 @@ static int parseInput(UserData data, int argc, char *argv[]){ |
|
|
case 'q': |
|
|
case 'q': |
|
|
data->TIgnFac=RCONST(atof(optarg)); |
|
|
data->TIgnFac=RCONST(atof(optarg)); |
|
|
break; |
|
|
break; |
|
|
|
|
|
case 'n': |
|
|
|
|
|
data->solveBVP=false; |
|
|
|
|
|
break; |
|
|
default: |
|
|
default: |
|
|
printInstructions(); |
|
|
printInstructions(); |
|
|
return(-1); |
|
|
return(-1); |
|
@@ -1484,6 +1492,7 @@ static void printInstructions(){ |
|
|
printf("\n-d :Enables manual dPdt entryspecies sensitivity output\n"); |
|
|
printf("\n-d :Enables manual dPdt entryspecies sensitivity output\n"); |
|
|
printf("\n-o :Enables sorting of ignition delay sensitivity output\n"); |
|
|
printf("\n-o :Enables sorting of ignition delay sensitivity output\n"); |
|
|
printf("\n-q :factor that determines ignition temperature (default: 0.2)\n"); |
|
|
printf("\n-q :factor that determines ignition temperature (default: 0.2)\n"); |
|
|
|
|
|
printf("\n-n :Disables the boundary value problem solver\n"); |
|
|
printf("\nexample: <executable> "); |
|
|
printf("\nexample: <executable> "); |
|
|
printf("-T 1200.0 -P 1.0 -m gri30.cti"); |
|
|
printf("-T 1200.0 -P 1.0 -m gri30.cti"); |
|
|
printf(" -c H2:1.0,N2:3.76,O2:1.0"); |
|
|
printf(" -c H2:1.0,N2:3.76,O2:1.0"); |
|
|