NGSolve
4.9
|
00001 #ifndef FILE_PDE 00002 #define FILE_PDE 00003 00004 /*********************************************************************/ 00005 /* File: pde.hpp */ 00006 /* Author: Joachim Schoeberl */ 00007 /* Date: 08. Jul. 2000 */ 00008 /*********************************************************************/ 00009 00010 namespace ngsolve 00011 { 00012 00013 00014 class EvalVariable : public NGS_Object 00015 { 00016 private: 00017 double * variable; 00018 EvalFunction evaluator; 00019 public: 00020 EvalVariable(const MeshAccess & ama, const string & aname) 00021 : NGS_Object(ama,aname), variable(NULL) 00022 { ; } 00023 00024 void SetVariable(double & avariable) 00025 { 00026 variable = &avariable; 00027 } 00028 00029 EvalFunction & GetEvaluator(void) 00030 { 00031 return evaluator; 00032 } 00033 00034 double Evaluate(void) 00035 { 00036 if(variable) 00037 return *variable = evaluator.Eval((double*)(0)); 00038 else 00039 return evaluator.Eval((double*)(0)); 00040 } 00041 00042 virtual string GetClassName () const 00043 { 00044 return "EvalVariable"; 00045 } 00046 00047 }; 00048 00049 00050 void BuildLineIntegratorCurvePoints ( const string filename, 00051 const MeshAccess & ma, 00052 Integrator & integrator, 00053 bool draw = true); 00054 void BuildLineIntegratorCurvePoints ( istream & infile, 00055 const MeshAccess & ma, 00056 Integrator & integrator, 00057 bool draw = true); 00058 00059 00060 00064 class NGS_DLL_HEADER PDE 00065 { 00067 MeshAccess ma; 00068 00070 string geometryfilename; 00072 string meshfilename; 00074 SymbolTable<double> constants; 00076 SymbolTable<string*> string_constants; 00078 SymbolTable<double*> variables; 00080 Array<EvalVariable*> evaluators; 00082 SymbolTable<CoefficientFunction*> coefficients; 00084 SymbolTable<FESpace*> spaces; 00086 SymbolTable<GridFunction*> gridfunctions; 00088 SymbolTable<BilinearForm*> bilinearforms; 00090 SymbolTable<LinearForm*> linearforms; 00092 SymbolTable<Preconditioner*> preconditioners; 00094 SymbolTable<NumProc*> numprocs; 00095 00097 bool isgood; 00099 Array<Integrator*> CurvePointIntegrators; 00100 Array<string*> CurvePointIntegratorFilenames; 00101 00103 Array<NGS_Object*> todo; 00104 00105 00107 int levelsolved; 00108 // 00109 string matfile; 00110 00111 string filename; 00112 string workingdirectory; 00113 00114 string evaluatefiles; 00115 00116 // #ifdef SOCKETS 00117 // ClientSocketAccess sa; 00118 // #endif 00119 00121 Tcl_Interp * tcl_interpreter; 00122 00123 public: 00125 PDE(); // MeshAccess & ama); 00127 ~PDE(); 00128 00130 void LoadPDE (const string & filename, const bool nomeshload = false, const bool nogeometryload = false); 00132 void LoadPDE (istream & input, const bool nomeshload = false, const bool nogeometryload = false); 00134 void SavePDE (const string & filename); 00135 00137 void SaveSolution (const string & filename, const bool ascii = false); 00139 void LoadSolution (const string & filename, const bool ascii = false); 00141 void Solve (); 00142 // void SolveBVP () { Solve(); } 00143 00145 void PrintReport (ostream & ost); 00146 00148 void PrintMemoryUsage (ostream & ost); 00149 00151 const MeshAccess & GetMeshAccess() const { return ma; } 00152 00154 MeshAccess & GetMeshAccess() { return ma; } 00155 00156 00157 00159 bool ConstantUsed (const string & aname) const; 00161 double GetConstant (const string & aname, bool opt = 0) const; 00163 bool StringConstantUsed (const string & aname) const; 00165 string GetStringConstant (const string & aname, bool opt = 0) const; 00167 bool VariableUsed (const string & aname) const; 00168 00170 double & GetVariable (const string & aname, bool opt = 0); 00172 CoefficientFunction * GetCoefficientFunction (const string & name, bool opt = 0); 00174 FESpace * GetFESpace (const string & name, bool opt = 0); 00176 GridFunction * GetGridFunction (const string & name, bool opt = 0); 00178 BilinearForm * GetBilinearForm (const string & name, bool opt = 0); 00180 LinearForm * GetLinearForm (const string & name, bool opt = 0); 00182 Preconditioner * GetPreconditioner (const string & name, bool opt = 0); 00184 NumProc * GetNumProc (const string & name, bool opt = 0); 00185 00186 00188 const CoefficientFunction * GetCoefficientFunction (const string & name, bool opt = 0) const; 00190 const FESpace * GetFESpace (const string & name, bool opt = 0) const; 00192 const GridFunction * GetGridFunction (const string & name, bool opt = 0) const; 00194 const BilinearForm * GetBilinearForm (const string & name, bool opt = 0) const; 00196 const LinearForm * GetLinearForm (const string & name, bool opt = 0) const; 00198 const Preconditioner * GetPreconditioner (const string & name, bool opt = 0) const; 00200 const NumProc * GetNumProc (const string & name, bool opt = 0) const; 00201 00202 00203 00204 00205 00207 void AddConstant (const string & name, double val); 00209 void AddStringConstant (const string & name, const string & val); 00211 void AddVariable (const string & name, double val, int im = 5); 00213 void AddVariable (const string & name, EvalVariable * eval); 00215 void AddCoefficientFunction (const string & name, CoefficientFunction* fun); 00217 FESpace * AddFESpace (const string & name, const Flags & flags); 00219 void AddFESpace (const string & name, FESpace * space); 00221 GridFunction * AddGridFunction (const string & name, Flags & flags); 00223 void AddGridFunction (const string & name, GridFunction * gf, bool addcf = false); 00225 BilinearForm * AddBilinearForm (const string & name, Flags & flags); 00227 LinearForm * AddLinearForm (const string & name, Flags & flags); 00229 Preconditioner * AddPreconditioner (const string & name, Flags & flags); 00231 void AddNumProc (const string & name, NumProc * np); 00232 00233 00235 void AddBilinearFormIntegrator (const string & name, BilinearFormIntegrator * part, 00236 const bool deletable = true); 00237 00238 /* 00239 void AddIndependentBilinearFormIntegrator (const string & name, BilinearFormIntegrator * part, 00240 const int master, const int slave, 00241 const bool deletable = true); 00242 */ 00243 00245 void AddLinearFormIntegrator (const string & name, LinearFormIntegrator * part); 00246 00248 void SetLineIntegratorCurvePointInfo(const string & filename, 00249 Integrator * integrator); 00250 00252 SymbolTable<double> & GetConstantTable () 00253 { return constants; } 00255 SymbolTable<string*> & GetStringConstantTable () 00256 { return string_constants; } 00258 SymbolTable<double*> & GetVariableTable () 00259 { return variables; } 00261 SymbolTable<CoefficientFunction*> & GetCoefficientTable () 00262 { return coefficients; } 00264 SymbolTable<FESpace*> & GetSpaceTable () 00265 { return spaces; } 00267 SymbolTable<GridFunction*> & GetGridFunctionTable() 00268 { return gridfunctions; } 00270 SymbolTable<BilinearForm*> & GetBilinearFormTable() 00271 { return bilinearforms; } 00273 SymbolTable<LinearForm*> & GetLinearFormTable() 00274 { return linearforms; } 00276 SymbolTable<Preconditioner*> & GetPreconditionerTable () 00277 { return preconditioners; } 00279 SymbolTable<NumProc*> & GetNumProcTable() 00280 { return numprocs; } 00281 00283 bool IsGood () { return isgood; } 00285 void SetGood (bool agood) { isgood = agood; } 00286 00287 00288 00289 00290 string GetMatfile() const 00291 { return matfile;} 00292 00293 void SetMatfile(const string str) 00294 { matfile = str; } 00295 00296 string GetDirectory() const 00297 { return workingdirectory; } 00298 00299 void SetDirectory(const string str) 00300 { workingdirectory = str; } 00301 00302 string GetFilename() const 00303 { return filename; } 00304 00305 void SetFilename(const string str); 00306 00307 00308 string & GetEvaluateFiles(void) { return evaluatefiles; } 00309 // #ifdef SOCKETS 00310 // ClientSocketAccess & GetClientSocketAccess(void); 00311 // void SolveBVPClientServer (LocalHeap & lh); 00312 // private: 00313 // void CallNumProcsClientServer(const int position, LocalHeap & lh); 00314 // #endif 00315 void SetMeshFileName(const string ameshfilename) { meshfilename = ameshfilename; } 00316 void SetGeoFileName ( const string ageofilename) { geometryfilename = ageofilename; } 00317 00318 string GetMeshFileName() const { return meshfilename; } 00319 string GetGeoFileName () const { return geometryfilename; } 00320 00321 void WritePDEFile ( string abspdefile, string geofile, 00322 string meshfile, string matfile, string oldpdefile ); 00323 00324 Tcl_Interp * GetTclInterpreter() const { return tcl_interpreter; } 00325 void SetTclInterpreter(Tcl_Interp * inter) { tcl_interpreter = inter; } 00326 00327 00328 void Tcl_Eval (string str); 00329 /* 00330 void Tcl_CreateCommand(Tcl_Interp *interp, 00331 CONST char *cmdName, Tcl_CmdProc *proc, 00332 ClientData clientData, 00333 Tcl_CmdDeleteProc *deleteProc); 00334 */ 00335 00336 #ifdef ASTRID 00337 void SaveZipSolution (const string & filename, const bool ascii = false); 00339 void LoadZipSolution (const string & filename, const bool ascii = false); 00340 #endif 00341 00342 }; 00343 00344 } 00345 00346 #endif