/[escript]/trunk/escript/py_src/linearPDEs.py
ViewVC logotype

Diff of /trunk/escript/py_src/linearPDEs.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 3283 by gross, Mon Oct 18 22:39:28 2010 UTC revision 3402 by gross, Tue Dec 7 07:36:12 2010 UTC
# Line 75  class SolverOptions(object): Line 75  class SolverOptions(object):
75      :cvar JACOBI: The Jacobi preconditioner      :cvar JACOBI: The Jacobi preconditioner
76      :cvar GMRES: The Gram-Schmidt minimum residual method      :cvar GMRES: The Gram-Schmidt minimum residual method
77      :cvar PRES20: Special GMRES with restart after 20 steps and truncation after 5 residuals      :cvar PRES20: Special GMRES with restart after 20 steps and truncation after 5 residuals
78      :cvar LUMPING: Matrix lumping      :cvar ROWSUM_LUMPING: Matrix lumping using row sum
79        :cvar HRZ_LUMPING: Matrix lumping using the HRZ approach
80      :cvar NO_REORDERING: No matrix reordering allowed      :cvar NO_REORDERING: No matrix reordering allowed
81      :cvar MINIMUM_FILL_IN: Reorder matrix to reduce fill-in during factorization      :cvar MINIMUM_FILL_IN: Reorder matrix to reduce fill-in during factorization
82      :cvar NESTED_DISSECTION: Reorder matrix to improve load balancing during factorization      :cvar NESTED_DISSECTION: Reorder matrix to improve load balancing during factorization
# Line 113  class SolverOptions(object): Line 114  class SolverOptions(object):
114      GMRES= 11      GMRES= 11
115      PRES20= 12      PRES20= 12
116      LUMPING= 13      LUMPING= 13
117        ROWSUM_LUMPING= 13
118        HRZ_LUMPING= 14
119      NO_REORDERING= 17      NO_REORDERING= 17
120      MINIMUM_FILL_IN= 18      MINIMUM_FILL_IN= 18
121      NESTED_DISSECTION= 19      NESTED_DISSECTION= 19
# Line 171  class SolverOptions(object): Line 174  class SolverOptions(object):
174          self.setMinCoarseMatrixSparsity()          self.setMinCoarseMatrixSparsity()
175          self.setNumRefinements()          self.setNumRefinements()
176          self.setNumCoarseMatrixRefinements()          self.setNumCoarseMatrixRefinements()
177            self.setUsePanel()
178            self.setUseDirectInterpolation()
179            self.setDiagonalDominanceThreshold()
180                    
181    
182      def __str__(self):      def __str__(self):
# Line 203  class SolverOptions(object): Line 209  class SolverOptions(object):
209              out+="\nApply preconditioner locally = %s"%self.useLocalPreconditioner()              out+="\nApply preconditioner locally = %s"%self.useLocalPreconditioner()
210              if self.getPreconditioner() == self.AMG:              if self.getPreconditioner() == self.AMG:
211                  out+="\nMaximum number of levels = %s"%self.getLevelMax()                  out+="\nMaximum number of levels = %s"%self.getLevelMax()
212                  out+="\nCoarsening method = %s"%self.getName(self.getCoarsening())                  out+="\nCoarsening threshold = %e"%self.getCoarseningThreshold()
213                  out+="\nCoarsening threshold = %e"%self.getMinCoarseMatrixSize()                  out+="\nMinimal sparsity on coarsest level = %e"%(self.getMinCoarseMatrixSparsity(),)
214                  out+="\nSmoother = %s"%self.getName(self.getSmoother())                  out+="\nSmoother = %s"%self.getName(self.getSmoother())
215                  out+="\nMinimum size of the coarsest level matrix = %e"%self.getCoarseningThreshold()                  out+="\nMinimum size of the coarsest level matrix = %e"%self.getCoarseningThreshold()
216                  out+="\nNumber of pre / post sweeps = %s / %s, %s"%(self.getNumPreSweeps(), self.getNumPostSweeps(), self.getNumSweeps())                  out+="\nNumber of pre / post sweeps = %s / %s, %s"%(self.getNumPreSweeps(), self.getNumPostSweeps(), self.getNumSweeps())
                 out+="\nMinimal sparsity on coarsest level = %e"%(self.getMinCoarseMatrixSparsity(),)  
