/[escript]/trunk/finley/test/python/rayleigh_taylor_instabilty.py
ViewVC logotype

Diff of /trunk/finley/test/python/rayleigh_taylor_instabilty.py

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

revision 876 by gross, Mon Oct 16 04:07:33 2006 UTC revision 877 by gross, Wed Oct 25 03:06:58 2006 UTC
# Line 20  import math Line 20  import math
20  ### DEFINITION OF THE DOMAIN ###  ### DEFINITION OF THE DOMAIN ###
21  l0=1.  l0=1.
22  l1=1.  l1=1.
23  n0=20  # IDEALLY 80...  n0=10  # IDEALLY 80...
24  n1=20  # IDEALLY 80...  n1=10  # IDEALLY 80...
25  mesh=esys.finley.Brick(l0=l0, l1=l1, l2=l0, order=2, n0=n0, n1=n1, n2=n0)  mesh=esys.finley.Brick(l0=l0, l1=l1, l2=l0, order=2, n0=n0, n1=n1, n2=n0)
26    
27  ### PARAMETERS OF THE SIMULATION ###  ### PARAMETERS OF THE SIMULATION ###
# Line 137  class StokesProblem(SaddlePointProblem): Line 137  class StokesProblem(SaddlePointProblem):
137        """        """
138        simple example of saddle point problem        simple example of saddle point problem
139        """        """
140        def __init__(self,domain):        def __init__(self,domain,debug=False):
141           super(StokesProblem, self).__init__(self)           super(StokesProblem, self).__init__(self,debug)
142           self.domain=domain           self.domain=domain
143           self.__pde_u=LinearPDE(domain,numEquations=self.domain.getDim(),numSolutions=self.domain.getDim())           self.__pde_u=LinearPDE(domain,numEquations=self.domain.getDim(),numSolutions=self.domain.getDim())
144           self.__pde_u.setSymmetryOn()           self.__pde_u.setSymmetryOn()
# Line 154  class StokesProblem(SaddlePointProblem): Line 154  class StokesProblem(SaddlePointProblem):
154             for j in range(self.domain.getDim()):             for j in range(self.domain.getDim()):
155               A[i,j,j,i] += self.eta               A[i,j,j,i] += self.eta
156               A[i,j,i,j] += self.eta               A[i,j,i,j] += self.eta
157           self.__pde_p.setValue(D=1./self.eta)           self.__pde_p.setValue(D=1/self.eta)
158           self.__pde_u.setValue(A=A,q=fixed_u_mask,Y=f)           self.__pde_u.setValue(A=A,q=fixed_u_mask,Y=f)
159    
160        def inner(self,p0,p1):        def inner(self,p0,p1):
# Line 172  class StokesProblem(SaddlePointProblem): Line 172  class StokesProblem(SaddlePointProblem):
172           dp=self.__pde_p.getSolution()           dp=self.__pde_p.getSolution()
173           return  dp           return  dp
174    
175    sol=StokesProblem(velocity.getDomain(),debug=True)
176  def solve_vel_uszawa(rho, eta, velocity, pressure):  def solve_vel_uszawa(rho, eta, velocity, pressure):
177  ### SOLVES THE VELOCITY PROBLEM USING A PENALTY METHOD FOR THE INCOMPRESSIBILITY ###  ### SOLVES THE VELOCITY PROBLEM USING A PENALTY METHOD FOR THE INCOMPRESSIBILITY ###
   sol=StokesProblem(velocity.getDomain())  
178    Y = Vector(0.0,Function(mesh))    Y = Vector(0.0,Function(mesh))
179    Y[1] -= rho*g    Y[1] -= rho*g
180    sol.initialize(fixed_u_mask=b_c,eta=eta,f=Y)    sol.initialize(fixed_u_mask=b_c,eta=eta,f=Y)
181    velocity,pressure=sol.solve(velocity,pressure,relaxation=1.,iter_max=50,tolerance=0.01)    velocity,pressure=sol.solve(velocity,pressure,iter_max=100,tolerance=0.01) #,accepted_reduction=None)
182    return velocity, pressure    return velocity, pressure
183        
184  def solve_vel_penalty(rho, eta, velocity, pressure):  def solve_vel_penalty(rho, eta, velocity, pressure):

Legend:
Removed from v.876  
changed lines
  Added in v.877

  ViewVC Help
Powered by ViewVC 1.1.26