westerberg/ivpNondimensional/ivpStepN.tcl

proc IVP.Integrator {qlfdid method relativeError maxNumberSolves initDeltaX maxDeltaX stopX} {

# qkfdid          stands for qualified id (this name occurs often in ASCEND
#                 tcl code) 
# method          must have value Am or Bdf and is used to set integration
#                 method
# relativeError   is the desired relative integration error, e.g., 1.0e-3.
# maxNumberSolves is max number of model solves the integrator is allowed
#                 to perform before forcibly stopping, e.g., 1000.
# initDeltaX      is the initial value to use for the stepsize for X, 
#                 e.g., 0.01 (unitless).
# maxDeltaX       is the maximum value for deltaX the integrator should use 
#                 when stepping, e.g., 2 (unitless).
# stopX           is the stopping value for X, e.g., 20 (unitless).
#
# Note: the last three items are dimensionless values (i.e., scaled by
# nominal values).

     global ascSolvStatVect;

# Solve initial point

     RUN $qlfdid.valuesForInitializing;
     RUN $qlfdid.specifyForInitializing;
     SOLVE $qlfdid.currentPt WITH QRSlv;
     if { !$ascSolvStatVect(converged) } {
         error "Initial point solve: Equations for model did not converge";
     }

# Prepare to take first integration step

     DELETE SYSTEM;
     RUN $qlfdid.valuesForStepping;
     RUN $qlfdid.specifyForStepping;

# Run one of the following to set the method

     switch $method {
         {Am}  -
         {Bdf} {}
         default {
             error "Expected Am or Bdf as method";
         }
     }
     RUN $qlfdid.setUseMethodTo$method;

# Set initial step size, integration error, stopping
# point

     ASSIGN $qlfdid.deltaX $initDeltaX;
     ASSIGN $qlfdid.maxDeltaX $maxDeltaX;
     ASSIGN $qlfdid.maxNominalSteppingError $relativeError;
     ASSIGN $qlfdid.stopX $stopX;
     
     set numberSolves          1;
     set ascStopCondHit        0;
     set ascThisIsTheFinalStep 0;

# Integrate

     while {$numberSolves < $maxNumberSolves} {

         incr numberSolves;
         set ascPolyOrder [u_getval $qlfdid.usePolyOrder];
         set ascPolyOrderValue [lindex $ascPolyOrder 0];
         set ivp_steps $ascPolyOrderValue;
         for {set ivp_i 1} {$ivp_i <= $ivp_steps} {incr ivp_i} {
#        puts "ivp_i for loop: $ivp_i $ivp_steps";
             if {[expr ($ascThisIsTheFinalStep == 0)]} {
#        puts "ascThisIsTheFinalStep is 0";
                 if {[expr ($ascStopCondHit == 0)]} {
#            puts "ascStepCondHit is 0.  Running stepX.";
                     RUN $qlfdid.stepX;
                 }
         SOLVE $qlfdid WITH QRSlv;
             if {!$ascSolvStatVect(converged)} {
                     error "Model solution $numberSolves: Equations for model did not converge.";
                 }
         }
             if  {[expr ($ascThisIsTheFinalStep == 1)]} {
                 error "Model solution $numberSolves: STOP condition reached";
         }

             RUN $qlfdid.setStopConditions;
             set ascStopCondHit [lindex [u_getval $qlfdid.stopCondHit] 0];
             set ascThisIsTheFinalStep [lindex [u_getval $qlfdid.thisIsTheFinalStep] 0];
#        puts "ascThisIsTheFinalStep $ascThisIsTheFinalStep";

         }

         RUN $qlfdid.computeMaxNominalStepsForEachVariable;

     }

     error "Model solution $numberSolves: Maximum solutions reached. Stopping.";

 }

proc getBoolean {varName} {
         # proc for getting value of boolean atoms
         qlfdid $varName
         set bval [inst atomvalue search]
         return $bval
}
proc getAttribute {varName attrName} {
         # proc for getting value of atom attributes
         # in the case of attributes with units, does not return units.
         qlfdid $varName
         set cnames "[inst child search]"
         set cindex [lsearch $cnames $attrName]
         if {$cindex < 0} {
                 error "No child named $attrName"
         }
         # children index from 1,not 0, so the list position is cindex++
         incr cindex
         set cdata [brow_child_list search $cindex VALUE]
         # modify the next lines if you want to return units also
         set cval [lindex [lindex $cdata 0] 2]
         return $cval
}

