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

revision 1486 by artak, Thu Apr 10 00:10:50 2008 UTC revision 1488 by artak, Fri Apr 11 00:22:31 2008 UTC
# Line 583  def GMRES(b, Aprod, Msolve, bilinearform Line 583  def GMRES(b, Aprod, Msolve, bilinearform
583     m=iter_restart     m=iter_restart
584     iter=0     iter=0
585     while True:     while True:
586        if iter  >= iter_max: raise MaxIterReached,"maximum number of %s steps reached."%iter_max        if iter  >= iter_max: raise MaxIterReached,"maximum number of %s steps reached"%iter_max
587        x,stopped=GMRESm(b, Aprod, Msolve, bilinearform, stoppingcriterium, x=x, iter_max=iter_max-iter, iter_restart=m)        x,stopped=GMRESm(b, Aprod, Msolve, bilinearform, stoppingcriterium, x=x, iter_max=iter_max-iter, iter_restart=m)
588        iter+=iter_restart        iter+=iter_restart
589        if stopped: break        if stopped: break
# Line 597  def GMRESm(b, Aprod, Msolve, bilinearfor Line 597  def GMRESm(b, Aprod, Msolve, bilinearfor
597     norm_b=math.sqrt(r_dot_r)     norm_b=math.sqrt(r_dot_r)
598
599     if x==None:     if x==None:
600        x=0*b        x=0
601     else:     else:
602        r=Msolve(b-Aprod(x))        r=Msolve(b-Aprod(x))
603        r_dot_r = bilinearform(r, r)        r_dot_r = bilinearform(r, r)
# Line 610  def GMRESm(b, Aprod, Msolve, bilinearfor Line 610  def GMRESm(b, Aprod, Msolve, bilinearfor
610     v=[]     v=[]
611
612     rho=math.sqrt(r_dot_r)     rho=math.sqrt(r_dot_r)
613
614     v.append(r/rho)     v.append(r/rho)
615     g[0]=rho     g[0]=rho
616
# Line 683  def GMRESm(b, Aprod, Msolve, bilinearfor Line 684  def GMRESm(b, Aprod, Msolve, bilinearfor
684     else : xhat=v[0]     else : xhat=v[0]
685
686     x += xhat     x += xhat
687     if iter!=iter_restart-1:     if iter<iter_restart-1:
688        stopped=True        stopped=True
689     else:     else:
690        stopped=False        stopped=False
1052        def getSubProblemTolerance(self):        def getSubProblemTolerance(self):
1053                return self.__reduction*self.getTolerance()                return self.__reduction*self.getTolerance()
1054
1055        def solve(self,v,p,max_iter=20, verbose=False, show_details=False, solver='PCG'):        def solve(self,v,p,max_iter=20, verbose=False, show_details=False, solver='PCG',iter_restart=10):
1056                """                """
1057                solves the saddle point problem using initial guesses v and p.                solves the saddle point problem using initial guesses v and p.
1058
1081                self.iter=0                self.iter=0
1082            if solver=='GMRES':                  if solver=='GMRES':
1083                  if self.verbose: print "enter GMRES method (iter_max=%s)"%max_iter                  if self.verbose: print "enter GMRES method (iter_max=%s)"%max_iter
1084                  p=GMRES(Bz,self.__Aprod_GMRES,self.__Msolve_GMRES,self.__inner_p,self.__stoppingcriterium_GMRES,iter_max=max_iter, x=p*1.)                  p=GMRES(Bz,self.__Aprod_GMRES,self.__Msolve_GMRES,self.__inner_p,self.__stoppingcriterium_GMRES,iter_max=max_iter, x=p*1.,iter_restart=iter_restart)
1085                  # solve Au=f-B^*p                  # solve Au=f-B^*p
1086                  #       A(u-v)=f-B^*p-Av                  #       A(u-v)=f-B^*p-Av
1087                  #       u=v+(u-v)                  #       u=v+(u-v)