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

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

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

revision 1413 by trankine, Fri Jan 11 07:45:58 2008 UTC revision 1414 by gross, Thu Feb 14 10:01:43 2008 UTC
# Line 621  class ArithmeticTuple(object): Line 621  class ArithmeticTuple(object):
621             self.__items[i]+=other[i]             self.__items[i]+=other[i]
622         return self         return self
623    
624    class HomogeneousSaddlePointProblem(object):
625          """
626          This provides a framwork for solving homogeneous saddle point problem of the form
627    
628                 Av+B^*p=f
629                 Bv    =0
630    
631          for the unknowns v and p and given operators A and B and given right hand side f.
632          B^* is the adjoint operator of B is the given inner product.
633    
634          """
635          def __init__(self,**kwargs):
636            self.setTolerance()
637            self.setToleranceReductionFactor()
638    
639          def initialize(self):
640            """
641            initialize the problem (overwrite)
642            """
643            pass
644          def B(self,v):
645             """
646             returns Bv (overwrite)
647             @rtype: equal to the type of p
648    
649             @note: boundary conditions on p should be zero!
650             """
651             pass
652    
653          def inner(self,p0,p1):
654             """
655             returns inner product of two element p0 and p1  (overwrite)
656            
657             @type p0: equal to the type of p
658             @type p1: equal to the type of p
659             @rtype: C{float}
660    
661             @rtype: equal to the type of p
662             """
663             pass
664    
665          def solve_A(self,u,p):
666             """
667             solves Av=f-Au-B^*p with accuracy self.getReducedTolerance() (overwrite)
668    
669             @rtype: equal to the type of v
670             @note: boundary conditions on v should be zero!
671             """
672             pass
673    
674          def solve_prec(self,p):
675             """
676             provides a preconditioner for BA^{-1}B^* with accuracy self.getReducedTolerance() (overwrite)
677    
678             @rtype: equal to the type of p
679             """
680             pass
681    
682          def stoppingcriterium(self,Bv,v,p):
683             """
684             returns a True if iteration is terminated. (overwrite)
685    
686             @rtype: C{bool}
687             """
688             pass
689                
690          def __inner(self,p,r):
691             return self.inner(p,r[1])
692    
693          def __stoppingcriterium(self,norm_r,r,p):
694              return self.stoppingcriterium(r[1],r[0],p)
695    
696          def setTolerance(self,tolerance=1.e-8):
697                  self.__tol=tolerance
698          def getTolerance(self):
699                  return self.__tol
700          def setToleranceReductionFactor(self,reduction=0.01):
701                  self.__reduction=reduction
702          def getSubProblemTolerance(self):
703                  return self.__reduction*self.getTolerance()
704    
705          def solve(self,v,p,max_iter=20, verbose=False, show_details=False):
706                  """
707                  solves the saddle point problem using initial guesses v and p.
708    
709                  @param max_iter: maximum number of iteration steps.
710                  """
711                  self.verbose=verbose
712                  self.show_details=show_details and self.verbose
713    
714              # Az=f is solved as A(z-v)=f-Av (z-v = 0 on fixed_u_mask)      
715    
716              self.__z=v+self.solve_A(v,p*0)
717                  Bz=self.B(self.__z)
718                  #
719              #   solve BA^-1B^*p = Bz
720                  #
721                  #   note that the residual r=Bz-BA^-1B^*p = B(z-A^-1B^*p) = Bv
722                  #
723                  #   with                    Av=Az-B^*p = f - B^*p (v=z on fixed_u_mask)
724                  #                           A(v-z)=Az-B^*p-Az = f -Az - B^*p (v-z=0 on fixed_u_mask)
725                  #
726                  self.iter=0
727                  if self.verbose: print "enter PCG method (iter_max=%s)"%max_iter
728                  p,r=PCG(ArithmeticTuple(self.__z*1.,Bz),self.__Aprod,self.__Msolve,self.__inner,self.__stoppingcriterium,iter_max=max_iter, x=p)
729              return r[0],p
730    
731          def __Msolve(self,r):
732              return self.solve_prec(r[1])
733    
734          def __Aprod(self,p):
735              # return BA^-1B*p
736              #solve Av =-B^*p as Av =f-Az-B^*p
737              v=self.solve_A(self.__z,p)
738              return ArithmeticTuple(v, self.B(v))
739    
740    
741  class SaddlePointProblem(object):  class SaddlePointProblem(object):
742     """     """
743     This implements a solver for a saddlepoint problem     This implements a solver for a saddlepoint problem
# Line 875  def MaskFromBoundaryTag(function_space,* Line 992  def MaskFromBoundaryTag(function_space,*
992     else:     else:
993        return util.whereNonZero(util.interpolate(out,function_space))        return util.whereNonZero(util.interpolate(out,function_space))
994    
995    
Legend:
Removed from v.1413  
changed lines
  Added in v.1414
  ViewVC Help
Powered by ViewVC 1.1.26