### module m_procBIE ''' Contains classes providing the solution of boundary integral equations (mainly using the fortran subroutines) Modules: bie (precompiled *.so module), numpy Calsses: solution ''' import bie from numpy import array,linspace,zeros from numpy.linalg import solve class solution(object): ''' Provides the solution of the boundary integral equations Inputs: Parameters from the module m_preBEM - iplane,el1,el2,pn1,pn2,elements1,elements2 Parameters: iplane,lements1,elements2 - see module m_prepBEM2 el1,el2,pn1,pn2,gl1,gl2 - [float] Young's modulus, Poissons' ratio and shear modulus el,pn,gl - [float - ndarray(2)] Young's modulus,Poissons'ratio and shear modulus recalculated for plane stress/strain nmax - [int] greater of the elements No. of the subdomains n - [int - ndarray(2)] No. of elements of subdomain1/2 x1,y1 - [float - ndarray(2,nmax)] coordinates of the first endpoint of the element x2,y2 - [float - ndarray(2,nmax)] coordinates of the second endpoint of the element xm,ym - [float - ndarray(2,nmax)] coordinates of the midpoint (node) of the element ub,vb - [float - ndarray(2,nmax)] coordinates of the boundary conditions at the node of the element pos - [int - ndarray(2,nmax)] position of the node in the G and H matrix of the subdomain1/2 kcode - [int - ndarray(2,nmax)] kind of the prescribed boundary condition at the node of the subdomain1/2 1 - uu, 2 - tt, 3 - ut, 4 - tu, 12 - 12 G,H - [float - ndarray(2,2*nmax,2*nmax)] G and H matrix calculated via the fortran subroutines gmatrel a hmatrel A - [float - ndarray(4*nmax,4*nmax)] matrix of the solved algebraic system composed of the matrices G and H r - [float - ndarray(4*nmax)] vector calculated from the matrices G and H and prescribed boundary conditions sol - [float - ndarray(4*nmax)] - solution of the algebraic system system A*sol=r u,v,tx,ty - [float - ndarray(2,nmax)] - displacement and traction coordinates at the nodes Methods: set_param() - reshapes and sets out the parameters, which enter to the fortran subroutines G_calc() - provides calculation of the matrix G via fortran subroutine gmatrel H_calc() - provides calculation of the matrix H via fortran subroutine hmatrel Ar_sol() - evaluates values of displacements and tranctions at the element nodes u,v,tx,ty, provides solution of the algebraic system A*u=r, where A is the matrix composed of the matrices G and H, u is vector of reordered unknown boundary conditions at the nodes of the elements and r is the vector calculated from the matrices G and H and prescribed boundary conditions inner_sol(sub,xin) - returns tuple of inner values of displacements and stresses (u,v,sx,sy,sxy), where u,v,sx,sy,sxy are ndarrays of rank 1 and length n, where n is the length of the input xin, which is the tuple of inner point coordinate tuple (x,y), sub=1 or 2 is No. of subdomain, where the xins lie ''' def __init__(self,iplane,el1,el2,pn1,pn2,elements1,elements2): self.iplane=iplane self.el1=el1 self.el2=el2 self.pn1=pn1 self.pn2=pn2 self.elements1=elements1 self.elements2=elements2 self.set_param() print 'parameters for fortran subroutines were established ...' self.G_calc() print 'G matrices were calculated ...' self.H_calc() print 'H matrices were calculated ...' self.Ar_sol() print 'algebraic system was solved ...' def set_param(self): elements=(self.elements1,self.elements2) nelem=(len(elements[0]),len(elements[1])) self.nmax=max(nelem) self.x1,self.y1=zeros((2,self.nmax)),zeros((2,self.nmax)) self.x2,self.y2=zeros((2,self.nmax)),zeros((2,self.nmax)) self.xm,self.ym=zeros((2,self.nmax)),zeros((2,self.nmax)) self.el,self.pn,self.gl=zeros(2),zeros(2),zeros(2) self.n=zeros(2) self.ub,self.vb=zeros((2,self.nmax)),zeros((2,self.nmax)) self.pos,self.kcode=zeros((2,self.nmax)),zeros((2,self.nmax)) for ii in range(2): for jj in range(nelem[ii]): self.xm[ii,jj]=elements[ii][jj][0][0] self.ym[ii,jj]=elements[ii][jj][0][1] self.x1[ii,jj]=elements[ii][jj][1][0] self.y1[ii,jj]=elements[ii][jj][1][1] self.x2[ii,jj]=elements[ii][jj][2][0] self.y2[ii,jj]=elements[ii][jj][2][1] cislo,cislob=0,0 for ii in range(2): kk=-1 for jj in elements[ii]: cislo+=1 kk+=+1 if jj['b$']!='12': self.ub[ii,kk]=jj['bxy'][0] self.vb[ii,kk]=jj['bxy'][1] self.pos[ii,kk]=cislo if jj['b$']=='12': self.ub[ii,kk]=0. self.vb[ii,kk]=0. if ii==0: cislob=cislo self.pos[ii,kk]=cislob cislo+=1 else: self.pos[ii,kk]=cislob cislo-=1 cislob-=2 for ii in range(2): for jj in range(nelem[ii]): if elements[ii][jj]['b$']=='uu': self.kcode[ii,jj]=1 elif elements[ii][jj]['b$']=='tt': self.kcode[ii,jj]=2 elif elements[ii][jj]['b$']=='ut': self.kcode[ii,jj]=3 elif elements[ii][jj]['b$']=='tu': self.kcode[ii,jj]=4 else: self.kcode[ii,jj]=12 self.n[0],self.n[1]=nelem[0],nelem[1] self.nsum=sum(self.n) if self.iplane==0: self.el[0]=self.el1/(1.-self.pn1**2.) self.el[1]=self.el2/(1.-self.pn2**2.) self.pn[0]=self.pn1/(1.-self.pn1) self.pn[1]=self.pn2/(1.-self.pn2) else: self.el[0]=self.el1 self.el[1]=self.el2 self.pn[0]=self.pn1 self.pn[1]=self.pn2 self.gl[0]=self.el[0]/(2.*(1.+self.pn[0])) self.gl[1]=self.el[1]/(2.*(1.+self.pn[1])) def G_calc(self): self.G=bie.bie.gmatrel(self.n,self.x1,self.y1,self.x2,self.y2,self.xm,self.ym,self.gl,self.pn) def H_calc(self): self.H=bie.bie.hmatrel(self.n,self.x1,self.y1,self.x2,self.y2,self.xm,self.ym,self.gl,self.pn) def Ar_sol(self): self.A,self.r=bie.bie.abmatrel(self.n,self.nsum,self.kcode,self.pos,self.ub,self.vb,self.G,self.H) self.sol=solve(self.A,self.r) self.u,self.v,self.tx,self.ty=bie.bie.reorderel(self.n,self.kcode,self.pos,self.ub,self.vb,self.sol) def inner_sol(self,sub,x): nx=len(x) xin,yin=zeros(nx),zeros(nx) for ii in x: ind=x.index(ii) xin[ind]=ii[0] yin[ind]=ii[1] uin,vin=bie.bie.uvinter(sub,xin,yin,self.n,self.x1,self.y1,self.x2,self.y2, \ self.gl,self.pn,self.u,self.v,self.tx,self.ty) sxin,syin,sxyin=bie.bie.stressin(sub,xin,yin,self.n,self.x1,self.y1,self.x2,self.y2, \ self.gl,self.pn,self.u,self.v,self.tx,self.ty) return uin,vin,sxin,syin,sxyin