1#ifndef PSCF_SWEEP_TMPL_TPP
2#define PSCF_SWEEP_TMPL_TPP
12#include <util/misc/Log.h>
21 template <
class State>
34 template <
class State>
41 template <
class State>
55 states_.allocate(historyCapacity_);
56 stateHistory_.allocate(historyCapacity_);
57 sHistory_.allocate(historyCapacity_);
58 c_.allocate(historyCapacity_);
61 template <
class State>
66 double ds = 1.0/double(
ns_);
70 Log::file() <<
"ds = " << ds << std::endl;
78 Log::file() <<
"===========================================\n";
79 Log::file() <<
"Attempt s = " << sNew << std::endl;
82 bool isContinuation =
false;
86 UTIL_THROW(
"Failure to converge initial state of sweep");
92 bool finished =
false;
102 Log::file() <<
"===========================================\n";
103 Log::file() <<
"Attempt s = " << sNew << std::endl;
114 isContinuation = reuseState_;
115 error = solve(isContinuation);
119 Log::file() <<
"Backtrack and halve sweep step size:"
128 UTIL_THROW(
"Sweep decreased ds too many times.");
138 if (sNew + ds > 1.0000001) {
142 Log::file() <<
"===========================================\n";
149 template <
class State>
161 paramType.name = name;
163 paramType.modifierPtr = &modifier;
167 template <
class State>
180 template <
class State>
183 for (
int i = 0; i < paramTypes.
size(); i++) {
198 template <
class State>
202 UTIL_CHECK(states_.capacity() == historyCapacity_);
203 UTIL_CHECK(sHistory_.capacity() == historyCapacity_);
204 UTIL_CHECK(stateHistory_.capacity() == historyCapacity_);
207 for (
int i = 0; i < historyCapacity_; ++i) {
214 for (
int i = 0; i < historyCapacity_; ++i) {
216 stateHistory_[i] = &states_[i];
220 template <
class State>
221 void SweepTmpl<State>::accept(
double sNew)
225 for (
int i = historyCapacity_ - 1; i > 0; --i) {
226 sHistory_[i] = sHistory_[i-1];
232 temp = stateHistory_[historyCapacity_-1];
233 for (
int i = historyCapacity_ - 1; i > 0; --i) {
234 stateHistory_[i] = stateHistory_[i-1];
236 stateHistory_[0] = temp;
240 if (historySize_ < historyCapacity_) {
251 template <
class State>
255 if (historySize_ == 1) {
260 for (i = 0; i < historySize_; ++i) {
263 for (j = 0; j < historySize_; ++j) {
265 num *= (sNew -
s(j));
266 den *= (
s(i) -
s(j));
286 template <
class State>
Base class allowing subclasses to define sweepable parameters.
void addParameterTypes(GArray< ParameterType > paramTypes)
Declare an array of specialized sweep parameter types.
virtual void checkAllocation(DArray< System::WField > &state)=0
std::string baseFileName_
Base name for output files.
virtual void cleanup()
Clean up operation at the end of a sweep.
void initialize()
Initialize variables that track history of solutions.
SweepTmpl(int historyCapacity)
Constructor (protected).
void setCoefficients(double sNew)
Compute coefficients of previous states for continuation.
void addParameterType(std::string name, int nId, ParameterModifier &modifier)
Declare a specialized parameter type.
virtual void sweep()
Iterate to solution.
virtual int solve(bool isContinuation)=0
virtual void readParameters(std::istream &in)
Read ns and baseFileName parameters.
GArray< ParameterType > parameterTypes_
Array of specialized parameter types.
double s(int i) const
Get the value of s for a stored solution, with i = 0 most recent.
int historyCapacity() const
Get the maximum number of stored previous states.
An automatically growable array, analogous to a std::vector.
int size() const
Return logical size of this array (i.e., current number of elements).
static std::ostream & file()
Get log ostream by reference.
ScalarParam< Type > & read(std::istream &in, const char *label, Type &value)
Add and read a new required ScalarParam < Type > object.
void setClassName(const char *className)
Set class name string.
ScalarParam< Type > & readOptional(std::istream &in, const char *label, Type &value)
Add and read a new optional ScalarParam < Type > object.
#define UTIL_CHECK(condition)
Assertion macro suitable for serial or parallel production code.
#define UTIL_THROW(msg)
Macro for throwing an Exception, reporting function, file and line number.
PSCF package top-level namespace.
Declaration of a specialized sweep parameter type.