⭐ 欢迎来到虫虫下载站! | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站

📄 libearingfilm.for

📁 边界远程序,加入油膜压力的接触问题求解的边界元程序
💻 FOR
📖 第 1 页 / 共 5 页
字号:
     1      +SUMZ(JY)*COSBEZ(I,L,3)
	A(JZ,KX)=SUMX(JZ)*COSBEX(I,L,1)+SUMY(JZ)*COSBEY(I,L,1)
     1      +SUMZ(JZ)*COSBEZ(I,L,1)
	A(JZ,KY)=SUMX(JZ)*COSBEX(I,L,2)+SUMY(JZ)*COSBEY(I,L,2)
     1      +SUMZ(JZ)*COSBEZ(I,L,2)
	A(JZ,KZ)=SUMX(JZ)*COSBEX(I,L,3)+SUMY(JZ)*COSBEY(I,L,3)
     1      +SUMZ(JZ)*COSBEZ(I,L,3)
600   CONTINUE
	print *, ' Coefficent matrix has been formed !' 
!       THIS PROGRAM IS USED TO REDUCE THE COEFFICIENT MATRIX
!              BY GAUSS'REDUCTION
	print *, ' Begin to reduce......'
      CALL Gauss_R(I,NUMBE,A,P,AC,PC,msi,nsi,ksi)
	IF(I.EQ.1)THEN
	CALL EXCHANGE(A,P,A1,P1,msi,nsi)
	ENDIF