1	proc IVP.Integrator {qlfdid method relativeError maxNumberSolves initDeltaX maxDeltaX stopX} {
2
3	# qkfdid stands for qualified id (this name occurs often in ASCEND
4	# tcl code)
5	# method must have value Am or Bdf and is used to set integration
6	# method
7	# relativeError is the desired relative integration error, e.g., 1.0e-3.
8	# maxNumberSolves is max number of model solves the integrator is allowed
9	# to perform before forcibly stopping, e.g., 1000.
10	# initDeltaX is the initial value to use for the stepsize for X,
11	# e.g., 0.01 (unitless).
12	# maxDeltaX is the maximum value for deltaX the integrator should use
13	# when stepping, e.g., 2 (unitless).
14	# stopX is the stopping value for X, e.g., 20 (unitless).
15	#
16	# Note: the last three items are dimensionless values (i.e., scaled by
17	# nominal values).
18
19	global ascSolvStatVect;
20
21	# Solve initial point
22
23	RUN $qlfdid.valuesForInitializing;
24	RUN $qlfdid.specifyForInitializing;
25	SOLVE $qlfdid.currentPt WITH QRSlv;
26	if { !$ascSolvStatVect(converged) } {
27	error "Initial point solve: Equations for model did not converge";
28	}
29
30	# Prepare to take first integration step
31
32	DELETE SYSTEM;
33	RUN $qlfdid.valuesForStepping;
34	RUN $qlfdid.specifyForStepping;
35
36	# Run one of the following to set the method
37
38	switch $method {
39	{Am} -
40	{Bdf} {}
41	default {
42	error "Expected Am or Bdf as method";
43	}
44	}
45	RUN $qlfdid.setUseMethodTo$method;
46
47	# Set initial step size, integration error, stopping
48	# point
49
50	ASSIGN $qlfdid.deltaX $initDeltaX;
51	ASSIGN $qlfdid.maxDeltaX $maxDeltaX;
52	ASSIGN $qlfdid.maxNominalSteppingError $relativeError;
53	ASSIGN $qlfdid.stopX $stopX;
54
55	set numberSolves 1;
56	set ascStopCondHit 0;
57	set ascThisIsTheFinalStep 0;
58
59	# Integrate
60
61	while {$numberSolves < $maxNumberSolves} {
62
63	incr numberSolves;
64	set ascPolyOrder [u_getval $qlfdid.usePolyOrder];
65	set ascPolyOrderValue [lindex $ascPolyOrder 0];
66	set ivp_steps $ascPolyOrderValue;
67	for {set ivp_i 1} {$ivp_i <= $ivp_steps} {incr ivp_i} {
68	# puts "ivp_i for loop: $ivp_i $ivp_steps";
69	if {[expr ($ascThisIsTheFinalStep == 0)]} {
70	# puts "ascThisIsTheFinalStep is 0";
71	if {[expr ($ascStopCondHit == 0)]} {
72	# puts "ascStepCondHit is 0. Running stepX.";
73	RUN $qlfdid.stepX;
74	}
75	SOLVE $qlfdid WITH QRSlv;
76	if {!$ascSolvStatVect(converged)} {
77	error "Model solution $numberSolves: Equations for model did not converge.";
78	}
79	}
80	if {[expr ($ascThisIsTheFinalStep == 1)]} {
81	error "Model solution $numberSolves: STOP condition reached";
82	}
83
84	RUN $qlfdid.setStopConditions;
85	set ascStopCondHit [lindex [u_getval $qlfdid.stopCondHit] 0];
86	set ascThisIsTheFinalStep [lindex [u_getval $qlfdid.thisIsTheFinalStep] 0];
87	# puts "ascThisIsTheFinalStep $ascThisIsTheFinalStep";
88
89	}
90
91	RUN $qlfdid.computeMaxNominalStepsForEachVariable;
92
93	}
94
95	error "Model solution $numberSolves: Maximum solutions reached. Stopping.";
96
97	}
98
99	proc getBoolean {varName} {
100	# proc for getting value of boolean atoms
101	qlfdid $varName
102	set bval [inst atomvalue search]
103	return $bval
104	}
105	proc getAttribute {varName attrName} {
106	# proc for getting value of atom attributes
107	# in the case of attributes with units, does not return units.
108	qlfdid $varName
109	set cnames "[inst child search]"
110	set cindex [lsearch $cnames $attrName]
111	if {$cindex < 0} {
112	error "No child named $attrName"
113	}
114	# children index from 1,not 0, so the list position is cindex++
115	incr cindex
116	set cdata [brow_child_list search $cindex VALUE]
117	# modify the next lines if you want to return units also
118	set cval [lindex [lindex $cdata 0] 2]
119	return $cval
120	}
121