45 std::string DriverArpack::arpackProblemTypeLookupIds[6] = {
46 LibUtilities::SessionReader::RegisterEnumValue(
"ArpackProblemType",
"LargestReal" ,0),
47 LibUtilities::SessionReader::RegisterEnumValue(
"ArpackProblemType",
"SmallestReal" ,1),
48 LibUtilities::SessionReader::RegisterEnumValue(
"ArpackProblemType",
"LargestImag" ,2),
49 LibUtilities::SessionReader::RegisterEnumValue(
"ArpackProblemType",
"SmallestImag" ,3),
50 LibUtilities::SessionReader::RegisterEnumValue(
"ArpackProblemType",
"LargestMag" ,4),
51 LibUtilities::SessionReader::RegisterEnumValue(
"ArpackProblemType",
"SmallestMag" ,5),
53 std::string DriverArpack::arpackProblemTypeDefault = LibUtilities::SessionReader::RegisterDefaultSolverInfo(
"ArpackProblemType",
"LargestMag");
54 std::string DriverArpack::driverLookupId = LibUtilities::SessionReader::RegisterEnumValue(
"Driver",
"Arpack",0);
58 std::string DriverArpack::ArpackProblemTypeTrans[6] =
59 {
"LR",
"SR",
"LI",
"SI",
"LM",
"SM" };
95 m_equ[0]->PrintSummary(out);
98 out <<
"\tArnoldi solver type : Arpack" << endl;
100 out <<
"\tArpack problem type : ";
117 int nq =
m_equ[0]->UpdateFields()[0]->GetNcoeffs();
140 if(
m_session->DefinesFunction(
"InitialConditions"))
142 out <<
"\tInital vector : input file " << endl;
148 out <<
"\tInital vector : random " << endl;
192 std::string name =
m_session->GetSessionName() +
".evl";
193 ofstream pFile(name.c_str());
203 &v[0], n, iparam, ipntr, &workd[0],
204 &workl[0], lworkl, info);
208 out <<
"\rIteration " << cycle <<
", output: " << info <<
", ido=" << ido <<
" " << std::flush;
210 if(!((cycle-1)%m_kdim)&&(cycle> m_kdim)&&(ido!=2))
212 pFile <<
"Krylov spectrum at iteration: " << cycle << endl;
216 pFile <<
"EV Magnitude Angle Growth Frequency Residual" << endl;
220 pFile <<
"EV Real Imaginary inverse real inverse imag Residual" << endl;
224 for(
int k = 0; k <
m_kdim; ++k)
227 WriteEvs(pFile,k, workl[ipntr[5]-1+k],workl[ipntr[6]-1+k]);
231 if (ido == 99)
break;
242 m_equ[0]->TransCoeffToPhys();
250 m_equ[1]->TransCoeffToPhys();
270 m_equ[0]->TransCoeffToPhys();
273 for (
int i = 0; i < fields.num_elements(); ++i)
275 fields[i]->IProductWRTBase(fields[i]->GetPhys(),
276 fields[i]->UpdateCoeffs());
284 ASSERTL0(
false,
"Unexpected reverse communication request.");
289 out<< endl <<
"Converged in " << iparam[8] <<
" iterations" << endl;
291 ASSERTL0(info >= 0,
" Error with Dnaupd");
315 z.get(), n, sigmar, sigmai, workev.get(), &B,
317 v.get(), n, iparam, ipntr, workd.get(),
318 workl.get(),lworkl,info);
320 ASSERTL0(info == 0,
" Error with Dneupd");
335 "Need to implement Ritz re-evaluation of"
336 "eigenvalue. Only one half of complex "
337 "value will be correct");
342 out <<
"Converged Eigenvalues: " << nconv << endl;
343 pFile <<
"Converged Eigenvalues: " << nconv << endl;
347 out <<
" Magnitude Angle Growth Frequency" << endl;
348 pFile <<
" Magnitude Angle Growth Frequency" << endl;
349 for(
int i= 0; i< nconv; ++i)
354 std::string file =
m_session->GetSessionName() +
"_eig_"
355 + boost::lexical_cast<std::string>(i)
362 out <<
" Real Imaginary " << endl;
363 pFile <<
" Real Imaginary " << endl;
364 for(
int i= 0; i< nconv; ++i)
371 std::string file =
m_session->GetSessionName() +
"_eig_"
372 + boost::lexical_cast<std::string>(i)
383 for(
int j = 0; j <
m_equ[0]->GetNvariables(); ++j)
387 if (
m_comm->GetRank() == 0)
389 out <<
"L 2 error (variable " <<
m_equ[0]->GetVariable(j) <<
") : " << vL2Error << endl;
390 out <<
"L inf error (variable " <<
m_equ[0]->GetVariable(j) <<
") : " << vLinfError << endl;
Array< OneD, NekDouble > m_imag_evl
static void Dneupd(const int &rvec, const char *howmny, const int *select, double *dr, double *di, double *z, const int &ldz, const double &sigmar, const double &sigmai, double *workev, const char *bmat, const int &n, const char *which, const int &nev, const double &tol, double *resid, const int &ncv, double *v, const int &ldv, int *iparam, int *ipntr, double *workd, double *workl, const int &lworkl, int &info)
Post-processing routine to computed eigenvector of computed eigenvalues in Dnaupd.
virtual ~DriverArpack()
Destructor.
#define ASSERTL0(condition, msg)
void CopyFwdToAdj()
Copy the forward field to the adjoint system in transient growth calculations.
NekDouble m_evtol
Maxmum number of iterations.
void CopyArnoldiArrayToField(Array< OneD, NekDouble > &array)
Copy Arnoldi storage to fields.
virtual void v_Execute(std::ostream &out=std::cout)
Virtual function for solve implementation.
SOLVER_UTILS_EXPORT void ArnoldiSummary(std::ostream &out)
void CopyFieldToArnoldiArray(Array< OneD, NekDouble > &array)
Copy fields to Arnoldi storage.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
void WriteFld(std::string file, std::vector< Array< OneD, NekDouble > > coeffs)
Write coefficients to file.
virtual void v_InitObject(std::ostream &out=std::cout)
Virtual function for initialisation implementation.
static const NekDouble kNekZeroTol
int m_nvec
Dimension of Krylov subspace.
Array< OneD, NekDouble > m_real_evl
Operator in solve call is negated.
void WriteEvs(std::ostream &evlout, const int k, const NekDouble real, const NekDouble imag, NekDouble resid=NekConstants::kNekUnsetDouble, bool DumpInverse=true)
#define WARNINGL0(condition, msg)
virtual void v_InitObject(std::ostream &out=std::cout)
static void Dnaupd(int &ido, const char *bmat, const int &n, const char *which, const int &nev, const double &tol, double *resid, const int &ncv, double *v, const int &ldv, int *iparam, int *ipntr, double *workd, double *workl, const int &lworkl, int &info)
Top level reverse communication interface to solve real double-precision non-symmetric problems...
enum EvolutionOperatorType m_EvolutionOperator
Evolution Operator.
static const NekDouble kNekUnsetDouble
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
int m_nfields
interval to dump information if required.
bool m_timeSteppingAlgorithm
Period of time stepping algorithm.
static std::string ArpackProblemTypeTrans[]
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
int m_infosteps
underlying operator is time stepping
DriverFactory & GetDriverFactory()
LibUtilities::CommSharedPtr m_comm
Communication object.
Base class for the development of solvers.
LibUtilities::SessionReaderSharedPtr m_session
Session reader object.
int m_nequ
number of equations
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
int m_nits
Number of vectors to test.