Python LinearPDE.getRightHandSide方法代码示例

# 需要导入模块: from esys.escript.linearPDEs import LinearPDE [as 别名]
# 或者: from esys.escript.linearPDEs.LinearPDE import getRightHandSide [as 别名]
class SmoothScalarDistributionFromTags(ParameterSet):
    """
    creates a smooth scalar distribution on a domain from region tags
            
    :ivar domain: domain
    :type domain: `esys.escript.Domain`
    :ivar default: default value 
    :ivar tag0: tag 0
    :type tag0: ``int``
    :ivar value0: value for tag 0
    :type value0: ``float``
    :ivar tag1: tag 1
    :type tag1: ``int``
    :ivar value1: value for tag 1
    :type value1: ``float``
    :ivar tag2: tag 2
    :type tag2: ``int``
    :ivar value2: value for tag 2
    :type value2: ``float``
    :ivar tag3: tag 3
    :type tag3: ``int``
    :ivar value3: value for tag 3
    :type value3: ``float``
    :ivar tag4: tag 4
    :type tag4: ``int``
    :ivar value4: value for tag 4
    :type value4: ``float``
    :ivar tag5: tag 5
    :type tag5: ``int``
    :ivar value5: value for tag 5
    :type value5: ``float``
    :ivar tag6: tag 6
    :type tag6: ``int``
    :ivar value6: value for tag 6
    :type value6: ``float``
    :ivar tag7: tag 7
    :type tag7: ``int``
    :ivar value7: value for tag 7
    :type value7: ``float``
    :ivar tag8: tag 8
    :type tag8: ``int``
    :ivar value8: value for tag 8
    :type value8: ``float``
    :ivar tag9: tag 9
    :type tag9: ``int``
    :ivar value9: value for tag 9
    :type value9: ``float``
    """
    def __init__(self,**kwargs):
        super(SmoothScalarDistributionFromTags, self).__init__(**kwargs)
        self.declareParameter(domain=None,
                              default=0.,
                              tag0=None,
                              value0=0.,
                              tag1=None,
                              value1=0.,
                              tag2=None,
                              value2=0.,
                              tag3=None,
                              value3=0.,
                              tag4=None,
                              value4=0.,
                              tag5=None,
                              value5=0.,
                              tag6=None,
                              value6=0.,
                              tag7=None,
                              value7=0.,
                              tag8=None,
                              value8=0.,
                              tag9=None,
                              value9=0.)


    def __update(self,tag,tag_value,value):
        if self.__pde==None:
           self.__pde=LinearPDE(self.domain,numSolutions=1)
        mask=Scalar(0.,Function(self.domain))
        mask.setTaggedValue(tag,1.)
        self.__pde.setValue(Y=mask)
        mask=wherePositive(abs(self.__pde.getRightHandSide()))
        value*=(1.-mask)
        value+=tag_value*mask
        return value

    def out(self):
        """
        returns a `esys.escript.Data` object
        Link against this method to get the output of this model.
        """
        d=Scalar(self.default,Solution(self.domain)) 
        self.__pde=None
        if not self.tag0 is None: d=self.__update(self.tag0,self.value0,d)
        if not self.tag1 is None: d=self.__update(self.tag1,self.value1,d)
        if not self.tag2 is None: d=self.__update(self.tag2,self.value2,d)
        if not self.tag3 is None: d=self.__update(self.tag3,self.value3,d)
        if not self.tag4 is None: d=self.__update(self.tag4,self.value4,d)
        if not self.tag5 is None: d=self.__update(self.tag5,self.value5,d)
        if not self.tag6 is None: d=self.__update(self.tag6,self.value6,d)
        if not self.tag7 is None: d=self.__update(self.tag7,self.value7,d)
#.........这里部分代码省略.........

# 需要导入模块: from esys.escript.linearPDEs import LinearPDE [as 别名]
# 或者: from esys.escript.linearPDEs.LinearPDE import getRightHandSide [as 别名]

#.........这里部分代码省略.........
      type Tensor4 on FunctionSpace
      """
      t=escript.Tensor4(0,escript.Function(self.__domain))
      st = self.__pool.map(getStressAndTangent2D,self.__scenes)
      for i in xrange(self.__numGaussPoints):
         t.setValueOfDataPoint(i,st[i][1])
      return t
      
   def getCurrentStrain(self):
      """
      return current strain
      type: Tensor on FunctionSpace
      """
      return self.__strain
   
   def exitSimulation(self):
      """finish the whole simulation, exit"""
      self.__pool.close()

## below is modified by Hongyang Cheng ##

   def setRHS(self,X,Y=escript.Data()):
      """
      set right hande side of PDE, including X_{ij} and Y_i
      X: stress tensor at (n) time step
      Y: vector, (equivalent) body force at (n) time step
      Note that boundary force, if any, is set inhere
      """
      # apply internal stress and equivalent body force
      self.__pde.setValue(X=X, Y=Y)
      # if exterior DE domain is given
      if self.__FEDENodeMap:
         FEf = self.getFEf()
         rhs = self.__pde.getRightHandSide()
      # !!!!!! apply boundary force to the right hande side of PDE
         for FEid in FEf.keys():
            rhs_i = rhs.getTupleForDataPoint(FEid)
            rhs_i_new = [sum(f) for f in zip(rhs_i,FEf[FEid])]
            rhs.setValueOfDataPoint(FEid,rhs_i_new)
      # consider Neumann boundary conditions
      rhs -= rhs*self.__pde.getCoefficient('q')

      """
      # !!!!!! apply boundary traction to the right hand side of PDE (deprecated)
         self.__pde.setValue(X=X, Y=Y, y=self.__Nbc)
      """
   
   def solve(self,damp=.2,dynRelax=False):
      """
      solve the equation of motion in centeral difference scheme
      """
      # get initial coordinate
      x = self.getDomain().getX()
      # density matrix
      kron = util.kronecker(self.getDomain().getDim())
      rho = self.__rho*kron
      # stress at (n) time step
      stress_last = self.getCurrentStress()
      # set coefficients for ODE
      self.__pde.setValue(D=rho)
      self.setRHS(X=-stress_last)
      # solve acceleration at (n) time step from ODE
      u_tt = self.__pde.getSolution()
      # update velocity at (n+1/2) time step
      self.__u_t = 1./(2.+damp*self.__dt)*((2.-damp*self.__dt)*self.__u_t + 2.*self.__dt*u_tt)
      # get displacement increment