217                  out+="\nNumber of refinement steps in coarsest level solver = %d"%(self.getNumCoarseMatrixRefinements(),)                  out+="\nNumber of refinement steps in coarsest level solver = %d"%(self.getNumCoarseMatrixRefinements(),)
218                                    out+="\nUse node panel = %s"%(self.usePanel())
219                    out+="\nUse direct interpolation only = %s"%(self.useDirectInterpolation())
220                    out+="\nThreshold for diagonal dominant rows = %s"%(setDiagonalDominanceThreshold(),)
221    
222              if self.getPreconditioner() == self.AMLI:              if self.getPreconditioner() == self.AMLI:
223                  out+="\nMaximum number of levels = %s"%self.getLevelMax()                  out+="\nMaximum number of levels = %s"%self.getLevelMax()
224                  out+="\nCoarsening method = %s"%self.getName(self.getCoarsening())                  out+="\nCoarsening method = %s"%self.getName(self.getCoarsening())
# Line 244  class SolverOptions(object): Line 252  class SolverOptions(object):
252          if key == self.JACOBI: return "JACOBI"          if key == self.JACOBI: return "JACOBI"
253          if key == self.GMRES: return "GMRES"          if key == self.GMRES: return "GMRES"
254          if key == self.PRES20: return "PRES20"          if key == self.PRES20: return "PRES20"
255          if key == self.LUMPING: return "LUMPING"          if key == self.ROWSUM_LUMPING: return "ROWSUM_LUMPING"
256            if key == self.HRZ_LUMPING: return "HRZ_LUMPING"
257          if key == self.NO_REORDERING: return "NO_REORDERING"          if key == self.NO_REORDERING: return "NO_REORDERING"
258          if key == self.MINIMUM_FILL_IN: return "MINIMUM_FILL_IN"          if key == self.MINIMUM_FILL_IN: return "MINIMUM_FILL_IN"
259          if key == self.NESTED_DISSECTION: return "NESTED_DISSECTION"          if key == self.NESTED_DISSECTION: return "NESTED_DISSECTION"
# Line 326  class SolverOptions(object): Line 335  class SolverOptions(object):
335              self.__converged = (value == True)              self.__converged = (value == True)
336          if name == "time_step_backtracking_used":          if name == "time_step_backtracking_used":
337              self.__time_step_backtracking_used = (value == True)              self.__time_step_backtracking_used = (value == True)
338            if name == "coarse_level_sparsity":
339               self.__coarse_level_sparsity = value
340        if name == "num_coarse_unknowns":
341               self.__num_coarse_unknowns= value
342      def getDiagnostics(self, name):      def getDiagnostics(self, name):
343          """          """
344          Returns the diagnostic information ``name``. Possible values are:          Returns the diagnostic information ``name``. Possible values are:
# Line 342  class SolverOptions(object): Line 355  class SolverOptions(object):
355          - "net_time": net execution time, excluding setup time for the solver and execution time for preconditioner          - "net_time": net execution time, excluding setup time for the solver and execution time for preconditioner
356          - "cum_net_time": cumulative net execution time          - "cum_net_time": cumulative net execution time
357          - "preconditioner_size": size of preconditioner [Bytes]          - "preconditioner_size": size of preconditioner [Bytes]
358          - "converged": return self.__converged          - "converged": return True if solution has converged.
359          - "time_step_backtracking_used": return self.__converged          - "time_step_backtracking_used": returns True if time step back tracking has been used.
360            - "coarse_level_sparsity": returns the sparsity of the matrix on the coarsest level
361            - "num_coarse_unknowns": returns the number of unknowns on the coarsest level
362            
363                    
364          :param name: name of diagnostic information to return          :param name: name of diagnostic information to return
# Line 366  class SolverOptions(object): Line 381  class SolverOptions(object):
381          elif name == "converged": return self.__converged                elif name == "converged": return self.__converged      
382          elif name == "preconditioner_size": return  self.__preconditioner_size          elif name == "preconditioner_size": return  self.__preconditioner_size
383          elif name == "time_step_backtracking_used": return  self.__time_step_backtracking_used          elif name == "time_step_backtracking_used": return  self.__time_step_backtracking_used
384            elif name == "coarse_level_sparsity": return self.__coarse_level_sparsity
385            elif name == "num_coarse_unknowns": return self.__num_coarse_unknowns
386          else:          else:
387              raise ValueError,"unknown diagnostic item %s"%name              raise ValueError,"unknown diagnostic item %s"%name
388      def hasConverged(self):      def hasConverged(self):
# Line 468  class SolverOptions(object): Line 485  class SolverOptions(object):
485          :param method: key of the solver method to be used.          :param method: key of the solver method to be used.
486          :type method: in `SolverOptions.DEFAULT`, `SolverOptions.DIRECT`, `SolverOptions.CHOLEVSKY`, `SolverOptions.PCG`,          :type method: in `SolverOptions.DEFAULT`, `SolverOptions.DIRECT`, `SolverOptions.CHOLEVSKY`, `SolverOptions.PCG`,
487                          `SolverOptions.CR`, `SolverOptions.CGS`, `SolverOptions.BICGSTAB`,                          `SolverOptions.CR`, `SolverOptions.CGS`, `SolverOptions.BICGSTAB`,
488                          `SolverOptions.GMRES`, `SolverOptions.PRES20`, `SolverOptions.LUMPING`, `SolverOptions.ITERATIVE`,                          `SolverOptions.GMRES`, `SolverOptions.PRES20`, `SolverOptions.ROWSUM_LUMPING`, `SolverOptions.HRZ_LUMPING`, `SolverOptions.ITERATIVE`,
489                          `SolverOptions.NONLINEAR_GMRES`, `SolverOptions.TFQMR`, `SolverOptions.MINRES`                          `SolverOptions.NONLINEAR_GMRES`, `SolverOptions.TFQMR`, `SolverOptions.MINRES`
490          :note: Not all packages support all solvers. It can be assumed that a package makes a reasonable choice if it encounters an unknown solver method.          :note: Not all packages support all solvers. It can be assumed that a package makes a reasonable choice if it encounters an unknown solver method.
491          """          """
492      if method==None: method=0      if method==None: method=0
493          if not method in [ SolverOptions.DEFAULT, SolverOptions.DIRECT, SolverOptions.CHOLEVSKY, SolverOptions.PCG,          if not method in [ SolverOptions.DEFAULT, SolverOptions.DIRECT, SolverOptions.CHOLEVSKY, SolverOptions.PCG,
494                             SolverOptions.CR, SolverOptions.CGS, SolverOptions.BICGSTAB,                             SolverOptions.CR, SolverOptions.CGS, SolverOptions.BICGSTAB,
495                             SolverOptions.GMRES, SolverOptions.PRES20, SolverOptions.LUMPING, SolverOptions.ITERATIVE,                             SolverOptions.GMRES, SolverOptions.PRES20, SolverOptions.ROWSUM_LUMPING, SolverOptions.HRZ_LUMPING,
496                               SolverOptions.ITERATIVE,
497                             SolverOptions.NONLINEAR_GMRES, SolverOptions.TFQMR, SolverOptions.MINRES ]:                             SolverOptions.NONLINEAR_GMRES, SolverOptions.TFQMR, SolverOptions.MINRES ]:
498               raise ValueError,"unknown solver method %s"%method               raise ValueError,"unknown solver method %s"%method
499          self.__method=method          self.__method=method
# Line 485  class SolverOptions(object): Line 503  class SolverOptions(object):
503    
504          :rtype: in the list `SolverOptions.DEFAULT`, `SolverOptions.DIRECT`, `SolverOptions.CHOLEVSKY`, `SolverOptions.PCG`,          :rtype: in the list `SolverOptions.DEFAULT`, `SolverOptions.DIRECT`, `SolverOptions.CHOLEVSKY`, `SolverOptions.PCG`,
505                          `SolverOptions.CR`, `SolverOptions.CGS`, `SolverOptions.BICGSTAB`,                          `SolverOptions.CR`, `SolverOptions.CGS`, `SolverOptions.BICGSTAB`,
506                          `SolverOptions.GMRES`, `SolverOptions.PRES20`, `SolverOptions.LUMPING`, `SolverOptions.ITERATIVE`,                          `SolverOptions.GMRES`, `SolverOptions.PRES20`, `SolverOptions.ROWSUM_LUMPING`, `SolverOptions.HRZ_LUMPING`, `SolverOptions.ITERATIVE`,
507                          `SolverOptions.NONLINEAR_GMRES`, `SolverOptions.TFQMR`, `SolverOptions.MINRES`                          `SolverOptions.NONLINEAR_GMRES`, `SolverOptions.TFQMR`, `SolverOptions.MINRES`
508          """          """
509          return self.__method          return self.__method
# Line 524  class SolverOptions(object): Line 542  class SolverOptions(object):
542          """          """
543          Returns the key of the reordering method to be applied if supported by the solver.          Returns the key of the reordering method to be applied if supported by the solver.
544    
545          :rtype ordering: in 'SolverOptions.NO_REORDERING', 'SolverOptions.MINIMUM_FILL_IN', 'SolverOptions.NESTED_DISSECTION', 'SolverOptions.DEFAULT_REORDERING'          :rtype: ordering in 'SolverOptions.NO_REORDERING', 'SolverOptions.MINIMUM_FILL_IN', 'SolverOptions.NESTED_DISSECTION', 'SolverOptions.DEFAULT_REORDERING'
546          """          """
547          return self.__reordering          return self.__reordering
548      def setRestart(self,restart=None):      def setRestart(self,restart=None):
# Line 559  class SolverOptions(object): Line 577  class SolverOptions(object):
577              return -1              return -1
578              else:              else:
579              return r              return r
580      
581        def setDiagonalDominanceThreshold(self,value=0.5):
582            """
583            Sets the threshold for diagonal dominant rows which are eliminated during AMG coarsening.
584    
585         :param value: threshold
586         :type value: ``float``
587            """
588            value=float(value)
589            if value<0 or value>1.:
590           raise ValueError,"Diagonal dominance threshold must be between 0 and 1."
591        self.__diagonal_dominance_threshold=value
592            
593        def getDiagonalDominanceThreshold(self):
594            """
595            Returns the threshold for diagonal dominant rows which are eliminated during AMG coarsening.
596    
597            :rtype: ``float``
598            """
599            return self.__diagonal_dominance_threshold
600            
601      def setTruncation(self,truncation=20):      def setTruncation(self,truncation=20):
602          """          """
603          Sets the number of residuals in GMRES to be stored for orthogonalization.  The more residuals are stored          Sets the number of residuals in GMRES to be stored for orthogonalization.  The more residuals are stored
# Line 900  class SolverOptions(object): Line 939  class SolverOptions(object):
939          self.__adapt_inner_tolerance=False          self.__adapt_inner_tolerance=False
940      def setInnerToleranceAdaption(self,adapt=True):      def setInnerToleranceAdaption(self,adapt=True):
941          """          """
942          Sets a flag to indicate automatic selection of the inner tolerance.          Sets the flag to indicate automatic selection of the inner tolerance.
943    
944          :param adapt: If ``True``, the inner tolerance is selected automatically.          :param adapt: If ``True``, the inner tolerance is selected automatically.
945          :type adapt: ``bool``          :type adapt: ``bool``
# Line 936  class SolverOptions(object): Line 975  class SolverOptions(object):
975          self.__accept_convergence_failure=False          self.__accept_convergence_failure=False
976      def setAcceptanceConvergenceFailure(self,accept=False):      def setAcceptanceConvergenceFailure(self,accept=False):
977          """          """
978          Sets a flag to indicate the acceptance of a failure of convergence.          Sets the flag to indicate the acceptance of a failure of convergence.
979    
980          :param accept: If ``True``, any failure to achieve convergence is accepted.          :param accept: If ``True``, any failure to achieve convergence is accepted.
981          :type accept: ``bool``          :type accept: ``bool``
# Line 959  class SolverOptions(object): Line 998  class SolverOptions(object):
998    
999      def setLocalPreconditionerOn(self):      def setLocalPreconditionerOn(self):
1000          """          """
1001          Sets a flag to use  local preconditioning to on          Sets the flag to use  local preconditioning to on
1002          """          """
1003          self.__use_local_preconditioner=True          self.__use_local_preconditioner=True
1004      def setLocalPreconditionerOff(self):      def setLocalPreconditionerOff(self):
1005          """          """
1006          Sets a flag to use  local preconditioning to off          Sets the flag to use  local preconditioning to off
1007          """          """
1008          self.__use_local_preconditioner=False          self.__use_local_preconditioner=False
1009            
1010      def setLocalPreconditioner(self,use=False):      def setLocalPreconditioner(self,use=False):
1011          """          """
1012          Sets a flag to use  local preconditioning          Sets the flag to use  local preconditioning
1013    
1014          :param accept: If ``True``, local proconditioning on each MPI rank is applied          :param use: If ``True``, local proconditioning on each MPI rank is applied
1015          :type accept: ``bool``          :type use: ``bool``
1016          """          """
1017          if use:          if use:
1018              self.setUseLocalPreconditionerOn()              self.setUseLocalPreconditionerOn()
# Line 1060  class SolverOptions(object): Line 1100  class SolverOptions(object):
1100        :rtype: non-negative ``int``        :rtype: non-negative ``int``
1101        """        """
1102        return self.__coarse_refinements        return self.__coarse_refinements
1103          
1104        def usePanel(self):
1105            """
1106            Returns ``True`` if a panel is used to serach for unknown in the AMG coarsening, The panel approach is normally faster
1107            but can lead to larger coarse level systems.
1108            
1109            :return: ``True`` if a panel is used to find unknowns in AMG coarsening
1110            :rtype: ``bool``
1111            """
1112            return self.__use_panel
1113    
1114        def setUsePanelOn(self):
1115            """
1116            Sets the flag to use a panel to find unknowns in AMG coarsening
1117            """
1118            self.__use_panel=True
1119            
1120        def setUsePanelOff(self):
1121            """
1122            Sets the flag to use a panel to find unknowns in AMG coarsening to off
1123            """
1124            self.__use_panel=False
1125            
1126        def setUsePanel(self,use=True):
1127            """
1128            Sets the flag to use  a panel to find unknowns in AMG coarsening
1129    
1130            :param use: If ``True``,a panel is used to find unknowns in AMG coarsening
1131            :type use: ``bool``
1132            """
1133            if use:
1134                self.setUsePanelOn()
1135            else:
1136                self.setUsePanelOff()
1137                
1138        def useDirectInterpolation(self):
1139            """
1140            Returns ``True`` if direct interpolation is used in AMG.
1141    
1142        :return: ``True`` if direct interpolation is used in AMG.
1143            :rtype: ``bool``
1144            """
1145            return self.__use_direct_interpolation
1146    
1147        def setUseDirectInterpolationOn(self):
1148            """
1149            Sets the flag to use direct interpolation in AMG
1150            """
1151            self.__use_direct_interpolation=True
1152            
1153        def setUseDirectInterpolationOff(self):
1154            """
1155            Sets the flag to use direct interpolation in AMG
1156            """
1157            self.__use_direct_interpolation=False
1158            
1159        def setUseDirectInterpolation(self,use=False):
1160            """
1161            Sets the flag to use direct interpolation in AMG
1162    
1163        :param use: If ``True``, direct interpolation in AMG
1164        :type use: ``bool``
1165            """
1166            if use:
1167                self.setUseDirectInterpolationOn()
1168            else:
1169                self.setUseDirectInterpolationOff()
1170    
1171    
1172  class IllegalCoefficient(ValueError):  class IllegalCoefficient(ValueError):
1173     """     """
1174     Exception that is raised if an illegal coefficient of the general or     Exception that is raised if an illegal coefficient of the general or
# Line 1664  class LinearProblem(object): Line 1772  class LinearProblem(object):
1772        :return: True if the current solver method is lumping        :return: True if the current solver method is lumping
1773        :rtype: ``bool``        :rtype: ``bool``
1774        """        """
1775        return self.getSolverOptions().getSolverMethod()==self.getSolverOptions().LUMPING        return self.getSolverOptions().getSolverMethod() in [ SolverOptions.ROWSUM_LUMPING, SolverOptions.HRZ_LUMPING ]
1776     # ==========================================================================     # ==========================================================================
1777     #    symmetry  flag:     #    symmetry  flag:
1778     # ==========================================================================     # ==========================================================================
# Line 2543  class LinearPDE(LinearProblem): Line 2651  class LinearPDE(LinearProblem):
2651        """        """
2652        Returns the system type which needs to be used by the current set up.        Returns the system type which needs to be used by the current set up.
2653        """        """
2654        solver_options=self.getSolverOptions()        if self.isUsingLumping():
2655        return self.getDomain().getSystemMatrixTypeId(solver_options.getSolverMethod(), solver_options.getPreconditioner(),solver_options.getPackage(), solver_options.isSymmetric())       return "__ESCRIPT_DATA"
2656          else:
2657         solver_options=self.getSolverOptions()
2658         return self.getDomain().getSystemMatrixTypeId(solver_options.getSolverMethod(), solver_options.getPreconditioner(),solver_options.getPackage(), solver_options.isSymmetric())
2659    
2660     def checkSymmetry(self,verbose=True):     def checkSymmetry(self,verbose=True):
2661        """        """
# Line 2595  class LinearPDE(LinearProblem): Line 2706  class LinearPDE(LinearProblem):
2706         if not self.isSolutionValid():         if not self.isSolutionValid():
2707            mat,f=self.getSystem()            mat,f=self.getSystem()
2708            if self.isUsingLumping():            if self.isUsingLumping():
2709             if not util.inf(abs(mat)) > 0.:
2710             raise ZeroDivisionError,"Lumped mass matrix as zero entry (try order 1 elements or HRZ lumping)."
2711               self.setSolution(f*1/mat)               self.setSolution(f*1/mat)
2712            else:            else:
2713               self.trace("PDE is resolved.")               self.trace("PDE is resolved.")
# Line 2668  class LinearPDE(LinearProblem): Line 2781  class LinearPDE(LinearProblem):
2781                   self.resetOperator()                   self.resetOperator()
2782                   operator=self.getCurrentOperator()                   operator=self.getCurrentOperator()
2783                   if hasattr(self.getDomain(), "addPDEToLumpedSystem") :                   if hasattr(self.getDomain(), "addPDEToLumpedSystem") :
2784                      self.getDomain().addPDEToLumpedSystem(operator, D_times_e, d_times_e)              hrz_lumping=( self.getSolverOptions().getSolverMethod() ==  SolverOptions.HRZ_LUMPING )
2785                      self.getDomain().addPDEToLumpedSystem(operator, D_reduced_times_e, d_reduced_times_e)              self.getDomain().addPDEToLumpedSystem(operator, D_times_e, d_times_e,  hrz_lumping )
2786                self.getDomain().addPDEToLumpedSystem(operator, D_reduced_times_e, d_reduced_times_e, hrz_lumping)
2787                   else:                   else:
2788                      self.getDomain().addPDEToRHS(operator, \                      self.getDomain().addPDEToRHS(operator, \
2789                                                   escript.Data(), \                                                   escript.Data(), \

Legend:
Removed from v.3283  
changed lines
  Added in v.3402

  ViewVC Help
Powered by ViewVC 1.1.26