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

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

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

revision 1467 by gross, Wed Apr 2 08:10:37 2008 UTC revision 1551 by artak, Wed May 7 23:11:44 2008 UTC
# Line 36  __date__="$Date:$" Line 36  __date__="$Date:$"
36  from escript import *  from escript import *
37  import util  import util
38  from linearPDEs import LinearPDE  from linearPDEs import LinearPDE
39  from pdetools import HomogeneousSaddlePointProblem  from pdetools import HomogeneousSaddlePointProblem,Projector
40    
41  class StokesProblemCartesian(HomogeneousSaddlePointProblem):  class StokesProblemCartesian(HomogeneousSaddlePointProblem):
42        """        """
# Line 82  class StokesProblemCartesian(Homogeneous Line 82  class StokesProblemCartesian(Homogeneous
82          for j in range(self.domain.getDim()):          for j in range(self.domain.getDim()):
83              A[i,j,j,i] += 1.              A[i,j,j,i] += 1.
84              A[i,j,i,j] += 1.              A[i,j,i,j] += 1.
85      self.__pde_prec.setValue(D=1./self.eta)      self.__pde_prec.setValue(D=1/eta) #1./self.eta
86          self.__pde_u.setValue(A=A*self.eta,q=fixed_u_mask,Y=f,y=surface_stress)          self.__pde_u.setValue(A=A*self.eta,q=fixed_u_mask,Y=f,y=surface_stress)
87    
88        def B(self,arg):        def B(self,arg):
# Line 96  class StokesProblemCartesian(Homogeneous Line 96  class StokesProblemCartesian(Homogeneous
96           s1=util.interpolate(p1,Function(self.domain))           s1=util.interpolate(p1,Function(self.domain))
97           return util.integrate(s0*s1)           return util.integrate(s0*s1)
98    
99          def inner_a(self,a0,a1):
100             p0=util.interpolate(a0[1],Function(self.domain))
101             p1=util.interpolate(a1[1],Function(self.domain))
102             alfa=(1/self.vol)*util.integrate(p0)
103             beta=(1/self.vol)*util.integrate(p1)
104         v0=util.grad(a0[0])
105         v1=util.grad(a1[0])
106         #print "NORM",alfa,beta,util.integrate((p0-alfa)*(p1-beta))+util.integrate(util.inner(v0,v1))
107             return util.integrate((p0-alfa)*(p1-beta)+((1/self.eta)**2)*util.inner(v0,v1))
108    
109    
110        def getStress(self,u):        def getStress(self,u):
111           mg=util.grad(u)           mg=util.grad(u)
112           return 2.*self.eta*util.symmetric(mg)           return 2.*self.eta*util.symmetric(mg)
# Line 105  class StokesProblemCartesian(Homogeneous Line 116  class StokesProblemCartesian(Homogeneous
116           solves Av=f-Au-B^*p (v=0 on fixed_u_mask)           solves Av=f-Au-B^*p (v=0 on fixed_u_mask)
117           """           """
118           self.__pde_u.setTolerance(self.getSubProblemTolerance())           self.__pde_u.setTolerance(self.getSubProblemTolerance())
119           self.__pde_u.setValue(X=-self.getStress(u)+p*util.kronecker(self.domain))           self.__pde_u.setValue(X=-self.getStress(u)-p*util.kronecker(self.domain))
120           return  self.__pde_u.getSolution(verbose=self.show_details)           return  self.__pde_u.getSolution(verbose=self.show_details)
121    
122    
123        def solve_prec(self,p):        def solve_prec(self,p):
124         #proj=Projector(domain=self.domain, reduce = True, fast=False)
125           self.__pde_prec.setTolerance(self.getSubProblemTolerance())           self.__pde_prec.setTolerance(self.getSubProblemTolerance())
126           self.__pde_prec.setValue(Y=p)           self.__pde_prec.setValue(Y=p)
127           q=self.__pde_prec.getSolution(verbose=self.show_details)           q=self.__pde_prec.getSolution(verbose=self.show_details)
128         q0=util.interpolate(q,Function(self.domain))
129             q-=(1/self.vol)*util.integrate(q0)
130           return q           return q
131    
132        def stoppingcriterium(self,Bv,v,p):        def stoppingcriterium(self,Bv,v,p):
133            n_r=util.sqrt(self.inner(Bv,Bv))            n_r=util.sqrt(self.inner(Bv,Bv))
134            n_v=util.Lsup(v)            n_v=util.Lsup(v)
# Line 123  class StokesProblemCartesian(Homogeneous Line 139  class StokesProblemCartesian(Homogeneous
139                return True                return True
140            else:            else:
141                return False                return False
142        def stoppingcriterium_GMRES(self,norm_r,norm_b):        def stoppingcriterium2(self,norm_r,norm_b,solver='GMRES',TOL=None):
143            if self.verbose: print "GMRES step %s: L2(r) = %s, L2(b)*TOL=%s"%(self.iter,norm_r,norm_b*self.getTolerance())        if TOL==None:
144                 TOL=self.getTolerance()
145              if self.verbose: print "%s step %s: L2(r) = %s, L2(b)*TOL=%s"%(solver,self.iter,norm_r,norm_b*TOL)
146            self.iter+=1            self.iter+=1
147            if norm_r <= norm_b*self.getTolerance():            
148                if self.verbose: print "GMRES terminated after %s steps."%self.iter            if norm_r <= norm_b*TOL:
149                  if self.verbose: print "%s terminated after %s steps."%(solver,self.iter)
150                return True                return True
151            else:            else:
152                return False                return False
153    
154    

Legend:
Removed from v.1467  
changed lines
  Added in v.1551

  ViewVC Help
Powered by ViewVC 1.1.26