988   CONTINUE
	CALL CONTACT(ITERL,MULTI,numbe,U0,nttyp,nutyp,A,P,A1,P1,AC,PC,
     #	msi,nsi,ksi,xn,yn,zn,cosbex,cosbey,cosbez,kod,nord,kodt,
     # kode,mors,morp,morc,merp,merc,mord,mort,node,MNB,film,TFORCE)
	return
	end

!       MAIN PROGRAM HAS BEEN FINISHED.
!       NEXT PROGRAMS ARE SUBROUTINE PROGRAMS.
	
      SUBROUTINE Gauss_R(I,NUMBE,A,P,AC,PC,msi,nsi,ksi)
      COMMON/S8/MS(2),MD(6),NS(2),ND(2)
 	dimension A(msi,nsi),P(msi),AC(2,3*numbe,ksi),PC(2,3*numbe)
	MD2=NS(I)-ND(I)
      NB=ND(I)
	DO 760 II=1,NB
	WRITE(*,*)I,II
      DO 750 JJ=II+1,MS(I)
      XM=A(JJ,II)/A(II,II)
      DO 740 L=II,NS(I)
       A(JJ,L)=A(JJ,L)-A(II,L)*XM
740	CONTINUE
       P(JJ)=P(JJ)-P(II)*XM
750	CONTINUE 
	 760	CONTINUE
!-----------------------------------
770	DO 781 I1=1,3*NUMBE
	I2=I1+ND(I)
	PC(i,I1)=P(I2)
	DO 780 J1=1,MD2
	J2=J1+ND(I)
	AC(i,I1,J1)=A(I2,J2)
780	continue
!	write(11,*) ac(i,i1,i1),pc(i,i1)
781	CONTINUE
      RETURN
	END
	SUBROUTINE EXCHANGE(A,P,A1,P1,msi,nsi)
	COMMON/S8/MS(2),MD(6),NS(2),ND(2)
	DIMENSION A(msi,nsi),P(msi)
	DIMENSION A1(msi,nsi),P1(msi)
	DO 114 I=1,msi
	P1(I)=P(I)
	do 114 J=1,nsi
	A1(I,J)=A(I,J)
114	CONTINUE 
	END
	SUBROUTINE  JACOB(NODAL,XW,YW,XK,YK,ZK,K)
	COMMON/JA/XD1,YD1,ZD1,XD2,YD2,ZD2
  	COMMON/S6/FJACOB,COSBX,COSBY,COSBZ
	DIMENSION NODAL(8),DFXW(8),DFYW(8),XK(8),YK(8),ZK(8)
	DO 10 J=5,8
	DFXW(J)=0.
10      DFYW(J)=0.
	IF(NODAL(5).NE.0) DFXW(5)=-XW*(1.-YW)
	IF(NODAL(6).NE.0) DFYW(6)=-YW*(1.+XW)
	IF(NODAL(7).NE.0) DFXW(7)=-XW*(1.+YW)
	IF(NODAL(8).NE.0) DFYW(8)=-YW*(1.-XW)
	IF(NODAL(5).NE.0) DFYW(5)=-0.5*(1.-XW*XW)
	IF(NODAL(6).NE.0) DFXW(6)= 0.5*(1.-YW*YW)
	IF(NODAL(7).NE.0) DFYW(7)= 0.5*(1.-XW*XW)
	IF(NODAL(8).NE.0) DFXW(8)=-0.5*(1.-YW*YW)
	DFXW(1)=-0.25*(1.-YW)-0.5*(DFXW(5)+DFXW(8))
	DFXW(2)= 0.25*(1.-YW)-0.5*(DFXW(5)+DFXW(6))
	DFXW(3)= 0.25*(1.+YW)-0.5*(DFXW(6)+DFXW(7))
	DFXW(4)=-0.25*(1.+YW)-0.5*(DFXW(7)+DFXW(8))
	DFYW(1)=-0.25*(1.-XW)-0.5*(DFYW(5)+DFYW(8))
	DFYW(2)=-0.25*(1.+XW)-0.5*(DFYW(5)+DFYW(6))
	DFYW(3)= 0.25*(1.+XW)-0.5*(DFYW(6)+DFYW(7))
	DFYW(4)= 0.25*(1.-XW)-0.5*(DFYW(7)+DFYW(8))
	XD1=DFXW(1)*XK(1)+DFXW(2)*XK(2)+DFXW(3)*XK(3)+DFXW(4)*XK(4)
     1   +DFXW(5)*XK(5)+DFXW(6)*XK(6)+DFXW(7)*XK(7)+DFXW(8)*XK(8)
	YD1=DFXW(1)*YK(1)+DFXW(2)*YK(2)+DFXW(3)*YK(3)+DFXW(4)*YK(4)
     1   +DFXW(5)*YK(5)+DFXW(6)*YK(6)+DFXW(7)*YK(7)+DFXW(8)*YK(8)
	ZD1=DFXW(1)*ZK(1)+DFXW(2)*ZK(2)+DFXW(3)*ZK(3)+DFXW(4)*ZK(4)
     1   +DFXW(5)*ZK(5)+DFXW(6)*ZK(6)+DFXW(7)*ZK(7)+DFXW(8)*ZK(8)
	XD2=DFYW(1)*XK(1)+DFYW(2)*XK(2)+DFYW(3)*XK(3)+DFYW(4)*XK(4)
     1   +DFYW(5)*XK(5)+DFYW(6)*XK(6)+DFYW(7)*XK(7)+DFYW(8)*XK(8)
	YD2=DFYW(1)*YK(1)+DFYW(2)*YK(2)+DFYW(3)*YK(3)+DFYW(4)*YK(4)
     1   +DFYW(5)*YK(5)+DFYW(6)*YK(6)+DFYW(7)*YK(7)+DFYW(8)*YK(8)
	ZD2=DFYW(1)*ZK(1)+DFYW(2)*ZK(2)+DFYW(3)*ZK(3)+DFYW(4)*ZK(4)
     1   +DFYW(5)*ZK(5)+DFYW(6)*ZK(6)+DFYW(7)*ZK(7)+DFYW(8)*ZK(8)
	 IF(K.NE.0) GOTO 75
	G1=YD1*ZD2-ZD1*YD2
	G2=ZD1*XD2-XD1*ZD2
	G3=XD1*YD2-YD1*XD2
	FJACOB=SQRT(G1*G1+G2*G2+G3*G3)
	COSBX=G1/FJACOB
	COSBY=G2/FJACOB
	COSBZ=G3/FJACOB
75	RETURN
	END

	SUBROUTINE ASSUM(I,NUMBE,MULTI,kod,nord,kodt,mors,morp,
     #           morc,merp,merc,mord,mort,node,msi,MNB)
	COMMON/S7/NODBS(2),NUMBS(2),NODPS(2),NODCS(2),NEDPS(2),NEDCS(2)
	dimension Nord(2,msi/3,8),KODT(2,msi/3)
	dimension kod(2,msi/3)
	dimension Mors(2,numbe),Morp(2,numbe),Morc(2,numbe),Merp(2,numbe),
     #          Merc(2,numbe)
      dimension Mord(2,msi/3),Mort(2,msi/3),Node(2,msi/3,8)
	DIMENSION MNB(2,2500)
	NUMBE=0
	NODDS=0
	NODPS(I)=0
	NODCS(I)=0
	DO 30 M=1,NODBS(I)
	GOTO (25,25,25,25,25,25,25,25,10,20),KOD(I,M)
10    NUMBE=NUMBE+1
	MORD(I,M)=NUMBE
	MORS(I,NUMBE)=M
	NODPS(I)=NODPS(I)+1
	MORT(I,M)=NODPS(I)
	MORP(I,NODPS(I))=M
	GOTO 30
20    NUMBE=NUMBE+1
	MORD(I,M)=NUMBE
	MORS(I,NUMBE)=M
	NODCS(I)=NODCS(I)+1
	MORT(I,M)=NODCS(I)
	MORC(I,NODCS(I))=M
	GOTO 30
25    NODDS=NODDS+1
	MORD(I,M)=NODDS
	MNB(I,NODDS)=M
30    CONTINUE
	DO 40 J=1,NUMBE
	L=0
	DO 35 K=1,MULTI
35      NODE(I,J,K)=0
	DO 40 K=1,NUMBS(I)
	DO 40 N=1,8
	IF(NORD(I,K,N).NE.MORS(I,J)) GOTO 40
	L=L+1
	NODE(I,J,L)=K
40      CONTINUE
	NEDPS(I)=0
	NEDCS(I)=0
	DO 90 K=1,NUMBS(I)
	GOTO (90,70,80),KODT(I,K)
70    NEDPS(I)=NEDPS(I)+1
	MERP(I,NEDPS(I))=K
	GOTO 90
80    NEDCS(I)=NEDCS(I)+1
	MERC(I,NEDCS(I))=K
90      CONTINUE
	RETURN
	END

	SUBROUTINE INITL(ID)
	COMMON/S2/AXX(8),AXY(8),AXZ(8),AYX(8),AYY(8),AYZ(8),
     1            AZX(8),AZY(8),AZZ(8),BXX(8),BXY(8),BXZ(8),
     2            BYX(8),BYY(8),BYZ(8),BZX(8),BZY(8),BZZ(8)
C
	IF(ID.NE.1) GOTO 20
	AXX(1)=-1.
	AXX(4)=-1.
	AXX(8)=-1.
	AXX(2)=1.
	AXX(3)=1.
	AXX(6)=1.
	AXX(5)=0.
	AXX(7)=0.
	AYY(1)=-1.
	AYY(2)=-1.
	AYY(5)=-1.
	AYY(3)=1.
	AYY(4)=1.
	AYY(7)=1.
	AYY(6)=0.
	AYY(8)=0.
	RETURN
20      DO 25 I=1,8
	AXX(I)=0.
	AXY(I)=0.
	AXZ(I)=0.
	AYX(I)=0.
	AYY(I)=0.
	AYZ(I)=0.
	AZX(I)=0.
	AZY(I)=0.
	AZZ(I)=0.
	BXX(I)=0.
	BXY(I)=0.
	BXZ(I)=0.
	BYX(I)=0.
	BYY(I)=0.
	BYZ(I)=0.
	BZX(I)=0.
	BZY(I)=0.
25    BZZ(I)=0.
	RETURN
	END

	SUBROUTINE COEFF(ID,NODAL,XP,YP,ZP,XK,YK,ZK)
	COMMON/S1/PI,PR(2),PR1(2),PR2(2),PR3(2),PR4(2),CON(2),DSMIN2(2),
     1  DSMAX2(2),E(2),G(2)
	COMMON/S2/AXX(8),AXY(8),AXZ(8),AYX(8),AYY(8),AYZ(8),
     1            AZX(8),AZY(8),AZZ(8),BXX(8),BXY(8),BXZ(8),
     1            BYX(8),BYY(8),BYZ(8),BZX(8),BZY(8),BZZ(8)
	COMMON/S5/VECTLC(6,3),WTFUN(6,3)
	COMMON/S6/FJACOB,COSBX,COSBY,COSBZ
	DIMENSION NODAL(8),FN(8),XK(8),YK(8),ZK(8)
     C
	XQP=0.125*(XK(1)+XK(2)+XK(3)+XK(4)+XK(5)+XK(6)+XK(7)+XK(8))-XP
	YQP=0.125*(YK(1)+YK(2)+YK(3)+YK(4)+YK(5)+YK(6)+YK(7)+YK(8))-YP
	ZQP=0.125*(ZK(1)+ZK(2)+ZK(3)+ZK(4)+ZK(5)+ZK(6)+ZK(7)+ZK(8))-ZP
	DSTAN2=XQP*XQP+YQP*YQP+ZQP*ZQP
	IF(DSTAN2.LE.DSMIN2(ID)) GOTO 10
	IF(DSTAN2.LE.DSMAX2(ID)) GOTO 20
	L=1
	M=2
 	N=1
	GOTO 30
10    L=1
	M=6
	N=2
	GOTO 30
20      L=3
	M=6
	N=1
30      DO 40 J=5,8
40      FN(J)=0.
	DO 50 IA=L,M
	XL=VECTLC(IA,N)
	WX=WTFUN(IA,N)
	DO 50 JA=L,M
	YL=VECTLC(JA,N)
	WY=WTFUN(JA,N)
	IF(NODAL(5).NE.0) FN(5)=0.5*(1.-XL*XL)*(1.-YL)
	IF(NODAL(7).NE.0) FN(7)=0.5*(1.-XL*XL)*(1.+YL)
	IF(NODAL(6).NE.0) FN(6)=0.5*(1.-YL*YL)*(1.+XL)
	IF(NODAL(8).NE.0) FN(8)=0.5*(1.-YL*YL)*(1.-XL)
	FN(1)=0.25*(1.-XL)*(1.-YL)-0.5*(FN(5)+FN(8))
	FN(2)=0.25*(1.+XL)*(1.-YL)-0.5*(FN(5)+FN(6))
	FN(3)=0.25*(1.+XL)*(1.+YL)-0.5*(FN(6)+FN(7))
	FN(4)=0.25*(1.-XL)*(1.+YL)-0.5*(FN(7)+FN(8))
	XQP=0.
	YQP=0.
	ZQP=0.
	DO 44 MM=1,8
	XQP=XQP+XK(MM)*FN(MM)
	YQP=YQP+YK(MM)*FN(MM)
44      ZQP=ZQP+ZK(MM)*FN(MM)
	XQP=XQP-XP
	YQP=YQP-YP
	ZQP=ZQP-ZP
	RQ2=XQP*XQP+YQP*YQP+ZQP*ZQP
	RQ1=SQRT(RQ2)
	RXX2=XQP*XQP/RQ2
	RYY2=YQP*YQP/RQ2
	RZZ2=ZQP*ZQP/RQ2
	RXY2=XQP*YQP/RQ2
	RYZ2=YQP*ZQP/RQ2
	RXZ2=XQP*ZQP/RQ2
	CALL JACOB(NODAL,XL,YL,XK,YK,ZK,0)
	FL1=CON(ID)*WX*WY*FJACOB/RQ1
	FL2=PR1(ID)*FL1/RQ2
	FL3=FL2*(XQP*COSBX+YQP*COSBY+ZQP*COSBZ)
	FL4=-3.*FL3/PR1(ID)
	AXXT=-FL3+FL4*RXX2
	AYYT=-FL3+FL4*RYY2
	AZZT=-FL3+FL4*RZZ2
	AXYT=FL4*RXY2+FL2*(XQP*COSBY-YQP*COSBX)
	AXZT=FL4*RXZ2+FL2*(XQP*COSBZ-ZQP*COSBX)
	AYXT=FL4*RXY2+FL2*(YQP*COSBX-XQP*COSBY)
	AYZT=FL4*RYZ2+FL2*(YQP*COSBZ-ZQP*COSBY)
	AZXT=FL4*RXZ2+FL2*(ZQP*COSBX-XQP*COSBZ)
	AZYT=FL4*RYZ2+FL2*(ZQP*COSBY-YQP*COSBZ)
	BXXT=FL1*(PR4(ID)+RXX2)
	BYYT=FL1*(PR4(ID)+RYY2)
	BZZT=FL1*(PR4(ID)+RZZ2)
	BXYT=FL1*RXY2
	BXZT=FL1*RXZ2
	BYZT=FL1*RYZ2
	BYXT=BXYT
	BZXT=BXZT
	BZYT=BYZT
	DO 50 J=1,8
	AXX(J)=AXX(J)+FN(J)*AXXT
	AXY(J)=AXY(J)+FN(J)*AXYT
	AXZ(J)=AXZ(J)+FN(J)*AXZT
	AYX(J)=AYX(J)+FN(J)*AYXT
	AYY(J)=AYY(J)+FN(J)*AYYT
	AYZ(J)=AYZ(J)+FN(J)*AYZT
	AZX(J)=AZX(J)+FN(J)*AZXT
	AZY(J)=AZY(J)+FN(J)*AZYT
	AZZ(J)=AZZ(J)+FN(J)*AZZT
c	IF(KK.GT.1) GOTO 50
	BXX(J)=BXX(J)+FN(J)*BXXT
	BXY(J)=BXY(J)+FN(J)*BXYT
	BXZ(J)=BXZ(J)+FN(J)*BXZT
	BYX(J)=BYX(J)+FN(J)*BYXT
	BYY(J)=BYY(J)+FN(J)*BYYT
	BYZ(J)=BYZ(J)+FN(J)*BYZT
	BZX(J)=BZX(J)+FN(J)*BZXT
	BZY(J)=BZY(J)+FN(J)*BZYT
	BZZ(J)=BZZ(J)+FN(J)*BZZT
50      CONTINUE
90      RETURN
	END

	SUBROUTINE IMPLE(ID,NODAL,I,XP,YP,ZP,XK,YK,ZK)
	COMMON/S1/PI,PR(2),PR1(2),PR2(2),PR3(2),PR4(2),CON(2),DSMIN2(2),
     1   DSMAX2(2),E(2),G(2)
	COMMON/S2/AXX(8),AXY(8),AXZ(8),AYX(8),AYY(8),AYZ(8),
     1            AZX(8),AZY(8),AZZ(8),BXX(8),BXY(8),BXZ(8),
     1            BYX(8),BYY(8),BYZ(8),BZX(8),BZY(8),BZZ(8)
	COMMON/S5/VECTLC(6,3),WTFUN(6,3)
	COMMON/S6/FJACOB,COSBX,COSBY,COSBZ
	COMMON/S13/COEF(3,10),COOR(3,10)
	DIMENSION NODAL(8),FN(8),XK(8),YK(8),ZK(8)
	DIMENSION XT(3),YT(3),FJ(3)
	M=2
	IF(I.GT.4) M=3
	DO 10 J=5,8
10      FN(J)=0.
	DO 90 IA=1,6
	XL=VECTLC(IA,2)
	WX=WTFUN(IA,2)
	DO 90 JA=1,6
	YL=VECTLC(JA,2)
	WY=WTFUN(JA,2)
	IF(I.GT.4) GOTO 40
	FL1=TAN(0.125*(1.+XL)*PI)
	FJ(1)=0.125*(1.+YL)*PI*(1.+FL1*FL1)
	FJ(2)=FJ(1)
	GOTO (15,20,25,30),I
15      XT(1)=YL
	YT(1)=(1.+YL)*FL1-1.
	XT(2)=YT(1)
	YT(2)=XT(1)
	GOTO 70
20      XT(1)=1.-(1.+YL)*FL1
	YT(1)=YL
	XT(2)=-YT(1)
	YT(2)=-XT(1)
	GOTO 70
25      XT(1)=-YL
	YT(1)=1.-(1.+YL)*FL1
	XT(2)=YT(1)
	YT(2)=XT(1)
	GOTO 70
30      XT(1)=(1.+YL)*FL1-1.
	YT(1)=-YL
	XT(2)=-YT(1)
	YT(2)=-XT(1)
	GOTO 70
40      FL1=TAN(0.5*(1.+XL)*ATAN(2.E0))
	FL2=TAN(XL*ATAN(0.5E0))
	FJ(1)=0.125*(1.+YL)*ATAN(2.E0)*(1.+FL1*FL1)
	FJ(2)=(1.+YL)*ATAN(0.5E0)*(1.+FL2*FL2)
	FJ(3)=FJ(1)
	GOTO (15,20,25,30,45,50,55,60),I
45      XT(1)=0.5*(1.+YL)
	YT(1)=0.5*(1.+YL)*FL1-1.
	XT(2)=-(1.+YL)*FL2
	YT(2)=YL
	XT(3)=-XT(1)
	YT(3)=YT(1)
	GOTO 70
50      XT(1)=1.-0.5*(1.+YL)*FL1
	YT(1)=0.5*(1.+YL)
	XT(2)=-YL
	YT(2)=-(1.+YL)*FL2
	XT(3)=XT(1)
	YT(3)=-YT(1)
	GOTO 70
55      XT(1)=-0.5*(1.+YL)
	YT(1)=1.-0.5*(1.+YL)*FL1
	YT(2)=-YL
	XT(2)=(1.+YL)*FL2
	XT(3)=-XT(1)
	YT(3)=YT(1)
	GOTO 70
60      XT(1)=0.5*(1.+YL)*FL1-1.
	YT(1)=-0.5*(1.+YL)
	XT(2)=YL
	YT(2)=(1.+YL)*FL2
	XT(3)=XT(1)
	YT(3)=-YT(1)
70      DO 90 K=1,M
	XW=XT(K)
	YW=YT(K)
	IF(NODAL(5).NE.0) FN(5)=0.5*(1.-XW*XW)*(1.-YW)
	IF(NODAL(7).NE.0) FN(7)=0.5*(1.-XW*XW)*(1.+YW)
	IF(NODAL(6).NE.0) FN(6)=0.5*(1.-YW*YW)*(1.+XW)
	IF(NODAL(8).NE.0) FN(8)=0.5*(1.-YW*YW)*(1.-XW)
	FN(1)=0.25*(1.-XW)*(1.-YW)-0.5*(FN(5)+FN(8))
	FN(2)=0.25*(1.+XW)*(1.-YW)-0.5*(FN(5)+FN(6))
	FN(3)=0.25*(1.+XW)*(1.+YW)-0.5*(FN(6)+FN(7))
	FN(4)=0.25*(1.-XW)*(1.+YW)-0.5*(FN(7)+FN(8))
	XQP=0.
	YQP=0.
	ZQP=0.

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -